椭偏仪Psi/Delta计算全解析:Fresnel公式与Python实现
2026/9/8 1:26:28 网站建设 项目流程

简介:椭偏测量技术作为一种非接触式光学手段,可用于分析材料表面和薄膜性质,在半导体检测、材料表征等领域应用广泛。基于Python的椭偏仪计算项目,面向材料科学、半导体制造、生物医学等方向,目标是从测量数据中计算关键参数ψ和δ,进而反演薄膜厚度、折射率、消光系数等光学常数。压缩包共包含4个文件,总大小462KB,其中Python脚本负责读取原始数据并执行数值计算;PDF文档讲解椭偏仪原理及代码用法;Markdown文档提供项目说明;另有配置文件便于管理依赖与环境。目前已有803人学习,许多初学者借助该资源入门了椭偏仪数据处理,说明其讲解与代码具有不错的实用性。通过研读源码和文档,使用者可以理解ψ和δ的计算流程,涉及复数运算、最小二乘拟合及反演算法等核心内容;同时也能掌握Python处理光学测量数据的完整思路,适用于科研数据分析、教学演示或进一步定制开发,有助于提升研究效率。 提到椭偏仪,很多刚接触薄膜测量的人第一反应是“这仪器能测膜厚”,但真把数据拿到手,看到一对叫psi和delta的参数,往往就卡住了。这两个值并不是直接可读的厚度,也不是折射率,而是反映了光在样品表面反射前后偏振状态的变化量。我们做光学薄膜表征、半导体工艺监控、甚至电化学原位分析的时候,天天跟这两个量打交道,但很多人只是把椭偏仪当成“黑盒子”,点一下测量就等结果出来,至于psi和delta到底是怎么算出来的,反演时为什么会出现多解,出了问题也不知道从哪里排查。

这篇文章我想把psi和delta的计算这件事彻底说透。从一个简单的三层膜模型出发,写清楚Fresnel公式怎么用、椭偏方程怎么来的、多层膜怎么递推、计算程序该怎么组织,再结合我自己实际做过的材料体系,把那些“公式里看不出来”的坑也一并列出来。无论你是要用椭偏仪测SiO2膜厚,还是做有机薄膜或者金属氧化物的光学常数拟合,这篇文章都能给你一个从头开始的可复现思路。

1. 椭偏测量到底在算什么:psi和delta的物理意义

1.1 一次反射就把光的偏振状态拆开了

椭偏测量的物理基础其实很简单:光是一种横波,电场矢量可以在垂直于传播方向的平面内任意取向。我们通常把电场在入射面内的分量称为P偏振(p-polarized),垂直于入射面的分量称为S偏振(s-polarized)。当一束线偏振光以一定角度斜入射到样品表面时,P分量和S分量的反射系数并不相同,反射后两个分量的振幅比会发生变化,同时两者之间的相位差也会被样品“调制”掉一部分。

椭偏仪做的事情,就是精确测量这种变化。它不直接测量反射光的绝对强度,而是测量反射前后P光和S光的振幅比以及相位差的变化,最终输出两个参数:

  • psi(Ψ):反射后P/S振幅比的反正切值,代表振幅的变化。
  • delta(Δ):反射后P光和S光之间的相位差,代表相位的变化。

这两个参数合在一起,就构成了椭偏测量最原始的数据单元。后续无论是算膜厚、算折射率、还是拟合多层膜结构,都建立在psi和delta的基础之上。

1.2 为什么必须同时看重和相位

有些人会问:我直接测反射率不也能判断膜厚吗?为什么非要搞这么复杂的偏振测量?这里的关键在于灵敏度。反射率测量实际上只取了光的“振幅平方”信息,相位信息全部丢掉了。而椭偏测量把振幅比和相位差同时保留下来,信息量翻了一倍,对亚纳米级别的薄膜变化尤其敏感。

举一个非常直观的例子:一束光垂直入射到样品上时,P和S是完全对称的,椭偏测量完全失效。但在Brewster角附近,P光的反射系数会经过一个从正到负的过零变化,相位会发生剧烈跳变,此时delta对膜厚和光学常数的变化极其敏感。这是我实际测量SiO2/Si样品时体会最深的一点:膜厚只有1到2纳米的自然氧化层,反射率差异几乎测不出来,但delta的变化可以达到好几度。这也是椭偏技术能够胜任原子级薄膜表征的根本原因。

2. 核心计算公式推导:从Fresnel公式到椭偏方程

2.1 三层膜模型和Fresnel反射系数

实际椭偏测量的样品结构千变万化,但最基本的模型是“空气/薄膜/基底”这样的三层结构。设入射介质为空气(n0=1),薄膜复折射率为n1,厚度为d,基底复折射率为n2。光以入射角θ0射入,经折射后在薄膜和基底界面处的折射角分别为θ1和θ2,由Snell定律确定。

在每一个界面上,P光和S光的Fresnel反射系数可以分别写出来。以薄膜-基底界面为例,P光和S光的反射系数分别为:

r_p12 = (n2·cosθ1 - n1·cosθ2) / (n2·cosθ1 + n1·cosθ2) r_s12 = (n1·cosθ1 - n2·cosθ2) / (n1·cosθ1 + n2·cosθ2)

空气-薄膜界面的r_p01和r_s01也以同样的形式写出,只是把下标换成对应的介质参数。

2.2 标准椭偏方程:rho = rp / rs = tan(psi)·exp(i·delta)

有了各层的折射角和界面反射系数,还需要考虑薄膜厚度带来的多次反射干涉效应。光在薄膜内部来回反射,最终从表面出射的总反射系数是无穷多项的叠加,这正好是等比级数求和,可以写成封闭形式:

Rp = (r_p01 + r_p12·exp(-i·2β)) / (1 + r_p01·r_p12·exp(-i·2β)) Rs = (r_s01 + r_s12·exp(-i·2β)) / (1 + r_s01·r_s12·exp(-i·2β))

其中β是薄膜的相位厚度,β = 2π·d·n1·cosθ1 / λ,λ是入射光在真空中的波长。

椭偏参数psi和delta就定义在这两个总反射系数的比值上:

ρ = Rp / Rs = tan(Ψ)·exp(i·Δ)

所以只要给定入射角、波长、各层复折射率和膜厚,psi和delta就可以精确计算出来。反过来说,实验中测得了psi和delta,就可以通过数值反演去反推膜厚和折射率。

2.3 多层膜情况下的递推实现

三层膜只是最简单的情况。实际样品经常是多层的,比如“空气/减反膜/硅氧化物/硅基底”,或者更复杂的膜系。此时每一层都可以用其特征矩阵来描述,或者用递推法从底层往上逐层合成等效反射系数。业界常用的算法是Abeles矩阵法,它的基本思路是把每一层的电磁场传播效果写成一个2x2矩阵,多层膜按顺序相乘,最终从矩阵元素直接读出总反射系数。

用Python实现的时候,我会直接对每一层计算相位厚度,从最底层的界面开始往上递推。这种方式不仅对三层模型有效,扩展任意层数也只需要增加一个循环,代码非常紧凑。后面的实操部分我会给出完整的示例。

3. 实操:从零写一个psi和delta计算程序

3.1 材料的光学常数从哪来

计算的第一步是拿到材料在各个波长下的光学常数n和k。有几个渠道:最常用的是使用商业椭偏仪自带的材料数据库,里面有硅、SiO2、Si3N4、各种金属的标准数据;如果你没有仪器配套的数据库,也可以用公开的折射率数据网站,比如refractiveindex.info,上面有海量材料的实测数据可以下载;再不然,就要自己用椭偏仪先对一个已知厚度的样品做反演拟合,把光学常数建立起来,这其实是很多新材料研究的日常工作。

关于数据,我要多提醒一句:同一个材料在不同制备条件下,光学常数可能有显著差异。比如PECVD生长的SiO2和热氧化生长的SiO2,虽然都叫二氧化硅,但密度、成分、缺陷态都有差别,折射率在可见光区可能差0.005左右,对精密拟合来说这个差异不可忽略。所以如果有条件,最好先用自己工艺条件下的样品建立数据库,而不是过分依赖“标准值”。

3.2 Python实现:三层模型的psi和delta计算

下面这段代码是我实际计算时常用的一个简化版本。输入是入射角、波长、薄膜和基底的复折射率以及膜厚,输出是psi和delta。代码只保留最核心的数学逻辑,方便你根据自己的需求扩展。

import numpy as np def snell_theta(n_i, n_t, theta_i): """根据Snell定律计算透射角(复角)""" ratio = n_i / n_t * np.sin(theta_i) # 超过临界角时ratio大于1,取arcsin的复数分支,保证物理一致性 return np.arcsin(ratio) def fresnel_amplitude(n0, n1, theta0, theta1): """P光和S光在单一界面的Fresnel反射系数""" cos0 = np.cos(theta0) cos1 = np.cos(theta1) r_p = (n1 * cos0 - n0 * cos1) / (n1 * cos0 + n0 * cos1) r_s = (n0 * cos0 - n1 * cos1) / (n0 * cos0 + n1 * cos1) return r_p, r_s def psi_delta_from_model(theta_deg, lam_nm, n_film, k_film, d_nm, n_sub, k_sub): """ 计算空气/薄膜/基底三层模型的psi和delta 参数: theta_deg: 入射角(度) lam_nm: 波长(纳米) n_film, k_film: 薄膜折射率与消光系数 d_nm: 薄膜厚度(纳米) n_sub, k_sub: 基底折射率与消光系数 返回: psi(度), delta(度) """ lam = lam_nm * 1e-9 d = d_nm * 1e-9 theta0 = np.deg2rad(theta_deg) n0 = 1.0 + 0j n1 = n_film + 1j * k_film n2 = n_sub + 1j * k_sub theta1 = snell_theta(n0, n1, theta0) theta2 = snell_theta(n1, n2, theta1) r_p01, r_s01 = fresnel_amplitude(n0, n1, theta0, theta1) r_p12, r_s12 = fresnel_amplitude(n1, n2, theta1, theta2) beta = 2 * np.pi * d * n1 * np.cos(theta1) / lam phase = np.exp(-2j * beta) r_p = (r_p01 + r_p12 * phase) / (1 + r_p01 * r_p12 * phase) r_s = (r_s01 + r_s12 * phase) / (1 + r_s01 * r_s12 * phase) rho = r_p / r_s psi = np.degrees(np.arctan(np.abs(rho))) delta = np.degrees(np.angle(rho)) return psi, delta # 示例:硅基底上100nm的SiO2,波长632.8nm,入射角70度 psi, delta = psi_delta_from_model( theta_deg=70, lam_nm=632.8, n_film=1.457, k_film=0.0, d_nm=100, n_sub=3.875, k_sub=0.016 ) print(f"psi = {psi:.3f} deg, delta = {delta:.3f} deg")

这段代码的写法有几个值得注意的地方。首先,在计算透射角时我没有把折射率强制限定为实数,因为在金属基底或吸收较强的材料中,折射率本身就是复数,透射角和cos值自然也是复数,用np.arcsin处理复数分支即可。其次,Fresnel系数公式里我使用了“n·cosθ”形式的组合方式,这种形式对复数折射率的处理更加鲁棒,不容易出现符号翻转的问题。

3.3 做一个厚度扫描算例

算了单点还不过瘾,我们把厚度从10nm扫到500nm,看看psi和delta会怎么变化。下面是扫描代码和结果分析:

import matplotlib.pyplot as plt d_list = np.linspace(10, 500, 200) psi_list = [] delta_list = [] for d in d_list: psi, delta = psi_delta_from_model( theta_deg=70, lam_nm=632.8, n_film=1.457, k_film=0.0, d_nm=d, n_sub=3.875, k_sub=0.016 ) psi_list.append(psi) delta_list.append(delta) fig, ax1 = plt.subplots() ax1.plot(d_list, psi_list, 'b', label='psi') ax1.set_xlabel('Thickness (nm)') ax1.set_ylabel('Psi (deg)') ax2 = ax1.twinx() ax2.plot(d_list, delta_list, 'r', label='delta') ax2.set_ylabel('Delta (deg)') plt.show()

从曲线可以看到,psi和delta随厚度的变化并不是单调的,而是呈现周期性振荡,振荡周期和薄膜的相位厚度β直接相关。这个特点带来一个很重要的推论:单波长椭偏测量存在“厚度周期模糊”问题,同一个psi和delta可能对应多个厚度值,相差一个周期。这也是为什么实际工作中我们会选择多波长椭偏仪,或者用多入射角测量来打破周期简并。如果你手里只有单波长仪器,那在拟合时务必给厚度一个合理的初始估计范围,否则非常容易陷入错误的局部最优解。

3.4 从psi和delta反推厚度:参数拟合的必要性

正向计算psi和delta是“由模型到数据”,反推厚度就是“由数据到模型”。实际拟合时,我们会定义一个评价函数,通常是实验测得的psi和delta与模型计算的psi和delta之间的均方误差,然后用优化算法不断调整膜厚和光学常数,让这个误差最小。

最简单的做法是scipy.optimize.curve_fit或者least_squares。这里给一个极简的反演示例,把之前的正向函数包一层就完成了:

from scipy.optimize import least_squares def model_residual(x, theta_deg, lam_nm, psi_meas, delta_meas): d = x[0] psi_calc, delta_calc = psi_delta_from_model( theta_deg=theta_deg, lam_nm=lam_nm, n_film=1.457, k_film=0.0, d_nm=d, n_sub=3.875, k_sub=0.016 ) return [psi_calc - psi_meas, delta_calc - delta_meas] # 假设实验测量得到一组psi和delta psi_meas, delta_meas = psi, delta result = least_squares( model_residual, x0=[50], args=(70, 632.8, psi_meas, delta_meas) ) print(f"拟合膜厚 = {result.x[0]:.2f} nm")

需要说明的是,这个反演只优化了厚度这一个变量。如果还要同时拟合折射率,变量数量会增多,这时就要特别小心参数之间的相关性。比如厚度和折射率在膜层很薄的时候会高度耦合,实测中你会发现厚度变大一点、折射率变小一点,拟合结果也能保持几乎同样好的误差。这个现象被称为“相关性问题”,在极薄膜(<10nm)的情况下尤其严重。破解的办法是:要么固定折射率只拟合厚度,要么用多波长数据一起拟合来约束参数,要么用已知厚度的样品先单独定出折射率,再反过来拟合其他样品。

4. 常见计算误差与排查技巧实录

4.1 典型问题对照表

我在实际使用中踩过不少坑,很多都让新手防不胜防。我把最常见的几类问题整理成了一张表,你在计算或拟合时如果发现结果不对劲,可以按图索骥:

现象可能原因排查方法
拟合误差很大,psi和delta系统性偏移入射角标定不准确用标准SiO2/Si样品校验角度,做角度偏移拟合
厚度反演变出奇异值模型层数不足,未考虑表面粗糙层或界面层增加粗糙层(EMA层)再拟合
膜厚结果周期性跳跃单波长数据简并改用多波长或增加入射角数据进行联合拟合
delta在180度附近剧烈变动相位包裹(wrapping)数据处理时做相位展开处理
同一组数据不同初始值得到不同结果存在多个局部最优解用全局优化算法(如差分进化)或手动网格扫描确定初值
波长变化时k值为负材料模型选择不当检查光学常数是否满足Kramers-Kronig一致性,更换振子模型

4.2 相位包裹问题

delta的取值范围在数学上通常定义在[-180°, 180°]或者[0°, 360°],但实际膜厚变化会让delta连续跨越边界,导致曲线看起来像“跳变”了。这个现象在厚度扫描时特别常见,前面那个算例如果厚度再大一些,你也会看到delta在180度附近“折断”。

处理办法有两种。一种是在反演时不做特殊处理,直接用角度制下的数值计算误差,但误差函数在边界附近会不平滑,优化算法容易卡住。另一种是先把测量数据和模型计算结果的delta都转换到连续相位,再做差值。更稳妥的做法是在计算评价函数时使用复数rho而不是直接用delta角度,因为rho是连续的,不存在角度包裹问题。

# 推荐:用复数rho的误差代替角度误差 def model_residual_rho(x, theta_deg, lam_nm, rho_meas): d = x[0] psi_calc, delta_calc = psi_delta_from_model(...) rho_calc = np.tan(np.deg2rad(psi_calc)) * np.exp(1j * np.deg2rad(delta_calc)) return [rho_calc.real - rho_meas.real, rho_calc.imag - rho_meas.imag]

我一直习惯在最终的反演代码里用复数rho作为评价对象。虽然初学时觉得直接用psi和delta更直观,但处理真实数据时就会发现,复数域拟合的稳定性和对边界情况的鲁棒性要好得多。

4.3 光学常数模型选择

在拟合膜厚的同时,很多人还想得到材料的光学常数,这时就不能只把n、k当作每个波长独立的变量来拟合了,那样自由度太多,必然导致过拟合。业界通行的做法是用物理模型来描述光学常数的色散关系。常用的有:

  • Cauchy模型:适用于透明薄膜或弱吸收薄膜,形式为n = A + B/λ² + C/λ⁴,参数少,拟合快。
  • Lorentz振子模型:适用于半导体和绝缘体在吸收带附近的色散,能同时给出n和k,满足Kramers-Kronig关系。
  • Drude模型:适用于金属中自由电子的响应,常用于金、银、铜等金属薄膜的光学常数建模。
  • Tauc-Lorentz模型:适合非晶半导体薄膜,在带隙附近表现良好。

选模型的原则是“够用就好”,参数越少越好。比如测量一个透明的PMMA薄膜,用Cauchy模型三个参数就能完美描述整个可见光范围内的色散,完全没必要上Lorentz振子。而测量一层非晶硅薄膜,用Tauc-Lorentz模型比用多个Lorentz振子系统更合理,因为非晶硅的带隙吸收形状和晶态硅不同,用振子模型硬拟合会得到不满足物理意义的参数组合。

4.4 膜层结构模型和初始值设置

最后还有一个容易被忽略但极其重要的点:拟合前一定要先把膜层结构模型定义清楚。椭偏拟合的本质是“先假设结构,再优化参数”,结构假设对了,拟合才会得到有意义的结果;结构假设错了,哪怕拟合误差很小,结果也是错的。典型的例子是:在硅基底上测一层很薄的氧化层,如果不加一个微小的表面粗糙层,拟合出来的氧化层厚度往往比真实值稍大一些,因为表面粗糙度的影响被“吸收”进了厚度里。我们的经验是,对于热氧化SiO2/Si这种界面非常平整的体系,三层模型就够用了;对于薄膜太阳能电池或者有机半导体器件,界面粗糙度不可忽略,通常会在模型中增加一个Bruggeman有效介质层来代表粗糙界面。

初始值方面,我的建议是不要完全依赖自动拟合。可以先做一个膜厚扫描的“伪图谱”:固定光学常数,在合理的厚度范围内扫描,看预测的psi和delta曲线与实验点的距离分布,大致能找到一个大致的厚度区间。把这个区间内的值作为优化初值,既快又稳,远比随便给一个初值硬拟合可靠。

结尾的一点经验分享

做椭偏数据分析和计算这些年,我自己最大的体会是:椭偏仪是个测量仪器,但psi和delta的计算本质上是一个光学建模问题。仪器只是帮你精确地测到了偏振变化,而真正决定结果可靠性的,是你对样品结构的理解和对模型的选择。刚开始接触椭圆偏振光谱的时候,我也曾试图把所有已知信息一股脑塞进模型,结果参数过多、拟合崩溃;后来学乖了,先做最简单假设,逐步增加复杂度,每一步都检查拟合误差是不是真的在减小。

如果你正要开始写自己的椭偏计算程序,我建议从三层模型开始,用已知样品验证代码正确性,再逐步扩展到多层膜、各向异性材料或者原位实时监测。另外一个小技巧:在编程实现时永远保留把模型计算值和实验数据画在同一张图上的习惯,肉眼观察比任何统计指标都更能快速发现问题。这也是为什么前面示例代码里要画厚度扫描曲线,数据一跑出来,周期性、边界跳变、灵敏度区间都一目了然,接下来怎么改模型、怎么设初值,心里就有底了。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询