clear all
thetaT=0;%平台波束斜视角
thetaT=thetaT*pi/180;%rad
c=3e8;%光速
fc=1.5e9;%载频
lambda=c/fc;%波长
%%测绘带区域
X0=200;%方位向[-X0,X0]
H=4000;%高度
Y=10000;
Rc=sqrt(H^2+Y^2)/cos(thetaT);
%Rc=5000;
R0=150;%距离向[Rc-R0,Rc+R0]
%%距离向(Range),r/t domain
Tr=5e-6;%LFM信号脉宽 5us (200m)
Br=50e6; %LFM信号带宽 50MHz
Kr=Br/Tr; %调频斜率
Nr=512;
r=Rc+linspace(-R0,R0,Nr);
t=2*r/c;%t域序列
dt=R0*4/c/Nr;%采样周期
f=linspace(-1/2/dt,1/2/dt,Nr);%f域序列
%%方位向(Azimuth,Cross-Range),x/u domain
v=100;%SAR 平台速度
D=4;
Lsar=lambda*Rc/D;%合成孔径长度
Na=1024;
x=linspace(-X0,X0,Na);%u域序列
u=x/v;
du=2*X0/v/Na;
fu=linspace(-1/2/du,1/2/du,Na);%fu域序列
fdc=2*v*sin(thetaT)/lambda;%Doppler调频中心频率
fdr=-2*(v*cos(thetaT))^2/lambda/Rc;%Doppler调频斜率
%%目标位置
Ptar=[Y,0,1];%距离向坐标,方位向坐标,sigma
%%产生回波
s_ut=zeros(Nr,Na);
U=ones(Nr,1)*u;%扩充为矩阵
T=t'*ones(1,Na);
rn=Ptar(1);xn=Ptar(2);sigma=Ptar(3);
R=sqrt(H^2+rn^2+(xn-v*U).^2);
DT=T-2*R/c;
phase=pi*Kr*DT.^2-4*pi/lambda*R;
s_ut=s_ut+sigma*exp(j*phase).*(abs(DT)<Tr/2).*(abs(v*U-xn)<Lsar/2);
%加海明窗
wr=hamming(Nr);
s_ut=s_ut.*(wr*ones(1,Na));%距离加窗
wa=hamming(Na);
s_ut=s_ut.*(ones(Nr,1)*wa');%方位加窗
%%距离压缩
p0_t=exp(j*pi*Kr*(t-2*Rc/c).^2).*(abs(t-2*Rc/c)<Tr/2);%距离向LFM信号
p0_f=fftshift(fft(fftshift(p0_t)));
s_uf=fftshift(fft(fftshift(s_ut)));%距离向FFT
src_uf=s_uf.*(conj(p0_f).'*ones(1,Na));%距离压缩
src_ut=fftshift(ifft(fftshift(src_uf)));%距离压缩后的信号
src_fut=fftshift(fft(fftshift(src_ut).')).';%距离多普勒域
%%距离徙动校正原理
F=f'*ones(1,Na);%扩充为矩阵
FU=ones(Nr,1)*fu;
Hrcc=exp(-j*pi/fc/fdr*FU.^2.*F-j*2*pi*Rc*sin(thetaT)/v*FU);%距离徙动校正函数
src_fuf=fftshift(fft(fftshift(src_uf).')).';%距离压缩后的二维频谱
s2rc_fuf=src_fuf.*Hrcc;
s2rc_fut=fftshift(ifft(fftshift(s2rc_fuf)));%距离多普勒域
%%方位压缩
p0_2fu=exp(j*pi/fdr*(FU-fdc).^2);%方位向压缩因子
s2rcac_fut=s2rc_fut.*p0_2fu;%方位压缩
s2rcac_fuf=fftshift(fft(fftshift(s2rcac_fut)));%距离方位压缩后的二维频谱
s2rcac_ut=fftshift(ifft(fftshift(s2rcac_fut).')).';%方位向IFFT
%
figure(1)
Sr=abs(src_ut(1:Nr,512));
Srmax=max(Sr);
Sr=20*log10(Sr/Srmax);
plot(1:Nr,Sr)
grid on
%最后成像
figure(2)
subplot(211)
G=(abs(src_ut));
gm=max(max(G));
gn=min(min(G));
G=255/(gm-gn)*(G-gn);
imagesc(x,sqrt((r*cos(thetaT)).^2-H^2),255-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('距离压缩后目标图象')
subplot(212)
G=(abs(s2rcac_ut));
gm=max(max(G));
gn=min(min(G));
G=255/(gm-gn)*(G-gn);
imagesc(x,sqrt((r*cos(thetaT)).^2-H^2),G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('目标图象')
figure(3)
Sr_ut=abs(src_ut).';
waterfall(Sr_ut((510:514),(1:Nr)));
axis tight
xlabel('Range')
ylabel('Azimuth')
title('距离压缩后目标图象')