%%%产生雷达碰基回波%%%%
echo_mobj_pulse=[zeros(1,100),1,1,zeros(1,N-102)]
s_pc=repmat(echo_mobj_pulse,1,16)
n=1:16*N
s_doppler=cos(n*2*pi*f_doppler/笑升谨fs)
s_pc=s_pc.*s_doppler
s_noise=0.1*rand(1,N*16)
s_pc=s_pc+s_noise
figure,
CS%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%产生Stripmap SAR的回波
clear all
clc
thetaT=0%T平台波束斜视角
thetaT=thetaT*pi/180%rad
thetaR=0%R平台波束斜视角
thetaR=thetaR*pi/180
c=3e8%光速
fc=1.5e9%载频
lambda=c/fc%波长
%%测绘带区域
X0=200%方位向[-X0,X0]
Rtc=3000
Rrc=3000
Rc=(Rtc+Rrc)/2
R0=150%距离向[Rc-R0,Rc+R0]
%%距离向(Range),r/t domain
Tr=1.33e-6%LFM信号脉宽1.33us (200m)
Br=150e6%LFM信号带宽 150MHz
Kr=Br/Tr%调频斜率
Nr=1024
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 平台速度
Lsar=300%合成孔径长度
Na=512
x=linspace(-X0,X0,Na)%u域序列
u=x/v
du=2*X0/v/Na
fu=linspace(-1/2/du,1/2/du,Na)%fu域序列
ftdc=v*sin(thetaT)
ftdr=-(v*cos(thetaT))^2/lambda/Rtc
frdc=v*sin(thetaR)
frdr=-(v*cos(thetaR))^2/lambda/Rrc
fdc=ftdc+frdc%Doppler调频中心频率
fdr=ftdr+frdr%Doppler调频斜率
%%目标位置
Ntar=3%目标个数
Ptar=[Rrc,0,1 %距离向坐标,方位向坐标,sigma
Rrc+50,-50,1
Rrc+50,50,1]
%%产生回波
s_ut=zeros(Nr,Na)
U=ones(Nr,1)*u%扩充为矩阵
T=t'*ones(1,Na)
for i=1:1:Ntar
rn=Ptar(i,1)xn=Ptar(i,2)sigma=Ptar(i,3)
rtn=rn+Rtc-Rrc
RT=sqrt(rtn^2+(rtn*tan(thetaT)+xn-v*U).^2)
RR=sqrt(rn^2+(rn*tan(thetaT)+xn-v*U).^2)
R=RT+RR
DT=T-R/c
phase=-pi*Kr*DT.^2-2*pi/lambda*R
s_ut=s_ut+sigma*exp(j*phase).*(abs(DT)<Tr/2).*(abs(v*U-xn)<Lsar/2)
end
%方位向fft
s_kt=fftshift(fft(fftshift(s_ut).')).'
%CS变换
kc=4*pi/lambda
kc=kc*ones(1,Na)
kx=fu/v
p_kx0=-sqrt(kc.^2-kx.^2)%相位项泰勒展开的系数函数
p_kx1=2*kc/c/p_kx0
p_kx2=-2.*kx.^2/早基c^2./p_kx0.^3
C_kx=-(c*p_kx1/2+1)
Ks_r=1-2*Kr*Rc.*p_kx2
Ks_kx_r=Kr/pi./Ks_r
r0=Rc
s2_ut=exp(j*pi*C_kx.*ones(Nr,1)*Ks_kx_r.*(t'*ones(1,Na)-2*r0*(1+C_kx)/c).^2)%设计的线性调频信号
S_cs=s_kt.*s2_ut
%距离向fft
S_kw=fftshift(fft(fftshift(S_cs)))
%距离向匹配滤波
w=2*pi*f
rmc_r=exp(j.*w*2*C_kx*r0/c).*exp(j.*w.^2/4/贺樱pi/Kr/(1+C_kx))
rmc_r=rmc_r'陆拍谨*ones(1,Na)
S_rmc=S_kw.*rmc_r
%距离向ifft
S_kt=fftshift(ifft(fftshift(S_rmc)))
d_kxr=4*pi/c^2*Kr*C_kx*(1+C_kx).*(Rc-r0).^2%CS变换带来的相位误差
S_kt=S_kt.*exp(-j*d_kxr)%消除相位误差
%方位向匹配滤波
FU=ones(Nr,1)*fu
H_kx=exp(j*pi/fdr*(FU-fdc).^2)%方位向压缩因子
I_ut=S_kt.*H_kx
I_ut=fftshift(ifft(fftshift(I_ut.'))).'
subplot(221)
G=20*log10(abs(s_ut)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(a)原始信号')
subplot(222)
G=20*log10(abs(S_rmc)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(b)距离向匹配滤波后频谱')
subplot(223)
G=20*log10(abs(S_kt)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(c)消除相位误差后频谱')
subplot(224)
G=20*log10(abs(I_ut)+1e-6)
gm=max(max(G))
gn=gm-60%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(d)目标图象')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
RD
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%产生Stripmap SAR的回波
clear all
thetaT=0%T平台波束斜视角
thetaT=thetaT*pi/180%rad
thetaR=0%R平台波束斜视角
thetaR=thetaR*pi/180
c=3e8%光速
fc=1.5e9%载频
lambda=c/fc%波长
%%测绘带区域
X0=200%方位向[-X0,X0]
Rtc=3000
Rrc=3000
Rc=(Rtc+Rrc)/2
R0=150%距离向[Rc-R0,Rc+R0]
%%距离向(Range),r/t domain
Tr=1.5e-6%LFM信号脉宽 1.5us (200m)
Br=150e6%LFM信号带宽 150MHz
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 平台速度
Lsar=300%合成孔径长度
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域序列
ftdc=v*sin(thetaT)
ftdr=-(v*cos(thetaT))^2/lambda/Rtc
frdc=v*sin(thetaR)
frdr=-(v*cos(thetaR))^2/lambda/Rrc
fdc=ftdc+frdc%Doppler调频中心频率
fdr=ftdr+frdr%Doppler调频斜率
%%目标位置
Ntar=3%目标个数
Ptar=[Rrc,0,1 %距离向坐标,方位向坐标,sigma
Rrc+50,-50,1
Rrc+50,50,1]
%%产生回波
s_ut=zeros(Nr,Na)
U=ones(Nr,1)*u%扩充为矩阵
T=t'*ones(1,Na)
for i=1:1:Ntar
rn=Ptar(i,1)xn=Ptar(i,2)sigma=Ptar(i,3)
rtn=rn+Rtc-Rrc
RT=sqrt(rtn^2+(rtn*tan(thetaT)+xn-v*U).^2)
RR=sqrt(rn^2+(rn*tan(thetaT)+xn-v*U).^2)
R=RT+RR
DT=T-R/c
phase=pi*Kr*DT.^2-2*pi/lambda*R
s_ut=s_ut+sigma*exp(j*phase).*(abs(DT)<Tr/2).*(abs(v*U-xn)<Lsar/2)
end
%%距离压缩
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).')).'%距离多普勒域
%%二次距离压缩,距离迁移校正原理仿真
src_fuf=fftshift(fft(fftshift(src_uf).')).'%距离压缩后的二维频谱
F=f'*ones(1,Na)%扩充为矩阵
FU=ones(Nr,1)*fu
p0_2f=exp(j*pi/fc^2/fdr*(FU.*F).^2+j*pi*fdc^2/fc/fdr*F-j*pi/fc/fdr*FU.^2.*F)
s2rc_fuf=src_fuf.*p0_2f
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
subplot(221)
G=20*log10(abs(s_ut)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(a)原始信号')
subplot(222)
G=20*log10(abs(src_fut)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(fu,r-Rc,-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(b)距离多普勒域频谱')
subplot(223)
G=20*log10(abs(s2rc_fut)+1e-6)
gm=max(max(G))
gn=gm-40%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(fu,r-Rc,-G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(c)RMC后的RD域频谱')
subplot(224)
G=20*log10(abs(s2rcac_ut)+1e-6)
gm=max(max(G))
gn=gm-60%显示动态范围40dB
G=255/(gm-gn)*(G-gn).*(G>gn)
imagesc(x,r-Rc,G),colormap(gray)
grid on,axis tight,
xlabel('Azimuth')
ylabel('Range')
title('(d)目标图象')
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)