-
相位恢复问题的目标是从信号强度或幅度测量值中恢复原始信号及其相位。使用一维离散向量
$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
Algorithm Pros Cons Alternating projections[29,54] Simple implementation
Wide applicabilityTime-consuming
Low precisionNon-
priorTIE[55] Deterministic solution
Simple computation
Require no phase unwrapping and complex optical systemsHarsh preconditions
Lower precision compared with interferometry
Limited applicationsKK relations[56] Reconstruction with high quality and high speed
Require less data collectionRequire 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 estimateSparsity-
basedGeneral[57-58] Reconstruction with high quality and stability
Low computational costApplicable 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对相位恢复问题解的唯一性进行了总结。Item Content Fourier measurements 1D No uniqueness ≥2D Uniqueness 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 probabilityReal 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的技术原理和最新研究进展。
最后需要强调的是,在求解相位恢复问题时,交替投影迭代算法虽然求解速度慢,求解精度不高,但是它具有无可比拟的优势:实现简单且适应性强。无需迭代的直接型相位恢复方法虽然求解速度快,但往往需要预设严苛的实验条件或者特殊的测量方案,不能灵活适用于任意场景下的相位恢复问题。相比之下,交替投影迭代算法只需要获取单个或多个平面的强度信息,并设定投影集的支撑域约束,无需对具体的成像方案做出限定。
-
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关系有望与光学衍射层析成像、非相干孔径扫描等技术等相结合,进一步提升光学系统无迭代无约束的性能,为三维成像、数字重聚焦、高通量数字病理学以及生物细胞和组织的大规模研究提供便捷高效的工具。
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) -
作为近年来备受关注的计算框架,深度学习在处理复杂数据集方面具有显著的优势,因而逐渐走入相位恢复领域研究者的视野,并且取得了一系列令人瞩目的开创性研究成果。实际上,深度学习同样包含稀疏性先验知识,它认为问题本身就是稀疏的,可以将输入信息从高维数据空间投影到低维认知空间进行处理。相比于传统基于稀疏性先验的相位恢复算法“先正向建模,再求解逆问题”的思路,深度学习方法规避了这种程序严格的数据表达过程,直接通过高维度特征拟合实现图像与信息的提取与重建[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],该网络框架结合基于物理模型算法的可解释性和深度神经网络的表达能力,在解决具有边界复杂性的非凸相位恢复问题方面显示出巨大的潜力。
Phase retrieval algorithms: principles, developments and applications (invited)
-
摘要: 研究表明,由于相位比振幅包含更多关于场的信息,因此相位测量在现代科学和工程的诸多分支中始终是研究的热点问题。在可见的电磁波范围内,相位信息很难通过现有的光电探测器直接采集获取。相位恢复技术提供了一种从捕获的强度信息中将相位信息“计算”出来的有效手段,并已成功应用于天文观测、生物医学成像和数字信号复原等多个科学领域。算法是相位恢复技术的核心,也是该技术发展和应用的关键。文中结合物理学原理和信号处理方法对相位恢复算法的基本原理进行阐述,综述了各类相位恢复算法的发展历程及其优缺点,并简单概述了相位恢复算法在光学领域的典型应用,最终指明其面临的挑战和未来的发展趋势:更优异的收敛性能和噪声鲁棒性、恢复更复杂物体相位信息的能力、多目标多任务集成的兼容性。Abstract: Because the phase contains more information about the field in contrast to the amplitude, phase measurement has always been a hot topic in many branches of modern science and engineering. Within the visible range of electromagnetic wave, it is quite difficult to directly obtain phase information by the existing photodetectors. Phase retrieval provides an effective method to “figure out” the phase information from the captured intensity information, and has achieved successful applications in several scientific fields including astronomical observation, biomedical imaging and digital signal restoration. Algorithm is not only the core of phase retrieval, but is also the key to its development and applications. This paper demonstrates the basic principles of phase retrieval algorithms in combination with physical principles and signal processing methods, summarizes the development of various kinds of algorithms as well as their advantages and disadvantages, and briefly lists some typical applications in the field of optics. Finally, the challenges are pointed out, and the future development directions are described as: better convergence performance and noise robustness, phase-retrieval ability for more complex objects, compatibility for integration of multiple objectives and tasks.
-
Key words:
- phase retrieval /
- computational imaging /
- signal processing /
- optimization theory
-
图 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
图 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
图 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
Algorithm Pros Cons Alternating projections[29,54] Simple implementation
Wide applicabilityTime-consuming
Low precisionNon-
priorTIE[55] Deterministic solution
Simple computation
Require no phase unwrapping and complex optical systemsHarsh preconditions
Lower precision compared with interferometry
Limited applicationsKK relations[56] Reconstruction with high quality and high speed
Require less data collectionRequire 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 estimateSparsity-
basedGeneral[57-58] Reconstruction with high quality and stability
Low computational costApplicable to objects with sparsity only Deep learning[59-62] Powerful expressive ability Lack of interpretability
Unable to deal with multiple models and applicationsItem Content Fourier measurements 1D No uniqueness ≥2D Uniqueness 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 probabilityReal signal (noisy) Stable uniqueness requires nlogn measurements Complex signal 4n–4 generic measurements are sufficient 表 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.表 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.表 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. -
[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