相干光成像的精确计算,说到底是波动光学里最绕不开、也最容易翻车的一块硬骨头。激光刚发明那会儿,全息术火了一阵,可算力跟不上,很多人只能靠解析公式做特殊场景的推导,碰到复杂光场就歇菜。我最早系统接触这套东西,是在一位长期深耕相干光计算成像领域的前辈所著的教材里,配合手推菲涅耳衍射积分和离散傅里叶变换,才慢慢把“数值计算”和“物理成像”这两条线扭到一起。做光学仿真和数字全息的人应该都有体会:公式写在纸上很漂亮,一旦写成程序,采样点数、传播距离、补零方式稍有不对,重建图像全是网格状伪影,连对焦都对不出来。这篇文章就是把相干光成像从物理原理到离散计算再到实际应用这条路完整过一遍,讲清楚每一步为什么这么做、哪些地方是隐藏的坑、真实项目中我是怎么取舍的。无论你是刚接触计算光学的研究生,还是已经在做数字全息、定量相位成像、衍射计算的工程师,这篇都能帮你把“算得出来”升级成“算得准”。
1. 为什么普通光学的“老一套”在相干光面前失灵:相位才是核心变量
1.1 相干成像与普通成像的底层差异
传统光学成像系统,比如单反镜头加可见光传感器,处理的基本是“光强分布”。对非相干或弱相干光源,像面上几何光学那套高斯成像公式就能用——物体上的每个点独立成像,点与点之间不干涉,最后叠加起来的就是肉眼看到的画面。这种模型只关心振幅平方,也就是强度,相位信息在成像过程中被彻底丢掉,而且丢了也就丢了,不影响你看清东西。
但相干光成像不一样。激光照射物体之后,物体表面的反射光或透射光携带完整的复振幅分布:既有振幅,又有相位。光场传播过程中,相位的变化恰恰决定了精细结构的走向。全息术能记录三维信息、干涉测量能达到亚波长精度、相干衍射成像能突破传统分辨率极限,全部依赖对相位分布的精确恢复和控制。
举个生活里的例子:你拍一张平面照片,能看出一个人长什么样,但看不出这个人离镜头多远,除非有阴影或透视关系。然而一张相位分布图,直接记录了光波走过的“路程差异”。一块透明的相衬样品,厚度差只有几百纳米,强度图几乎全透明、什么都看不见,相位图却能清楚显示每一个细胞边缘。这就是为什么相干光成像中,相位不是可有可无的“配角”,它才是真正的主角。
1.2 记录介质只认强度:相位信息的丢失与恢复难题
物理上,人眼和探测器都只能响应光强(振幅平方),相位本身不能被直接“拍下来”。光波频率太高,电磁场的振动周期是飞秒量级,没有任何传感器能实时跟踪场矢量。我们看到的强度估算式是I = |U|²,这个平方运算把相位全部消掉了。
那要怎么恢复相位?两条路。一条是干涉法,让携带信息的物光和一束已知的参考光叠加,把相位差转变成强度上的干涉条纹,这就是全息记录的原理。另一条是迭代计算法,在已知强度分布和边界条件的前提下,用约束条件反复迭代求解相位,这就是相干衍射成像(CDI)和分片相位恢复(PIE)的路子。
但要注意,干涉法获得的“全息图”依然是强度分布,只是这个强度被参考光调制成光栅结构,相位信息“藏”在条纹的分布里,想还原出来还是得经过一次数值计算。所以无论是干涉还是非干涉,相干光成像最终都逃不过一个问题:如何从离散采样的强度数据里精确重建复振幅。这一关过不了,后面的应用全是空中楼阁。
2. 把衍射积分“算准”的关键:采样、边界与传播距离的三方博弈
2.1 从基尔霍夫到菲涅耳:工程计算用哪一版公式
光场从物平面传播到像平面,最严格的是求解麦克斯韦方程组得到的基尔霍夫衍射积分。理论上它描述了任意距离下的衍射场,但积分形式复杂,数值计算代价很高。工程上绝大多数场景都在近轴近似下做简化,最常用的就是标量衍射框架下的菲涅耳衍射积分:
U(x,y,z) = (exp(ikz) / (iλz)) ∫∫ U₀(ξ,η) exp[ i(k/2z) ((x-ξ)² + (y-η)²) ] dξ dη
这个公式可以理解为:光在自由空间传播距离z之后,物面上每个点发出的球面波在观察面上产生的叠加。这里关键的一个近似是 z 远大于孔径尺寸和观察区域尺寸,即菲涅耳近似,它把球面波的准确径向距离近似成抛物面展开,角度误差对相位的影响在可接受范围内。
如果传播距离更远,观察面只关心远场分布,菲涅耳积分可以再退化为夫琅禾费衍射。这时候积分内部的二次相位因子消失,计算变成了单纯的傅里叶变换。但工程上要注意:很多“远场测量”实验,实际距离未必满足严格夫琅禾费条件,用过快的傅里叶变换近似,结果会跟实验对不上。我个人的习惯是,除明确测量远场光斑的情况外,一律用菲涅耳积分做离散版本,多算一点开销,换来鲁棒性。
2.2 离散采样的硬约束:奈奎斯特条件在光场计算中的变形
在计算机里做衍射计算,第一步就是把连续的物平面光场采样成N×N的网格。采样间隔的选择是第一个大坑。参考信号处理里的奈奎斯特采样定理,光场离散化同样要求采样率大于最高空间频率的两倍。
关键问题是:这个“最高空间频率”哪儿来的?两个来源。一是物体本身的精细结构,最小周期记为Δx_min,采样间隔必须满足 Δx < Δx_min/2,否则细节直接混叠;二是菲涅耳积分里的二次相位因子 exp(iπ(x-ξ)²/(λz)) 引入了相位随位置快速变化的分量。当传播距离z很小、横向范围很大时,这个相位变化非常剧烈,要求采样间隔足够小,直观上讲就是 Δξ < λz / L,其中L是观察面或物面的横向尺寸。
所以采样间隔实际上被两方夹住:一方是物体的空间频率,另一方是传播距离和横向范围共同决定的空间频谱带宽。很多初学者只照顾了物体边界,忽略了二次相位因子的变化,算出来的光场全是卷绕伪影。这类问题在计算全息做“近距离投影”时尤其严重——距离只有几毫米,像素尺寸又大,不满足约束时无论如何加窗函数都救不回来。
提示:判断采样是否合理的快速方法是计算相位差分。把衍射积分里二次相位部分离散后,相邻采样点的相位差绝对值如果超过π,必然出现卷绕误差;小于π也只是勉强合格,工程上建议留三到五倍余量。
2.3 边界效应与补零:一个常被忽略的误差来源
离散傅里叶变换有一个周期性假设:时域是周期延拓的,空域也一样。做衍射计算时,物面和观察面的有限尺寸会被当成周期重复的无限阵列来处理,于是边缘部分会产生相互干扰,这叫“环状卷绕效应”。表现为重建像的四个边出现另外一边的亮度“泄漏”,严重的时候图像正中央会叠上复制出来的鬼影。
怎么处理?最直接的手段是补零。把输入光场矩阵周围的空白区域用零填充,扩大有效计算窗口,让边界的相互影响减弱。实际操作中我一般至少补两倍尺寸,也就是把512×512的物面补成2048×2048再算衍射。补得越大,边界效应越小,但计算量按平方上升,四倍补零再叠加多次傅里叶变换,耗时很容易翻十几倍,必须平衡。
还有一个不太直观的边界效应来自光场本身的相位倾斜。如果入射光不是垂直照明,而是带有一个横向波矢分量,光场整体有一个线性相位梯度。离散采样时这个梯度也必须满足奈奎斯特条件,否则倾斜部分直接变成混叠,最终重建像会移动甚至扭曲。处理办法是去掉线性相位因子,先算垂直入射的分布,最后再把平移量换算回来。特别是做离轴全息重建时,这个操作几乎是必须的。
3. 数值衍射计算的三种主流路线与我的工程选型心得
3.1 SR法(单次傅里叶变换):适合远距离传播的快速方案
把菲涅耳衍射积分离散化后,可以直接改写成单个傅里叶变换形式,这就是SR法(Single Fourier Transform)。它的操作是:输入光场乘以二次相位因子,做一次FFT,再乘以另一个二次相位因子和系数,就得到传播z距离后的光场。这个算法最大的优点是快,只做一次FFT,特别适合实时性要求高的场景,比如全息视频的快速重建。
但SR法的局限性也很明显。它的采样间隔必须满足条件:输出面采样间隔 Δu = λz / (N Δx),直接由输入面的尺寸N Δx决定。这相当于说,输入窗有多宽,输出面的采样分辨率就被钉死了。当传播距离很短、物体窗口又很大时,Δu会非常小,输出面的视野范围 N·Δu 变得极其有限,只看得到中心一小块区域,周边信息全部丢失。反过来,传播距离很长时输出视野太大、分辨率降低,物体细节又糊在一起。
所以SR法最适合的场景是:物面尺寸和观察面尺寸差别不大、传播距离适中或偏远、只需要观察中心区域的分布。举个例子,研究光学涡旋束从出射端到焦平面的演变,或者计算经过一段自由空间后的远场光斑,SR法都很好用。
3.2 卷积法(IR法):聚焦到场采样细分时的可靠选择
另一种思路是把菲涅耳衍射看成卷积运算:输出是输入与固定脉冲函数 h(z; x) 的卷积。根据卷积定理,时域做卷积等价于频域做乘积,于是可以先把输入和脉冲函数分别做FFT,相乘后做逆FFT得到输出。这就是IR法(Impulse Response Method),也叫卷积法。
跟SR法比,卷积法最大的优势是输出面的采样间隔可以独立选择,不再受输入窗口尺寸约束。你可以设定输出面采样间隔等于输入面采样间隔,保持像素对应关系;也可以根据实际需要重新定义观察区域。这意味着,即使物面和观察面都很小,只要补零补得够大,也能得到合理的结果。
但代价是计算量。卷积法至少要做三次FFT,只算一次结果时比SR法慢很多。而且卷积法的隐含条件是光场在采样区域内必须是循环卷积,同样要处理边界周期性问题,补零的需求比SR法更强烈。实际操作中我通常在两种场景下用卷积法:一是模拟通过光学透镜聚焦后的焦平面场分布,因为此时输出面边界非常窄、需要精细采样;二是在做三维体光场的多层传播时,确保每层之间像素格点对齐,方便后续三维重建算法。
3.3 三种方法的效果对比与适用范围
实际工程中不只有SR和卷积两个极端。在菲涅耳近似的框架下,还可以按传播距离分区域选择算法。比如角谱法(Angular Spectrum Method),本质上也属于卷积法的一种,区别在于它在频域直接乘上传播算符而不是脉冲函数,对短距离传播特别稳定。
为了直观说明差异,我整理了一张对比表,列出不同传播距离下应该优先选择哪种方法:
| 传播距离 | 推荐方法 | 核心原因 | 主要风险 |
|---|---|---|---|
| 毫米级(近场) | 角谱法 / 卷积法 | 需要在极近距离保持采样精度 | 二次相位变化快,必须极密采样,计算量大 |
| 厘米级到米级(菲涅耳区) | SR法 / 卷积法均可 | SR法快,卷积法灵活 | SR法输出视野受限;卷积法需注意补零 |
| 远场(夫琅禾费区) | SR法(一次FFT) | 二次相位消失,直接对应傅里叶变换 | 需确认距离满足远场条件 |
| 多层传播(三维体场) | 卷积法 | 保持层间网格对齐 | 每层都要补零,内存消耗大 |
注意,这只是经验性的第一判断。真实项目里还需要结合物体本身的频谱分布,如果物体含有超高空间频率成分,即使传播距离“理论上”适合SR法,也仍然可能混叠。保守做法是先算一次只包含物体中心小区域的快速试算,对比几种方法的结果差异,差异超过5%就说明偏离近似条件了。
3.4 我踩过的坑:采样间隔算错导致重建像“虚胖”
说起来,我在做数字全息重建时踩过一个非常经典的坑。当时用CCD记录了一张离轴全息图,物面大概512×512像素,传感器感光面尺寸是6.4mm×4.8mm,像素尺寸7.4μm。按SR法重建时,我还没转换单位,直接把传播距离设成30cm,用像素数代入公式,得到的结果图像是花掉的、完全对不上焦。一开始我怀疑是参考光角度没滤干净,去零级像处理也做了,问题依旧。
后来我把物面尺寸从像素单位换算成毫米单位再做一次,重建图像立刻就清晰了。问题就出在FFT的结果索引跟物理尺寸量纲之间的映射上。SR法公式里,FFT输出阵列对应的空间采样间隔是 λz/(N Δx),这个 Δx必然是物理长度单位,比如米或毫米。如果拿像素数当长度代入,整个坐标全部错位,重建像会“虚胖”——高频细节集中在中心、低频外移,看起来就像原地放大了一团糊影。
这个错误虽然低级,但非常容易漏。因为很多示例代码直接拿像素数算,看起来跑通了,实际上是对某个特殊距离下像素数与毫米数相等的假象所误导。建议所有人在写衍射计算程序时,从一开始就把单位统一成毫米或米,坐标网格用meshgrid生成物理坐标,不要用整数索引直接当空间坐标。这个习惯能省掉大量排查时间。
4. 相位恢复:探测器只有强度图像时,如何把相位“找回来”
4.1 为什么相位丢失会让成像“缺了一半信息”
上一节讲的衍射计算,前提是输入面的复振幅已知——振幅和相位都摆在桌上。但真实实验不是这样。用CCD或CMOS记录图像,得到的永远是强度分布,相位信息在记录瞬间就被丢掉了。
丢失相位的问题有多严重?一个标准透明的样品,比如生物细胞涂片,在透射光照明下几乎所有位置的光强都差不多,普通显微镜拍出来是透明的、看不清结构的。但相位分布里记录着细胞内部各区域的折射率差异和厚度变化,这才是区分细胞结构的核心信息。如果把相位扔掉,等于把一幅图像的纹理信息全擦掉,只留下一张纯色底图。
这里的核心难点是:能否仅从一张或多张强度图反推出相位?数学上这是相位恢复问题。因为强度图的数量通常小于未知数(每个像素一个振幅、一个相位),从单张强度图严格求解是高度欠定问题,很难一步到位。工程上往往采用迭代方法,利用已知约束逐步逼近。
4.2 迭代相位恢复的基本思路:先猜一个、反复修正
最经典的迭代方法是GS算法(Gerchberg-Saxton algorithm,后又有Fienup改进),最初用于电子显微成像。基本原理可以用一句话概括:在物平面和频域平面之间来回迭代,每一次在已知约束的位置强制满足约束条件,其余位置保持自由,经过几十上百次循环,相位分布会逐步收敛到某个合理解。
具体流程大概是这样的:已知观察面强度|g|²和照明光场振幅|U₀|,先给物面相位一个随机初猜;接着做正向衍射计算得到观察面的复振幅,把振幅替换为实测强度,保留计算得到的相位;再做逆向衍射计算回到物面,把物面振幅替换为已知的照明振幅,保留计算得到的相位作为新的估计。重复这个过程,直到观察面的更新强度与实测差异足够小。
听起来简单,但GS算法有两个致命弱点:一是对初猜敏感,随机初猜可能导致陷入局部极值,收敛到完全错误的相位;二是存在“孪生像”问题,即真假解互为镜像,迭代只能收敛到其中一个,而两个都可能满足强度约束。实际使用中需要仔细设计约束条件和光场的多样性,比如通过平移照明光斑或改变入射角度制造多个独立强度测量。
4.3 从GS算法到PIE:工程上更实用的改进思路
GS算法的收敛问题催生了结构化的相位恢复方法,最有代表性的就是PIE(Ptychographic Iterative Engine)。PIE的思路是:不用固定光斑照明整个样品,而是让一个已知的照明光斑(比如汇聚激光束产生的微型光斑)在样品表面扫描,每扫一个位置记录一张衍射图。位置之间必须有重叠区域,比如文档通常建议相邻光斑平移步长不超过光斑直径的40%到60%。
每一张衍射图、每一个扫描位置,都是一个独立的强度约束。多个约束交叉作用,大幅减少了迭代解空间,相位恢复的收敛速度和稳定性远超GS单图迭代。这也是为什么近些年X射线显微镜、电子显微镜领域,PIE类方法成了绝对主流。
从工程实现角度,PIE的核心参数有三个:光斑重叠率、初始相位猜测的位置数量、迭代轮数。重叠率太低,扫描位置之间没有交叉约束,等效于一堆独立的GS;重叠率太高,数据量增大但信息增益趋缓,计算成本浪费。我习惯先用25%到35%重叠率做快速原理验证,确认光强数据质量没问题后再提高到50%左右增加鲁棒性。
提示:相位恢复算法对噪声非常敏感。如果实验强度图周围有大量暗电流噪声或读出噪声,建议先做背景扣除和预处理,否则迭代时这些噪声会被当成真实强度约束,导致最终相位分布出现横向扫描条纹状伪影。这比单纯提高迭代次数有效得多。
5. 精确计算的实际落地:数字全息的数值重建与自聚焦
5.1 全息图到重建像:一条完整的计算链路
数字全息的流程是:激光照明物体后,物光和参考光干涉,CCD记录干涉条纹全息图。关键就是离轴全息:参考光与物光之间设置一个微小角度,让干涉条纹的频率移到频谱中便于分离的位置。这时候全息图在频谱上三个分量:零级像(直流项)和两个共轭孪生像。用FFT可以直观看到,中心亮斑是零级像,旁边对称分布的亮斑是正负一级分量。
数值重建第一步是频域滤波:把全息图做FFT,找到一级像频谱周围的一个区域,用矩形窗或圆窗提取出来,再做逆FFT回到空域。滤波窗口大小很讲究,窗口太小会丢掉高频细节,重建像分辨率下降;窗口太大又会混入零级像边缘的残留,使重建像叠加一个低频背景。实践中可以先大致看一下频谱图,用灰度阈值自动估计一级像中心位置和半径,再做精细调节。
第二步是相位补偿。离轴全息里的倾斜参考光给全息图叠加了一个线性相位因子,即使做完频谱滤波,这个倾斜项也不会自动消失。通常用数值方法构造一个和实验参考光匹配的虚拟参考光,乘以提取的频谱分量,消除倾斜产生的位移。具体操作时,可以在全息图边缘选一个无物体区域测量条纹周期,反推参考光的角度参数。
最后,经过补偿后的复振幅做一次衍射传播,就等于重建了物面上各点的高度和相位信息。这个传播距离通常等于记录距离。整个过程环环相扣,任何一步的参数不准确,最终重建像的横向分辨率和相位精度都会打折。
5.2 重建距离不知道怎么办:自聚焦指标与一次扫描计算
全息重建最让人头疼的问题之一是光学平台的记录距离未必精确,尤其是手持式或流水线在线的全息测量系统,CCD到物体之间的距离是未知的或实时变化的。距离错一毫米,重建像的模糊程度立刻体现出来。这时候就需要自动对焦——用一系列重建距离扫描,找到最清晰位置对应的距离。
自聚焦的核心是清晰度评价指标,常见的有:
| 评价函数 | 原理 | 适用场景 |
|---|---|---|
| 拉普拉斯方差 | 对图像求二阶微分后计算方差 | 对比度一般的生物样品 |
| 梯度能量(平方梯度) | 计算相邻像素差的平方累加 | 具有锐利边缘的结构 |
| 频谱能量(高频能量) | 计算FFT后高频分量的能量占比 | 纹理细节丰富的样品 |
| 修正的倒数方差 | 像素强度倒数变化幅度 | 低信噪比条件下的稳健判断 |
实际操作时,不需要每个距离都从全息图重新做完整传播。可以先把全息图频谱滤波后的复振幅固定下来,只在这个复振幅上做不同距离的衍射传播。因为传播计算可以预先缓存输入面的FFT,每次改变距离只需要更新传播函数的频域表达式,做一次乘法再加一个逆FFT,速度可以大幅提升。一个典型的512×512全息图,扫描50个距离做一次自动对焦,在普通台式机上几秒内就能完成。
关于评价指标,我个人的体会是拉普拉斯方差往往比频谱能量更稳定,因为频谱能量容易被周期性伪影欺骗,遇到条纹噪声时容易把噪声当成“高频清晰度”而错误对焦。如果样品本身含有轻微周期结构,建议至少同时计算两三个指标,取一致对焦距离,减少偶然情况。
5.3 衍射计算在其他场景的扩展:相干衍射成像与超分辨
数字全息只是相干光成像计算应用的一个分支。把衍射计算和相位恢复结合起来,还能做传统方法做不到的事。相干衍射成像(CDI)的思路是,不用参考光,只用一束照明光照射孤立样品,在远场记录衍射图,然后用相位恢复迭代反推出样品复振幅。因为没有参考光,系统结构简单、对机械稳定性要求降低,在极紫外和X射线波段尤其受关注。
基于CDI又衍生出布拉格相干衍射成像(Bragg CDI),用来测量纳米晶体的应变场;还有前面提到的PIE,把CDI对孤立样品的要求放松到任意延展样品,配合扫描机制做高分辨成像。这些方法的核心,都是衍射计算加上相位恢复迭代,区别只在数据采集策略和约束条件的设计上。可以说,掌握了离散衍射计算的本质,再学会相位恢复的工程实现,就等于解锁了计算成像的一整个工具箱。
超分辨方向上也有类似的思路。常规光学系统受衍射极限限制,分辨率大概在照明波长的一半量级。结合结构光照明显微(SIM)、频谱边缘外推或叠层成像中的频谱拼接技术,可以把有效频谱范围向外扩展,获得超越经典衍射极限的分辨能力。这些技术最终都落到“算”上——如何精确建模光学系统的结构函数、如何在被噪声破坏的频谱数据中提取真实信息。学透相干光成像的数值计算,看这些前沿方法时会轻松很多。
最后再分享一个我用了很久的工作习惯:每次搭建衍射计算流程,先写一个严格解析解可检验的简单场景,比如圆孔衍射或高斯光束传播,用解析结果跟数值结果做逐点对比。只有测试用例吻合到足够精度,才敢把代码用于未知样品的处理。这个习惯帮我挡掉了太多因单位、补零、距离换算搞错而导致的无效实验结果。相干光成像的计算没有太多玄学,把每一步的物理假设和数值离散条件搞清楚,结果自然就对了。