信息光学教程李俊昌第三章.pdf

1、第3章暋衍射积分的快速傅里叶变换计算在第2章,我们基于标量衍射理论推导出了严格满足亥姆霍兹方程的角谱衍射公式及其傍轴近似 菲涅耳衍射积分;并且指出,理论上还存在严格满足亥姆霍兹方程的基尔霍夫公式及瑞利索末菲公式 1 .为便于叙述,以上4种衍射计算公式统一称为经典衍射计算公式 2 .对经典的衍射公式研究表明 2 - 4 ,它们均能用傅里叶变换表示,从而通过傅里叶变换求解.但是,对于实际给定的衍射问题,能够直接从傅里叶变换求出解析表达式的函数非常有限,在研究实际问题时,不得不将函数按一定规律在二维空间进行取样及延拓,变成该函数的周期离散分布作离散傅里叶变换.然而,离散傅里叶变换的计算仍然十分繁杂,

2、如果没有计算机,事实上很难完成一个可以解决实际问题的计算工作. 1 9 6 5年,由库利图基( C o o l e y 灢 T u k e y ) 5 提出的快速傅里叶变换技术( t h e f a s t F o u r i e r t r a n s f o r m , F F T )彻底改变了这种状况,计算机的普及应用为这种快速计算方法的推广创造了良好的条件.因此,利用快速傅里叶变换技术计算衍射的方法逐渐被广泛采用.本章首先介绍离散傅里叶变换与傅里叶变换的关系,然后,基于离散傅里叶变换的基本理论及取样定理,对经典衍射积分的快速傅里叶变换计算方法进行研究.由于柯林斯公式也可以表示成傅里叶变

3、换的形式,其计算方法可以推广到柯林斯公式的计算.本章的所有计算均能通过附录A中的M A T L A B程序L I M 1 . m及L I M 2 . m实现.3 . 1 暋离散傅里叶变换与傅里叶变换的关系3.1.1暋空域连续函数的离散及延拓函数作二维离散傅里叶变换时,要求被变换函数是二维空间的周期离散函数 6 - 8 .由于实际需要作傅里叶变换的函数通常是在空域无限大平面上均有定义的连续函数,必须将函数截断在有限的区域进行取样及延拓.图3 灢 1 灢 1给出二维空域连续函数的离散及延拓示意图.图3 灢 1 灢 1 暋空域连续函数的离散及延拓暋 暋图中,左上方灰色图像给出一连续函数的分布区域.对

4、连续函数通常的取样方法是,先将函数的主要部分通过坐标变换放在第一象限,并沿平行于坐标轴的方向将函数截断在一个Lx暳Ly的矩形区域内;然后,取样周期为Tx=Lx/Nx,Ty=Ly/Ny,从坐标原点开始将函数离散为Nx暳Ny个点的二维离散分布值.图3 灢 1 灢 1 ( a ) 、 ( b )描述了上述过程(图中用黑点标注取样点落在函数定义区域上的位置,用小圆圈表示取样为零的位置) ,图3 灢 1 灢 1 ( c )是二维周期延拓结果.3.1.2暋离散傅里叶变换与傅里叶变换的关系很明显,函数经截断及离散处理后无论在空域及频域均会引入误差.现以x方向的傅里叶变换为例进行研究,以后再将结果推广到二维空

5、间.图3 灢 1 灢 2表示对于某一给定的y,函数沿x方向进行离散傅里叶变换的过程.图中,左边为一列空域的原函数图像,右边一列图像是它们的频谱的模,符号表示它们为傅里叶变换对.例如,图3 灢 1 灢 2 ( a 1 )为空域的原函数g(x,y) ,图3 灢 1 灢 2 ( a 2 )为它的频谱G(fx,y)的模G(fx,y) .对未经截断函数的取样,等于用图3 灢 1 灢 2 ( b 1 )的梳状函数毮Tx(x)乘以图3 灢 1 灢 2 ( a 1 )的原函数,数学表达式为 7 , 8 gTx(x,y)=g(x,y)毮Tx(x)=g(x,y) 暺曓n=-曓毮(x-nTx) ( 3灢1灢1 )暋

6、 暋由于梳状函数毮Tx(x)为周期Tx的毮函数,可以表示为傅里叶级数毮Tx(x)=暺曓n=-曓毮(t-nTx)=暺曓k=-曓Ake x p jk2 毿Txx式中, j=-1 ,Ak=1Tx曇Tx/ 2-Tx/ 2毮Tx(x) e x p-jk2 毿Txxdx=1Tx.于是gTx(x,y)=g(x,y) 1Tx暺曓k=-曓e x p jk2 毿Txx暋 暋上式表明,取样信号已经不是原信号,而是无穷多个截波信号1Tx暺曓k=-曓e x p jk2 毿Txx被信号gx,( )y调制的结果图3 灢 1 灢 2 ( c 1 ) .现在,通过傅里叶变换来考察信号经取样后的频谱与原信号频谱的关系.对上式作傅

7、里叶变换得GTxfx,( )y=曇曓-曓gTxx,( )ye x p-j 2 毿fx( )xdx=曇曓-曓g(x,y) 1Tx暺曓k=-曓e x p jk2 毿Txxe x p-j 2 毿fx( )xdx=1Tx暺曓k=-曓曇曓-曓gx,( )ye x p-j 2 毿fx-kT xtdx=1Tx暺曓k=-曓Gfx-kTx, y( 3灢1灢2 )结果表明,在取样信号频谱GTxfx,( )y中除了包含原信号频谱Gfx,( )y外,还包含了无穷多个被延拓的频谱,延拓的周期为1 /Tx图3 灢 1 灢 2 ( c 2 ) .并且,由于原函数的频谱宽度大于延拓的周期1 /Tx,相邻的频谱曲线产生了混叠.

8、64信息光学教程图3 灢 1 灢 2 暋函数的离散傅里叶变换过程根据傅里叶变换中频域的卷积定律,图3 灢 1 灢 2 ( c 2 )也可以通过原函数的频谱函数Gfx,( )y图3 灢 1 灢 2 ( a 2 ) 与梳状函数的频谱函数殼Txf( )x图3 灢 1 灢 2 ( b 2 ) 的卷积求出GTxfx,( )y=Gfx,( )y*殼Txf( )x( 3灢1灢3 )为强调这个关系,图3 灢 1 灢 2 ( c 2 )的纵坐标由这个卷积表达式标注.由此可见,连续函数经过周期为Tx的无穷毮序列取样离散后,其频谱与原函数频谱相比有如下两点区别:( 1 )频谱发生了周期为1 /Tx的周期延拓.如果原

9、函数的频谱宽度大于1 /Tx时,则产生频谱混叠,引入失真.( 2 )离散信号频谱GTxfx,( )y的幅度是原函数频谱Gfx,( )y的1 /Tx倍.74第3章暋衍射积分的快速傅里叶变换计算然而,上面对连续函数被无穷毮序列取样离散后的频谱研究只是一个理论结果,因为实际上不可能作取样点为无限多的数值计算.并且,由于离散傅里叶变换事实上讨论的是在空域及频域均是周期离散函数的傅里叶变换问题,还要将离散函数截断及延拓才能满足要求.因此,将空域非周期的离散函数图3 灢 1 灢 2 ( c 1 ) 先通过下述矩形窗函数图3 灢 1 灢 2 ( d 1 ) 截断:rLx(x)=1, 暋-Tx/ 2 xLx-

10、Tx/ 20( 3灢1灢4 )得到具有Nx个点的离散分布图3 灢 1 灢 2 ( e 1 ) gTxr(x,y)=g(x,y)毮Tx(x)rLx(x) ( 3灢1灢5 )然后,再将截断后的部分进行周期为Lx的延拓,形成如图3 灢 1 灢 2 ( g 1 )所示的周期离散序列gTxrk(x,y)=gTxr(x+kLx,y) , 暋k=0 ,暲1 ,暲2 , ( 3灢1灢6 )按照傅里叶变换理论,空域中矩形窗函数图3 灢 1 灢 2 ( d 1 ) 与离散序列图3 灢 1 灢 2 ( c 1 ) 的乘积的频谱函数,可表示为矩形函数的频谱函数RLxf( )x图3 灢 1 灢 2 ( d 2 ) 与图

11、3 灢 1 灢 2 ( c 1 )的频谱函数图3 灢 1 灢 2 ( c 2 ) 的卷积GTxrfx,( )y=Gfx,( )y*殼Txf( ) x*RLxf( )x( 3灢1灢7 )对应的频谱函数曲线如图3 灢 1 灢 2 ( e 2 )所示.由图可见,由于矩形窗函数的频谱RLxf( )x具有较大的起伏变化的旁瓣,卷积运算的结果使图3 灢 1 灢 2 ( e 2 )的频谱曲线形状产生了失真(为说明问题,图中略有夸大) .将图3 灢 1 灢 2 ( e 2 )与图3 灢 1 灢 2 ( a 2 )比较不难发现,现在得到的是带有畸变的原函数频谱的周期延拓曲线,延拓周期为1 /Tx.由于离散傅里叶

12、变换是对空域及频域均为周期离散函数的变换,因此,图3 灢 1 灢 2 ( e 2 )的曲线还将被周期为1 /Lx的梳状函数图3 灢 1 灢 2 ( f 2 ) 取样,其结果是一个周期为Nx的频域的离散函数图3 灢 1 灢 2 ( g 2 ) .在频域进行上面频谱函数与梳状函数的乘积取样时,就对应着它们在空域原函数的卷积运算.图3 灢 1 灢 2 ( e 1 )与图3 灢 1 灢 2 ( f 1 )的函数在空域卷积运算的结果成为一周期为Nx的空域离散函数图3 灢 1 灢 2 ( g 1 ) .空域及频域离散函数均以Nx为周期,我们只要分别知道一个周期内的离散值或样本点便可以了解离散函数全貌.离散

13、傅里叶变换或其快速算法F F T便是完成从空域到频域,以及从频域到空域的这Nx个样本点的计算方法.至此,我们已经知道,离散傅里叶变换是傅里叶变换的一种近似计算.只要能够将衍射的计算表示为卷积的形式,并了解离散傅里叶变换与傅里叶变换间的定量关系,采取合适的措施抑制畸变,便能对衍射问题求解.关于抑制矩形窗函数的旁瓣引起的畸变问题,可以选择其他形式的窗函数对离散序列进行截断,但这样做将会在空域引入额外的能量损失,对于特别关心衍射场能量问题的衍射计算未必可行,详细的讨论可以参看有关数字信号处理的专著(例如文献 8 ) .3 . 2 暋菲涅耳衍射积分的快速傅里叶变换计算菲涅耳衍射积分是应用研究中最广泛使

14、用的公式,因此,首先对它的计算进行讨论.鉴于菲涅耳衍射积分可以表示为傅里叶变换及卷积两种形式,存在一次快速傅里叶变换及快速卷84信息光学教程积算法.快速卷积算法需要作一次快速傅里叶变换及一次快速傅里叶反变换运算,为表述方便,以下称第一种算法为S 灢 F F T算法,第二种算法为D 灢 F F T算法.我们将看到,合理选择两种方法才能正确处理实际衍射问题.3.2.1暋菲涅耳衍射积分的S灢FFT算法第2章式( 2 灢 2 灢 3 4 b )中设U0x0 ,y( )0为物平面光波复振幅,令Ux,y,( )d=Ux,( )y及Ux,y,( )0=U0x,( )y,经距离d的衍射到达观测平面的光波复振幅

15、Ux,( )y可由以下形式的菲涅耳衍射积分表示:Ux,( )y=e x p j( )kdj毸de x p jk2dx2+y( ) 2暳曇 曇曓-曓U0x0 ,y( )0 e x p jk2dx20+y( ) 20暋 暳 e x p-j 2 毿x0x毸d+y0y毸 ddx0 dy0 ( 3灢2灢1 )暋 暋式( 3 灢 2 灢 1 )与熟知的二维傅里叶变换式进行比较可以看出,式( 3 灢 2 灢 1 )的主要计算是进行物函数与指数相位因子乘积U0x0 ,y( )0 e x p jk2dx20+y( ) 20的傅里叶变换,但变换结果还要再乘以一个二次相位因子.若利用快速傅里叶变换F F T计算上式

16、,物平面取样宽度为L0 ,取样数为N暳N,取样间距殼x0=殼y0=L0 /N,式( 3 灢 2 灢 1 )可写为Up殼x,q殼( )y=e x p j( )kdj毸de x p jk2dp殼( )x2+q殼( )y( ) 2暋暳F F TU0m殼x0 ,n殼y( )0 e x p jk2dm殼x( )0 2+n殼y( )0( ) 2p殼x毸d,q殼y毸d暋 暋暋 ( 3灢2灢2 )(p,q,m,n=-N2 ,-N2+1 , ,N2-1 )式中, 殼x=殼y是离散傅里叶变换后对应的空域取样间距.为确定这个数值,根据本章离散傅里叶变换与傅里叶变换关系的讨论,上式的计算结果是取值范围1 / 殼x0的

17、N暳N个点的离散值,即L毸d=1殼x0=NL0 ( 3灢2灢3 )或者L=毸dNL0( 3灢2灢4 )因此殼x=殼y=LN=毸dL0( 3灢2灢5 )暋 暋以上结果表明,如果保持物平面取样间隔及观测区域不变,离散傅里叶变换计算结果在观测平面上衍射图样的取样范围L不但是光波长毸及取样数N的函数,并且将随衍射距离d增加而增加.对于给定的光波长毸及物平面取样范围L0 ,当衍射距离d很小时,如果取样数N保持不变,计算结果只对应于观测平面上邻近光轴的很小区域的衍射图像.因此, S 灢 F F T算法主要适应于衍射距离d较大的情况.然而,只有满足取样定律的计算才能获得正确的计算结果.分析式( 3 灢 2

18、灢 1 )可知,被变换函94第3章暋衍射积分的快速傅里叶变换计算数由物函数与指数相位因子的乘积组成,理论分析容易证明,指数相位因子e x p jk2dx20+y( ) 20的傅里叶变换为毸dj e x p-j毸d毿x毸 d2+y毸 d 2,是在整个频域都有取值的非带限函数.按照频域卷积定理,U0x0 ,y( )0 e x p jk2dx20+y( ) 20的频谱是指数相位因子频谱与物函数U0x0 ,y( )0频谱的卷积,无论物函数是否是带限函数,卷积运算结果都是非带限函数.因此,U0x0 ,y( )0 e x p jk2dx20+y( ) 20也将是整个频域都有取值的非带限函数,要通过式( 3

19、 灢 2 灢 2 )的D F T计算准确对式( 3 灢 2 灢 1 )求解是不可能的.然而,在形式上奈奎斯特取样定理可以视为空域取样间距的倒数大于或等于函数最高频谱的两倍,即在最高频谱所对应的空间周期上至少要有两个取样点.在实际衍射计算中,通常基于下面的分析来让计算近似地满足取样定律 2 - 4 .通常情况下,物函数相对于指数相位因子的空间变化率不高,如果L0确定的范围中指数相位e x p jk2dm殼x( )0 2+n殼y( )0 2每变化一个周期2 毿时至少有两个取样点,则认为F F T计算近似满足取样定理.由于二次相位因子空间频率最高点对应于m及n等于暲N/ 2时的取样值,因此,求解以下

20、不等式可得到近似满足奈奎斯特取样定理的条件:毠毠mk2dm殼x( )02+n殼y( )0 2m,n=N/ 2曑 毿 ( 3灢2灢6 )即毿毸d暳2 殼x20N2 曑 毿由此可得殼x20 曑毸dN( 3灢2灢7 )如果只考虑衍射场的强度分布,式( 3 灢 2 灢 7 )可以作为近似满足奈奎斯特取样定理的条件.但是,一次F F T算法对式( 3 灢 2 灢 2 )的最终结果是F F T计算结果与前方二次指数相位因子的乘积.如果我们期望计算结果是满足奈奎斯特取样定理的取样,还应考虑式( 3 灢 2 灢 2 )中F F T前方二次指数相位因子的取样问题.将获得式( 3 灢 2 灢 7 )的讨论方法应用

21、于F F T前方的相位因子,可得殼x2 曑毸dN( 3灢2灢8 )但根据式( 3 灢 2 灢 4 ) ,有N殼x=毸dNN殼x0, 暋或暋 殼x=毸dN殼x0代入式( 3 灢 2 灢 8 )得毸dN殼x 02曑毸dN,即殼x20 曒毸dN.与式( 3 灢 2 灢 7 )比较可以看出,这是一组基本相互矛盾的条件.于是,只有两不等式取等号殼x0=殼x=毸dN暋或者暋L0=L=毸dN( 3灢2灢9 )才可以通过一次离散傅里叶变换计算获得满足奈奎斯特取样定理的菲涅耳衍射场离散分布.综上所述,若利用S 灢 F F T方法计算菲涅耳衍射积分,当光波长毸给定而物平面取样宽度殼L0和取样数N是可变参数时,有三

22、条主要结论:05信息光学教程其一,根据式( 3 灢 2 灢 4 ) ,计算的衍射场宽度为L=毸dN/L0 .由于衍射距离d趋近于0时,L将趋近于0 .当观测平面临近物平面时,对于给定的L0 ,必须使用庞大的取样数N才能得到期望宽度L的解.因此,使用S 灢 F F T方法将无法计算距离d趋近于0的衍射图样.其二,根据式( 3 灢 2 灢 7 ) ,物平面取样间隔满足殼x0 毸d/N或者L0 毸dN时,可以较好地计算菲涅耳衍射场的强度分布.其三,根据式( 3 灢 2 灢 9 ) ,如果要使离散计算结果近似满足奈奎斯特取样定理的衍射场,物平面及衍射场平面的取样宽度必须相等,并满足L0=L=毸dN.由

23、此可见,不可能使用S 灢 F F T变换法解决整个菲涅耳衍射区域的衍射计算问题.附录A中L I M 1 . m是用M A T L A B语言编写的S 灢 F F T变换法计算衍射的程序.通过程序的阅读及执行,可以进一步加深对S 灢 F F T变换法计算衍射的理解.3.2.2暋菲涅耳衍射的S灢FFT计算与实际测量的比较现在,以图3 灢 2 灢 1激光通过光阑的衍射实验为例,来验证菲涅耳衍射的S 灢 F F T计算及取样条件讨论的可行性.在该实验研究中,激光波长为5 3 2 n m ,经扩束及准直后照射透光孔为正十二边形孔的光阑,孔的直径(两对角顶点距离)约4 . 8 m m .用1 5 3 6

24、暳 2 0 4 8像素,像素宽度3 . 2 毺 m的C C D直接探测衍射场的强度分布.由于C C D面阵尺寸约为4 . 9 m m 暳 6 . 6 m m ,经透光孔衍射的光波的主要能量能够被C C D接收. C C D探测到的光波场能量分布将为衍射计算的可行性提供比较依据.图3 灢 2 灢 1 暋衍射实验光路使用衍射的S 灢 F F T算法时,物平面与观测平面宽度L0 ,L通常不相同.为便于比较,令L0=L=毸dN= 1 0 m m ,选择N= 2 5 6 , 5 1 2 , 1 0 2 4 ,分别求得d7 3 4 m m , 3 6 7 m m , 1 8 4 m m .将上述参数代入式

25、( 3 灢 2 灢 3 ) ,通过计算得到的光斑图像与实际测量的图像如图3 灢 2 灢 2所示.图中, C C D测量结果通过周边补零形成1 0 m m 暳 1 0 m m ( 3 1 0 3 暳 3 1 0 3像素)的图像.可以看出,理论计算与实验测量吻合很好.比较不同距离计算时使用的取样点数还可看出,观测屏离光阑越近,N越大.当使用菲涅耳衍射的S 灢 F F T计算十分邻近光阑平面的光波场时,必须使用庞大的取样数才能完成计算.利用附录A中的L I M 1 . m程序,读者很容易验证上述结果.3.2.3暋菲涅耳衍射的D灢FFT算法第2章式( 2 灢 2 灢 3 4 a )中设U0x0 ,y(

26、 )0为物平面光波复振幅,令Ux,y,( )d=Ux,( )y及Ux,y,( )0=U0x,( )y,经距离d的衍射到达观测平面的光波复振幅Ux,( )y可由以下形式15第3章暋衍射积分的快速傅里叶变换计算图3 灢 2 灢 2 暋理论模拟光斑与实际测量图像的比较( 1 0 m m 暳 1 0 m m )的菲涅耳衍射积分表示:暋Ux,( )y=e x p j( )kdj毸d曇曓-曓曇曓-曓U0x0 ,y( )0 e x p jk2dx-x( )0 2+y-y( )0 2 dx0 dy0 ( 3灢2灢1 0 )式中, j=-1 ,毸为光波长,k=2 毿 /毸.可以看出,观测平面的光波复振幅Ux,(

27、 )y被表示为物平面光波复振幅U0x0 ,y( )0与一复函数的卷积.对式( 3 灢 2 灢 1 0 )两边作傅里叶变换得FUx,( ) y=FU0x0 ,y( ) 0Fe x p j( )kdj毸de x p jk2dx2+y( ) 2 ( 3灢2灢1 1 )令fx,fy是频域坐标,定义菲涅耳衍射传递函数HFfx,f( )y=Fe x p j( )kdj毸de x p jk2dx2+y( ) 2 ( 3灢2灢1 2 )容易证明,式( 3 灢 2 灢 1 2 )存在解析解HFfx,f( )y=e x p jkd1-毸22f2x+f2( ) y( 3灢2灢1 3 )暋 暋由于傅里叶变换可以通过离

28、散傅里叶变换作近似计算,通过离散傅里叶变换求解卷积形式的菲涅耳衍射积分时,理论上使用式( 3 灢 2 灢 1 2 )或式( 3 灢 2 灢 1 3 )表示传递函数是等价的.然而,为获得满足取样定理的计算结果,解析形式的式( 3 灢 2 灢 1 3 )能够用较少的取样数得到足够准确的解.因此,在进行菲涅耳衍射积分的卷积计算时,实际上均使用解析形式的式( 3 灢 2 灢 1 2 )作为传递函数.关于使用式( 3 灢 2 灢 1 2 )计算时遇到的问题我们将在后面讨论基尔霍夫衍射积分及瑞利索末菲衍射积分的计算时再作分析.这里,仅就解析形式的传递函数进行研究.对式( 3 灢 2 灢 1 1 )两边作傅

29、里叶逆变换得到用傅里叶变换及逆变换表述的菲涅耳衍射表达式Ux,( )y=F-1FU0x0 ,y( ) 0HFfx,f( ) y( 3灢2灢1 4 )暋 暋可以看出,菲涅耳衍射过程相当于将物面光波场通过一个线性空间不变系统的过程,观测25信息光学教程平面的光波场的频谱是物平面光波场的频谱与一个菲涅耳衍射传递函数HF的乘积.设衍射场的计算宽度是L0 ,取样间隔为殼x0 ,根据离散傅里叶变换与傅里叶变换关系的讨论,物函数的F F T完成后,其取值范围是1 / 殼x0 =N/L0.为实现在同一坐标尺度下与传递函数的乘积运算,传递函数在频域的取样单位必须满足殼fx=殼fy=1 /L0.于是,当乘积运算完

30、成并进行快速逆傅里叶变换I F F T回到空域时,空域宽度还原为L=1 / 殼fx=L0.因此,利用传递函数法计算衍射时,物平面及衍射观测平面保持相同的取样宽度.此外,由于传递函数不改变物函数的频谱宽度,当物平面的取样满足奈奎斯特取样定理时,计算结果也必然是满足奈奎斯特取样定理的衍射场.现在考虑如何让离散取样满足取样定理的问题.由于式( 3 灢 2 灢 1 4 )中逆变换函数由物函数的频谱与菲涅耳衍射传递函数HFfx,f( )y的乘积组成.HFfx,f( )y及其原函数分别是在整个频域及空域都有取值的函数,基于离散傅里叶变换与傅里叶变换关系的讨论,无论物函数是否是带限函数,衍射运算结果都是非带

31、限函数.因此,要使式( 3 灢 2 灢 1 4 )的离散计算严格地满足奈奎斯特取样定理是不可能的.在实际衍射计算中,可用下面的方法来确定近似满足取样定理的条件.由于菲涅耳衍射传递函数是解析函数,对于给定频率值能准确地得到函数值,它与FU0x0 ,y( ) 0的乘积不改变FU0x0 ,y( ) 0频谱的宽度.因此,只要U0x0 ,y( )0的取样满足取样定理,便能让衍射的D 灢 F F T计算满足取样定理.由于U0x0 ,y( )0的频谱在计算前未知,现对如何正确取样进行讨论.设衍射场总能量为E,根据能量守恒定理及3 . 1节讨论的离散傅里叶变换与傅里叶变换的量值关系,若取样数为N,计算宽度为L

32、0 ,能量值为连续变换的N4 /L40倍.由于离散傅里叶变换后在频域的取样点对应面积为1 /L20 ,离散变换后频域的取样点p/L0 ,q/L( )0对应的能量是D F TU0mL0N,nL0 Np,( )q2L0N4 1L20 ,正确的U0x0 ,y( )0离散傅里叶变换计算必然满足E=L20N4 暺N/ 2-1p=-N/ 2暺N/ 2-1q=-N/ 2D F TU0mL0N,nL0 Np,( )q2曋 c o n s t a n t ( 3灢2灢1 5 )当U0x0 ,y( )0取样合适时,增加取样数将不改变总能量E.因此,在应用研究中可以首先按需要的空间分辨率给定某取样数N,并利用式(

33、3灢2灢1 5 )计算E.此后,在同一物理尺度下,将取样数减小为N/ 2或增加为2N,再计算总能量.若计算结果无本质区别,则认为取样数N满足要求.以下通过理论计算与实验测量的比较证明上述结论.3.2.4暋使用解析形式的传递函数计算衍射的实验证明沿用图3 灢 2 灢 1的衍射实验系统,选择宽度L0 = 1 0 m m的包含透光孔的方形区域为取样区,分别用N= 1 2 8 , 2 5 6 , 5 1 2对U0 (x0 ,y0 )取样、对U0的频谱强度用0 2 5 5的灰度归一化.三种取样计算获得的频谱图像如图3 灢 2 灢 3所示.由图可见,频谱主要能量均集中于中部.此外,根据计算得三幅图像的总能

34、量比为6 4 6 9 暶 6 4 6 6 暶 6 4 6 7 .按照上面对式( 3 灢 2 灢 1 5 )的讨论,三种取样均能较好地满足取样条件.图3 灢 2 灢 4 ( a )给出d= 1 8 4 m m时C C D探测的图像.图3 灢 2 灢 4 ( b ) 、图3 灢 2 灢 4 ( c )及图3 灢 2 灢 4 ( d )依次给出N= 1 2 8 , 2 5 6 , 5 1 2时数值计算的光波场强度图像.35第3章暋衍射积分的快速傅里叶变换计算图3 灢 2 灢 3 暋不同取样数对应的初始场频谱强度图像(N/ 1 0 m m - 1 暳N/ 1 0 m m - 1 )图3 灢 2 灢 4

35、 暋衍射距离d= 1 8 4 m m时实验测量图像与不同取样数的衍射计算图像的比较比较图3 灢 2 灢 4中各图像可以看出,理论计算与实验测量吻合较好.但对于取样数N= 1 2 8的情况,因分辨率相对较低,图像的细节 衍射斑内部衍射条纹变得模糊.因此,当采用解析形式的传递函数计算菲涅耳衍射积分时,通常可以参照分辨率的需要确定取样数.但是,按照衍射的角谱理论 1 , 7 ,光波的衍射是光波场各角谱衍射的叠加.随着衍射距离的增加,衍射场的范围将线性展宽.由于D 灢 F F T算法中物平面及衍射观测平面保持相同的取样宽度,当衍射距离较大时, D 灢 F F T算法不能完整地给出衍射场.因此,这种算法

36、主要用于物光场的高频角谱分量较小及衍射距离较短的情况.附录A中用M A T L A B语言编写的L I M 2 灡 m程序较好地体现了D 灢 F F T计算衍射的功能,读者可以利用该程序验证上述结论,并利用该程序解决许多实际的衍射问题.回顾菲涅耳衍射的S 灢 F F T算法可知,距离较长的菲涅耳衍射计算问题可以使用S 灢 F F T解决.因此,根据实际情况对S 灢 F F T及D 灢 F F T这两种算法作合理选择,才是有效求得衍射结果的途径.45信息光学教程* 3 . 3 暋菲涅耳衍射变换及其在二元光学设计中的应用作为衍射运算的一个重要应用,以下基于二元光学设计的G e r c h b e

37、r g 灢 S a x t o n ( G S )算法 9 - 1 1 ,给出可以用于产品表面进行激光标记的衍射元件设计实例.3.3.1暋菲涅耳衍射变换根据式( 3 灢 2 灢 1 3 )及式( 3 灢 2 灢 1 4 ) ,可以将菲涅耳衍射积分表示为Ux,( )y=F-1FU0x0 ,y( ) 0 e x p jkd1-毸22f2x+f2( ) y( 3灢3灢1 a )对式( 3 灢 3 灢 1 a )进行一次正变换及一次逆变换后容易得到U0x0 ,y( )0=F-1FUx,( ) ye x p-jkd1-毸22f2x+f2( ) y( 3灢3灢1 b )暋 暋将式( 3 灢 3 灢 1 a

38、 )表征的菲涅耳衍射定义为菲涅耳衍射正变换,用符号F(d) 表示.相应地,式( 3 灢 3 灢 1 b )定义为菲涅耳衍射逆变换,简写为F-1(d) ,每个符号的下标d为菲涅耳衍射距离.于是菲涅耳衍射变换对可以简单地表示为Ux,( )y=F(d)U0x0 ,y( ) 0 ( 3灢3灢2 a )U0x0 ,y( )0=F-1(d)Ux,( ) y( 3灢3灢2 b )暋 暋由于上述正变换及逆变换均表示为卷积形式,按照上面D 灢 F F T算法的讨论,作适当的离散处理后就能通过F F T进行卷积的快速计算求解.并且,为得到逆变换满足取样定理的条件,只要将正变换条件的d换为-d就可以了.3.3.2暋

39、二元光学元件图3 灢 3 灢 1 暋二元光学元件及坐标定义在光束整形的应用研究中,通常期望光波通过二元光学系统后成为一个给定强度分布并沿某预定方向传播的平行光.这种既变换振幅又变换波面的二元光学系统通常可以由两个二元光学元件组成 9 - 1 1 ,图3 灢 3 灢 1为所研究问题的示意图.图中,平面x0y0上的第一个元件作振幅变换,使到达xy平面的第二个元件表面的光波强度分布满足设计要求;第二个元件则作波面整形,使透过元件的光波成为沿光轴传播的平面波.以下通过菲涅耳衍射变换讨论光学元件的设计问题.设第一个元件的复振幅变换函数为T0x0 ,y( )0 ,复振幅为Uix0 ,y( )0的光束自左向

40、右传播,期望通过距离d的传播到达xy平面时形成复振幅为Ux,( )y的光波场.根据图3 灢 3 灢 1有Ux,( )y=F(d)Uix0 ,y( )0T0x0 ,y( ) 0 ( 3灢3灢3 )暋 暋对式( 3 灢 3 灢 3 )两边作菲涅耳衍射逆变换后容易得到T0x0 ,y( )0=F-1(d)Ux,( ) yUix0 ,y( )0 ( 3灢3灢4 )令TA(x0 ,y0 ) ,p0 (x0 ,y0 )及e x p j毤(x0 ,y0 ) 分别为二元光学元件的振幅透过率、光瞳及位相变换因子,可将该元件的复振幅透过函数表示为55第3章暋衍射积分的快速傅里叶变换计算T0x0 ,y( )0=TAx

41、0 ,y( )0p0x0 ,y( )0 e x p j毤x0 ,y( ) 0 ( 3灢3灢5 )于是得p0x0 ,y( )0 e x p j毤x0 ,y( ) 0=F-1(d)Ux,( ) yTAx0 ,y( )0Uix0 ,y( )0 ( 3灢3灢6 )暋 暋由于光瞳内p0 (x0 ,y0 ) = 1是实函数,理想的纯位相型元件应满足TAx0 ,y( )0=1 ( 3灢3灢7 )给定入射到元件表面的光波场Uix0 ,y( )0及期望通过光学系统形成的光波场强度分布I(x,y)后,二元光学元件的设计主要任务是求出满足制作工艺要求的位相变换因子,由于第二个光学元件的设计比较简单,以下主要对第一个

42、元件的设计进行讨论.3.3.3暋二元光学元件设计的Gerchberg灢Saxton(GS)算法基于G S算法,现提出一种用菲涅耳衍射变换进行上述衍射元件设计的方法.( 1 )令Q(x,y)为0 2 毿满足给定约束条件的随机数,观测平面的初始振幅可设为U1x,( )y=Ux,( )ye x p jQx,( )( )y( 3灢3灢8 )式中,Ux,( )y=Ix,( )y,约束条件是Q(x,y)所确定的波面法线方向是来自第一个元件光瞳并指向I(x,y)的非零区域的方向.( 2 )二元光学元件的理论尝试解为氞T0 1x0 ,y( )0=F-1(d)U1x,( ) yUix0 ,y( )0 ( 3灢3

43、灢9 )暋 暋 ( 3 )对上面得到的相位分布作量化处理,得到符合光刻要求的尝试解T0 1x0 ,y( )0T0 1x0 ,y( )0=氞T0 1x0 ,y( )0 ( 3灢3灢1 0 )暋 暋例如,利用二值化掩膜处理工艺的尝试解的幅角可按下式作量化:a r gT0 1x0 ,y( ) 0=I N T 2La r g氞T0 1x0 ,y( ) 02 毿2 毿2L( 3灢3灢1 1 )式中,L为正整数, I N T 表示对内的数据进行取整操作.当设计完成后,L次掩膜融刻处理便能生成具有2L级不同位相调制的衍射元件 9 .( 4 )将尝试解归一化,重新表示出观测平面的复振幅U曚1x,( )y=F(

44、d)Uix0 ,y( )0T0 1x0 ,y( )0T0 1x0 ,y( ) 0( 3灢3灢1 2 )暋 暋 ( 5 )式( 3 灢 3 灢 1 2 )规一化并将观测平面复振幅重新设为U1x,( )y=Ux,( )yU曚1x,( )yU曚1x,( )y( 3灢3灢1 3 )暋 暋 ( 6 )将以上结果作为新的迭代计算初始值,反复进行( 2 ) ( 5 )的操作,直到获得满足误差要求或达到设定迭代次数的复振幅变化函数T0 1x0 ,y( )0 .如果应用研究中只需要在xy平面形成期待的强度分布,可以不使用第二个光学元件.反之,如果期望经过xy平面后的光波变为具有期待强度分布的平面波,在上面的设计

45、已经达到要求后,将最后一次迭代时到达观测平面的复振幅U曚1 (x,y)作类同于式( 3 灢 3 灢 1 1 )的量化处理,第二个光学元件的复振幅透过函数应满足65信息光学教程T1x,( )y=e x p j毤1x,( ) y( 3灢3灢1 4 )毤1x,( )y=-I N T 2La r gU曚1x,( ) y2 毿 2 毿2L( 3灢3灢1 5 )暋 暋至此,基本完成了将光束强度分布进行变换并准直的二元光学设计.3.3.4暋二元光学标记元件设计实例在激光对材料表面改性处理应用研究中,利用激光在材料表面烧融成特殊的图案或文字的激光标记技术获得了重要应用.只要将给定衍射距离的光波场强度分布设计成

46、与被标记的图案相对应的形式,二元光学技术就可以方便地用于激光标记元件的设计,获得能量利用率高、标记时间短、标记图案丰富多彩的标记元件.以下,以波长毸= 1 0 . 6 毺 m的千瓦级C O 2激光及2 0 0 8年北京奥运会一吉祥物作为标记图案图3 灢 3 灢 2 ( a ) 为例,给出二元光学元件的一个设计实例 7 .图3 灢 3 灢 2 暋标记图案( a )及标记光束强度图案( b )比较(图面尺寸3 0 m m 暳 3 0 m m )到达二元光学元件的C O 2激光强度分布如图3 灢 3 灢 3 ( a )所示,图3 灢 3 灢 3 ( b )给出经过1 0次迭代运算后在平面xy上光斑的

47、强度分布.图3 灢 3 灢 3 暋入射光束( a )及变换后光束( b )的强度分布比较(图面尺寸3 0 m m 暳 3 0 m m )为简明地显示变换后光束在材料表面的标记效果,令材料表面对作用光束强度分布响应的阈值为变换后光束强度分布极大值的5 % ,即大于该阈值的光束作用于材料表面后,将在材料表面留下热化学作用的印迹.图3 灢 3 灢 2 ( b )给出大于阈值后的光束强度分布图案.与原图比较容易看出, “标记图案暠是原图较忠实的复现.75第3章暋衍射积分的快速傅里叶变换计算* 3 . 4 暋经典衍射公式及其快速傅里叶变换计算根据标量衍射理论,基尔霍夫公式、瑞利索末菲公式及衍射的角谱传播

48、公式是亥姆霍兹方程的准确解 1 , 7 .本节将衍射问题视为光波场通过线性空间不变系统的变换问题,对经典衍射公式:基尔霍夫公式、瑞利索末菲公式、衍射的角谱传播公式及菲涅耳衍射积分的计算进行综合研究,导出每一个公式的衍射传递函数及让计算能够较好地满足取样定理的条件 2 .最后,通过衍射计算实例验证所得结论.3.4.1暋基尔霍夫公式及瑞利索末菲公式的卷积形式与角谱衍射公式一样,基尔霍夫公式与瑞利索末菲公式的物理意义仍然是用某空间平面的光波复振幅U0x0 ,y( )0表示沿光传播方向上经距离d传播后在新的空间平面上的光波复振幅Ux,( )y.图3 灢 4 灢 1给出衍射计算的初始平面x0y0与观测平面xy的关系.图3 灢 4 灢 1 暋衍射计算的初始面与观测面的关系根据式( 2 灢 2 灢 2 8 ) ,观测平面上的光波复振幅可由基尔霍夫公式表示为 1 U(x,y)=1j毸犽毑0U0x0 ,y( )0 e x p j( )krr暳c o s毴+12 dx0 dy0 ( 3灢4灢1 )式中, j=-1 ,毸为光波长;k=2 毿 /毸;r为初始平面任意点(x0 ,y0 )到观察点(x,y)的距离,r表示(x0 ,y0 )到(x,y)的矢径;毴表示r和z轴的夹角.在同一坐标定义下,