# 【雷达通信】基于matlab NCP算法SAR回波生成和成像【含Matlab源码 1185期】

## 二、部分源代码

``````%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%对应于HOP的NCS成像算法
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
clear all;
clc;
c=299792458;                        %light speed
lamda = c/fc;
K=B/Tp;                        %调频率
dt=1/fs;                       %采样间隔
delaymin=2*Rmin/c-Tp/2;             %距离波门起始时刻
Na=rawdatawidth;
Nr=slicelen;
T=Tp;
lammda=lamda;
fx0=round(Na/2);
ft0=round(Nr/2);
cj=sqrt(-1);

%方位向傅立叶变换
sfa=fftshift(fft(sr,[],2),2);

%距离向FFt，三次相位滤波，距离向逆FFT
fx=([1:Na]-fx0)/(Na-1)*PRF;
fy=([1:Nr]-ft0)/(Nr-1)*fs;             %距离向频率
indexfx=find(abs(fx)<2*V/lammda&abs(fx)<4*V*sin(theta/2)/(c/(fc+B/2)));%有效多普勒位置
indexlen=length(indexfx);
fxm=2*V/lammda;
costheta=sqrt(1-(fx(indexfx)/fxm).^2);
sintheta=sqrt((fx(indexfx)/fxm).^2);
% Cs=1./sqrt(1-(lammda*fx(indexfx)/2).^2)-1; %CS因子
Cs=1./costheta-1;
% ksrc=2*lammda/c.^2*(lammda*fx(indexfx)/2).^2./(1-(lammda*fx(indexfx)/2).^2).^1.5;%距离弯曲因子
ksrc=(2*lammda/c^2)*sintheta.^2./costheta.^3;
% R=(delaymin+[0:Nr-1]'*dt)*c/2;
R=(Src+[1:slicelen]*dy-Tp/2*c/2)';
Ref=Src;                  %参考距离
Kref=K./(1-K*Ref*ksrc); %距离调频率
Yfa=(2-costheta).*(1+costheta)*lammda./costheta.^2./Kref/c/2;  %三次相位滤波系数

sfar1=fftshift(fft(sfa,[],1),1);%距离向FFT
H_TP_Filter=exp(cj*2/3*pi*(Yfa'*ones(1,slicelen)).*(ones(indexlen,1)*fy.^3));%三次相位滤波参考函数

phase_ref1=(4*pi*Ref*sqrt((ones(indexlen,1)*(fc+fy)/c).^2-(fx(indexfx).'/2/V).^2*ones(1,slicelen)))';
phase_ref2=ones(slicelen,1)*(-2*pi*Ref/V*sqrt(fxm.^2-fx(indexfx).^2));
phase_ref3=(-4*pi*Ref/c*(ones(indexlen,1)*fy).*((1+Cs)'*ones(1,slicelen)))';
phase_ref4=pi*Ref*2*lammda/c^2*(fy.^2).'*(((fx(indexfx).*lammda/2/V).^2).*(1+Cs).^2);
H_PhaseC=exp(cj*(phase_ref1+phase_ref2+phase_ref3+phase_ref4));
sfar1(:,indexfx)=sfar1(:,indexfx).*H_PhaseC.*H_TP_Filter';%三次相位滤波
sfa1=ifft(sfar1);%距离向IFFT

%NCS操作
q2=Kref.*(1-costheta)./costheta;       %NCS运算参考函数系数
q3=(lammda/2/c)*Kref.^2.*(1-costheta.^2)./costheta.^2.*(1-costheta)./costheta;%NCS运算参考函数系数
T_t=2/c*((R-Ref)*(1+Cs));
H_NCS=exp(cj*pi*(ones(slicelen,1)*q2.*T_t.^2)+...
cj*2/3*pi*(ones(slicelen,1)*q3.*T_t.^3));%NCS参考函数
sfa1(:,indexfx)=sfa1(:,indexfx).*H_NCS;
%距离向傅立叶变换
sfar2=fft(sfa1,[],1);
%sfar=(fftshift(fft(sfa,[],1),1));
% sfar=fft(sfa);
%距离压缩和二次距离压缩，距离迁移校正
phase1=pi*(fy.^2)'*(costheta./Kref);
phase2=4*pi*Ref/c*fy.'*((1-costheta)./costheta);
phase3=pi*lammda/3/c*(fy.^3).'*((((1-costheta).*(1-costheta.^2))+...
((1+costheta).*(2-costheta).*costheta))./Kref);
H_Rng_Cmp=exp(cj*(phase1+phase2+phase3));%距离压缩、迁移校正参考函数
sfar2(:,indexfx)=sfar2(:,indexfx).*H_Rng_Cmp;

% sfar(:,indexfx)=sfar(:,indexfx).*exp(cj*pi*(ft.^2*ones(1,indexlen))...
%                                 ./(ones(Nr,1)*(Kref.*(1+Cs))));%.*(hamming(Nr)*ones(1,indexlen));%距离压缩和二次距离压缩
% sfar(:,indexfx)=sfar(:,indexfx).*exp(cj*4*pi*Ref/c*ft*Cs);      %距离迁移校正

%距离逆傅立叶变换
sfr2=ifft(sfar2,[],1);
sfr=sfr2;
figure;imagesc(abs(sfr));

%方位向压缩
indext=[1:Nr]';
% indext=[1:ImgRng]';
R1=(indext*dt+delaymin)*c/2;
pres=4*pi/c.^2*(R1-Ref).^2*(Kref.*Cs.*(1+Cs));%剩余相位
sfr(indext,indexfx)=sfr(indext,indexfx).*exp(cj*pres);%剩余相位补偿
sfr(indext,indexfx)=sfr(indext,indexfx).*exp(cj*2*pi*R1*2/lammda...
*(sqrt(1-(lammda*fx(indexfx)/2/V).^2)-1));%.*(hamming(Nr)*ones(1,indexlen));%方位压缩

%方位向逆傅立叶变换
sfra1=ifft(sfr,[],2);
% sfra=sfra1(:,Azm_lower:Azm_upper);
figure; imagesc(abs(sfra1));
[pm,pn] = plotsection(sfra1.');
% Nt =pn;
kp_a = 2;
deltA = V/PRF*kp_a;
deltR = c/2/fs;
indexa =[1:Na].';
Mr = 64;
Ma = round(Mr*kp_a);
%幅度量化dB图
% my_plot(indexa(pm-Ma:kp_a:pm+Ma),indext(pn-Mr:pn+Mr),sfra1(pn-Mr:pn+Mr,pm-Ma:kp_a:pm+Ma),50);
% xlabel('方位向(点数)');ylabel('距离向(点数)');title('NCS成像结果幅度量化dB图');
figure(8);%彩色等高线图
contour(indext(pn-Mr:pn+Mr),indexa(pm-Ma:kp_a:pm+Ma),abs(sfra1(pn-Mr:pn+Mr,pm-Ma:kp_a:pm+Ma).'),40);colorbar;
function [m,n]= plotsection(y)
[x,b] = max(abs(y));
[z,n] = max(x);
[A,B] = size(y);
m = b(n);
figure;
plot(20*log10(abs(y(m,:))/z)),grid,axis([1,B,-40,0]);
title('距离向剖面图');
figure;
plot(20*log10(abs(y(:,n))/z)),grid,axis([1,A,-40,0]);
title('方位向剖面图');