​微信公众号：智能电磁频谱算法 QQ交流群：949444104

主要内容

MATLAB代码

% Parameters
M = 16;
N = 4; % Number of circles for CQAM
SNR_dB = 0:2:25; % Extended SNR range to reach higher values
num_symbols = 1e5; % Total number of symbols for simulation
​
% Define the constellation points for 16-PAM, 16-PSK, and 16-QAM
constellation_pam = pammod(0:(M-1), M, 0, 'gray');  % 16-PAM
constellation_psk = pskmod(0:(M-1), M, 0, 'gray');  % 16-PSK
constellation_qam = qammod(0:(M-1), M, 'gray');     % 16-QAM
constellation_cqam = generate_16CQAM(M, N);
​
​
% Calculate SEP vs SNR for each modulation scheme
SEP_vs_SNR = zeros(length(SNR_dB), 4);
​
for i = 1:length(SNR_dB)
    SNR_linear = 10^(SNR_dB(i)/10);
    SEP_vs_SNR(i, 1) = sep_pam(SNR_linear);
    SEP_vs_SNR(i, 2) = sep_psk(SNR_linear);
    SEP_vs_SNR(i, 3) = sep_qam(SNR_linear);
    SEP_vs_SNR(i, 4) = sep_cqam(SNR_linear);
end
​
​
​
% Plot SEP vs SNR
figure;
semilogy(SNR_dB, SEP_vs_SNR(:, 1), 'r-', 'LineWidth', 2);
hold on;
semilogy(SNR_dB, SEP_vs_SNR(:, 2), 'g--', 'LineWidth', 2);
semilogy(SNR_dB, SEP_vs_SNR(:, 3), 'b-.', 'LineWidth', 2);
semilogy(SNR_dB, SEP_vs_SNR(:, 4), 'm:', 'LineWidth', 2);
xlabel('E_b/N_0 (dB)');
ylabel('SEP');
legend('16-PAM', '16-PSK', '16-QAM', '16-CQAM');
title('SEP vs SNR for 16-PAM, 16-PSK, 16-QAM, and 16-CQAM');
grid on;
​
​
​
​
% Function to generate 16-CQAM constellation
function constellation = generate_16CQAM(M, N)
    n = M / N;
    R1 = 1; % Smallest radius
    constellation = [];
    for i = 1:N
        Ri = R1 * i; % Increasing radius for each circle
        base_angle = pi/4 * (i-1); % Rotate each circle by 45 degrees more than the previous one
        angles = (0:n-1) * (2 * pi / n) + base_angle;
        constellation = [constellation, Ri * exp(1i * angles)]; % Append points
    end
end
​
% Function to calculate theoretical SEP for 16-PAM
function sep = sep_pam(SNR)
    M = 16;
    sep = 2 * (M - 1) / M * qfunc(sqrt(6 * SNR / (M^2 - 1)));
end
​
% Function to calculate theoretical SEP for 16-PSK
function sep = sep_psk(SNR)
    M = 16;
    sep = 2 * qfunc(sqrt(2 * SNR) * sin(pi / M));
end
​
% Function to calculate theoretical SEP for 16-QAM
function sep = sep_qam(SNR)
    M = 16;
    Ps_M = 2 * (1 - 1/sqrt(M)) * qfunc(sqrt(3 / (M - 1) * SNR));
    sep = 1 - (1 - Ps_M)^2;
end
​
% Function to calculate theoretical SEP for 16-CQAM
function sep = sep_cqam(SNR)
    M = 16;
    R = 1; % Smallest radius
    N = 4;
    sum_sq = sum((1:N).^2);
    Es_avg = (4 / M) * R^2 * sum_sq;
    SNR_eff = SNR * Es_avg;
    sep = 2 * (M - 1) / M * qfunc(sqrt(6 * SNR_eff / (M^2 - 1)));
end

% Define the constellation points for 16-PAM, 16-PSK, and 16-QAM
M = 16;
constellation_pam = pammod(0:(M-1), M, 0, 'gray');  % 16-PAM
constellation_psk = pskmod(0:(M-1), M, 0, 'gray');  % 16-PSK
constellation_qam = qammod(0:(M-1), M, 'gray');  % 16-QAM
​
% Normalize constellations so that average symbol energy Es = 1
constellation_pam = normalize_energy(constellation_pam);
constellation_psk = normalize_energy(constellation_psk);
constellation_qam = normalize_energy(constellation_qam);
​
% Generate 16-CQAM constellation
N = 4; % Change N to other values as needed
constellation_cqam = generate_16CQAM(M, N);
constellation_cqam = normalize_energy(constellation_cqam);
​
% Calculate Euclidean distances for each modulation
dmin_pam = calculate_dmin(constellation_pam);
dmin_psk = calculate_dmin(constellation_psk);  % Calculate dmin for 16-PSK
dmin_qam = calculate_dmin(constellation_qam);
dmin_cqam = calculate_dmin(constellation_cqam);
​
% Calculate PAPR for each modulation
calculatePAPR = @(signal) max(abs(signal).^2) / mean(abs(signal).^2);
papr_pam = calculatePAPR(constellation_pam);
papr_psk = calculatePAPR(constellation_psk);
papr_qam = calculatePAPR(constellation_qam);
papr_cqam = calculatePAPR(constellation_cqam);
​
% Plot PAPR vs dmin
figure;
hold on;
scatter(dmin_pam, papr_pam, 'filled', 'DisplayName', '16-PAM');
scatter(dmin_psk, papr_psk, 'filled', 'DisplayName', '16-PSK');
scatter(dmin_qam, papr_qam, 'filled', 'DisplayName', '16-QAM');
scatter(dmin_cqam, papr_cqam, 'filled', 'DisplayName', '16-CQAM');
xlabel('d_{min}');
ylabel('PAPR');
title('Comparison of PAPR and Euclidean Distance / d_{min} in 16-PAM, 16-PSK, 16-QAM, and 16-CQAM');
legend('Location', 'best');
grid on;
hold off;
​
% Plot constellations
figure;
scatter(real(constellation_pam), imag(constellation_pam), 'filled', 'o');
title('16-PAM Constellation');
grid on;
axis equal;
for i = 1:length(constellation_pam)
    text(real(constellation_pam(i)), imag(constellation_pam(i)), ['s', num2str(i)], 'VerticalAlignment','bottom', 'HorizontalAlignment','right')
end
​
figure;
scatter(real(constellation_psk), imag(constellation_psk), 'filled', 'o');
title('16-PSK Constellation');
grid on;
axis equal;
for i = 1:length(constellation_psk)
    text(real(constellation_psk(i)), imag(constellation_psk(i)), ['s', num2str(i)], 'VerticalAlignment','bottom', 'HorizontalAlignment','right')
end
​
figure;
scatter(real(constellation_qam), imag(constellation_qam), 'filled', 'o');
title('16-QAM Constellation');
grid on;
axis equal;
for i = 1:length(constellation_qam)
    text(real(constellation_qam(i)), imag(constellation_qam(i)), ['s', num2str(i)], 'VerticalAlignment','bottom', 'HorizontalAlignment','right')
end
​
figure;
scatter(real(constellation_cqam), imag(constellation_cqam), 'filled', 'o');
title('16-CQAM Constellation');
grid on;
axis equal;
for i = 1:length(constellation_cqam)
    text(real(constellation_cqam(i)), imag(constellation_cqam(i)), ['s', num2str(i)], 'VerticalAlignment','bottom', 'HorizontalAlignment','right')
end
​
% Function to generate 16-CQAM constellation
function constellation = generate_16CQAM(M, N)
    n = M / N;
    R1 = 1; % Smallest radius
    constellation = [];
    for i = 1:N
        Ri = R1 * i; % Increasing radius for each circle
        base_angle = pi/4 * (i-1); % Rotate each circle by 45 degrees more than the previous one
        angles = (0:n-1) * (2 * pi / n) + base_angle;
        constellation = [constellation, Ri * exp(1i * angles)]; % Append points
    end
end
​
% Function to calculate the minimum distance (dmin) in a constellation
function dmin = calculate_dmin(constellation)
    num_points = length(constellation);
    distances = inf(num_points*(num_points-1)/2, 1);
    k = 1;
    for i = 1:num_points
        for j = i+1:num_points
            distances(k) = abs(constellation(i) - constellation(j));
            k = k + 1;
        end
    end
    dmin = min(distances);
end
​
% Function to normalize the average symbol energy to 1
function normalized_constellation = normalize_energy(constellation)
    avg_energy = mean(abs(constellation).^2);
    normalized_constellation = constellation / sqrt(avg_energy);
end
​

MATLAB仿真结果