留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

相位恢复算法:原理、发展与应用(特邀)

王爱业 潘安 马彩文 姚保利

王爱业, 潘安, 马彩文, 姚保利. 相位恢复算法:原理、发展与应用(特邀)[J]. 红外与激光工程, 2022, 51(11): 20220402. doi: 10.3788/IRLA20220402
引用本文: 王爱业, 潘安, 马彩文, 姚保利. 相位恢复算法:原理、发展与应用(特邀)[J]. 红外与激光工程, 2022, 51(11): 20220402. doi: 10.3788/IRLA20220402
Wang Aiye, Pan An, Ma Caiwen, Yao Baoli. Phase retrieval algorithms: principles, developments and applications (invited)[J]. Infrared and Laser Engineering, 2022, 51(11): 20220402. doi: 10.3788/IRLA20220402
Citation: Wang Aiye, Pan An, Ma Caiwen, Yao Baoli. Phase retrieval algorithms: principles, developments and applications (invited)[J]. Infrared and Laser Engineering, 2022, 51(11): 20220402. doi: 10.3788/IRLA20220402

相位恢复算法:原理、发展与应用(特邀)

doi: 10.3788/IRLA20220402
基金项目: 国家自然科学基金(12104500)
详细信息
    作者简介:

    王爱业,男,硕士生,主要从事傅里叶叠层成像技术方面的研究

    潘安,男,副研究员,博士生导师,博士,主要从事计算成像技术方面的研究

    马彩文,男,研究员,博士生导师,博士,主要从事光电测控技术方面的研究

  • 中图分类号: O438

Phase retrieval algorithms: principles, developments and applications (invited)

  • 摘要: 研究表明,由于相位比振幅包含更多关于场的信息,因此相位测量在现代科学和工程的诸多分支中始终是研究的热点问题。在可见的电磁波范围内,相位信息很难通过现有的光电探测器直接采集获取。相位恢复技术提供了一种从捕获的强度信息中将相位信息“计算”出来的有效手段,并已成功应用于天文观测、生物医学成像和数字信号复原等多个科学领域。算法是相位恢复技术的核心,也是该技术发展和应用的关键。文中结合物理学原理和信号处理方法对相位恢复算法的基本原理进行阐述,综述了各类相位恢复算法的发展历程及其优缺点,并简单概述了相位恢复算法在光学领域的典型应用,最终指明其面临的挑战和未来的发展趋势:更优异的收敛性能和噪声鲁棒性、恢复更复杂物体相位信息的能力、多目标多任务集成的兼容性。
  • 图  1  相位恢复算法发展时间线

    Figure  1.  Development timeline of phase retrieval algorithm

    图  2  单强度交替投影迭代算法流程图与原理图。(a) GS算法流程图;(b) HIO算法流程图;(c1)~(c3) ER、HIO、CHIO算法下一次迭代输入与当前迭代输出之间的关系曲线

    Figure  2.  Flow diagrams and schematic diagrams of single-intensity alternating projection iterative algorithms. (a) Flow diagram of GS algorithm; (b) Flow diagram of HIO algorithms; (c1)-(c3) The relationship between the input of the next iteration and the output of the current iteration for ER, HIO, CHIO respectively

    图  3  径向多强度信息相位恢复原理示意图。(a) 叠层成像原理示意图;(b) ePIE算法流程图;(c) 傅里叶叠层显微术的成像原理与系统配置[95]

    Figure  3.  Schematic diagrams of radial multi-intensity phase retrieval. (a) Schematic diagram of ptychography; (b) Flow diagram of ePIE algorithm; (c) Principles and experimental setup of FPM[95]

    图  4  TIE成像系统设置与成像结果[126]。(a)~(b) 基于4f系统和无透镜成像系统的TIE成像配置;(c1)~(c4) 使用TIE方法采集的三张MCF-7细胞离焦强度图像以及最终恢复的相位图像

    Figure  4.  Systematic setup and imaging result of TIE[126]. (a)-(b) TIE systematic setup based on 4f system and lensless imaging; (c1)-(c4) Three defocused intensity images of an MCF-7 cell acquired by TIE and its reconstructed phase image

    图  5  KKSAI的基本原理与可视化流程图

    Figure  5.  Principles and visual flowchart of KKSAI

    图  6  基于优化理论的相位恢复算法的性能比较。(a1)~(a2) PhaseLift和PhaseCut算法在高斯噪声以及图像扫描线存在的情况下相位恢复的相对误差曲线[35];(b1) 使用谱方法和截断谱方法两种初始化方案对一维高斯信号进行相位重构的相对误差曲线[37];(b2) 使用WF和TWF算法对实值信号进行相位恢复的经验成功率曲线[37];(c1)~(c2) 使用基于强度损失函数和基于幅度损失函数的算法对实值和复值信号进行相位恢复的经验成功率曲线[38]

    Figure  6.  Performance comparison of phase retrieval algorithms based on optimization theory. (a1)-(a2) Relative error curves of phase retrieval for PhaseLift and PhaseCut under Gaussian noise and image scan-lines[35]; (b1) Relative error curves of spectral method and truncated spectral method for the phase retrieval of 1D Gaussian signal[37]; (b2) Empirical success rate curves of WF and TWF for the phase retrieval of real-valued signal[37]; (c1)-(c2) Empirical success rate curves of algorithms based on intensity/amplitude loss function for the phase retrieval of real-valued and complex-valued signals[38]

    图  7  GESPAR算法与半正定规划和稀疏性Fienup算法的比较[58]。(a1)~(a2) 对不同矢量长度信号进行相位恢复的算法复杂度对比;(b)相位恢复成功率对比曲线

    Figure  7.  Comparison between GESPAR, SDP and Sparse-Fienup algorithm[58]. (a1)-(a2) Comparison of algorithm complexity for the phase recovery of signals with different vector lengths; (b) Comparison curves of phase recovery successful rate

    图  8  基于深度学习的相位恢复方法。(a) 使用深度学习解决相位恢复问题的物理模型;(b)使用深度学习进行单帧无透镜相位恢复[59];(c)使用深度学习仅采集少量图像实现快速FPM成像[61];(d) prDeep算法中的去噪卷积神经网络的框架 [62]

    Figure  8.  Phase retrieval based on deep learning. (a) Physical model of deep learning solving phase retrieval problems; (b) Single-frame lensless phase retrieval using deep learning[59]; (c) Fast FPM imaging with a few images using deep learning[61]; (d) Basic framework of DnCNN in prDeep algorithm[62]

    图  9  几种典型的相位恢复应用场景。(a) 使用X射线对晶体微观结构进行衍射成像;(b) 采用4f系统配置的基于相位恢复的光学加密系统;(c) 天文观测中的自适应光学系统;(d) 大口径光学镜面的面型检测

    Figure  9.  Several typical application scenarios of phase retrieval. (a) X-ray diffraction imaging of crystal microstructure; (b) Optical encryption system based on phase retrieval with a 4f setup; (c) Adaptive optics system in astronomical observation; (d) Surface shape detection of large-aperture optical mirrors

    表  1  相位恢复算法的分类

    Table  1.   Classification of phase retrieval algorithms

    AlgorithmProsCons
    Alternating projections[29,54] Simple implementation
    Wide applicability
    Time-consuming
    Low precision
    Non-
    prior
    TIE[55] Deterministic solution
    Simple computation
    Require no phase unwrapping and complex optical systems
    Harsh preconditions
    Lower precision compared with interferometry
    Limited applications
    KK relations[56]Reconstruction with high quality and high speed
    Require less data collection
    Require high precision of NA matching
    SDP[34-35] Guarantee for convergence and reconstruction accuracy Not applicable to large-scale signals
    Non-convex
    optimization[36-38]
    Struggle with parameter
    optimization
    Rely on initial estimate
    Sparsity-
    based
    General[57-58] Reconstruction with high quality and stability
    Low computational cost
    Applicable to objects with sparsity only
    Deep learning[59-62] Powerful expressive ability Lack of interpretability
    Unable to deal with multiple models and applications
    下载: 导出CSV

    表  2  相位恢复问题解的唯一性[66]

    Table  2.   Uniqueness of phase retrieval problem[66]

    ItemContent
    Fourier measurements1DNo uniqueness
    ≥2DUniqueness for real nonreducible signals with a twice oversampling rate
    General measurements Real signal Satisfying the complement property is necessary
    2N–1 random measurements guarantee uniqueness with high probability
    Real signal (noisy) Stable uniqueness requires nlogn measurements
    Complex signal 4n–4 generic measurements are sufficient
    下载: 导出CSV

    表  3  Misell算法流程表

    Table  3.   Flow chart for Misell algorithm

    Algorithm 1:Misell
    Input: Intensity measurements in the frequency domain without modulation ${\left| {{F_0}} \right|^2}$ and phase modulation factors ${t_m}$ with their corresponding intensity measurements ${\left| {{F_m}} \right|^2}\left( {1 \leqslant m \leqslant n} \right)$.
    Output:Spatial complex distribution f.
    S1. Select the initial Guess of spatial complex distribution without modulation f0 at random.
    S2. Add factor t1 to the phase of f0 , and obtain the spatial image after the first phase modulation f1.
    S3. Fourier transform f1 to obtain the complex distribution in the frequency domain G1.
    S4. Replace the amplitude of G1 with the measurement ${\left| {{F_1}} \right|^{}}$, and obtain the new complex distribution in the frequency domain $G_1'$.
    S5. Inverse Fourier transform $G_1'$ and substract the phase factor t1 , and obtain the new spatial image without modulation $f_1'$.
    S6. Add the factor $ {t}_{2} $ to the phase of $f_1'$, and obtain the spatial image after the second phase modulation f2 .
    S7. Repeat S2-S6 to cover all modulated images.
    S8. Repeat S2-S7 until the convergence conditions are satisfied, and output the final result.
    下载: 导出CSV

    表  4  KKSAI技术流程表

    Table  4.   Flow chart for KKSAI technique

    Algorithm 2:KKSAI
    Input:Intensity measurements Ii (x, y) at each sub-aperture and reference plane wave ri (x, y). Their definitions are given respectively as:
    $ I_{i}(x, y)=\left|\mathcal{F}^{-1}\left\{S_{i}(u, v)\right\}\right|^{2}$
    $v_{1}(x, y)={\rm{e}}^{-{y\left(\omega_{1}-x+y_{1}-y\right)} }$
    where $ S_{i}(u, v)=S(u, v) D\left(u-u_{i}, v-v_{i}\right)$ denotes the subregion spectrum at the imaging plane cropped by the sub-aperture at position (ui, vi).
    Output: Pupil-limited object function S(u, v).
    S1. For the starting intensity I1, introduce a Hilbert kernel $ H_{1}(u, v)=-j \operatorname{sgn}\left(v_{1}\right) \cdot \operatorname{sgn}(v)$ and define an intermediate variable as:
    $ X_{1}=\ln \left\{\left[\tilde{s}_{1}(x, y)+r_{1}(x, y)\right] / r_{1}(x, y)\right\}$
    where $ \tilde{s}_{1}(x, y)$ denotes the shifted scattered field exiting from the object plane.
    S2. Impose KK relations to calculate the real part and imaginary part of the intermediate variable:
    $\operatorname{Re}\left\{X_{1}\right\}=\dfrac{1}{2} \ln \left[I_{1}(x, y) /\left|r_{1}(x, y)\right|^{2}\right]$
    $ \operatorname{Im}\left\{X_{1}\right\}=\mathcal{F}^{-1}\left\{\mathcal{F}\left\{\operatorname{Re}\left\{X_{1}\right\}\right\} \cdot H_{1}(u, v)\right\}$
    S3. Obtain the shifted form of the subregion:
    $S_{1}\left(u+u_{1}, v+v_{1}\right)=F\left\{{\rm{e}}^{\operatorname{Re}\left\{X_{1}\right\}+j \ln \left\{X_{1}\right\} } \cdot v_{1}(x, y)\right\}$
    move it back to the correct position to get Si (u, v).
    S4. Repeat above operations for other subregions.
    S5. Subsequently stitch all subregional complex fields into a full pupil plane spectrum S(u, v):
    $S(u, v)=\displaystyle\sum_{i=1}^{4} S_{i}(u, v) /\left[\displaystyle\sum_{i=1}^{4} D\left(u-u_{i}, v-v_{i}\right)+\varepsilon\right]$
    where $\varepsilon$ is a small constant to keep numerical stability in the zero-valued region.
    下载: 导出CSV

    表  5  截断维尔丁格流算法流程表

    Table  5.   Flow chart for TWF algorithm

    Algorithm 3:Truncated Wirtinger Flow (TWF)
    Input:Intensity measurements yi and sampling vectors αi; trimming thresholds αy and truncation criteria sets ε1, ε2.
    Output:Iteration value zT.
    S1. Initialize estimated value ${{\textit{z}}_0} = \lambda \tilde {\textit{z}}$,where λ is defined by Eq. (23) and $\tilde {\textit{z}}$ is the leading eigenvector of the following matrix:
    $Y = \dfrac{1}{m}\displaystyle\sum\limits_{i = 1}^m { {y_i} } {a_i} a_i^{\rm{H}}{1_{\left\{ {\left| { {y_i} } \right| \leqslant \alpha _y^2\left( {\dfrac{1}{m}\displaystyle\sum\nolimits_{i = 1}^m { {y_i} } } \right)} \right\} } }$
    S2. Impose iterations on the initial estimated value by gradient update, and the update formula is:
    ${ {\textit{z} }_{k + 1} } = { {\textit{z} }_k} + \dfrac{ {2{\mu _k} } }{m}\displaystyle\sum\limits_{i = 1}^m {\dfrac{ { {y_i} - \left| {a_i^{\rm{H}}{ {\textit{z} }_k} } \right|} }{ {a_i^{\rm{H}}{ {\textit{z} }_k} } } } {\text{ } }{a_i}{1_{\varepsilon _1^i \cap \varepsilon _2^i} }$
    S3. Output the final result after completing iterations of predetermined number.
    下载: 导出CSV

    表  6  运行时间对比表[58]

    Table  6.   Comparison of runtime performance[58]

    SDP/sSparse-Fienup/sGESPAR/s
    s=31.320.090.04
    s=51.780.120.05
    s=83.850.500.06
    下载: 导出CSV
  • [1] Lax M, Louisell W H, Mcknight W B. From Maxwell to paraxial wave optics [J]. Phys Rev A, 1975, 11(4): 1365-1370. doi:  10.1103/PhysRevA.11.1365
    [2] Cowley J M. Diffraction Physics[M]. Amsterdam: Elsevier, 1995.
    [3] Stratton J A. Electromagnetic Theory [M]. New Jersey: John Wiley & Sons, 2007.
    [4] Hecht E. Optics [M]. 4th ed. San Francisco: Addison Wesley, 2001.
    [5] Oppenhein A V, Lim J S. The importance of phase in signals [C]//Proceedings of IEEE, 1981, 69(5): 529-541.
    [6] Giloh H, Sedat J W. Fluorescence microscopy: Reduced photobleaching of rhodamine and fluorescein protein conjugates by n-propyl gallate [J]. Science, 1982, 217(4566): 1252-1255. doi:  10.1126/science.7112126
    [7] Stephens D J, Allan V J. Light microscopy techniques for live cell imaging [J]. Science, 2003, 300(5616): 82-86. doi:  10.1126/science.1082160
    [8] Zernike F. Phase contrast, a new method for the microscopic observation of transparent objects part II [J]. Physica, 1942, 9(10): 974-986. doi:  10.1016/S0031-8914(42)80079-8
    [9] Nomarski G. Differential microinterferometer with polarized waves [J]. J Phys Radium, 1955, 16: 9-13.
    [10] Park Y, Depeursinge C, Popescu G. Quantitative phase imaging in biomedicine [J]. Nat Photonics, 2018, 12(10): 578-589. doi:  10.1038/s41566-018-0253-x
    [11] Jo Y, Cho H, Lee S Y, et al. Quantitative phase imaging and artificial intelligence: A review [J]. IEEE J Sel Top Quantum Electron, 2018, 25(1): 1-14.
    [12] Cuche E, Bevilacqua F, Depeursinge C. Digital holography for quantitative phase-contrast imaging [J]. Opt Lett, 1999, 24(5): 291-293. doi:  10.1364/OL.24.000291
    [13] Schnars U, Jptner W P O. Digital recording and numerical reconstruction of holograms [J]. Meas Sci Technol, 2002, 13(9): R85. doi:  10.1088/0957-0233/13/9/201
    [14] Tahara T, Quan X, Otani R, et al. Digital holography and its multidimensional imaging applications: A review [J]. Microscopy, 2018, 67(2): 55-67. doi:  10.1093/jmicro/dfy007
    [15] Hartmann J. Bemerkungen uber den bau und die justirung von spektrographen [J]. Zt Instrumentenkd, 1990, 20(47): 17-27.
    [16] Shack R V, Platt B C. Production and use of a lecticular hartmann screen [J]. J Opt Soc Am, 1971, 61: 656-661.
    [17] Platt B C, Shack R V. History and principles of shack-hartmann wavefront sensing [J]. J Cataract Refr Surg, 2001, 17(5): S573-577.
    [18] Esposito S, Riccardi A. Pyramid wavefront sensor behavior in partial correction adaptive optic systems [J]. Astron Astrophys, 2001, 369(2): L9-L12. doi:  10.1051/0004-6361:20010219
    [19] Ragazzoni R, Diolaiti E, Vernet E. A pyramid wavefront sensor with no dynamic modulation [J]. Opt Commun, 2002, 208(1-3): 51-60. doi:  10.1016/S0030-4018(02)01580-8
    [20] Neil M A A, Booth M J, Wilson T. New modal wave-front sensor: A theoretical analysis [J]. J Opt Soc Am A, 2000, 17(6): 1098-1107. doi:  10.1364/JOSAA.17.001098
    [21] Booth M J. Wave front sensor-less adaptive optics: A model-based approach using sphere packings [J]. Opt Express, 2006, 14(4): 1339-1352. doi:  10.1364/OE.14.001339
    [22] Booth M J. Adaptive optics in microscopy [J]. Phil Trans R Soc A, 2007, 365: 2829-2843. doi:  10.1098/rsta.2007.0013
    [23] Cha J W, Ballesta J, So P T C. Shack-hartmann wavefront-sensor-based adaptive optics system for multiphoton microscopy [J]. J Biomed Opt, 2010, 15(4): 046022. doi:  10.1117/1.3475954
    [24] Dayton D, Gonglewski J, Pierson B, et al. Atmospheric structure function measurements with a Shack-Hartmann wave-front sensor [J]. Opt Lett, 1992, 17(24): 1737-1739. doi:  10.1364/OL.17.001737
    [25] Ricklin J C, Davidson F M. Atmospheric turbulence effects on a partially coherent Gaussian beam: Implications for free-space laser communication [J]. J Opt Soc Am A, 2002, 19(9): 1794-1802. doi:  10.1364/JOSAA.19.001794
    [26] Liang J, Grimm B, Goelz S, et al. Objective measurement of wave aberrations of the human eye with the use of a hartmann-shack wave-front sensor [J]. J Opt Soc Am A, 1994, 11(7): 1949-1957. doi:  10.1364/JOSAA.11.001949
    [27] Moreno-Barriuso E, Navarro R. Laser ray tracing versus hartmann-shack sensor for measuring optical aberrations in the human eye [J]. J Opt Soc Am A, 2000, 17(6): 974-985. doi:  10.1364/JOSAA.17.000974
    [28] Sayre D. Some implications of a theorem due to Shannon [J]. Acta Crystallogr, 1952, 5(6): 843-843.
    [29] Gerchberg R W. A practical algorithm for the determination of phase from image and diffraction plane pictures [J]. Optik, 1972, 35(2): 237-246.
    [30] Gerchberg R W. Phase determination for image and diffraction plane pictures in the electron microscope [J]. Optik, 1971, 34: 275-284.
    [31] Miao J, Charalambous P, Kirz J, et al. Extending the methodology of X-ray crystallography to allow imaging of micrometre-sized non-crystalline specimens [J]. Nature, 1999, 400(6742): 342-344. doi:  10.1038/22498
    [32] Miao J, Sayre D, Chapman H N, et al. Phase retrieval from the magnitude of the Fourier transforms of nonperiodic objects [J]. J Opt Soc Am A, 1998, 15: 1662. doi:  10.1364/JOSAA.15.001662
    [33] Howells M R, Beetz T, Chapman H N, et al. An assessment of the resolution limitation due to radiation-damage in X-ray diffraction microscopy [J]. J Electron Spectrosc Relat Phenom, 2009, 170(1-3): 4-12. doi:  10.1016/j.elspec.2008.10.008
    [34] Candès E J, Thomas S, Vladislav V. PhaseLift: Exact and stable signal recovery from magnitude measurements via convex programming [J]. Commun Pure Appl Math, 2013, 66(8): 1241-1274. doi:  10.1002/cpa.21432
    [35] Waldspurger I, D’aspremont A, Mallat S. Phase recovery, max-cut and complex semidefinite programming [J]. Math Program, 2015, 149(1-2): 47-81. doi:  10.1007/s10107-013-0738-9
    [36] Candès E J, Li X, Soltanolkotabi M. Phase retrieval via wirtinger flow: Theory and algorithms [J]. IEEE Trans Inf Theory, 2015, 61(4): 1985-2007. doi:  10.1109/TIT.2015.2399924
    [37] Chen Y, Candès E J. Solving random quadratic systems of equations is nearly as easy as solving linear systems [J]. Commun Pure Appl Math, 2017, 70(5): 822-883. doi:  10.1002/cpa.21638
    [38] Wang G, Giannakis G B, Eldar Y C. Solving systems of random quadratic equations via truncated amplitude flow [J]. IEEE Trans Inf Theory, 2018, 64(2): 773-794. doi:  10.1109/TIT.2017.2756858
    [39] Eldar Y C, Kutyniok G. Compressed Sensing: Theory and Applications [M]. Cambridge: Cambridge University Press, 2012.
    [40] Sanz J L C. Mathematical considerations for the problem of Fourier transform phase retrieval from magnitude [J]. SIAM J Appl Math, 1985, 45(4): 651-664. doi:  10.1137/0145038
    [41] Fannjiang A. Absolute uniqueness of phase retrieval with random illumination [J]. Inverse Probl, 2012, 28(7): 075008. doi:  10.1088/0266-5611/28/7/075008
    [42] Hofstetter E. Construction of time-limited functions with specified auto-correlation functions [J]. IEEE Trans Inf Theory, 1964, 10(2): 119-126. doi:  10.1109/TIT.1964.1053648
    [43] Bruck Y M, Sodin L. On the ambiguity of the image reconstruction problem [J]. Opt Commun, 1979, 30(3): 304-308. doi:  10.1016/0030-4018(79)90358-4
    [44] Hayes M. The reconstruction of a multidimensional sequence from the phase or magnitude of its Fourier transform [J]. IEEE Trans Acoust, Speech, Signal Process, 1982, 30(2): 140-154. doi:  10.1109/TASSP.1982.1163863
    [45] Bates R H T. Fourier phase problems are uniquely solvable in more than one dimension. I: Underlying theory [J]. Optik, 1982, 61(3): 247-262.
    [46] Van H, Hayes M, Lim J, et al. Signal reconstruction from signed fourier transform magnitude [J]. IEEE Trans Acoust, Speech, Signal Process, 1983, 31(5): 1286-1293.
    [47] Beinert R. Non-negativity constraints in the one-dimensional discrete-time phase retrieval problem [J]. Information and Inference: A Journal of the IMA, 2017, 6(2): 213-224.
    [48] Beinert R, Plonka G. Ambiguities in one-dimensional discrete phase retrieval from Fourier magnitudes [J]. J Fourier Anal Appl, 2015, 21(6): 1169-1198. doi:  10.1007/s00041-015-9405-2
    [49] Elad M. Sparse And Redundant Representations: From Theory To Applications In Signal And Image Processing [M]. New York: Springer-Verlag, 2010.
    [50] Ranieri J, Chebira A, Lu Y M, et al. Phase retrieval for sparse signals: Uniqueness conditions [EB/OL]. (2013-08-14) [2022-10-11]. https://arxiv.org/abs/1308.3058.
    [51] Ohlsson H, Eldar Y C. On conditions for uniqueness in sparse phase retrieval [C]//2014 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2014: 1841-1845.
    [52] Jaganathan K, Oymak S, Hassibi B. Sparse phase retrieval: Uniqueness guarantees and recovery algorithms [J]. IEEE Trans Signal Process, 2017, 65(9): 2402-2410. doi:  10.1109/TSP.2017.2656844
    [53] Candès E J, Romberg J, Tao T. Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information [J]. IEEE Trans Inf Theory, 2006, 52(2): 489-509. doi:  10.1109/TIT.2005.862083
    [54] Fienup J R. Phase retrieval algorithms: A comparison [J]. Appl Optics, 1982, 21(15): 2758-2769. doi:  10.1364/AO.21.002758
    [55] Teague M R. Irradiance moments: Their propagation and use for unique retrieval of phase [J]. J Opt Soc Am, 1982, 72(9): 1199-1209. doi:  10.1364/JOSA.72.001199
    [56] Shen C, Liang M, Pan A, et al. Non-iterative complex wave-field reconstruction based on Kramers–Kronig relations [J]. Photonics Res, 2021, 9(6): 1003-1012. doi:  10.1364/PRJ.419886
    [57] Mukherjee S, Seelamantula C S. An iterative algorithm for phase retrieval with sparsity constraints: Application to frequency domain optical coherence tomography [C]//2012 IEEE International Conference on Acoustics, Speech and Signal Processing Processing (ICASSP), 2012: 553-556.
    [58] Shechtman Y, Beck A, Eldar Y C. GESPAR: Efficient phase retrieval of sparse signals [J]. IEEE Trans Signal Process, 2014, 62(1-4): 928-938.
    [59] Rivenson Y, Zhang Y, Günaydın H, et al. Phase recovery and holographic image reconstruction using deep learning in neural networks [J]. Light: Sci Appl, 2017, 7(2): 17141.
    [60] Sinha A, Lee J, Li S, et al. Lensless computational imaging through deep learning [J]. Optica, 2017, 4(9): 1117. doi:  10.1364/OPTICA.4.001117
    [61] Nguyen T, Xue Y, Li Y, et al. Deep learning approach for Fourier ptychography microscopy [J]. Opt Express, 2018, 26(20): 26470-26484. doi:  10.1364/OE.26.026470
    [62] Metzler C, Schniter P, Veeraraghavan A, et al. prDeep: Robust phase retrieval with a flexible deep network [C]//Proceedings of the 35th International Conference on Machine Learning, 2018, 80: 3501-3510.
    [63] Balan R, Casazza P, Edidin D. On signal reconstruction without phase [J]. Appl Comput Harmon Anal, 2006, 20(3): 345-356. doi:  10.1016/j.acha.2005.07.001
    [64] Conca A, Edidin D, Hering M, et al. An algebraic characterization of injectivity in phase retrieval [J]. Appl Comput Harmon Anal, 2015, 38(2): 346-356. doi:  10.1016/j.acha.2014.06.005
    [65] Eldar Y C, Mendelson S. Phase retrieval: Stability and recovery guarantees [J]. Appl Comput Harmon Anal, 2014, 36(3): 473-494. doi:  10.1016/j.acha.2013.08.003
    [66] Shechtman Y, Eldar Y C, Cohen O, et al. Phase retrieval with application to optical imaging: A contemporary overview [J]. IEEE Signal Proc Mag, 2015, 32(3): 87-109. doi:  10.1109/MSP.2014.2352673
    [67] Wackerman C C, Yagle A E. Use of fourier domain real-plane zeros to overcome a phase retrieval stagnation [J]. J Opt Soc Am A, 1991, 8(12): 1898-1904. doi:  10.1364/JOSAA.8.001898
    [68] Lu G, Zhang Z, Yu F T S, et al. Pendulum iterative algorithm for phase retrieval from modulus data [J]. Opt Eng, 1994, 33(2): 548-555. doi:  10.1117/12.153152
    [69] Takajo H, Takahashi T, Kawanami H, et al. Numerical investigation of the iterative phase-retrieval stagnation problem: Territories of convergence objects and holes in their boundaries [J]. J Opt Soc Am A, 1997, 14(12): 3175-3187. doi:  10.1364/JOSAA.14.003175
    [70] Soifer V A, Kotlar V, Doskolovich L. Iteractive Methods For Diffractive Optical Elements Computation [M]. London: Taylor & Francis Group, 1997.
    [71] Fienup J R. Reconstruction of an object from the modulus of its Fourier transform [J]. Opt Letters, 1978, 3(1): 27-29. doi:  10.1364/OL.3.000027
    [72] Fienup J R, Crimmins T, Holsztynski W. Reconstruction of the support of an object from the support of its autocorrelation [J]. J Opt Soc Am A, 1982, 72(5): 610-624. doi:  10.1364/JOSA.72.000610
    [73] Fienup J R. Phase retrieval with continuous version of hybrid input-output [C]//Frontiers in Optics, OSA Technical Digest (CD), 2003: ThI3.
    [74] Elser V. Phase retrieval by iterated projections [J]. J Opt Soc Am A, 2003, 20(1): 40-55. doi:  10.1364/JOSAA.20.000040
    [75] Luke D R. Relaxed averaged alternating reflections for diffraction imaging [J]. Inverse Probl, 2004, 21(1): 37.
    [76] Fienup J R, Wackerman C C. Phase-retrieval stagnation problems and solutions [J]. J Opt Soc Am A, 1986, 3: 1897-1907. doi:  10.1364/JOSAA.3.001897
    [77] Williams G J, Pfeifer M A, Vartanyants I A, et al. Three-dimensional imaging of microstructure in Au nanocrystals [J]. Phys Rev Lett, 2003, 90(17): 175501. doi:  10.1103/PhysRevLett.90.175501
    [78] Robinson I K, Vartanysnts I A, Williams G J, et al. Recon-struction of the shapes of gold nanocrystals using coherent X-Ray diffraction [J]. Phys Rev Lett, 2001, 87(19): 195505. doi:  10.1103/PhysRevLett.87.195505
    [79] Loh N T D, Eisebitt S, Flewett S, et al. Recovering magne-tization distributions from their noisy diffraction data [J]. Phys Rev E, 2010, 82(6): 061128. doi:  10.1103/PhysRevE.82.061128
    [80] Misell D L. A method for the solution of the phase problem in electron microscopy [J]. J Phys D: Appl Phys, 1973, 6(1): L6. doi:  10.1088/0022-3727/6/1/102
    [81] Bao P, Zhang F, Pedrini G, et al. Phase retrieval using multiple illumination wavelengths [J]. Opt Lett, 2008, 33(4): 309-311. doi:  10.1364/OL.33.000309
    [82] Allen L J, Oxley M P. Phase retrieval from series of images obtained by defocus variation [J]. Opt Commun, 2001, 1999(1-4): 65-75. doi:  10.1016/S0030-4018(01)01556-5
    [83] Zhang Y, Pedrini G, Osten W, et al. Whole optical wave field reconstruction from double or multi in-line holograms by phase retrieval algorithm [J]. Opt Express, 2003, 11(24): 3234-3241. doi:  10.1364/OE.11.003234
    [84] Almoro P, Pedrini G, Osten W. Complete wavefront reconstruction using sequential intensity measurements of a volume speckle field [J]. Appl Optics, 2006, 45(34): 8596-8605. doi:  10.1364/AO.45.008596
    [85] Gao P, Pedrini G, Zuo C, et al. Phase retrieval using spatially modulated illumination [J]. Opt Lett, 2014, 39(12): 3615. doi:  10.1364/OL.39.003615
    [86] Zhang F, Pedrini G, Osten W. Phase retrieval of arbitrary complex-valued fields through aperture- plane modulation [J]. Phys Rev A, 2007, 75(4): 043805. doi:  10.1103/PhysRevA.75.043805
    [87] Morris D. Simulated annealing applied to the Misell algorithm for phase retrieval [J]. Microwaves, Antennas and Propagation, 1996, 143(4): 298-303. doi:  10.1049/ip-map:19960446
    [88] Xu Y, Ye Q, Meng G. Hybrid phase retrieval algorithm based on modified very fast simulated annealing [J]. Int J Microw Wirel Technol, 2018, 10(9): 1072-1080. doi:  10.1017/S1759078718000971
    [89] Faulkner H M L, Rodenburg J M. Movable aperture lensless transmission microscopy: A novel phase retrieval algorithm [J]. Phys Rev Lett, 2004, 93(2): 023903. doi:  10.1103/PhysRevLett.93.023903
    [90] Rodenburg J M, Faulkner H M L. A phase retrieval algorithm for shifting illumination [J]. Appl Phys Lett, 2004, 85(20): 4795-4797. doi:  10.1063/1.1823034
    [91] Maiden A M, Rodenburg J M. An improved ptychographical phase retrieval algorithm for diffractive imaging [J]. Ultramicroscopy, 2009, 109(10): 1256-1262. doi:  10.1016/j.ultramic.2009.05.012
    [92] Maiden A M, Humphry M J, Rodenburg J. Ptychographic transmission microscopy in three dimensions using a multi-slice approach [J]. J Opt Soc Am A, 2012, 29(8): 1606-1614. doi:  10.1364/JOSAA.29.001606
    [93] Maiden A, Johnson D, Li P. Further improvements to the ptychographical iterative engine [J]. Optica, 2017, 4(7): 736-745. doi:  10.1364/OPTICA.4.000736
    [94] Zheng G, Horstmeyer R, Yang C. Wide-field, high-resolution Fourier ptychographic microscopy [J]. Nat Photonics, 2013, 7(9): 739-745. doi:  10.1038/nphoton.2013.187
    [95] Ou X, Horstmeyer R, Yang C, et al. Quantitative phase imaging via Fourier ptychographic microscopy [J]. Opt Lett, 2013, 38(22): 4845. doi:  10.1364/OL.38.004845
    [96] Pan A, Zhang Y, Zhao T, et al. System calibration method for Fourier ptychographic microscopy [J]. J Biomed Optics, 2017, 22(9): 096005.
    [97] Zhang Y, Pan A, Lei M, et al. Data preprocessing methods for robust Fourier ptychographic microscopy [J]. Opt Eng, 2017, 56(12): 123107.
    [98] Pan A, Zuo C, Xie Y, et al. Vignetting effect in Fourier ptychographic microscopy [J]. Opt Laser Eng, 2019, 120: 40-48. doi:  10.1016/j.optlaseng.2019.02.015
    [99] Pan A, Shen C, Yao B, et al. Single-shot Fourier ptychographic microscopy via annular monochrome LED array [C]//Frontiers in Optics + Laser Science APS/DLS, 2019: FTh3 F.4.
    [100] Pan A, Zhang Y, Wen K, et al. Subwavelength resolution Fourier ptychography with hemispherical digital condensers [J]. Opt Express, 2018, 26(18): 23119-23131. doi:  10.1364/OE.26.023119
    [101] Pan A, Yao B. Three-dimensional space optimization for near-field ptychography [J]. Opt Express, 2019, 27(4): 5433-5446. doi:  10.1364/OE.27.005433
    [102] Chan A, Kim J, Pan A, et al. Parallel Fourier ptychographic microscopy for high-throughput screening with 96 cameras (96 Eyes) [J]. Sci Rep, 2019, 9: 11114. doi:  10.1038/s41598-019-47146-z
    [103] Pan A, Wen K, Yao B. Linear space-variant optical cryptosystem via Fourier ptychography [J]. Opt Lett, 2019, 44(8): 2032-2035. doi:  10.1364/OL.44.002032
    [104] Xiang M, Pan A, Zhao Y, et al. Coherent synthetic aperture imaging for visible remote sensing via reflective Fourier ptychography [J]. Opt Lett, 2021, 46(1): 29-32. doi:  10.1364/OL.409258
    [105] Gao Y, Chen J, Wang A, et al. High-throughput fast full-color digital pathology based on Fourier ptychographic microscopy via color transfer [J]. Sci China-Phys Mech, 2021, 64: 114211. doi:  10.1007/s11433-021-1730-x
    [106] Wang A, Zhang Z, Wang S, et al. Fourier ptychographic microscopy via alternating direction method of multipliers [J]. Cells, 2022, 11(9): 1512. doi:  10.3390/cells11091512
    [107] Sun Jiasong, Zhang Yuzhen, Chen Qian, et al. Fourier ptychographic microscopy: Theory, advances, and applications [J]. Acta Optica Sinica, 2016, 36(10): 1011005. (in Chinese)
    [108] Pan An, Yao Baoli. High-throughput and fast-speed Fourier ptychographic microscopy [J]. Infrared and Laser Engi-neering, 2019, 48(6): 0603012. (in Chinese) doi:  10.3788/IRLA201948.0603012
    [109] Konda P C, Loetgering L, Zhou K C, et al. Fourier ptycho-graphy: Current applications and future promises [J]. Opt Express, 2020, 28(7): 9603-9630. doi:  10.1364/OE.386168
    [110] Pan A, Zuo C, Yao B. High-resolution and large field-of-view Fourier ptychographic microscopy and its applications in biomedicine [J]. Rep Prog Phys, 2020, 83(9): 096101. doi:  10.1088/1361-6633/aba6f0
    [111] Zheng G, Shen C, Jiang S, et al. Concept, implementations and applications of Fourier ptychography [J]. Nat Rev Phys, 2021, 3(3): 207-223. doi:  10.1038/s42254-021-00280-y
    [112] Zhang Shaohui, Zhou Guocheng, Cui Baiqi, et al. Review of Fourier ptychographic microscopy: Models, algorithms, and systems [J]. Laser & Optoelectronics Progress, 2021, 58(14): 1400001. (in Chinese)
    [113] Teague M R. Deterministic phase retrieval: A green’s function solution [J]. J Opt Soc Am A, 1983, 73(11): 1434-1441. doi:  10.1364/JOSA.73.001434
    [114] Guigay J P. Fourier transform analysis of fresnel diffraction patterns and in-line holograms [J]. Optik, 1977, 49: 121-125.
    [115] Paganin D, Nugent K A. Noninterferometric phase imaging with partially coherent light [J]. Phys Rev Lett, 1998, 80(12): 2586.
    [116] Gureyev T E, Raven C, Snigirev A, et al. Hard X-ray quantitative non-interferometric phase- contrast microscopy [J]. J Phys Appl Phys, 1999, 32(5): 563. doi:  10.1088/0022-3727/32/5/010
    [117] Pinhasi S V, Alimi R, Perelmutter L, et al. Topography retrieval using different solutions of the transport intensity equation [J]. J Opt Soc Am A, 2010, 27(10): 2285-2292. doi:  10.1364/JOSAA.27.002285
    [118] Xue B, Zheng S. Phase retrieval using the transport of intensity equation solved by the FMG-CG method [J]. Opt-Int J Light Electron Opt, 2011, 122(23): 2101-2106. doi:  10.1016/j.ijleo.2011.01.004
    [119] Voitsekhovich V V. Phase-retrieval problem and orthogonal expansions: Curvature sensing [J]. J Opt Soc Am A, 1995, 12(10): 2194-2202. doi:  10.1364/JOSAA.12.002194
    [120] Ros S, Acosta E, Bar S. Modal phase estimation from wavefront curvature sensing [J]. Opt Commun, 1996, 123(4): 453-456.
    [121] Volkov V V, Zhu Y, De Graef M. A new symmetrized solution for phase retrieval using the transport of intensity equation [J]. Micron, 2002, 33(5): 411-416. doi:  10.1016/S0968-4328(02)00017-3
    [122] Frank J, Altmeyer S, Wernicke G. Non-interferometric, non-iterative phase retrieval by green’s functions [J]. J Opt Soc Am A, 2010, 27(10): 2244-2251. doi:  10.1364/JOSAA.27.002244
    [123] Zuo C, Chen Q, Asundi A. Boundary-artifact-free phase retrieval with the transport of intensity equation: fast solution with use of discrete cosine transform [J]. Opt Express, 2014, 22(8): 9220. doi:  10.1364/OE.22.009220
    [124] Zuo C, Chen Q, Li H, et al. Boundary-artifact-free phase retrieval with the transport of intensity equation II: Applications to microlens characterization [J]. Opt Express, 2014, 22(15): 18310. doi:  10.1364/OE.22.018310
    [125] Huang L, Zuo C, Idir M, et al. Phase retrieval with the transport-of-intensity equation in an arbitrarily shaped aperture by iterative discrete cosine transforms [J]. Opt Lett, 2015, 40(9): 1976. doi:  10.1364/OL.40.001976
    [126] Zuo C, Li J, Sun J, et al. Transport of intensity equation: A tutorial [J]. Opt Laser Eng, 2020, 152: 106187.
    [127] Hall S H, Heck H L. Advanced Signal Integrity for High-speed Digital Designs [M]. Hoboken, NJ, USA: Wiley, 2009.
    [128] Graf U. Introduction to Hyperfunctions and Their Integral Transforms: An Applied and Computational Approach [M]. Basel: Birkhauser, 2010.
    [129] Ikeda T, Popescu G, Dasari R R, et al. Hilbert phase microscopy for investigating fast dynamics in transparent systems [J]. Opt Lett, 2005, 30(10): 1165-1167. doi:  10.1364/OL.30.001165
    [130] Beak Y, Lee K, Shin S, et al. Kramers–Kronig holographic imaging for high-space-bandwidth product [J]. Optica, 2019, 6(1): 45-51. doi:  10.1364/OPTICA.6.000045
    [131] Huang Z, Cao L. High bandwidth-utilization digital holo-graphic multiplexing: An approach using Kramers-Kronig relations [J]. Adv Photonics Res, 2022, 3: 2100273. doi:  10.1002/adpr.202100273
    [132] Gillis N, Glineur F. Low-rank matrix approximation with weights or missing data is NP-hard [J]. SIAM J Matrix Anal Appl, 2011, 32(4): 1149-1165. doi:  10.1137/110820361
    [133] Lu Y M, Li G. Phase transitions of spectral initialization for high-dimensional non-convex estimation [J]. Information and Inference: A Journal of the IMA, 2017, 9(3): 507-541.
    [134] Netrapalli P, Jain P, Sanghavi S. Phase retrieval using alternating minimization [J]. IEEE Trans Signal Process, 2015, 63(18): 4814-4826.
    [135] Cai T, Fan J, Jiang T. Distributions of angles in random packingon spheres [J]. J Mach Learn Res, 2013, 14(1): 1837-1864.
    [136] Wang G, Giannakis G B, Saad Y. Phase retrieval via reweighted amplitude flow [J]. IEEE Trans Signal Process, 2018, 66: 2818-2833.
    [137] Yuan Z, Wang H, Wang Q. Phase retrieval via sparse Wirtinger flow[J]. J Comput Appl Math, 2019, 355: 162-173.
    [138] Kolte R, Özgür A. Phase retrieval via incremental truncated Wirtinger flow [EB/OL]. (2016-06-10) [2022-08-03]. https://arxiv.org/abs/1606.03196.
    [139] Zhang H, Chi Y, Liang Y. Median-truncated nonconvex approach for phase retrieval with outliers [J]. IEEE Trans Inf Theory, 2018, 64(11): 7287-7310.
    [140] Wang G, Giannakis G B, Saad Y, et al. Scalable solvers of random quadratic equations via stochastic truncated amplitude flow [J]. IEEE Trans Signal Process, 2017, 65(8): 1961-1974.
    [141] Pinilla S, Bacca J, Arguello H. Phase retrieval algorithm via nonconvex minimization using a smoothing function [J]. IEEE Trans Signal Process, 2018, 66(17): 4574-4584. doi:  10.1109/TSP.2018.2855667
    [142] Luo Q, Wang H, Lin S. Phase retrieval via smoothed amplitude flow [J]. Signal Process, 2020, 177: 107719. doi:  10.1016/j.sigpro.2020.107719
    [143] Loock S, Plonka G. Phase retrieval for Fresnel measurements using a Shearlet sparsity constraint [J]. Inverse Probl, 2014, 30(5): 055005. doi:  10.1088/0266-5611/30/5/055005
    [144] Yang Zhenya, Zheng Chujun. Phase retrieval of pure phase object based on compressed sensing [J]. Acta Physica Sinica, 2013, 62(10): 104203. (in Chinese) doi:  10.7498/aps.62.104203
    [145] Shechtman Y, Eldar Y C, Szameit A, et al. Sparsity based sub-wavelength imaging with partially incoherent light via quadratic compressed sensing [J]. Opt Express, 2011, 19(16): 14807-14822. doi:  10.1364/OE.19.014807
    [146] Ohlsson H, Yang A Y, Dong R, et al. Compressive phase retrieval from squared output measurements via semidefinite programming [J]. IFCA Proceedings Volumes, 2012, 45(16): 89-94. doi:  https://doi.org/10.3182/20120711-3-BE-2027.00415
    [147] Zuo Chao, Feng Shijie, Zhang Xiangyu, et al. Deep learning based computational imaging: Status, challenges, and future [J]. Acta Optica Sinica, 2020, 40(1): 0111003. (in Chinese) doi:  10.3788/AOS202040.0111003
    [148] Nishizaki Y, Horisaki R, Kitaguchi K, et al. Analysis of non-iterative phase retrieval based on machine learning [J]. Opt Rev, 2020, 27(1): 136-141. doi:  10.1007/s10043-019-00574-8
    [149] Wang F, Bian Y, Wang H, et al. Phase imaging with an untrained neural network [J]. Light: Sci Appl, 2020, 9: 77. doi:  10.1038/s41377-020-0302-3
    [150] Naimipour N, Khobahi S, Soltanalian M. UPR: A model-driven architecture for deep phase retrieval [C]//54th Asilomar Conference on Signals, Systems, and Computers, 2020: 205-209.
    [151] Boutet S, Robinson I K. Coherent X-ray diffractive imaging of protein crystals [J]. J Synchrotron Radiat, 2008, 15(6): 576-583. doi:  10.1107/S0909049508029439
    [152] Parker M W. Protein structure from X-ray diffraction [J]. Journal of Biological Physics, 2003, 29(4): 341-362. doi:  10.1023/A:1027310719146
    [153] Miao J, Ohsuna T, Terasaki O, et al. Atomic resolution three-dimensional electron diffraction microscopy [J]. Phys Rev Lett, 2002, 89(15): 155502. doi:  10.1103/PhysRevLett.89.155502
    [154] Kamimuraet O, Kawahara K, Doi T, et al. Diffraction microscopy using 20 kV electron beam for multiwall carbon nanotubes [J]. Appl Phys Lett, 2008, 92(2): 024106. doi:  10.1063/1.2834372
    [155] Maiden A M, Humphry M J, Zhang F, et al. Superresolution imaging via ptychography [J]. J Opt Soc Am A, 2011, 28(4): 604-612. doi:  10.1364/JOSAA.28.000604
    [156] Maiden A M, Rodenburg J M, Humphry M J. Optical ptychography: A practical implementation with useful resolu-tion [J]. Opt Lett, 2010, 35(15): 2585-2587. doi:  10.1364/OL.35.002585
    [157] Chapman H N, Nugent K A. Coherent lensless X-ray imaging [J]. Nat Photonics, 2010, 4(12): 833-839. doi:  10.1038/nphoton.2010.240
    [158] Wu J, Zhang H, Zhang W, et al. Single-shot lensless imaging with fresnel zone aperture and incoherent illumination [J]. Light: Sci Appl, 2020, 9(1): 53. doi:  10.1038/s41377-020-0289-9
    [159] Alfalou A, Brosseau C. Optical image compression and encryption methods [J]. Adv Opt Photonics, 2009, 1(3): 589. doi:  10.1364/AOP.1.000589
    [160] Wang R K, Watson I A, Chatwin C R. Random phase encoding for optical security [J]. Opt Eng, 1996, 35(9): 2464-2469. doi:  10.1117/1.600849
    [161] Situ Guohai, Zhang Jingjuan, Zhang Yan, et al. A cascaded-phases retrieval algorithm for optical image encryption [J]. Journal of Optoelectronics · Laser, 2004, 15(3): 341. (in Chinese)
    [162] Shi Y, Situ G, Zhang J. Multiple-image hiding in the Fresnel domain [J]. Opt Lett, 2007, 32(13): 1914-1916. doi:  10.1364/OL.32.001914
    [163] Huang J J, Hwang, H E, Chen C Y, et al. Lensless multiple-image optical encryption based on improved phase retrieval algorithm [J]. Appl Opt, 2012, 51(13): 2388. doi:  10.1364/AO.51.002388
    [164] Guo C, Liu S, Sheridan J T. Iterative phase retrieval algorithms. Part I: Optimization [J]. Appl Opt, 2015, 54(15): 4698-4708. doi:  10.1364/AO.54.004698
    [165] Guo C, Liu S, Sheridan J T. Iterative phase retrieval algorithms. Part II: Attacking optical encryption systems [J]. Appl Opt, 2015, 54(15): 4709-4718. doi:  10.1364/AO.54.004709
    [166] Guo C, Wei C, Tan J, et al. A review of iterative phase retrieval for measurement and encryption [J]. Opt Lasers Eng, 2017, 89: 2-12. doi:  10.1016/j.optlaseng.2016.03.021
    [167] Wang S, Meng X, Wang Y, et al. Phase retrieval algorithm for optical information security [J]. Chin Phys B, 2019, 28(8): 084203. doi:  10.1088/1674-1056/28/8/084203
    [168] Hu L, Liu C, Shen W, et al. Advancement of adaptive optics in astronomical observation [J]. Sci China-Phys Mech, 2017, 47(8): 084202. doi:  https://doi.org/10.1360/SSPMA2016-00425
    [169] Sandler D G, Stahl S, Angel J R P, et al. Adaptive optics for diffraction-limited infrared imaging with 8-m telescopes [J]. J Opt Soc Am A, 1994, 11(2): 925-945. doi:  10.1364/JOSAA.11.000925
    [170] Rao C, Zhu L, Rao X, et al. 37-element solar adaptive optics for 26-cm solar fine structure telescope at Yunnan Astronomical Observatory [J]. Chin Opt Lett, 2010, 8(10): 966-968. doi:  10.3788/COL20100810.0966
    [171] Baranova N B, Mamaev A V, Pilipetsky N F, et al. Wave-front dislocations: Topological limitations for adaptive systems with phase conjugation [J]. J Opt Soc Am, 1983, 73(5): 525-528. doi:  10.1364/JOSA.73.000525
    [172] Ping Y, Ming A, Yuan L, et al. Intracavity transverse modes controlled by a genetic algorithm based on Zernike mode coefficients [J]. Opt Express, 2007, 15(25): 17051. doi:  10.1364/OE.15.017051
    [173] Zommer S, Ribak E N, Lipson S G, et al. Simulated annealing in ocular adaptive optics [J]. Opt Letters, 2006, 31(7): 939-941. doi:  10.1364/OL.31.000939
    [174] El-Agmy R, Bulte H, Greenaway A H, et al. Adaptive beam profile control using a simulated annealing algorithm [J]. Opt Express, 2005, 13(16): 6085. doi:  10.1364/OPEX.13.006085
    [175] Zakynthinaki M S, Saridakis Y G. Stochastic optimization for adaptive real-time wavefront correction [J]. Numerical Algorithms, 2003, 33(1-4): 509-520.
    [176] Feng L, Zeng Z, Wu Y. Phase retrieval hybrid algorithm for optical surface testing of the high dynamic range error [C]//Proceedings of SPIE, 2014, 9282: 92822Y.
    [177] Fienup J R, Marron J C, Schulz T J, et al. Hubble space telescope characterized by using phase retrieval algorithms [J]. Appl Optics, 1993, 32(10): 1747-1767.
    [178] Dean B H, Aronstein D L, Smith J S, et al. Phase retrieval algorithm for JWST flight and testbed telescope [C]//Procee-dings of SPIE, 2006, 6265: 626511.
    [179] Li Shengyi, Hu Xiaojun, Wu Yulie. Phase retrieval on site testing for large mirrors [J]. Acta Photonica Sinica, 2009, 38(2): 365. (in Chinese)
    [180] Wu Yulie, Hu Xiaojun, Dai Yifan, et al. In-situ surface measurement for large aperture optical mirror based on phase retrieval technology [J]. Journal of Mechanical Engineering, 2009, 45(2): 157-163. (in Chinese) doi:  10.3901/JME.2009.02.157
  • [1] 胡以华, 赵禄达.  激光成像处理技术研究现状及展望(特邀) . 红外与激光工程, 2023, 52(6): 20230169-1-20230169-21. doi: 10.3788/IRLA20230169
    [2] 朱硕, 郭恩来, 柏连发, 韩静.  高效学习的透过未知散射介质的相位恢复方法 . 红外与激光工程, 2022, 51(2): 20210889-1-20210889-9. doi: 10.3788/IRLA20210889
    [3] 郭恩来, 师瑛杰, 朱硕, 程倩倩, 韦一, 苗金烨, 韩静.  深度学习下的散射成像:物理与数据联合建模优化(特邀) . 红外与激光工程, 2022, 51(8): 20220563-1-20220563-13. doi: 10.3788/IRLA20220563
    [4] 和河向, 黎永耀, 黄锦圣.  基于相位共轭实现散射成像及光学幻像的双功能散射光调控方法(特邀) . 红外与激光工程, 2022, 51(8): 20220266-1-20220266-12. doi: 10.3788/IRLA20220266
    [5] 李丽艳, 范松涛, 周燕.  消光强扰动的相位载波解调算法 . 红外与激光工程, 2021, 50(9): 20210485-1-20210485-8. doi: 10.3788/IRLA20210485
    [6] 程鸿, 熊帮玲, 王金成, 马慧敏, 张芬, 韦穗.  透镜模型下基于色散和强度传输方程的相位恢复技术 . 红外与激光工程, 2019, 48(6): 603018-0603018(6). doi: 10.3788/IRLA201948.0603018
    [7] 韩庆阳, 陈赟, 张红胜, 高胜英, 张晰.  航天级反射式光电编码器信号处理系统设计 . 红外与激光工程, 2019, 48(S1): 135-141. doi: 10.3788/IRLA201948.S117010
    [8] 张佳琳, 陈钱, 张翔宇, 孙佳嵩, 左超.  无透镜片上显微成像技术:理论、发展与应用 . 红外与激光工程, 2019, 48(6): 603009-0603009(33). doi: 10.3788/IRLA201948.0603009
    [9] 赵楠翔, 胡以华.  激光反射层析成像相位恢复算法研究 . 红外与激光工程, 2019, 48(10): 1005005-1005005(7). doi: 10.3788/IRLA201948.1005005
    [10] 潘安, 姚保利.  高通量快速傅里叶叠层显微成像技术研究进展 . 红外与激光工程, 2019, 48(6): 603012-0603012(19). doi: 10.3788/IRLA201948.0603012
    [11] 刘正君, 耿勇, 谭久彬.  基于柱透镜多旋转测量的计算成像 . 红外与激光工程, 2019, 48(6): 603016-0603016(5). doi: 10.3788/IRLA201948.0603016
    [12] 程鸿, 邓会龙, 沈川, 王金成, 韦穗.  光强传输方程与图像插值融合的相位恢复 . 红外与激光工程, 2018, 47(10): 1026003-1026003(7). doi: 10.3788/IRLA201847.1026003
    [13] 刘正君, 郭澄, 谭久彬.  基于多距离相位恢复的无透镜计算成像技术 . 红外与激光工程, 2018, 47(10): 1002002-1002002(16). doi: 10.3788/IRLA201847.1002002
    [14] 刘子伟, 许廷发, 王洪庆, 申子宜, 饶志涛.  深度成像理论与实现 . 红外与激光工程, 2016, 45(7): 726001-0726001(5). doi: 10.3788/IRLA201645.0726001
    [15] 马鑫雪, 王建立, 王斌.  利用相位恢复波前传感技术检测球面镜面形 . 红外与激光工程, 2014, 43(10): 3428-3433.
    [16] 曾文雯, 钟小品, 李景镇.  从单幅干涉图中恢复相位的区间反转方法 . 红外与激光工程, 2014, 43(9): 3151-3156.
    [17] 郭腾霄, 丁学全, 董晓强, 穆宁, 黄启斌, 李翠萍, 温红宇.  基于EMD的红外遥测光谱信号预处理新方法 . 红外与激光工程, 2013, 42(12): 3196-3200.
    [18] 曹文勇, 马明, 赵彬, 何幸锴, 侯天晋, 周鼎富, 邓华荣.  滤波器组在相干激光测风雷达信号处理中的应用 . 红外与激光工程, 2013, 42(5): 1179-1183.
    [19] 万俊, 张晓晖, 饶炯辉, 梁瑞涛.  基于独立成分分析的舰船气泡尾流后向散射光信号处理 . 红外与激光工程, 2013, 42(1): 244-250.
    [20] 俞建杰, 马晶, 谭立英, 韩琦琦.  基于相位混合算法的衍射光学元件优化设计方法 . 红外与激光工程, 2013, 42(9): 2472-2477.
  • 加载中
图(9) / 表(6)
计量
  • 文章访问数:  1453
  • HTML全文浏览量:  261
  • PDF下载量:  563
  • 被引次数: 0
出版历程
  • 收稿日期:  2022-06-13
  • 修回日期:  2022-08-03
  • 刊出日期:  2022-11-30

相位恢复算法:原理、发展与应用(特邀)

doi: 10.3788/IRLA20220402
    作者简介:

    王爱业,男,硕士生,主要从事傅里叶叠层成像技术方面的研究

    潘安,男,副研究员,博士生导师,博士,主要从事计算成像技术方面的研究

    马彩文,男,研究员,博士生导师,博士,主要从事光电测控技术方面的研究

基金项目:  国家自然科学基金(12104500)
  • 中图分类号: O438

摘要: 研究表明,由于相位比振幅包含更多关于场的信息,因此相位测量在现代科学和工程的诸多分支中始终是研究的热点问题。在可见的电磁波范围内,相位信息很难通过现有的光电探测器直接采集获取。相位恢复技术提供了一种从捕获的强度信息中将相位信息“计算”出来的有效手段,并已成功应用于天文观测、生物医学成像和数字信号复原等多个科学领域。算法是相位恢复技术的核心,也是该技术发展和应用的关键。文中结合物理学原理和信号处理方法对相位恢复算法的基本原理进行阐述,综述了各类相位恢复算法的发展历程及其优缺点,并简单概述了相位恢复算法在光学领域的典型应用,最终指明其面临的挑战和未来的发展趋势:更优异的收敛性能和噪声鲁棒性、恢复更复杂物体相位信息的能力、多目标多任务集成的兼容性。

English Abstract

    • 在波动光学中,单色相干光场可以使用二维复振幅函数来描述,包含幅值和相位两个组成部分[1-2]。可见光的高时间频率和现有成像探测器的有限带宽使得光学相位信息无法被直接记录[3],通常只能探测到光场的强度(幅值的平方)信息[4]。基于强度测量的光学成像系统,包括照相相机和人眼等,足以满足日常生活中对宏观世界观察的需求。然而,当人们试图在光学显微镜下观察微观世界中的生命时,相位信息的重要性[5]便表现得尤为突出。由于大多数生物细胞和组织是透明的,对光的吸收很少,因此使用传统的亮场显微镜来观察会存在成像对比度低的问题。为了克服这一缺陷,研究人员设计了各种外源性标记策略,其中荧光标记手段利用细胞不同部分对于特定染料的亲和力发射出区分明显的荧光光谱,来提高观测的对比度[6]。然而,荧光剂的光漂白特性以及光毒性可能会干扰内源性生物过程,同时也对活细胞的长时间成像提出了严峻挑战[7]。此后,人们开始探究无标记成像的新思路,即通过测量光波透过样品前后的相位变化来清楚地“看见”样品的内部结构。

      无标记高对比度的活细胞成像最早是通过泽尼克发明的相衬显微镜[8]来实现的,该方法基于4f系统的空间滤波将相位差转化为可视化的强度变化信息。尽管相衬显微镜及其变体(如差分干涉相衬技术[9])实现了对透明生物细胞的活体观察,但是测量强度与样品相位之间并非严格的线性关系,因此相位信息的可视化只是定性层面的。随着数字图像传感器的出现以及信息光学理论的完善,定量相位成像(QPI)[10-11]成为了成像技术领域的一大重要课题。以数字全息术[12-14]为代表的干涉测量方法将目标光场的相位添加到额外引入的参考光场中发生干涉,然后从获取的干涉条纹中将可见的相位信息解调出来。干涉测量方法的缺点是对于外部干扰非常敏感,因此往往对实验环境的稳定性有很高的要求,这在天文观测、晶体学、大孔径光学测试等应用场景中是很难严格满足的。

      波前传感技术[15-21]是一类实现简单、无需参考光的非干涉相位测量手段。典型的夏克-哈特曼波前传感器[15-17]使用微透镜阵列将入射光场聚焦到传感器的点阵上,通过测量该点阵相对于理想无像差情况下光斑阵列的偏移量重构波前相位分布,目前已在显微像差校正[22-23]、大气激光通信[24-25]以及眼科学[26-27]等领域取得广泛的应用。另一类基于衍射成像的定量相位测量方法称为相位恢复,即从测量的衍射场强度模式中反演出相位分布信息。尽管图像传感器无法直接记录光场的相位信息,但是光场具备的某些特性使得相位分布可以借助于算法恢复出来,即光波远场的信息对应于近场信息的傅里叶变换。相位恢复算法提供了在已知远场强度测量值以及某些先验知识的前提下重构相位的可能性。

      相位恢复的最初设想来自于晶体学研究领域,具有周期性结构的晶体在远场信息中通常会包含傅里叶变换峰值。Sayre[28]猜想,在一些简单的情形下可以根据傅里叶强度测量的某些精细特征(峰值信息)将散射波的相位信息恢复出来。20世纪70年代,以GS算法[29-30]为代表的交替投影迭代算法首先将这一设想付诸实践,从二维图像的傅里叶模值以及约束中重构相位信息。20世纪末,Miao等人首次使用X射线记录和重建非晶体结构的连续衍射图案,引发了人们对光学相位恢复新的兴趣[31]。实际上,人们的初衷是提高光学成像系统的空间分辨率。根据阿贝衍射公式,分子尺度的特征无法在显微镜中使用可见光观察,X射线则可以提供与原子尺度相当的衍射分辨率。遗憾的是,简单地采用波长更短的电磁波是不切实际的,因为在X射线波段范围内工作的成像光学元件和高相干光源的制造难度极高[31-32]。此外,随着分辨率的提高,必须适当地增加曝光时间以提供合理的成像信噪比(SNR)[33],高能量的X射线会对样品造成不可逆的损伤,因而远场强度测量值无法真实反映样品的微观结构。

      在研究X射线成像早期遭遇的难题迫使研究人员选择相位恢复的替代手段来重建微观尺度样品的结构信息。X射线衍射、过采样和相位恢复的结合催生出目前非常活跃的研究领域,称为相干衍射成像(CDI)。在距离样品足够远处放置一定范围的探测器来记录远场衍射强度信息(与样品傅里叶变换的幅值平方成正比),通过对样品施加合适的支持性约束,能够实现纳米分辨率的相位恢复。

      近十年来,光学相位恢复问题的研究热潮对信号处理领域也产生了深远的影响。人们跳出光学领域从衍射场傅里叶幅度测量值中恢复二维图像相位信息的思维枷锁,构建起针对多维度信号的广义相位恢复问题模型,利用最优化理论开发出了一系列性能强大的相位恢复算法。基于凸松弛策略的优化算法[34-35]可以提供相位信息精确恢复的保证,但需要将问题提升到更高的维度空间,解决大规模半正定规划;非凸优化算法[36-38]的恢复保证依赖于精确的初始值设定,能够在给定足够多测量值的情况下获得成功。

      随着压缩感知技术的发展与普及,研究人员已经开始探索相位恢复与稀疏性信息处理的关联。通过引入稀疏性先验知识,相位恢复算法的重建性能和稳定性得到改善。值得注意的是,由于光学系统对信号的采集过程是非线性的,在相位恢复中使用基于稀疏性的概念通常需要修改传统压缩感知理论中的线性采样模型[39],反过来也对非线性压缩感知的理论研究起到了推动作用。深度学习作为近年来备受关注的计算框架,将图像稀疏性通过正则化项融合到相位恢复中来提升重建质量,为相位恢复算法的发展注入了崭新的活力。

      相位恢复提供了一种从捕获的强度信息中将相位信息“计算”出来的有效手段,并已成功应用于天文观测、生物医学成像和数字信号复原等多个科学领域。算法是相位恢复的核心,也是该技术发展和应用的关键。文中将相关的物理原理与信号处理方法和算法联系起来,回顾了相位恢复算法的发展历程。图1简要绘制了相位恢复算法发展的时间线。文中的目标是描述相位恢复领域当前的技术状态和应用场景,确定其面临的挑战,并提出未来的发展趋势。相位恢复问题在诸多科学分支中具有关键性意义,不断探索与开发性能更加优异的相位恢复算法会对这些领域的成功实践与良性发展做出重要的贡献。

      图  1  相位恢复算法发展时间线

      Figure 1.  Development timeline of phase retrieval algorithm

    • 相位恢复问题的目标是从信号强度或幅度测量值中恢复原始信号及其相位。使用一维离散向量$x \in {C^n}$来表示待恢复信号,通过已知的采样向量${a_i} \in {C^n}$得到信号的测量值$\left\langle {{a_i},x} \right\rangle $。信号的测量值不包含相位,因此用模值的平方来表示强度测量值,即有:

      $$ {y_i} = {\left| {\left\langle {{a_i},x} \right\rangle } \right|^2},i = 1,2,\cdots m $$ (1)

      相位恢复的目标是求解上述二次方程组来获取原始信号x丢失的相位信息。如果定义采样矩阵${A^{\rm{T}}} = \left[ {{a_1},{a_2},\cdots {a_m}} \right] \in {C^{n \times m}}$以及强度向量${y^{\rm{T}}} = \left[ {y_1},{y_2},\cdots \right. \left.{y_m} \right] \in {C^{1 \times m}}$,则相位恢复问题可以进一步地归纳为:

      $$ \begin{split} &{\rm{find}}\;\; x \in {C^n}\\ &{\rm{s.t.}}\;\; y = {\left| {Ax} \right|^2} \end{split} $$ (2)

      以上给出了相位恢复问题的数学模型,该模型具有广泛的适用性,可以非常直观地向更高维度的信号领域拓展。例如,在光学相位恢复问题中常使用二维连续函数来表示感兴趣的物体以及光场的分布。实际上,光学系统中数据的采集以及后续的数据处理通常是离散化的,因此上述模型同样适用于光学相位恢复问题的研究。

    • 对于傅里叶相位恢复问题,采样向量$ {a_i} $对应于适当维度的离散傅里叶变换,即${\alpha _i}\left[ k \right] = {{\rm{e}}^{j2\pi \left( {i - 1} \right)\left( {k - 1} \right)/m}}$,采样矩阵则为m行的离散傅里叶矩阵取前n列组成的m×n矩阵${F_{mn}}$,故傅里叶相位恢复问题可以表述为:

      $$ \begin{split} &{\rm{find}}\;\;x \in C^n\\ & {\rm{s.t.}} \;\;y = {\left| {{F_{mn}}x} \right|^2} \end{split} $$ (3)

      一维傅里叶相位恢复问题是多解的[40]。许多完全不同的信号具有相同的傅里叶变换强度,这意味着从傅里叶幅度之中恢复相位信息的解并不是唯一的。以下三种变换以及它们的线性组合具有相同的傅里叶变换幅度:

      (1) 全局相移:$x\left[ k \right] \to x\left[ k \right]{{\rm{e}}^{j\theta }}$

      (2) 空间移位:$x\left[ k \right] \to x\left[ {k + {k_0}} \right]$

      (3) 共轭转置:$x\left[ k \right] \to x\overline {\left[ { - k} \right]} $

      这一类解并不会改变相位恢复解的整体轮廓,只是对原图像进行位移或者旋转,由此产生的解模糊性称为平凡模糊[41]。此外,任何一对具有相同自相关函数的一维信号也会产生相同的傅里叶幅度,这类解的存在则会导致非平凡模糊[41-42]。有关二维及以上信号的解的唯一性的探讨,可以参考有关文献[43-46]

      消除解模糊性的一种方式是对信号进行过采样以获取更多关于信号的信息。通常而言,在多维信号的离散傅里叶相位恢复中,当过采样倍数为2时,可以从傅里叶测量值中唯一地恢复出原始信号的相位[44]。然而,过采样意味着增加测量值的数量,导致观测效率降低,不利于数据的处理、储存与传输。此外,为了尽可能得到傅里叶相位恢复问题的唯一解,还可以添加关于原始信号的先验信息,对信号的重建施加约束。例如,如果真实信号表示某种强度或者概率分布,那么信号本身必须是实值且非负的[47-48]。目前比较热门的一种先验信息基于信号的稀疏性[49-53],即信号向量中仅包含少量的非零元素。尽管大部分自然信号本身不具有稀疏性,但是经过适当的线性变换之后便成为稀疏信号。

      不同的先验知识催生了大量相位恢复算法,表1将现有的相位恢复算法根据先验信息的不同进行了分类,具体内容将从下文开始逐一进行介绍。

      表 1  相位恢复算法的分类

      Table 1.  Classification of phase retrieval algorithms

      AlgorithmProsCons
      Alternating projections[29,54] Simple implementation
      Wide applicability
      Time-consuming
      Low precision
      Non-
      prior
      TIE[55] Deterministic solution
      Simple computation
      Require no phase unwrapping and complex optical systems
      Harsh preconditions
      Lower precision compared with interferometry
      Limited applications
      KK relations[56]Reconstruction with high quality and high speed
      Require less data collection
      Require high precision of NA matching
      SDP[34-35] Guarantee for convergence and reconstruction accuracy Not applicable to large-scale signals
      Non-convex
      optimization[36-38]
      Struggle with parameter
      optimization
      Rely on initial estimate
      Sparsity-
      based
      General[57-58] Reconstruction with high quality and stability
      Low computational cost
      Applicable to objects with sparsity only
      Deep learning[59-62] Powerful expressive ability Lack of interpretability
      Unable to deal with multiple models and applications
    • 考虑采样向量为一般测量向量的广义相位恢复问题能够更深入地分析问题的本质,从而可以更简单地推导理论保证。关于确保解的唯一性所需要的测量值数量和性质,已有的数理结论主要涉及独立同分布的随机测量向量。具体而言,解的唯一性探讨分为实值和复值两种情形:

      (1) 当模型采用实值测量,即$ {a_i} $x均为实数时,$m \geqslant (2 n - 1)$个随机测量值以大概率确保解的唯一性[63]

      (2) 当模型采用复值测量时,$m \geqslant (4 n - 4)$个随机测量值足以确保解的唯一性[64]

      从稳定性的角度来看,当实值测量中存在噪声干扰时,可以从$m \geqslant O\left( {n\log n} \right)$个随机测量值中实现稳定唯一性[65]表2对相位恢复问题解的唯一性进行了总结。

      表 2  相位恢复问题解的唯一性[66]

      Table 2.  Uniqueness of phase retrieval problem[66]

      ItemContent
      Fourier measurements1DNo uniqueness
      ≥2DUniqueness for real nonreducible signals with a twice oversampling rate
      General measurements Real signal Satisfying the complement property is necessary
      2N–1 random measurements guarantee uniqueness with high probability
      Real signal (noisy) Stable uniqueness requires nlogn measurements
      Complex signal 4n–4 generic measurements are sufficient
    • 有关相位恢复算法最早的研究可以追溯到1972年Gerchberg和Saxton提出的GS算法[29-30],该算法开创了交替投影迭代算法的先河。GS算法的原理如图2(a)所示,简单来说,初始相位的估计值通过正逆傅里叶变换在空域约束集和频域约束集之间交替映射,每次映射后计算的振幅分布都被替换为预先测量的振幅数据。随着迭代的进行,空域相位信息逐步从初始估计值往真实值方向收敛,最终可实现对 “丢失”相位信息的重建。尽管GS算法能够简单有效地实现相位恢复,但其本身也存在着明显的不足。交替投影迭代算法属于典型的贪婪型算法,由于在迭代过程中缺乏全局搜索机制,在求解非凸非线性的相位恢复问题时不可避免地会出现以下问题[67-69]

      图  2  单强度交替投影迭代算法流程图与原理图。(a) GS算法流程图;(b) HIO算法流程图;(c1)~(c3) ER、HIO、CHIO算法下一次迭代输入与当前迭代输出之间的关系曲线

      Figure 2.  Flow diagrams and schematic diagrams of single-intensity alternating projection iterative algorithms. (a) Flow diagram of GS algorithm; (b) Flow diagram of HIO algorithms; (c1)-(c3) The relationship between the input of the next iteration and the output of the current iteration for ER, HIO, CHIO respectively

      (1) 初值敏感:不同的迭代初始值将得到不同的收敛结果,最终恢复的收敛解依赖于迭代算法的初始值;

      (2) 局部最优:迭代易陷入停滞而无法达到精确的全局解,最终的收敛解为局部最优值;

      (3) 噪声敏感:当所测强度数据有较强的噪声干扰时,相位重建质量显著下降。

      为了改善GS算法的收敛性能,确保算法能够收敛到全局最优解,学者们通过优化算法的迭代机制提出了一系列改进算法。Soifer等人基于线性加权思想,在频域约束条件中添加线性加权因子,提出自适应附加(AA)算法[70],其频域约束条件表示为:

      $$ {F'} = \left[ {\beta \left| {\sqrt {{I_0}} } \right| + \left( {1 - \beta } \right)\left| F \right|} \right]{{\rm{e}}^{j \theta }} $$ (4)

      式中:β是取值区间为[0, 1]的加权常数,当β=1时AA算法退化为GS算法。AA算法有效避免了迭代停滞的问题,但是对于收敛速度的提升非常有限,收敛精度也有待提高。Fienup在GS算法的基础上,将空域的幅值约束条件替换为对迭代估计值进行修正,提出了误差下降(ER)算法[71-72],其空域估计值的更新方式为:

      $$ {f_{k + 1}}\left( {x,y} \right) = \left\{ {\begin{array}{*{20}{l}} {f_k'\left( {x,y} \right),}&{\left( {x,y} \right) \notin \gamma } \\ {0,}&{\left( {x,y} \right) \in \gamma } \end{array}} \right. $$ (5)

      式中:γ表示设定的空域约束,这里的空域约束可以是支持域约束(比如信号在应该为零的地方是非零的)或者非负性约束。ER算法具有更加广泛的适用性,对于仅已知频域幅度测量值的信号也能在空域内实现重构。随后,Fienup又在ER算法的基础上进一步加强空域约束的支撑性,引入负反馈机制,提出了混合输入输出(HIO)算法[54](原理如图2(b)所示),其空域估计值的更新方式可以描述为:

      $$ {f_{k + 1}}\left( {x,y} \right) = \left\{ {\begin{array}{*{20}{l}} {f_k'\left( {x,y} \right),}&{\left( {x,y} \right) \notin \gamma } \\ {{f_k}\left( {x,y} \right) - \beta f_k'\left( {x,y} \right){\text{,}}}&{\left( {x,y} \right) \in \gamma } \end{array}} \right. $$ (6)

      式中:β是一个取值区间为[0.5, 1]的常数,通常在取0.7时算法效果最佳。HIO算法对于GS算法的改进侧重于将上一次迭代的输入值与本次迭代的输出值相结合,有效降低了计算载荷,不变的是HIO算法依旧使用测量数据代替计算数据来实现迭代更新。HIO算法的负反馈机制进一步加快了算法的收敛速度,但是随着迭代的进行图像会出现振荡现象,因此仍然存在较大的收敛误差。为此,Fienup对HIO算法进行改进,引入分段评估函数解决了HIO算法估计值更新不连续的问题,该算法被称为连续混合输入输出(CHIO)算法[73],其空域估计值的更新方式表示为:

      $$ {f_{k + 1}}\left( {x,y} \right) = \left\{ {\begin{array}{*{20}{l}} {f_k'\left( {x,y} \right),}&{\alpha {f_k} \leqslant f_k'} \\ {{f_k}\left( {x,y} \right) - \dfrac{{1 - \alpha }}{\alpha }f_k'\left( {x,y} \right),}&{0 \leqslant f_k' \leqslant \alpha {f_k}} \\ {{f_k}\left( {x,y} \right) - \delta f_k'\left( {x,y} \right),}&{{\text{otherwise}}} \end{array}} \right. $$ (7)

      式中:αδ是取值区间为[0, 1]的常数。图2(c1)~(c3) 分别展示了ER、HIO以及CHIO算法的下一次迭代输入与当前次迭代输出之间的关系曲线。研究表明迭代过程中空域分布连续性的增强使得CHIO算法以更快的收敛速度实现了更高的收敛精度。Fienup的一系列改进继承了GS算法的框架结构,在改善算法收敛性能的同时拓展了算法的适用范围,为后续迭代相位恢复算法的发展奠定了理论基础。Elser和Luke等人在Fienup算法的基础上建立了更加完整的空域约束机制,分别提出了差分映射(DM)算法[74]和松弛平均交替映射(RAAR)算法[75],针对特定的问题类型选取合适的参数,算法的收敛效果得到了进一步提升。

      针对交替投影迭代算法的噪声敏感性,在实际应用中,人们发现将ER与HIO算法相结合可以有效提高算法的抗噪声性能,尤其对于高斯噪声具有良好的抑制作用[76-78]。2010年,Loh等人通过改变空域与频域约束类型,提出了改进的DM算法[79],并验证了该算法在泊松噪声的环境下能够实现鲁棒的相位重构。

      GS算法及其衍生的一系列改进算法都属于单强度类型的迭代算法,即仅借助单个频域面的强度测量实现空域迭代相位恢复。总体来说,这类算法解决相位恢复不适定性问题的能力比较薄弱,因此无法最大程度地提升算法的收敛性能。

    • 对于任何线性系统,相位恢复问题都可以视为求解复数域内一个多元高维病态方程组的过程,这也是导致解模糊性与低收敛性能的根本原因。可以获取输出平面内的多个强度测量图像,增加问题的先验信息,从而弱化相位恢复问题的不适定性,改善迭代相位恢复算法的收敛效果。

      1976年,Misell在研究电子显微成像时提出了Misell算法[80],通过两个离焦面图像的强度分布实现了波函数相位分布的重构。简单来说,光学系统记录两个远场平面(例如一组距离已知的离焦面和对焦面)的强度信息,两个远场平面之间的传播通过可引入适当离焦相位变化的孔径平面来实现。算法首先对孔径分布进行初始猜测,然后通过正逆傅里叶变换在两个远场平面之间交替迭代,用振幅的测量值代替计算出的振幅值,相位分布保持不变,直到满足收敛标准。表3直观展示了Misell算法的具体步骤。Misell算法的优点在于能够快速得到收敛解,因为在交替投影的迭代过程中误差是持续减小的。相比于GS算法,Misell算法更加灵活实用,也为改进迭代相位恢复算法开辟了新思路。后续基于Misell算法的完善通过多波长照明[81]、多重传播距离[82-83]、调整照明模式[84-86]等多种手段来实现。

      表 3  Misell算法流程表

      Table 3.  Flow chart for Misell algorithm

      Algorithm 1:Misell
      Input: Intensity measurements in the frequency domain without modulation ${\left| {{F_0}} \right|^2}$ and phase modulation factors ${t_m}$ with their corresponding intensity measurements ${\left| {{F_m}} \right|^2}\left( {1 \leqslant m \leqslant n} \right)$.
      Output:Spatial complex distribution f.
      S1. Select the initial Guess of spatial complex distribution without modulation f0 at random.
      S2. Add factor t1 to the phase of f0 , and obtain the spatial image after the first phase modulation f1.
      S3. Fourier transform f1 to obtain the complex distribution in the frequency domain G1.
      S4. Replace the amplitude of G1 with the measurement ${\left| {{F_1}} \right|^{}}$, and obtain the new complex distribution in the frequency domain $G_1'$.
      S5. Inverse Fourier transform $G_1'$ and substract the phase factor t1 , and obtain the new spatial image without modulation $f_1'$.
      S6. Add the factor $ {t}_{2} $ to the phase of $f_1'$, and obtain the spatial image after the second phase modulation f2 .
      S7. Repeat S2-S6 to cover all modulated images.
      S8. Repeat S2-S7 until the convergence conditions are satisfied, and output the final result.

      尽管Misell算法通过多强度测量实现了收敛速度和精度的提高,但并没有从根本上解决算法的局部收敛特性以及初值敏感性。当孔径分布的初始估计值极大地偏离全局最优解,或者远场平面的强度测量中出现噪声以及测量误差时,算法会在局部最优解处陷入停滞。Morris利用模拟退火算法的全局搜索机制获取Misell算法的最优化初值,有效降低了算法陷入局部最优解的可能[87]。但是模拟退火算法的引入提高了算法的复杂度,大大增加了计算成本。Xu等人对传统的模拟退火算法加以改进,在迭代时加入新的扰动模型,进一步增强了全局搜索的能力,该算法被称为改进的极快模拟退火(MVFSA)算法[88]。由于MVFSA-Misell算法综合了MVFSA强大的全局搜索能力和Misell算法优越的局部收敛速度,因此该算法能够在短时间内从不同的初始相位猜测中精确收敛到全局最优解。

      Misell算法通过轴向多图像之间的关联信息完成信号的重构。2004年,Rodenburg等人提出的叠层迭代引擎(PIE)算法[89-90]实现了径向平面上多强度图像的获取,其原理该算法的核心在于:对被测物体的不同区域逐一扫描获取衍射场强度,扫描过程中保证相邻区域间存在一定程度的交叠,从而造成衍射信息的数据冗余。空域孔径的重叠为相位恢复问题提供先验信息,有助于消除相位求解的模糊性,提升收敛速度和求解精度。图3(a)给出了PIE算法的成像示意图。实际应用中PIE算法的照明函数不易精确获取,为此Maiden等人提出了扩展叠层迭代引擎(ePIE)算法[91](具体步骤见图3(b)),将照明函数也作为未知量,与样本的复振幅分布一起迭代求解,降低了对照明函数先验信息的需求,提高了PIE算法的稳定性和适用性。2012年,mPIE算法被进一步地被推广到三维物体重建领域,实现了大尺寸薄样品的三维折射率成像[92]。近年来,基于PIE算法改进的rPIE和mPIE[93]算法相继被提出,算法的收敛性能和鲁棒性得到进一步提升。

      PIE算法中多强度约束信息的获取是在空域中进行的,在频域中也可以通过多角度照明的方式引入多强度信息的交叠,从而实现高分辨率的相位重建。2013年,Zheng等人提出了傅里叶叠层显微术(FPM)[94-95],使用LED阵列对样品提供多角度照明,同时保证不同角度的入射光照明下产生的相邻频域间满足一定的交叠率,就能够采用类似PIE的相位恢复方法对样品的频谱进行拼接,最终实现大视场、高空间分辨率的定量相位成像。FPM的实验平台配置和成像原理如图3(c)所示。笔者课题组围绕FPM开展了一系列研究工作,主要关注提升该技术的稳定性[96-98]和成像效率[99],探索成像通量的极限[100],并在此基础上推进其在三维成像[101]、药物筛选[102]、光学加密[103]以及遥感技术[104]等领域的应用。最近,笔者成功实践基于色彩传递的FPM彩色化方案[105],在几乎不牺牲重建精度的前提下实现了快速FPM全彩色成像,提高了病理切片筛查的效率。还将在深度学习和统计学领域广泛应用的交替方向乘子法(ADMM)引入到FPM技术中,进一步提高了成像的噪声鲁棒性[106]。可以参考近年来的多篇综述文章[107-112]详细了解FPM的技术原理和最新研究进展。

      图  3  径向多强度信息相位恢复原理示意图。(a) 叠层成像原理示意图;(b) ePIE算法流程图;(c) 傅里叶叠层显微术的成像原理与系统配置[95]

      Figure 3.  Schematic diagrams of radial multi-intensity phase retrieval. (a) Schematic diagram of ptychography; (b) Flow diagram of ePIE algorithm; (c) Principles and experimental setup of FPM[95]

      最后需要强调的是,在求解相位恢复问题时,交替投影迭代算法虽然求解速度慢,求解精度不高,但是它具有无可比拟的优势:实现简单且适应性强。无需迭代的直接型相位恢复方法虽然求解速度快,但往往需要预设严苛的实验条件或者特殊的测量方案,不能灵活适用于任意场景下的相位恢复问题。相比之下,交替投影迭代算法只需要获取单个或多个平面的强度信息,并设定投影集的支撑域约束,无需对具体的成像方案做出限定。

    • 1983年,Teague使用一个二阶椭圆偏微分方程揭示了相干光束纵向强度变化与相位之间的定量关系,即强度传输方程(TIE)[55,113],其数学形式表述如下:

      $$ - k\frac{{\partial I}}{{\partial {\textit{z}}}} = \nabla \cdot \left( {I\nabla \varphi } \right) = \nabla I \cdot \nabla \varphi + I{\nabla ^2}\varphi $$ (8)

      式中:$k = 2\pi /\lambda $表示波数;z为轴向离焦位置;Iφ分别为光波在z处平面上的强度和相位分布。在最初的论述中,TIE的推导是在近轴近似的条件下,通过将光波的复振幅表达式引入到亥姆霍兹方程中,然后将实部和虚部分离来实现的。后来的研究表明,TIE的推导方式并不是单一的,在近轴近似和小距离传播的限制下也可以从菲涅耳衍射公式[114]和坡印廷定理[115]中推导出TIE。根据TIE方法,需要采集至少两个离焦面的光强进行有限差分估计来获取方程左侧的光强轴向微分,通过直接测量可以获得待测平面的光强分布,之后就可以使用数值求解的方式得到待测平面的相位分布信息,不需要任何迭代过程。图4(a)、(b)展示了TIE定量相位成像技术的两种典型应用光路,图4(c1)~(c4)展示了使用TIE方法采集的三张离焦强度图像以及最终恢复的相位伪彩图像。

      TIE相位恢复的本质是一个边值问题,即在给定的边界条件下求解二阶偏微分方程,因此如何精确地获取TIE相位恢复的数值解是极具挑战性的研究内容。Teague的方法是引入辅助函数将原方程简化为两个泊松方程,通过格林函数得到TIE的解析解。此后,学者们又围绕这一问题提出了各种数值解法,包括多重网格法[116-118]、泽尼克多项式展开法[119-120]以及基于傅里叶变换的方法[121-122]等。然而,与这些方法相关的边界条件在实际中往往难以测量或作为先验信息来获取,因此TIE 方法受限于求解精确度始终没有得到广泛的应用。

      2014年,Zuo等人提出基于离散余弦变换(DCT)的TIE解法[123-124],其解的形式为:

      $$ \varphi = - k\nabla _D^{ - 2}{\nabla _D} \cdot \left[ {{I^{ - 1}}{\nabla _D}\nabla _D^{ - 2}\frac{{\partial I}}{{\partial {\textit{z}}}}} \right] $$ (9)

      式中:$\nabla _D^{ - 2}$${\nabla _D}$分别表示逆拉普拉斯算子和梯度算子。该解法将所有数据和运算都严格定义在一个与孔径函数有关的闭合矩形区域内,所有边界信号都被封闭在其中,因而不需要额外的检测方案来获取边界测量值。Huang等人通过迭代补偿机制将DCT解法进一步地拓展到任意形状的孔径,在提高适用性的同时获得更加精确的数值解[125]。TIE的相关数学基础的确立与其数值解法的成功实践成为了其后续应用快速发展的强劲推动力。

      TIE的最大亮点在于通过限制成像范围与近似条件,使得原本只能通过迭代算法求解的相位恢复问题能以“解析”的形式直接求解。相比于其他方法,具有计算简单、无需相位解包裹、无需复杂的光学系统等优势,在光学显微成像以及计量学等领域展现出广阔的应用前景。需要强调的是,文中将TIE归为无需先验知识的相位恢复算法的范畴,仅仅意味着该方法无需添加有关样品额外的先验信息来消除解模糊性,而上文提及的近似条件则是TIE实现的必要前提,不要对两个概念产生混淆。南京理工大学左超课题组在这一领域开展了系统性的研究工作,可以参考其撰写的教程论文[126]全面地归纳与总结该领域的理论基础与应用进展。此外,TIE方法的缺点同样十分突出,主要包括以下几方面:

      (1) 实现的前提条件严苛:须满足一定的近似条件(近轴近似,弱离焦)才能建立强度与相位信息之间的线性化关系,从而实现对相位恢复问题的直接求解;

      (2) 相位测量精度有待提高:与干涉测量方法对比,测量精度还存在一定的差距。在实际工程应用中实现TIE高精度相位测量的制约因素包括:对轴向离焦距离的精确校准,对边界条件的准确获取,对数值解的正确实现;

      (3) 应用场景有限:只在傍轴近似的菲涅耳衍射场成立,不适用于夫琅禾费衍射及大角度菲涅耳衍射,因此在X射线和电镜的相干衍射成像、大型望远镜的面形检测等领域无法发挥作用。

      图  4  TIE成像系统设置与成像结果[126]。(a)~(b) 基于4f系统和无透镜成像系统的TIE成像配置;(c1)~(c4) 使用TIE方法采集的三张MCF-7细胞离焦强度图像以及最终恢复的相位图像

      Figure 4.  Systematic setup and imaging result of TIE[126]. (a)-(b) TIE systematic setup based on 4f system and lensless imaging; (c1)-(c4) Three defocused intensity images of an MCF-7 cell acquired by TIE and its reconstructed phase image

    • 考虑一个线性系统,其时域响应函数满足因果律,即$f\left( t \right) = 0\left( {t < 0} \right)$,因此该系统的频域响应函数$F\left( \omega \right)$在上半平面应当是解析的。Kramers-Kronig(KK)关系基于因果性描述了线性系统频域响应函数的实部与虚部之间的相互约束[127],可以表示如下:

      $$ {Re} \left[ {F\left( \omega \right)} \right] = \frac{1}{\pi }{\text{p}}{\text{.v}}{\text{.}}\int_{ - \infty }^\infty {\frac{{{Im} \left[ {F\left( {\omega '} \right)} \right]}}{{\omega ' - \omega }}{\text{d}}\omega '} $$
      $$ {Re} \left[ {F\left( \omega \right)} \right] = \frac{1}{\pi }{\text{p}}{\text{.v}}{\text{.}}\int_{ - \infty }^\infty {\frac{{{Im} \left[ {F\left( {\omega '} \right)} \right]}}{{\omega ' - \omega }}{\text{d}}\omega '} $$ (10)

      式中:p.v.表示柯西主值积分。KK关系适用于很多系统的响应函数,目前在光学领域中的应用主要集中于折射率、消光系数等物理量的测量。KK关系是希尔伯特变换[128]的一个特例,对于具有因果关系的平方可积函数,希尔伯特变换便演化为KK关系。实际上,早在2005年,Gabriel等人就提出使用希尔伯特变换从全息图中直接求解透明物体的定量相位信息[129]。不过那时并没有涉及KK关系的概念,并且数理结论的推导也没有足够细致。

      传统离轴全息系统在参考光和目标光之间引入一个小角度来实现全息图的单次曝光采集,这导致系统的空间带宽积(SBP)无法得到有效利用。2019年,Park等人提出基于KK关系的高SBP离轴全息成像技术[130],该方法无需对样本施加假设性约束,通过引入解析性条件准确获取样本的复振幅信息,相比于传统方法获得四倍左右的SBP提升。利用KK关系直接求解全息图相位的实现揭示了这样一个核心思想:对于衍射受限的光学系统,尽管物函数的实部与虚部没有关联,但是在因果律的约束下,已知其中之一便能根据KK关系得到另外相应的部分。随后,Huang等人将KK关系与全息复用技术结合,突破了传统技术带宽利用率的极限[131]。这种显著的性能提升来自于对传感器空间带宽的充分利用:允许自相关谱与波前谱重叠,以便在单张全息图中记录多通道信息。

      Yang等人提出基于KK关系的合成孔径成像(KKSAI)技术[56],采用带限信号分析与瞳孔调制的方案从强度测量中计算相位信息。如图5所示,该方法使用四个在光瞳中心处相切的子孔径实现光瞳调制,并采集相对应的强度图像,然后基于KK关系对频谱进行移位拼接,恢复样本的全部频谱信息,最终通过傅里叶逆变换实现相位重构。与瞳孔调制FPM(PM-FPM)技术相比,该方法无需数据冗余,无需迭代,只需采集少量强度图像就能够获取更高精度的定量相位成像结果,在重建时间上也存在显著的优势。表4对KKSAI技术的具体流程进行了描述。然而该方法获取最佳的相位重建结果对于光瞳扫描的精度有着极高的要求,即覆盖整个光瞳范围的各孔径区域应当在光瞳中心处精确相切,不存在任何交叠部分。未来,KK关系有望与光学衍射层析成像、非相干孔径扫描等技术等相结合,进一步提升光学系统无迭代无约束的性能,为三维成像、数字重聚焦、高通量数字病理学以及生物细胞和组织的大规模研究提供便捷高效的工具。

      图  5  KKSAI的基本原理与可视化流程图

      Figure 5.  Principles and visual flowchart of KKSAI

      TIE和KK关系这两种相位恢复的直接解法均无需迭代操作,因此求解速度快是它们的共同优势。此外,两者在实现过程中均需要满足一定的限制条件,不同之处在于TIE的近似条件需要严格满足且往往难以精确校准,这也是导致TIE相位测量精度较低的部分原因;KK关系对数值孔径(NA)匹配的要求很高,但是实现简单,并且在某些情形下少量误差是可以接受的,最终的相位重建精度优于TIE。从分辨率角度来看,TIE对分辨率没有提升但在环形照明下可以提升至两倍,KK关系可以实现两倍分辨率但是不存在继续提升的空间。

      表 4  KKSAI技术流程表

      Table 4.  Flow chart for KKSAI technique

      Algorithm 2:KKSAI
      Input:Intensity measurements Ii (x, y) at each sub-aperture and reference plane wave ri (x, y). Their definitions are given respectively as:
      $ I_{i}(x, y)=\left|\mathcal{F}^{-1}\left\{S_{i}(u, v)\right\}\right|^{2}$
      $v_{1}(x, y)={\rm{e}}^{-{y\left(\omega_{1}-x+y_{1}-y\right)} }$
      where $ S_{i}(u, v)=S(u, v) D\left(u-u_{i}, v-v_{i}\right)$ denotes the subregion spectrum at the imaging plane cropped by the sub-aperture at position (ui, vi).
      Output: Pupil-limited object function S(u, v).
      S1. For the starting intensity I1, introduce a Hilbert kernel $ H_{1}(u, v)=-j \operatorname{sgn}\left(v_{1}\right) \cdot \operatorname{sgn}(v)$ and define an intermediate variable as:
      $ X_{1}=\ln \left\{\left[\tilde{s}_{1}(x, y)+r_{1}(x, y)\right] / r_{1}(x, y)\right\}$
      where $ \tilde{s}_{1}(x, y)$ denotes the shifted scattered field exiting from the object plane.
      S2. Impose KK relations to calculate the real part and imaginary part of the intermediate variable:
      $\operatorname{Re}\left\{X_{1}\right\}=\dfrac{1}{2} \ln \left[I_{1}(x, y) /\left|r_{1}(x, y)\right|^{2}\right]$
      $ \operatorname{Im}\left\{X_{1}\right\}=\mathcal{F}^{-1}\left\{\mathcal{F}\left\{\operatorname{Re}\left\{X_{1}\right\}\right\} \cdot H_{1}(u, v)\right\}$
      S3. Obtain the shifted form of the subregion:
      $S_{1}\left(u+u_{1}, v+v_{1}\right)=F\left\{{\rm{e}}^{\operatorname{Re}\left\{X_{1}\right\}+j \ln \left\{X_{1}\right\} } \cdot v_{1}(x, y)\right\}$
      move it back to the correct position to get Si (u, v).
      S4. Repeat above operations for other subregions.
      S5. Subsequently stitch all subregional complex fields into a full pupil plane spectrum S(u, v):
      $S(u, v)=\displaystyle\sum_{i=1}^{4} S_{i}(u, v) /\left[\displaystyle\sum_{i=1}^{4} D\left(u-u_{i}, v-v_{i}\right)+\varepsilon\right]$
      where $\varepsilon$ is a small constant to keep numerical stability in the zero-valued region.
    • 从理论角度来看,相位恢复在大多数情况下缺乏唯一的解决方案。即使该解决方案存在,也不能确保能够精确地找到它。近年来,在相位恢复算法领域掀起的研究热潮很大程度得益于光学领域各种新型成像技术的蓬勃发展,这一趋势已经逐步影响信号处理领域。该领域的研究人员利用现代优化理论的工具开发出一系列性能强大相位恢复方法,展现出了广阔的应用前景。

      其中一类解决相位恢复问题的方法是基于半正定规划的优化算法,通过将公式(1)中的二次方程组在更高的维度上重写为线性方程,将相位恢复问题转化为凸优化问题,然后使用梯度下降法等传统最优化算法直接迭代求解。2011年,Candes等人基于矩阵半正定规划开发出PhaseLift算法[34],使用向量提升的方法将原始相位恢复问题转化一个的低秩矩阵优化问题。具体来说,在公式(1)的基础上引入矩阵$X = xx^*$并做适当的变形,从而得到:

      $$ {y_i} = {\left| {\left\langle {{a_i},x} \right\rangle } \right|^2} = {x^ * }{a_i}a_i^ * x = {x^ * }{A_i}x = {\text{Tr}}\left( {{A_i}X} \right) $$ (11)

      式中:矩阵$ {A_i} = {a_i}a_i^ * $X均为秩1矩阵;${\text{Tr}}\left( \cdot \right)$表示矩阵的迹范数,定义为矩阵主对角线上的元素之和。显然矩阵X是半正定的,其所有特征值非负,记作$X \geqslant 0$。因此,相位恢复问题可以表述为:

      $$ \begin{split} &{\rm{find}} \;\;X \\ &{\rm{s.t.}}\;\; {y_i} = {\rm{Tr}}\left( {{A_i}X} \right),i = 1,2,\cdots m\\ &{\text{rank}}\left( X \right) = 1,{\text{ }}X \geqslant {\text{0}}{} \end{split} $$ (12)

      上式等价于求解以下的矩阵秩的最小化问题:

      $$ \begin{split} &\text{minimize rank} (X)\\ &{\rm{s.t.}} \;\;{y_i} = {\rm{Tr}}\left( {{A_i}X} \right),i = 1,2,\cdots m\\ &X \geqslant 0 \end{split} $$ (13)

      低秩优化问题在计算上一般是NP难的[132],其复杂度随矩阵维度的增加呈平方指数关系增长,采用全局优化方法去寻求最优解难以实现。针对这类问题,凸松弛策略是最为简单有效的解决方案,将原本难以求解的非凸问题转变为求解单个或者多个易于处理的凸问题,最终得到概率性的全局最优解。这里使用迹范数将相位恢复问题松弛为:

      $$ \begin{split} &\text{minimize Tr} (X) \\ &{\rm{s.t.}} \;\;{y_i} = {\rm{Tr}}\left( {{A_i}X} \right),i = 1,2,\cdots m\\ &X \geqslant 0 \end{split} $$ (14)

      PhaseLift算法则通过求解上述半正定规划问题从幅度测量值中实现信号的重构。研究结果表明,当信号的测量值服从高斯分布且测量次数满足一定的条件时,该方法能以较大概率准确重构信号,并且对于测量中出现的噪声具有鲁棒性。

      Waldspurger等人在PhaseLift算法的基础上加以改进提出PhaseCut算法[35],该算法首先将待恢复信号的幅值与相位分离,然后把相位信息作为未知分量进行优化,最后同样使用凸松弛策略将相位恢复问题转化为半正定规划问题求解。首先令$Ax = diag\left( y \right)u$,则相位向量$u \in {C^m}$中的每个元素都满足$\left| {{u_i}} \right| = 1,{\text{ }}i = 1,2,\cdots,m$。于是公式(2)转化为:

      $$ \begin{split} &{\rm{minimize}}\;\;\left\| {Ax - diag\left( y \right)u} \right\|_2^2 \\ &{\rm{s.t.}}\;\; \left| {{u_i}} \right| = 1,{\text{ }}i{\text{ = 1,2,}}\cdots ,m \end{split} $$ (15)

      式中:$diag\left( y \right)$表示由向量y的元素构成的对角矩阵。根据最小二乘法,使得目标函数最小的必要条件为$x = {A^ + }diag\left( y \right)u$,其中${A^ + }$为矩阵A的伪逆矩阵,因此上述问题可以改写为:

      $$ \begin{split} &{\rm{minimize}} \;\;\left\| {\left( {A{A^ + } - I} \right)diag\left( y \right)u} \right\|_2^2 \\ &{\rm{s.t.}}\;\; \left| {{u_i}} \right| = 1,{\text{ }}i= 1,2,\cdots ,m \end{split} $$ (16)

      对上式中的二次表达式进行变换,可以得到:

      $$ \begin{split} &{\rm{minimize}}\;\; {u^ * }Mu\\ &{\rm{s.t.}}\;\; \left| {{u_i}} \right| = 1,{\text{ }}i = 1,2,\cdots,m \end{split} $$ (17)

      式中:矩阵$M = diag\left( y \right)\left( {I - A{A^ + }} \right)diag\left( y \right)$是一个正定埃米尔特矩阵。令$U = u{u^ * }$,使用向量提升方法将问题转化为关于半正定矩阵U的优化问题:

      $$ \begin{split} &{\rm{minimize}}\;\; {\text{Tr}}\left( {UM} \right) \\ &{\rm{s.t.}} \;\;diag\left( U \right) = 1,{\text{ rank}}\left( U \right) = 1,{\text{ }}U \geqslant 0 \end{split} $$ (18)

      引入凸松弛策略,将秩约束条件忽略,得到PhaseCut算法最终需要解决的半正定规划问题:

      $$ \begin{split} &{\rm{minimize}}\;\; {\text{Tr}}\left( {UM} \right) \\ &{\rm{s.t.}}\;\; diag\left( U \right) = 1,{\text{ }}U \geqslant 0 \end{split} $$ (19)

      与PhaseLift算法相比,PhaseCut算法具有更高的恢复效率,且在噪声环境中更加稳定。图6(a1)、(a2)比较了两种算法在高斯噪声和图像扫描线等外部干扰存在的情形下的重建误差。相比于传统的迭代算法,基于半正定规划的凸松弛优化算法可以精确地恢复信号并提供收敛性保证,但由于其使用向量提升的方法将向量维度的问题提升到矩阵维度中加以解决,繁琐地表示未知的高维度提升矩阵大大增加了算法复杂度,因此不适用于大规模信号(如一维长信号或大尺度高分辨率图像)的相位恢复问题。

      图  6  基于优化理论的相位恢复算法的性能比较。(a1)~(a2) PhaseLift和PhaseCut算法在高斯噪声以及图像扫描线存在的情况下相位恢复的相对误差曲线[35];(b1) 使用谱方法和截断谱方法两种初始化方案对一维高斯信号进行相位重构的相对误差曲线[37];(b2) 使用WF和TWF算法对实值信号进行相位恢复的经验成功率曲线[37];(c1)~(c2) 使用基于强度损失函数和基于幅度损失函数的算法对实值和复值信号进行相位恢复的经验成功率曲线[38]

      Figure 6.  Performance comparison of phase retrieval algorithms based on optimization theory. (a1)-(a2) Relative error curves of phase retrieval for PhaseLift and PhaseCut under Gaussian noise and image scan-lines[35]; (b1) Relative error curves of spectral method and truncated spectral method for the phase retrieval of 1D Gaussian signal[37]; (b2) Empirical success rate curves of WF and TWF for the phase retrieval of real-valued signal[37]; (c1)-(c2) Empirical success rate curves of algorithms based on intensity/amplitude loss function for the phase retrieval of real-valued and complex-valued signals[38]

    • 最近提出的非凸优化相位恢复算法直接在原始信号的维度上进行迭代优化,避免了大规模数据处理的过程,显示出了以高效和准确的方式解决相位恢复问题的潜力。此外,相比于交替投影迭代算法,这类算法无需任何关于信号的先验信息,因而激发了人们对于相位恢复新的兴趣。非凸优化相位恢复算法通常可以分为两个步骤:

      (1) 初始化阶段:对于非凸的相位恢复问题,寻找到一个足够接近全局最优解的良好初值是优化算法所必需的。目前常见的初值获取方法包括谱方法[133-134]、正交法[135]、最大相关法[136]等。

      (2) 迭代优化阶段:该一阶段通常设计合理的梯度下降算法对初始估计值进行更新。为了提供更好的梯度下降方向以提高算法的收敛性能,对梯度施加截断、加权等额外操作也是必要的。常见的优化损失函数可根据相位恢复模型的不同分为以下两种:

      基于强度的最小二乘损失函数[36]

      $$ {\text{minimize }}f\left( {\textit{z}} \right): = \frac{1}{{2m}}{\sum\limits_{i = 1}^m {\left( {{y_i} - {{\left| {a_i^{\rm{H}}{\textit{z}}} \right|}^2}} \right)} ^2} $$ (20)

      基于幅度的最小二乘损失函数[38]

      $$ {\text{minimize }}f\left( {\textit{z}} \right): = \frac{1}{{2m}}{\sum\limits_{i = 1}^m {\left( {\sqrt {{y_i}} - \left| {a_i^{\rm{H}}{\textit{z}}} \right|} \right)} ^2} $$ (21)

      综上所述,使用非凸优化算法进行相位恢复的准确性依赖于良好的初始估计值和梯度下降方向,因此初始化策略的选择以及迭代阶段步长等参数的优化往往是繁琐的,算法的改进也集中在这两个维度进行。

    • 2015年,Candès使用一种基于Wirtinger导数的梯度下降算法从有限测量次数中实现了原始信号的精确重构,该算法被命名为维尔丁格流(WF)算法[36],为后续非凸优化相位恢复问题的研究开创了理论基础。WF算法求解相位恢复问题的具体步骤如下:

      (1)初始化阶段:使用谱方法进行初始化,根据采样向量$ {\alpha _i} $和观测值${y_i}$共同构造半正定埃米尔特矩阵:

      $$ Y = \frac{1}{m}\sum\limits_{i = 1}^m {{y_i}} {a_i}a_i^{\rm{H}} $$ (22)

      求出该矩阵的主特征向量$\tilde {\textit{z}}$,即最大特征值所对应的特征向量。最终的初始估计值在主特征向量的基础上进行适当缩放得到,即${{\textit{z}}_0} = \lambda \tilde {\textit{z}}$,其中缩放尺度设定为:

      $$ \lambda = \sqrt {n\frac{{\displaystyle\sum\nolimits_i {{y_i}} }}{{\displaystyle\sum\nolimits_i {{{\left\| {{a_i}} \right\|}^2}} }}} $$ (23)

      (2)迭代优化阶段:根据公式(20)定义强度损失函数的Wirtinger梯度为:

      $$ \nabla {f_{WF}}\left( {\textit{z}} \right) = \frac{1}{m}\sum\limits_{i = 1}^m {\left( {\left| {a_i^{\rm{H}}{\textit{z}}} \right| - {y_i}} \right)\left( {{a_i}a_i^{\rm{H}}} \right){\textit{z}}} $$ (24)

      将初始化阶段得到的良好估计值${{\textit{z}}_0}$作为迭代初值,使用如下梯度下降法进行优化更新:

      $$ {{\textit{z}}_{k + 1}} = {{\textit{z}}_k} - \dfrac{{{\mu _{k + 1}}}}{{{{\left\| {{{\textit{z}}_0}} \right\|}^2}}}\nabla {f_{WF}}\left( {{{\textit{z}}_k}} \right) $$ (25)

      式中:k表示迭代次数;µ为算法迭代的步长,其大小随着迭代过程的推进而发生变化。实验表明,步长的取值为${\mu _k} = \min \left( {1 - {{\rm{e}}^{ - k/{k_{_0}}}},{\mu _{\max }}} \right)$时算法的收敛效果最佳,其中${k_0}$${\mu _{\max }}$分别取330和0.4。WF算法在无噪声情况下可以确保收敛到全局最优解,并且在高斯噪声存在的情况下也能够稳定收敛到最大似然估计。

      在此基础上,截断维尔丁格流(TWF)[37]算法通过设计合理的截断规则,提供更精确的初始猜测和更好的梯度下降方向,从而降低了算法复杂度,进一步提高了算法的收敛性能。表5给出了TWF算法的具体流程。该算法在初始化阶段根据预先设定好的截断参数将数值过大的观测值丢弃,在迭代优化阶段引入截断规则集合有效滤除那些诱导更新朝着远离最优解方向进行的梯度分量。图6 (b1)展示了使用谱方法和截断谱方法两种初始化方案对一维高斯信号进行相位重构的相对经验误差,图6(b2)比较了WF和TWF算法的相位恢复成功率,不难发现TWF算法需要较少的测量值便能确保相位信息的精确重构。近年来,在WF和TWF算法的基础上又涌现出许多改进算法,例如可以快速实现稀疏信号重构的SWF算法[137],将截断规则分别与增量梯度下降方法和样本中值相结合提出的ITWF算法[138]和median-TWF算法[139]等。

      表 5  截断维尔丁格流算法流程表

      Table 5.  Flow chart for TWF algorithm

      Algorithm 3:Truncated Wirtinger Flow (TWF)
      Input:Intensity measurements yi and sampling vectors αi; trimming thresholds αy and truncation criteria sets ε1, ε2.
      Output:Iteration value zT.
      S1. Initialize estimated value ${{\textit{z}}_0} = \lambda \tilde {\textit{z}}$,where λ is defined by Eq. (23) and $\tilde {\textit{z}}$ is the leading eigenvector of the following matrix:
      $Y = \dfrac{1}{m}\displaystyle\sum\limits_{i = 1}^m { {y_i} } {a_i} a_i^{\rm{H}}{1_{\left\{ {\left| { {y_i} } \right| \leqslant \alpha _y^2\left( {\dfrac{1}{m}\displaystyle\sum\nolimits_{i = 1}^m { {y_i} } } \right)} \right\} } }$
      S2. Impose iterations on the initial estimated value by gradient update, and the update formula is:
      ${ {\textit{z} }_{k + 1} } = { {\textit{z} }_k} + \dfrac{ {2{\mu _k} } }{m}\displaystyle\sum\limits_{i = 1}^m {\dfrac{ { {y_i} - \left| {a_i^{\rm{H}}{ {\textit{z} }_k} } \right|} }{ {a_i^{\rm{H}}{ {\textit{z} }_k} } } } {\text{ } }{a_i}{1_{\varepsilon _1^i \cap \varepsilon _2^i} }$
      S3. Output the final result after completing iterations of predetermined number.
    • 以上都是基于强度损失函数延伸出来的一系列非凸优化算法。实际上对基于强度的损失函数做一定的泰勒逼近就可以得到基于幅度的损失函数,并且大量实验验证与数值结果表明,基于后者构建的优化算法具有更加优异的收敛性能。

      截断幅度流(TAF)算法[38]是一种典型的基于幅度损失函数的非凸优化算法,该算法的实现同样包含初始化和迭代优化两个阶段:

      (1) 初始化阶段:利用正交初始化方法得到较为精确的初始估计值。该方法的理论基础是高维随机向量之间的正交性,其核心思想是寻找一个与采样向量$\left\{ {{a_i}} \right\}$的子集最为正交的向量${{\textit{z}}_0}$来逼近真实解。具体方法是求解矩阵

      $$ Y = \frac{1}{{\left| {{{\bar I}_0}} \right|}}\sum\limits_{i \in {{\bar I}_0}} {\frac{{{a_i}a_i^{\rm{H}}}}{{{{\left\| {{a_i}} \right\|}^2}}}} $$ (26)

      的主特征向量$\tilde {\textit{z}}$,其中${\bar I_0}$表示${I_0}$$\left\{ {1,2,\cdots,m} \right\}$中的补集。初始估计值${{\textit{z}}_0}$基于主特征值的尺度放缩,即${{\textit{z}}_0} = \lambda \tilde {\textit{z}}$,其中缩放尺度定义为:

      $$ \lambda = \sqrt {\sum\nolimits_{i = 1}^m {{y_i}/m} } $$ (27)

      (2) 迭代优化阶段:类似于TWF算法,TAF算法在迭代过程中也添加了截断规则对初始估计值${{\textit{z}}_0}$进行更新。迭代更新方式表示为:

      $$ {{\textit{z}}_{k + 1}} = {{\textit{z}}_k} - \frac{\mu }{m}\sum\limits_{i \in {I_{k + 1}}} {\left( {a_i^{\rm{H}}{{\textit{z}}_t} - {\psi _i}\frac{{a_i^{\rm{H}}{{\textit{z}}_t}}}{{\left| {a_i^{\rm{H}}{{\textit{z}}_t}} \right|}}} \right)} {a_i} $$ (28)

      其中截断规则定义为:

      $$ {I_{k + 1}} = \left\{ {1 \leqslant i \leqslant m|\left| {a_i^{\rm{H}}{{\textit{z}}_k}} \right| \geqslant \left( {1/1 + \gamma } \right){\psi _i}} \right\} $$ (29)

      式中:$\gamma $为预先设定好的截断参数。该截断规则丢弃了模值过长的梯度分量,在一定程度上消除了错误的估计值带来的不良影响,因此初始估计值的选择更加精确,从而显著改善算法的性能。图6 (c1)、(c2)比较了TAF算法与基于强度损失函数的非凸优化算法的相位恢复成功率,结果表明,TAF算法所需的迭代次数和测量值数量均明显低于WF和TWF算法。

      在TAF算法框架的基础上,Wang等人提出更加适用于大规模数据信号恢复的随机截断幅度流(STAF)算法[140],该算法具有简单高效且扩展性强的特点,并且在有附加噪声存在时具有鲁棒性;重加权幅度流(RAF)算法[136]在初始化和迭代优化阶段以重加权的方式对不同梯度分量进行权重分配,不需要任何截断规则和对待恢复信号的额外假设,只需最小次数的观测就能精确求解全局最优值。值得注意的是,基于幅度的相位恢复目标函数由于模值的存在呈现出非光滑的特性,为了解决这一问题,通常会引入光滑性函数来替代原损失函数中的非光滑函数,这一类算法包括PG-SCG[141]、SSPR[141]、SAF[142]等。

    • 稀疏性可以作为一种特殊的先验知识添加到相位恢复问题中,通过施加约束或使用正则化手段将对原始信号的搜索限制在稀疏向量集上,从而使用少量的数据就能够近乎完美地恢复原始信号,有效提高了算法的性能与稳定性。事实上,所感兴趣的对象在某些已知的表达中是稀疏的,可以简单地如下表示为:

      $$ x = \psi \alpha $$ (30)

      式中:x表示所感兴趣的对象;$\psi $$\alpha $分别表示表达矩阵(稀疏基)和稀疏向量,“稀疏”的含义可以理解为向量中仅包含少量非零元素。稀疏性先验知识与相位恢复算法的融合有多种形式。

    • 理论上,任何一种交替投影算法都可以在空域施加稀疏性约束,从而提高相位重建的性能。例如在Fienup型迭代算法的基础上,Mukherjee等人通过求解稀疏编码问题并保留K个最大的稀疏系数来施加稀疏性约束,用来替换原有算法中的空域估计值修正,有效抑制了相位重建过程中出现的自相关伪影和背景噪声,这一算法被命名为Max-K算法[57]。Stefan等人基于剪切波变换稀疏约束,使用RAAR算法实现了菲涅耳域的相位恢复[143]。Yang等人提出将梯度稀疏性融合到HIO算法中,该方法能够从部分测量数据中有效重建纯相位物体[144]。大量实验结果表明,基于稀疏性的交替投影算法无论是重建质量还是收敛速度均明显优于其对应的传统算法。

      另一种引入稀疏性先验的方式是贪婪策略,其中较为经典的相位恢复方案是贪婪相位检索(GESPAR)算法[58]。该算法在迭代中动态更新信号的支撑信息,同时在支撑域和非支撑域之间交替执行快速的局部搜索,使用阻尼高斯-牛顿算法最小化幅值误差函数,寻找在当前支撑域约束下对应于测量值的稀疏性向量。GESPAR算法的计算复杂度较小,运行效率更高,适用于大规模信号重建的场景,此外,对于噪声和其他不确定干扰因素具有鲁棒性。图7比较了GESPAR算法与其他相位恢复算法的计算复杂度以及相位重构成功率,表6提供了算法运行时间的对比结果。

      基于半正定规划的优化算法也可以考虑进行适当修改与信号的稀疏性先验知识合并。稀疏性半正定规划相位恢复的最初工作来自于亚波长成像领域,这项工作中提出的二次压缩感知(QCS)方法[145]可以在图像稀疏的前提下将部分相干光照明的亚波长图像从其远场模糊图像中恢复出来,并且可以进一步地推广,用来解决从部分相关测量中恢复信号的其他问题。该方法在公式(12)所示的矩阵秩最小化问题的基础上施加稀疏性约束。当信号向量x稀疏时,矩阵X = xx* 也应当是包含少量非零行的稀疏矩阵,矩阵X行向量的l2范数构成的向量同样具有稀疏性,对该向量的l1范数施加一定的约束条件,待解决的最优化问题转化为:

      图  7  GESPAR算法与半正定规划和稀疏性Fienup算法的比较[58]。(a1)~(a2) 对不同矢量长度信号进行相位恢复的算法复杂度对比;(b)相位恢复成功率对比曲线

      Figure 7.  Comparison between GESPAR, SDP and Sparse-Fienup algorithm[58]. (a1)-(a2) Comparison of algorithm complexity for the phase recovery of signals with different vector lengths; (b) Comparison curves of phase recovery successful rate

      $$ \begin{split} &\text{minimize rank} (X) \\ &{\rm{s.t.}} \left| {{\text{Tr}}\left( {{A_i}X} \right) - {y_i}} \right| < \varepsilon ,{\text{ }}i = 1,2,\cdots m\\ &\quad X \geqslant {\text{0}}\\ &\quad \sum\nolimits_j {{{\left( {\sum\nolimits_i {{{\left| {{X_{ji}}} \right|}^2}} } \right)}^{1/2}}} \leqslant \eta \end{split} $$ (31)

      式中:εη分别为噪声参数和稀疏性约束参数。随后采用迭代算法搜索与测量值和稀疏性先验信息一致的低秩矩阵$\hat X$,最后使用奇异值分解获得$\hat X$的最佳秩1矩阵近似值,待恢复信号x便可以简单求出。

      此后,许多类似的方法被提出,例如Ohlsson等人在PhaseLift算法中加入额外的l1范数正则化项来诱导稀疏性,开发出基于提升压缩(CPRL)的相位恢复算法[146],该算法需解决如下凸优化问题:

      $$ \begin{split} &{\rm{minimize}}\,\, {\text{Tr}}\left( X \right) + \lambda {\left\| X \right\|_1}\\ &{\rm{s.t.}} \left| {{\text{Tr}}\left( {{A_i}X} \right) - {y_i}} \right| < \varepsilon ,{\text{ }}i = 1,2,\cdots m\\ &\quad X \geqslant 0 \end{split} $$ (32)

      表 6  运行时间对比表[58]

      Table 6.  Comparison of runtime performance[58]

      SDP/sSparse-Fienup/sGESPAR/s
      s=31.320.090.04
      s=51.780.120.05
      s=83.850.500.06
    • 作为近年来备受关注的计算框架,深度学习在处理复杂数据集方面具有显著的优势,因而逐渐走入相位恢复领域研究者的视野,并且取得了一系列令人瞩目的开创性研究成果。实际上,深度学习同样包含稀疏性先验知识,它认为问题本身就是稀疏的,可以将输入信息从高维数据空间投影到低维认知空间进行处理。相比于传统基于稀疏性先验的相位恢复算法“先正向建模,再求解逆问题”的思路,深度学习方法规避了这种程序严格的数据表达过程,直接通过高维度特征拟合实现图像与信息的提取与重建[147]

      图8(a)展示了基于深度学习的相位恢复的基本模型。假设待恢复的对象为 f,使用照明操作符Hi描述光与物体的相互作用,即从光源出射的光线经过光学元件的调制后照射在物体上,使用采集操作符Hc表示光线在光学系统剩余部分的传播,即经过物体的散射光在探测器上成像。探测器采集输出的原始强度图像为g,考虑到光学测量中信号的泊松统计特性以及加性噪声的影响,成像模型可以描述为$g = \mathcal{P}\left( {Hf} \right) + \varGamma$,其中$ \mathcal{P} $表示信号光子到达的泊松随机过程,$\varGamma$为加性噪声。进一步地,相位恢复可以描述为最小化Tikhonov泛函的数学问题,表达式如下:

      $$ f' = \mathop {\arg \min }\limits_f \left\{ {\left\| {Hf - g} \right\|_2^2 + \alpha \phi \left( f \right)} \right\} $$ (33)

      其中的二次范数项使用最小二乘法将测量值与假设对象的正向模型相匹配,正则化项通过引入有关测量对象的先验信息弱化正向算子的不适定性。目前应用比较广泛的正则化项基于稀疏性先验知识,一方面对于大部分信号具有适用性,另一方面在应对原始图像被噪声破坏的情形时更有优势。

      图  8  基于深度学习的相位恢复方法。(a) 使用深度学习解决相位恢复问题的物理模型;(b)使用深度学习进行单帧无透镜相位恢复[59];(c)使用深度学习仅采集少量图像实现快速FPM成像[61];(d) prDeep算法中的去噪卷积神经网络的框架 [62]

      Figure 8.  Phase retrieval based on deep learning. (a) Physical model of deep learning solving phase retrieval problems; (b) Single-frame lensless phase retrieval using deep learning[59]; (c) Fast FPM imaging with a few images using deep learning[61]; (d) Basic framework of DnCNN in prDeep algorithm[62]

      现有的基于深度学习的相位恢复方法大致可以分为两类[148],一类是使用深度神经网络(DNN)来计算相位恢复问题的反函数,即直接从可用测量值或通过简单的基于模型的反演方法获得的初始估计值来重建未知图像,可以实现比传统方法更快速的非迭代相位恢复。这一类方法已经在同轴全息[59]图8(b))、无透镜成像[60]以及FPM[61]图8(c))等场景中得到广泛应用;另一类是将DNN作为迭代相位恢复过程中的去噪器,可以提高算法运行的稳定性。2018年,Metzler等人将所示的去噪卷积神经网络(DnCNN)(图8(d))嵌入到RED正则化框架中,提出prDeep算法[62]来解决散射介质成像中的相位恢复问题。研究表明,该算法适用于各种测量矩阵模型,并且对泊松噪声具有鲁棒性。

      深度学习为相位恢复算法的进一步发展注入了新的活力,但是其面临的挑战也是艰巨的。现有的方法主要专注于为有限的算法和应用场景设计神经网络,这些工作由于无法处理多种类型的模型以及更为复杂的相位恢复算法而受到限制。此外,深度学习在各领域应用普遍存在的问题包括训练数据的获取成本高、训练后的网络泛化能力不足以及可解释性不明晰等。针对这些迫切需要解决的问题,研究人员们也进行了深刻的探索与讨论。物理模型驱动是当前深度学习发展的一个重要方向,即在深度学习中嵌入特征规则先验,代替单一的纯数据驱动式学习[147]。Wang等人将传统的深度神经网络与描述图像生成过程的物理模型相结合,提出了一种物理增强型的神经网络架构PhysenNet[149],无需获取并标记大规模训练数据,通过神经网络和物理模型之间的相互作用就能自动优化网络并最终生成待测对象的相位信息。Naimipour等人在物理模型驱动数据的深度展开技术基础上提出UPR混合框架[150],该网络框架结合基于物理模型算法的可解释性和深度神经网络的表达能力,在解决具有边界复杂性的非凸相位恢复问题方面显示出巨大的潜力。

    • 光学成像中的相位恢复起源于从实验测量中恢复图像的尝试。需要强调的是,与数值模拟仿真不同,实验条件下的相位恢复还面临着一些额外的挑战。最为显著的差别在于实际应用中的强度测量值被各种类型的噪声破坏,例如光子散粒噪声、探测器的读出噪声以及光学系统内部反射产生的杂散光,相位重构的准确度无法得到保证。尽管这些问题增加了相位恢复算法的复杂度,但是光学成像中成功的相位恢复已得到广泛实现。由于篇幅所限,这里只简单介绍其中四项具有代表性的应用,不再尽数罗列。

      图  9  几种典型的相位恢复应用场景。(a) 使用X射线对晶体微观结构进行衍射成像;(b) 采用4f系统配置的基于相位恢复的光学加密系统;(c) 天文观测中的自适应光学系统;(d) 大口径光学镜面的面型检测

      Figure 9.  Several typical application scenarios of phase retrieval. (a) X-ray diffraction imaging of crystal microstructure; (b) Optical encryption system based on phase retrieval with a 4f setup; (c) Adaptive optics system in astronomical observation; (d) Surface shape detection of large-aperture optical mirrors

    • X射线的波长远小于可见光,具有非常高的分辨率,因此可以对微小的晶体结构等进行纳米级成像(图9(a))[151-152]。尽管微观物体本身无法被直接测量,但是借助衍射的“放大”效应可以将其频谱拓展到易测量的宏观尺度。待测物体与其夫琅禾费衍射场之间存在傅里叶变换关系,衍射场强度即对应物体的频谱。衍射场强度同样缺失相位信息,若要从宏观频谱中将物体原本的微观结构重构出来,求解丢失的相位信息是不可或缺的,因此X射线衍射成像也可以归结为典型的相位恢复问题。此外,相干衍射成像(CDI)的媒介涵盖电磁波的多个波段,衍射电子显微镜[153-154]以及可见光范围内的叠层成像[155-156]、傅里叶叠层显微成像[94-95]、无透镜成像[157-158]等新型成像技术均采用了相似的原理。

    • 与传统的数字加密相比,光学加密系统将信息编码到多个光学参量(振幅或相位)之中,使其通过一系列光学元器件(透镜、偏振片等),实现诸如菲涅耳衍射、傅里叶变换等数学变换,从而达到数据加密或隐藏的目的[159]。光学加密过程通常是快速的,这得益于光学信息处理对数据的并行计算能力。在普通光学系统中,不同平面上的相位通常是相互关联的,根据某一已知平面可以推导出另一个平面的相位分布。然而在光学加密系统中,由于随机相位调制,任意两个平面之间的相位关系被破坏,因此获取各平面的相位信息变得至关重要。最初的基于相位恢复的光学加密方案[160]采用4f系统配置,如图9(b)所示,其中输出平面的图像通过投影算法被编码为傅里叶平面中的相位掩码,位于输入平面的另一个相位掩码保持固定。此后,基于级联相位恢复算法的光学加密方案[161]被提出,在每次迭代中两个相位掩码同时调整,提高了算法运行速度,有效增强了加密系统的安全性。基于相位恢复的光学加密系统也由通常使用的4f系统扩展到菲涅耳域[162-163]甚至线性正则域[164-165]。此外,多图加密、彩色图像加密等方法的研究进一步挖掘了相位恢复在光学加密领域的应用前景,使得该领域始终焕发着蓬勃的生机,这一部分的内容可以参考综述文章[166-167]

    • 自适应光学系统在天文观测(图9(c))、空间探测以及大气激光传输等领域具有广泛的应用价值[168-170]。在这些应用场景中,大气湍流在光学系统的观测光路上引入随机的相位扰动,使得成像质量和稳定性大大降低。自适应光学通过实时探测畸变的波前并引入像差补偿机制来克服大气湍流对成像的不良影响。目前应用较为成熟的自适应光学技术按照有无波前传感可以分为两类:传统的自适应光学系统使用波前传感器实现畸变波前的测量,往往造价昂贵、系统复杂,并且无法解决振幅大范围波动情形下光束波前产生的不连续性[171],因此难以大规模推广应用;无波前传感的自适应光学系统无法直接探测波前相位信息,因而借助一系列优化算法[172-175]不断搜索使目标函数达到最优的控制向量,驱动波前校正的实现。就相位畸变的测量而言,两者本质上都是应用非干涉的定量相位测量方法来实现的,分别对应于波前传感和相位恢复两类技术。

    • 随着空间探测技术的发展,人们对大口径高精度光学镜面的需求也日益增长。因此在加工过程中需要对镜面面形进行检测,以便根据测量结果及时修正加工方案,确保达到最终所要求的面形精度。传统的干涉测量方法对振动等外界干扰敏感,且无法解决动态范围误差过大引发的相位缠绕问题[176]。哈特曼波前传感器虽然不易受振动影响,但横向测量分辨率较低,无法识别镜面上的局部误差。相位恢复已在大型空间望远镜的装配检测方面取得成功的应用,比较典型的是哈勃望远镜的像差诊断与校正[177],美国航空航天局(NASA)拼接式大型天文望远镜JWST的在轨装调[178]等。如图9(d)所示,基于相位恢复的光学表面检测方法[179-180]利用同一束测量光产生的衍射强度图像,通过迭代优化计算来恢复入射波面相位,从而获取镜面面形信息。该方法具有光路简单、抗振动和干扰能力强、精度较高等特点,可以较好地适应大型光学镜面面形在位测量的需求。

    • 文中综述了相位恢复算法的原理、发展历程以及在光学领域的应用。所述的相位恢复算法各有所长,从求解精度来看,基于优化理论的算法能够提供相位信息精确恢复的保证;从求解速度来看,TIE、KK关系等直接解法无需迭代操作,只涉及代数计算,因此速度上存在优势;从实验条件来看,稀疏信号处理方法只适用于具备稀疏性的实验对象,TIE方法的直接求解需满足若干近似条件,而传统的迭代算法不需要预设严苛的实验条件和特殊的测量方案,因而能灵活适用于任意场景下的相位恢复问题。

      尽管各类相位恢复技术已经成功应用于许多领域,但是它们在科学和工程角度仍有巨大的研究价值和发展空间,未来的研究应当针对各类算法的不足加以改进。相位恢复算法未来的发展趋势可以简单地概括如下。

      (1) 更优异的收敛性能和噪声鲁棒性:相位恢复问题的病态性以及实际应用场景中的噪声干扰使得算法易陷入局部最优解,甚至导致相位重构的失败。对于迭代型以及基于最优化理论的相位恢复算法,可以通过开发合适的初始化策略以及优化算法结构的设计来弱化相位恢复问题的不适定性,确保算法在存在噪声干扰时仍能稳定收敛到全局最优解,同时又不过多地增加算法的复杂度。

      (2) 恢复更复杂物体相位信息的能力:最初相位恢复技术的发展是为了满足人们观察透明生物细胞的需求。如今,“看见”的目标早已能够轻易实现,但是从更长远的角度来看,人们希望对高度复杂的生物分子结构直接成像,从而可以从分子尺度把握和理解生命科学,破译生物分子间相互作用的微观机理,甚至通过实时成像来追踪它们随时间推进产生的结构演变。要实现这样的愿景,显然还有很长的路要走,不过多种相位恢复方法的同时使用或有机结合也许能够为更复杂物体的相位重构,尤其是低频相位信息的重构提供新的思路。例如,将TIE快速计算得出的相位信息作为迭代算法的初始估计值,替代原本预设的随机猜测,既可以加速迭代算法的收敛,同时克服TIE在精度方面的不足;深度学习技术能够借助大量训练数据直接建立TIE输入(离焦强度图像)与理想输出(定量相位分布)之间复杂的非线性关系,从而巧妙地规避传统关系的限制。

      (3) 多目标多任务集成的兼容性:复杂的相位重构场景不仅体现在物体本身的高度复杂性,更重要的在于人们期待获取“一石二鸟”甚至“多鸟”的额外价值。综合解决多问题融合的复杂情形需要进一步挖掘样品的隐含信息,从而多维度地恢复更多内容。例如,FPM技术的优雅便体现在一次性实现三个维度的目标:相位恢复、合成孔径(实现更高的空间分辨率)和上采样(实现更高的像素分辨率)。此外,在相位恢复的同时迭代更新光学系统的光瞳函数实现了系统像差的数字矫正。

      最后需要指明的是,相位恢复并非是完全独立的技术领域,未来新型光学器件、数学与算法工具以及计算性能的提升有助于突破现有算法的固有瓶颈,必将极快地推动相位恢复技术的巨大进步,以上目标的实现便不再遥远。

参考文献 (180)

目录

    /

    返回文章
    返回