复现不出来求助！！！matlab仿真复现有偿求助：基于信息论的密钥容量分析及验证-仿真实现图2-3和图2-4

尝试复现很久没有成功，关键无线信道矢量的理解、均值矢量的计算、协方差矩阵的计算、信噪比关系还有相干时间、发送功率等等很乱，感觉自己生成的信道矢量不对，以及相对速度怎么莱斯信道和瑞利信道中体现以及这两个信道是如何=根据上述ITU-VehA生成的h。（2-21）函数已经i邪恶好了。请各位大佬帮助，可有偿。

（2-21）

function C_key = calculate_key_capacity(T1, P, L, sigma_w1, sigma_w2, gamma, m_eq)
    % 计算密钥容量C_key
    product = 1;
    for l = 1:L
        % product = product * (1 - (T1^2 * P^2 * (gamma(l) - (m_eq(l))^2)^2) / ((T1 * P * gamma(l) + sigma_w1) * (T1 * P * gamma(l) + sigma_w2)));
        product = product *(1-((T1^2 * P^2 * (gamma(l) - (m_eq(l))^2)^2)/((T1 * P * gamma(l) + sigma_w1) * (T1 * P * gamma(l) + sigma_w2))));
    end
    C_key = (-1/(4 * T1) )* (log(product));
end

信道建模：代码有误、求大佬建议 尤其是相干时间T和T1、瑞利和莱斯信道产生无线信道矢量、均值矢量的计算等等

f1 = 2e9; % 载频1，单位：赫兹
f2 = 3e9; % 载频2，单位：赫兹
tau=[0,310,710,1090,1730,2510]/(10^9); % 相对时延 单位ns
pdb=[0,-1,-9,-10,-15,-20];             % 每径增益 单位dB
N = length(tau);%N个并行的子流、循环矩阵的大小
L = length(tau); % 路径数
sigma = 10.^(pdb/10); % 各径的方差
h_Rayleigh = zeros(L, 1); % 初始化h_Rician向量
h_Rician = zeros(L, 1); % 初始化h_Rician向量
h_Rayleigh_FDD = zeros(2*L-1, 1);
h_Rician_FDD = zeros(2*L-1, 1);
k_dB = 3; % 莱斯因子，3dB
K = 10^(k_dB/10);
Ps = 1; % 传输功率为1
T = 0.016;
T1 = T/2;
SNR_dB = 0:4:40;
C_key_Rician = zeros(1,length(SNR_dB));
C_key_Rayleigh = zeros(1,length(SNR_dB));
C_key_Rician_FDD = zeros(1,length(SNR_dB));
C_key_Rayleigh_FDD = zeros(1,length(SNR_dB));

v = 120000/3600; % 定义相对移动速度，单位为m/s
c = 3e8; % 光速，单位为m/s
% 计算多普勒频移
fd1 = v/c * f1;
fd2 = v/c * f2;



simu_times = 10000;

for kk = 1:1:simu_times

    for l = 1:L

        doppler_shift_phase1 = exp(1i * 2 * pi * fd1 * tau(l));
        doppler_shift_phase2 = exp(1i * 2 * pi * fd2 * tau(l));

        h_Rayleigh(l) = sqrt(sigma(l)) * randn + 1i * sqrt(sigma(l)) * randn ; % 零均值高斯变量Rayleigh信道
        h_Rician(l) = sqrt(K/(K+1))  + sqrt(1/(K+1)) * h_Rayleigh(l); % Rician信道
    
        h_Rayleigh(l) = h_Rayleigh(l) * doppler_shift_phase1;
        h_Rician(l) = h_Rician(l) * doppler_shift_phase1;

        h_Rayleigh_2(l) = sqrt(sigma(l)) * randn + 1i * sqrt(sigma(l)) * randn ; % 零均值高斯变量Rayleigh信道
        h_Rician_2(l) = sqrt(K/(K+1))  + sqrt(1/(K+1)) * h_Rayleigh_2(l); % Rician信道
        delay_2 = exp(1i * 2 * pi *f1* tau(l)); % 考虑信道时延
        h_Rayleigh_2(l) = h_Rayleigh_2(l) * delay_2;
        h_Rician_2(l) = h_Rician_2(l) * delay_2;

    end



    h_Rayleigh_FDD = conv(h_Rayleigh,h_Rayleigh_2);
    h_Rician_FDD = conv(h_Rician,h_Rician_2);

    V_L_Rayleigh = diag(h_Rayleigh*h_Rayleigh');
    V_L_Rician = diag(h_Rician*h_Rician');
    % 计算均值矢量
    m_eq_Rayleigh = zeros(1,L);
    % mu_Rician = abs(sum(h_Rician)/L);
    mu_Rician = sqrt(K/(K+1))*sigma(1);

    m_eq_Rician = zeros(1,L);
    m_eq_Rician(1) = mu_Rician;

    V_L_Rayleigh_FDD = diag(h_Rayleigh_FDD'*h_Rayleigh_FDD);
    V_L_Rician_FDD = diag(h_Rician_FDD'*h_Rician_FDD);
    m_eq_Rayleigh_FDD = zeros(1,2*L-1);
    mu_Rician_FDD = sqrt(sigma(1));
    m_eq_Rician_FDD = zeros(1,2*L-1);
    m_eq_Rician_FDD(1) = mu_Rician_FDD;


    for k = 1:length(SNR_dB)

        Pn = Ps/SNR_dB(k); % 计算传输信号功率
        C_key_temp_Rician(k) = calculate_key_capacity(T1, Ps, L, Pn, Pn, V_L_Rician, m_eq_Rician);
        C_key_temp_Rayleigh(k) = calculate_key_capacity(T1, Ps, L, Pn, Pn, V_L_Rayleigh, m_eq_Rayleigh);
        C_key_temp_Rician_FDD(k) = calculate_key_capacity(2*T1, Ps, 2*L-1, Pn, Pn, V_L_Rician_FDD, m_eq_Rician_FDD);
        C_key_temp_Rayleigh_FDD(k) = calculate_key_capacity(2*T1, Ps, 2*L-1, Pn, Pn, V_L_Rayleigh_FDD, m_eq_Rayleigh_FDD);
    end
    C_key_Rician = C_key_Rician + C_key_temp_Rician;
    C_key_Rayleigh = C_key_Rayleigh + C_key_temp_Rayleigh;
    C_key_Rician_FDD = C_key_Rician_FDD + C_key_temp_Rician_FDD;
    C_key_Rayleigh_FDD = C_key_Rayleigh_FDD + C_key_temp_Rayleigh_FDD;
end

C_key_Rician_mean = C_key_Rician./simu_times;
C_key_Rayleigh_mean = C_key_Rayleigh./simu_times;
C_key_Rician_mean_FDD = C_key_Rician_FDD./simu_times;
C_key_Rayleigh_mean_FDD = C_key_Rayleigh_FDD./simu_times;

disp('仿真次数:')
disp(simu_times)
disp('TDD_Rician密钥容量:')
disp(C_key_Rician_mean)
disp('TDD_Rayleigh密钥容量:')
disp(C_key_Rayleigh_mean)
disp('FDD_Rician密钥容量:')
disp(C_key_Rician_mean_FDD)
disp('FDD_Rayleigh密钥容量:')
disp(C_key_Rayleigh_mean_FDD)