function main(Macro_Bkhaul,Pico_Bkhaul)

%M = 28;
num_user = 30;
K = num_user;

debug_mode = 0;

if debug_mode == 0
    
    current_dir = cd;
    
    %cd('C:\Users\Admin\Desktop\cvx-w64\cvx');
    cd('/fs2/comm/wygrp/ppatil/cvx-a64/cvx');
    cvx_setup
    
    cd(current_dir);
end

%simulation parameters
Sim.snr_gap_dB = 9;
Sim.snr_gap = dB2dec(Sim.snr_gap_dB);
Sim.sqr_gap = 2.7;
Sim.sqr_gap_dB = dec2dB(Sim.sqr_gap);
noise_dBm = [-130 -140 -150 -160];
Sim.noise_dBm = noise_dBm(3);
Sim.noise = dB2dec(Sim.noise_dBm);
Sim.Num_Cell = 7;%total num of cells
Sim.Num_BS = 4; %num base-stations per cell
Sim.Num_Macro = 1; %num macro base-stations
Sim.Num_Pico = 3; %num pico base-stations
Sim.Num_User = K; %num users per cell
Sim.Num_Ant_Macro = 4; %num ant's per macro BS
Sim.Num_Ant_Pico = 2; %num ant's per pico BS
Sim.Num_Rx_Ant = 2;
Sim.B2Bdist = 0.8; %km
Sim.BW = 10e6;%bandwidth
Sim.P_max_Macro_dBm = 43 - dec2dB(Sim.BW); %dBm/Hz
Sim.P_max_Pico_dBm = 30 - dec2dB(Sim.BW); %dBm/Hz
Sim.MinWMMSEIter = 15; %minimum WMMSE iteration
Sim.snr_bias = 0;%snr_bias;
Sim.rcv_power_gap = 0;%rcv_power_gap;

Power_Mode = 3;
intf_2nd = zeros(Sim.Num_Rx_Ant,Sim.Num_Rx_Ant);
%[Chn, path_loss_dB, dist] = GenHetNetChnSingleCell(Sim);

%Chn_file_name = sprintf('Chn%d_K%d.mat',Chn_idx,K);
%load(Chn_file_name);
load Chn.mat;

M= 14;

%M = 2; %maximum number of serving BSs per User
if debug_mode == 0
    save_file_name = sprintf('/fs2/comm/wygrp/ppatil/Results_20151014_1400_FullCSI/DataSharing_MIMO_Macro%d_Pico%d_NumUser%d.mat',Macro_Bkhaul,Pico_Bkhaul,num_user);
    save(save_file_name,'Sim');
end
% 
% Bkhaul_ratio = 4; %the ratio of the amount of backhaul between macro BS and pico BS
% 
% Pico_Bkhaul_Limit = [2 4 6 8];
% Macro_Bkhaul_Limit = Pico_Bkhaul_Limit*Bkhaul_ratio;


Cells = struct([]);
User = struct([]);
BS = struct([]);

for c = 1:Sim.Num_Cell
    for b = 1:Sim.Num_BS
        if b <= Sim.Num_Macro
            BS(c,b).total_power_budget = dB2dec(Sim.P_max_Macro_dBm);
            BS(c,b).total_bkhaul_budget = Macro_Bkhaul;%Macro_Bkhaul_Limit(Bkhaul_idx);
            BS(c,b).Per_Ant_Power = dB2dec(Sim.P_max_Macro_dBm);
        else
            BS(c,b).total_power_budget = dB2dec(Sim.P_max_Pico_dBm);
            BS(c,b).total_bkhaul_budget = Pico_Bkhaul;%Pico_Bkhaul_Limit(Bkhaul_idx);
            BS(c,b).Per_Ant_Power = dB2dec(Sim.P_max_Pico_dBm);
        end
        
        BS(c,b).compression_power = 0;
        BS(c,b).current_compression_bkhaul = 1e10; %some large value
        
    end
end


%initialization
for iCell = 1:Sim.Num_Cell
    for k = 1:Sim.Num_User
        User(iCell,k).long_avg_rate = 1e-3;
        User(iCell,k).avg_rate = 0;
    end
    Cells(iCell).UnScheduled_User = 1:Sim.Num_User;
end

[User Asso_User BS] = Neighbour(Sim, Sim.Num_Cell, Sim.Num_User, path_loss_dB, BS,User,M,Sim.noise_dBm);

out_iter = 1;
previous_utility = 1;
current_utility = 2;
Log_Data = struct([]);
Num_Scheduled_User = Sim.Num_User;
%%
while out_iter <= 70
    
    %abs((current_utility - previous_utility)/current_utility) > 1e-2 || out_iter <=70
    
    previous_utility = current_utility;
    fprintf('$$$$$$ User Scheduling iter = %d $$$$$$$ \n', out_iter);
    tic;
    %User scheduler
    for iCell = 1:Sim.Num_Cell
        [Cells(iCell).Scheduled_User Cells(iCell).UnScheduled_User]= Round_Robin(Cells(iCell).UnScheduled_User,Num_Scheduled_User);

        %         Log_Data(out_iter).Cells(iCell).Scheduled_User = Cells(iCell).Scheduled_User;
        %         Log_Data(out_iter).Cells(iCell).UnScheduled_User = Cells(iCell).UnScheduled_User; % log user scheduling
        if isempty(Cells(iCell).UnScheduled_User)
            Cells(iCell).UnScheduled_User = 1:Sim.Num_User;
        end
    end
    
    Num_Serving_User = zeros(Sim.Num_Cell,Sim.Num_BS); %num of users each BS supports at the begining
    for iCell = 1:Sim.Num_Cell
        for ik = 1:length(Cells(iCell).Scheduled_User)
            k = Cells(iCell).Scheduled_User(ik);
            User(iCell,k).ServingCluster = User(iCell,k).NeighbourCluster;
            User(iCell,k).ServingAnt = 0;
            for iBS = 1:size(User(iCell,k).ServingCluster,1)
                
                curr_cell = User(iCell,k).ServingCluster(iBS,1);
                curr_BS = User(iCell,k).ServingCluster(iBS,2);
                User(iCell,k).bkhaul_wgt(curr_cell,curr_BS) = 0.01;
                Num_Serving_User(curr_cell,curr_BS) = Num_Serving_User(curr_cell,curr_BS) + 1;
                
                if curr_BS <= Sim.Num_Macro %Macro BS
                    User(iCell,k).ServingAnt = User(iCell,k).ServingAnt + Sim.Num_Ant_Macro;
                else
                    User(iCell,k).ServingAnt = User(iCell,k).ServingAnt + Sim.Num_Ant_Pico;
                end
            end
        end
    end
    
    %Initialize wgts
    max_wgt = 0;
    for iCell = 1:Sim.Num_Cell
        for ik = 1:length(Cells(iCell).Scheduled_User)
            k = Cells(iCell).Scheduled_User(ik);
            User(iCell,k).wgt_wsr = 1/User(iCell,k).long_avg_rate;
            User(iCell,k).prev_inst_rate = 0.01;
            if User(iCell,k).wgt_wsr > max_wgt
                max_wgt = User(iCell,k).wgt_wsr;
            end
        end
    end
    %normalize user wgts by max_wgt
    for iCell = 1:Sim.Num_Cell
        for ik = 1:length(Cells(iCell).Scheduled_User)
            k = Cells(iCell).Scheduled_User(ik);
            User(iCell,k).wgt_wsr = User(iCell,k).wgt_wsr/max_wgt;
            
        end
    end
    
    % initialize TX beamformers and sparse Matrix
    %load Beam_Initial.mat;
    for l = 1:Sim.Num_Cell
        for kk = 1:length(Cells(l).Scheduled_User)
            k = Cells(l).Scheduled_User(kk); %scheduled user in cell l
            Num_ServingAnt = User(l,k).ServingAnt;
            temp_tx = randn(Num_ServingAnt,1) + 1i*randn(Num_ServingAnt,1);
            
            ant_head = 1;
            temp_tx2 = [];
            for m = 1:size(User(l,k).ServingCluster,1)
                
                curr_cell = User(l,k).ServingCluster(m,1);
                curr_BS = User(l,k).ServingCluster(m,2);
                if curr_BS <= Sim.Num_Macro
                    Num_Tx_Ant = Sim.Num_Ant_Macro;
                else
                    Num_Tx_Ant = Sim.Num_Ant_Pico;
                end
                ant_tail = ant_head + Num_Tx_Ant - 1;
                norm_tx = temp_tx(ant_head:ant_tail)/norm(temp_tx(ant_head:ant_tail));
                temp_tx2 = [temp_tx2;sqrt(BS(curr_cell,curr_BS).total_power_budget/Num_Serving_User(curr_cell,curr_BS))*norm_tx];
                ant_head = ant_tail+1;
          %      User(l,k).BS_Power(m) = BS(curr_cell,curr_BS).total_power_budget/Num_Serving_User(curr_cell,curr_BS);
            end
            
            User(l,k).beam_tx = temp_tx2;
           % User(l,k).beam_tx = User_initial(l,k).beam_tx;
            if length(User(l,k).beam_tx) ~= User(l,k).ServingAnt
                fprintf('$$$ Tx Initialization Not Right $$$ \n');
            end
        end
    end
    
    %%% Fix BF's, optimize rate and backhaul for compression using iterior-point method
    
    % %%%%% parameters define %%%%%%%%%
    % %[User,Cells,BS,Beta] = parameter_def_rate_backhaul(User,Cells,BS,Sim,Chn);
    %
    % %%% initialize the rate and the backhaul for compression %%%
    % for c = 1:Sim.Num_Cell
    %
    %     for b = 1:Sim.Num_BS
    %         BS(c,b).current_compression_bkhaul = BS(c,b).current_min_compression_bkhaul;
    %     end
    %
    %     for ik = 1:length(Cell(c).Scheduled_User)
    %         k = Cell(c).Scheduled_User(ik);
    %         User(c,k).current_rate = 0;
    %     end
    %
    % end
    %
    % t=1;
    % while t <=1e10
    %
    %     %for each t, optimize the log-barrier function
    %     while 1 >0
    %         %%% find the Newton's direction %%%
    %         [User, BS, decreasement, curr_obj] = Find_Newton_Direction(User,Cells,BS,Sim,Beta,t);
    %
    %         %%% backtrack line search %%%
    %         [step_size, new_curr_obj, User, BS] = Backtrack_Line_Search(User,Cells,BS,Sim,Beta,t,decreasement,curr_obj);
    %
    %         %%% update variables
    %         for c = 1:Sim.Num_Cell
    %             for ik = 1:length(Cells(c).Scheduled_User)
    %                 k = Cells(c).Scheduled_User(ik);
    %                 User(c,k).current_rate = User(c,k).new_rate;
    %
    %             end
    %         end
    %         for c = 1:Sim.Num_Cell
    %             for b = 1:Sim.Num_BS
    %
    %                 BS(c,b).current_compression_bkhaul = BS(c,b).new_compression_bkhaul;
    %             end
    %         end
    %
    %         %stop criteria
    %         if sqrt(-decreasement) <= 1e-3
    %             break;
    %         end
    %     end
    %
    %     %update t
    %     t = 10*t;
    %
    % end
    %
    %
    % %% Fix backhaul for compression and quantization noise, design BF's using WMMSE approach
    
    %define the constant parameters to be used in optimizing the bf's
    [User, BS] = parameter_def_bf_optimize(User,Cells,BS,Sim,Chn);
    
%     previous_wsr = 1;
%     while 1>0
        %WMMSE approach to optimize bf's
        [User,sumrate, wgt_sumrate, Per_BS_Bkhaul, Cells, BS] = WSR_Solver_WMMSE(Cells,User,Chn,Sim,intf_2nd,Power_Mode,BS);
        
        
%         [User, max_diff] = PowerMinAltMthd(User,Cells,Chn,BS,Sim,intf_2nd);
%         
%         max_rate_diff = ...
%             rate_DL_after_power_min(User,Cells,Chn,Sim,intf_2nd);
%         
%         %      diff = (wgt_sumrate_after_power_min - wgt_sumrate)/wgt_sumrate;
%         fprintf('$$$$$  max rate diff to log(1+SINR) after power minimization = %f $$$$$\n', max_rate_diff);
%         
%         diff = (wgt_sumrate - previous_wsr)/wgt_sumrate;
%         fprintf('$$$$$ wsr diff with previous iteration = %f  $$$$$\n',diff);
%         
%         if abs(diff) < 1e-2
%             break;
%         end
%         
%         previous_wsr = wgt_sumrate;
%         
%     end
    
    
    Log_Data(out_iter).Per_BS_Bkhaul = Per_BS_Bkhaul;
    
    [Log_Utility,User,Log_Data] = Compute_Utility(Sim.Num_Cell,Sim.Num_User,User,Cells,BS,out_iter,Sim.Num_BS,Log_Data);
    
    fprintf('$$$$ Current Log Utility = %f $$$$ \n\n\n\n',Log_Utility);
    
    current_utility = Log_Utility;
    
    toc;
    Log_Data(out_iter).elapsetime = toc;
    
    
    out_iter = out_iter + 1;
    
    if debug_mode == 0
    save(save_file_name, 'Log_Data', 'User', 'BS', 'Cells','Sim');
    end
    
end
