rd.m

rd算法仿真

clear all;
%%%(I) parameters' definition
c=3e+8;										% speed of light
pi=3.1415926	 					   	% pi
j00=sqrt(-1);								% square root of -1

res_a=1;										% required azimuth resolution
res_r=1;										% required range resolution
k_a=1.2;										% azimuth factor
k_r=1.2;										% range factor

Ra=4000.;									% radar working distance
va=70.; 										% radar/platform forward velocity
Tp=1.e-6; 									% transmitted pulse width      
fc=1.e+9;	 						   	% carrier frequency 
FsFactor = 1.0

lamda=c/fc;									% wavelength
B=k_r*c/2./res_r;							% required transmitted bandwidth
Fs=B*FsFactor;								% A/D sampling rate
bin_r=c/2./Fs;								% range bin
gama=B/Tp;									% range chirp rate	

delta_theta=k_a*lamda/2/res_a;		% required azimuth opening angle
												% corresponding to processed bandwith
												% (should be less than azimuth 3dB opening angle)
La=Ra*delta_theta;						% required synthetic aperture
Ta=La/va;									% required synthetic aperture time
thetas_c=pi/2;								% zero squint
fdc=2*va*cos(thetas_c)/lamda;  		% doppler centriod
fdr=-2*va*va*sin(thetas_c)^2/lamda/Ra;					% doppler rate
Bd=abs(fdr)*Ta;							% doppler bandwidth
prf=100.;									% PRF	

%%%(II) echo return modelling(point target)
na=fix(Ta*prf/2);							% azimuth sampling number
ta=-na:na;											
ta=ta/prf;									% slow time along azimuth
xa=va*ta;									% azimuth location along flight track
na=2*fix(na);

%x0=[0  0 0 0 0 0   -60 -30 30 60]; 
%R0=[-10 0 10 20 30 40 0 0 0 0];
x0=[0 ]; 
R0=[0 ];

Npt_num = length(x0);
												% single point target's azimuth location 
												% positive towards forward velocity
												% single point target's slant range location
ra=zeros(Npt_num, length(xa));      % positive towards far range
for i=1:Npt_num                                    
	ra(i,:)=sqrt((Ra+R0(i)).^2+(xa+x0(i)).^2);		% single point target's slant range history 
end   
rmax=max(max(ra));						% max. slant range
rmin=min(min(ra));						% min. slant range
rmc=fix((rmax-rmin)/bin_r);			% range migration	
rg=0*ra;
rg=fix((ra-rmin)/bin_r+1);				% range gate beause of range migration
rgmax=max(max(rg));
rgmin=min(min(rg));

nr=round(Tp*Fs);							% min range samples
tr=1:fix(nr)+1;									
tr=tr/Fs-Tp/2;								% fast time along range
                                               
[m,Na]=size(ta);
Nr=nr+rgmax;
Na
Nr

sig=zeros(Na,Nr);
for i=1:Na			
	for k=1:Npt_num
  		sig(i,rg(k,i):rg(k,i)+nr)=sig(i,rg(k,i):rg(k,i)+nr)+exp(-j00*4*pi/lamda*ra(k,i))*exp(-j00*pi*gama*(tr).^2);
   end  
end
disp('end of echo return modelling')
%figure
%contour(abs(sig))


%%%===============================================================================

for i=1:Na
   sig(i,:)=fftshift(fft(sig(i,:)));		
end
disp('end of range fft')

dfr=Fs/Nr;
fr=((0:Nr-1)-Nr/2)*dfr;  				
phase1=exp(-j00*pi*(fr).^2/gama);
for i=1:Na
   sig(i,:)=ifft(fftshift( sig(i,:).*phase1) );		
end   
disp('end of range compression')

%figure
%contour(abs(sig))

for i=1:Nr
  sig(:,i)=fftshift(fft(sig(:,i)));
end
disp('end of azimuth fft')

fdr=-2*va*va/lamda/Ra;			
dfa=prf/Na;
fa=((0:Na-1)-Na/2)*dfa;  			

figure
%contour(abs(sig(Na,:)))
plot(abs(sig(Na,:)))

R=4;
Up00=zeros(size(1:Nr*R+2000));
kkk=1:Nr;
for i=1:Na
  Up00=0*Up00;
  Up00(1:Nr*R)=interp(sig(i,:),R);
  dRa=Ra/sqrt(1-(lamda*fa(i)/2/va)^2)-Ra;
  kpos(i)=round(R*dRa/bin_r);
  sig(i,:)=Up00((kkk-1)*R+kpos(i)+1);
end

%figure
%contour(abs(sig))
figure
%contour(abs(sig(Na,:)))
plot(abs(Up00))


disp('Range motion correction')

rr=rmin+((1:Nr)-1)*bin_r;
delt_beta = sqrt(1-(lamda/2/va*fa).^2)-1;
phase2 = exp((j00*4*pi/lamda)*(delt_beta'*rr));
sig = sig.*phase2;
for i=1:Nr
   sig(:,i)=ifft(sig(:,i));
end
figure
contour(abs(sig));