function [User Asso_User BS] = Neighbour(Sim, Num_Cell, K, path_loss_dB, BS,User,M,noise_dBm)

%%% -----------------------------------------------------------------------
%%% This function finds the nearest M neighbouring BS's for every user
%%% accoring to the pathloss and BS trasnmit power budget
%%% -----------------------------------------------------------------------
Total_BS = Sim.Num_BS;
%User = struct([]);
% Tx_Power = [P_max_Macro_dBm*ones(1,L_Macro) P_max_Pico_dBm*ones(1,L_Pico)];
% Tx_Power_Matrix = kron(Tx_Power,ones(Num_Cell,1));
Asso_User = zeros(M,Num_Cell,Total_BS);
%num of users claiming each BS as the strongest BS

per_BS_power_constraint_matrix = zeros(Num_Cell,Total_BS);
snr_bias_mtx = zeros(Num_Cell,Total_BS);

for c=1:Num_Cell
    for b = 1:Total_BS
        if b <= Sim.Num_Macro
            per_BS_power_constraint_matrix(c,b) = dB2dec(Sim.P_max_Macro_dBm);
        else
            per_BS_power_constraint_matrix(c,b) = dB2dec(Sim.P_max_Pico_dBm);
            snr_bias_mtx(c,b) = Sim.snr_bias;
        end
        BS(c,b).asso_user = [];
    end
end


for iCell = 1:Num_Cell
    for k = 1:K
        
        path_loss = squeeze(path_loss_dB(:,:,iCell,k));
        Rcv_Power = dec2dB(per_BS_power_constraint_matrix) - path_loss + snr_bias_mtx;
        Rcv_Power_vec =reshape(Rcv_Power',[],1);
        
        %         sum_rcv_power = sum(dB2dec(Rcv_Power_vec));
        %         intf = sum_rcv_power - Rcv_Power_vec;
        %         intf_noist = intf + dB2de
        
        [SortedPower,SortedBS] = sort(Rcv_Power_vec,'descend');
        NeighbourCluster = SortedBS(1:M);
        %User(iCell,k).NeighbourCluster = SortedBS(find(SortedPower(1) - SortedPower <= Sim.rcv_power_gap));
        
        User(iCell,k).NeighbourCluster = zeros(length(NeighbourCluster),2);
        
        for c = 1:length(NeighbourCluster)
            BS_idx = NeighbourCluster(c);
            cell_idx = ceil(BS_idx / Total_BS);
            bs_idx = BS_idx - (cell_idx-1)*Total_BS;
            User(iCell,k).NeighbourCluster(c,:) = [cell_idx bs_idx];
            
            %User(iCell,k).ServingCluster(c,1) = ceil(BS_idx / Total_BS);
            %User(iCell,k).ServingCluster(c,2) = BS_idx - (User(iCell,k).ServingCluster(c,1)-1)*Total_BS;
            
            BS(cell_idx,bs_idx).asso_user = [BS(cell_idx,bs_idx).asso_user; iCell k];
            
        end
        
        User(iCell,k).data_sharing_BSs = User(iCell,k).NeighbourCluster;
        
        for m = 1:M
            Cell_index = ceil(SortedBS(m)/Total_BS);
            BS_index = mod(SortedBS(m)-1,Total_BS)+1;
            User(iCell,k).Neighbour(m).Cell_BS = [Cell_index BS_index];
            
            Asso_User(m, Cell_index, BS_index) = Asso_User(m, Cell_index, BS_index) + 1;
            
        end
        
    end
end
