Initial ARSS. Need to improve and add functions.

This commit is contained in:
2026-03-02 14:48:21 +09:00
commit 775668afdb
20 changed files with 1610 additions and 0 deletions
+31
View File
@@ -0,0 +1,31 @@
% chirpIdx , t와 fs
function tx_matrix = apply_mimo_mode(x_amp, t, Timing, NumTx, mimoMode, fs)
N_samples = length(x_amp);
tx_matrix = zeros(NumTx, N_samples);
T_chirp = Timing.IdleTime + Timing.RampEndTime;
SamplesPerChirp = round(T_chirp * fs);
% []: (1, 2, 3...)
chirp_indices = floor((0:N_samples-1) / SamplesPerChirp) + 1;
switch upper(mimoMode)
case 'TDM'
% 1,2,3,1,2,3
active_tx_seq = mod(chirp_indices - 1, NumTx) + 1;
for tx = 1:NumTx
%
tx_matrix(tx, :) = x_amp .* (active_tx_seq == tx);
end
case 'DDMA'
%
for tx = 1:NumTx
phase_shift_seq = 2 * pi * (tx - 1) * (chirp_indices - 1) / NumTx;
tx_matrix(tx, :) = x_amp .* exp(1j * phase_shift_seq);
end
otherwise
error(' MIMO . TDM DDMA를 .');
end
end
+26
View File
@@ -0,0 +1,26 @@
function [t, tx_mask] = generate_waveform_timing(RadarParams)
% RadarParams TX
% TX .
%
% :
% - RadarParams:
%
% :
% - t: ()
% - tx_mask: TX (0 1)
Timing = RadarParams.Waveform.Timing;
NumChirps = RadarParams.Waveform.NumChirps;
T_chirp = Timing.IdleTime + Timing.RampEndTime;
T_frame = T_chirp * NumChirps; %
fs = RadarParams.Waveform.fs_waveform;
N_samples = round(T_frame * fs);
% ( )
t = linspace(0, T_frame, N_samples);
% (t_mod) ( )
t_mod = mod(t, T_chirp);
% TX
tx_mask = zeros(1, N_samples);
tx_start_abs = Timing.IdleTime + Timing.TxStartTime;
tx_mask(t_mod >= tx_start_abs) = 1;
end
+25
View File
@@ -0,0 +1,25 @@
function TxOut = radiate_antenna(Target, TxPatternArray)
% : , TX
% : TxOut.G_tx_amp (, TX별 [NumTargets x NumTx])
NumTx = length(TxPatternArray);
NumTargets = Target.NumTargets;
TxOut.G_tx_amp = zeros(NumTargets, NumTx);
for k = 1:NumTargets
for tx = 1:NumTx
pat = TxPatternArray(tx);
%
if strcmpi(pat.Type, '2D')
g_dBi = interp2(pat.az_angles, pat.el_angles, pat.gain_dBi, Target.az(k), Target.el(k), 'linear', -20);
else
g_az = interp1(pat.az_angles, pat.gain_az_dBi, Target.az(k), 'linear', -20);
g_el = interp1(pat.el_angles, pat.gain_el_dBi, Target.el(k), 'linear', -20);
g_dBi = g_az + g_el - pat.max_gain_dBi;
end
% (sqrt(10^(G_dBi/10)))
TxOut.G_tx_amp(k, tx) = sqrt(10^(g_dBi/10));
end
end
end
+145
View File
@@ -0,0 +1,145 @@
function fig = visualize_multi_tx_waveform(t, Timing, fc, f_start, Slope, tx_mask, TotalNumChirps, mimoMode, NumTx)
fig = figure('Name', 'Multi-Chirp & MIMO Modulation', 'Position', [150, 150, 1100, 750]);
T_chirp = Timing.IdleTime + Timing.RampEndTime;
NumChirpsToPlot = min(TotalNumChirps, 4);
plot_idx = (t <= NumChirpsToPlot * T_chirp);
t_plot = t(plot_idx);
tx_mask_plot = tx_mask(plot_idx);
t_mod = mod(t_plot, T_chirp);
t_ramp = t_mod - Timing.IdleTime;
inst_freq_theoretical = f_start + Slope * t_ramp;
inst_freq_masked = inst_freq_theoretical;
inst_freq_masked(tx_mask_plot == 0) = NaN;
% =========================================================
% [Subplot 1] :
% =========================================================
subplot(2, 1, 1);
y_min = (f_start / 1e9) - 0.25;
y_max = (f_start / 1e9) + (Slope * Timing.RampEndTime / 1e9) + 0.1;
plot(t_plot * 1e6, inst_freq_masked / 1e9, 'b', 'LineWidth', 2);
grid on; hold on;
title(sprintf('Continuous Multi-Chirp Sequence (Showing %d of %d Chirps)', NumChirpsToPlot, TotalNumChirps), 'FontSize', 12);
ylabel('Absolute Frequency (GHz)', 'FontSize', 11);
ylim([y_min, y_max]);
if ~isempty(t_plot)
xlim([0, max(t_plot)*1e6]);
else
xlim([0, 1]); % fallback range when no data
end
f_valid_start = f_start + (Slope * Timing.AdcStartTime);
f_valid_end = f_start + (Slope * (Timing.AdcStartTime + Timing.AdcSampTime));
valid_bandwidth = f_valid_end - f_valid_start;
info_str = {
sprintf(' Valid Start Freq : %.4f GHz', f_valid_start / 1e9), ...
sprintf(' Center Freq : %.4f GHz', fc / 1e9), ...
sprintf(' Valid End Freq : %.4f GHz', f_valid_end / 1e9), ...
sprintf(' Transmit Bandwidth : %.4f GHz', valid_bandwidth / 1e9), ...
sprintf(' Total Gen Chirps : %d', TotalNumChirps)
};
text(0.02, 0.96, info_str, 'Units', 'normalized', 'FontSize', 10, 'FontWeight', 'bold', 'BackgroundColor', [1 1 1 0.85], 'EdgeColor', 'k', 'VerticalAlignment', 'top', 'Margin', 5);
guide_line_args = {'Color', [0 0 0.5], 'LineStyle', '--', 'LineWidth', 1};
for i = 0:NumChirpsToPlot
plot([i * T_chirp * 1e6, i * T_chirp * 1e6], [y_min, y_max], guide_line_args{:});
end
draw_dim_arrow = @(x1, x2, y, label_str) ...
[plot([x1, x2], [y, y], 'k-', 'LineWidth', 1.2), ...
fill([x1, x1 + min(0.6, (x2-x1)*0.35), x1 + min(0.6, (x2-x1)*0.35)], [y, y + 0.015, y - 0.015], 'k', 'EdgeColor', 'none'), ...
fill([x2, x2 - min(0.6, (x2-x1)*0.35), x2 - min(0.6, (x2-x1)*0.35)], [y, y + 0.015, y - 0.015], 'k', 'EdgeColor', 'none'), ...
text((x1+x2)/2, y + 0.03, label_str, 'HorizontalAlignment', 'center', 'VerticalAlignment', 'bottom', 'FontSize', 9, 'FontWeight', 'bold', 'BackgroundColor', 'w', 'EdgeColor', 'k')];
y_pri = y_min + 0.08;
pri_us = T_chirp * 1e6;
prf_khz = (1 / T_chirp) / 1e3;
pri_prf_label = sprintf('PRI: %.1f \\mus\nPRF: %.1f kHz', pri_us, prf_khz);
draw_dim_arrow(0, pri_us, y_pri, pri_prf_label);
hold off;
% =========================================================
% [Subplot 2] : MIMO Active TX & Phase
% =========================================================
subplot(2, 1, 2);
hold on; grid on;
% --- [ ]: TX ---
% lines() .
tx_colors = lines(NumTx);
for tx = 1:NumTx
plot([0, NumChirpsToPlot * T_chirp * 1e6], [tx, tx], ':', 'Color', [0.8 0.8 0.8], 'HandleVisibility', 'off');
end
for i = 0 : NumChirpsToPlot - 1
t_center_us = (i + 0.5) * T_chirp * 1e6;
if strcmpi(mimoMode, 'TDM')
active_tx_list = mod(i, NumTx) + 1;
elseif strcmpi(mimoMode, 'DDMA')
active_tx_list = 1:NumTx;
else
active_tx_list = 1;
end
for tx = 1:NumTx
% TX
current_color = tx_colors(tx, :);
if ismember(tx, active_tx_list)
if strcmpi(mimoMode, 'DDMA')
phase_rad = 2 * pi * (tx - 1) * i / NumTx;
phase_deg = mod(rad2deg(phase_rad), 360);
else
phase_deg = 0;
end
% :
plot(t_center_us, tx, 'o', 'MarkerSize', 12, 'MarkerFaceColor', current_color, 'MarkerEdgeColor', current_color, 'HandleVisibility', 'off');
text_offset_us = T_chirp * 1e6 * 0.08;
text(t_center_us + text_offset_us, tx, sprintf('%d^\\circ', round(phase_deg)), ...
'VerticalAlignment', 'middle', 'HorizontalAlignment', 'left', ...
'FontSize', 9, 'FontWeight', 'bold', 'Color', current_color);
else
% :
plot(t_center_us, tx, 'o', 'MarkerSize', 12, 'MarkerFaceColor', 'w', 'MarkerEdgeColor', current_color, 'HandleVisibility', 'off');
end
end
end
title(sprintf('%s Active TX Antenna & Phase Map', upper(mimoMode)), 'FontSize', 12);
ylabel('TX Antenna', 'FontSize', 11);
yticks(1:NumTx);
% Y축 ( )
% yticklabels(arrayfun(@(x) sprintf('\\color[rgb]{%f,%f,%f}TX %d', tx_colors(x,1), tx_colors(x,2), tx_colors(x,3), x), 1:NumTx, 'UniformOutput', false));
yticklabels(arrayfun(@(x) sprintf('TX %d', x), 1:NumTx, 'UniformOutput', false));
ylim([0.5, NumTx + 0.5]);
% --- [ ]: ( ) ---
neutral_color = [0.3 0.3 0.3];
h_on = plot(NaN, NaN, 'o', 'MarkerSize', 10, 'MarkerFaceColor', neutral_color, 'MarkerEdgeColor', neutral_color, 'DisplayName', 'Active');
h_off = plot(NaN, NaN, 'o', 'MarkerSize', 10, 'MarkerFaceColor', 'w', 'MarkerEdgeColor', neutral_color, 'DisplayName', 'Inactive');
legend([h_on, h_off], 'Location', 'best', 'FontSize', 10, 'Box', 'on');
xlabel('Time (\mus)', 'FontSize', 11);
xlim([0, max(t_plot)*1e6]);
y_lims_bottom = ylim;
for i = 0:NumChirpsToPlot
plot([i * T_chirp * 1e6, i * T_chirp * 1e6], y_lims_bottom, guide_line_args{:}, 'HandleVisibility', 'off');
end
hold off;
end
+78
View File
@@ -0,0 +1,78 @@
function fig = visualize_tx_waveform(t, Timing, fc, f_start, Slope, tx_mask, peak_phase_error, f_ripple)
fig = figure('Name', 'FMCW Radar Timing & True RF Frequency', 'Position', [100, 100, 1100, 700]);
t_ramp = t - Timing.IdleTime;
%
inst_freq_theoretical = f_start + Slope * t_ramp;
freq_nonlin = peak_phase_error * f_ripple * cos(2 * pi * f_ripple * t_ramp);
inst_freq_theoretical = inst_freq_theoretical + freq_nonlin;
inst_freq_masked = inst_freq_theoretical;
inst_freq_masked(tx_mask == 0) = NaN;
% Y축
y_min = (f_start / 1e9) - 0.2;
y_max = (f_start / 1e9) + (Slope * Timing.RampEndTime / 1e9) + 0.1;
% --- [ ] ---
plot(t * 1e6, inst_freq_masked / 1e9, 'b', 'LineWidth', 2.5);
title(sprintf('Mathematical Instantaneous RF Frequency & HW Timing (Center fc = %.2f GHz)', fc/1e9), 'FontSize', 12);
xlabel('Time (\mus)', 'FontSize', 11); ylabel('Absolute Frequency (GHz)', 'FontSize', 11);
ylim([y_min, y_max]);
grid on; hold on;
% --- [ : ] ---
% ADC /
f_valid_start = f_start + (Slope * Timing.AdcStartTime);
f_valid_end = f_start + (Slope * (Timing.AdcStartTime + Timing.AdcSampTime));
valid_bandwidth = f_valid_end - f_valid_start; %
info_str = {
sprintf(' Valid Start Freq : %.4f GHz', f_valid_start / 1e9), ...
sprintf(' Center Freq : %.4f GHz', fc / 1e9), ...
sprintf(' Valid End Freq : %.4f GHz', f_valid_end / 1e9), ...
sprintf(' Valid Bandwidth : %.4f GHz', valid_bandwidth / 1e9) %
};
text(0.02, 0.96, info_str, 'Units', 'normalized', ...
'FontSize', 10, 'FontWeight', 'bold', 'BackgroundColor', [1 1 1 0.85], ...
'EdgeColor', 'k', 'VerticalAlignment', 'top', 'Margin', 5);
% --- [ ] ---
y_lims = ylim;
guide_line_args = {'Color', [0 0 0.5], 'LineStyle', '--', 'LineWidth', 1};
x_coords = [0, Timing.IdleTime, ...
(Timing.IdleTime + Timing.TxStartTime), ...
(Timing.IdleTime + Timing.AdcStartTime), ...
(Timing.IdleTime + Timing.AdcStartTime + Timing.AdcSampTime), ...
(Timing.IdleTime + Timing.RampEndTime)];
x_coords_us = x_coords * 1e6;
for i = 1:length(x_coords_us)
plot([x_coords_us(i), x_coords_us(i)], y_lims, guide_line_args{:});
end
% --- [ ] ---
draw_dim_arrow = @(x1, x2, y, name, val_us) ...
[plot([x1, x2], [y, y], 'k-', 'LineWidth', 1.2), ...
fill([x1, x1 + min(0.6, (x2-x1)*0.35), x1 + min(0.6, (x2-x1)*0.35)], [y, y + 0.015, y - 0.015], 'k', 'EdgeColor', 'none'), ...
fill([x2, x2 - min(0.6, (x2-x1)*0.35), x2 - min(0.6, (x2-x1)*0.35)], [y, y + 0.015, y - 0.015], 'k', 'EdgeColor', 'none'), ...
text((x1+x2)/2, y + 0.03, sprintf('%s\n(%.1f \\mus)', name, val_us), 'HorizontalAlignment', 'center', 'VerticalAlignment', 'bottom', 'FontSize', 9, 'FontWeight', 'bold', 'BackgroundColor', 'w', 'EdgeColor', 'k')];
%
y1 = y_min + 0.05;
y2 = y_min + 0.18;
y3 = y_min + 0.31;
draw_dim_arrow(x_coords_us(1), x_coords_us(2), y1, 'Idle Time', Timing.IdleTime * 1e6);
draw_dim_arrow(x_coords_us(2), x_coords_us(3), y2, 'TX Start', Timing.TxStartTime * 1e6);
draw_dim_arrow(x_coords_us(2), x_coords_us(4), y3, 'ADC Delay', Timing.AdcStartTime * 1e6);
draw_dim_arrow(x_coords_us(4), x_coords_us(5), y2, 'ADC Sampling (Valid)', Timing.AdcSampTime * 1e6);
draw_dim_arrow(x_coords_us(5), x_coords_us(6), y1, 'Excess', Timing.ExcessTime * 1e6);
hold off;
end
+59
View File
@@ -0,0 +1,59 @@
function ChannelOut = apply_channel_effects(RadarParams)
Target = RadarParams.Target;
TxPos = RadarParams.Antenna.TxPos;
RxPos = RadarParams.Antenna.RxPos;
NumChirps = RadarParams.Waveform.NumChirps;
Timing = RadarParams.Waveform.Timing;
fs_adc = RadarParams.Waveform.fs_adc;
f_start = RadarParams.Waveform.f_start;
Slope = RadarParams.Waveform.Slope;
c = RadarParams.Basic.c;
NumTx = size(TxPos, 2);
NumRx = size(RxPos, 2);
NumTargets = Target.NumTargets;
T_pri = RadarParams.Waveform.PRI;
t_adc = 0 : 1/fs_adc : (Timing.AdcSampTime - 1/fs_adc);
N_adc = length(t_adc);
% []
% t_adc에 (lambda)
f_inst = f_start + Slope * t_adc;
lambda_t = c ./ f_inst;
ChannelOut.tau = zeros(NumTargets, NumRx, NumTx, NumChirps, N_adc);
ChannelOut.space_loss_amp = zeros(NumTargets, NumRx, NumTx, NumChirps, N_adc);
ref_point = (TxPos(:,1) + RxPos(:,1)) / 2;
for k = 1:NumTargets
az_rad = deg2rad(Target.az(k));
el_rad = deg2rad(Target.el(k));
rel_pos0 = Target.R(k) * [cos(el_rad)*sin(az_rad); cos(el_rad)*cos(az_rad); sin(el_rad)];
target_world_pos0 = ref_point + rel_pos0;
v_vec = Target.v(k) * (rel_pos0 / norm(rel_pos0));
sigma = 10^(Target.rcs(k)/10);
for m = 0:NumChirps-1
%
target_pos_at_chirp = target_world_pos0 + v_vec * (m * T_pri);
for tx = 1:NumTx
for rx = 1:NumRx
%
curr_target_pos = target_pos_at_chirp + v_vec * t_adc;
d_tx = sqrt(sum((curr_target_pos - TxPos(:,tx)).^2, 1));
d_rx = sqrt(sum((curr_target_pos - RxPos(:,rx)).^2, 1));
% TTD
ChannelOut.tau(k, rx, tx, m+1, :) = (d_tx + d_rx) / c;
% [ ]
denom = (4*pi)^3 * (d_tx.^2 .* d_rx.^2);
ChannelOut.space_loss_amp(k, rx, tx, m+1, :) = sqrt((lambda_t.^2 * sigma) ./ denom);
end
end
end
end
end
+39
View File
@@ -0,0 +1,39 @@
function quantized_data = apply_adc_quantization(adc_raw_data, n_bits, v_full_scale, mode)
% :
% - adc_raw_data: HPF/LPF가 (Voltage Scale)
% - n_bits: ADC
% - v_full_scale: ADC Peak-to-Peak (e.g., 2.0V)
% - mode: 'IQ' (Complex ) 'Real' (Real )
v_max = v_full_scale / 2;
lsb = v_full_scale / (2^n_bits);
if strcmpi(mode, 'IQ')
% --- Case 1: I/Q Demodulation (Complex) ---
% I채널과 Q채널 Clipping Quantization
r_part = real(adc_raw_data);
i_part = imag(adc_raw_data);
r_part = max(min(r_part, v_max), -v_max);
i_part = max(min(i_part, v_max), -v_max);
q_r = round(r_part / lsb) * lsb;
q_i = round(i_part / lsb) * lsb;
quantized_data = q_r + 1j * q_i;
elseif strcmpi(mode, 'Real')
% --- Case 2: Real Demodulation ---
% (I채널) ADC
%
r_part = real(adc_raw_data);
r_part = max(min(r_part, v_max), -v_max);
q_r = round(r_part / lsb) * lsb;
% (Real)
quantized_data = q_r;
else
error(' IQ Real .');
end
end
+39
View File
@@ -0,0 +1,39 @@
function filtered_data = apply_analog_hpf(adc_raw_data, fs_adc, fc_hpf_Hz)
% :
% - adc_raw_data: [NumRx, NumTx, NumChirps, N_samples]
% - fs_adc: ADC (Hz)
% - fc_hpf_Hz: HPF (Cut-off Frequency, Hz)
[NumRx, NumTx, NumChirps, N_samples] = size(adc_raw_data);
filtered_data = zeros(size(adc_raw_data));
% 1. 1 HPF (S-plane -> Z-plane )
% Analog Transfer Function: H(s) = s / (s + omega_c)
omega_c = 2 * pi * fc_hpf_Hz;
% (Bilinear Transform)
% [b, a] = butter(1, fc_hpf_Hz / (fs_adc/2), 'high'); %
T = 1 / fs_adc;
alpha = 2 / T;
b0 = alpha / (alpha + omega_c);
b1 = -alpha / (alpha + omega_c);
a1 = (omega_c - alpha) / (alpha + omega_c);
b = [b0, b1];
a = [1, a1];
% 2.
for rx = 1:NumRx
for tx = 1:NumTx
for m = 1:NumChirps
% 퀀
raw_sig = squeeze(adc_raw_data(rx, tx, m, :));
% ( 'filtfilt' 'filter' )
% (Causal) filter
filtered_data(rx, tx, m, :) = filter(b, a, raw_sig);
end
end
end
end
+203
View File
@@ -0,0 +1,203 @@
function [adc_raw_data, adc_combined] = apply_lna_and_mixer(RxOut, TxOut, RadarParams)
% : adc_raw_data [NumRx, NumTx, NumChirps, N_adc_samples] -
% adc_combined [NumRx, NumChirps, N_adc_samples] - ADC ( TX )
% : RxOut, TxOut, RadarParams
% MIMO TX .
Timing = RadarParams.Waveform.Timing;
fs_adc = RadarParams.Waveform.fs_adc;
f_start = RadarParams.Waveform.f_start;
Slope = RadarParams.Waveform.Slope;
pn_level = RadarParams.Waveform.nonideal.pn_level;
f_ripple = RadarParams.Waveform.nonideal.f_ripple;
peak_phase_error = RadarParams.Waveform.nonideal.peak_phase_error;
enable_phase_noise = true;
enable_nonlinearity = true;
if isfield(RadarParams.Waveform.nonideal, 'enable_phase_noise')
enable_phase_noise = logical(RadarParams.Waveform.nonideal.enable_phase_noise);
end
if isfield(RadarParams.Waveform.nonideal, 'enable_nonlinearity')
enable_nonlinearity = logical(RadarParams.Waveform.nonideal.enable_nonlinearity);
end
use_datasheet_phase_noise = false;
phase_noise_cfg = struct();
if isfield(RadarParams.Waveform.nonideal, 'phaseNoise')
phase_noise_cfg = RadarParams.Waveform.nonideal.phaseNoise;
if isfield(phase_noise_cfg,'enabled') && phase_noise_cfg.enabled && ...
isfield(phase_noise_cfg,'offset_Hz') && isfield(phase_noise_cfg,'level_dBc_Hz')
use_datasheet_phase_noise = true;
end
end
mimoMode = RadarParams.Waveform.mimoMode;
NumTx = RadarParams.Antenna.NumTx;
rxPathGain_dB = RadarParams.Rxpath.rxPathGain_dB;
system_NF_dB = RadarParams.Rxpath.system_NF_dB;
PA_Profile = RadarParams.RFOutput.PA_Profile;
SpurParams = RadarParams.SpurParams;
enable_spur = true;
if isstruct(SpurParams) && isfield(SpurParams, 'enabled')
enable_spur = logical(SpurParams.enabled);
end
[NumTargets, NumRx, ~, NumChirps, N_adc_samples] = size(RxOut.space_loss_amp);
% TDM : TX별
if strcmpi(mimoMode, 'TDM') && mod(NumChirps, NumTx) ~= 0
error('TDM mode requires NumChirps to be a multiple of NumTx. Current NumChirps=%d, NumTx=%d.', NumChirps, NumTx);
end
% t=0 ( ) , ADC
t_adc = Timing.AdcStartTime : 1/fs_adc : (Timing.AdcStartTime + Timing.AdcSampTime - 1/fs_adc);
f_inst = f_start + Slope * t_adc;
% --- 1. PA_Profile을 (Vtx) ---
% PA_Profile의 ( PA_Profile(:,1), PA_Profile(:,2) )
if isstruct(PA_Profile)
P_inst_dBm = interp1(PA_Profile.freqs, PA_Profile.power_dBm, f_inst, 'linear', 'extrap');
else
P_inst_dBm = interp1(PA_Profile(:,1), PA_Profile(:,2), f_inst, 'linear', 'extrap');
end
% dBm -> Watt -> (V) (50 )
Vtx_inst = sqrt(10.^((P_inst_dBm - 30) / 10) * 50);
% --- 2. ---
T_ref = RadarParams.Basic.T0;
P_noise_floor_W = RadarParams.Basic.kb * T_ref * fs_adc;
rxPathGain_lin = 10^(rxPathGain_dB/10);
system_NF_lin = 10^(system_NF_dB/10);
% ADC
P_noise_total_W = P_noise_floor_W * system_NF_lin * rxPathGain_lin;
sigma_n = sqrt(P_noise_total_W * 50);
% --- MIMO TX ---
% TDM: TX만
% DDMA: TX +
% legacy phase-noise model ( )
window_size = 50;
ma_filter = ones(1,window_size)/window_size;
% --- 3. ADC Raw Data ---
adc_raw_data = zeros(NumRx, NumTx, NumChirps, N_adc_samples);
for m = 1:NumChirps
% generate phase noise & nonlinearity for this chirp
t_rel = t_adc; % time within chirp
if enable_phase_noise && use_datasheet_phase_noise
phi_noise = generate_phase_noise_from_datasheet(t_rel, phase_noise_cfg);
elseif enable_phase_noise
phi_noise = pn_level * filter(ma_filter, 1, randn(1, N_adc_samples));
else
phi_noise = zeros(1, N_adc_samples);
end
phi_nonlin = zeros(1, N_adc_samples);
if enable_nonlinearity
ramp_idx2 = (t_rel >= Timing.IdleTime);
phi_nonlin(ramp_idx2) = peak_phase_error * sin(2 * pi * f_ripple * (t_rel(ramp_idx2) - Timing.IdleTime));
end
% MIMO TX
if strcmpi(mimoMode, 'TDM')
% TDM: TX (1,2,3,1,2,3,...)
active_tx_list = mod(m - 1, NumTx) + 1; % : m=1TX1, m=2TX2, m=3TX1...
elseif strcmpi(mimoMode, 'DDMA')
% DDMA: TX
active_tx_list = 1:NumTx;
else
% : TX
active_tx_list = 1:NumTx;
end
for tx = 1:NumTx
% TDM TX는
if strcmpi(mimoMode, 'TDM') && ~ismember(tx, active_tx_list)
adc_raw_data(:, tx, m, :) = 0; % TX는 0
continue;
end
for rx = 1:NumRx
signal_mix = zeros(1, N_adc_samples);
for kt = 1:NumTargets
% (1) : Vtx(f) * G_tx * Space_Loss(f) * G_rx * rxPathGain
A_path = squeeze(RxOut.space_loss_amp(kt, rx, tx, m, :))';
% (Vtx_inst)
A_total = Vtx_inst .* (TxOut.G_tx_amp(kt, tx) * A_path * sqrt(rxPathGain_lin));
% (2) TTD Beat Phase ( t_adc )
tau = squeeze(RxOut.tau(kt, rx, tx, m, :))';
beat_phase = 2*pi * (f_start * tau + Slope * tau .* t_adc - 0.5 * Slope * tau.^2);
% add transmit impairments (phase noise + nonlinearity)
beat_phase = beat_phase + phi_noise + phi_nonlin;
% (3) DDMA
if strcmpi(mimoMode, 'DDMA')
phase_shift = 2 * pi * (tx - 1) * (m - 1) / NumTx;
beat_phase = beat_phase + phase_shift;
end
% (4)
base_sig = A_total .* exp(1j * beat_phase);
if enable_spur && exist('SpurParams','var') && ~isempty(SpurParams)
spur_vec = generate_spur(SpurParams, t_adc, base_sig);
signal_mix = signal_mix + base_sig + spur_vec;
else
signal_mix = signal_mix + base_sig;
end
end
% (5)
noise = (sigma_n/sqrt(2)) * (randn(1, N_adc_samples) + 1j*randn(1, N_adc_samples));
adc_raw_data(rx, tx, m, :) = signal_mix + noise;
end
end
end
% --- 4. ADC ( TX RX별로 ) ---
% RX가 TX
adc_combined = zeros(NumRx, NumChirps, N_adc_samples);
for m = 1:NumChirps
for rx = 1:NumRx
% TX의
combined_signal = squeeze(sum(adc_raw_data(rx, :, m, :), 2));
% ADC spur이 (orientation )
if enable_spur && exist('SpurParams','var') && isfield(SpurParams,'adc')
spur_adc = generate_spur(SpurParams, t_adc);
% ensure same shape as combined_signal (Nx1 vs 1xN)
spur_adc = reshape(spur_adc, size(combined_signal));
combined_signal = combined_signal + spur_adc;
end
adc_combined(rx, m, :) = combined_signal;
end
end
% --- 5. : DC offset/LO leakage clipping ---
if enable_spur && exist('SpurParams','var')
% LO leakage: /DC
if isfield(SpurParams,'lo_leak')
dc_amp = SpurParams.lo_leak.amp;
adc_combined = adc_combined + dc_amp;
end
% Clipping spur:
if isfield(SpurParams,'clip')
target_range = SpurParams.clip.range;
beat_freq = 2 * Slope * target_range / RadarParams.Basic.c;
clip_amp = SpurParams.clip.amp;
for m = 1:NumChirps
% ()
clip_tone = clip_amp * cos(2*pi*beat_freq*t_adc);
for rx = 1:NumRx
adc_combined(rx, m, :) = adc_combined(rx, m, :) + reshape(clip_tone, [1, 1, length(clip_tone)]);
end
end
end
end
end
@@ -0,0 +1,58 @@
function phi_noise = generate_phase_noise_from_datasheet(t, phaseNoiseCfg)
% DATASHEET SSB phase noise(dBc/Hz) (rad)
% :
% t : [1xN] [Nx1] ()
% phaseNoiseCfg.offset_Hz : (Hz)
% phaseNoiseCfg.level_dBc_Hz : SSB phase noise (dBc/Hz)
% :
% phi_noise : [1xN] (rad)
t = t(:);
N = numel(t);
phi_noise = zeros(1, N);
if N < 2
return;
end
dt = mean(diff(t));
fs = 1 / dt;
df = fs / N;
offsets = phaseNoiseCfg.offset_Hz(:);
levels_dBc = phaseNoiseCfg.level_dBc_Hz(:);
valid = isfinite(offsets) & isfinite(levels_dBc) & offsets > 0;
offsets = offsets(valid);
levels_dBc = levels_dBc(valid);
if isempty(offsets)
return;
end
[offsets, order] = sort(offsets, 'ascend');
levels_dBc = levels_dBc(order);
f_pos = (1:floor(N/2))' * df;
if isempty(f_pos)
return;
end
if numel(offsets) == 1
L_dBc = levels_dBc(1) * ones(size(f_pos));
else
L_dBc = interp1(log10(offsets), levels_dBc, log10(f_pos), 'linear', 'extrap');
end
% SSB phase noise L(f) PSD : L(f) ~= 0.5 * S_phi(f)
% => S_phi(f) ~= 2 * 10^(L(f)/10) [rad^2/Hz]
S_phi = 2 * 10.^(L_dBc / 10);
%
amp = sqrt(2 * S_phi * df); % (rad)
rand_phase = 2 * pi * rand(size(f_pos));
t_rel = t - t(1);
% NxK
phase_matrix = 2*pi*(t_rel * f_pos.') + rand_phase.';
phi_noise = (cos(phase_matrix) * amp).';
phi_noise = phi_noise(:).';
end
+83
View File
@@ -0,0 +1,83 @@
function spur = generate_spur(SpurParams, t, baseSignal)
%GENERATE_SPUR create spur signals based on a parameter structure
% spur = GENERATE_SPUR(SpurParams, t, baseSignal)
% t : time vector (row or column)
% baseSignal : (optional) complex baseband signal used for some spur
% mechanisms such as mixer nonlinearity. If omitted the
% function will only generate independent tones.
%
% SpurParams is a struct that may contain any of the following fields
% .lo : struct with fields amp, freq, phase (optional) -
% LO-related spur added to phase of signal
% .mixer : struct with fields alpha2, alpha3 - coefficients for
% 2nd/3rd order nonlinearity. Requires baseSignal input.
% .adc : struct with fields amp, freq, phase - tone added after
% ADC (complex). Does not require baseSignal.
% .switch : struct with fields amp, freq - square/pulse spur
% .pulse : struct with fields amp, freq, duty (01), phase -
% periodic gating/pulse train (lock) added to output
%
% Example:
% SpurParams = struct();
% SpurParams.lo = struct('amp',0.01,'freq',1e6,'phase',0);
% SpurParams.mixer= struct('alpha2',1e-4,'alpha3',1e-6);
% SpurParams.adc = struct('amp',1e-3,'freq',2e6,'phase',0);
% SpurParams.switch = struct('amp',5e-4,'freq',5e5);
% SpurParams.pulse = struct('amp',0.02,'freq',500e3,'duty',0.1,'phase',0);
% LO , , ADC ,
% 10% duty pulse train .
%
% The returned spur vector has the same dimensions as t. If baseSignal
% is supplied the mixer nonlinearity is computed elementwise using that
% signal; otherwise only independent spur terms are returned.
if nargin < 3
baseSignal = [];
end
% ensure t is row for consistent operations
t = t(:)';
spur = zeros(size(t));
if isfield(SpurParams, 'lo')
p = SpurParams.lo;
phase = 0;
if isfield(p, 'phase'); phase = p.phase; end
spur = spur + p.amp .* sin(2*pi*p.freq .* t + phase);
end
if isfield(SpurParams, 'mixer') && ~isempty(baseSignal)
m = SpurParams.mixer;
% apply polynomial nonlinearity to the provided base signal
if isfield(m, 'alpha2')
spur = spur + m.alpha2 .* (baseSignal.^2);
end
if isfield(m, 'alpha3')
spur = spur + m.alpha3 .* (baseSignal.^3);
end
end
if isfield(SpurParams, 'adc')
a = SpurParams.adc;
phase = 0;
if isfield(a, 'phase'); phase = a.phase; end
spur = spur + a.amp .* exp(1j*(2*pi*a.freq .* t + phase));
end
if isfield(SpurParams, 'switch')
s = SpurParams.switch;
spur = spur + s.amp .* square(2*pi*s.freq .* t);
end
if isfield(SpurParams, 'pulse')
p = SpurParams.pulse;
% duty default 50%
duty = 0.5;
if isfield(p,'duty'); duty = p.duty; end
phase = 0;
if isfield(p,'phase'); phase = p.phase; end
% generate pulse train (0/1) scaled by amp
% square with duty cycle multiplied then shifted to [0,1]
pulse_wave = (square(2*pi*p.freq .* t + phase, duty*100)+1)/2;
spur = spur + p.amp .* pulse_wave;
end
end
+21
View File
@@ -0,0 +1,21 @@
function RxOut = receive_antenna(ChannelOut, Target, RxPattern)
[NumTargets, NumRx, NumTx, NumChirps, N_adc] = size(ChannelOut.space_loss_amp);
RxOut = ChannelOut;
for k = 1:NumTargets
for rx = 1:NumRx
pat = RxPattern(rx);
if strcmpi(pat.Type, '2D')
g_dBi = interp2(pat.az_angles, pat.el_angles, pat.gain_dBi, Target.az(k), Target.el(k), 'linear', -20);
else
g_az = interp1(pat.az_angles, pat.gain_az_dBi, Target.az(k), 'linear', -20);
g_el = interp1(pat.el_angles, pat.gain_el_dBi, Target.el(k), 'linear', -20);
g_dBi = g_az + g_el - pat.max_gain_dBi;
end
%
rx_gain_amp = sqrt(10^(g_dBi/10));
RxOut.space_loss_amp(k, rx, :, :, :) = RxOut.space_loss_amp(k, rx, :, :, :) * rx_gain_amp;
end
end
end
+140
View File
@@ -0,0 +1,140 @@
function [det_mask, threshold_map, detections] = detect_targets_cfar(rd_map, RadarParams)
% CFAR target detection for range-doppler map
% - method: 'CA' or 'OS'
% - dimension: '1D' or '2D'
% - axis (for 1D): 'range' or 'doppler'
%
% Input:
% rd_map : [Ndoppler x Nrange] complex or real map
% RadarParams : struct containing SP.CFAR options
%
% Output:
% det_mask : logical detection mask [Ndoppler x Nrange]
% threshold_map : threshold map [Ndoppler x Nrange]
% detections : [Ndet x 2] = [doppler_bin, range_bin]
cfg = RadarParams.SP.CFAR;
method = upper(string(cfg.method));
dim_mode = upper(string(cfg.dimension));
axis_mode = lower(string(cfg.axis));
pfa = cfg.pfa;
train = cfg.train;
guard = cfg.guard;
if numel(train) == 1
train = [train, train];
end
if numel(guard) == 1
guard = [guard, guard];
end
os_rank_ratio = cfg.rank;
os_scale = cfg.os_scale;
rd_power = abs(rd_map).^2;
[n_dop, n_rng] = size(rd_power);
det_mask = false(n_dop, n_rng);
threshold_map = nan(n_dop, n_rng);
switch dim_mode
case "2D"
td = train(1); tr = train(2);
gd = guard(1); gr = guard(2);
for d = (td+gd+1):(n_dop-(td+gd))
for r = (tr+gr+1):(n_rng-(tr+gr))
d_idx = (d-(td+gd)):(d+(td+gd));
r_idx = (r-(tr+gr)):(r+(tr+gr));
win = rd_power(d_idx, r_idx);
cut_d = td+gd+1;
cut_r = tr+gr+1;
guard_mask = false(size(win));
guard_mask((cut_d-gd):(cut_d+gd), (cut_r-gr):(cut_r+gr)) = true;
train_cells = win(~guard_mask);
th = local_cfar_threshold(train_cells, method, pfa, os_rank_ratio, os_scale);
threshold_map(d, r) = th;
det_mask(d, r) = rd_power(d, r) > th;
end
end
case "1D"
switch axis_mode
case "range"
tr = train(2); gr = guard(2);
for d = 1:n_dop
[det_row, th_row] = cfar_1d_line(rd_power(d, :), tr, gr, method, pfa, os_rank_ratio, os_scale);
det_mask(d, :) = det_row;
threshold_map(d, :) = th_row;
end
case "doppler"
td = train(1); gd = guard(1);
for r = 1:n_rng
[det_col, th_col] = cfar_1d_line(rd_power(:, r).', td, gd, method, pfa, os_rank_ratio, os_scale);
det_mask(:, r) = det_col.';
threshold_map(:, r) = th_col.';
end
otherwise
error('CFAR axis must be ''range'' or ''doppler'' when dimension is 1D.');
end
otherwise
error('CFAR dimension must be ''1D'' or ''2D''.');
end
[d_idx, r_idx] = find(det_mask);
detections = [d_idx, r_idx];
end
function [det_line, th_line] = cfar_1d_line(x, t, g, method, pfa, rank_ratio, os_scale)
n = numel(x);
det_line = false(1, n);
th_line = nan(1, n);
left = t + g;
right = t + g;
for i = (left+1):(n-right)
l_train = x((i-g-t):(i-g-1));
r_train = x((i+g+1):(i+g+t));
train_cells = [l_train, r_train];
th = local_cfar_threshold(train_cells, method, pfa, rank_ratio, os_scale);
th_line(i) = th;
det_line(i) = x(i) > th;
end
end
function th = local_cfar_threshold(train_cells, method, pfa, rank_ratio, os_scale)
train_cells = train_cells(:);
n_train = numel(train_cells);
if n_train == 0
th = inf;
return;
end
switch method
case "CA"
noise_hat = mean(train_cells);
alpha = n_train * (pfa^(-1/n_train) - 1);
th = alpha * noise_hat;
case "OS"
sorted_cells = sort(train_cells, 'ascend');
k = max(1, min(n_train, round(rank_ratio * n_train)));
noise_hat = sorted_cells(k);
th = os_scale * noise_hat;
otherwise
error('CFAR method must be ''CA'' or ''OS''.');
end
end
+18
View File
@@ -0,0 +1,18 @@
function target_rd_map = integrate_nci_rdm(rd_cube, RadarParams)
% NCI (Noncoherent Integration) 2D RDM
% :
% rd_cube: [RX, TX, Doppler, Range]
% RadarParams.Waveform.mimoMode: 'TDM'
% :
% target_rd_map: [Doppler, Range]
%
% :
% - TDM : TX-RX (|.|^2)
% - : TX=1 RX
if strcmpi(RadarParams.Waveform.mimoMode, 'TDM')
target_rd_map = squeeze(sum(sum(abs(rd_cube).^2, 1), 2));
else
target_rd_map = squeeze(sum(abs(rd_cube(:, 1, :, :)).^2, 1));
end
end
@@ -0,0 +1,45 @@
function [rd_map, doppler_axis] = process_doppler_fft(range_profile, RadarParams)
% :
% - range_profile: [NumRx, NumTx, NumChirps, NumRangeBins]
% - RadarParams:
NumChirps = RadarParams.Waveform.NumChirps;
window_type = RadarParams.SP.RDM.window_type_doppler;
[~, ~, ~, ~] = size(range_profile);
%
lambda = RadarParams.Waveform.lambda_c; % = c/fc
% 1. (PRI, Pulse Repetition Interval)
T_pri = RadarParams.Waveform.PRI;
% 2. Doppler-FFT용 ( )
% (3 ) .
if strcmpi(window_type, 'none')
win_doppler = ones(1, NumChirps);
elseif strcmpi(window_type, 'hamming')
win_doppler = hamming(NumChirps)';
elseif strcmpi(window_type, 'blackman')
win_doppler = blackman(NumChirps)';
elseif strcmpi(window_type, 'hann')
win_doppler = hann(NumChirps)';
else % default: 'chebwin'
win_doppler = chebwin(NumChirps, 60)'; % 60dB
end
win_data = range_profile .* reshape(win_doppler, [1, 1, NumChirps, 1]);
% 3. Doppler-FFT (3 : Chirp)
% NFFT를 NumChirps보다 .
rd_fft = fft(win_data, NumChirps, 3);
% 4. fftshift ( 0 )
% [-, 0, +] .
rd_map = fftshift(rd_fft, 3);
% 5. (Velocity Axis)
% Vmax = lambda / (4 * T_pri)
v_max = lambda / (4 * T_pri);
% dv = lambda / (2 * NumChirps * T_pri)
doppler_axis = linspace(-v_max, v_max, NumChirps);
end
@@ -0,0 +1,48 @@
function [range_profile, range_axis] = process_range_fft_lpf(adc_raw_data, RadarParams)
% :
% - adc_raw_data: [NumRx, NumTx, NumChirps, N_samples]
% - RadarParams:
fs_adc = RadarParams.Waveform.fs_adc;
Slope = RadarParams.Waveform.Slope;
fc_lpf_Hz = RadarParams.Rxpath.fc_lpf;
window_type = RadarParams.SP.RDM.window_type_range;
[~, ~, ~, N_samples] = size(adc_raw_data);
c = RadarParams.Basic.c;
% 1.
if strcmpi(window_type, 'none')
win = ones(1, N_samples);
elseif strcmpi(window_type, 'hamming')
win = hamming(N_samples)';
elseif strcmpi(window_type, 'blackman')
win = blackman(N_samples)';
elseif strcmpi(window_type, 'chebwin')
win = chebwin(N_samples, 60)';
else % default: 'hann'
win = hann(N_samples)';
end
win_data = adc_raw_data .* reshape(win, [1, 1, 1, N_samples]);
% 2. Range-FFT
range_fft = fft(win_data, N_samples, 4);
% 3.
df = fs_adc / N_samples;
freq_axis = (0 : N_samples-1) * df;
full_range_axis = (freq_axis * c) / (2 * Slope);
% 4. Ideal LPF
cutoff_idx = floor(fc_lpf_Hz / df);
%
half_idx = floor(N_samples/2);
final_idx = min(cutoff_idx, half_idx);
% 5. [ ] (Truncation)
% 1 final_idx까지만
range_profile = range_fft(:, :, :, 1:final_idx);
range_axis = full_range_axis(1:final_idx);
end
@@ -0,0 +1,159 @@
function fig = visualize_rd_map_with_spurs(target_rd_map, r_axis, v_axis, Target, SpurParams, RadarParams, Slope)
% 2D Range-Doppler .
%
% :
% - target_rd_map: 2D RD [Doppler x Range]
% - r_axis:
% - v_axis:
% - Target: (R, v, NumTargets )
% - SpurParams: (lo, adc, mixer, pulse, lo_leak, clip )
% - RadarParams:
% - Slope: (Hz/s)
%
% :
% - fig: figure
fig = figure('Name', '2D Range-Doppler Map');
imagesc(r_axis, v_axis, 20*log10(abs(target_rd_map)));
axis xy; % y축 ()
colorbar;
xlabel('Range (m)');
ylabel('Velocity (m/s)');
title('Range-Doppler Map (Single Channel)');
colormap(jet);
%
hLines = gobjects(0);
hLabels = {};
% --- 1. ---
hold on;
hTarget = gobjects(0);
if exist('Target','var') && isstruct(Target)
for k = 1:Target.NumTargets
% /
if Target.R(k) < min(r_axis) || Target.R(k) > max(r_axis) || ...
Target.v(k) < min(v_axis) || Target.v(k) > max(v_axis)
continue; %
end
[~, ix] = min(abs(r_axis - Target.R(k)));
[~, iy] = min(abs(v_axis - Target.v(k)));
hTarget(end+1) = plot(r_axis(ix), v_axis(iy), 'ro', 'MarkerSize',8, 'LineWidth',1.5);
text(r_axis(ix), v_axis(iy), sprintf(' T%d', k), 'Color','r','FontSize',8);
end
end
hold off;
%
if exist('hTarget','var') && ~isempty(hTarget)
hLines(end+1) = hTarget(1);
hLabels{end+1} = 'actual target';
end
% --- nonideal/spur ---
show_nonideal_overlay = true;
if exist('SpurParams','var') && isstruct(SpurParams) && isfield(SpurParams,'enabled')
show_nonideal_overlay = logical(SpurParams.enabled);
end
if exist('RadarParams','var') && isstruct(RadarParams) && ...
isfield(RadarParams,'Waveform') && isfield(RadarParams.Waveform,'nonideal')
ni = RadarParams.Waveform.nonideal;
enable_phase_noise = false;
enable_nonlinearity = false;
enable_datasheet_phase_noise = false;
if isfield(ni,'enable_phase_noise')
enable_phase_noise = logical(ni.enable_phase_noise);
end
if isfield(ni,'enable_nonlinearity')
enable_nonlinearity = logical(ni.enable_nonlinearity);
end
if isfield(ni,'phaseNoise') && isstruct(ni.phaseNoise) && isfield(ni.phaseNoise,'enabled')
enable_datasheet_phase_noise = logical(ni.phaseNoise.enabled);
end
show_nonideal_overlay = show_nonideal_overlay && ...
(enable_phase_noise || enable_nonlinearity || enable_datasheet_phase_noise);
end
% --- 2. LO/ADC ---
if show_nonideal_overlay && exist('SpurParams','var') && isstruct(SpurParams)
spur_freqs = [];
spur_names = {};
if isfield(SpurParams,'lo') && ~isempty(SpurParams.lo)
spur_freqs(end+1) = SpurParams.lo.freq;
spur_names{end+1} = 'LO spur';
end
if isfield(SpurParams,'adc') && ~isempty(SpurParams.adc)
spur_freqs(end+1) = SpurParams.adc.freq;
spur_names{end+1} = 'ADC spur';
end
hold on;
for idx = 1:length(spur_freqs)
f = spur_freqs(idx);
name = sprintf('%s @ %.2f MHz', spur_names{idx}, f/1e6);
r_spur = (RadarParams.Basic.c * f) / (2 * Slope);
hLines(end+1) = xline(r_spur, '--m', name);
hLabels{end+1} = name;
end
% --- 3. Pulse comb ---
if isfield(SpurParams,'pulse') && ~isempty(SpurParams.pulse)
p = SpurParams.pulse;
num_harm = 5;
for n = 1:num_harm
f_n = n * p.freq;
r_n = (RadarParams.Basic.c * f_n) / (2 * Slope);
if r_n <= max(r_axis)
lbl = sprintf('pulse comb %d', n);
hLines(end+1) = xline(r_n, '--g', lbl);
hLabels{end+1} = lbl;
end
end
end
% --- 4. Mixer ---
if isfield(SpurParams,'mixer') && ~isempty(SpurParams.mixer) && exist('Target','var')
base_ranges = Target.R;
harmonics = [2,3];
colors = {'--c','--y'};
for hi = 1:length(harmonics)
mul = harmonics(hi);
for r0 = base_ranges
r_spur2 = r0 * mul;
if r_spur2 <= max(r_axis)
lbl = sprintf('mixer x%d', mul);
hLines(end+1) = xline(r_spur2, colors{hi}, lbl);
hLabels{end+1} = lbl;
end
end
end
end
% --- 5. DC offset / LO leakage ---
if isfield(SpurParams,'lo_leak') && ~isempty(SpurParams.lo_leak)
hLines(end+1) = line([min(r_axis), max(r_axis)], [0, 0], 'Color','r','LineStyle',':', 'LineWidth',1.5);
hLabels{end+1} = 'DC offset/LO leak';
end
% --- 6. Clipping spur ---
if isfield(SpurParams,'clip') && ~isempty(SpurParams.clip)
r_clip = SpurParams.clip.range;
if r_clip >= min(r_axis) && r_clip <= max(r_axis)
hLines(end+1) = line([r_clip, r_clip], [min(v_axis), max(v_axis)], 'Color', 'r', 'LineStyle', ':', 'LineWidth', 1.2);
hLabels{end+1} = sprintf('clipping @ %.1f m', r_clip);
end
end
hold off;
end
% --- 7. ---
if ~isempty(hLines)
legend(hLines, hLabels, 'Location', 'northeastoutside');
end
end
@@ -0,0 +1,68 @@
function PatternArray = build_antenna_patterns(NumAnt, PatternType, MaxGain, AzSquintStep, AzBeamwidth, ElBeamwidth, FileList)
% build_antenna_patterns: .
%
% []
% - FileList (): Cell Array.
% : {'tx1.mat', 'tx2.mat', ''} ( )
% ( )
az_angles_default = -90:1:90;
el_angles_default = -90:1:90;
% nargin을 FileList가
if nargin < 7
FileList = {};
end
for n = 1:NumAnt
PatternArray(n).Type = PatternType;
% -------------------------------------------------------------
% Case A: , (Load)
% -------------------------------------------------------------
if length(FileList) >= n && ~isempty(FileList{n}) && isfile(FileList{n})
filePath = FileList{n};
loadedData = load(filePath); % .mat
% (az_angles, el_angles ) .
PatternArray(n).az_angles = loadedData.az_angles;
PatternArray(n).el_angles = loadedData.el_angles;
if strcmpi(PatternType, '2D')
PatternArray(n).gain_dBi = loadedData.gain_dBi;
elseif strcmpi(PatternType, '1D')
PatternArray(n).max_gain_dBi = loadedData.max_gain_dBi;
PatternArray(n).gain_az_dBi = loadedData.gain_az_dBi;
PatternArray(n).gain_el_dBi = loadedData.gain_el_dBi;
end
%
continue; %
end
% -------------------------------------------------------------
% Case B: ( )
% -------------------------------------------------------------
PatternArray(n).az_angles = az_angles_default;
PatternArray(n).el_angles = el_angles_default;
% (Squint Angle)
% (3): (1-2)*5 = -5, (2-2)*5 = 0, (3-2)*5 = +5
squint_az = (n - (NumAnt + 1)/2) * AzSquintStep;
if strcmpi(PatternType, '2D')
[AZ, EL] = meshgrid(PatternArray(n).az_angles, PatternArray(n).el_angles);
PatternArray(n).gain_dBi = MaxGain - 3*((AZ - squint_az)/AzBeamwidth).^2 - 3*(EL/ElBeamwidth).^2;
PatternArray(n).gain_dBi(PatternArray(n).gain_dBi < -20) = -20;
elseif strcmpi(PatternType, '1D')
PatternArray(n).max_gain_dBi = MaxGain;
PatternArray(n).gain_az_dBi = MaxGain - 3*((PatternArray(n).az_angles - squint_az)/AzBeamwidth).^2;
PatternArray(n).gain_az_dBi(PatternArray(n).gain_az_dBi < -20) = -20;
PatternArray(n).gain_el_dBi = MaxGain - 3*(PatternArray(n).el_angles/ElBeamwidth).^2;
PatternArray(n).gain_el_dBi(PatternArray(n).gain_el_dBi < -20) = -20;
end
end
end
@@ -0,0 +1,81 @@
function fig = visualize_antenna_pattern(PatternArray, AntennaName)
% visualize_antenna_pattern: .
% AntennaName (: 'TX', 'RX') .
if nargin < 2
AntennaName = 'Antenna'; %
end
fig_name = sprintf('%s Radiation Pattern Cuts', AntennaName);
fig = figure('Name', fig_name, 'Position', [200, 200, 1000, 700]);
NumAnt = length(PatternArray);
colors = lines(NumAnt);
% [ Y축 ]
y_min_total = inf;
y_max_total = -inf;
for i = 1:NumAnt
pat = PatternArray(i);
if strcmpi(pat.Type, '2D')
y_min_total = min(y_min_total, min(pat.gain_dBi(:)));
y_max_total = max(y_max_total, max(pat.gain_dBi(:)));
else
y_min_total = min(y_min_total, min([pat.gain_az_dBi(:); pat.gain_el_dBi(:)]));
y_max_total = max(y_max_total, max([pat.gain_az_dBi(:); pat.gain_el_dBi(:)]));
end
end
% ---------------------------------------------------------
% [Subplot 1] : Azimuth Pattern Cut
% ---------------------------------------------------------
subplot(2, 1, 1);
hold on; grid on;
for i = 1:NumAnt
pat = PatternArray(i);
if strcmpi(pat.Type, '2D')
el0_idx = find(abs(pat.el_angles - 0) < 1e-6, 1);
az_gain_cut = pat.gain_dBi(el0_idx, :);
else
az_gain_cut = pat.gain_az_dBi;
end
plot(pat.az_angles, az_gain_cut, '-', 'Color', colors(i,:), 'LineWidth', 2.5, 'DisplayName', sprintf('%s %d', AntennaName, i));
end
title(sprintf('[%s] Azimuth Radiation Pattern Cut (at Elevation = 0^\\circ)', AntennaName), 'FontSize', 12);
xlabel('Azimuth Angle (deg)', 'FontSize', 11);
ylabel('Antenna Gain (dBi)', 'FontSize', 11);
xlim([min(PatternArray(1).az_angles), max(PatternArray(1).az_angles)]);
ylim([y_min_total - 2, y_max_total + 2]);
legend('show', 'Location', 'south');
hold off;
% ---------------------------------------------------------
% [Subplot 2] : Elevation Pattern Cut
% ---------------------------------------------------------
subplot(2, 1, 2);
hold on; grid on;
for i = 1:NumAnt
pat = PatternArray(i);
if strcmpi(pat.Type, '2D')
az0_idx = find(abs(pat.az_angles - 0) < 1e-6, 1);
el_gain_cut = pat.gain_dBi(:, az0_idx);
else
el_gain_cut = pat.gain_el_dBi;
end
plot(pat.el_angles, el_gain_cut, '-', 'Color', colors(i,:), 'LineWidth', 2.5, 'DisplayName', sprintf('%s %d', AntennaName, i));
end
title(sprintf('[%s] Elevation Radiation Pattern Cut (at Azimuth = 0^\\circ)', AntennaName), 'FontSize', 12);
xlabel('Elevation Angle (deg)', 'FontSize', 11);
ylabel('Antenna Gain (dBi)', 'FontSize', 11);
xlim([min(PatternArray(1).el_angles), max(PatternArray(1).el_angles)]);
ylim([y_min_total - 2, y_max_total + 2]);
legend('show', 'Location', 'south');
hold off;
end
+244
View File
@@ -0,0 +1,244 @@
% =========================================================================
% Tx_Main.m
% Pure MATLAB FMCW TX Simulator (Multi-Chirp Frame Generation)
% =========================================================================
clear; clc; close all;
% =============== ===============
currentFilePath = mfilename('fullpath');
currentFolder = fileparts(currentFilePath);
addpath(genpath(currentFolder));
%% =================== ===================
% 1)
RadarParams.Basic.c = 3e8; % (m/s)
RadarParams.Basic.kb = physconst('Boltzmann'); %
RadarParams.Basic.T0 = 290; % (K)
% 2)
RadarParams.Waveform.fc = 77e9; % (Hz)
RadarParams.Waveform.lambda_c = RadarParams.Basic.c / RadarParams.Waveform.fc; % (m)
RadarParams.Waveform.NumChirps = 128; %
RadarParams.Waveform.B_valid = 1e9; % (Hz)
RadarParams.Waveform.fs_adc = 10e6; % ADC (Hz)
RadarParams.Waveform.fs_waveform = 1e6; % (Hz)
% 2-1)
RadarParams.Waveform.Timing.IdleTime = 7e-6; % (s)
RadarParams.Waveform.Timing.TxStartTime = 1e-6; % Tx on (s)
RadarParams.Waveform.Timing.AdcStartTime = 6e-6; % ADC on (s)
RadarParams.Waveform.Timing.AdcSampTime = 50e-6; % ADC (s)
RadarParams.Waveform.Timing.ExcessTime = 1e-6; % (s)
RadarParams.Waveform.Timing.RampEndTime = RadarParams.Waveform.Timing.IdleTime + RadarParams.Waveform.Timing.TxStartTime + RadarParams.Waveform.Timing.AdcSampTime + RadarParams.Waveform.Timing.ExcessTime; % (s)
RadarParams.Waveform.PRI = RadarParams.Waveform.Timing.IdleTime + RadarParams.Waveform.Timing.RampEndTime; % Pulse Repetition Interval (s)
RadarParams.Waveform.PRF = 1 / RadarParams.Waveform.PRI; % Pulse Repetition Frequency (Hz)
RadarParams.Waveform.Slope = RadarParams.Waveform.B_valid / RadarParams.Waveform.Timing.AdcSampTime; % (Hz/s)
RadarParams.Waveform.f_start = RadarParams.Waveform.fc - RadarParams.Waveform.B_valid/2; % (Hz)
% 2-2) Phase noise & nonlinearity parameters
RadarParams.Waveform.nonideal.pn_level = 0.05; %
RadarParams.Waveform.nonideal.f_ripple = 300e3; %
RadarParams.Waveform.nonideal.peak_phase_error = 0.1; %
RadarParams.Waveform.nonideal.power_drop_edge = 0.8; %
RadarParams.Waveform.nonideal.enable_phase_noise = false; % phase noise on/off
RadarParams.Waveform.nonideal.enable_nonlinearity = false; % nonlinearity on/off
% Datasheet phase noise (: -89 dBc/Hz @ 1 MHz)
RadarParams.Waveform.nonideal.phaseNoise.offset_Hz = [1e6];
RadarParams.Waveform.nonideal.phaseNoise.level_dBc_Hz = [-89];
RadarParams.Waveform.nonideal.phaseNoise.enabled = false;
% 2-3) MIMO
RadarParams.Waveform.mimoMode = 'TDM'; % 'TDM' 'DDMA'
% 3) RF - [dBm]
RadarParams.RFOutput.PA_Profile.freqs = [76e9, 76.5e9, 77e9, 77.5e9, 78e9, 79e9];
RadarParams.RFOutput.PA_Profile.power_dBm = [10.5, 11.5, 12.0, 11.5, 10.5, 8.0];
% 4)
% 4-1) (3xN , : m)
RadarParams.Antenna.NumTx = 2;
RadarParams.Antenna.NumRx = 4;
RadarParams.Antenna.lambda = RadarParams.Waveform.lambda_c;
RadarParams.Antenna.TxPos = [ (0:RadarParams.Antenna.NumTx-1)*2*RadarParams.Antenna.lambda; zeros(1,RadarParams.Antenna.NumTx); zeros(1,RadarParams.Antenna.NumTx) ]; % 3xNumTx
RadarParams.Antenna.RxPos = [ (0:RadarParams.Antenna.NumRx-1)*0.5*RadarParams.Antenna.lambda; zeros(1,RadarParams.Antenna.NumRx); zeros(1,RadarParams.Antenna.NumRx) ]; % 3xNumRx
% 4-2) ( NumTx, NumRx )
RadarParams.Antenna.tx_files = {}; % ( )
RadarParams.Antenna.rx_files = {};
RadarParams.Antenna.TxPattern = build_antenna_patterns(RadarParams.Antenna.NumTx, '1D', 12, 5, 30, 10, RadarParams.Antenna.tx_files);
RadarParams.Antenna.RxPattern = build_antenna_patterns(RadarParams.Antenna.NumRx, '2D', 14, 2, 40, 15, RadarParams.Antenna.rx_files);
% 5)
RadarParams.Target.R = [12, 80, 120]; % (m)
RadarParams.Target.v = [7, -5, 0]; % (m/s)
RadarParams.Target.rcs = [10, 5, 20]; % RCS (dBsm)
RadarParams.Target.az = [10, -10, 0]; % (deg)
RadarParams.Target.el = [0, 0, 5]; % (deg)
RadarParams.Target.NumTargets = length(RadarParams.Target.R);
% 6)
RadarParams.Rxpath.fc_hpf = 1400; % HPF (Hz)
RadarParams.Rxpath.fc_lpf = 0.8*RadarParams.Waveform.fs_adc/2; % LPF (Hz)
RadarParams.Rxpath.adc_bits = 12; % ADC
RadarParams.Rxpath.adc_v_full_scale = 2.0; % full-scale (Vp-p)
RadarParams.Rxpath.receiver_mode = 'Real'; % 'IQ' 'Real'
RadarParams.Rxpath.rxPathGain_dB = 50; % RX (dB)
RadarParams.Rxpath.system_NF_dB = 15; % (dB)
% 6-5) Spur parameters ( )
% Spur의 , LO leakage, , ADC , .
% SpurParams .
% spur의 , spur apply_lna_and_mixer SpurParams를 spur adc_combined에 .
RadarParams.SpurParams = struct();
RadarParams.SpurParams.lo = struct('amp',0.01,'freq',1e6); %
RadarParams.SpurParams.mixer= struct('alpha2',1e-4,'alpha3',1e-6); % 2/3
RadarParams.SpurParams.adc = struct('amp',1e-3,'freq',2e6); % ADC
RadarParams.SpurParams.switch = struct('amp',5e-4,'freq',500e3); %
RadarParams.SpurParams.pulse = struct('amp',0.0,'freq',200e3,'duty',0.1,'phase',0);%
RadarParams.SpurParams.lo_leak = struct('amp',0.05); % DC/LO leakage
RadarParams.SpurParams.clip = struct('amp',0.03,'range',20); % Clipping spur
RadarParams.SpurParams.enabled = false; % spur on/off
% 7)
RadarParams.SP.RDM.window_type_range = 'hann'; % Range FFT용
RadarParams.SP.RDM.window_type_doppler = 'chebwin'; % Doppler FFT용
RadarParams.SP.CFAR.method = 'OS'; % 'CA' 'OS'
RadarParams.SP.CFAR.dimension = '2D'; % '1D' '2D'
RadarParams.SP.CFAR.axis = 'doppler'; % 1D일 : 'range' 'doppler'
RadarParams.SP.CFAR.pfa = 1e-6; % false alarm
RadarParams.SP.CFAR.train = [8, 8]; % [doppler, range] training cell (1D면 )
RadarParams.SP.CFAR.guard = [2, 2]; % [doppler, range] guard cell (1D면 )
RadarParams.SP.CFAR.rank = 0.75; % OS-CFAR rank (0~1)
RadarParams.SP.CFAR.os_scale = 15.0; % OS-CFAR
%% 2. (TX )
% step 0. TX
Target = RadarParams.Target;
TxPattern = RadarParams.Antenna.TxPattern;
RxPattern = RadarParams.Antenna.RxPattern;
NumTx = RadarParams.Antenna.NumTx;
NumChirps = RadarParams.Waveform.NumChirps;
T_chirp = RadarParams.Waveform.Timing.IdleTime + RadarParams.Waveform.Timing.RampEndTime;
T_frame = T_chirp * NumChirps;
[t, tx_mask] = generate_waveform_timing(RadarParams);
% ================= ADC Raw Data (step 1 ~ 6) =================
% step 1. (Ptx )
TxOut = radiate_antenna(Target, TxPattern);
% step 2. (TTD )
ChannelOut = apply_channel_effects(RadarParams);
% step 3.
RxOut = receive_antenna(ChannelOut, Target, RxPattern);
% step 4-1. (TX_Gain * Space_Loss * RX_Gain)
% RF (V_tx) "시스템 전달 함수" .
SystemAmp = zeros(size(RxOut.space_loss_amp));
for k = 1:Target.NumTargets
for tx = 1:NumTx
SystemAmp(k, :, tx, :, :) = RxOut.space_loss_amp(k, :, tx, :, :) * TxOut.G_tx_amp(k, tx);
end
end
% step 4-2. RF/IF ADC
% - : /, PA , , , MIMO , non-ideal, spur를
% RX .
% - 1) adc_raw_data : [NumRx x NumTx x NumChirps x N_adc_samples]
% TX별로 ( TX )
% - 2) adc_combined : [NumRx x NumChirps x N_adc_samples]
% ADC ( RX에서 TX )
[adc_raw_data, adc_combined] = apply_lna_and_mixer(RxOut, TxOut, RadarParams);
% step 5. ADC Combined ( ADC ) HPF
adc_combined_hpf = apply_analog_hpf(adc_combined, RadarParams.Waveform.fs_adc, RadarParams.Rxpath.fc_hpf);
% step 6. ADC (Quantization & Clipping)
%
v_peak = max(abs(adc_combined_hpf(:)));
lsb_val = RadarParams.Rxpath.adc_v_full_scale / (2^RadarParams.Rxpath.adc_bits);
% ADC
adc_digital = apply_adc_quantization(adc_combined_hpf, RadarParams.Rxpath.adc_bits, RadarParams.Rxpath.adc_v_full_scale, RadarParams.Rxpath.receiver_mode);
% ================= (step 7 ~ ) =================
% step 7. Range-FFT Ideal LPF
% LPF
% (RadarParams.Rxpath.fc_lpf, RadarParams.SP.RDM.*)
% adc_digital을 [NumRx, 1, NumChirps, N_samples] (process_range_fft_lpf )
adc_digital_expanded = reshape(adc_digital, [size(adc_digital,1), 1, size(adc_digital,2), size(adc_digital,3)]);
[range_data, r_axis] = process_range_fft_lpf(adc_digital_expanded, RadarParams);
% step 8. Range Profile (1 , 1 )
figure('Name', 'Range Profile with Ideal LPF');
plot(r_axis, 20*log10(abs(squeeze(range_data(1,1,1,:)))));
grid on; hold on;
xlabel('Range (m)');
ylabel('Magnitude (dB)');
title(['Range Profile (LPF Cut-off: ', num2str(RadarParams.Rxpath.fc_lpf/1e6), ' MHz)']);
% LPF
xline((RadarParams.Rxpath.fc_lpf * RadarParams.Basic.c)/(2*RadarParams.Waveform.Slope), '--r', 'LPF Cut-off');
% step 9. Doppler-FFT
[rd_cube, v_axis] = process_doppler_fft(range_data, RadarParams);
% step 10. NumRx (Noncoherent Integration) RDM
target_rd_map = integrate_nci_rdm(rd_cube, RadarParams); % [Doppler x Range]
% step 11. CFAR (OS/CA, 1D/2D, range/doppler )
[cfar_mask, cfar_threshold, cfar_detections] = detect_targets_cfar(target_rd_map, RadarParams);
%% 3. (Single Chirp & Multi-Chirp)
% [Figure 1]
T_chirp = RadarParams.Waveform.Timing.IdleTime + RadarParams.Waveform.Timing.RampEndTime;
idx_single = (t <= T_chirp);
t_single = t(idx_single);
tx_mask_single = tx_mask(idx_single);
fig_single = visualize_tx_waveform(t_single, RadarParams.Waveform.Timing, RadarParams.Waveform.fc, RadarParams.Waveform.f_start, RadarParams.Waveform.Slope, tx_mask_single, RadarParams.Waveform.nonideal.peak_phase_error, RadarParams.Waveform.nonideal.f_ripple);
% [Figure 2] 퀀 MIMO
% (: mimoMode와 NumTx )
fig_multi = visualize_multi_tx_waveform(t, RadarParams.Waveform.Timing, RadarParams.Waveform.fc, RadarParams.Waveform.f_start, RadarParams.Waveform.Slope, tx_mask, NumChirps, RadarParams.Waveform.mimoMode, NumTx);
% --- [ Figure 3, 4: ] ---
% 'TX' 'RX' .
fig_tx_ant = visualize_antenna_pattern(TxPattern, 'TX');
fig_rx_ant = visualize_antenna_pattern(RxPattern, 'RX');
% [Figure 5] Range-Doppler Map (NumRx NCI RDM)
fig_rd_map = visualize_rd_map_with_spurs(target_rd_map, r_axis, v_axis, Target, RadarParams.SpurParams, RadarParams, RadarParams.Waveform.Slope);
% [Figure 6] CFAR (2D MAP )
fig_cfar = figure('Name', 'CFAR Detections on NCI RDM');
imagesc(r_axis, v_axis, 20*log10(abs(target_rd_map)));
axis xy;
colormap(jet);
colorbar;
xlabel('Range (m)');
ylabel('Velocity (m/s)');
title('CFAR Detections (Circle Markers)');
hold on;
if ~isempty(cfar_detections)
det_d_idx = cfar_detections(:,1); % Doppler bin index
det_r_idx = cfar_detections(:,2); % Range bin index
det_r = r_axis(det_r_idx);
det_v = v_axis(det_d_idx);
plot(det_r, det_v, 'wo', 'MarkerSize', 7, 'LineWidth', 1.5);
end
hold off;