%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   %   PULSED SPOTLIGHT SAR SIMULATION AND RECONSTRUCTION   %   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%colormap(gray(256))cj=sqrt(-1);pi2=2*pi;%c=3e8;                   % propagation speedf0=50e6;                 % baseband bandwidth is 2*f0w0=pi2*f0;fc=200e6;                % carrier frequencywc=pi2*fc;lambda_min=c/(fc+f0);    % Wavelength at highest frequencylambda_max=c/(fc-f0);    % Wavelength at lowest frequencykc=(pi2*fc)/c;           % wavenumber at carrier frequencykmin=(pi2*(fc-f0))/c;    % wavenumber at lowest frequencykmax=(pi2*(fc+f0))/c;    % wavenumber at highest frequency%Xc=1000;                 % Range distance to center of target areaX0=20;                   % target area in range is within [Xc-X0,Xc+X0]Yc=300;                  % Cross-range distance to center of target areaY0=60;                  % target area in cross-range is within                         % [Yc-Y0,Yc+Y0]% Case 1: L < Y0; requires zero-padding of SAR signal in synthetic% aperture domain%  L=100;                 % synthetic aperture is 2*L% Case 2: L > Y0; slow-time Doppler subsampling of SAR signal spectrum% reduces computation%% L=400;                 % synthetic aperture is 2*Ltheta_c=atan(Yc/Xc);     % Squint angleRc=sqrt(Xc^2+Yc^2);      % Squint radial rangeL_min=max(Y0,L);         % Zero-padded aperture is 2*L_min%Xcc=Xc/(cos(theta_c)^2); % redefine Xc by Xcc for squint processing%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% u domain parameters and arrays for compressed SAR signal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%duc=(Xcc*lambda_min)/(4*Y0);      % sample spacing in aperture domain                                  % for compressed SAR signalduc=duc/1.2;                      % 10 percent guard band; this guard band                                  % would not be sufficient for targets                                  % outside digital spotlight filter (use                                  % a larger guard band, i.e., PRF)mc=2*ceil(L_min/duc);             % number of samples on apertureuc=duc*(-mc/2:mc/2-1);            % synthetic aperture arraydkuc=pi2/(mc*duc);                % sample spacing in ku domainkuc=dkuc*(-mc/2:mc/2-1);          % kuc array%dku=dkuc;                         % sample spacing in ku domain%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    u domain parameters and arrays for SAR signal     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%if Yc-Y0-L < 0,                            % minimum aspect angle theta_min=atan((Yc-Y0-L)/(Xc-X0));else, theta_min=atan((Yc-Y0-L)/(Xc+X0));end;theta_max=atan((Yc+Y0+L)/(Xc-X0));         % maximum aspect angle%du=pi/(kmax*(sin(theta_max)- ...                     sin(theta_min))); % sample spacing in aperture                                       % domain for SAR signaldu=du/1.4;                        % 20 percent guard bandm=2*ceil(pi/(du*dku));            % number of samples on aperturedu=pi2/(m*dku);                   % readjust duu=du*(-m/2:m/2-1);                % synthetic aperture arrayku=dku*(-m/2:m/2-1);              % ku array%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%       Fast-time domain parmeters and arrays          %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Tp=2.5e-7;                     % Chirp pulse durationalpha=w0/Tp;                   % Chirp ratewcm=wc-alpha*Tp;               % Modified chirp carrier%if Yc-Y0-L < 0, Rmin=Xc-X0;else, Rmin=sqrt((Xc-X0)^2+(Yc-Y0-L)^2);end;Ts=(2/c)*Rmin;                 % start time of samplingRmax=sqrt((Xc+X0)^2+(Yc+Y0+L)^2);Tf=(2/c)*Rmax+Tp;              % end time of samplingT=Tf-Ts;                       % fast-time interval of measurementTs=Ts-.1*T;                    % start slightly earlier (10% guard band)Tf=Tf+.1*T;                    % end slightly later (10% guard band)T=Tf-Ts;Tmin=max(T,(4*X0)/(c*cos(theta_max)));  % Minimum required T%dt=1/(4*f0);                 % Time domain sampling (guard band factor 2)n=2*ceil((.5*Tmin)/dt);      % number of time samplest=Ts+(0:n-1)*dt;             % time array for data acquisitiondw=pi2/(n*dt);               % Frequency domain samplingw=wc+dw*(-n/2:n/2-1);        % Frequency array (centered at carrier)k=w/c;                       % Wavenumber array%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Resolution for Broadside: (x,y) domain rotated by theta_c %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%DX=c/(4*f0);                      % range resolution (broadside)DY=(Xcc*lambda_max)/(4*L);         % cross-range resolution (broadside)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%           Parameters of Targets                 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%ntarget=9;                        % number of targets% Set ntarget=1 to see "clean" PSF of target at origin% Try this with other targets% xn: range;            yn= cross-range;    fn: reflectivity  xn=zeros(1,ntarget);  yn=xn;              fn=xn;% Targets within digital spotlight filter%  xn(1)=0;              yn(1)=0;            fn(1)=1;  xn(2)=.7*X0;          yn(2)=-.6*Y0;       fn(2)=1.4;  xn(3)=0;              yn(3)=-.85*Y0;      fn(3)=.8;  xn(4)=-.5*X0;         yn(4)=.75*Y0;       fn(4)=1.;  xn(5)=-.5*X0+DX;      yn(5)=.75*Y0+DY;    fn(5)=1.;% Targets outside digital spotlight filter% (Run the code with and without these targets)%    xn(6)=-1.2*X0;        yn(6)=.75*Y0;       fn(6)=1.;  xn(7)=.5*X0;          yn(7)=1.25*Y0;      fn(7)=1.;  xn(8)=1.1*X0;         yn(8)=-1.1*Y0;      fn(8)=1.;  xn(9)=-1.2*X0;        yn(9)=-1.75*Y0;     fn(9)=1.;  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                   SIMULATION                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%s=zeros(n,mc);     % SAR signal array%for i=1:ntarget;   % Loop for each target td=t(:)*ones(1,mc)-2*ones(n,1)*sqrt((Xc+xn(i)).^2+(Yc+yn(i)-uc).^2)/c; s=s+fn(i)*exp(cj*wcm*td+cj*alpha*(td.^2)).*(td >= 0 & td <= Tp & ...   ones(n,1)*abs(uc) <= L & t(:)*ones(1,mc) < Tf);end;%s=s.*exp(-cj*wc*t(:)*ones(1,mc));      % Fast-time baseband conversion% User may apply a slow-time domain window, e.g., power window, on% simulated SAR signal array "s" here.G=abs(s)';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);image(t,uc,256-cg*(G-ng));axis('square');axis('xy')xlabel('Fast-time t, sec')ylabel('Synthetic Aperture (Slow-time) U, meters')title('Measured Spotlight SAR Signal')print P5.1.pspause(1)%td0=t(:)-2*sqrt(Xc^2+Yc^2)/c;s0=exp(cj*wcm*td0+cj*alpha*(td0.^2)).*(td0 >= 0 & td0 <= Tp);s0=s0.*exp(-cj*wc*t(:));            % Baseband reference fast-time signals=ftx(s).*(conj(ftx(s0))*ones(1,mc));  % Fast-time matched filtering%G=abs(iftx(s))';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);tm=(2*Rc/c)+dt*(-n/2:n/2-1);    % fast-time array after matched filteringimage(tm,uc,256-cg*(G-ng));axis('square');axis('xy')xlabel('Fast-time t, sec')ylabel('Synthetic Aperture (Slow-time) U, meters')title('SAR Signal after Fast-time Matched Filtering')print P5.2.pspause(1)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  Slow-time baseband conversion for squint %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%kus=2*kc*sin(theta_c)*ones(1,n);     % Doppler frequency shift in ku                                     % domain due to squint%s=s.*exp(-cj*kus(:)*uc);             % slow-time baseband conversionfs=fty(s);% Display aliased SAR spectrum%G=abs(fs)';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);image(k*c/pi2,kuc,256-cg*(G-ng));axis('square');axis('xy')xlabel('Fast-time Frequency, Hertz')ylabel('Synthetic Aperture (Slow-time) Frequency Ku, rad/m')title('Aliased Spotlight SAR Signal Spectrum')print P5.3.pspause(1)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  Digital Spotlighting and Bandwidth Expansion in ku Domain  %%%%          via Slow-time Compression and Decompression        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%s=s.*exp(cj*kus(:)*uc);      % Original signal before baseband                             % conversion for squintcs=s.*exp(cj*2*(k(:)*ones(1,mc)).* ...       (ones(n,1)*sqrt(Xc^2+(Yc-uc).^2))-cj*2*k(:)*Rc*ones(1,mc));% compressionfcs=fty(cs);            % F.T. of compressed signal w.r.t. u%G=abs(fcs)';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);image(k*c/pi2,kuc,256-cg*(G-ng));axis('square');axis('xy')xlabel('Fast-time Frequency, Hertz')ylabel('Synthetic Aperture (Slow-time) Frequency Ku, rad/m')title('Compressed Spotlight SAR Signal Spectrum')print P5.4.pspause(1)%fp=iftx(fty(cs));      % Narrow-bandwidth Polar Format Processed                       % reconstruction%PH=asin(kuc/(2*kc));   % angular Doppler domainR=(c*tm)/2;            % range domain mapped from reference                       % fast-time domain%% Full Aperture Digital-Spotlight Filter%W_d=((abs(R(:)*cos(PH+theta_c)-Xc) < X0).* ...    (abs(R(:)*sin(PH+theta_c)-Yc) < Y0));%G=(abs(fp)/max(max(abs(fp)))+.1*W_d)';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);image((Rc/Xc)*(.5*c*tm-Rc),(kuc*Rc)/(2*kc),256-cg*(G-ng));xlabel('Range x, m')ylabel('Cross-range y, m')title('Polar Format SAR Reconstruction with Digital Spotlight Filter')axis image; axis xy;print P5.5.pspause(1)fd=fp.*W_d;                % Digital Spotlight Filteringfcs=ftx(fd);               % Transform to (omega,ku) domain% Zero-padding in ku domain for slow-time upsampling%mz=m-mc;        % number is zerosfcs=(m/mc)*[zeros(n,mz/2),fcs,zeros(n,mz/2)];%cs=ifty(fcs);              % Transform to (omega,u) domains=cs.*exp(-cj*2*(k(:)*ones(1,m)).* ...       (ones(n,1)*sqrt(Xc^2+(Yc-u).^2))+cj*2*k(:)*Rc*ones(1,m));% decompression%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                           CAUTION                             %% For TDC or backprojection, do not subsample in Doppler domain %% and do not perform slow-time baseband conversion               %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%s_ds=s;                    % Save s(omega,u) array for TDC and                           % backprojection algorithms%s=s.*exp(-cj*kus(:)*u);    % Slow-time baseband conversion for squintfs=fty(s);                 % Digitally-spotlighted SAR signal spectrum%G=abs(fs)';xg=max(max(G)); ng=min(min(G)); cg=255/(xg-ng);image(k*c/pi2,ku,256-cg*(G-ng));axis('square');axis('xy')xlabel('Fast-time Frequency, Hertz')ylabel('Synthetic Aperture (Slow-time) Frequency Ku, rad/m')title('Spotlight SAR Signal Spectrum after DS & Upsampling')print P5.6.pspause(1)