matlab-光学薄膜讲解学习

安徽工业大学光信息科学与技术专业《光学软件课程设计》
课程报告
年级：光082
姓名：
指导教师：黄仙山（2900）
日期：2011年6月28日
报告正文：
一、 设计题目
一维周期性光学薄膜结构中光学电磁传输特性的理论研究
二、 设计理论
薄膜特征矩阵
在界面1
2
'
2
2111111'2
1111cos cos cos cos i r i t i r i i r t r i H H H H H E E E E E θθθθ-=-=+=+=
由《电磁场与电磁波》H 与E
的关系=得：E H μ
ε=，
取0μμ=则21'210
1011001cos )(cos )(i r t i r i n E E n E E H θμεθμε-=-= 在界面2
2
20
21220022
222cos cos )(t G t i r i t r i n E n E E H E E E E θμεθμε=-==+= 在不考虑薄膜对光能的吸收时，
1
2E t i E 和，2'
2E r r E 和的关系如下：
)
12'2112exp()exp(δδi E E i E E r r t i ==
平面波通过薄膜一次A,B 两点的位相变化2111cos 2i h n θλ
π
δ=
⎪⎭
⎪
⎬⎫
⎪⎩⎪⎨⎧--=-+=2100'
21121'2112cos ]exp()exp([)exp()exp(1
i r t r t n i E i E H i E i E E θμεδδδδ
⎪⎪⎪⎪⎭⎪
⎪⎪⎪
⎬⎫
⎪⎪⎪⎪
⎩⎪⎪⎪
⎪⎨⎧=-=
+-=2
1001222'212211cos )(2)exp()(2)exp(i r t n H E i E H E i E θμεηηδηδ ⎪
⎭⎪⎬⎫⎪⎩⎪⎨⎧+-=-=)cos(sin )sin (cos 121121
112
121δδηηδδi H i E H i H E E t
写成矩阵形式：
波）
（对波）
对P cos S (cos cos 2cos sin sin cos 2
1
0012100
12
1112211111
111i i i n n h n H E i i
H E θμεηθμεηθλ
πδδδηδηδ=
==⎥⎦⎤⎢⎣⎡⎥⎥⎦⎤
⎢⎢⎣⎡
--
=⎥⎦
⎤⎢⎣⎡
当膜系包含N 层膜时，则有
⎥⎦
⎤⎢⎣⎡=⎥
⎦⎤
⎢⎣⎡++112111...N N N H E M M M H E 整个膜系的特征矩阵N M M M M
...21=
膜系反射率的计算
令M 的矩阵元为A,B,C,D 即⎥
⎦
⎤
⎢⎣⎡=D C B A M
则⎥⎦⎤
⎢⎣⎡⎥⎦⎤⎢⎣
⎡=⎥⎦⎤⎢⎣⎡++11
11N N H E D C B A H E 膜系反射系数G
G G
G i r D C B A D C B A E E r ηηηηηηηη+++--+=
=
000011
膜系反射系数
G
G D C B A t ηηηηη+++=
000
2
反射率：*⋅=r r R 投射率：*⋅=t t T
由两种不同介电常数（,A B εε）和厚度（,A B d d ）的电磁材料所组成的一维周期性光学薄膜结构。

利用传输矩阵计算光在介质薄膜。

本文中所选取的色散材料，如[11,12]，LiF （模型中的A 层）的折射率在
1.25~
2.35m μ范围内的关系为
2241.387610.0017960.000410.00230450.00000557A n L L λλ=+---
Si （模型中的B 层）的折射率在同一波段为
2243.416960.1384970.0139240.00002090.000000148B n L L λλ=++-+，
其中21/0.028L λ=-。

在本课程设计中取m μλ5.1=，则==b n 1.3831,a n 3.4770
1.4 1.6
1.82
2.2 2.4 2.6
1.3741.376
1.378
1.38
1.3821.384
1.386
波长d/μm
折射率n
介质A 折射率随波长变化关系
1.4 1.6
1.82
2.2 2.4 2.6
波长d/ m
折射率n
介质B 折射率随波长变化关系
三、 计算程序
（运行环境Matlab R2009a ） clear all;
n0=1;na=1.3831;nb=3.4770;%折射率 n=input('输入光学薄膜的周期数:'); theta1=input('输入入射角:'); da=input('输入介质A 的实际厚度:'); db=input('输入介质B 的实际厚度:'); for lambda=1250:2350;
theta2=asin(n0*sin(theta1)/na);%介质A中的折射角theta3=asin(n0*sin(theta1)/nb);%介质B中的折射角deltaa=2*pi*na*da*cos(theta2)/lambda*1e9; deltab=2*pi*nb*db*cos(theta3)/lambda*1e9;
%由于波长取值不是SI单位制，此处化为SI单位制epsilon0=1e-9/(36*pi);%真空介电常数
mu0=4*pi*1e-7;%真空磁导率
m=sqrt(mu0/epsilon0);%本征阻抗
eta0=n0*cos(theta1)/m;
etaa=na*cos(theta2)/m;
etab=nb*cos(theta3)/m;
%传输矩阵特征参数
MA1=cos(deltaa);
MA2=-1i*sin(deltaa)/etaa;
MA3=-1i*etaa*sin(deltaa);
MA4=cos(deltaa);
MA=[MA1,MA2;MA3,MA4];%介质A中传输矩阵
MB1=cos(deltab);
MB2=-1i*sin(deltab)/etab;
MB3=-1i*etab*sin(deltab);
MB4=cos(deltab);
MB=[MB1,MB2;MB3,MB4];%介质B中传输矩阵
O=MA*MB;%介质A和介质B等效层的传输矩阵
M=O^n; % n层光学薄膜总的传输矩阵
A=M(1,1);
B=M(1,2);
C=M(2,1);
D=M(2,2);
r=(A*eta0+B*eta0^2-C-D*eta0)/(A*eta0+B*eta0^2+C+D*eta0) ;%反射系数
R=r*r';%反射率
t=2*eta0/(A*eta0+B*eta0^2+C+D*eta0); %透射率
T=t*t’;%透射率
s=lambda-1250+1;%对选取的数据点进行编号，以方便下一步存储
k(1,s)=R;%建立存储数据的向量，以方便下一步作图
end
lambda=1250:2350;
plot(lambda,k);
xlabel('\lambda/nm');
ylabel('R');
四、设计结果
1、同周期相同结构光学薄膜无吸收时反射率R和透射率T的关
系（T+R=1）
λ/nm
（n=8 ,da=m db m μμ2.0,8.0=）
2、同结构不同周期一维光学薄膜的反射率
λ/nm
R
R
λ/nm
R
λ/nm
R
λ/nm
λ/nm
R
3、同周期（n=8）不同结构单元对一维光学薄膜的反射率的影响
（入射角0
=θ）
λ/nm
R
（图1：da=0.5m μ，db=0.35m μ）
λ/nm
R
（图2：da=0.5m μ，db=0.36m μ）
1200
1400
1600
18002000
2200
2400
00.10.20.30.40.50.60.70.80.91λ/nm
R
（图3：da=0.5m μ，db=0.37m μ）
4、同周期（n=8）不同入射角对一维光学薄膜的反射率的影响
（da=m db m μμ2.0,8.0=）
λ/nm
R
)
6
(0πθ=
1200
1400
1600
18002000
2200
2400
00.10.20.30.40.50.60.70.80.91λ/nm
R
)
4
(0πθ=
1200
1400
1600
1800200022002400
00.10.20.30.40.50.60.70.80.9
1λ/nm
R
)
3(0πθ=
五、
六、 参考文献
Matlab 实用教程 马莉 清华大学出版社 物理光学（第3版） 梁铨廷 电子工业出版社
Matlab 环境下多层光学薄膜的数值计算 广西师范学院学报
传输矩阵法分析一维光子晶体的传光特性。