%%聚束SAR仿真
%%pfa算法，摘選自sap工具箱
%+-----------------------------------------------------------------------+
%|?This?is?the?main?program?that?controls?the?SAP?toolbox.???????????????|
%|?It?brings?up?the?menu?bar?and?sets?up?appropriate?pulldown?menus??????|
%|?corresponding?to?the?menu?bar.????????????????????????????????????????|
%|???????????????????????????????????????????????????????????????????????|
%|?Copyright?（c）?David?Kil?（1999）????????????????????????????????????????|
%|???????????????????????????????????????????????????????????????????????|
%|?Permission?to?use?the?toolbox?for?mine?countermeasure?and/or?other????|
%|?related?research?program?is?granted?as?long?as?the?user?acknowledges??|?
%|?the?use?of?the?toolbox?in?any?publication?and?no?commercial?use?of????|
%|?the?toolbox?is?planned?or?executed.?Please?acknowledge?the?Office?of??|
%|?Naval?Research?Code?321TS?and?Rockwell?Science?Center.?This?work????|
%|?sponsored?in?part?by?ONR?under?Grant?N00014-99-1-0320.????????????????|?
%+-----------------------------------------------------------------------+

clcclearclose?all;
Rs=1000;%天線(xiàn)相位中心到成像中心的距離
L=102;%天線(xiàn)沿方位向移動(dòng)的距離為
deltaX=0.4;%方位向分辨率
Tp=1e-6;%脈沖時(shí)寬
fc=3804.4e6;%中心頻率
B=133.5e6;%信號(hào)帶寬
Fs=185e6;%采樣頻率
Nfast=256;
C=3e8;%光速
K=B/Tp;%調(diào)頻斜率
Nslow?=?L/deltaX+1;
if?mod（Nslow2）
Nslow?=?Nslow+1;
end
H=200;%平臺(tái)高度
Li=8;%FFT升采樣率
Rc=sqrt（Rs^2+H^2）;
Tc=2*Rc/C;
t=linspace（Tc-Nfast/（2*Fs）Tc+Nfast/（2*Fs）?...
????Nfast）;
Yr=t*C/2-Rc;
t_m=ones（Nslow1）*t;
Ffast=（-Nfast/2:Nfast/2-1）*Fs/Nfast;
Ffast_m=ones（Nslow1）*Ffast;
Xa?=?-L/2+deltaX*（0:Nslow-1）;%天線(xiàn)位置
Xa_m=Xa‘*ones（1Nfast）;
Ntarget=2;
Ptarget=[0?Rs+0
????20?Rs+10];
signal=zeros（NslowNfast）;
for?i=1:Ntarget
????Rt=sqrt（（Xa_m-Ptarget（i1））.^2+H.^2+（Ptarget（i2））.^2）;
????tau=2*Rt/C;
????Dfast=t_m-tau;
????phase=pi*K*Dfast.^2-2*pi*fc*tau;
????signal=signal+exp（1i*phase）.*（abs（Dfast）end
figure;
imagesc（abs（signal））;
%%pfa算法仿真
%dechirp
Ra=sqrt（Rs.^2+H.^2+Xa_m.^2）;
tau=2*Ra/C;
Dfast=t_m-tau;
phase=pi*K*Dfast.^2-2*pi*fc*tau;
Sref=exp（1i*phase）;
S1=signal.*conj（Sref）;
%去除RVP（視頻相位）
Srvp=exp（-1i*pi*Ffast_m.^2/K）;
S2=ifty（fty（S1）.*Srvp）;
%極坐標(biāo)轉(zhuǎn)化為直角坐標(biāo)參數(shù)設(shè)定
phia=pi/2-asin（H./Ra）;
KR=4*pi*K/C*（fc/K+（t_m-2*Ra/C））;
Kp=KR.*（sin（phia））;
theta_p=asin（Xa_m./sqrt（Xa_m.^2+Rs^2））;
V=Kp.*cos（theta_p）;
U=Kp.*sin（theta_p）;
Uregion=linspace（U（11）U（Nslow1）Nslow）;
Vregion=linspace（V（11）V（1Nfast）Nfast）;
%距離插值
tli=linspace（Tc-Nfast/（2*Fs）Tc+Nfast/（2*Fs）?...
????Nfast*Li）;
tli_m=ones（Nslow1）*tli;
Temp=zeros（NslowNfast）;
min_value=zeros（NslowNfast）;
V1d=zeros（1Nfast）;
for?k?=?1:Nslow
????sig1d?=?（S2（k:））‘;
????temp?=?pfa_via_FFT（sig1dNfastLi）;
????KR1d?=?4*pi*K/C*（fc/K+tli-2*Ra（k1）/C）;
????Kp1d?=?KR1d*sin（phia（k））;
????V1d?=?Kp1d*cos（theta_p（k1））;
????[Temp（k:）?min_value（k:）]?=?findInterpolated（tempV1dVregion）;
end
%方位插值
msin_value?=?zeros（NslowNfast）;
for?k?=?1:Nfast
????sig1d?=?Temp（:k）;
????temp?=?pfa_via_