clear; % 清除工作区变量
clc; % 清除命令行窗口
close all; % 关闭所有图形窗口
%% 参数设置(保持不变)
fc=10e9;
lightspeed=299792458;
lambda=lightspeed/fc;
pd = 0.9;
pfa = 1e-6;
trcs = 1.5;
%% 创建阵列(保持不变)
Array = phased.URA('Size',[25 25],...
'Lattice','Rectangular','ArrayNormal','x');
Array.ElementSpacing = [0.5 0.5]*lambda;
% rwind = hamming(25);
% cwind = hamming(25);
% taper = rwind*cwind.';
% Array.Taper = taper.';
arraygain = phased.ArrayGain('SensorArray',Array,'PropagationSpeed',physconst('LightSpeed'));
ag = arraygain(fc,[0;0]);% 天线法向增益
Elem = phased.IsotropicAntennaElement;
Elem.BackBaffled = true;
Elem.FrequencyRange = [9.99e9 10.01e9];
Array.Element = Elem;
figure;
pattern(Array.Element,fc,'PropagationSpeed',physconst('LightSpeed'),...
'Type','powerdb');% 单个天线基本阵子方向图
figure;
pattern(Array,fc,'PropagationSpeed',physconst('LightSpeed'),...
'Type','powerdb');% 阵列天线方向图
%% 获取方向图数据
az = -60:0.01:60; % 方位角范围
el = 0; % 固定俯仰角
P = pattern(Array, fc, az, el, 'Type', 'powerdb', ...
'PropagationSpeed', physconst('LightSpeed'));
P = P(:); % 确保为列向量
%% 核心:寻找所有峰值(主瓣+旁瓣)
% 1. 找到主瓣峰值
[main_peak_val, main_peak_idx] = max(P);
main_peak_az = az(main_peak_idx);
% 2. 确定主瓣区域(根据3dB波束宽度)
half_power = main_peak_val - 3;
left_idx = find(P(1:main_peak_idx) <= half_power, 1, 'last');
right_idx = main_peak_idx + find(P(main_peak_idx:end) <= half_power, 1, 'first') - 1;
beamwidth = az(right_idx) - az(left_idx);
mainlobe_start_idx2 = max(1, main_peak_idx - round(beamwidth/0.01/2));
mainlobe_end_idx2 = min(length(az), main_peak_idx + round(beamwidth/0.01/2));
% 3. 排除主瓣区域
mainlobe_start_idx = max(1, main_peak_idx - round(3.4*beamwidth/0.01/2));
mainlobe_end_idx = min(length(az), main_peak_idx + round(3.4*beamwidth/0.01/2));
% 4. 标记非主瓣区域
sidelobe_mask = true(size(P));
sidelobe_mask(mainlobe_start_idx:mainlobe_end_idx) = false;
% 5. 在非主瓣区域寻找所有旁瓣峰值
sidelobe_data = P(sidelobe_mask);
sidelobe_az = az(sidelobe_mask);
% 使用findpeaks函数寻找旁瓣峰值
[peak_values, peak_indices] = findpeaks(sidelobe_data, ...
'MinPeakDistance', 50,... % 至少50个点(0.5度)的间隔
'MinPeakHeight', main_peak_val - 60); % 只找高于主瓣-60dB的峰值
% 转换回原始索引
sidelobe_peak_indices = find(sidelobe_mask);
sidelobe_peak_indices = sidelobe_peak_indices(peak_indices);
sidelobe_peak_az = az(sidelobe_peak_indices);
% 6. 确定关键旁瓣参数
if ~isempty(peak_values)
% 最高旁瓣(SLL)
[max_sll_val, max_sll_idx] = max(peak_values);
max_sll_az = sidelobe_peak_az(max_sll_idx);
max_sll_db = max_sll_val - main_peak_val; % 相对于主瓣的dB值
% 第一旁瓣(最靠近主瓣的旁瓣)
[~, first_sll_idx] = min(abs(sidelobe_peak_az - main_peak_az));
first_sll_val = peak_values(first_sll_idx);
first_sll_az = sidelobe_peak_az(first_sll_idx);
first_sll_db = first_sll_val - main_peak_val;
else
warning('未找到有效旁瓣峰值,请调整角度范围或参数');
max_sll_db = -inf;
first_sll_db = -inf;
end
%% 绘图显示阵列方向图结果
figure;
plot(az, P, 'b', 'LineWidth', 1.2);
hold on;
% 添加-3dB水平线(红色虚线)
half_power_line = repmat(half_power, size(az));
plot(az, half_power_line, 'r--', 'LineWidth', 1.5);
% 标记主瓣峰值
plot(main_peak_az, main_peak_val, 'ro', 'MarkerSize', 2, 'MarkerFaceColor', 'r');
text(main_peak_az+1, main_peak_val, '主瓣峰值', 'Color', 'r');
% 标记主瓣区域
fill([az(mainlobe_start_idx2) az(mainlobe_end_idx2) az(mainlobe_end_idx2) az(mainlobe_start_idx2)], ...
[min(P) min(P) max(P) max(P)], 'y', 'FaceAlpha', 0.2);
text(az(main_peak_idx), min(P)+5, '主瓣区域', 'HorizontalAlignment', 'center');
% 标记所有旁瓣峰值
if ~isempty(peak_values)
plot(sidelobe_peak_az, peak_values, 'go', 'MarkerSize', 2, 'MarkerFaceColor', 'g');
% 标记最高旁瓣
plot(max_sll_az, max_sll_val, 'mo', 'MarkerSize', 2, 'MarkerFaceColor', 'm');
text(max_sll_az+1, max_sll_val, sprintf('最高旁瓣: %.1fdB', max_sll_db), 'Color', 'm');
end
%% 添加-3dB点和波束宽度标注(优化布局)
% 标记-3dB点(红色五角星)
plot(az(left_idx), half_power, 'rp', 'MarkerSize', 2, 'MarkerFaceColor', 'r');
plot(az(right_idx), half_power, 'rp', 'MarkerSize', 2, 'MarkerFaceColor', 'r');
% 添加-3dB线文本(移到主瓣中央上方)
text(main_peak_az-50, half_power, '-3dB', 'Color', 'r', ...
'VerticalAlignment', 'bottom', 'HorizontalAlignment', 'center');
% 绘制垂线(红色点线)
line([az(left_idx), az(left_idx)], [half_power, main_peak_val-75], ...
'LineStyle', ':', 'Color', 'r', 'LineWidth', 1.5);
line([az(right_idx), az(right_idx)], [half_power, main_peak_val-75], ...
'LineStyle', ':', 'Color', 'r', 'LineWidth', 1.5);
% 标记垂足(红色圆点)
plot(az(left_idx), main_peak_val-70, 'ro', 'MarkerSize', 2, 'MarkerFaceColor', 'r');
plot(az(right_idx), main_peak_val-70, 'ro', 'MarkerSize', 2, 'MarkerFaceColor', 'r');
% 添加垂足数值标注(红色) - 调整位置避免重叠
text(az(left_idx)-2, main_peak_val-71, sprintf('%.1f°', az(left_idx)), ...
'HorizontalAlignment', 'center', 'Color', 'b', 'FontSize', 10); % 向左移动
text(az(right_idx)+2, main_peak_val-71, sprintf('%.1f°', az(right_idx)), ...
'HorizontalAlignment', 'center', 'Color', 'b', 'FontSize', 10); % 向右移动
% 添加波束宽度标注(红色双箭头+文字) - 下移避免重叠
arrow_y = main_peak_val - 18; % 从-15下移到-18
text(mean([az(left_idx), az(right_idx)])-5, arrow_y+0.5, sprintf('BW=%.2f°', beamwidth), ...
'HorizontalAlignment', 'center', 'Color', 'g', 'FontWeight', 'bold');
% 调整Y轴范围保证标注可见
ylim([main_peak_val-80 main_peak_val+5]);
grid on;
xlabel('方位角(度)');
ylabel('功率(dB)');
title('阵列方向图及旁瓣峰值标记');
xlim([min(az) max(az)]);
%% 输出阵列参数结果
fprintf('主瓣峰值: %.2f dB\n', main_peak_val);
fprintf('波束宽度: %.2f 度\n', beamwidth);
fprintf('最高旁瓣电平(SLL): %.2f dBc\n', max_sll_db); % dBc表示相对于载波(主瓣)
if ~isempty(peak_values)
fprintf('第一旁瓣电平: %.2f dBc\n', first_sll_db);
end
%% 设计线性调频信号
Waveform = phased.LinearFMWaveform('SampleRate',20e+06,...
'PulseWidth',5e-6,'PRF',5000,'SweepBandwidth',10e6,...
'SweepDirection','Up','SweepInterval','Positive'...
,'Envelope','Rectangular','NumPulses',1,'FrequencyOffset',10e9);
Fs = Waveform.SampleRate;
prf = Waveform.PRF;
NumOfPulses = 6; % 用于多普勒分析的脉冲数
tb_product = Waveform.SweepBandwidth * Waveform.PulseWidth;% 时宽带宽积
delta_R = lightspeed / (2*Waveform.SweepBandwidth);% 距离分辨率
delta_v = (lambda * prf) / (2*NumOfPulses);% 速度分辨率
compression_gain = 10*log10(tb_product); % 脉冲压缩增益(dB)
snr_min_raw = albersheim(pd, pfa, NumOfPulses); % 未考虑压缩增益的SNR
snr_min_actual = snr_min_raw - compression_gain; % 实际所需SNR
% 生成基带信号
x = Waveform();
% 匹配滤波
coeff = getMatchedFilter(Waveform);
Compression = phased.MatchedFilter('Coefficients',coeff(:,1), ...
'SpectrumWindow','None');
y = Compression(x);
%% 设计发射机和接收机
transmitter=phased.Transmitter('PeakPower',20e3,...
'Gain',20,'InUseOutputPort',true,...
'CoherentOnTransmit',true);
receiver=phased.ReceiverPreamp(...
"Gain",20,'NoiseFigure',7.0,...
'SampleRate',Fs,'EnableInputPort',true);
% 噪声参数计算
F_dB = 7; % 噪声系数
T0 = 290; % 参考温度 (K)
L_dB = 10; % 系统损失
k = 1.38e-23; % 玻尔兹曼常数
max_range = 20e3; % 最大探测距离
P_noise = k * T0 * Waveform.SweepBandwidth * db2pow(F_dB); % 噪声功率
% 峰值功率和最大作用距离计算
peak_power = ((4*pi)^3*k*T0*Waveform.SweepBandwidth*db2pow(F_dB)*max_range^4*...
db2pow(snr_min_raw))/(db2pow(2*(ag+transmitter.Gain))*trcs*lambda^2);
Rmax=((transmitter.PeakPower * db2pow(ag+transmitter.Gain) * db2pow(ag+receiver.Gain)...
* lambda^2 * trcs) ./ (P_noise ...
* (4 * pi)^3 * db2pow(snr_min_raw))).^(1/4);
%% 目标和传感器运动设置
initialAz = 60;
endAz = -60;
scanstep = -2; % 扫描步进角度
scangrid = initialAz+scanstep/2:scanstep:endAz;
numscans = length(scangrid);
pulsenum = NumOfPulses*numscans;
% 目标初始位置和速度
tgtpos = [[8000; 6000; 0],[8020; 6000; 0],[8000; 6020; 0],...
[12000; 5000; 0],[12030; 5000; 0],[12000; 5030; 0]];
tgtvel = [[10; 0; 0],[10; 0; 0],[10; 0; 0],...
[-30; 0; 0],[-30; 0; 0],[-30; 0; 0]];% 目标真实径向速度
tgtmotion = phased.Platform('InitialPosition',tgtpos,'Velocity',tgtvel);
tgtrcs = [1.5 1.5 1.5 1.5 1.5 1.5];% 目标RCS
target = phased.RadarTarget('MeanRCS',tgtrcs,'OperatingFrequency',fc);
% 传感器位置(静止)
sensorpos = [0; 0; 0];
sensorvel = [0; 0; 0];
sensormotion = phased.Platform('InitialPosition',sensorpos,'Velocity',sensorvel);
% 初始目标距离和角度
[tgtrng,tgtang] = rangeangle(tgtmotion.InitialPosition,...
sensormotion.InitialPosition);
numtargets = length(target.MeanRCS);
% 导向矢量和波束形成器
steeringvec = phased.SteeringVector('SensorArray',Array,...
'PropagationSpeed',lightspeed);
beamformer = phased.PhaseShiftBeamformer('SensorArray',Array,...
'OperatingFrequency',fc,'PropagationSpeed',lightspeed,...
'DirectionSource','Input port');
% 传播信道
channel = phased.FreeSpace(...
'SampleRate',Fs,...
'TwoWayPropagation',true,...
'OperatingFrequency',fc);
% 快时间网格和距离轴
fast_time_grid = unigrid(0, 1/Fs, 1/prf, '[)');
range_axis = lightspeed * fast_time_grid / 2; % 距离轴
% 接收脉冲存储
waveform_samples = length(Waveform());
rxpulses = zeros(waveform_samples, pulsenum, 'double');
% 辐射器和收集器
radiator = phased.Radiator(...
'Sensor',Array,...
'PropagationSpeed', lightspeed,...
'OperatingFrequency', fc,...
'CombineRadiatedSignals', true,...
'WeightsInputPort',true);
collector = phased.Collector(...
'Sensor',Array,...
'PropagationSpeed', lightspeed,...
'OperatingFrequency', fc,...
'Wavefront', 'Plane',...
'WeightsInputPort',false);
%% 回波模拟
for m = 1:pulsenum
fprintf('脉冲数: %d \n', m);
% 更新位置
[sensorpos,sensorvel] = sensormotion(1/prf);
[tgtpos,tgtvel] = tgtmotion(1/prf);
% 目标角度计算
[tgtrng,tgtang] = rangeangle(tgtpos,sensorpos);
% 扫描角度和导向矢量
scanid = floor((m-1)/NumOfPulses) + 1;
sv = steeringvec(fc,scangrid(scanid));
w = conj(sv);
% 发射信号和传播
pulse = Waveform();
%% ===== 计算基带信号功率 =====
% 1. 计算瞬时功率(每个采样点)
instant_power = abs(pulse).^2; % |x|²
% 2. 计算平均功率(整个脉冲)
avg_power_baseband = mean(instant_power);
% 3. 计算峰值功率(脉冲内最大瞬时功率)
peak_power_baseband = max(instant_power);
% 4. 计算脉冲能量
pulse_energy = sum(instant_power) / Waveform.SampleRate;
%% ===== 显示结果 =====
fprintf('基带信号功率分析:\n');
fprintf(' 平均功率: %.6f W\n', avg_power_baseband);
fprintf(' 峰值功率: %.6f W\n', peak_power_baseband);
fprintf(' 脉冲能量: %.6f J\n', pulse_energy);
fprintf(' 信号长度: %d 采样点\n', length(pulse));
fprintf(' 采样率: %.0f Hz\n', Waveform.SampleRate);
fprintf(' 脉冲宽度: %.6f s\n', Waveform.PulseWidth);
figure;
subplot(2,1,1);
plot(real(pulse), 'b');
title('实部 (Re)');
xlabel('样本索引');
ylabel('幅度');
subplot(2,1,2);
plot(imag(pulse), 'r');
title('虚部 (Im)');
xlabel('样本索引');
ylabel('幅度');
[txsig,txstatus] = transmitter(pulse);
% 1. 理论功率计算(包含发射机增益)
pulse_peak_power = transmitter.PeakPower; % 20 kW
duty_cycle = Waveform.PulseWidth * Waveform.PRF;
avg_power_theoretical = pulse_peak_power * duty_cycle * db2pow(transmitter.Gain);
% 2. 仿真功率计算(使用发射后的信号)
% 计算信号的实际功率(考虑发射机放大)
pulse_power_instantaneous = abs(txsig).^2;
avg_power_simulated = mean(pulse_power_instantaneous);
% 3. 峰值功率计算
% 理论峰值功率(包含发射机增益)
peak_power_theoretical = transmitter.PeakPower * db2pow(transmitter.Gain);
% 仿真峰值功率(取瞬时功率最大值)
peak_power_simulated = max(pulse_power_instantaneous);
% 4. 结果对比
fprintf('脉冲 %d 功率分析:\n', m);
fprintf(' 理论峰值功率 = %.2f W (%.2f kW)\n', ...
peak_power_theoretical, peak_power_theoretical/1000);
fprintf(' 仿真峰值功率 = %.2f W (%.2f kW)\n', ...
peak_power_simulated, peak_power_simulated/1000);
fprintf(' 峰值功率相对误差 = %.4f%%\n', ...
100*abs(peak_power_theoretical-peak_power_simulated)/peak_power_theoretical);
fprintf(' 理论平均功率 = %.2f W\n', avg_power_theoretical);
fprintf(' 仿真平均功率 = %.2f W\n', avg_power_simulated);
fprintf(' 平均功率相对误差 = %.4f%%\n\n', ...
100*abs(avg_power_theoretical-avg_power_simulated)/avg_power_theoretical);
% 绘制发射信号 txsig 的波形
figure(11);
set(gcf, 'Position', [100, 100, 1000, 600]);
% 绘制实部和虚部
subplot(3,1,1);
plot(real(txsig), 'b');
hold on;
plot(imag(txsig), 'r');
hold off;
title('发射信号:实部(蓝)与虚部(红)');
xlabel('样本索引');
ylabel('幅度');
grid on;
legend('实部(I)', '虚部(Q)');
% 绘制幅度
subplot(3,1,2);
plot(abs(txsig), 'g', 'LineWidth', 1.5);
title('发射信号幅度');
xlabel('样本索引');
ylabel('幅度');
grid on;
% 绘制相位
subplot(3,1,3);
plot(angle(txsig), 'm', 'LineWidth', 1.2);
title('发射信号相位');
xlabel('样本索引');
ylabel('相位(弧度)');
grid on;
% 添加总标题
sgtitle(sprintf('发射信号波形 - 脉冲 %d (扫描角度 %.1f°)', m, scangrid(scanid)));
txsig = radiator(txsig,tgtang,w);
txsig = channel(txsig,sensorpos,tgtpos,sensorvel,tgtvel);
% 目标反射和接收处理
tgtsig = target(txsig);
rxsig = collector(tgtsig,tgtang);
rxsig = receiver(rxsig,~(txstatus>0));
rxpulses(:,m) = beamformer(rxsig,[scangrid(scanid);0]);
end
最新发布