工程光学综合练习2(改后)

工程光学综合练习一 干涉仿真

一、要求

3-4人组成小组，对下面给出的各题目利用Matlab 等工具进行仿真。练习结束时每组提交一份报告及仿真程序。在报告中应注明各仿真结果所对应的参数，如相干光源间距、光入射倾角等。

二、仿真题目

一）、对于杨氏双缝干涉，改变双缝的缝宽和缝间距，观察干涉图样变化 （1）Matlab 程序及注释 建立如下图所示坐标系

两束平面波干涉采用的接收屏是x-y 平面 在接收屏上坐标(xs,ys)点处，易求得 ：

光程

2

22

2)2

(D ys d xs r +++=

相位差

光强2

)2

cos(

4phi B =

程序： clear

lam=500e-9; %设定波长为500nm

d=input('两缝之间宽度d=')/1000; %输入初始双缝宽d

D=1; %设定光源中心到接收屏的距离 A=input('光波的振幅A=');

xm=5*lam*D/d;ym=xm; %设定光屛的范围

n=101;xs=linspace(-xm,xm,n); %把光屛的x 方向分成101点 ys=linspace(-ym,ym,n); %把光屛的y 方向分成101点

for i=1:n for j=1:n

r1=sqrt((xs(i)-d/2)^2+ys(j)^2+D^2); %光程r1 r2=sqrt((xs(i)+d/2)^2+ys(j)^2+D^2); %光程r2

xiangwei=2*pi*(r2-r1)/lam ; %屏上各点的相位差 I(i,j)=4*A^2*cos(xiangwei/2)^2; %屏上各点光强 end

end

N=255; %确定用的灰度等级为255

级

Br=(I/4.0)*N; %使最大光强对应于最大灰

度级（白色）

figure(1) %绘制图形

image(xs,ys,Br); %画干涉条纹

xlabel('光屏x方向');

ylabel('光屏y方向');

title('杨氏双缝干涉条纹');

运行结果如图1—1

图1—1

将程序中振幅A初始设定为1，然后在改变振幅A为2，观察运行之后观察干涉图样；程序中通过改变振幅来改变强度，而改变强度即改变了缝宽，运行结果如图1—2

将d的初始缝宽设定为2mm，即输入数值2，然后改变缝间距，设定缝宽为3mm，输入数值3再次运行，运行结果如图1—3

图1—2 图1—3

二）、对于杨氏双孔干涉，改变双孔的直径和孔间距，观察干涉图样变化

D=1; %设定光源中心到接收屏的距离

d=input('两空间距(d)=')/1000;

%设定两光源间距为0.002mm

R1=input('空2的半径(R2)=');

R2=input('空2的半径(R2)=');

I0=1;

ny=101;

lam=5e-7; %设定波长为500nm

I1=I0*R1*R1*R1*R1;

I2=I0*R2*R2*R2*R2;

ymax=5*lam*D/d; %屏幕上y的最大范围

y=linspace(-ymax,ymax,ny) %设定光屛的范围

x=y;

[x,y]=meshgrid(x,y); %在坐标系中绘制网格图

r1=sqrt((x-d/2).^2+y.^2+(D-d/2).^2); %光源1到接收屏的光程r1

r2=sqrt((x+d/2).^2+y.^2+(D+d/2).^2); %光源2到接收屏的光程r2 xiangwei=(r1-r2)*2*pi/lam; %光屏上各点相位差

I=I1+I2+2*sqrt(I1*I2)*cos(xiangwei); %光屏上各点光强

N=255; %设置灰度等级为255

Br=(I/4.0)*N; %使最大光强对应于最大灰度级

figure(1);

image(Br);

colormap(gray(N)); %根据光强生成图像

设置两光源的初始距离d为0.002mm，即输入的数值为0.002，之后变更d为0.02mm，即输入的数值为0.02，达到增大光源距离的目的，生成图像如图2-1所示，孔的直径大小R1和R2初始大小都设置为1，即R1=1，R2=1，之后变更为R1=1，R2=3观察运行结果，所得图像如图2-2所示

图2-1 图2-2

三）、改变下列光波场分布，观察干涉图样变化

图1 图2 图3 1、如图1所示，两平面光波叠加，改变光波振幅比、两光波夹角，观察在接收屏上的干涉图样变化；

（1）Matlab程序及注释

A1=input('第一束光的振幅(m)=');

A2=input('第二束光的振幅(m)=');

th=input('夹角=');

lam=input('光波长（nm）=')/1000000000;

ymax=0.00000005;

ny=161;

y=linspace(-ymax,ymax,ny);

z=y;

j=1:ny

I1=A1^2; %计算第一束光的光强 I2=A2^2; %计算第二束光的光强 phi=pi*th*y(i)/lam; %计算相位差

I(i,j)=I1+I2+2*sqrt(I1*I2)*cos(phi); %根据公式计算两束光叠加后的光强 end

nclevels=255;

br=I.*255/max(max(I)); figure(1); image(y,z,br); xlabel('y(n)'); ylabel('z(n)'); title('干涉条纹');

axis([-ymax,ymax,-ymax,ymax]); colormap(copper(nclevels));

（2） 仿真结果截图

设波长λ=500，改变参数：

振幅1A =10，2A =10，夹角θ=30 如图3-1-1 振幅1A =10，2A =30，夹角θ=30，如图3-1-2 振幅1A =10，2A =10，夹角θ=20，如图3-1-3

图3-1-1 图3-1-2

图3-1-3

结论：振幅比增大，条纹间距增加；夹角减小，条纹间距增加

2、如图2所示，两点光源前后放置，改变其间距，观察在接收屏上的干涉图样变化；

（1）Matlab程序及注释

（2）仿真结果截图

3、如图3所示，两点光源并排放置，改变其聚散性（会聚球面波、发散球面波）和间距，观察在接收屏上的干涉图样变化。

（1）Matlab程序及注释

d=input('两点光源之间间距(m)=');

D=input('右侧点距屏距离(m)=')+(d/2);

lanbda=input('光波长(nm)=')/1000000000;

ymax=0.0055;

ny=161;

y=linspace(-ymax,ymax,ny);

z=y;

for i=1:ny

for j=1:ny

l1=sqrt(y(i)^2+(D-d/2)^2+z(j)^2); %计算r1

l2=sqrt(y(i)^2+(D+d/2)^2+z(j)^2); %计算r2

phi=2*pi*(l2-l1)/lanbda; %计算相位差

I(i,j)=exp(2)/l1+exp(2)/l2+2*exp(2)/l1/l2*cos(phi); %根据公式计算叠加的光强end

end

nclevels=255;

br=I.*255/max(max(I));

figure(2);

image(y,z,br);

xlabel('y(n)');

ylabel('z(n)');

title('单色光双缝干涉条纹');

axis([-ymax,ymax,-ymax,ymax]);

colormap(copper(nclevels));

（2）仿真结果截图

设波长λ=500，与屏幕距离D=5，改变参数：

两光源间距离d=5 ，如图3-2-1

两光源间距离d=2 ，如图3-2-2

图3-2-1 图3-2-2

结论：两光源距离增加，条纹间距减小

3．设波长λ=500，与屏幕距离D=10，且屏幕与两光源中心连线平行，改变参数（1）Matlab程序及注释

d=input('两点光源之间间距(m)=');

D=input('与屏距离(m)=');

lanbda=input('光波长(nm)=')/1000000000;

ymax=0.55;

ny=1001;

y=linspace(-ymax,ymax,ny);

z=y;

i=(ny-1)/2;

for i=1:ny

for j=1:ny

l1=sqrt(y(i)^2+D^2); %第一个点光源到屏幕某点的距离

l2=sqrt((d-y(i))^2+D^2); %第二个点光源到屏幕某点的距离

phi=2*pi*(l2-l1)/lanbda; %计算相位差

I=I1+I2+2*sqrt(I1*I2)*cos(phi); %根据公式计算叠加后的光强

end

end

nclevels=255;

br=I.*255/max(max(I));

figure(2);

image(y,z,br);

xlabel('y(n)');

ylabel('z(n)');

title('单色光双缝干涉条纹'); axis([-ymax,ymax,-ymax,ymax]); colormap(gray(nclevels));

（2）仿真结果截图

两光源间距离d=2，如图3-3-1 两光源间距离d=10，如图

3-3-2

图3-3-1 图3-3-2

结论：光源距离增大，条纹距离减小；当其中一个光源聚散性改变后，条纹由椭圆形变为双曲线。

四）、如图4-6所示，改变平面光波场分布，观察干涉图样变化

图4 图5 图6 1.如图4，平行光垂直于y轴入射

（1）Matlab程序及注释

clear

lamda=500e-9; %波长

z=5 ; %点光源距离theta=0 ; %平行光波前与x

轴夹角

ymax=0.005;

xmax=ymax;

yn=101;ys=linspace(-ymax,ymax,yn); %把光屛的x方向分成101点

xn=101;xs=linspace(-xmax,xmax,xn); %把光屛的y方向分成101点

for i=1:yn

for j=1:xn

r=sqrt(xs(j).^2+ys(i).^2+z^2)-z*sin(theta)-xs(j)*cos(theta); %计算光程差Phi=2*pi*r/lamda; %计算相位差b(i,j)=4*cos(Phi/2).^2; %光强

end

end

clf;figure(gcf);

nclevels=255; %确定用的灰度等级为225级

br=(b/4)*nclevels;

image(xs,ys,br); %显示colormap(gray(nclevels));

（2）仿真结果如图4-1

2.如图5所示，平行光与y轴正半轴成45度角入射

（1）Matlab程序及注释

clear

lamda=500e-9; %波长

z=5; %点光源距离

theta=45; %平行光波前与x轴夹角

ymax=0.005;

xmax=ymax;

yn=101;ys=linspace(-ymax,ymax,yn) ; %取点数

xn=101;xs=linspace(-xmax,xmax,xn);

for i=1:yn

for j=1:xn

r=sqrt(xs(j).^2+ys(i).^2+z^2)-z*sin(theta)-xs(j)*cos(theta); %计算光程差Phi=2*pi*r/lamda; %计算相位差b(i,j)=4*cos(Phi/2).^2; %光强

end

end

clf;figure(gcf);

nclevels=255; %确定用的灰度等级为225级

br=(b/4)*nclevels;

image(xs,ys,br); %显示colormap(gray(nclevels));

（2）仿真结果如图4-2

3. 如图6所示，平行光与y轴负半轴成45度角入射

（1）Matlab程序及注释

clear

lamda=500e-9; %波长

z=5; %点光源距离theta=45; %平行光波前与x轴夹角

ymax=0.005;

xmax=ymax;

yn=101;ys=linspace(-ymax,ymax,yn); %取点数

xn=101;xs=linspace(-xmax,xmax,xn);

for i=1:yn

for j=1:xn

r=sqrt(xs(j).^2+ys(i).^2+z^2)-z*sin(theta)+xs(j)*cos(theta); %计算光程差Phi=2*pi*r/lamda ; %计算相位差

b(i,j)=4*cos(Phi/2).^2; %光强

end

end

clf;figure(gcf);

nclevels=255; %确定用的灰度等级为225级

br=(b/4)*nclevels;

image(xs,ys,br); %显示colormap(gray(nclevels));

（2）仿真结果如图4-3

图4-1 图4-2

图4-3

五）、用平行光（点光源+准直镜）不同形状楔板，观察干涉图样（1）

图5-1-1

图5-1-2 图5-1-3 如图5-1-2，设置将要显示的楔板部分垂直的高度差为设0.0005m，设置楔板的夹角为10-6rad，，在垂直方向上将楔板分割取微元，计算出水平方向上

各点对应的光程差。

程序代码如下：

clear %清除原有变量

Lambda=500*(1e-9); %设定波长

theta=1*(1e-6); %设定倾角

ni=500; %微元个数

ds=linspace(0,0.0005,ni); %竖直方向分割

for k=1:ni

y(k)=ds(k)/sin(theta); %水平方向的对应坐标

Delta=2*ds(k)+Lambda/2; %对应点的光程差

Phi=2*pi*Delta/Lambda; %对应点的相位差

B(k,:)=4*cos(Phi/2).^2; %光强与相位差相联系

end

figure(gcf); %以下为形成图像过程，将光强与图像灰度相对应

NCLevels=250;

Br=(B/4.0)*NCLevels;

image(0,y,Br);

colormap(gray(NCLevels));

title('二维强度分布');

运行程序，仿真结果如图5-1-3所示。

（2）

图5-2-1

图5-2-2 图5-2-3 如图5-2-2，设置将要显示的介质圆弧部分半径为10cm，设置介质最大厚

度为5cm，计算出水平方向上各点对应的光程差，显示中心两侧各1mm的范围内干涉图样。

程序代码如下：

clear %清除原有变量

Lambda=500*(1e-9); %设置波长

r=10*(1e-2); %设置介质的尺寸

t=5*(1e-2);

ni=1000; %分割份数

y=linspace(-0.001,0.001,ni); %y方向上进行分割

for k=1:ni

h(k)=t-(r-sqrt(r^2-y(k)^2)); %对应点的介质厚度

Delta=2*h(k)+Lambda/2; %对应点的光程差

Phi=2*pi*Delta/Lambda; %对应点的相位差

B(k,:)=4*cos(Phi/2).^2; %对应点的光强

end

figure(gcf); %以下为形成图像过程，将光强与图像灰度相对应

NCLevels=250;

Br=(B/4.0)*NCLevels;

image(0,y,Br);

colormap(gray(NCLevels));

title('二维强度分布');

运行程序，仿真结果如图5-2-3所示。

（3）

图5-3-1

图5-3-2 图

5-3-3

如图5-3-2，设置将要显示的介质圆弧部分半径为5m，设置介质上半部分柱体厚度为2m，计算出水平方向上各点对应的光程差，显示位于中心的半径为5mm的圆的范围内的干涉图样。

程序代码如下：

clear %清除原有变量lambda=500*1e-9;R=5;H=2; %设置波长和介质尺寸ni=2000; %设置分割份数

x=linspace(-0.005,0.005,ni); %在x，y方向上进行分割

y=linspace(-0.005,0.005,ni);

for i=1:ni

for j=1:ni

r=x(i).^2+y(j).^2; %点距介质中心的直线距离

delta=2*H+2*(R-sqrt(R^2-r))+lambda/2; %相应点的光程差

phi=2*pi*delta/lambda; %相应点的相位差

B(i,j)=4*cos(phi/2).^2; %相应点的光强

end

end

figure(gcf); %以下为形成图像过程，将光强与图

像辉度相对应

NCLevels=250;

Br=(B/4.0)*NCLevels;

image(x,y,Br);

colormap(gray(NCLevels));

title('二维强度分布');

运行程序，仿真结果如图5-3-3所示。

（4）

图5-4-1

图5-4-2 图5-4-3

如图5-4-2，拟选用圆柱体为介质的基础，下方横截面的下边界为正弦曲线。设置将要介质圆柱高度为10cm，将厚度变化限制在波长量级（正弦函数前乘以10-6作为系数），设置正弦函数周期为5mm，计算出水平方向上各点对应的光程差，显示中心左右一个正弦周期范围内的干涉图样。

程序代码如下：

lambda=500*1e-9; %设置波长

H=0.1; %介质的基础厚度

ni=2000; %分割份数

T=0.005; %显示范围，也是正弦曲线的一个周期omega=2*pi/T; %正弦函数的角频率

x=linspace(-T,T,ni); %在x和y方向上进行分割

y=linspace(-T,T,ni);

for i=1:ni

for j=1:ni

r=sqrt(x(i).^2+y(j).^2); %点距中心的直线距离

delta=2*H+sin(omega*r)*1e-6; %相应点的光程差

phi=2*pi*delta/lambda; %相应点的相位差

B(i,j)=4*cos(phi/2).^2; %相应点的光强

end

end

figure(gcf); %生成图像

NCLevels=250;

Br=(B/4.0)*NCLevels;

image(x,y,Br);

colormap(gray(NCLevels));

title('二维强度分布');

运行程序，仿真结果如图5-4-3所示。

六）、改变平行平板折射率，观察干涉图样

（1）Matlab程序及注释

clear

lamda=[600 600]*1e-5; %波长

RGB=[1,0,0;0,0,0;]; %设置红光RGB值h=0.05; %点光源距离

pc=zeros(500,500,3); %设置光屏

r=0.54; %反射系数

n=1.5; %折射率（可变）d=1.7; %薄膜厚度

[x,y]=meshgrid(linspace(-0.005,0.005,500));

z=x+1i*y;

theta=atan(abs(z)/h); %设置入射角

for k=1:2

phi=4*pi*n*d*cos(theta)/lamda(k); %相位差

It=1./(1+4*r^2*sin(phi/2).^2/(1-r^2)^2); %光强

p(:,:,1)=It*RGB(k,1);

p(:,:,2)=It*RGB(k,2);

p(:,:,3)=It*RGB(k,3); %将包含颜色信息的光强用矩阵保存下来

pc=pc+p; p=[];

end; %显示

Br=1/max(max(max(pc))); %调整矩阵元素的最大值为1的系数pcl=pc*Br; %调节

imshow(pcl,[]) %显示仿真结果

title('平行平板等倾干涉结果')

（2）仿真结果截图

n=1.33 n=1.5

n=1.65 n=1.7