SAR信号处理.rar_SAR_SAR信号处理_factorygbk

%%================================================================ %%Filename: stripmapSAR.m %%Help file: stripmapSAR.doc %%Project: Stripmap SAR Simulation using point targets and Reconstrction %% with Range-Doppler Algorithm %%Author: Zhihua He ,National University of Defence Tecnology ,2005/3 %%E-mail: skynismile@yahoo.com.cn %%================================================================ clear;clc;close all; %%================================================================ %%Parameter--constant 恒定参数 C=3e8; %propagation speed %%Parameter--radar characteristics Fc= 1.5000e+010; %carrier frequency 1GHz 载波 lambda=C/Fc; %wavelength 波长 %%Parameter--target area Xmin=0; %target area in azimuth is within[Xmin,Xmax] Xmax=50; Yc=9000; %center of imaged area Y0=200; %target area in range is within[Yc-Y0,Yc+Y0] r=0.02;wr=32*pi;beta0=30*pi/180;fai0=0;alpha0=30*pi/180; %imaged width 2*Y0 %%Parameter--orbital information 轨道信息 V=200; %SAR velosity 100 m/s H=4500; %height 5000 m R0 =sqrt(Yc^2+H^2); %%Parameter--antenna 天线参数 D=0.5; %antenna length in azimuth direction 雷达天线孔径长度 Lsar=lambda*R0/D; %SAR integration length 合成孔径长度 Tsar=Lsar/V; %SAR integration time 合成孔径时间 %%Parameter--slow-time domain 慢时间域参数设置 Ka=-2*V^2/lambda/R0; %doppler frequency modulation rate 慢时间调频率 Ba=abs(Ka*Tsar); %doppler frequency modulation bandwidth 多普勒宽度 PRF=2048; %pulse repitition frequency % 发射脉冲重复频率 PRT=1/PRF; %pulse repitition time 脉冲重复间隔 ds=PRT; %sample spacing in slow-time domain 慢时间域的采样间隔,就等于脉冲重复周期 Nslow=ceil((Xmax-Xmin+Lsar)/V/ds); %sample number in slow-time domain 采样数 %Nslow=2^nextpow2(Nslow); %for fft 为快速傅里叶变换做准备 tm=linspace((Xmin-Lsar/2)/V,(Xmax+Lsar/2)/V,Nslow);%discrete time array in slow-time domain 慢时间离散化 PRT=(Xmax-Xmin+Lsar)/V/Nslow; %refresh 修正以后的PRT 更新 PRF=1/PRT; ds=PRT; %%Parameter--fast-time domain 快时间域参数设置 Tr=10e-7; %pulse duration 10us 发射的脉冲宽度 Br=150e6; %chirp frequency modulation bandwidth 30MHz 线性调频信号的带宽 Kr=Br/Tr; %chirp slope 快时间调频斜率 Fsr=2*Br; %sampling frequency in fast-time domain 快时间域的采样频率 dt=1/Fsr; %sample spacing in fast-time domain 采样间隔 Rmin=sqrt((Yc-Y0)^2+H^2); %最小作用距离 Rmax=sqrt((Yc+Y0)^2+H^2+(Lsar/2)^2); %最大作用距离为什么要加(Lsar/2)^2？？？？？ Nfast=ceil(2*(Rmax-Rmin)/C/dt+Tr/dt); %sample number in fast-time domain %采样点数 向上取整 % Rwid = Rmax-Rmin; %接收距离窗 下面三句和上面一句功能一样 % Twid = 2*Rwid/C; %时间窗 % Nfast = ceil(Twid/dt); %为什么加Tr/dt 因为最小距离从零开始的 %Nfast=2^nextpow2(Nfast); %for fft tfast=linspace(2*Rmin/C,2*Rmax/C+Tr,Nfast); %discrete time array in fast-time domain 离散化 dt=(2*Rmax/C+Tr-2*Rmin/C)/Nfast; %refresh Fsr=1/dt; %%Parameter--resolution 距离维参数 DY=C/2/Br; %range resolution 距离分辨率 DX=0.25; %cross-range resolution 横向距离分辨力 也就是方位分辨力 %%Parameter--point targets 点目标 Ntarget=1; %number of targets 点目标数目 %format [x, y, reflectivity] Ptarget=[Xmin,Yc,1; %position of targets 点目标的位置和RCS Xmin,Yc+10*DY,1; Xmin+20*DX,Yc+50*DY,1] disp('Parameters:') disp('Sampling Rate in fast-time domain');disp(Fsr/Br) disp('Sampling Number in fast-time domain');disp(Nfast) disp('Sampling Rate in slow-time domain');disp(PRF/Ba) disp('Sampling Number in slow-time domain');disp(Nslow) disp('Range Resolution');disp(DY) disp('Cross-range Resolution');disp(DX) disp('SAR integration length');disp(Lsar) disp('Position of targets');disp(Ptarget) %%================================================================ %%Generate the raw signal data 产生需要处理的数据矩阵 K=Ntarget; %number of targets N=Nslow; %number of vector in slow-time domain M=Nfast; %number of vector in fast-time domain T=Ptarget; %position of targets Srnm=zeros(N,M); %初始化回波数据矩阵 for k=1:1:K sigma=T(k,3); %取出目标的反射系数 Dslow=tm*V-T(k,1); %V*tm-x R=sqrt(Dslow.^2+T(k,2)^2+H^2)+r*cos(beta0)*cos(wr*tm+fai0+alpha0);%斜距计算公式 tau=2*R/C; %计算得到每个目标的延迟 Dfast=ones(N,1)*tfast-tau'*ones(1,M); %计算得到快时间变量 phase=pi*Kr*Dfast.^2-(4*pi/lambda)*(R'*ones(1,M)); Srnm=Srnm+sigma*exp(j*phase).*(0<Dfast&Dfast<Tr).*((-Tsar/2<(tm)&(tm)<Tsar/2)'*ones(1,M)); end %% %%================================================================ %%Range compression 距离维匹配滤波 tr=tfast-2*Rmin/C; % Refr=exp(j*pi*Kr*tr.^2).*(0<tr&tr<Tr); Sr=ifty(fty(Srnm).*(ones(N,1)*conj(fty(Refr))));%复共轭 Gr=abs(Sr); %插值法完成距离徙动校正 %方位向FFT 变为距离多普域进行距离弯曲校正 srSa_RD = fftshift(fft(fftshift(Sr))); % %第二种方法进行截断sinc插值进行距离徒动校正 % h = waitbar(0,'Sinc插值'); % P=4; %4点sinc核插值 % N=Nslow; M=Nfast; % RMCmaxtix = zeros(N,M); % for n=1:N % for m=P:M % delta_R = (1/8)*(lambda/V)^2*(R0+(m-M/2)*C/2/Fsr)*((n-N/2)*PRF/N)^2;%首先计算距离迁移量 计算方法就是把斜距变换到距离多普勒域就知道了 % RMC=2*delta_R*Fsr/C; %距离徒动了几个距离单元 % delta_RMC = RMC-round(RMC); %距离徒动量的小数部分 % for i = -P/2:P/2-1 % if m+RMC+i>M %判断是否超出边界 % RMCmaxtix(n,m)=RMCmaxtix(n,m)+srSa_RD(n,M)*sinc((-i+RMC)); % else % RMCmaxtix(n,m)=RMCmaxtix(n,m)+srSa_RD(n,m+round(RMC)+i)*sinc((-i+delta_RMC)); % end % end % end % waitbar(n/N) % end % close(h) % %======================== % Sr_rmc=fftshift(ifft(fftshift(RMCmaxtix))); %%距离徙动校正后还原到时域 % Gar = abs(Sr_rmc); %测试距离压缩效果 colormap(jet); figure('Name','距离压缩'); imagesc((Gr)); axis on; title('SAR图像'); xlabel('距离向采样点点'); ylabel('方位向采样点'); %Azimuth compression 方位向匹配滤波 ta=[1:Nslow]*PRT; Sa_ref=rectpuls(ta-Tsar/2,Tsar).*exp(j*pi*Ka*(ta-Tsar/2).^2); %Sr=ifty(fty(Srnm).*(ones(N,1)*conj(fty(Refr))));%复共轭 Sa=iftx(ftx(Sr).*(conj(ftx(Sa_ref)).'*ones(1,M))); Ga=abs(Sa); %测试方位压缩效果 colormap(jet); figure('Name','方位压缩'); imagesc((Ga)); axis on; title('SAR图像'); xlabel('距离向采样点'); ylabel('方位向采样点'); % %等高线图 % colormap(jet) % figure('Name','SAR图像'); % contour(Ga,30); % axis on; % title('SAR图像'); % xlabel('距离向采样点'); % ylabel('方位向采样点'); % %%================================================================ % %%graw the intensity image of signal % colormap(gray); % figure(1) % subplot(211); % row=tfast*C/2-2008;col=tm*V-26; %