SAR成像学习(三)距离方向成像matlab代码解析 1
来源:互联网 发布:疯狂联盟座狼升级数据 编辑:程序博客网 时间:2024/05/17 20:33
本文将结合matlab代码讲解SAR距离向成像问题。
本文只研究距离向,且是正侧视情况。
文中以同一方位向坐标上四个目标点的成像为例,这四个目标的关系如下:
目标的相关信息:
% 关于目标%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Xc=2.e3; % Range distance to center of target areaX0=50; % target area in range is within [Xc-X0,Xc+X0]ntarget=4; % number of targets%%%%%%%%%%%%% Targets' parameters %%%%%%%%%%%%%%%%%%%% xn: range; fn: reflectivity 发射系数%xn(1)=0; fn(1)=1;xn(2)=.7*X0; fn(2)=.8;xn(3)=xn(2)+2*dx; fn(3)=1.;xn(4)=-.5*X0; fn(4)=.8;%注意这里的xn是相对于中间的Xc的位置%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
发射信号是线性调频信号:% 关于发射信号%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%colormap(gray(256))cj=sqrt(-1);pi2=2*pi;%c=3e8; % Propagation speedB0=100e6; % Baseband bandwidth is plus/minus B0 注意这里的带宽是2B0w0=pi2*B0;fc=1e9; % Carrier frequencywc=pi2*fc;Tp=.1e-6; % Chirp pulse durationalpha=w0/Tp; % Chirp ratewcm=wc-alpha*Tp; % Modified chirp carrier 即是式子中的beta%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1 问题的描述
成像的过程只要分为两步: 发射信号到接收信号;后处理,即接收信号到影像如图:
第一步是个正问题,主要由硬件完成,第二步是个逆问题,主要由软件完成。
第一个问题的输入
第二个问题的输入是
理想情况下
2 问题的分析和解决
其实问题的解决办法在上一个图中已经暗示了,他就是匹配滤波技术,这个技术我们已经在“急救箱系列”分析过了。但是面对实际问题,还是有许多问题需要说明,尤其是像我这样的新手。
2.1 再谈匹配滤波
匹配滤波的实施公式:
其中
接着推导匹配滤波的实施公式得到:
其中
这样就得到了我们的问题的答案:
这里的转换时用到了时间和距离的比例关系:
利用这个关系同样可以得到距离域的点扩散函数。
2.2 采样
这个问题是必然会碰见的问题,因为我们要进行的是“数字处理”。这里的关键问题有:
-
采样间隔,这个由信号带宽2B 0 决定:δt<=1/2B 0 =π/ω 0 -
采样起止时刻:T s =2(X c −X 0 )c T f =2(X c +X 0 )c + T p
其中,T p 表示信号持续的时间,最后加上这一项是为了避免出现回波信号中含有部分波的情况。 -
采样的时间段:t∈[T s ,T f ] T f −T s =4X 0 c + T p -
采样点数:n=2∗ceil((.5∗(Tf−Ts))/dt)
取这么多采样点是考虑到离散傅里叶变换的问题。 -
空间域和频率域采样序列:x i =ct i 2 =(X c −X 0 )+(i−1)cδt2 ω i =ω c +(i−n2 −1)δω δω=2πnδt
该部分的代码是:
% 采样问题%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%dt=pi/(2*alpha*Tp); % Time domain sampling (gurad band plus minus % 50 per) or use dt=1/(2*B0) for a general % radar signal%Ts=(2*(Xc-X0))/c; % Start time of samplingTf=(2*(Xc+X0))/c+Tp; % End time of sampling% Measurement parametersn=2*ceil((.5*(Tf-Ts))/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)x=Xc+.5*c*dt*(-n/2:n/2-1); % range bins (array); reference signal is % for target at x=Xc.%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2.3 结果
四个目标的回波信号对应的代码是:% 获得回波信号%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%s=zeros(1,n); % Initialize echoed signal array%暂时可以忽略na=8; % Number of harmonics in random phase ar=rand(1,na); % Amplitude of harmonicster=2*pi*rand(1,na); % Phase of harmonicsfor i=1:ntarget; td=t-2*(Xc+xn(i))/c; %注意这里的时间是对采样时间序列的一个延迟,另外有一个平移 %这个要跟后面的参考信号也就是匹配滤波器对应上,而且他会影响恢复出来的f(x)的解译 pha=wcm*td+alpha*(td.^2); % Chirp (LFM) phase for j=1:na; % Loop for CPM harmonics 相位调制!可以暂时忽略 pha=pha+ar(j)*cos((w(n/2+1+j)-wc)*td+ter(j)); end; s=s+fn(i)*exp(cj*pha).*(td >= 0 & td <= Tp); %回波信号 注意是各个目标的回波信号的叠加 %特别注意这里的Tp,是脉冲持续的时间end;
匹配滤波用到的参考信号对应的代码是:% Reference echoed signal%td0=t-2*(Xc+0)/c; %注意这里的时间参考点%这个要跟前面的回波信号的时间对应上,而且他会影响恢复出来的f(x)的解译pha0=wcm*td0+alpha*(td0.^2); % Chirp (LFM) phasefor j=1:na; % Loop for CPM harmonics pha0=pha0+ar(j)*cos((w(n/2+1+j)-wc)*td0+ter(j));end;s0=exp(cj*pha0).*(td0 >= 0 & td0 <= Tp);
解调及显示:% Baseband conversion 解调 % 这里和该系列前面讲的正交解调还是不一样的,当时讨论过% 为什么会出来了复数图像,现在看来,这个问题有点幼稚。% 线性调频信号本身就可以表示成复数,回波信号也是。% 任何一个数只要我们愿意,都可以在形式上表示成复数(无非实数的虚数部分为零嘛)。% 不应该觉得复数是什么奇怪的事情,而应该觉得它实实在在,就在那里,这是见得少了罢了。%sb=s.*exp(-cj*wc*t);sb0=s0.*exp(-cj*wc*t);plot(t,real(sb))%只是显示了振幅xlabel('Time, sec')ylabel('Real Part')title('Baseband Echoed Signal')axis('square')axis([Ts Tf 1.1*min(real(sb)) 1.1*max(real(sb))])print P1.1a.pspause(1)%plot(t,real(sb0))xlabel('Time, sec')ylabel('Real Part')title('Baseband Reference Echoed Signal')axis('square')axis([Ts Tf 1.1*min(real(sb0)) 1.1*max(real(sb0))])print P1.2a.pspause(1)%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
![这里写图片描述]()
注意图中只是画出了振幅信息,而忽略了相位信息,实际上相位信息是非常有用的,特别是在有的应用中。
匹配滤波及结果:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 目标函数重建%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Fourier Transform%fsb=fty(sb);fsb0=fty(sb0);% Power equalization% 这里是对参考信号做了一个处理,暂时可以先忽略mag=abs(fsb0);amp_max=1/sqrt(2); % Maximum amplitude for equalizationafsb0=abs(fsb0);P_max=max(afsb0);I=find(afsb0 > amp_max*P_max);nI=length(I);fsb0(I)=((amp_max*(P_max^2)*ones(1,nI))./afsb0(I)).* ... exp(cj*angle(fsb0(I)));%% Apply a window (e.g., power window) on fsb0 here%E=sum(mag.*abs(fsb0));%% Matched Filtering%fsmb=fsb.*conj(fsb0);%%% Inverse Fourier Transform%smb=ifty(fsmb); % Matched filtered signal (range reconstruction)%% Display%plot(x,(n/E)*abs(smb))xlabel('Range, meters')ylabel('Magnitude')title('Range Reconstruction Via Matched Filtering')axis([Xc-X0 Xc+X0 0 1.1]); axis('square')print P1.6a.pspause(1)
结果:
3 如果不采用匹配滤波技术
如果不采用匹配滤波技术,可以用一种叫pulse compression的方法。但是它对信号的形式有一定的限制,比如线性调频信号。上面介绍的匹配滤波方法并不局限于线性调频信号。
代码:
% 另一种方法Time domain Compression%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%解调并压缩td0=t-2*(Xc+0)/c;scb=conj(s).* ... exp(cj*wcm*td0+cj*alpha*(td0.^2)); % Baseband compressed signal%plot(t,real(scb))xlabel('Time, sec')ylabel('Real Part')title('Time Domain Compressed Signal')axis('square')print TimeDomainCompressedSignal.pngpause(1)%因为上面的信号包含了与目标的位置相关的频率成分,所以只要进行一个傅里叶变换即可得到问题的解。fscb=fty(scb);X=(c*(w-wc))/(4*alpha); % Range array for time domain compression 注意参考点plot(X+Xc,(dt/Tp)*abs(fscb)) %?xlabel('Range, meters')ylabel('Magnitude')title('Range Reconstruction Via Time Domain Compression')axis([Xc-X0 Xc+X0 0 1.1]); axis('square')%axis([Xc-X0 Xc+X0 0 1.1*max(abs(fscb))]); axis('square')pause(1)
结果:
说明:
观察上面的结果和匹配滤波得到的结果可以发现明显的不同,这里的原因我暂时也不是特别清楚,也可能是这段代码有问题。但是必须注意到:这肯定是由于信号的压缩造成的。这个后面再回来说。
以上所有的分析中有一个问题被我掩盖了,那就是幅宽问题,如果幅宽很大,以上的分析没有问题,但是如果幅宽相对脉冲持续时间比较小,会带来一个问题,这个问题的解决需要在匹配滤波的过程中添加一些步骤。
未完待续。
- SAR成像学习(三)距离方向成像matlab代码解析 1
- SAR成像学习(四)距离方向成像matlab代码解析 2
- SAR成像学习(五)方位向成像及matlab代码解析
- SAR成像学习(二) “Hello World”版本matlab代码
- 【matlab】雷达成像系列 之 RM(Range Migration,距离迁徙)成像算法
- 成像
- SAR成像学习(一)信号到原始数据&原始数据到图像
- SAR成像基础知识急救箱(二)关于离散傅里叶变换
- 摄像机成像、畸变模型(三)
- 点光源小孔成像matlab
- 摄像头成像 (待续)
- 基于二维傅里叶变换法的MRI成像原理的Matlab仿真(1)
- 小孔成像的解析-20151205
- SAR成像基础知识急救箱(一)卷积 相关 滤波器那些事儿
- SAR成像基础知识急救箱(零)关于傅里叶变换的几个小困惑
- 超声弹性成像代码公开
- PIV(粒子成像测速)
- 成像圈(Image Circle)
- C++第4次实验
- IOS学习 网络HTTP Get和Post请求与登录界面加密 涉及按钮边框圆角与颜色设置
- 电话监听器
- lightoj 1198 - Karate Competition 【贪心】
- IOS App常用界面结构解析,让开发更简单
- SAR成像学习(三)距离方向成像matlab代码解析 1
- 四 单例模式
- leetcode刷题8. String to Integer (atoi)
- ThreadLocal类源码解析
- java---序列六(合并流)SequenceInputStream ——对多个流进行合并
- 致自己
- 开发工具
- js中call与apply用法
- 深入理解JVM(九)——类加载的过程