基于RD算法的点目标SAR成像，用于单点或多点目标的距离徙动校正，距离压缩和方位压缩matlab.zip

%校正线性距离走动的RD算法 clear;clc;close all; %% Parameter C=3e8; %光速 Fc=1e9; %载波频率 1GHz lambda=C/Fc; %波长 Xmin=0; %目标方位向位置 Xmax=100; Yc=10000; %目标距离向中心位置 Y0=500; %目标距离向位置[Yc-Y0,Yc+Y0] %imaged width 2*Y0 V=100; %机体速度 H=5000; %雷达高度 5000 m R0=sqrt(Yc^2+H^2); %雷达轨迹到场景中心线垂直距离 D=8; %天线真实尺寸 theta0=12/180*pi; %斜视角为4° Ls=lambda*R0/D/(cos(theta0))^2; %合成孔径长度 Ts=Ls/V; %合成孔径时间 %% Parameter--慢时间域 Ka=-2*V^2*(cos(theta0))^2/lambda/R0; %多普勒调频斜率 Bd=abs(Ka*Ts); %多普勒带宽 PRF=Bd; %脉冲重复频率 PRT=1/PRF; %脉冲重复时间 ds=PRT; %慢时间域采样间隔 Ns=ceil((Xmax-Xmin+Ls)/V/ds); %慢时间域采样数 Ns=2^nextpow2(Ns); %采样数表示为2的指数 ts=linspace((Xmin-Ls/2)/V,(Xmax+Ls/2)/V,Ns); %离散时间序列 PRT=(Xmax-Xmin+Ls)/V/Ns; PRF=1/PRT; ds=PRT; theta=atan(abs(V*ts)./R0); fa=2*V/lambda*sin(theta)*cos(theta0); fam=2*V/lambda*cos(theta0); %% Parameter--快时间域 Tr=5e-6; %脉宽 5us Br=30e6; %chrip信号调频带宽 30MHz Kr=Br/Tr; %chirp线性调频斜率 Fsr=3*Br; %采样频率 dt=1/Fsr; %采样间隔 Rmin=sqrt((Yc-Y0)^2+H^2)/cos(theta0); Rmax=sqrt((Yc+Y0)^2+H^2)/cos(theta0); Nf=ceil(2*(Rmax-Rmin)/C/dt+Tr/dt); %快时间域采样数 Nf=2^nextpow2(Nf); %采样数表示为2的指数 tf=linspace(2*Rmin/C,2*Rmax/C+Tr,Nf); %离散时间序列 dt=(2*Rmax/C+Tr-2*Rmin/C)/Nf; Fsr=1/dt; fr=linspace(-Fsr/2,Fsr/2,Nf); %% Parameter--分辨率 DY=C/2/Br; %距离向分辨率 DX=D/2; %角分辨率 %% Parameter--点目标 Ntarget=2; %目标数 Ptarget=[Xmin,Yc,1 Xmin+10*DX,Yc+40*DY,1 ]; %目标位置(注意合成孔径长度) disp('Parameters:') disp('快时间域采样频率');disp(Fsr/Br) disp('快时间域采样数');disp(Nf) disp('慢时间域采样频率');disp(PRF/Bd) disp('慢时间域采样数');disp(Ns) disp('距离向分辨率');disp(DY) disp('方位向分辨率');disp(DX) disp('SAR合成孔径长度');disp(Ls) disp('点目标位置');disp(Ptarget) %% 原始信号 K=Ntarget; %目标数 N=Ns; %慢时间域向量数 M=Nf; %快时间域向量数 T=Ptarget; %目标位置 Snm=zeros(N,M); %采样空间 for k=1:1:2 sigma=T(k,3); xn=T(k,1)-abs(tan(theta0)*sqrt(H^2+T(k,2)^2)); Ds=ts*V-T(k,1)+abs(tan(theta0)*sqrt(H^2+T(k,2)^2)); Rb=sqrt(H^2+T(k,2)^2); %点目标到航迹垂直距离 Rb0=Rb/cos(theta0); %点目标到雷达波束中心斜距 % R=sqrt((Ds).^2+T(k,2)^2+H^2); % R=Rb+1/8*(lambda/V)^2*Rb*fa.^2; % R=sqrt(Rb0^2+V^2*cos(theta0)^2*(ts-xn/V).^2)-(ts-xn/V)*V*sin(theta0); R=Rb-sin(theta0)*Ds+Ds.^2*cos(theta0)^2/2/Rb; tau=2*R/C; Df=ones(N,1)*tf-tau'*ones(1,M); phase=pi*Kr*Df.^2-(4*pi/lambda)*(R'*ones(1,M)); %时间域相位 Snm=Snm+sigma*exp(1i*phase).*(0<Df&Df<Tr).*((abs(Ds)<(Ls/2+abs(tan(theta0)*sqrt(H.^2+T(k,2).^2))))'*ones(1,M)); end %% RCM补偿 fdc=-2*V*sin(theta0)/lambda; for i=1:N S1(i,:)=fftshift(fft(fftshift(Snm(i,:)))); %对信号每一行即距离向fft end tr=tf-2*Rmin/C; Refr=exp(1i*pi*Kr*tr.^2).*(0<tr&tr<Tr); %时域匹配 S2=ifty(S1.*(ones(N,1)*conj(fty(Refr))).*exp(-1i*1/4*pi/C*sin(theta0).*Ds.'*fr));%%%补偿距离的线性相位因子 Gr=abs(S2); %% 方位向压缩 for i=1:M S3(:,i)=fftshift(fft(fftshift(S2(:,i)))); %对信号每一列即方位向fft end ta=ts-Xmin/V; Refa=exp(1i*pi*Ka*ta.^2+1i*2*pi*fdc*ta).*(abs(ta)<Ts/2); S4=iftx(S3.*(conj(ftx(Refa)).'*ones(1,M))); Ga=abs(S4); figure(3) subplot(321); imagesc(abs(Snm)); subplot(322); imagesc(abs(S1)); subplot(323); imagesc(abs(S2)); subplot(324); imagesc(abs(S3)); subplot(325); imagesc(Gr); subplot(326); imagesc(Ga); %% 信号强度图像 figure(1) subplot(211); colormap(gray); row=tf*C/2-1800; col=ts*V-50; imagesc(row,col,255-Gr); %intensity image of Sr axis([Yc-Y0,Yc+Y0,Xmin-Ls/2,Xmax+Ls/2]); xlabel('距离向');ylabel('方位向'); title('距离压缩后'); subplot(212); imagesc(row,col,255-Ga); %intensity image of Sa axis([Yc-Y0,Yc+Y0,Xmin-Ls/2,Xmax+Ls/2]); xlabel('距离向');ylabel('方位向'); title('方位压缩后'); %% 点目标3D图 figure(2) mesh(Ga((100:180),(400:700)));axis tight xlabel('距离');ylabel('方位'); title('点目标成像');