Implementing Doppler Shift in FMCW Radar simulation

Hi, Im trying to implement Doppler Shift in my FMCW Radar simulation code.

Problem is that i can't find the code which transmits the FMCW signal of vradar(victom radar)

As far as i know, the FMCW signal is transmitted by the function 'step'

ex) sig = step(fmcwwav1);

But there is no code line like this for vrader in my simulation

Oddly, When i run the entire code, the results appear to be correct.

So my question,

1. what part of the code is responsible for the signal transmission of vradar
2. what part of the code should i modify to add some doppler shift to the vradar signal

I'll attach my code at the end

Also my reffernce code isright below my code

Thank you for the help, in advance

%% My Code
clear;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Define Victim FMCW Radar Waveform

rng(2017);             % Set random number generator for repeatable results
fc = 24e9;             % Center frequency (Hz)
c = physconst('LightSpeed');             % Speed of light (m/s)
lambda = freq2wavelen(fc,c);             % Wavelength (m)
rangeMax = 40.05;           % Maximum range (m)
rangeRes = 0.15;            % Desired range resolution (m)
vMax = 230*1000/3600;       % Maximum velocity of the Radar (m/s)

tSw = 6*range2time(rangeMax,c);                 % Sweep Time [s]
bw = rangeres2bw(rangeRes,c);                   % Bandwith [Hz]
slope = bw/tSw;                                 % Sweep slope 
fr_max = range2beat(rangeMax,slope,c);          % Maximum frequency corresponding to the maximum range [Hz]
fd_max = speed2dop(2*vMax,lambda);              % Maximum frequency corresponding to the Doppler Shift [Hz]
fb_max = fr_max+fd_max;                         % Maximum beat frequency
% fs = max(2*fb_max,bw);                        % Sampling rate
fs = 2*fb_max;
fs = round(fs*tSw) / tSw;                       % sampling rate * sweep time must be integer

frequnecy_range = [fc fc+bw];                   % Frequency range [Hz]

VtargetMax = 20*1000/3600;
Max_doppler = (2*VtargetMax*fc)/c                % Maximum Speed []: ??

fmcwwav1 = phased.FMCWWaveform('SweepTime',tSw,'SweepBandwidth',bw,'SampleRate',fs);
sig = fmcwwav1();
Ns = numel(sig);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Define Interfering FMCW Radar Waveform

fcRdr2 = 24e9;                              % Center frequency (Hz)
lambdaRdr2 = freq2wavelen(fcRdr2,c);        % Wavelength (m)
rangeMaxRdr2 = 40.05;                       % Maximum range (m)
rangeResRdr2 = 0.15;                        % Desired range resolution (m)
vMaxRdr2 = vMax;                            % Maximum Velocity of cars (m/s)

tSwRdr2 = 5*range2time(rangeMaxRdr2,c);                 % Sweep Time [s]
bwRdr2 = rangeres2bw(rangeResRdr2,c);                   % Bandwith [Hz]
slopeRdr2 = bwRdr2/tSwRdr2;                             % Sweep slope 
fr_maxRdr2 = range2beat(rangeMaxRdr2,slopeRdr2,c);      % Maximum frequency corresponding to the maximum range [Hz]
fd_maxRdr2 = speed2dop(2*vMaxRdr2,lambdaRdr2);          % Maximum frequency corresponding to the Doppler Shift [Hz]
fb_maxRdr2 = fr_maxRdr2+fd_maxRdr2;                     % Maximum beat frequency
% fsRdr2 = max(2*fb_maxRdr2,bwRdr2);                    % Sampling rate
fsRdr2 = 2*fb_maxRdr2;
fsRdr2 = round(fsRdr2*tSwRdr2) / tSwRdr2;                       % sampling rate * sweep time must be integer

frequnecy_rangeRdr2 = [fcRdr2 fcRdr2+bwRdr2];                   % Frequency range [Hz]

Max_dopplerRdr2 = (c*fsRdr2)/(4*fcRdr2);                        % Maximum Speed []: ??


fmcwwav2 = phased.FMCWWaveform('SweepTime',tSwRdr2,'SweepBandwidth',bwRdr2, 'SampleRate',fsRdr2);
sigRdr2 = fmcwwav2();

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Plot the victim radar and interfering radar waveforms
% figure(1)
% pspectrum(repmat(sig,3,1),fmcwwav1.SampleRate,'spectrogram', ...
%     'Reassign',true,'FrequencyResolution',10e6)
% axis([0 15 -80 80]); title('Victim Radar'); 
% figure(2)
% pspectrum(repmat(sigRdr2,3,1),fmcwwav1.SampleRate,'spectrogram',...
%     'Reassign',true,'FrequencyResolution',10e6,'MinThreshold',-13)
% axis([0 15 -80 80]); title('Interfering Radar'); 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Model Victim and Interfering Radar Transceivers

% Construct the array for Victim Radar
Nve = 2*3;                                      % Number of virtual antenna with 2Tx, 3Rx
antElmnt = phased.IsotropicAntennaElement('BackBaffled',false);
rxVArray = phased.ULA('Element',antElmnt,'NumElements',Nve,'ElementSpacing',freq2wavelen(fc)/2,'Taper',taylorwin(Nve));

antAperture = 6.06e-4;                          % Antenna aperture (m^2)
antGain = aperture2gain(antAperture,lambda);    % Antenna gain (dB)
txPkPower = db2pow(13)*1e-3;                    % Tx peak power (W)
rxNF = 4.5;                                     % Receiver noise figure (dB)
% Model Victim Radar Transceiver
vradar = radarTransceiver("Waveform",fmcwwav1);
vradar.Transmitter = phased.Transmitter('PeakPower',txPkPower,'Gain',antGain);
vradar.TransmitAntenna = phased.Radiator('Sensor',antElmnt,'OperatingFrequency',fc);
vradar.ReceiveAntenna = phased.Collector('Sensor',rxVArray,'OperatingFrequency',fc);
vradar.Receiver = phased.ReceiverPreamp('Gain',antGain,'NoiseFigure',rxNF,'SampleRate',fmcwwav1.SampleRate);



% Define Victim Radar
VictimRadar = radarTransceiver('Waveform',fmcwwav1,'Transmitter',vradar.Transmitter, ...
    'TransmitAntenna',vradar.TransmitAntenna,'ReceiveAntenna',vradar.ReceiveAntenna ,'Receiver',vradar.Receiver);

% Model Interfering Radar Transceiver
iradar = clone(vradar);
iradar.Waveform = fmcwwav2;
vradar.Receiver = phased.ReceiverPreamp('Gain',antGain,'NoiseMethod','Noise power','NoisePower',0);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Model Free Space Propagation Channel for Interference

% Create a multipath channel for interfering radar
channelIradar = helperCreateMultiPathChannel('SampleRate',fmcwwav2.SampleRate, ...
    'OperatingFrequency',fcRdr2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Driving scenario
% Create driving scenario
[scenario, egoCar] = helperAutoDrivingScenario;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Radar simulation

Nft = fmcwwav1.SweepTime*fmcwwav1.SampleRate; % Number of fast-time samples
Nft = round(Nft);
Nsweep = 256;                                 % Number of slow-time samples
% Nr = 2^nextpow2(Nft);                         % Number of range samples
% Nd = 2^nextpow2(Nsweep); 
Nr = 512;
Nd = 256;

snr = 0.001;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Domain Calculate
% GPT Code

% Define the range and Doppler axes
rangeAxis = linspace(0, rangeMax, Nr/2);


dopplerAxis = linspace(-VtargetMax, VtargetMax, Nd);

% Doppler axis in m/s
dopplerAxis_mps = dopplerAxis * lambda / 2;

% For visualization in km/h, convert m/s to km/h
dopplerAxis_kmph = dopplerAxis_mps * 3.6;

%%
% Run the simulation loop
sweepTime = fmcwwav1.SweepTime;

while advance(scenario)

    % Get the current scenario time
    time = scenario.SimulationTime;

    % Get current target poses in ego vehicle's reference frame
    tgtPoses = targetPoses(egoCar);
    actProf = actorProfiles(scenario);

    % Assemble data cube at current scenario time
    Xcube = zeros(Nft,Nve,Nsweep);
    XcubeVradar = zeros(Nft,Nve,Nsweep);

    % Generate waveform signal to transmit for interfering radar
    sIsig = step(iradar.Waveform);

    % Match the length of interfering signal with victim sweep
    sIsig = repmat(sIsig,ceil(Nft/length(sIsig)),1);

    for m = 1:Nsweep
        sIcar = sIsig(1:Nft);
        sIsig = circshift(sIsig,length(sIsig)-Nft);

        % Generate all paths to victim radar
        [vpaths, ipaths] = helperGeneratePaths(tgtPoses,actProf,[lambda lambdaRdr2]);  % Global to local coord system(s) inside

        % Send signal from interfering radar to power amplifier for transmit
        IcarSig = step(iradar.Transmitter,sIcar);

        % Radiate signal from antenna elements of interfering radar
        xtIcar = step(iradar.TransmitAntenna,IcarSig,[ ipaths(1).AngleOfDeparture ipaths(2).AngleOfDeparture]);

        % Propagate signal along paths
        [xrxIradar,angIradar ]= step(channelIradar, xtIcar, ipaths,time);

        % Pass the signals through the receive antenna array elements
        xcIradar = step(vradar.ReceiveAntenna, xrxIradar, angIradar);

        % Pass the signals at the output of the receive antenna element to the receiver amplifier
        rxInt = step(vradar.Receiver,xcIradar);

        % Get the received signal without interference
        rxVictim = VictimRadar(vpaths,time);
        % Dechirp the received signal
        rxVsig = dechirp(rxVictim,sig);
        % Save sweep to data cube - create ADC data
        XcubeVradar(:,:,m) = awgn(rxVsig, snr, 'measured');

        % Dechirp the received signal + Interference and save sweep to data cube
        rx = dechirp(rxVictim+rxInt,sig);
        Xcube(:,:,m) = awgn(rx, snr, 'measured');

        % Get the current scenario time
        time = time + fmcwwav1.SweepTime;

        % Move targets forward in time for next sweep
        tgtPoses(1).Position=[tgtPoses(1).Position]+[tgtPoses(1).Velocity]*fmcwwav1.SweepTime;
        tgtPoses(2).Position=[tgtPoses(2).Position]+[tgtPoses(2).Velocity]*fmcwwav1.SweepTime;
    end
%     break;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Sig Process for interfered signal
% Perform 2D FFT on the interfered signal cube
Xcube_tmp = Xcube(:,1,:);
Xcube_tmp = squeeze(Xcube_tmp);

Xcube_tmp2 = Xcube_tmp(1:Nr,1:Nd);

%fft, fftshift
Xcube_fft = fftshift(fft2(Xcube_tmp2, Nr, Nd), 2);
% fft calibration
Xcube_fft = Xcube_fft(1:Nr/2, :);

%% Sig Process for Original signal
% Perform 2D FFT on the non-interfered signal cube
XcubeVradar_tmp = XcubeVradar(:,1,:);
XcubeVradar_tmp = squeeze(XcubeVradar_tmp);

%Nr, Nd만큼만 남기고 나머지는 버림
XcubeVradar_tmp2 = XcubeVradar_tmp(1:Nr,1:Nd);

%fft, fft shift
XcubeVradar_fft = fftshift(fft2(XcubeVradar_tmp2, Nr, Nd), 2);
% fft calibration
XcubeVradar_fft = XcubeVradar_fft(1:Nr/2, :);

%% Both signal Process
% Compute the magnitude of the FFT results
Xcube_fft_mag = abs(Xcube_fft);
XcubeVradar_fft_mag = abs(XcubeVradar_fft);

% Scaling for better visualization
Xcube_fft_mag_db = 10*log10(Xcube_fft_mag);
XcubeVradar_fft_mag_db = 10*log10(XcubeVradar_fft_mag);

%% Visualization with surf
% Interfered Signal Range-Doppler Map
figure;
surf(dopplerAxis, rangeAxis, Xcube_fft_mag_db);
shading interp;
xlabel('Doppler (Hz)');
ylabel('Range (m)');
zlabel('Magnitude (dB)');
title('Interfered Signal Range-Doppler Map');
colorbar;

% Non-interfered Signal Range-Doppler Map
figure;
% surf(XcubeVradar_fft_mag_db);
surf(dopplerAxis, rangeAxis, XcubeVradar_fft_mag_db);
shading interp;
xlabel('Doppler (m/s)');
ylabel('Range (m)');
zlabel('Magnitude (dB)');
title('Non-interfered Signal Range-Doppler Map');
colorbar;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end

https://kr.mathworks.com/help/radar/ug/simulate-fmcw-interference-between-automotive-radars.html

The example code seems to be changed.

so i leave the original code that i downloaded before

https://drive.google.com/file/d/1M4ahfFvZHQKs1uqF2biRSOqrbrbKo5Yk/view?usp=sharing

Answer 1

[akshatsood]

Dear @YongBin Lee,

I understand that you are trying to implement Doppler Shift in your FMCW Radar simulation code and are having difficulties understanding the example code. Please find below the response as per my understanding:

What part of the code is responsible for the signal transmission of "vradar"?

The signal transmission for "vradar" is handled by the "radarTransceiver" object and its associated methods. Specifically, the transmission occurs within the while loop where you call

rxVictim = VictimRadar(vpaths, time);

Here, "VictimRadar" is the "radarTransceiver" object for vradar, and this line of code essentially simulates the reception of the transmitted FMCW signal after it has interacted with the environment.

What part of the code should I modify to add some Doppler shift to the vradar signal?

To add Doppler shift to the "vradar" signal, you need to modify the signal processing part where you handle the received signal. Doppler shift is a frequency shift that occurs due to relative motion between the radar and the target.You can introduce Doppler shift by modifying the "rxVictim" signal before the dechirping process.

relativeVelocity = 10; % RELATIVE VELOCITY
fd = (2 * relativeVelocity * fc) / c; % Doppler Shift FREQ
t = (0:length(rxVictim)-1) / fs;

rxVictim = rxVictim .* exp(1j * 2 * pi * fd * t).';

Insert this code snippet right after the following line in your simulation loop.

rxVictim = VictimRadar(vpaths, time);

I hope this helps!