function [Psi_k,theta_k,sum_rate] = H_MIMO_precoding(K,NumdS,N_F,dS,G_k,Psi_k,theta_k,f_base,P_T,sigma_2,L_x,L_y)

Q=4;
sum_rate_temp=zeros(Q,1);
for q=1:Q

[sum_rate_temp(q),~] = Capacity_calculate(K,NumdS,dS,G_k,Psi_k,theta_k,sigma_2);

rho_k = rho_k_update(K,NumdS,dS,G_k,Psi_k,theta_k);

Psi_k = Psi_update(K,NumdS,dS,rho_k,G_k,theta_k,sigma_2);

h_k = h_k_update(K,NumdS,G_k,Psi_k);

alpha_k_uv = alpha_k_uv_update(K,NumdS,N_F,dS,h_k,f_base,L_x,L_y);
alpha_k = alpha_k_uv_2_alpha_k(K,N_F,alpha_k_uv);

w_k = w_k_optimize(K,dS,NumdS,N_F,f_base,rho_k,alpha_k,P_T); %二分法优化w_k

w_k_uv = w_k_2_w_k_uv(K,N_F,w_k);
theta_k = w_k_2_theta_k(K,NumdS,N_F,w_k_uv,f_base);    %转换为发射机电流密度

temp=0;
for n=1:NumdS
    for k=1:K         
        temp=temp+norm(theta_k(:,k,n))^2*dS;
    end
end
theta_k=theta_k*sqrt(P_T)/sqrt(temp);

% [sum_rate_temp(q),SINR_k] = Capacity_calculate(K,NumdS,dS,G_k,Psi_k,theta_k,sigma_2);
end

% plot(sum_rate_temp);

sum_rate=max(sum_rate_temp);

end

