|
一维ca-cfar(均值类恒虚警检测)与二维ca-cfar算法matlab实现
- % Implement 1D CFAR using lagging cells on the given noise and target scenario.
- % Close and delete all currently open figures
- close all;
- % Data_points
- Ns = 1000;
- % Generate random noise
- s=abs(randn(Ns,1));
- %Targets location. Assigning bin 100, 200, 300 and 700 as Targets with the amplitudes of 8, 9, 4, 11.
- s([100 ,200, 350, 700])=[8 15 7 13];
- %plot the output
- figure ('Name','original_signal'),plot(s);
- % TODO: Apply CFAR to detect the targets by filtering the noise.
- % 1. Define the following:
- % 1a. Training Cells
- T = 12;
- % 1b. Guard Cells
- G = 4;
- % Offset : Adding room above noise threshold for desired SNR
- offset=5;
- % Vector to hold threshold values
- threshold_cfar = []; % 存储阈值的数组
- %Vector to hold final signal after thresholding
- signal_cfar = []; % 是否有目标?
- % 2. Slide window across the signal length
- for i = 1:(Ns-(G+T+1))
- % 2. - 5. Determine the noise threshold by measuring it within the training cells
-
- noise_level =sum(s(i:i+T-1));
-
- % 6. Measuring the signal within the CUT
-
- threshold = (noise_level/T)*offset; % 这里是乘以offset
- threshold_cfar=[threshold_cfar,{threshold}];
-
- signal=s(i+T+G);
-
- % 8. Filter the signal above the threshold
- if (signal<threshold)
- signal=0; % 不是目标
- end
- signal_cfar = [signal_cfar, {signal}];
- end
- % plot the filtered signal
- % signal_cfar2 = cell2mat(signal_cfar);
- figure,plot (cell2mat(signal_cfar),'g--');
- % plot original sig, threshold and filtered signal within the same figure.
- figure,plot(s);
- hold on,plot(cell2mat(circshift(threshold_cfar,G)),'r--','LineWidth',2)
- hold on, plot (cell2mat(circshift(signal_cfar,(T+G))),'g--','LineWidth',4);
- legend('Ori_Signal','CFAR Threshold','detection')
复制代码
- clear all
- close all;
- clc;
- %% Radar Specifications
- %%%%%%%%%%%%%%%%%%%%%%%%%%%
- % Frequency of operation = 77GHz
- % Max Range = 200m
- % Range Resolution = 1 m
- % Max Velocity = 100 m/s
- %%%%%%%%%%%%%%%%%%%%%%%%%%%
- max_range=200;
- c = 3e8;
- range_resolution= 1;
- %Operating carrier frequency of Radar
- fc= 77e9; %carrier freq
- wavelength = c/fc;
- %% User Defined Range and Velocity of target
- % *%TODO* :
- % define the target's initial position and velocity. Note : Velocity
- % remains contant
- target_pos=100;
- target_speed=30;
- %% FMCW Waveform Generation
- % *%TODO* :
- %Design the FMCW waveform by giving the specs of each of its parameters.
- % Calculate the Bandwidth (B), Chirp Time (Tchirp) and Slope (slope) of the FMCW
- % chirp using the requirements above.
- B_sweep = c/(2*range_resolution); %Calculate the Bandwidth (B)
- T_chirp = 5.5*2*max_range/c; % 为何啊?
- slope=B_sweep/T_chirp;
-
- %The number of chirps in one sequence. Its ideal to have 2^ value for the ease of running the FFT
- %for Doppler Estimation.
- Nd=128; % #of doppler cells OR #of sent periods % number of chirps
- %The number of samples on each chirp.
- Nr=1024; %for length of time OR # of range cells
- % Timestamp for running the displacement scenario for every sample on each
- % chirp
- t=linspace(0,Nd*T_chirp,Nr*Nd); %total time for samples
- %Creating the vectors for Tx, Rx and Mix based on the total samples input.
- Tx=zeros(1,length(t)); %transmitted signal
- Rx=zeros(1,length(t)); %received signal
- Mix = zeros(1,length(t)); %beat signal
- %Similar vectors for range_covered and time delay.
- r_t=zeros(1,length(t)); % 每时每刻的距离
- td=zeros(1,length(t)); % 距离算出来的时延
- %% Signal generation and Moving Target simulation
- % Running the radar scenario over the time.
- doppler_resolution = 1 / (Nd * T_chirp);
- v_resolution = (wavelength*doppler_resolution)/2; % 速度分辨率为2.07 为何?
- for i=1:length(t)
-
-
- % *%TODO* :
- %For each time stamp update the Range of the Target for constant velocity.
-
- r_t(i) = target_pos+(target_speed*t(i));
- td(i) = 2*r_t(i)/c; % Time delay
-
- % *%TODO* :
- %For each time sample we need update the transmitted and
- %received signal.
- Tx(i) = cos(2 * pi * (fc * t(i) + slope * (t(i)^2)/2));
- Rx(i) = cos(2 * pi * (fc * (t(i) - td(i)) + slope * ((t(i)-td(i))^2)/2));
-
- % *%TODO* :
- %Now by mixing the Transmit and Receive generate the beat signal
- %This is done by element wise matrix multiplication of Transmit and
- %Receiver Signal
- Mix(i) = Tx(i) .* Rx(i); % 不是发射信号的复共轭 不是复信号啊
-
- end
- %% RANGE MEASUREMENT
- % *%TODO* :
- %reshape the vector into Nr*Nd array. Nr and Nd here would also define the size of
- %Range and Doppler FFT respectively.
- Mix_reshape = reshape(Mix,[Nr,Nd]);
- % *%TODO* :
- %run the FFT on the beat signal along the range bins dimension (Nr) and
- %normalize.
- sig_fft1_0=fft(Mix_reshape,Nr)./Nr;
- % *%TODO* :
- % Take the absolute value of FFT output
- sig_fft1_1=abs(sig_fft1_0);
- % *%TODO* :
- % Output of FFT is double sided signal, but we are interested in only one side of the spectrum.
- % Hence we throw out half of the samples.
- sig_fft1=sig_fft1_1(1:Nr/2, 1); % 这个切片从矩阵(1024*128)变向量(1*512)了 只取一列?
- %plotting the range
- figure ('Name','Range from First FFT')
- subplot(1,1,1)
- % *%TODO* :
- % plot FFT output
- plot(sig_fft1);
- %plot(range_resolution*Nr/2, sig_fft1);
- axis ([0 200 0 1]);
- %% RANGE DOPPLER RESPONSE
- % The 2D FFT implementation is already provided here. This will run a 2DFFT
- % on the mixed signal (beat signal) output and generate a range doppler
- % map.You will implement CFAR on the generated RDM
- % Range Doppler Map Generation.
- % The output of the 2D FFT is an image that has reponse in the range and
- % doppler FFT bins. So, it is important to convert the axis from bin sizes
- % to range and doppler based on their Max values.
- Mix=reshape(Mix,[Nr,Nd]);
- % 2D FFT using the FFT size for both dimensions.
- sig_fft2 = fft2(Mix,Nr,Nd);
- % Taking just one side of signal from Range dimension.
- sig_fft2 = sig_fft2(1:Nr/2,1:Nd);
- sig_fft2 = fftshift (sig_fft2);
- RDM = abs(sig_fft2);
- RDM = 10*log10(RDM) ;
- %use the surf function to plot the output of 2DFFT and to show axis in both
- %dimensions
- doppler_axis = linspace(-100,100,Nd);
- range_axis = linspace(-200,200,Nr/2)*((Nr/2)/400);
- figure ('Name','2D FFT'),surf(doppler_axis,range_axis,RDM);
- xlabel('Speed');
- ylabel('Range');
- zlabel('db of rdm');
- %% CFAR implementation
- %Slide Window through the complete Range Doppler Map
- % *%TODO* :
- %Select the number of Training Cells in both the dimensions.
- Tr=10;
- Td=8;
- % *%TODO* :
- %Select the number of Guard Cells in both dimensions around the Cell under
- %test (CUT) for accurate estimation
- Gr=4;
- Gd=4;
- % *%TODO* :
- % offset the threshold by SNR value in dB
- off_set=1.4;
- % *%TODO* :
- %design a loop such that it slides the CUT across range doppler map by
- %giving margins at the edges for Training and Guard Cells.
- %For every iteration sum the signal level within all the training
- %cells. To sum convert the value from logarithmic to linear using db2pow
- %function. Average the summed values for all of the training
- %cells used. After averaging convert it back to logarithimic using pow2db.
- %Further add the offset to it to determine the threshold. Next, compare the
- %signal under CUT with this threshold. If the CUT level > threshold assign
- %it a value of 1, else equate it to 0.
- % Use RDM[x,y] as the matrix from the output of 2D FFT for implementing
- % CFAR
- RDM = RDM/max(max(RDM)); % Normalizing
- % *%TODO* :
- % The process above will generate a thresholded block, which is smaller
- %than the Range Doppler Map as the CUT cannot be located at the edges of
- %matrix. Hence,few cells will not be thresholded. To keep the map size same
- % set those values to 0.
- %Slide the cell under test across the complete martix,to note: start point
- %Tr+Td+1 and Td+Gd+1
- for i = Tr+Gr+1:(Nr/2)-(Tr+Gr) % 这里只需要一半
- for j = Td+Gd+1:(Nd)-(Td+Gd) % trip个数 慢时间维 多普勒维
- %Create a vector to store noise_level for each iteration on training cells
- noise_level = zeros(1,1);
- %Step through each of bins and the surroundings of the CUT
- for p = i-(Tr+Gr) : i+(Tr+Gr)
- for q = j-(Td+Gd) : j+(Td+Gd)
- %Exclude the Guard cells and CUT cells
- if (abs(i-p) > Gr || abs(j-q) > Gd)
- %Convert db to power 噪声水平求和
- noise_level = noise_level + db2pow(RDM(p,q));
- end
- end
- end
-
- %Calculate threshould from noise average then add the offset
- threshold = pow2db(noise_level/(2*(Td+Gd+1)*2*(Tr+Gr+1)-(Gr*Gd)-1));
- %Add the SNR to the threshold db值就是加offset
- threshold = threshold + off_set;
- %Measure the signal in Cell Under Test(CUT) and compare against
- CUT = RDM(i,j);
-
- if (CUT < threshold)
- RDM(i,j) = 0;
- else
- RDM(i,j) = 1;
- end
-
- end
- end
- RDM(RDM~=0 & RDM~=1) = 0;
- % *%TODO* :
- %display the CFAR output using the Surf function like we did for Range
- %Doppler Response output.
- % figure,surf(doppler_axis,range_axis,'replace this with output');
- % colorbar;
- figure('Name', 'CA-CFAR Filtered RDM')
- surf(doppler_axis,range_axis,RDM);
- colorbar;
- title( 'CA-CFAR Filtered RDM surface plot');
- xlabel('Speed');
- ylabel('Range');
- zlabel('Normalized Amplitude');
-
-
复制代码
|
|