人类能看懂的衍射光学（含基尔霍夫衍射，瑞利--索末菲衍射，夫琅禾费衍射，角谱衍射，菲涅尔衍射积分，菲涅尔衍射的S-FFT算法,T-FFT算法,D-FFT算法）_瑞利索末菲衍射模型

考虑一些有基础的小伙伴的需求，开头先给结论，感兴趣的可以慢慢看后面的推导过程。

结论

基尔霍夫，瑞利–索末菲公式：$U(x,y,d)=\frac{1}{j\lambda }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)\frac{exp(jkr)}{r}K(\theta )d{x}_{0}d{y}_0$‘

1.$K(\theta )=\frac{cos(\theta )+1}{2}$ 称为基尔霍夫公式;

2.$K(\theta )=cos(\theta )$ 称为第一种瑞利–索末菲公式;

3.$K(\theta )=1$ 称为第二种瑞利–索末菲公式;

傍轴近似：积分函数$\frac{exp(jkr)}{r}$中的分母r可以由d来代替

角谱：$G_d(f_x,f_y)=G_0(f_x,f_y)exp[j\frac{2\pi }{\lambda }z\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}]$
$U(x,y,d)={ϝ}^{-}[ϝ[U(x,y,0)]H(f_x,f_y)]$

菲涅尔近似：积分函数$\frac{exp(jkr)}{r}$中的分子指数部分中的$r\approx d+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d}$

菲涅尔衍射积分：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$

夫琅禾费近似：积分函数$\frac{exp(jkr)}{r}$中的分子指数部分中的$r\approx d+\frac{x^2+y^2}{2d}-\frac{x{x}_0+y{y}_0}{d}$

夫琅禾费衍射：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{-j2\pi }{\lambda d}(x{x}_0+y{y}_0)]d{x}_{0}d{y}_0$

S-FFT：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]ϝ[U_0(x_0,y_0,0)exp[\frac{jk}{2d}(x_0^2+y_0^2)]]$

T-FFT：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{ϝ}^{-}[ϝ[U_0(x_0,y_0,0)]ϝ[exp[\frac{jk}{2d}(x^2+y^2)]]]$

D-FFT：$U(x,y,d)={ϝ}^{-}[ϝ[U_0(x_0,y_0,0)]exp[jkd(1-\frac{{\lambda }^2}{2}(f_x^2+f_y^2))]]$








接下来是上述衍射公式的推导过程：

首先，对一些基本概念和衍射的意义进行说明。

基本概念和衍射的意义

衍射（diffraction）是指波遇到障碍物时偏离原来直线传播的物理现象。

波面的任何变形（通过相位物体）或者说波面（波前）上光场的复振幅分布受到任何空间调制，都将导致衍射现象的发生，而使通过障碍物以后的光场的复振幅重新分布。

衍射屏：导致衍射发生的障碍物

如果不涉及光传播与变换过程中障碍物或光学元件结构尺寸接近于光波长的情况，对衍射问题的研究不邻接衍射平面。
(这句话我TM想了好久没想懂邻接什么意思，查阅的时候一个个引用这句话引用得飞起但是都没解释，以我浅薄的理解就是在纳米级别的情况下衍射作用可能不会发生，也就是本来应该导致衍射发生的障碍物不再认为是障碍物。)

标量衍射理论：将电场强度E视为标量，描述光传播的物理过程和求解波动方程的方法。

傍轴近似：当光束发散较小且观测区域的宽度远小于光传播距离。

实际应用中需要计算的往往是垂直于光学系统光轴或光束传播方向的空间平面上的光波场，并且通常将z轴视为光学系统的光轴。也就是说，需要表述的平面是与光轴垂直的平面。
所以我们先以如下三维图建立一个坐标系方便后续公式推导。设衍射屏与观察屏的距离为d，矢径$r=\sqrt{(x-x_0{)}^2+(y-y_0{)}^2+d^2}$,$\theta$ 代表点$(x_0,y_0,0)$ 到点$(x,y,d)$ 的矢径 r 与 $(x_0,y_0,0)$点法线$n$的夹角，波矢$k=2\pi /\lambda$ ，$U_0(x_0,y_0,0)$以及$U(x,y,d)$分别为衍射屏及观察屏上的光波的复振幅，$G_0(f_x,f_y)$,$G_d(f_x,f_y)$分别为其对应的频谱函数。

基尔霍夫公式，瑞利–索末菲公式，以及角谱传播公式分别是亥姆霍兹方程的三个解（光是一种电磁波，故其传播特性可由经典电磁场理论的麦克斯韦方程推出，这里不具体展开，其实是笔者对于波动方程真的深恶痛绝）。

在标量衍射理论框架下，光传播的物理过程可以严格地由基尔霍夫公式和瑞利–索末菲公式以及衍射的角谱的传播公式表示。

如果知道空间中垂直于光传播方向的一个平面上的光波场，便能计算该平面前后的空间光波场分布尤其是在观察屏上的分布。
（xd后面没有图片了，全是密密麻麻的公式和文字）

基尔霍夫，瑞利–索末菲公式

如果只是告诉大家基尔霍夫，瑞利–索末菲公式是亥姆霍兹方程的解而来可能难以理解其物理含义，但如果大家通过惠更斯–菲涅尔原理，以及场论中的格林定理电磁场的边值条件等一步步推导就能逐渐窥见那些大佬的逻辑所在。这里就不展开来了，移步工程光学一书吧。

基尔霍夫，瑞利–索末菲公式两者在数学上有统一的公式：$U(x,y,d)=\frac{1}{j\lambda }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)\frac{exp(jkr)}{r}K(\theta )d{x}_{0}d{y}_0$

（其中$r=\sqrt{(x-x_0{)}^2+(y-y_0{)}^2+d^2}$,$\theta$ 代表点$(x_0,y_0,0)$ 到点$(x,y,d)$ 的矢径 r 与 $(x_0,y_0,0)$点法线$n$的夹角，,$K(\theta )$为倾斜因子（由惠更斯–菲涅尔原理可知，$K(\theta )$表示子波的振幅随面元法线与子波所在点到传播点的夹角$\theta$的变化），而不同的倾斜因子对应不同的公式：

1.$K(\theta )=\frac{cos(\theta )+1}{2}$ 称为基尔霍夫公式;

2.$K(\theta )=cos(\theta )$ 称为第一种瑞利–索末菲公式;

3.$K(\theta )=1$ 称为第二种瑞利–索末菲公式;

在实际问题中，$\theta$通常较小，三个公式的倾斜因子均接近1。

观察屏上的光波的复振幅 $U(x,y,d)=\frac{1}{j\lambda }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)\frac{exp(jkr)}{r}d{x}_{0}d{y}_0$

在傍轴情况（当光束发散较小且观测区域的宽度远小于光传播距离）下，积分函数$\frac{exp(jkr)}{r}$中的分母r可以由d来代替，因为这对于光波振幅无大的影响。$r=\sqrt{(x-x_0{)}^2+(y-y_0{)}^2+d^2}$，这是由于傍轴情况下光束与z轴平行传播时不会发散，也就是在x轴y轴上的分量不会有移动故$(x-x_0{)}^2+(y-y_0{)}^2=0$，$r\approx d$。
但是对于指数部分的r不能进行这样简单的处理，其原因是指数部分r的轻微变化会引起甚至大于$2\pi$的相位变化。

由$r=\sqrt{(x-x_0{)}^2+(y-y_0{)}^2+d^2}$(起初建立坐标系时就有的)

对$exp(jkr)$中的 r 进行泰勒展开，$r\approx d[1+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d^2}-\frac{[(x-x_0{)}^2+(y-y_0{)}^2{]}^2}{8d^4}]$……,第两项之后的值过小，只保留泰勒展开的前两项，故

$exp(jkr)$中的$r\approx d+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d}$
综上，$\frac{exp(jkr)}{r}\approx \frac{exp[jk[d+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d}]]}{d}$

化简后观察屏上的光波的复振幅 $U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$，

上述那条公式也就是菲涅尔衍射积分。

而指数部分$r\approx d+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d}$一般被称为菲涅尔近似，满足菲涅尔近似条件的衍射称为菲涅尔衍射。

菲涅耳衍射是光源或观察屏或两者距离衍射屏有限远时产生的衍射，故也称为近场衍射。

角谱衍射

亥姆霍兹方程的第三个解便是角谱公式，在角谱理论中，光波沿z方向传播的结果，在频域内表现为将衍射屏上的光波场的频谱 $G_0（f_x,f_y）$乘以一个与z有关的相位延迟因子$exp[j\frac{2\pi }{\lambda }z\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}]$，等效于通过一个半径为1/$\lambda$的理想低通滤波器。在线性系统理论中，该相位延迟因子即为衍射在频域的传递函数$H(f_x,f_y)$，表明衍射问题可以视为是光波场通过一个线性空间不变系统的变换过程。（$f_x$ ，$f_y$为频域谱上对应的坐标系，也就是空间频率，$f_x=\frac{x}{\lambda d}$,$f_y=\frac{y}{\lambda d}$。）

叨叨叨一堆话是不是不理解为什么呢？其实角谱理论就是将平面上的复振幅分布分解为多个不同传播方向的平面波的线性叠加，类似于傅里叶变换把一个时域分量分解成多个频率分量的总和那样。

总结成公式就更为简洁易懂了，就是

$G_d(f_x,f_y)=G_0(f_x,f_y)exp[j\frac{2\pi }{\lambda }z\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}]$

（$G_0(f_x,f_y)$,$G_d(f_x,f_y)$分别为$U(x,y,0)$和$U(x,y,d)$的频谱函数。）

由上式可知若$1-(\lambda f_x{)}^2-(\lambda f_y{)}^2<0$，则角谱分量将随z的增大按指数规律急剧衰减，光波只存在于邻近衍射屏的一个非常薄的区域，故只有$1-(\lambda f_x{)}^2-(\lambda f_y{)}^2>0$，即$f_x^2+f_y^2<\frac{1}{{\lambda }^2}$，角谱分量才会到达观察屏。

（以$f_x=f_y$，$\lambda$=632nm为例，$f_x=f_y<1.1188e+06$ m）

简写下上面的公式，设衍射距离为d，定义角谱衍射的传递函数为

$H(f_x,f_y)=exp[j\frac{2\pi }{\lambda }d\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}$

则角谱衍射公式可简写为

$U(x,y,d)={ϝ}^{-}[ϝ[U_0(x_0,y_0,0)]H(f_x,f_y)]$

当傍轴近似（当光束发散较小且观测区域的宽度远小于光传播距离）时，即$f_x,f_y$较小时，$\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}$泰勒系数展开，同样只保留前两项

（$\sqrt{(1-(\lambda f_x{)}^2-(\lambda f_y{)}^2}=1-\frac{1}{2}{\lambda }^2(f_x^2+f_y^2)+\frac{1}{8}{\lambda }^4(f_x^2+f_y^2)$……）,故

角谱衍射的传递函数$H(f_x,f_y)=exp[jkd(1-\frac{{\lambda }^2}{2}(f_x^2+f_y^2))]$

这只是对$H(f_x,f_y)$的一个小近似，然后根据傅里叶变换性质，频域分量的相乘等于空域分量的卷积。

$G_d(f_x,f_y)=G_0(f_x,f_y)H(f_x,f_y)$便变换为

$U(x,y,d)=U_0(x_0,y_0,0)\ast {ϝ}^{-}[H(f_x,f_y)]$

而${ϝ}^{-}[H(f_x,f_y)]$有解析解，即${ϝ}^{-}[H(f_x,f_y)]=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]$

角谱衍射卷积形式的公式即可改写为菲涅尔衍射积分

$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$

兜兜转转，最后还是你。尽管角谱衍射公式，基尔霍夫公式以及瑞利–索末菲公式有不同的形式，但它们的傍轴近似具有相同的形式，也就是菲涅尔积分形式。

夫琅禾费衍射

若菲涅尔衍射积分的傅里叶变换中，如果$d>>\frac{k(x_0^2+y_0^2{)}_{max}}{2}$，即衍射屏距离观察屏无限远，则$\frac{jk(x_0^2+y_0^2)}{2d}\approx 0$，$exp[\frac{jk(x_0^2+y_0^2)}{2d}]\approx 1$，故夫琅禾费衍射场可简单地变为$U_0(x_0,y_0,0)$的傅里叶变换。

$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{-j2\pi }{d\lambda }(x{x}_0+y{y}_0)]d{x}_{0}d{y}_0$

将上式对比于基尔霍夫–瑞利索末菲公式$U(x,y,d)=\frac{1}{j\lambda }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)\frac{exp(jkr)}{r}d{x}_{0}d{y}_0$，$exp[\frac{jk(x_0^2+y_0^2)}{2d}]\approx 1$也就是对指数部分r进行近似。

$r\approx d+\frac{x^2+y^2}{2d}-\frac{x{x}_0+y{y}_0}{d}$

也就是在菲涅尔近似 $r\approx d+\frac{(x-x_0{)}^2+(y-y_0{)}^2}{2d}$ 的基础上$x_0^2+y_0^2=0$ ,可以理解为对比于在无限远的观察屏，衍射屏上的坐标分量已经小到可以忽略不计了。

上述近似即被称为夫琅禾费近似，满足夫琅禾费近似的衍射即被称为夫琅禾费衍射。

夫琅禾费近似成立所要求的条件较为苛刻，但如果物平面光波场为一个会聚的球面波，则夫琅禾费近似很容易满足，而会聚球面波应用上可以很容易通过透镜实现，这也是观察夫琅禾费衍射情况下在观察屏前加透镜的原因。








接下来是对菲涅尔衍射积分的三种常用算法的推导：

S-FFT（single fast Fourier transform algorithm）:单次傅里叶变换算法

菲涅尔衍射积分： $U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$

对上式二次相位因子$exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]$进行展开，波矢 $k=2\pi /\lambda$ ，并将与积分变量无关的项提到积分号前则得到

$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x,y,0)exp[\frac{jk}{2d}[(x_0^2+y_0^2]]exp[-j2\pi (\frac{x}{\lambda d}{x}_0+\frac{y}{\lambda d}{y}_0)]d{x}_{0}d{y}_0$

对比于傅里叶变换的公式

$F(f_x,f_y)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }f(x,y)exp[-j2\pi (f_{x}x+f_{y}y)]dxdy$

便可以将二重积分项等同为$U_0(x,y,0)exp[\frac{jk}{2d}(x_0^2+y_0^2)]$的傅里叶变换，设$G_1(f_x,f_y)$是$U_0(x,y,0)exp[\frac{jk}{2d}(x_0^2+y_0^2)]$的傅里叶变换

且$f_x=\frac{x}{\lambda d}$ ， $f_y=\frac{y}{\lambda d}$

故菲涅尔衍射积分的S-FFT算法为：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]G_1(f_x,f_y)$

T-FFT(triple fast Fourier transform algorithm):三次傅里叶变换算法

（两次傅里叶变换一次逆傅里叶变换），这种算法常被称为菲涅尔衍射的卷积算法。

菲涅尔衍射积分： $U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$

类比于卷积公式 $y(t)={\int }_{-\infty }^{\infty }x(p)h(t-p)dp=x(t)⨂$h(t)

定义脉冲响应

$h(x,y)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]$

便可将菲涅尔衍射积分可转化为$U_0(x_0,y_0,0)$和$h(x,y)$的卷积。

$U(x,y,d)=U_0(x_0,y_0,0)⨂h(x,y)$

把脉冲项的无关项提前，则式子变为

$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{U}_0(x_0,y_0,0)⨂exp[\frac{jk}{2d}(x^2+y^2)]$

但由于在空域完成卷积运算比较花机器时间，所以通常要转换到频域进行计算。根据频谱特性，空域的卷积等于频域的乘积。故需要对$ U_0(x_0,y_0,0)$和$exp[\frac{jk}{2d}(x^2+y2)]$各做一次傅里叶变换，再对它们的乘积做一次逆傅里叶变换。

故菲涅尔衍射积分的T-FFT算法为：$U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{ϝ}^{-}[ϝ[U_0(x_0,y_0,0)]ϝ[exp[\frac{jk}{2d}(x^2+y^2)]]]$

与S-FFT算法不同的是，观察面的尺寸与衍射面的尺寸是一样的，是因为在计算中在频域进行相乘时，要求频谱中相同频率分量对应相乘，故最高频率也就必须相同，抽样数一样对应的几何尺寸也就相等了。

D-FFT（double fast Fourier transform algorithm）:两次傅里叶变换算法

与T-FFT算法一样，该算法也常被称为卷积算法。

菲涅尔衍射积分： $U(x,y,d)=\frac{exp(jkd)}{j\lambda d}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{U}_0(x_0,y_0,0)exp[\frac{jk}{2d}[(x-x_0{)}^2+(y-y_0{)}^2]]d{x}_{0}d{y}_0$

把D-FFT放在T-FFT后面讲是因为这两者有相似之处。

同样设脉冲响应 $h(x,y)=\frac{exp(jkd)}{j\lambda d}exp[\frac{jk}{2d}(x^2+y^2)]$

与T-FFT算法类似，菲涅尔积分变换为$U(x,y,d)=U_0(x_0,y_0,0)⨂h(x,y)$

从这里开始，D-FFT与T-FFT的算法的推导就不一样了。

设$H(f_x,f_y)$为$h(x,y)$对应的频谱函数（$H(f_x,f_y)$也被称为菲涅尔传递函数。），由于$h(x,y)$ 的傅里叶变换有解析解，即$H(f_x,f_y)=exp[jkd(1-\frac{{\lambda }^2}{2}(f_x^2+f_y^2))]$

则菲涅尔衍射积分的D-FFT算法：$U(x,y,d)={ϝ}^{-}[ϝ[U_0(x_0,y_0,0)]exp[jkd(1-\frac{{\lambda }^2}{2}(f_x^2+f_y^2))]]$

与T-FFT算法不同，D-FFT，已经知道传递函数$H(f_x,f_y)$的解析表达式，相当于把T-FFT对$h(x,y)$ 的FFT变换先行计算了出来，只需要将空间频率$f_x$ ，$f_y$代入解析式即可得到传递函数的值而不用做FFT变换，但同时也要求必须通过实验参数给出正确的$f_x,$$f_y$

对比于角谱函数，可以发现角谱衍射的传递函数在傍轴近似条件下等于菲涅尔衍射的传递函数，也就是角谱衍射的公式在傍轴条件下等于菲涅尔积分的D-FFT算法。

总结

到此为止，已对常见的衍射公式进行了大致的介绍和简单的推导，建议你再到开头去看下结论，希望你能对这些经典的衍射公式有更深的理解。（然后还是要吐槽下国内的很多光学教材都编得不咋地，很多语句都是抄来抄去，引用起来也没有相关解释，很多浅显的公式道理也讲得繁琐复杂(╯‵□′)╯︵┻━┻）

（附注：因为是手敲公式所以可能会有些许疏漏还望见谅QAQ，此外，也欢迎感兴趣有疑问的小伙伴在评论区留言讨论）

参考文献：
《信息光学数字实验室（Matlab版）》——钱晓凡
《衍射计算及数字全息》—— 李俊昌
《计算机制全息图》——金国藩，张浩，苏萍，蒋强，曹良才
《工程光学》——郁道银，谈恒英
《信息光学》——苏显渝