Untitled7.zip_wk_wk 合成孔径_wk成像_合成孔径成像

%回波信号仿真 %parameter c=3e8; %光速 fc=5.3e9; lambda=c/fc; %SAR基本参数 %SAR平台 v=100; %雷达平台水平速度 H=5000;%高度5000m %测绘带数据 xmin=-100; xmax=100; yc=1500; y0=200; thita=0; R0=sqrt(yc^2+H^2);%目标与SAR垂直斜距 D=3; %天线等效孔径 La=3; Lsar=lambda*R0/D/cos(thita)^2;%SAR长度 Tsar=Lsar/v; %方位向参数 Ka=2*v^2*cos(thita)^3/lambda/R0;%调频率 BW=Ka*Tsar;%带宽 PRF=2.0*BW;%脉冲频率 PRT=1/PRF;%脉冲时间 ds=PRT; Nslow=ceil((xmax-xmin+Lsar)/v/ds); Nslow=2^nextpow2(Nslow); sn=linspace((xmin-Lsar/2)/v,(xmax+Lsar/2)/v,Nslow); PRT=(xmax-xmin+Lsar)/v/Nslow; PRF=1/PRT; ds=PRT; %距离向参数 Tr=2e-6; Br=100e6; %脉冲信号带宽 Kr=Br/Tr; Fsr=1.2*Br; %距离向采样频率 dt=1/Fsr; Rmin=sqrt((yc-y0)^2+H^2); Rmax=sqrt((yc+y0)^2+H^2); Nfast=ceil(2*(Rmax-Rmin)/c/dt+Tr/dt); tm=linspace(2*Rmin/c-Tr/2,2*Rmax/c+Tr/2,Nfast); dt=(2*Rmax/c+Tr-2*Rmin/c)/Nfast; f0=0; Fsr=1/dt; %分辨率 DY=c/2/Br;%距离向 DX=D/2;%方位向 %目标点 Ntarget=1; Ptarget=[(xmin+xmax)/2, yc, 1 %目标位置 xmin, yc-y0, 1 xmin yc+y0 1 xmax yc-y0 1 xmax yc+y0 1 1/3*xmin yc-y0 1 1/3*xmax yc-y0 1 2/3*xmin yc-y0 1 2/3*xmax yc+y0 1 1/3*xmin yc-1/3*y0 1 1/3*xmin yc-2/3*y0 1 1/3*xmin yc+1/3*y0 1 1/3*xmin yc+2/3*y0 1 1/3*xmax yc-1/3*y0 1 1/3*xmax yc-2/3*y0 1 1/3*xmax yc+1/3*y0 1 1/3*xmax yc+2/3*y0 1 2/3*xmin yc-1/3*y0 1 2/3*xmin yc-2/3*y0 1 2/3*xmin yc+1/3*y0 1 2/3*xmin yc+2/3*y0 1 2/3*xmax yc-1/3*y0 1 2/3*xmax yc-2/3*y0 1 2/3*xmax yc+1/3*y0 1 2/3*xmax yc+2/3*y0 1]; M=Nfast; N=Nslow; K=Ntarget; T=Ptarget; Srnm=zeros(N,M); for k=1:1:K sigma=T(k,3); Dslow=sn*v-T(k,1); R=sqrt(Dslow.^2+T(k,2)^2+H^2);%斜距 tau=2*R/c; Dfast=ones(N,1)*tm-tau'*ones(1,M); phase=pi*Br/Tr*Dfast.^2-(4*pi/lambda)*(R'*ones(1,M));%相位 ssssssss=-Tr/2<Dfast&Dfast<Tr/2; Srnm=Srnm+sigma*exp(1i*phase).*ssssssss.*((abs(Dslow)<Lsar/2)'*ones(1,M)); end tr=tm-2*Rmin/c; Refr=exp(j*pi*Kr*tr.^2).*(-Tr/2<tr&tr<Tr/2); Sr=ifft((fft(Srnm')'.*(ones(N,1)*conj(fft(Refr')')))'); Gr=abs(Sr); ta=sn-xmin/v; Refa=exp(j*pi*Ka*ta.^2).*(abs(ta)<Tsar/2); Sa=ifft(fft(Sr).*((conj(fft(Refa)).'*ones(1,M)))'); Ga=abs(Sa); %w_fast=zeros(1,512); %for ii=1:512; % ss=tm(1,ii); %w_fast(1,ii)=rectpuls(ss/Tr); %end %W_fast=fft(w_fast); %pa_thita=sinc(La*thita/lambda); %yitac=0; %w_slow=(sinc(atan(v*(sn-yitac)/R0))).^2; %W_slow=fft(w_slow); uuu=(atan(v*sn/R0)*La/lambda);%时域方位向天线方向图中间变量 pa=sinc(uuu).^2; Srnm=Srnm.*(pa'*ones(1,Nfast)); figure; title('点目标仿真回波信号'); imagesc(abs(Srnm)); colormap(gray); Srnm1=fft2(Srnm); imagesc(abs(Srnm1)); colormap(gray); Srnm1=fftshift(Srnm1); Ka=2*v^2*(cos(thita))^3/lambda/R0; fyita=linspace(-PRF/2,PRF/2,Nslow); %获得慢时间频率变量;在(-Fslow/2,Fslow/2)范围; ftao= linspace(-Fsr/2,Fsr/2,Nfast); %获得快时间频率变量;在(-Ffast/2,Ffast/2)范围; Theta_range_comp = ones(Nslow,1)*(pi/Kr*ftao.^2); %用于距离压缩的相位 Theta_sqrt_1 = ones(Nslow,1) * ((fc+ftao).^2); %根号内第一项 Theta_sqrt_2 = (c^2*fyita.^2/4/v^2).'*ones(1,Nfast); %根号内第二项 Theta_sqrt = 4*pi*R0/c*sqrt(Theta_sqrt_1-Theta_sqrt_2); %根号项 Theta_2d_ref = Theta_range_comp + Theta_sqrt; %获得总相位 Ref_2d= exp(1i*Theta_2d_ref); %获得参考函数 Srnm2= Ref_2d.*Srnm1; %参考函数相乘 imagesc(abs(Srnm2)); Delta_ffast = sqrt(ones(Nslow,1)*(fc + ftao).^2 + (c^2*fyita.^2/4/v^2).'*ones(1,Nfast)) - ones(Nslow,1) * (ftao+fc); Delta_ffast = Delta_ffast * Nfast/Fsr;%归一化 Nd = max(max(fix(Delta_ffast))); temp = zeros(Nslow,Nfast); %临时数据矩阵temp; Ni = 16/2; %sinc插值核长度Ni=8; for i = 1 : 1 : Nslow data = Srnm2(i,:); for m = Ni+Nd : 1 : Nfast-Nd-Ni Insert_AMT = Delta_ffast(i,m); AMT = Insert_AMT; %获得插值量(位移量); AMT_N = fix(AMT); %取插值位移量的整数部分 AMT_F = AMT - AMT_N; %取插值位移量的小数部分 NN = AMT_N + (-Ni+1:Ni); %所取的插值位置 NN1 = AMT_F - (-Ni+1:Ni); %所取的sinc核位置 %% 进行sinc核数据插值计算 sinc_kernal = sinc(NN1); %sinc核函数的值; temp(i,m) = sinc_kernal * data(1,m+NN).';%重建方程 end end Srnm3 = temp; Wfast=kaiser(Nfast,2.5); Wfast=ones(Nslow,1)*(Wfast.'); Srnm3=Srnm3.*Wfast; imagesc(abs(Srnm3)); Srnm4 = ifft2(Srnm3); colormap(gray); imagesc(255-abs(Srnm4)); figure; title('Image after 2-d ifft'); xlabel('距离向');ylabel('方位向'); [a,b]=max(Srnm4); [x,y]=max(a); z=b(y); Nslice=32; times=16; slice_2D=zeros(Nslice,Nslice); slice_2D(1:Nslice,1:Nslice)=Srnm4(z-Nslice/2:1:z+Nslice/2-1,y-Nslice/2:1:y+Nslice/2-1); slice_2D_fft=fft2(slice_2D); imagesc(abs(slice_2D_fft)); slice_2D_interp=zeros(Nslice*times,Nslice*times); slice_2D_interp(1:Nslice,1:Nslice)=slice_2D_fft; slice_2D_ifft=ifft2(fftshift(slice_2D_interp)); slice=abs(slice_2D_ifft); figure(3) [a1,b1]=max(slice); [x1,y1]=max(a1); z1=b1(y1); plot(20*log10(slice(:,y1)/max(slice(:,y1)))); xlabel('方位向'); ylabel('幅度（db）'); title('方位向剖面图'); figure(4) plot(20*log10(slice(z1,:)/max(slice(z1,:)))); xlabel('距离向'); ylabel('幅值（db）'); title('距离向剖面图');