基于压缩感知的SAR成像算法，有详细说明，结合压缩感知技术应用于SAR成像.zip

%分析SAR成像时的PSF %编写者：林月冠 %时间：2010-7-12 21:26:26 %%================================================================ clear;clc;close all; %%================================================================ %%雷达参数 lambda=1.245e-2; %载波波长 C=3e8; %光速 Fc=C/lambda; %载波频率 %%参数--目标区域 Xmin=-100; %方位向目标区域[Xmin,Xmax] Xmax=100; Yc=850e3; %图像中心 Y0=100; %地距目标区域[Yc-Y0,Yc+Y0] %图像宽度 2*Y0 %%参数--平台 Vr=7090; %平台速度100 m/s H=0; %平台高度 5000 m R0=sqrt(Yc^2+H^2); %图像中心斜距 %%参数--天线 D=10; %方位向天线长度 Lsar=0.886*lambda*R0/D; %SAR 合成孔径宽度 Tsar=2*Lsar/Vr; %SAR 合成孔径时间 %参数--慢时间域 Ka=2*Vr^2/lambda/R0; %多普勒调频率 Ba=abs(Ka*Tsar); %多普勒带宽 PRF=1.2*Ba/3; %脉冲重复频率 PRT=1/PRF; %脉冲重复时间 dat=PRT; %方位向采样时间间隔 Na=ceil((Xmax-Xmin+Lsar)/Vr/dat); %方位向采样数 sn=linspace((Xmin-Lsar/2)/Vr,(Xmax+Lsar/2)/Vr,Na);%方位向离散时间采样序列 df=-PRF/2:PRF/Na:(PRF/2-PRF/Na); %方位向频谱采样处 dat=PRT; %%%%参数--快时间域 Br=30e6; %距离向带宽 Fsr=1.1*Br; %距离采样率 Tau=1e-6; %脉冲持续时间 Kr=-Br/Tau; %距离向调频率 drt=1/Fsr; %距离向采样时间间隔 Rmin=sqrt((Yc-Y0)^2+H^2); %距离采样最近斜距 Rmax=sqrt((Yc+Y0)^2+H^2+(Lsar/2)^2); %距离采样的最远斜距 Nr=ceil(2*(Rmax-Rmin)/C/drt+Tau/drt); %距离向的采样数 tm=linspace(2*Rmin/C,2*Rmax/C+Tau,Nr); %距离向离散时间采样序列 fu=-Fsr/2:(Fsr/Nr):(Fsr/2 - Fsr/Nr); %距离向频率采样处 %%参数--分辨率 DY=C/2/Br; %距离向分辨率 DX=D/2; %方位向分辨率 Ntarget=25; %目标数 %%%%现实参数信息和目标属性 disp('参数显示 单位m');disp(' '); disp('距离向分辨率：');disp(DY) disp('方位向分辨率：');disp(DX) disp('SAR合成孔径长度：');disp(Lsar) %%================================================================ %%基信号 Ia=31;Ir=41; S0=zeros(Na,Nr); nn=ceil(numel(S0)*0.1); Basis=zeros(nn,Ia*Ir); for k_Ia=1:1:Ia for k_Ir=1:1:Ir Dslow=sn*Vr-(k_Ia-16)*PRT*Vr; %各个方位时刻平台相对于目标的方位距离序列 R=sqrt(Dslow.^2+(Yc-(k_Ir-21)*drt*C/2)^2+H^2); %各个方位时刻平台相对于目标的斜距 tau=2*R/C; %各个方位时刻平台到目标的来回时间序列 Dfast=ones(Na,1)*tm-tau'*ones(1,Nr)-Tau/2; %调制时间 phase=pi*Kr*Dfast.^2-(4*pi/lambda)*(R'*ones(1,Nr)); %见书本P156 S0=exp(j*phase).*(abs(Dfast)<Tau/2).*((abs(Dslow)<Lsar/2)'*ones(1,Nr)); inde=randperm(numel(S0)); ind=sort(inde(1:nn)); Basis(:,(k_Ia-1)*Ir+k_Ir)=S0(ind); end end %%=============================================================== figure,imagesc(abs(Basis)); %%% 计算PSF PSF=zeros(1,Ia*Ir); for k=1:Ia*Ir PSF(k)=abs(Basis(:,636)'*Basis(:,k))/norm(Basis(:,636))/norm(Basis(:,k)); end hh=figure; semilogy(PSF,'k');axis tight;grid on;hold on dd=PSF;dd((dd==max(dd)))=0;dd=max(dd)*ones(size(PSF)); semilogy(dd,'r','LineWidth',2); title('PSF');ylim([10^-2.5 1]) set(hh,'units','centimeters'); set(hh,'position',[5 5 11 7]); %%===========================CS重建=============================== Ns=Ia*Ir; %场景网格的点数 Nt=10; %场景中的目标点数 x=zeros(Ns,1); inde=randperm(Ns); ind=sort(inde(1:Nt)); x(ind)=1; y=Basis*x; %回波数据，Basis为观测矩阵 % 重建观测场景 tic Theta=[real(Basis) -imag(Basis);imag(Basis) real(Basis)]; Y=[real(y);imag(y)]; sigma2 = std(Y)^2/1e2; eta = 1e-20; lambda_init = []; adaptive = 0; optimal = 1; [weights,used]=FastLaplace(Theta,Y,sigma2,eta,lambda_init, adaptive,optimal); out=zeros(1,2*Ns);out(used)=weights; xr=out(1:Ns)+1i*out(Ns+1:end); scene=reshape(xr,31,41); %按31*41排列原来的1271*1的向量 toc hh=figure; imagesc(abs(scene));axis tight;axis off colormap(copper) set(hh,'units','centimeters'); set(hh,'position',[5 5 11 7]);