西安邮电大学光学光波偏振态的仿真报告 本文关键词：偏振，光波，西安，仿真，光学
西安邮电大学光学光波偏振态的仿真报告 本文简介：电子工程学院光学课程设计实验报告姓名：系部：光电子技术系专业：年级：学号：指导教师：地点：2号实验楼234时间：2015/12/21--2015/12/31光波偏振态的仿真一、实验目的通过对两相互垂直偏振态的合成1．掌握圆偏振、椭圆偏振及线偏振的概念及基本特性；2．掌握偏振态的分析方法。任务与要求：
西安邮电大学光学光波偏振态的仿真报告 本文内容：
电子工程学院
光学课程设计
实验报告
姓
名：
系
部：
光电子技术系
专
业：
年
级：
学
号：
指导教师：
地
点：
2号实验楼234
时
间：
2015/12/21--2015/12/31
光波偏振态的仿真
一、实验目的
通过对两相互垂直偏振态的合成
1．掌握圆偏振、椭圆偏振及线偏振的概念及基本特性；
2．掌握偏振态的分析方法。
任务与要求：
对两相互垂直偏振态的合成进行计算，绘出电场的轨迹。要求计算在j=0、j=p/4、j=p/2、j=3p/4、j=p、j=5p/4、j=3p/2、j=7p/4时，在Ex=Ey及Ex=2Ey情况下的偏振态曲线并总结规律
2、
实验原理
平面光波是横电磁波，其光场矢量的振动方向与光波传播方向垂直。一般情况下，在垂直平面光波传播方向的平面内，光场振动方向相对光传播方向是不对称的，光波性质随光场振动方向的不同而发生变化。将这种光振动方向相对光传播方向不对称的性质，称为光波的偏振特性。它是横波区别于纵波的最明显标志。
1)
光波的偏振态
根据空间任一点光电场E的矢量末端在不同时刻的轨迹不同，其偏振态可分为线偏振、圆偏振和椭圆偏振。
设光波沿z方向传播，电场矢量为
为表征该光波的偏振特性，可将其表示为沿x、y方向振动的两个独立分量的线性组合，即
其中
将上二式中的变量t消去，经过运算可得
式中，φ=φy-φx。这个二元二次方程在一般情况下表示的几何图形是椭圆，如图1-1所示。
图1-1
椭圆偏振诸参量
在上式中，相位差φ和振幅比Ey/Ex的不同，决定了椭圆形状和空间取向的不同，从而也就决定了光的不同偏振状态。图1-2画出了几种不同φ值相应的椭圆偏振态。实际上，线偏振态和圆偏振态都可以被认为是椭圆偏振态的特殊情况。
图1-2
不同j值相应的椭圆偏振
(1)
线偏振光
当Ex、Ey二分量的相位差φ=mπ(m=0,±1,±2，…)时，椭圆退化为一条直线，称为线偏振光。此时有
当m为零或偶数时，光振动方向在Ⅰ、Ⅲ象限内；当m为奇数时，光振动方向在Ⅱ、Ⅳ象限内。
由于在同一时刻，线偏振光传播方向上各点的光矢量都在同一平面内，因此又叫做平面偏振光。通常将包含光矢量和传播方向的平面称为振动面。
(2)
圆偏振光
当Ex、Ey的振幅相等(E0x=E0y=E0)，相位差φ=mπ/2(m=±1,±3,±5…)时，椭圆方程退化为圆方程
该光称为圆偏振光。用复数形式表示时，有
式中，正负号分别对应右旋和左旋圆偏振光。所谓右旋或左旋与观察的方向有关，通常规定逆着光传播的方向看，E为顺时针方向旋转时，称为右旋圆偏振光，反之，称为左旋圆偏振光。
(3)
椭圆偏振光
在一般情况下，光场矢量在垂直传播方向的平面内大小和方向都改变，它的末端轨迹是椭圆，故称为椭圆偏振光。在某一时刻，传播方向上各点对应的光矢量末端分布在具有椭圆截面的螺线上(图1-3)。椭圆的长、短半轴和取向与二分量Ex、Ey的振幅和相位差有关。其旋向取决于相位差φ：当2mπ＜φ＜(2m+1)π时，为右旋椭圆偏振光；当(2m-1)π＜φ＜2mπ时，为左旋椭圆偏振光。
图1-3
椭圆偏振光
3、
程序流程图
开始
显示图像，结束程序
画出二维图像：
subplot(4,4,n);
n=n+1;
plot(Ex,Ey);
循环计算;
Fy=0:pi/4:7*pi/4
Ex=Eox*cos(w*t-k*z);
Ey=Eoy*cos(w*t-k*z+Fy);
定义c、lamd、w、k、Eox=5、
Eoy=10、t、z、i=1、n=9
画出三维图像：
subplot(4,4,i);
i=i+1;
plot3(Ex,Ey,z);
循环计算：
Fy=0:pi/4:7*pi/4;
Ex=Eox*cos(w*t-k*z);
Ey=Eoy*cos(w*t-k*z+Fy);
4、
结果分析
由理论知识可以知道光的偏振态的合成与振幅和相位差有关，即相位差ψ和振幅比Ey/Ex的不同，决定了椭圆形状和空间取向的不同，从而决定了光的不同偏振状态。
如上图取得是2Ex=Ey
的不同相位时的偏振合成，当二者的相位差ψ=mπ（m=0，±1，±2，…）时合成为线偏振光，即第一幅图和第五幅图为线偏振光的图像，可以看出合成图为一条线。而椭圆的长、短半轴和取向与二分量Ex、Ey的振幅和相位差有关，其旋向取决于相位差ψ：当2mπ<ψ<(2m+1)π时，为右旋椭圆偏振光；(2m-1)π<ψ<2mπ时，为左旋椭圆偏振光。第二种方法：迎着光的传播方向看，若光矢量沿顺时针方向转动，称为右旋椭圆偏振光，反之称为左旋的，这个方法也可以判断圆偏振光的旋向。如果把振幅改为Ex=Ey进行仿真会发现只要相位差ψ=mπ/2(m=±1,±3,±5,…)时，偏振合成为圆偏振光。此时ψ值仿真结果会出现线偏振，圆偏振和椭圆偏振的合成图像。
思考题
1.说明偏振的定义；
答：光场的振动方向相对光的传播方向的不对称性叫光的偏振。
为什么圆偏振2.椭圆和线偏振是完全偏振光？
答：应为它们
3.如何确定光的左右旋？
答：规定逆着光传播方向看，E为顺时针方向旋转时，称为右旋圆偏振光，反之，称为左旋圆偏振光。
2.
如何区分圆偏振和自然光？
答：通过1/4波片，再通过偏振片，然后旋转偏振片，若光强不变化，为自然光；若光强有变化，出现两次消光，为圆偏光。
3.
如何区分椭圆偏振和部分偏振光？
答：通过1/4波片，并且最大或最小方向与波片光轴方向一致或垂直，再通过偏振片，并旋转偏振片
有消光现象为椭圆偏振，无消光的为部分偏振光。
6.根据仿真结果总结左右旋的规律。
答：?=mπ时候为线偏光，m=0/偶数时，在一、三象限；m=奇数时，在二、四象限；?=mπ/2时，为圆偏振光；其它为椭圆偏振光。
五、仿真小结
这是仿真的第一个题目，而且我也不是第一次接触matlab，因此也很快的仿真出结果。但这是我头一次使用matlab来仿真物理现象，这让我对matlab有了新的认识。在仿真过程中还学了不少实用的语法以及指令，总之仿真实习不仅巩固了我光学的基础，还帮助我提高了matlab的编程能力真是一举两得。
附录：
clear
all;
c=3e+8;
%光速
lamd=5e-7;
%波长
T=lamd/c;
t=linspace(0,T,1000);
z=linspace(0,5,1000);
w=2*pi/T;
k=2*pi/lamd;
%波数
Eox=5;
Eoy=10;
i=1;
for
Fy=0:pi/4:7*pi/4;
Ex=Eox*cos(w*t-k*z);
Ey=Eoy*cos(w*t-k*z+Fy);
subplot(4,4,i);
i=i+1;
plot3(Ex,Ey,z);
axis
equal;
axis
normal;
zlabel(
z
);
xlabel(
x
);
ylabel(
y
);
end
n=9;
for
Fy=0:pi/4:7*pi/4
Ex=Eox*cos(w*t-k*z);
Ey=Eoy*cos(w*t-k*z+Fy);
subplot(4,4,n);
n=n+1;
plot(Ex,Ey);
ylabel(
y
);
xlabel(
x
);
axis
equal;
end
光波场的时域频谱
1、
实验目的
1．掌握单色光、复色光的概念；
2．掌握准单色光的概念及光波频谱宽窄的影响因素。
任务与要求：
对常见光波
①无限长等幅振荡
②持续有限时间的等幅振荡，持续时间为1ns、1ms、1s、10s、100s
③指数衰减振荡E(t)=e-bte-i2pn0，(t≥0)，b=0、1、5、10、100
进行傅里叶变换计算并绘出频谱图，总结影响频谱宽窄的因素。等时间进行计算，
2、
实验原理
实际上，严格的单色光波是不存在的，我们所能得到的各种光波均为复色波。所谓复色波是指某光波由若干单色光波组合而成，或者说它包含有多种频率成分，它在时间上是有限的波列。
复色波的电场是所含各个单色光波电场的叠加，
即
在一般情况下，若只考虑光波场在时间域内的变化，可以表示为时间的函数E(t)。通过傅里叶变换，它可以展成如下形式：
即一个随时间变化的光波场振动E(t)，可以视为许多单频成分简谐振荡的叠加，
各成分相应的振幅E(ν)，并且E(ν)按下式计算：
|E(ν)|2表征了ν频率分量的功率，称|E(ν)|2为光波场的功率谱。
一个时域光波场E(t)可以在频率域内通过它的频谱描述。下面，给出几种经常运用的光波场E(t)的频谱分布。
(1)
无限长时间的等幅振荡
它的频谱为
表明，等幅振荡光场对应的频谱只含有一个频率成分ν0，称其为理想单色振动。
图3-1
等幅振荡及其频谱图
(2)
持续有限时间的等幅振荡
图3-2
有限正弦波及其频谱图
(3)
衰减振荡
图3-3
衰减振荡及其频谱图
显示函数
plot(t，E1)
plot(f,abs(F1))
衰减震荡E3=Eo*exp(-B*t3)*exp(-2i*pi*f.*t3)*(heaviside(t3))
F3=fourier(E3)
持续有限时间等幅震荡
E2=Eo*exp(-2i*pi*f*t2)*(heaviside(t2+tao)-heaviside(t2-tao))
F2=fourier(E2)
无限长时间的等幅震荡
E1=Eo*exp(-2i*
pi*f*t1)
F1=fft(E1)
定义变量:t,w,f,Eo,B,tao
开始
三、程序流程图
显示函数
ezplot(t,E3,[-1,10])
ezplot(w,abs(F3),[-10,-2])
显示函数
ezplot(t，E2,[-2,2])
ezplot(w,abs(F2),[-22,10])
结束
4、
结果分析
从上面的仿真结果可以看出，当光波为无限长等振幅时它的频域为一冲击函数，表明该光波为单色波只包含一种频率。而有限长等振幅光波场的频域包含多种频率。最后的衰减振荡的频域有一个中心频率v0并且具有一定谱宽，随着衰减因子β的减小其频谱宽度越来越小，逐渐趋于单色波。
实际上第二种光波与单色波的不同是，单色波是无限延伸的，而第二种波只是单色波的一段，通常称为波列。根据公式:
表明波列长度2L和波列所包含的单色分波的波长范围成反比关系，波列越短，波列所包含的单色波的波长范围就越宽；相反，波列越长，波列所包含的单色分波的波长范围就越窄。当波列长度等于无穷大时，等于零，这就是单色波。
思考题
1.
如何获得准单色光？
答：对于一个实际的表观频率为υ0的振荡，若其振幅随时间的变化比振荡本身缓慢很多，则这种振荡的平率就集中于υ0附近的一个很窄的频段内，可认为是中心频率为υ0的准单色光。
4.
影响光的单色性的因素有哪些？
答：β和频率，振幅。
5.
衰减震荡中β的含义？
答：衰减因子。
五、仿真小结
本次实验虽然看上去很简单，但是在编写完后无论如何也调试不出来，检查了好几遍没没发现究竟什么地方有错误，感觉是matlab里面的傅里叶变换和阶跃函数之间存在bug，最后用了fft函数才解决这个问题。本次实验不仅锻炼我们的书本知识，也磨练了我们分析问题，解决问题的能力，合作的能力，而且这次试验也告诉我结束们往往在你想放弃的时候，也许就在成功路上的90%，再坚持一下就能成功了。
附录：
clear
all;
syms
t
;
Eo=1;
f=1;
T=2;
b=0.5;
E=Eo*exp(-2i*pi*f*t);
subplot(3,2,1);
ezplot(t,E,[-10,10]);
F=fourier(E);
subplot(3,2,2);
ezplot(abs(F),[-10,10]);
E2=Eo*exp(-2i*pi*f*t)*(heaviside(t+T)-heaviside(t-T));
F2=fourier(E2);
subplot(3,2,3);
ezplot(t,E2,[-6,6]);
subplot(3,2,4);
ezplot(abs(F2),[-15,5]);
axis
equal;
E3=Eo*exp(-b*t)*exp(-2i*pi*f*t)*(heaviside(t));
F3=fourier(E3);
subplot(3,2,5);
ezplot(t,E3,[-2,10]);
subplot(3,2,6);
ezplot(abs(F3),[-15,0]);
axis
equal;
双光束干涉
一、实验目的
1．掌握光的相干条件；
2．掌握分波阵面双光束干涉的特点。
任务与要求：
对双缝干涉进行计算，分别绘出单色光和复色光（白光）的干涉条纹，总结双缝干涉的特点。
二、实验原理
1.
两束光的干涉现象
光的干涉是指两束或多束光在空间相遇时，
在重叠区内形成稳定的强弱强度分布的现象。
例如，图5-1所示的两列单色线偏振光
图5-1
两列光波在空间重叠
在空间P点相遇，E1与E2振动方向间的夹角为θ，则在P点处的总光强为
式中，I1、I2是二光束的光强；φ是二光束的相位差，且有
由此可见，二光束叠加后的总强度并不等于这两列波的强度和，而是多了一项交叉项I12，它反映了这两束光的干涉效应，通常称为干涉项。干涉现象就是指这两束光在重叠区内形成的稳定的光强分布。所谓稳定是指，用肉眼或记录仪器能观察到或记录到条纹分布，即在一定时间内存在着相对稳定的条纹分布。显然，如果干涉项I12远小于两光束光强中较小的一个，就不易观察到干涉现象；如果两束光的相位差随时间变化，使光强度条纹图样产生移动，且当条纹移动的速度快到肉眼或记录仪器分辨不出条纹图样时，就观察不到干涉现象了。
在能观察到稳定的光强分布的情况下，满足
m=0,±1,±2，…
的空间位置为光强极大值处，且光强极大值IM为
满足
φ=(2m+1)π
m=0,±1,±2，
的空间位置为光强极小值处，且光强极小值Im为
当两束光强相等，即I1=I2=I0时，相应的极大值和极小值分别为
IM=2I0(1+cosθ)
Im=2I0(1-cos
θ)
2.
产生干涉的条件
首先引入一个表征干涉效应程度的参量——干涉条纹可见度，由此深入分析产生干涉的条件。
1)
干涉条纹可见度(对比度)
干涉条纹可见度定义为
当干涉光强的极小值Im=0时，V=1，二光束完全相干，条纹最清晰；当IM=Im时，V=0，二光束完全不相干，无干涉条纹；当IM≠Im≠0时，0＜V＜1，二光束部分相干，条纹清晰度介于上面两种情况之间。
2)
产生干涉的条件
由上述二光束叠加的光强分布关系可见，影响光强条纹稳定分布的主要因素是：二光束频率；二光束振动方向夹角和二光束的相位差。
(1)
对干涉光束的频率要求
由二干涉光束相位差的关系式可以看出，当二光束频率相等，Δω=0时，干涉光强不随时间变化，可以得到稳定的干涉条纹分布。当二光束的频率不相等，Δω≠0时，干涉条纹将随着时间产生移动，且Δω愈大，条纹移动速度愈快，当Δω大到一定程度时，肉眼或探测仪器就将观察不到稳定的条纹分布。因此，为了产生干涉现象，要求二干涉光束的频率尽量相等。
(2)
对二干涉光束振动方向的要求
当二光束光强相等时
V=cosθ
因此，当θ=0、二光束的振动方向相同时，V=1，干涉条纹最清晰；当θ=π/2、二光束正交振动时，V=0，不发生干涉；当0＜θ＜π/2时，0＜V＜1，干涉条纹清晰度介于上面两种情况之间。所以，为了产生明显的干涉现象，要求二光束的振动方向相同。
(3)
对二干涉光束相位差的要求
由式可见，为了获得稳定的干涉图形，二干涉光束的相位差必须固定不变，即要求二等频单色光波的初相位差恒定。实际上，考虑到光源的发光特点，这是最关键的要求。
可见，要获得稳定的干涉条纹，则：
①
两束光波的频率应当相同；
②
两束光波在相遇处的振动方向应当相同；
③
两束光波在相遇处应有固定不变的相位差。
这三个条件就是两束光波发生干涉的必要条件，通常称为相干条件。
三、程序流程图
开始程序
定义波长l，狭缝的间隔d，狭缝到屏的距离D，设置光屏范围等参数。
定义等间距的矢量矩阵，即仿真光屏y方向分成n个点
第二次调用for函数，对各采样点进行计算，
实现对复色光的计算。
调用for函数，据光强式，对各采样点进行计算。
调用imagesc(x,y,I)绘制图像,调用plot(I,y)绘制光强分布曲线。
结束程序
四、结果分析
从仿真结果可以得知，单色光的相干性非常好，在无限远处仍可以看见明暗相间的干涉条纹，而复色波随着光程差的增大其条纹对比度逐渐下降，最后降为零，完全看不清条纹。这是由于复色波有一定的光谱宽度，这实际上是限制了所产生清晰条纹的光程差。
对于光谱宽度为的光源，能产生干涉条纹的最大光程差称为相干长度。假定在某一光程差下，波长为的第m级条纹和波长为的第m+1级条纹重合，即这两种波长条纹的相对移动量达到一个条纹，那么波长为的第m级和第m-1级条纹之间便充满范围内其他波长的条纹，因而该处各点强度相等，条纹对比度降为零，无法看到条纹。故可以求得相干长度为
表明能够发生干涉的最大光程差与光源的光谱宽度成反比。另外相干长度实际上等于波列长度。这说明利用波列长度和光谱宽度的概念来讨论问题完全等效。
光波在一定光程差下能够发生干涉的事实表明了光波的时间相干性。我们把光通过相干长度所需的时间称为相干时间。
思考题
1光的相干条件
？
答：在相遇的地方，频率相同，振动方向相同，相位相同或有恒定的相位差。
2试讨论光源分波面法和分振幅法的相干性并说明如何用非相干光源获得相干？
答：分波面法是将一个波列的波面分成几部分，由这每一部分发出的波再相遇时，必然是相干的；分振幅法是利用透明薄板的第一，第二表面对入射光的依次反射，将入射光的振幅分解成若干部分，将这些不同部分的光波相遇时将产生干涉；要获得相干光，要把一个波列的光分成两束或几束光波，然后令其重合而产生稳定的干涉效应，这样的方法可以使相干光束初相位差保持恒定。
3为什么双光束干涉是分波阵面法
？
答：一束光透过两个缝，分成两束光在观察屏上叠加，有恒定的相位差。
4解释干涉的时间相干概念并用复色光的仿真进行解释
？
答：实际光源都包含有一定的光谱宽度，在干涉试验中，范围内的每一种波长的光都生成各自的一组干涉条纹，因此，光源的光谱宽度限制了干涉条纹的可见度。复色光在宽度内各光谱分量产生的总光强为
对于一定的，可见度V随着增大而下降；当=0时，光源为单色光，V=1；
当0＜＜时，0＜V＜1；当=时，V=0。
5假如利用光的干涉现象进行长度的测量，试分析光源用宽谱还是窄谱的精度高
？
答：用窄谱近似于单色光，单色性更好；用宽谱时，干涉的光强分布集中，精度更高。
五、仿真小结
通过本次光学仿真，使我对书本的知识有了更深的理解。本来在光学实验室已经做了关于干涉的实验，如果说那个是宏观的话，那么这次仿真就是很好的微观教学，本来书本上的东西时间久了容易混淆，这次实验那些仿真图十分生动形象，给我留下了很深的印象，作为仿真的第三个实验，刚开始接触觉得还是很有难度，但随着理解和小伙伴们一起研究，最终我们还是出色完成了这个实验，给人很大的成就感。
附录：
clear
all;
lamd=5e-7;
d=0.005;
D=1;
x=1;
k=1e-3;
%干涉场长度
y=linspace(-k,k,1000)
for
n=1:1000;
r1=sqrt((y(n)-d/2)^2+D^2);
r2=sqrt((y(n)+d/2)^2+D^2);
phase=2*pi*(r2-r1)/lamd;
I(n,:)=4*cos(phase/2)^2;
end
colormap(gray);
subplot(1,4,1)
imagesc(x,y,I);
subplot(1,4,2);
plot(I(:),y)
for
n=1:1000;
s=0;
r1=sqrt((y(n)-d/2)^2+D^2);
r2=sqrt((y(n)+d/2)^2+D^2);
dl=linspace(0,0.2,5);
for
N=1:5;
%各个频点在干涉场上的光强叠加
lamd1=lamd*(1+dl(N));
phase2=2*pi*(r2-r1)/lamd1;
s=s+4*cos(phase2/2)^2;
end
I(n,:)=s;
end
subplot(1,4,3)
imagesc(x,y,I);
subplot(1,4,4);
plot(I(:),y)
光的圆孔衍射
一、实验目的
1．掌握近场和远场的概念；
2．掌握夫琅禾费圆孔衍射特点及艾里斑的概念；
3．掌握菲涅尔圆孔衍射的特点。
任务与要求：
利用教材3.1-15式对圆孔衍射进行计算，其中入射波长为632.8nm，圆孔半径为1mm，光源位于系统的轴线上。改变光源位置及观察屏位置，观察远场衍射图案及艾里斑，近场观察距离改变衍射图案的变化；对仿真结果进行总结分析。
二、实验原理
光的衍射是指光波在传播过程中遇到障碍物时，所发生的偏离直线传播的现象。光的衍射，也可以叫光的绕射，即光可绕过障碍物，传播到障碍物的几何阴影区域中，并在障碍物后的观察屏上呈现出光强的不均匀分布。通常将观察屏上的不均匀光强分布称为衍射图样。
如图9-1所示，让一个足够亮的点光源S发出的光透过一个圆孔Σ，照射到屏幕K上，并且逐渐改变圆孔的大小，就会发现：当圆孔足够大时，在屏幕上看到一个均匀光斑，光斑的大小就是圆孔的几何投影(图3-1(a))；随着圆孔逐渐减小，起初光斑也相应地变小，而后光斑开始模糊，并且在圆斑外面产生若干围绕圆斑的同心圆环(图3-1(b))，当使用单色光源时，这是一组明暗相间的同心环带，当使用白色光源时，这是一组色彩相间的彩色环带；此后再使圆孔变小，光斑及圆环不但不跟着变小，反而会增大起来。这就是光的衍射现象。
图9-1
光的衍射现象
由于光学仪器的光瞳通常是圆形的，因而讨论圆孔衍射现象对光学仪器的应用，具有重要的实际意义。
夫朗和费圆孔衍射的讨论方法与矩形孔衍射的讨论方法相同，只是由于圆孔结构的几何对称性，采用极坐标处理更加方便。如图9-2所示，设圆孔半径为a,圆孔中心O1位于光轴上，则圆孔上任一点Q的位置坐标为ρ1、φ1，与相应的直角坐标x1、y1的关系为
图9-2
夫朗和费圆孔衍射光路
x1=ρ1cosφ1
y1=ρ1
sinφ1
类似地，观察屏上任一点P的位置坐标ρ、φ与相应的直角坐标的关系为
由此，P点的光场复振幅在经过坐标变换后为
式中
是衍射方向与光轴的夹角，称为衍射角。在这里，已利用了sinθ≈θ的近似关系。
根据零阶贝塞尔函数的积分表示式
可将P点的光场复振幅变换为
可得
式中，J1(x)为一阶贝塞尔函数。因此，P点的光强度为
由上式，可以得到夫朗和费圆孔衍射的如下特点：
(1)
衍射图样
由于Φ=kaθ，夫朗和费圆孔衍射的光强度分布仅与衍射角θ有关(或者，由于θ=ρ/f，仅与ρ有关)，而与方位角φ坐标无关。这说明，夫朗和费圆孔衍射图样是圆形条纹(图9-3)。
图9-3
圆孔夫朗和费衍射图样图
9-3
夫朗和费圆孔衍射光强度分布
(2)
衍射图样的极值特性
由贝塞尔函数的级数定义，可将P点的光强度表示为
该强度分布曲线如图9-3所示。当Φ=0时，即对应光轴上的P0点，有I=I0，它是衍射光强的主极大值。当Φ满足J1(Φ)=0
时，I=0，这些Φ值决定了衍射暗环的位置。在相邻两个暗环之间存在一个衍射次极大值，其位置由满足下式的Φ值决定：
这些次极大值位置即为衍射亮环的位置。上式中，J2(Φ)为二阶贝塞尔函数。表9-1列出了中央的几个亮环和暗环的Φ值及相对光强大小。
表9
–
1
圆孔衍射的光强分布
(3)
爱里斑
由表9-1
可见，中央亮斑集中了入射在圆孔上能量的83.78%，这个亮斑叫爱里斑。爱里斑的半径ρ0由第一光强极小值处的Φ值决定，即
因此
或以角半径θ0表示
爱里斑的面积为
式中，S为圆孔面积。可见，圆孔面积愈小，爱里斑面积愈大，衍射现象愈明显。
只有在S=0.61πfλ
时，S0=S。
三、程序流程图
开始
结束
plot(abs(Br));
画出接收屏上对应光谱
imshow(abs(U));
画衍射图样
h=exp(j*k*z)*exp((j*k*(x.^2+y.^2))/(2*z))/(j*lamda*z);
H
=fftshift(fft2(h));
B=fftshift(fft2(I));
G=H.*B;
U=
fftshift(ifft2(G));
菲涅尔衍射的复振幅分布是孔径平面的复振幅分布和一个二次相位因子乘积的傅里叶变换
Br=(U/max(U));
归一化
i=find(D<=r);
I(i)=1;
画出圆孔，衍射屏图像
定义lamda，k，z，I，圆孔半径r，以及屏幕取样点数N
四、结果分析
光的衍射是指光波在传播过程中遇到障碍物时，所发生的偏离直线传播的现象。光的衍射，也叫光的绕射，即光可绕过障碍物，传播到障碍物的几何阴影区域中，并在障碍物后的观察屏上呈现出光强的不均匀分布。上图是菲涅尔圆孔衍射图样，根据圆孔的大小我们可以看到衍射中心有明暗变换。本次程序的原理是由傅里叶光学内容推出。菲涅尔衍射的傅里叶变换表达式表明除了积分号前面的一个与x1、y1无关的振幅和相位因子外，菲涅尔衍射的复振幅分布是孔径平面的复振幅分布和一个二次相位因子乘积的傅里叶变换。
根据采用的距离近似的不同，衍射区还有另一种划分方法：衍射效应可以忽略的几何投影区，衍射效应不能忽略的菲涅尔衍射区（包括在几何投影区以后的所有区域），以及衍射图样基本形状保持不变的夫琅禾费区。这种衍射区的划分
方法认为，夫琅禾费衍射只是菲涅尔衍射的特殊情况。
菲涅尔衍射和夫琅禾费衍射是傍轴近似下的两种衍射情况，二者的区别条件是观察屏到衍射屏的距离z1与衍射孔的线度
x1，y1
之间的相对大小。
思考题
1什么是近场，远场，远场和近场衍射各有什么特点？
答：近场是菲涅尔衍射，发散光入射，会聚光出射；远场是夫琅禾费衍射，平行光入射，平行光出射，当衍射孔径平面到观察平面距离r无限时时远场衍射。
2根据仿真结果，总结圆孔的近场及远场衍射的特点
？
答：观测时，会看到菲涅耳衍射所产生的圆孔成像，大小与形状会与原来的圆孔不一样，即是说边缘多少会有一些锯齿在，但是夫琅和费衍射的成像则只有大小的改变，这是因为远场的波动比较接近平行光束及平面波的性质在远场衍射中，如果观测屏幕在圆孔不动时往後移动，则产生的条纹会一致地改变大小。但近场衍射则不会这样，衍射条纹的大小与影状都会改变。
3结合菲涅尔半波带法说明近场和远场的距离划分？
答：菲涅耳圆孔衍射露出的波带数N与圆孔半径的关系如下
当在一定范围内变化时，波带数会随着变化，此区域内衍射称之为菲涅耳衍射；若大到一定程度时，，波带数将是一个定值，此区域内衍射称之为夫朗和费衍射。
五、仿真小结
这是本次实习的最后一个仿真，也是最难的一个。最早设计的程序虽然能够运行，但是由于程序设计之初所运用的原理导致仿真精度和运算时间上无法很好的协调。每次运行都要十几分钟，因此我从原理上出发利用傅里叶光学的内容改进了程序算法直接提升了运行效率。总之本次仿真作为收官之作，结果还是令我满意的，不仅成功改良程序还学到了新的知识。
附录：
lamda=632.8e-6;
k=2*pi/lamda;
z=1000000;
%先确定衍射屏
N=300;
%圆屏采样点数
r=50;
a=0;
b=0;
I=zeros(N,N);
[m,n]=meshgrid(linspace(-N/2,N/2-1,N));
D=((m-a).^2+(n-b).^2).^(1/2);
i=find(D<=r);
I(i)=1;
%孔半径范围内透射系数为1
subplot(2,2,1);
%圆孔图像画在2行2列的第一个位置
imagesc(I)
%画衍射屏的形状
axis
image
;
title(
衍射屏形状
)
;
L=300;
M=300;
%取相同点数用于矩阵运算
[x,y]=meshgrid(linspace(-L/2,L/2,M));
h=exp(j*k*z)*exp((j*k*(x.^2+y.^2))/(2*z))/(j*lamda*z);%接收屏
H
=fftshift(fft2(h));
B=fftshift(fft2(I));
%圆孔频谱
G=H.*B;
%公式中为卷积，空间域中相卷相当于频域相乘
U=
fftshift(ifft2(G));
%求逆变换，得到复振幅分布矩阵
Br=(U/max(U));
%归一化
subplot(2,2,2);
imshow(abs(U));
axis
image;
colormap(gray);
title(
衍射后的图样
);
subplot(2,2,3);
mesh(x,y,abs(U));
%画三维图形
subplot(2,2,4);
plot(abs(Br));
篇2：南京大学C波段双偏振雷达_垂直90度观测_分析总结报告_第1部分_20XX年8月5日_7日_21日补充
南京大学C波段双偏振雷达_垂直90度观测_分析总结报告_第1部分_2014年8月5日_7日_21日补充 本文关键词：南京大学，偏振，波段，观测，垂直
南京大学C波段双偏振雷达_垂直90度观测_分析总结报告_第1部分_2014年8月5日_7日_21日补充 本文简介：南京大学C波段双偏振雷达垂直90度观测分析总结报告（第1部分）南京大学2014年8月目录1概述12垂直90度观测时雷达的相关参数13Matlab分析程序24ZDR和PHIDP偏差与方差的分析结果34.1ZDR偏差与方差的分析结果34.1.1IOP144.1.2IOP264.1.3IOP364.1.4
南京大学C波段双偏振雷达_垂直90度观测_分析总结报告_第1部分_2014年8月5日_7日_21日补充 本文内容：
南京大学C波段双偏振雷达
垂直90度观测
分析总结报告
（第1部分）
南京大学
2014年8月
目
录
1
概述1
2
垂直90度观测时雷达的相关参数1
3
Matlab分析程序2
4
ZDR和PHIDP偏差与方差的分析结果3
4.1
ZDR偏差与方差的分析结果3
4.1.1
IOP14
4.1.2
IOP26
4.1.3
IOP36
4.1.4
IOP47
4.1.5
IOP59
4.1.6
IOP611
4.1.7
IOP713
4.1.8
IOP815
4.1.9
IOP917
4.1.10
IOP1019
4.1.11
IOP1125
4.2
PHIDP偏差与方差的分析结果27
5
对整个观测试验进行分析（分析偏差的中期变化）32
5.1
ZDR系统偏差、方差的变化33
5.2
PHIDP系统偏差、方差的变化38
5.3
列表汇总43
5.4
其它几个参数的方差44
6
结论47
6.1
通过垂直90度观测，说明这部雷达的技术指标如何？47
6.2
对观测数据应该如何将偏差订正47
7
下一步工作48
7.1
其它几个剩余的IOP继续分析48
7.2
推广到别的双偏振雷达，也对垂直90度的数据进行同样的分析48
7.3
将各参数的方差和理论值进行对比分析，看是否符合理论。48
8
基于IQ数据重新计算的基数据的部分分析结果48
8.1
ZDR偏差与方差的分析结果48
8.2
PHIDP偏差与方差的分析结果52
1
概述
2014年6、7月份，由北京敏视达公司制造的南京大学C波段双偏振雷达在安徽长丰站进行了持续的降水观测，总计有11个IOP观测。在每次体扫观测中，都进行了垂直90度扫描，录取了宝贵的数据。
所有IOP的基数据都保存在了“2014年长丰站所有IOP的基数据:/NJU_CPol/BaseData”
目录下。.bz2
是雷达保存的压缩数据格式；.
AR2是将bz2解压之后的数据；.mat是matlab能直接读取的数据格式（参见Fun_Read_MSD_Radar_BaseData.m程序）。
本文对垂直90度观测的数据进行了系统的分析，总共有5687个基数据文件，从中发现有1068个有合适降水的垂直90度的回波。然后绘制了各种统计分析的图片和曲线。
从这些垂直90度的数据得出以下初步结论：
l
6月14日的数据偏差有个跳变（这是由于6月10日更换了旋转关节造成的）
l
各个时刻的PHIDP的系统偏差非常稳定；各个时刻的ZDR的系统偏差基本稳定；
l
ZDR的随机差（方差）为0.25~0.4dB（和谱宽有关系，方位上采用64点积分时）；
l
PHIDP的随机差在2~3度（方位上采用64点积分时）
l
但是，ZDR的偏差还是存在“小时”量级的波动的，波动范围约为0.3~0.4dB。原因还需进一步分析IQ数据才能确认。
l
第9个IOP，有11个VCP的ZDR系统偏差为0.3dB；第10个IOP，前7个VCP的ZDR的系统偏差为0.2dB。这些偏差结果和其它IOP下、其它VCP的结果有0.3~0.4dB的差异。需要进一步分析IQ数据，找出突然发生变化的原因。
l
第10个IOP（台风）中，有几个体扫的ZDR的随方位的起伏较大，可能是风和雨的共同作用导致的
8月7日，补充了基于IQ数据重新计算的基数据的部分分析结果（只有IOP8和IOP9的部分结果，其余的IOP的IQ数据还没有计算出来），见第8
章。
8月21日，将所有的IOP的数据都补充完整了。
2
垂直90度观测时雷达的相关参数
在5月27日进行了指标测试，相关的指标如下：（详见《南大C波段雷达扫描表_2014年5月30日.doc》和《南大C波段雷达标校情况_2014年5月30日.doc》）
l
天线波束宽度：1.19度（H）、
1.18度（V）
l
天线增益：41.5dB
l
馈线损耗：H:
3.6dB、V:3.3dB（含波导、旋转关节、TR管等，不含大功率的功分器）
l
发射载频：5625MHz
l
发射脉冲宽度：0.5μs
l
发射峰值功率：254kW
l
发射采用初始相位随机编码的方式
l
重复周期：2000Hz（奇怪的事情是：从6月20日开始，重复周期竟然变为了1000Hz）
l
不模糊速度：26.66m/s
l
接收机噪声系数：2.75dB（H）、2.78dB（V）
l
接收机带宽：2MHz
l
接收机灵敏度：-108.24dBm（水平通道）
-108.21dBm（垂直通道）
l
库长：150m（从6月11日，变成了75m，即距离上不再进行2点积分了）
l
库数：133（从6月11日，变成了266个）
l
方位平均数（即脉冲累积数）：128（奇怪的事情是：从6月20日开始，脉冲累积数竟然变为了64）
l
天线扫描速度：15.625deg/s
l
雷达常数C(dB)：79.93（水平通道）
79.63（垂直通道）
l
雷达系统灵敏度：[email protected]（水平通道）
[email protected]（垂直通道）
3
Matlab分析程序
我们编写了能进行批量分析的Matlab程序，自动对基数据进行读取、解析、判断90度有无合适的降水，然后绘制PPI图、随高度变化的图、随方位变化的图，并进行统计，绘制直方图。
同时，以.txt文本形式将该体扫下的分析结果记录下来，保存在Result子目录中，以便后续的分析。同时记录了ZDR、std(ZDR)
、SNR、W、dBZ、std(dBZ)
、V、std(V)
、CC、std(CC)
、PHIDP、std(PHIDP)
等参数随高度（距离）的变化。
一个典型的txt文件如下所示：
从上面的txt文件中，可以很容易的看出某个基数据中，各个参数的偏差和方差的情况。
同时，txt文件中的内容，也保存为mat格式，以方便进一步的分析。
那么如何判断90度有合适的降水呢？
我们选取零度层以下，且不含地物的相关系数CC的数据来进行判断。
首先选取距离（即高度）从1000m~3000m的CC值，总共有约2500个数据。如果在这些数据中，有超过5%的CC小于0.95（包括没有回波的数据），则说明该文件不是降水回波，该文件就不再进行后续的分析。
经过对5月31日、6月1日的数据的验证，这种判断方法还是非常准确可靠的（通过人工识别出有合适降水的情况，与通过上述准则自动判断的结果是一致的）。
4
ZDR和PHIDP偏差与方差的分析结果
4.1
ZDR偏差与方差的分析结果
所有图片在“L:/NJU_CPol/Analysis/IOPx/90deg/ZDR”目录下。
4.1.1
IOP1
其中一个5月31日的基数据的图片如下：
注意：上面这张图中，是选取了高度在1000~2000m之间的数据，然后先进行距离上的平均（平均的目的是为了减少ZDR的随机起伏，从而能更好的发现ZDR偏差随方位的变化），然后再按照横轴是方位来绘制的。
注意：上面这张图中，是选取了高度在1000~2000m之间的数据，然后将两维的数据用
reshape函数变为一维的数组（注意：这里因为要计算方差，因此绝对不能进行距离上的平均），然后用hist函数统计得到直方图，最后用mean、median、std函数得到平均值、中间值和方差的。
4.1.2
IOP2
没有合适的垂直90度的降水，因此这个IOP不分析。
4.1.3
IOP3
由于
6月10日更换过关节，因此这里再绘制一个6月14日的降水数据的ZDR图：
对比关节更换前后的结果可以看出，ZDR的系统偏差从原来的0.16dB，改为了0.70dB，但方差没有变。
另外，从ZDR随方位变化的图可以看出，不管是国产还是进口关节，都不存在ZDR系统偏差随方位变化的现象（但在别的雷达中，这是一个普遍的现象）。
4.1.4
IOP4
4.1.5
IOP5
4.1.6
IOP6
4.1.7
IOP7
4.1.8
IOP8
4.1.9
IOP9
4.1.10
IOP10
第10个IOP，前7个VCP的ZDR的系统偏差为0.2dB。这些偏差结果和其它IOP下、其它VCP的结果有0.3~0.4dB的差异。
从UTC时间08:07:25开始，ZDR的系统偏差就正常了：
但是，其中有几个体扫的ZDR的随方位的起伏较大，可能是风和雨的共同作用导致的。
4.1.11
IOP11
4.2
PHIDP偏差与方差的分析结果
所有图片在“L:/NJU_CPol/Analysis/IOPx/90deg/PHIDP”目录下。
其中一个5月31日的基数据的图片如下：
注意：
6月10日更换过关节，因此这里再绘制一个6月14日的降水数据的PHIDP图：
对比关节更换前后的结果可以看出，PHIDP的系统偏差从原来的25.79度，变为了-97.5度，但方差没有变。另外，都存在PHIDP系统偏差随方位而波动的现象，波动范围为±3度，但两者波动的曲线不同。
我们再看看其它日子的PHIDP随方位的变化：
从上面几张不同日子的PHIDP系统偏差随方位而波动的图可以看出，这个波动是稳定的，因此可以很容易的修正掉。
5
对整个观测试验进行分析（分析偏差的中期变化）
读取前面分析得到的mat文件，然后选取特定高度的数据（这里取2500m的高度，经过观察回波发现，在这个距离处，地物的回波比较弱）的数据，绘制成横轴是时间、数轴是数值的曲线的图，以便可以清晰的看出在整个观测试验中，ZDR、PHIDP等参数系统偏差的变化情况。
生成的图片保存在
Analysis/IOPx/90deg/Result/Pic目录下。
5.1
ZDR系统偏差、方差的变化
下面，将各个IOP观测试验下，ZDR的系统偏差、方差，随各个VCP扫描的变化曲线绘制如下：
ZDR存在“小时”量级的变化
此处ZDR变化，是因为时间上相隔了12个小时
此处ZDR变化，原因未明？
备注：
从上面几张图中可以看出，ZDR的偏差还是存在“小时”量级的波动的。但是，这个波动的原因，是由于雷达的硬件的确发生了增益、灵敏度的变化呢？还是由于RVP9在计算的时候，其雷达常数等参数发生了改变造成的呢（就像RVP9输出的强度值会在VCP之间突然发生变化，但IQ数据重新计算的强度值就没有变化）？
需要进一步分析IQ数据才能得出结论。
另外，如果ZDR的偏差的确有这么大的变化（也就是说，雷达的硬件的确有这么大的变化），则今后偏差修正，以及基于ZDR的目标识别的时候，就麻烦了，必须想办法如何降低这个偏差的波动造成的影响。
可能也跟硬件中，接收机的温度没有稳定有关。这个可以通过对IQ数据进行茅草的统计来分析是否是这个原因。
另外，由于雷达处于自动标定有效的状态（尽管我们要求不要将标定的结果传入信号处理器，但实际上是传入的），也就是说，每次VCP体扫结束之后，都会执行自动标定，并将标定的结果传入信号处理器，作为计算的基准。因此，也有很大的可能是自动标定的问题。
5.2
PHIDP系统偏差、方差的变化
下面，将各个IOP观测试验下，PHIDP的系统偏差、方差，随各个VCP扫描的变化曲线绘制如下：
5.3
列表汇总
将所有的IOP的结果下，ZDR的偏差、方差，PHIDP的偏差和方差，列表显示出来
表1
ZDR和PHIDP的偏差、方差的结果
IOP
起止时间
90度观测次数
方位*距离积分次数
ZDR偏差
(dB)
ZDR方差(dB)
PHIDP偏差(度)
PHIDP方差(度)
谱宽(m/s)
1
5.31~6.1
26
128*2
0.12
0.21
27.33
2.13
1.09
2
6.11
0
128*1
3
6.14~6.17
249
128*1
0.98
0.29
-93.29
2.38
0.93
4
6.20~6.21
26
64*1
0.86
0.35
-95.76
2.74
0.75
5
6.24~6.27
253
64*1
1.02
0.31
-95.49
2.53
0.84
6
6.30~7.3
130
64*1
0.85
0.27
-95.16
2.30
0.92
7
7.4~7.5
93
64*1
0.91
0.28
-94.17
2.38
0.99
8
7.11~7.13
130
64*1
0.97
0.30
-95.17
2.49
1.02
9
7.15~7.18
73
64*1
0.54
0.25
-94.78
2.21
0.92
10
7.24~7.25
55
64*1
0.61
0.30
-93.33
2.51
1.01
11
7.30
33
64*1
0.70
0.29
-93.13
2.42
1.14
备注：6月14日的数据偏差有个跳变（这是由于6月10日更换了旋转关节造成的）
从上表，我们终于可以得出一点结论了：
l
各个时刻的PHIDP的系统偏差非常稳定；各个时刻的ZDR的系统偏差基本稳定；
l
ZDR的随机差（方差）为0.25~0.4dB（和谱宽有关系，方位上采用64点积分时）；
l
PHIDP的随机差在2~3度（方位上采用64点积分时）
l
但是，ZDR的偏差还是存在“小时”量级的波动的，波动范围约为0.3~0.4dB。原因还需进一步分析IQ数据才能确认。
l
第9个IOP，有11个VCP的ZDR系统偏差为0.3dB；第10个IOP，前7个VCP的ZDR的系统偏差为0.2dB。这些偏差结果和其它IOP下、其它VCP的结果有0.3~0.4dB的差异。需要进一步分析IQ数据，找出突然发生变化的原因。
l
第10个IOP（台风）中，有几个体扫的ZDR的随方位的起伏较大，可能是风和雨的共同作用导致的
5.4
其它几个参数的方差
下面，绘制一下其它几个参数随VCP的变化（这里仅将IOP1的结果展示出来）。所有图片在“L:/
NJU_CPol/Analysis/IOPx/90deg/Result/Pic/”目录下。
6
结论
6.1
通过垂直90度观测，说明这部雷达的技术指标如何？
l
6月14日的数据偏差有个跳变（这是由于6月10日更换了旋转关节造成的）
l
基本上，各个时刻的ZDR和PHIDP的系统偏差都非常稳定；
l
ZDR的随机差（方差）为0.25~0.4dB（和谱宽有关系，方位上采用64点积分时）；
l
PHIDP的随机差在2~3度（方位上采用64点积分时）
l
但是，ZDR的偏差还是存在“小时”量级的波动的，波动范围约为0.3~0.4dB。原因还需进一步分析IQ数据才能确认。
l
第9个IOP，有11个VCP的ZDR系统偏差为0.3dB；第10个IOP，前7个VCP的ZDR的系统偏差为0.2dB。这些偏差结果和其它IOP下、其它VCP的结果有0.3~0.4dB的差异。需要进一步分析IQ数据，找出突然发生变化的原因。
l
第10个IOP（台风）中，有几个体扫的ZDR的随方位的起伏较大，可能是风和雨的共同作用导致的
6.2
对观测数据应该如何将偏差订正
根据基数据的观测的时间，查找第5.3
章
表1中的ZDR和PHIDP的偏差值，然后进行订正。
注意：第2个IOP
6月11日的数据，由于没有合适的垂直90度降水，因此只能按照第3个IOP
的偏差进行订正。（因为6月10日更换了旋转关节）
7
下一步工作
7.1
其它几个剩余的IOP继续分析
7.2
推广到别的双偏振雷达，也对垂直90度的数据进行同样的分析
7.3
将各参数的方差和理论值进行对比分析，看是否符合理论。
8
基于IQ数据重新计算的基数据的部分分析结果
8.1
ZDR偏差与方差的分析结果
所有图片在“L:/NJU_CPol/Analysis_DSP/IOPx/90deg/ZDR”目录下。
其中一个7月11日的基数据的图片如下：
然后，我们分析偏差的中期变化。
读取前面分析得到的mat文件，然后选取特定高度的数据（这里取2500m的高度，经过观察回波发现，在这个距离处，地物的回波比较弱）的数据，绘制成横轴是时间、数轴是数值的曲线的图，以便可以清晰的看出在整个观测试验中，ZDR、PHIDP等参数系统偏差的变化情况。
生成的图片保存在
Analysis_DSP/IOPx/90deg/Result/Pic目录下。
下面，将各个IOP观测试验下，ZDR的系统偏差、方差，随各个VCP扫描的变化曲线绘制如下（目前只分析了IOP8和IOP9的数据，而且IOP9只分析了部分数据）：
备注：
从IQ数据重新计算的结果，其ZDR的系统偏差和RVP9给出的有一点差别。这个原因是因为：基于IQ数据重新计算，其计算中的雷达常数、接收机灵敏度等值，是严格按照第2章的参数设定的。而RVP9计算中的参数，其具体的值不清楚（黑盒子）。
从IOP8和IOP9的两张图中可以看出，ZDR的偏差还是存在“小时”量级的波动的。
由于IQ计算中，雷达常数、接收机灵敏度等值是固定的。因此也就说明此时雷达的确发生了雷达常数、增益、灵敏度等硬件上的变化。（是否是由于雷达刚开机，接收机的温度没有达到稳定的状态？）
（当然，还有一种可能性：此时外部降水的ZDR的特性就是发生了变化。）
8.2
PHIDP偏差与方差的分析结果
所有图片在“L:/NJU_CPol/Analysis_DSP/IOPx/90deg/PHIDP”目录下。
其中一个7月11日的基数据的图片如下：
同样的，我们分析偏差的中期变化。
下面，将各个IOP观测试验下，PHIDP的系统偏差、方差，随各个VCP扫描的变化曲线绘制如下（目前只分析了IOP8和IOP9的数据）：
54