1.   引　言

针对光场强度二阶关联成像受大气湍流影响的研究, 2009年华南科技大学Cheng曾在文献[1]中提出了二阶关联成像在大气湍流中成像的理论和解析式. 理论结果表明, 关联成像视为原始物体与光学点扩散函数卷积, 湍流强度增大、传播距离增加和光源尺寸降低均将使得点扩散函数变宽, 导致图像质量下降. 随后华侨大学课题组^[2]推导了任意尺寸和相干度的关联成像与大气湍流相互作用理论, 与文献[1]的结论一致, 关联成像的像质将因湍流强度的增大而受到影响. 2010年, 中国科学院上海光学精密机械研究所研究团队在文献[3]中报道了利用随机相位屏数值仿真湍流下关联成像, 结论是相比于传统成像, 关联成像能够一定程度抵御湍流的影响, 获得相对较高的成像分辨率. 至此, 学术界达成的共识是: 关联成像在大气湍流中成像质量受到湍流影响, 只是影响程度较弱. 然而, 2011年美国陆军实验室Meyers等^[4,5]连续发表两篇文章均指出: 关联成像不受大气湍流影响, 将其称为“turbulence-free ghost imaging”, 并且利用理论和实验验证了该现象. 这一惊人的反直觉结论一出, 立即引起学术界争论. 2012年1月, 美国MIT教授Shapiro^[6]指出关联成像并非不受大气湍流影响, 文献[4, 5]给出的抗湍流成像具有适用条件限制. 文献[6]认为, Meyers等^[4]的实验因为选取了特殊物理参数才导致关联成像对大气湍流的免疫, 主要原因是: 照明光场在物体表面横向相干尺寸$l_{\rm c}$小于湍流的相干尺寸$\rho$, 仅在这个区域内关联成像对湍流免疫. 尽管文献[4, 6]的观点存在巨大差异, 关联成像并不是所有情况下都不受大气湍流影响, 但此争论的结果却给抗大气湍流关联成像研究提出了解决思路和具体实施方案. 本文通过对比研究文献[1—6], 总结出以下关联成像方案不受大气湍流影响需满足的条件: 1)照明光场横向相干尺寸$l_{\rm c}$需小于链路湍流平均相干尺寸$\rho$; 2)关联成像方案中的参考光场传输路径和通过物体的光场传输路径需经历同样的湍流扰动; 3)关联成像分辨率取决于照射在物体上的光斑平均相干尺寸$l_{\rm c}$, 远小于湍流尺寸$\rho$; 4)关联成像信噪比与测量次数密切相关, 理论测量次数为无穷大, 实际应用中测量次数应远大于待重建物体的有效像素数.

尽管关于关联成像能够抗湍流已有明确结论, 但其在实际应用中仍然受到限制. 受限的因素主要包括如下三方面. 1)赝热光源的亮度受限, 无法实现远距离成像. 在光学领域, 主要采用的方案仍是赝热光方案, 对强激光进行随机调制, 受激光亮度影响无法做到较远距离的大视场成像. 例如: 文献[7—9]等报道的成像距离最远在千米级. 2)高亮度赝热光源调制速度受限, 无法实现快速成像. 文献[7, 8]中调制速度在2 kHz左右, 若采用数字微镜器件调制技术, 速度可提升至3.2 kHz, 但这种调制方式对光能损失较大, 器件本身损伤阈值低, 仍无法应用于远距离照明光源. 3)太阳光作为真热光源同时满足上述两个条件, 高亮度和高的涨落速度($10^{10}$ Hz), 但现有探测器无法在一个相干时间皮秒量级($10^{-12}$ s)内完成对光场的连续直接记录, 原有理论框架下, 无法利用现有技术将太阳光关联成像实用化.

2012年, 美国罗彻斯特大学Zerom等^[10]发表的研究成果使得关联成像应用出现转机. 文章指出了关联成像中不必只利用光场在同一相干时间内的信息, 利用光场在一定时间内的强度涨落累加, 同时保证空间上光场的关联, 同样可实现关联成像, 文章对该结论进行了详细理论推导和实验验证. 简言之, 利用响应速度比光场涨落慢的探测器探测光场信息, 不影响关联成像结果. 这使得同时具有高亮度、高速涨落和完全非相干的太阳光场用于实现关联成像成为可能.

出于上述问题, 本文提出一种具体的消除大气湍流的自关联成像理论和实验, 并验证了在远距离条件下水平大气消湍流成像的可行性, 本文提出的点对点自关联成像算法建立在关联成像理论基础上, 但算法更加简单易实现, 成像质量相比于传统方法提升明显. 本文主要内容分为四部分: 1)介绍实验原理及系统; 2)介绍相关理论对比及分析; 3)介绍实验结果及定性分析; 4)定义了图像质量评价方法, 对图像质量进行了定量分析.

2.   实验系统

2.1   实验方案

经典关联成像方案有两种: 第一种方案, 光路需要利用分束器(beam spliter, BS)分束, 在理论模型中考虑两路光束的独立传播, 经分束的光场, 一束光直接被探测器探测, 称为参考探测器, 另一束光经由物体透射或反射被桶探测器收集, 桶探测器只收集光强, 无空间分辨能力; 第二种方案, 物体在光源与分束器之间, 通过物体透射或反射的光场被BS分成两束, 两束光均包含物体信息, 但仍被参考探测器和桶探测器探测. 对于第二种探测方案的研究, 也被称为光场的自相关(或自关联)成像.

本文采用第二种实验方案, 但也略有不同, 实验原理如图1所示. 本文只采用一个参考探测器探测目标, 传统方法对参考探测器求和作为桶探测器信号强度, 本文则直接利用探测器各点与自身进行关联运算. 表面上看, 形式上等同于经典成像, 处理方法上采用了关联算法. Goodman^[11]曾研究了光场高阶关联和光场通过随机非均匀介质时的成像, 与其研究的切入点不同, 本文在二阶关联成像理论框架下研究太阳光强度涨落二阶自关联成像受大气湍流影响模型, 对比了理论上细微差别导致的实验结果差别, 并从实用性角度给出了明确的物理约束条件和实现方法.

图1中, 太阳光场经过大气湍流照射到物体表面, 物体反射光经过长度L (m)的大气湍流被望远镜接收进行关联成像探测; 望远镜系统入射光束经50 : 50的BS分束, 透射光被黑白电荷耦合相机(CCD, MER-030-120UC)记录, 反射光束由光电探测器(PD)进行强度探测. 实际实验中未采用单独的桶探测器, 由相机像素值求和代替桶探测器. 望远镜的口径为280 mm, 焦距F = 2800 mm. 成像相机最大帧频120帧/秒 (frame per second, FPS), 分辨率$656\times492$像素, 像元尺寸5.6 μm, 感光芯片1/3英寸. 在相机与望远镜之间采用了焦距$f$ = 50 mm的镜头, 镜头内置光阑可限制通光孔径和光强. 实验中, 相机的曝光时间$t =$100 μs, 相机与控制单元之间图像传输$50\ \ \rm{FPS}$, 相机的增益在不同实验条件下, 根据光强进行调节后固定.

2.2   实验环境

大气湍流随季节和地面温度变化较大, 因此首先对开展实验的时间和地点进行说明. 实验开展时间为2017年11 月16日上午, 空气轻度污染, 实验系统距离地面高度约20 m, 目标物体为北京中央电视塔, 利用商用地图GPS定位, 可获得目标与实验地点距离为$L = L_{2} = 17.2$ km (图2), GPS误差最大不超过$\pm$ 20 m, 相比于实验距离L可忽略. 根据公开资料显示, 北京电视塔高度386.5 m, 成像实验的区域为塔尖最高段, 按文献[12]中给出的参数, 塔尖最后三节自下而上直径分别为2, 1和0.75 m. 由于大气湍流将随地面高度H和温度T而变化, 可按实验高度H和温度参数T查看大致的大气湍流结构常数的范围, 根据实验数据20—386.5 m高度, 温度$T = [-4\;^{\circ}{\rm C}—8\;^{\circ}{\rm C}]$, 大气湍流强度系数数据的变化范围在文献[13]中可查.

3.   理论模型及分析

3.1   太阳光关联成像理论

如前文所述, 文献[10]利用理论和实验得出了重要的结论: 对于热光关联成像, 利用光场在一定时间内的强度涨落叠加, 同时保证空间上光场的关联, 同样可实现关联成像. 简言之, 利用响应速度比光强涨落慢的探测器探测光强信息, 不影响关联成像结果, 只要光场强度涨落高于探测器散粒噪声水平即可. 文献[10]还指出, 在同一相干时间内的光场对比度与关联成像对比度无关, 此时的关联成像对比度受光场的峰度$\gamma$影响. 传统二阶强度涨落关联鬼成像函数$G({x})$可写为

$$G({x}) \!=\! \frac{1}{N}\!\sum\limits^{N}_{n = 1}I^{(n)}_{\rm B}({x})I^{(n)}({x})\!-\!\frac{1}{N^{2}}\!\sum\limits^{N}_{n = 1}I^{(n)}_{\rm B}\sum\limits^{N}_{n = 1}I^{(n)}({x}),$$

(1)

式中,

$$I^{(n)}_{\rm B} = \sum\limits_{x}I^{(n)}({x})T({x}),$$

(2)

其中$T({x})$为待成像物体反射率函数, $I^{(n)}({x})$是照明光场, N为测量次数. 在利用太阳光场高阶关联成像时, 由于探测器响应速度远比太阳光场涨落慢, CCD曝光时间$\tau = 100$ μs, 太阳光场涨落时间$\tau_{0} = 100$ ps, $\tau\gg\tau_{0}$. 因此CCD相机获得的光场需描述为

$$I^{(n)}_{{\rm B}, M} = \sum\limits_{m = 1}^{M}I_{\rm B}^{(n, m)}$$

(3)

和

$$I^{(n)}_{M}({x}) = \sum\limits_{m}I^{(n, m)}({x}),$$

(4)

其中$I^{(n)}_{M}$下标$M$表示光强度在探测时间内的$M$次叠加. 根据赝热光理论, 光强的对比度定义为$V =$$\sigma(I({x}))/\left\langle I({x})\right\rangle$, 其中强度方差定义为$\sigma(I({x})) =$$\left\langle I^{2}({x})\right\rangle-\left\langle I({x})\right\rangle^{2}$, $\left\langle\cdot\right\rangle$代表系综平均. 结合(4)式可得到$M$次叠加后光强对比度为

$$\begin{split}V_{M} & = \sigma(I_{M}({x}))/\left\langle I_{M}({x})\right\rangle \\ & = \frac{1}{\sqrt{M}}\left\langle I_{M}({x})\right\rangle/\left\langle I_{M}({x})\right\rangle = \frac{1}{\sqrt{M}}. \end{split}$$

(5)

对于热光场, 电磁场按高斯随机分布, 光强按$\rm e$的负指数规律分布, 满足条件$\sigma(I_{M}({x})) = M \left\langle I\right\rangle^{2} =$$\left\langle I_{M}\right\rangle^{2}/M$, 该条件是(5)式中等号成立的条件. 由(5)式可见, 对于高速涨落的热光, $M\rightarrow \infty$, 对比度$V_{M}\rightarrow 0$, 按传统理论$G_{M}({x}) = 0$, 无法实现关联成像. 然而文献[10]证明关联成像在上述情况下仍可成像, 且成像对比度与所测量光场多次叠加后的对比度$V_{M}$无关, 而是与光场峰度$\gamma$有关. 峰度也称为四阶标准矩, 其定义为

$$\gamma_{M} = \left\langle (I_{M}({x})-\left\langle I_{M}({x})\right\rangle)^{4}\right\rangle/\sigma^{4}(I_{M}({x}))-3.$$

(6)

最终, 关联成像在慢探测器条件下的对比度噪声比(contrast to noise ratio, CNR)满足如下条件^[10]:

$${ CNR} = \left[\frac{N-1}{(2T-2+3/N)+(1-1/N)(\gamma_{M}+3)}\right]^{\frac{1}{2}},$$

(7)

式中T代表图像透光面积与光场平均相干面积比值. (7)式说明基于平均光场关联成像的对比度噪声比与图像有效面积T、测量次数N和光强峰度值$\gamma_{M}$有关. (7)式还表明, 利用慢探测器探测太阳光也可实现关联成像.

3.2   消大气湍流关联成像理论

根据3.1节结论, 太阳光下消大气湍流关联成像具有可行性. 关于湍流成像的讨论^[1-6,11,14], 前文提到, 文献[4, 5]提出的完全消大气湍流关联成像理论及实验, 一度在学术界引发争议, 因为在2010 年之前的理论表明关联成像也受大气湍流的影响, 对此, 文献[6]还针对文献[4]提出修正理论. 本文研究了大气湍流对关联成像的影响, 为消湍流成像提供理论支撑. 光场从太阳到目标物体, 再由目标物体到CCD相机, 其光学传递过程可用下式描述：

$$E_{\rm {ccd}}(x_{\rm c}) = \iint {\rm d}x{\rm d}\xi E_{\rm s}(x)h_{1}(\xi, x)T(x)h_{2}(x, x_{\rm c}),$$

(8)

其中$E_{\rm {ccd}}(x_{c})$为CCD表面光场分布, $T(x)$为物体表面反射率函数, $h_{1}(\xi, x)$和$h_{2}(x, x_{\rm c})$分别是从光源到目标和从目标到探测器CCD的光学点扩散函数. 根据扩展的Huygens-Fresnel积分, $h_{1}(\xi, x)$和$h_{2}(x, x_{\rm c})$有如下形式:

$$h_{1}(\xi, x) = \frac{1}{\sqrt{{\rm j}\lambda L_{1}}} {\rm e}^{\frac{{\rm j}{\text{π}}}{\lambda L_{1}}(x-\xi)^{2}+\phi_{1}(\xi, x)}$$

(9)

和

$$h_{2}(x_{c}, \xi) = \frac{1}{{\rm j}\lambda\sqrt{F L_{2}}}\int_{A_{\rm c}} {\rm d}\eta {\rm e}^{-\frac{{\rm j}{\text{π}}}{\lambda L_{2}}(\xi-x_{\rm c}/\beta)\eta+\phi_{2}(\xi, \eta)},$$

(10)

其中, $L_{1}$是太阳到目标的距离; $\phi_{1}(\xi, x)$和$\phi_{2}(\xi, \eta)$分别对应光源太阳到目标物体、目标物体到探测器CCD之间大气湍流引起的光场相位扰动; $\beta = F/L_{2}$是图像的放大率; $A_{\rm c}$为望远镜口径. 由于望远镜焦距$F\gg z_{1}$且$F\gg f$, (10)式中忽略了$z_{1}$和相机镜头焦距$f$.

参考二阶强度涨落关联公式, 本文提出利用光场空间强度自相关算法, 利用空间各点自关联, 有如下公式:

$$\begin{split} & G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c}) \\ = & \left\langle I_{M}(x_{\rm c})I_{M}(x'_{\rm c})\right\rangle-\left\langle I_{M}(x_{\rm c})\right\rangle\left\langle I_{M}(x'_{\rm c})\right\rangle \\ =& \int {\rm d}x_{1}{\rm d}x_{2}{\rm d} x_{3}{\rm d}x_{4}{\rm d}\xi {\rm d}\xi'\left\langle E_{\rm s}(x_{1})E^{\ast}_{\rm s}(x_{2})\right\rangle \\ &\times\left\langle E_{\rm s}(x_{3})E^{\ast}_{\rm s}(x_{4})\right\rangle\left\langle h_{1}(\xi, x_{1})h^{\ast}_{1}(\xi', x_{2})\right\rangle\\ &\times\left\langle h_{1}(\xi, x_{3})h^{\ast}_{1}(\xi', x_{4})\right\rangle \left\langle h_{2}(x_{\rm c}, \xi)h^{\ast}_{2}(x_{\rm c}, \xi')\right\rangle\\ &\times\left\langle h_{2}(x'_{\rm c}, \xi)h^{\ast}_{2}(x'_{\rm c}, \xi')\right\rangle T^{2}(\xi){T^{\ast}}^{2}(\xi'), \end{split}$$

(11)

式中太阳光作为光源, 在可见光波段$\lambda = 0.55$ μm, 太阳在地面相干面积的直径为0.08 mm, 在本实验中可认为是完全非相干光. 考虑到太阳直径, 在源平面处光源为无限大, $\left\langle E_{\rm s}(x_{i})E^{\ast}_{\rm s}(x_{j})\right\rangle =$$I_{0}\delta(x_{i}-x_{j})$, 代入到(11)式并结合(8)—(10)式得

$$\begin{split}& G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c}) \\= &\;\frac{I^{2}_{0}}{\lambda^{3}F L_{2} L_{1}}\int {\rm d}x_{1}{\rm d}x_{2}{\rm d}x_{3}{\rm d}x_{4}{\rm d}\xi {\rm d}\xi'{\rm d}\eta {\rm d}\eta \\ &\times T^{2}(\xi){T^{\ast}}^{2}(\xi') \delta(x_{1}-x_{2})\delta(x_{3}-x_{4}) \\ & \times \big\langle{\rm e}^{\phi_{1}(\xi, x_{3} )+\phi^{\ast}_{1}(\xi', x_{4})}\big\rangle \big\langle {\rm e}^{\phi_{1}(\xi, x_{1} )+\phi^{\ast}_{1}(\xi', x_{2})}\big\rangle \\ &\times {\rm e}^{\frac{{\rm j}{\text{π}}}{\lambda L_{1}}\left[(x_{1}-\xi)^{2}-(x_{2}-\xi')^{2}\right]} \\ &\times {\rm e}^{\frac{{\rm j}{\text{π}}}{\lambda L_{1}}\left[(x_{3}-\xi)^{2}-(x_{4}-\xi')^{2}\right]} \\ &\times {\rm e}^{-\frac{{\rm j}2{\text{π}}}{\lambda L_{2}}\left[(\xi-x_{\rm c}/\beta)\eta-(\xi'-x'_{\rm c}/M)\eta'\right]} \\ &\times \big\langle \rm e^{2\phi_{2}(\xi, \eta)+2\phi^{\ast}_{2}(\xi', \eta')}\big\rangle, \end{split}$$

(12)

式中光场涨落发生在光源面上, 光源涨落和湍流涨落随机过程互不相关, 统计上相互独立. 探测器的像面上光场空间分布受到湍流的调制, 包含在系统传函中. 其中, 湍流项利用文献[15, 16]中的湍流结构畸变函数做正交近似:

$$\begin{split}& \big\langle \rm e^{\phi_{1}(\xi, x_{i})+\phi^{\ast}_{1}(\xi', x_{j})}\big\rangle\\\cong &\; {\rm e}^{-{\frac{1}{\rho^{2}_{1}}\left[(\xi-\xi')^{2}+(\xi-\xi')(x_{i}-x_{j})+(x_{i}-x_{j})^{2}\right]}}\end{split}$$

(13)

和

$$\begin{split}& \big\langle {\rm e}^{\phi_{2}(\xi, \eta)+\phi^{\ast}_{2}(\xi', \eta')}\big \rangle\\\cong & \; \rm e^{-\frac{1}{\rho^{2}_{2}}\left[(\xi-\xi')^{2}+(\xi-\xi')(\eta-\eta')+(\eta-\eta')^{2}\right]} \end{split},$$

(14)

其中$\rho_{1}$是从太阳进入大气层到目标的湍流空间相干尺寸, $\rho_{2}$是目标到探测器路径上湍流空间相干尺寸, 两者均采用了平面波近似, ${(i, j)}\!\subseteq\{(1, 2),$$(3, 4)\}$. 平面波大气湍流相干尺寸

$\rho_{i} = (2.91k^{2}\displaystyle\int^{L}_{0}$${\rm d}z C^{2}_{n}(z)(1-z_{i}/L_{i})^{5/3})^{-3/5}$,

$C^{2}_{n}(z)$是大气湍流强度系数. 将(13)和(14)式代入(12)式, 经过整理得

$$\begin{split} G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c}) & \propto\int {\rm d}\eta \left|T(\xi)\right|^{4}{\rm e}^{-\frac{\rm j2{\text{π}}}{\lambda L}(x_{\rm c}/\beta-x'_{\rm c}/\beta)\eta} \\ & \propto\left|T(x_{\rm c}/\beta)\right|^{4}\delta(x_{\rm c}/\beta-x'_{\rm c}/\beta) \\ & = G^{(2)}_{\rm {ccd}}(x'_{\rm c}, x'_{\rm c}) = \sigma^{2}(I(x'_{\rm c}/\beta)). \end{split}$$

(15)

(15)式表明, 太阳到物体之间的大气湍流相干尺寸$\rho_{1}$和物体到探测器之间的大气湍流$\rho_{2}$对自关联成像均无影响, 因关联效应公式中湍流相关项可完全消除, 即理想情况下自相关成像完全不受大气湍流影响.

值得注意的是, 此处的推导未利用传统关联成像中的桶探测器, 即空间一个点和区域内所有点求和的关联, 而是空间各点和各点自身的关联. 若利用桶探测器和CCD采集图像自关联, 则需修改(11)式对$I_{M}(x_{\rm c})$进行积分:

$$\begin{split}&G^{2}(x'_{\rm c}) = \int {\rm d}x_{\rm c}G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c}) \\ = & \bigg\langle\!\int I_{M}(x_{\rm c}){\rm d}x_{\rm c}I_{M}(x'_{\rm c})\bigg\rangle\!-\!\bigg\langle \!\int I_{M}(x_{\rm c}){\rm d}x_{\rm c}\bigg\rangle\big\langle I_{M}(x'_{\rm c})\big\rangle \\ = &\bigg\langle I_{M}(x'_{\rm c})I_{M}(x'_{\rm c})\bigg\rangle-\bigg\langle I_{M}(x'_{\rm c})\bigg\rangle\big\langle I_{M}(x'_{\rm c})\big\rangle \\ & +\bigg\langle \int_{x_{\rm c}\neq x'_{\rm c}} I_{M}(x_{\rm c}){\rm d}x_{\rm c}I_{M}(x'_{\rm c})\bigg\rangle\\ &-\bigg\langle \int_{x_{\rm c}\neq x'_{\rm c}} I_{M}(x_{\rm c}){\rm d}x_{\rm c}\bigg\rangle\big\langle I_{M}(x'_{\rm c})\big\rangle \\ = &\;G^{(2)}_{\rm {ccd}}(x'_{\rm c}, x'_{\rm c})+G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c})|_{x_{\rm c}\neq x'_{\rm c}}, \end{split}$$

(16)

式中最后一行等式右侧第一项与(15)式得到的结论相同, 第二项为空间点交叉关联项, 在理想情况下$G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c}) = 0$, 理论上(16)式和(15)式一致, 前提是经过无限次测量的系综平均. 实际应用中往往难以满足无限次测量的条件, $G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c})\neq0$, 则(15)式与(16)式不相等.

根据(15)式第一行右侧, 当$x_{\rm c} = x'_{\rm c}$时, 成像结果是狄拉克函数$\delta(x)$与物体$T(x_{\rm c}/\beta)$的卷积, 图像分辨率等于经典成像的分辨率, 即望远镜与相机分辨率综合结果. 而根据(16)式最后一行, 当$x_{\rm c}\neq x'_{\rm c}$时, $G^{(2)}_{\rm {ccd}}(x_{\rm c}, x'_{\rm c})\neq0$, 此时(15)式第一行右侧不能写为狄拉克函数$\delta(x)$, 此时物体$T(\xi)$和${\rm e}^{-\frac{{\rm j}2{\text{π}}}{\lambda L}(x_{\rm c}/\beta-x'_{\rm c}/\beta)\eta}$指数项的卷积与湍流平均相干尺寸有关, 因此将导致分辨率降低. 类似分辨率的讨论在文献[1—3]中均有详细论述. (15)式与(16)式还表明基于桶探测方案的传统关联成像与本文提到的空间点对点自关联, 即强度方差, 成像结果将有所不同, 相比于传统强度关联成像, 本文提出的空间强度点对点自关联在少量测量情况下, 可很好地消除大气湍流, 获得高质量图像. 相比于文献[17]报道的其他关联成像方案, 例如单像素成像^{[18, 19]}, 本文方案真正利用了光场关联特性, 在消除大气湍流方面更具优势.

4.   实验结果

在实验中, 采用的相机工作在可见光波段, 实验过程中观测连续拍摄的视频可明显看出大气湍流对目标图像的影响. 尽管CCD曝光时间$\tau =$100 μs, 由于成像距离$L_{2} = 17.2$ km, 且距离地面较近, 在每次测量的曝光时间内除太阳光场强度涨落进行了多次累加, 大气湍流引起的光场变化也进行了累加. 根据本文2.1节中的实验系统参数, 光学系统的缩放倍数$V = -F/L = 1/6143$, 图像像素数为$656\times492$, 像元尺寸5.6 μm, 塔尖最后三节自下而上直径分别为2, 1和 0.75 m, 其中最后一节0.75 m, 在像面线度122 μm, 对应20 个像素左右, $\lambda = 0.550$ μm时, 光学系统衍射极限$1.22 \lambda F/D =$6.7 μm $>$ 5.6 μm, 相机的分辨率高于光学系统衍射极限, 故CCD相机的物方分辨率极限0.04 m. 图3中, 四幅图像来源于同一时间段内, 不同时刻连续拍摄, 仔细观察可发现大气湍流影响了图像质量, 包括模糊、信噪比低、对比度低, 且发生了几何畸变. 图3中所有图片的背景出现了不均匀亮度的环形区域, 这是由相机光圈和卡式望远镜特殊结构引起的强度不均匀造成的, 调节相机光圈可降低不均匀性, 但同时导致图像对比度下降, 因此实验保留了图中拍摄效果. 图3所有图像都无法直接看到塔尖最下部分的明暗(红白色)对比, 塔尖上部的明暗(红白色)对比也并不清晰.

利用上述连续的测量数据, 本文测量次数$M = 500$, 并将多次测量的结果进行累加求平均, 得到累加图像, 如图4所示. 图4(a)为某一时刻CCD记录图像, 其中塔尖部分的明暗对比无法直接观测到, 图4(b)为500次测量累加求平均后的图像, 与图4(a)对比可看出累加后效果提升不明显, 多次平均只是导致图像变得平滑, 清晰度并未得到改善. 图4(c)和图4(d)分别是传统有桶探测器的光强涨落二阶关联成像和光强点对点二阶强度自关联成像, 后者即图像强度方差. 比较图4(c)、图4(d)两幅图和图4(a)、图4(b)两幅图, 图像质量均有显著改善, 其中图4(c)中的光圈背景得到了很好地去除, 塔尖整体上突显出来, 图4(d)在所有成像结果中最为清晰, 塔尖最下部分的明暗(红白)对比得到了显著提高, 清晰度也有较大改善, 如明暗(红白色)分界限处. 理论上, 图4(c) 因消除了部分湍流影响, 较直接成像效果略有改善, 但根据(16)式其等号右侧第二项仍为噪声项, 理论上只有当测量次数无穷大, 图4(c)和图4(d)才能趋于一致. 至此, 本文用光强点对点二阶强度自关联实现了成像, 且实验验证了在17.2 km距离的消湍流成像, 实验结果与理论预测的结果一致, 研究成果具有较强实用性.

5.   图像质量评价

由于图像峰值信噪比公式^[20]对图像客观评价具有一定局限性, 为了将图像质量进行量化比较, 定量分析图像质量提升的幅度, 本文定义一种图像的评价方法, 选取图像部分灰度值的方差$S_{i}$作为信号, $S_{i}$一定程度上反映了重建图像灰度偏离平均值的大小, 偏离越大表明对比度越大; 选取图像以外对应的均匀背景部分方差$\sigma_{i}$作为噪声, 图像均匀部分可通过先验知识来选取, 例如天空、白墙成像后的方差可视为噪声, 因为理论上均匀背景的方差应为0, 此时对应信噪比将为无限大. 根据上述分析, 重新定义图像信噪比公式SNR (signal-to-noise ratio, SNR):

$$SNR = 10{\rm {log}}_{10}(S_{i}/\sigma_{i}),$$

(17)

其中$S_{i} = \sum_{x} (T_{\rm {obj}}(x)-\overline{T_{\rm {obj}}})^{2}/P$, $T_{\rm {obj}}(x)$为目标区域灰度值, $\overline{T_{\rm {obj}}}$为目标灰度平均值, $P$为目标区域像素数; $\sigma_{i}$为图像背景部分方差, $\sigma_{i} = \sum_{x} (T_{\rm b}(x)-$$\overline{T_{\rm b}})^{2}/Q$, $T_{\rm b}(x)$为背景区域灰度值, $\overline{T_{\rm b}}$为背景区域灰度平均值, $Q$为背景区域像素数.

根据(17)式, 首先对原始图像进行处理, 将电视塔与天空背景进行分离以用于分别提取各部分方差$S_{i}$和$\sigma_{i}$. 用于分离目标区域与背景区域的掩模图像如图(5)所示, 白色区域为计算目标物体方差$S_{i}$, 黑色区域用于计算$\sigma_{i}$. 分别将图4(a)—(d)结合图5掩模, 按(17)式进行信噪比的计算, 得到的结果为: 图4(a)—(d)的信噪比分别为SNR = 0.60, 0.58, 3.71, 7.11 dB, 其数值结果与人眼判断结果一致. 定量分析的结果表明, 关联成像相比于传统直接成像SNR提升 3—6 dB, 点对点的自关联图像SNR相对于传统“桶”对点自关联成像可提升3.40 dB.

6.   结　论

本文提出了太阳光强度涨落点对点二阶关联消大气湍流成像理论, 在深入分析消大气湍流关联成像理论及相关文献的分岐后, 针对高速涨落的太阳光场是否能用于实现关联成像进行了理论推导和实验测试, 实验上获得了图像信噪比的提升, 与理论预测相符. 在大气湍流影响下, 传统图像去除湍流影响的方法需要复杂的计算^[21], 以及利用图像先验知识来进行图像处理和增强. 本文提出的强度涨落点对点二阶关联算法, 经理论证明和17 km距离外场实验验证, 可显著提高经典成像的图像信噪比, 提升效果优于6 dB. 该算法可直接移植到经典成像领域, 无需改变经典成像系统的光机结构, 直接利用图像数据采用在线或离线方式处理, 达到消大气湍流成像效果, 算法简单, 实用性强.