【雷達通信】基于matlab的雷達數字信號處理【含Matlab原始碼 281期】

一、簡介

利用MATLAB設計經典的雷達數字信號處理，該系統具備對雷達目標回波的處理能力，能夠從噪聲中將目標檢測出來，并提取目標的距離、速度、角度資訊，教程分五節完成，主要包括：
第一節，雷達LFM信號分析；
第二節，脈沖壓縮處理；
第三節，相參積累處理；
第四節，恒虛警CFAR處理；
第五節，目標資訊提取處理，
1.雷達LFM信號分析
脈沖壓縮雷達最常見的調制信號是線性調頻（Linear Frequency Modulation）信號,接收時采用匹配濾波器（Matched Filter）壓縮脈沖，脈沖壓縮雷達能同時提高雷達的作用距離和距離解析度，這種體制采用寬脈沖發射以提高發射的平均功率，保證足夠大的作用距離；而接受時采用相應的脈沖壓縮演算法獲得窄脈沖，以提高距離解析度，較好的解決雷達作用距離與距離解析度之間的矛盾，
LFM信號的數學運算式為：


2.脈沖壓縮處理
脈沖壓縮指雷達在發射時采用寬脈沖信號，接收和處理回波后輸出窄脈沖，脈沖壓縮技術是匹配濾波理論和相關接收理論的一個很好的實際應用，很好地解決了這樣的一個問題：在發射端發射大時寬、帶寬信號，以提高信號的發射能量，而在接收端，將寬脈沖信號壓縮為窄脈沖，以提高雷達對目標的距離分辨精度和距離分辨力，該技術解決了雷達遠距離探測與高精度測距性能不可兼顧的問題，是現代雷達中不可缺少的關鍵技術，
脈沖壓縮的DSP處理方法有時域相關或頻域相乘，對于點數較多的回波信號，采用頻域相乘方法可以獲得較快的運算速度，頻域脈沖壓縮的原圖如下圖所示，

DSP對采樣后的資料進行FFT變換，變換至頻域后，與其匹配濾波器頻率資料進行復數相乘，相乘后，再與復數補償因子進行相乘解決脈沖間距離走動問題，最后將結果做IFFT，重新變換回時域，其中，FFT點數、復數相乘點數、IFFT點數均為1024點，
信號的匹配濾波器的時域脈沖回應為：

3.相參積累處理
在信號理論中，相參又稱為相干，定義為脈沖之間存在確定的相位關系，簡單來說，脈沖間的相位可以互相對照，知道其中一個相位就有辦法知道另外一個，相參處理的意義在于脈沖積累時提高信噪比，提高多普勒頻率的準確度，由于雷達回波信號不但有微弱的信號，還會有很強的噪聲，雷達的主要目的就是要把微弱的目標信號從噪聲中分離出來，設法提高信噪比，要想把信號提取出來，必須要將信號放大，但放大的同時噪聲也被放大，因為它們總是同時存在的，并且放大電路自己本身也有噪聲，放大后信號與噪聲的比值反而變小了，更不利于提取有用的回波信號，解決的方法是進行相參積累，可以對n個回波進行累加，由于噪聲是隨機的，累加的結果是信號變強，而噪聲因是隨機的，強度反而變小，這樣信號與噪聲比就提高了，相參積累中多個脈沖之間相位關系固定且明確，所以理論上積累后信噪比可提高到n倍，
下圖為相參積累原理示意圖：


4.恒虛警CFAR處理
恒虛警率CFAR是ConstantFalse-AlarmRate的縮寫，在雷達信號檢測中，當外界干擾強度變化時，雷達能自動調整其靈敏度，使雷達的虛警概率保持不變，這種特性稱為恒虛警率特性，恒虛警率檢測是雷達目標自動檢測的一個重要組成部分，作為從雷達中提取目標的第一步，是進一步識別目標的基礎，虛警率是指偵察設備在單位時間內將噪聲或其他干擾信號誤判為威脅輻射源信號的概率，而恒虛警率檢測則證明了檢測演算法的穩定性和可靠性，
雷達信號恒虛警率檢測就是要求虛警概率保持恒定，這主要是因為在雷達信號檢測中，信號的最佳檢測通常采用奈曼-皮爾遜準則，即在保持恒定的虛警概率的條件下，使正確檢測的概率達到最大值，

精細的恒虛警檢測，
5.目標資訊提取處理
本文介紹的雷達系統采用單脈沖體制，具備精密跟蹤的能力，每發射一個脈沖，天線能同時形成若干個波束，將各波束回波信號的振幅和相位進行比較，當目標位于天線軸線上時，各波束回波信號的振幅和相位相等，信號差為零；當目標不在天線軸線上時，各波束回波信號的振幅和相位不等，產生信號差，驅動天線轉向目標直至天線軸線對準目標，這樣可測出目標的方位角與俯仰角，從各波束接收的信號之和，可測出目標的距離，從而實作對目標的測量和跟蹤功能，單脈沖雷達已經廣泛應用，在軍事上主要用于目標識別、靶場精密跟蹤測量、導彈預警和跟蹤、導彈再入彈道測量、火箭和衛星跟蹤、武器火力控制、炮位偵查、地形跟隨、導航、地圖測繪等，在民用上主要用于交通管制，
單脈沖自動測角屬于同時波瓣測角法，在一個角平面內，兩個相同的波束部分重疊，交疊方向即為等信號軸的方向，將這兩個波束接收到的回波信號進行比較，就可取得目標在這個平面上的角誤差信號，然后將此誤差電壓放大變換后加到驅動電動機控制天線向減小誤差的方向運動，因為兩個波束同時接收到回波，故單脈沖測角獲得目標角誤差資訊的時間可以很短，理論上只要分析一個回波脈沖就可以確定角誤差，所以叫“單脈沖”，這種方法可以獲得很高的測角精度，故精密跟蹤雷達通常采用它，
由于取出角度誤差信號的具體方法不同，單脈沖雷達的種類很多，應用最廣的是振幅和差式單脈沖雷達，該方法的實質實際上是利用兩個偏置天線方向圖的和差波束，

二、源代碼

%% 線性調頻信號的產生
clc;close all;clear all;

T = 10e-6;                          %脈寬
B = 25e6;                           %信號帶寬
K = B/T;                            %調頻斜率
Fs = 2*B;Ts = 1/Fs;                 %采樣頻率及采樣間隔
N = T/Ts;
t = linspace(-T/2,T/2,N);
St = exp(j*pi*K*t.^2);            %生成線性調頻信號

figure;
subplot(211)
plot(t*1e6,St);
xlabel('Time in u sec');
title('線性調頻信號');
grid on;axis tight;
subplot(212)
freq = linspace(-Fs/2,Fs/2,N);
plot(freq*1e-6,fftshift(abs(fft(St))));
xlabel('Frequency in MHz');
title('線性調頻信號的幅頻特性');
grid on;axis tight;
%% CFAR

% 先運行 processCoherent.m檔案

%%
%%%%%%%%%二維篩選%%%%%%%%%
N_ref_2D = M/4;  %參考區域大小
PC_data_ifft_CA_abs = abs(PC_data_ifft_CA);

N_point = length(PC_data_ifft_CA_abs);
%計算區域1的和
sum1 = sum(PC_data_ifft_CA_abs(1:N_ref_2D,1:N_ref_2D)); %按行相加
sum_ref_2D(1,1) = sum(sum1(1,:));            %按列相加
%計算區域2的和
sum2 = sum(PC_data_ifft_CA_abs(1:N_ref_2D,(N_point-N_ref_2D+1):N_point));
sum_ref_2D(1,2) = sum(sum2(1,:));
%計算區域3的和
sum3 = sum(PC_data_ifft_CA_abs((M-N_ref_2D+1):M,1:N_ref_2D));
sum_ref_2D(1,3) = sum(sum3(1,:));
%計算區域4的和
sum4 = sum(PC_data_ifft_CA_abs((M-N_ref_2D+1):M,(N_point-N_ref_2D+1):N_point));
sum_ref_2D(1,4) = sum(sum4(1,:));

%%%產生門限Threshold%%%
SNR_Threshold = 15;   % 檢測門限 dB
Threshold = ((sum(sum_ref_2D(:,1:4)) - max(sum_ref_2D))/(N_ref_2D^2*3))*10^(SNR_Threshold/10);

%%%尋找 PC_data_ifft_CA 中大于門限的點并記錄其位置%%%
location = zeros(M,N_point);  %填1表示該點大于門限
for i=1:1:M
    for j=1:1:N_point
        if(PC_data_ifft_CA_abs(i,j) >= (1*Threshold))
            location(i,j) = 1;
        else
            location(i,j) = 0;
        end
    end
end
% figure(4),mesh(location),title('二維篩選后的結果/location');axis tight;

%%%篩選結果顯示%%%
selected = zeros(M,N_point);
for i=1:1:M
    for j=1:1:N_point
        if(location(i,j) == 1)
            selected(i,j) = PC_data_ifft_CA_abs(i,j);
        else
            selected(i,j) = 0;
        end
    end
end
figure(5),mesh(t_x_ca_cut,f_x_ca,selected),title('二維CFAR后的結果(信號1的SNR=-5,信號2的SNR=2)/selected');axis tight;
xlabel('時間/s','FontSize',12);ylabel('頻率/Hz','FontSize',12);zlabel('檢測結果(尖峰處為目標)','FontSize',12);
%%%%%%%%%%二維篩選結束%%%%%%%%%%%

%%%%%%%%%%一維篩選%%%%%%%%%%
%%%%距離維篩選%%%%
N_prot_R = 30; %保護單元覆寫第二旁瓣，次大點與次次大點
N_ref_R = 135;  %參考單元
threshold_R = zeros(M,N_point); %記錄距離維每個觀察點的門限值
for i=1:1:M
    for j=1:1:N_point
        if(location(i,j) == 1) %只篩選經過二維篩選之后的點
            if(j <= (N_prot_R + N_ref_R) ) %觀察點過于偏左，左側點數不足
                threshold_R(i,j) = (sum(PC_data_ifft_CA_abs(i,(j+N_prot_R+1):(j+N_prot_R+N_ref_R))))/N_ref_R; %觀察點右側參考點相加求平均
            elseif( j >= (N_point - N_prot_R - N_ref_R + 1)) %觀察點過于偏右，右側點數不足
                threshold_R(i,j) = (sum(PC_data_ifft_CA_abs(i,(j-N_prot_R-N_ref_R):(j-N_prot_R-1))))/N_ref_R; %觀察點左側參考點求平均
            else %觀察點居中，左右點數均足夠
                sum_R_left = sum(PC_data_ifft_CA_abs(i,(j-N_prot_R-N_ref_R):(j-N_prot_R-1)));
                sum_R_right = sum(PC_data_ifft_CA_abs(i,(j+N_prot_R+1):(j+N_prot_R+N_ref_R)));
                threshold_R(i,j) = max(sum_R_left,sum_R_right)/(N_ref_R); %觀察點兩側取大的一側求均值
            end
        else
            threshold_R(i,j) = 0;
        end
    end
end
%利用距離維門限進行一維篩選
location_R = zeros(M,N_point);
for i=1:1:M
    for j=1:1:N_point
        if(location(i,j) == 1) %只篩選經過二維篩選之后的點
            if(PC_data_ifft_CA_abs(i,j) >= (1*threshold_R(i,j)))
                location_R(i,j) = 1;
            else
                location_R(i,j) = 0;
            end
        else
            location_R(i,j) = 0;
        end
    end
end

三、運行結果

四、備注

完整代碼或者代寫添加QQ1564658423
往期回顧>>>>>>
【信號處理】基于matlab HMM的睡眠狀態檢測【含Matlab原始碼 050期】
【信號處理】基于matlab CDR噪聲和混響抑制【含Matlab原始碼 051期】
【信號處理】基于matlab最小二乘法解決稀疏信號恢復問題【含Matlab原始碼 052期】
【信號處理】基于matlab小波變換的音頻水印嵌入提取【含Matlab原始碼 053期】
【信號處理】基于matlab ICA演算法信號分離【含Matlab原始碼 054期】
【信號處理】基于matlab GUI界面的脈搏信號之脈率存檔【含Matlab原始碼 237期】
【信號處理】基于matlab GUI界面的虛擬信號發生器（各種波形）【含Matlab原始碼 271期】
【信號處理】基于matlab GUI界面信號發生器之電子琴【含Matlab原始碼 272期】
【信號處理】基于matlab的數字電子琴設計與實作【含Matlab原始碼 273期】