CN110807822A - 基于Wirtinger Flow算法的散斑相关成像方法及装置 - Google Patents

基于Wirtinger Flow算法的散斑相关成像方法及装置 Download PDF

Info

Publication number
CN110807822A
CN110807822A CN201910972577.1A CN201910972577A CN110807822A CN 110807822 A CN110807822 A CN 110807822A CN 201910972577 A CN201910972577 A CN 201910972577A CN 110807822 A CN110807822 A CN 110807822A
Authority
CN
China
Prior art keywords
target
speckle
cost function
power spectrum
image
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201910972577.1A
Other languages
English (en)
Other versions
CN110807822B (zh
Inventor
边丽蘅
李一鸣
张军
曹先彬
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Institute of Technology BIT
Original Assignee
Beijing Institute of Technology BIT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Institute of Technology BIT filed Critical Beijing Institute of Technology BIT
Priority to CN201910972577.1A priority Critical patent/CN110807822B/zh
Publication of CN110807822A publication Critical patent/CN110807822A/zh
Application granted granted Critical
Publication of CN110807822B publication Critical patent/CN110807822B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20048Transform domain processing
    • G06T2207/20056Discrete and fast Fourier transform, [DFT, FFT]

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明公开了一种基于Wirtinger Flow算法的散斑相关成像方法及装置,其中,该方法包括:S1,获取目标散斑图像,根据维纳‑辛钦定理对目标散斑图像进行自相关的傅里叶变换,得到目标功率谱;S2,通过Wirtinger Flow算法建立目标图像和目标功率谱的代价函数,通过优化算法对代价函数进行优化并求解代价函数的最优解,根据最优解重建目标图像。该方法可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。

Description

基于Wirtinger Flow算法的散斑相关成像方法及装置
技术领域
本发明涉及抗散射成像技术领域,特别涉及一种基于Wirtinger Flow算法的散斑相关成像方法及装置。
背景技术
携带目标信息的光波在透射过强散射介质如生物组织、烟雾、鸡蛋薄膜等会发生散射,使得光波所携带的目标信息被重新“编码”,探测器所接收到的是随机散斑图像,无法分辨出目标的轮廓和细节。
散斑相关成像技术是基于强散射介质固有的光学记忆效应性质所发展而来的一种可以在传输路径上存在强散射介质的情况下进行成像的技术。光学记忆效应最早是由以色列科学家I.Freund等人首次提出的概念,即光波透过散射介质后,当小范围的改变光波的入射角时,不同入射角得到的散斑场之间存在较强的相关性,可以近似看作散斑场随着入射角度的变化而移动的一种现象。基于光学记忆效应的散射成像研究就此展开。2012年,意大利科学家J.Bertolotti等人在Nature上首次提出利用散斑相关法实现透过散射介质成像,借助在光学记忆效应范围内散斑场之间的相关性,结合交替投影的相位恢复算法实现了随机散射介质的成像。该方法成功实现了非侵入式散射成像,但是由于要在光学记忆效应范围内进行扫描,整个过程需要耗费几十分钟,因此无法应用在实际情况中。为解决这个问题,2014年,以色列科学家O.Katz等人在J.Bertolotti研究的基础上提出了一种基于单帧散斑相关的无透镜成像方法,该方法不仅完全避免了原有成像系统需要进行扫描的不足,而且由于系统无透镜,也避免了由透镜带来的像差对成像质量的影响。然而,该方法依然存在需要改进的地方,即该方法在相位恢复的算法中应用的是Fienup算法,该算法需要目标的先验信息,即目标的尺寸。2015年,国内学者代伟佳在他的研究中,对于散斑图像的相位恢复使用的是通用近似信息传递相位恢复算法,该算法优势在于不需要目标尺寸的先验信息,转而需要较易实现的估计目标稀疏度的先验信息。虽然有所改进,但仍旧没有避免需要先验信息这一不足,在重建过程中对测量噪声和系统畸变的鲁棒性上呈现的效果一般。
如何提高透过强散射介质实现高分辨率成像,是光学成像领域亟待解决的一项重要问题,其在生物医学成像、海洋环境探测、公共安全等领域具有重要的应用价值。
发明内容
本发明旨在至少在一定程度上解决相关技术中的技术问题之一。
为此,本发明的一个目的在于提出一种基于Wirtinger Flow算法的散斑相关成像方法,该方法可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
本发明的另一个目的在于提出一种基于Wirtinger Flow算法的散斑相关成像装置。
为达到上述目的,本发明一方面实施例提出了一种基于Wirtinger Flow算法的散斑相关成像方法,包括:
S1,获取目标散斑图像,根据维纳-辛钦定理对所述目标散斑图像进行自相关的傅里叶变换,得到目标功率谱;
S2,通过Wirtinger Flow算法建立目标图像和所述目标功率谱的代价函数,通过优化算法对所述代价函数进行优化并求解所述代价函数的最优解,根据所述最优解重建所述目标图像。
本发明实施例的基于Wirtinger Flow算法的散斑相关成像方法,通过搭建散斑相关成像系统,在光学记忆效应范围内照射目标,获取目标的散斑图像;对获取的散斑图像做相关计算;根据维纳-辛钦定理,目标的功率谱是散斑图像的自相关的傅里叶变换,通过对散斑图像进行自相关运算并进行傅里叶变换,可得到目标的功率谱;利用Wirtinger Flow相关算法建立目标功率谱与目标图像相关的代价函数,并求取函数对目标图像的类梯度并利用类似梯度下降的方法进行优化,最终得到最优解。由此,可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
另外,根据本发明上述实施例的基于Wirtinger Flow算法的散斑相关成像方法还可以具有以下附加的技术特征:
进一步地,在本发明的一个实施例中,所述获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取所述目标散斑图像。
进一步地,在本发明的一个实施例中,还包括:
通过窗口函数对所述目标功率谱进行处理,其中,所述窗函数包括矩形窗口函数或塔基窗口函数。
进一步地,在本发明的一个实施例中,所述S2进一步包括:
通过Wirtinger Flow算法建立所述目标图像和所述目标功率谱的最小二乘代价函数,通过梯度优化算法对所述最小二乘代价函数进行迭代优化,求解所述最小二乘代价函数的最优解。
进一步地,在本发明的一个实施例中,所述S2进一步包括:
S211,建立所述目标图像和所述目标功率谱的所述最小二乘代价函数
Figure BDA0002232579460000031
其中,ai为所述目标功率谱,bi为傅里叶变换向量,
Figure BDA0002232579460000032
是bi的共轭向量,o为所述目标图像;
S212,设定常数λ,其中,
Figure BDA0002232579460000033
n为所述目标图像的维数,根据所述目标功率谱与傅里叶变换向量建立矩阵
Figure BDA0002232579460000034
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S213,根据初始向量o0值和梯度优化算法对所述最小二乘代价函数进行迭代优化使得f(o)取最小值,其中,
Figure BDA0002232579460000035
优化公式为:
其中,μτ+1是随迭代次数变化的迭代步长,τ0和μmax为经验值,
Figure BDA0002232579460000038
为所述最小二乘代价函数对所述目标图像的类梯度。
进一步地,在本发明的一个实施例中,所述S2进一步包括:
通过Truncated Wirtinger Flow算法建立所述目标图像和所述目标功率谱的log-似然代价函数,通过梯度优化算法对所述log-似然代价函数进行迭代优化,求解所述log-似然代价函数的最优解。
进一步地,在本发明的一个实施例中,所述S2进一步包括:
S221,建立所述目标图像和所述目标功率谱的所述log-似然代价函数
Figure BDA0002232579460000039
i=1、2、3...m,其中,ai为所述目标功率谱,bi为傅里叶变换向量,
Figure BDA00022325794600000310
是bi的共轭向量,o为所述目标图像;
S222,设定常数λ,其中,
Figure BDA00022325794600000311
n为所述目标图像的维数,根据所述目标功率谱与傅里叶变换向量建立矩阵
Figure BDA00022325794600000312
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S223,根据初始向量o0值和梯度优化算法对所述log-似然代价函数进行迭代优化使得l(o;ai)取最小值,其中,
Figure BDA0002232579460000041
优化公式为:
Figure BDA0002232579460000042
Figure BDA0002232579460000043
其中,
Figure BDA0002232579460000044
为所述log-似然代价函数对所述目标图像的类梯度,μτ是随迭代次数变化的迭代步长,
Figure BDA0002232579460000045
表示在集合
Figure BDA0002232579460000047
的交集中取1,不在时取0,
Figure BDA0002232579460000048
分别表示为:
Figure BDA00022325794600000410
Figure BDA00022325794600000411
Figure BDA00022325794600000412
αh,αy为设定的阈值。
为达到上述目的,本发明另一方面实施例提出了一种基于Wirtinger Flow算法的散斑相关成像装置,包括:
变换模块,用于获取目标散斑图像,根据维纳-辛钦定理对所述目标散斑图像进行自相关的傅里叶变换,得到目标功率谱;
成像模块,用于通过Wirtinger Flow算法建立目标图像和所述目标功率谱的代价函数,通过优化算法对所述代价函数进行优化并求解所述代价函数的最优解,根据所述最优解重建所述目标图像。
本发明实施例的基于Wirtinger Flow算法的散斑相关成像装置,通过搭建散斑相关成像系统,在光学记忆效应范围内照射目标,获取目标的散斑图像;对获取的散斑图像做相关计算;根据维纳-辛钦定理,目标的功率谱是散斑图像的自相关的傅里叶变换,通过对散斑图像进行自相关运算并进行傅里叶变换,可得到目标的功率谱;利用Wirtinger Flow相关算法建立目标功率谱与目标图像相关的代价函数,并求取函数对目标图像的类梯度并利用类似梯度下降的方法进行优化,最终得到最优解。由此,可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
另外,根据本发明上述实施例的基于Wirtinger Flow算法的散斑相关成像装置还可以具有以下附加的技术特征:
进一步地,在本发明的一个实施例中,所述获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取所述目标散斑图像。
进一步地,在本发明的一个实施例中,还包括:处理模块;
所述处理模块,用于通过窗口函数对所述目标功率谱进行处理,其中,所述窗函数包括矩形窗口函数或塔基窗口函数。
本发明附加的方面和优点将在下面的描述中部分给出,部分将从下面的描述中变得明显,或通过本发明的实践了解到。
附图说明
本发明上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解,其中:
图1为根据本发明一个实施例的基于Wirtinger Flow算法的散斑相关成像方法流程图;
图2为根据本发明一个实施例的散斑相关成像光学系统示意图;
图3为根据本发明一个实施例的散斑相关技术成像理论模型示意图;
图4为根据本发明一个实施例的基于Wirtinger Flow算法的散斑相关成像装置结构示意图。
具体实施方式
下面详细描述本发明的实施例,所述实施例的示例在附图中示出,其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的,旨在用于解释本发明,而不能理解为对本发明的限制。
下面参照附图描述根据本发明实施例提出的基于Wirtinger Flow算法的散斑相关成像方法及装置。
首先将参照附图描述根据本发明实施例提出的基于Wirtinger Flow算法的散斑相关成像方法。
图1为根据本发明一个实施例的基于Wirtinger Flow算法的散斑相关成像方法流程图。
如图1所示,该基于Wirtinger Flow算法的散斑相关成像方法包括以下步骤:
在步骤S1中,获取目标散斑图像,根据维纳-辛钦定理对目标散斑图像进行自相关的傅里叶变换,得到目标功率谱。
进一步地,在本发明的一个实施例中,获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取目标散斑图像。
具体地,如图2所示,高性能窄带单频激光器射出的激光经过旋转毛玻璃产生非相干赝热光,赝热光经过孔径光阑过滤掉一部分光,以合适的范围照射目标,使得被照目标范围在光学记忆效应范围内。目标后方u处放置宽度为L的强散射介质,光透过强散射介质后经过孔径光阑,在散射介质后方v处被探测器接收。携带有目标信息的光束经过距离u入射到强散射介质中,并在其内部发生多重散射,出射的散射光已无入射光的场分布,而是产生了新的光场,携带的目标信息被重新“编码”,出射光在距离v处形成的像即是散斑图像。在光学记忆效应范围内,目标被照射的部分的各个点产生的散斑场是基本不变的,只是发生了位移变化,该范围可表示为:Δx<<u·λ/πL,其中,Δx是目标被照射的范围,u是目标到散射介质的距离,L是散射介质的厚度。根据光学记忆效应的特点,可以将成像系统看成位移不变性的光学成像系统。
如图3所示,探测器接收到的散斑图像是目标图像与系统的点扩散函数PSF的卷积。假设探测器上r处的测量值为I(r),那么探测器上测得的测量值可表示为:
Figure BDA0002232579460000061
其中,O(r)是目标图像,S(r)是系统的PSF。上述表达式可简写为:
I=O*S (2)
对I计算自相关:
IΘI=(O*S)Θ(O*S)=(OΘO)*(SΘS)=OΘO+C (3)
其中,Θ为自相关运算,*为卷积运算,(SΘS)为峰值函数,用常数C表示,其本质上是由噪声引入的背景项。由上述公式可以看出,散斑图像的自相关运算等于目标的自相关运算加上一个背景项。除去背景项,再利用相位恢复算法可以有效地恢复出目标图像。
进一步地,在得到目标散斑图像后,对获取的目标散斑图像进行数学建模,并根据维纳-辛钦定理,对目标散斑图像进行相关运算后再对其傅里叶变换,得到目标的功率谱。
具体地,假设R表示散斑场的自相关,根据维纳-辛钦定理:目标的功率谱可通过其自相关的傅里叶变换得到。由于光学记忆效应的有效范围限制,散斑图像只有中心部分是携带有效目标信息的,因此可以设置一个窗口函数W(x,y)求取有效功率谱:
A(kx,ky)=|FT{W(x,y)R(x,y)}| (4)
其中,W(x,y)可以选择矩形窗口或者塔基窗口,大小一般设置为100*100-320*320个像素之间。
由公式(3)和公式(4)可知:
A(kx,ky)=|FT{O(x,y)}|2 (5)
将公式(5)转成下面这种形式:
ai=|<bi,o>|2 (6)
其中,
Figure BDA0002232579460000071
为功率谱,
Figure BDA0002232579460000072
为傅里叶变换向量,
Figure BDA0002232579460000073
为目标向量。从ai恢复目标o属于典型的非凸二次规划问题。
在步骤S2中,通过Wirtinger Flow算法建立目标图像和目标功率谱的代价函数,通过优化算法对代价函数进行优化并求解代价函数的最优解,根据最优解重建目标图像。
进一步地,步骤S2包括:
通过Wirtinger Flow算法建立目标图像和目标功率谱的最小二乘代价函数,通过梯度优化算法对最小二乘代价函数进行迭代优化,求解最小二乘代价函数的最优解。
进一步地,步骤S2进一步包括:
S211,建立目标图像和目标功率谱的最小二乘代价函数
Figure BDA0002232579460000074
i=1、2、3...m,其中,ai为目标功率谱,bi为傅里叶变换向量,
Figure BDA0002232579460000075
是bi的共轭向量,o为目标图像;
S212,设定常数λ,其中,
Figure BDA0002232579460000076
n为目标图像的维数,根据目标功率谱与傅里叶变换向量建立矩阵
Figure BDA0002232579460000077
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S213,根据初始向量o0值和梯度优化算法对最小二乘代价函数进行迭代优化使得f(o)取最小值,其中,优化公式为:
Figure BDA0002232579460000079
其中,μτ+1是随迭代次数变化的迭代步长,
Figure BDA00022325794600000710
τ0和μmax为经验值,为最小二乘代价函数对目标图像的类梯度。
为了解决非凸二次规划问题,建立一个代价函数,那么公式(6)的求解就相当于下面公式的求解:
对于公式(7),求得代价函数对目标图像的类梯度:
从而利用类似梯度下降的方法求得目标o。
可以理解的是,首先求取一个接近真实解的初始向量,其初始向量为目标功率谱与傅里叶变换向量共同构建的矩阵
Figure BDA0002232579460000083
的最大特征值对应的特征向量,令初始向量的欧式距离等于一个跟目标功率谱和傅里叶变换向量有关的常数
Figure BDA0002232579460000084
得到初始向量值。然后利用获得的初始向量和类梯度并设置好合适步长进行迭代下降。
具体地,Wirtinger Flow算法主要包括最下面步骤:
(1)通过特征值获得一个接近正确解的初始向量o0。首先,假设一个常数λ,其值与测量值ai和测量向量bi有关:
Figure BDA0002232579460000085
然后,建立一个矩阵:
Figure BDA0002232579460000086
令o0为矩阵Y的最大特征值所对应的特征向量,并令:
||o0||=λ
得到初始估计值o0
(2)通过一个类似梯下降的方法,通过由第(1)步骤得到的初始值o0按照以下方式进行迭代:
Figure BDA0002232579460000087
其中,μτ+1是随迭代次数变化的迭代步长,它的值定义为:
Figure BDA0002232579460000088
其中,τ0和μmax为经验值,根据实际情况设置数值。通过不断迭代,直到使f(o)最小,此时得到的o即为最优解。算法流程如表1所示,表1为散斑相关成像的Wirtinger Flow算法流程。
表1
Figure BDA0002232579460000091
进一步地,为了进一步提高重建过程的鲁棒性并减少计算复杂度,应用改进Wirtinger Flow算法——Truncated Wirtinger Flow算法。
在本发明的一个实施例中,步骤S2包括:
通过Truncated Wirtinger Flow算法建立目标图像和目标功率谱的log-似然代价函数,通过梯度优化算法对log-似然代价函数进行迭代优化,求解log-似然代价函数的最优解。
进一步地,S2具体包括:
S221,建立目标图像和目标功率谱的log-似然代价函数
Figure BDA0002232579460000092
i=1、2、3...m,其中,ai为目标功率谱,bi为傅里叶变换向量,
Figure BDA0002232579460000093
是bi的共轭向量,o为目标图像;
S222,设定常数λ,其中,n为目标图像的维数,根据目标功率谱与傅里叶变换向量建立矩阵
Figure BDA0002232579460000101
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S223,根据初始向量o0值和梯度优化算法对log-似然代价函数进行迭代优化使得l(o;ai)取最小值,其中,
Figure BDA0002232579460000102
优化公式为:
Figure BDA0002232579460000104
其中,
Figure BDA0002232579460000105
为log-似然代价函数对目标图像的类梯度,μτ是随迭代次数变化的迭代步长,
Figure BDA0002232579460000106
表示在集合
Figure BDA0002232579460000107
Figure BDA0002232579460000108
的交集中取1,不在时取0,分别表示为:
Figure BDA00022325794600001011
Figure BDA00022325794600001012
αh,αy为设定的阈值。
具体地,对于公式(6)的求解可以用最大似然估计的方法将其转换为log-似然函数:
Figure BDA00022325794600001014
然后,求l(o;ai)对o的类梯度并正则化:
Figure BDA00022325794600001015
最后,用一个类似梯度下降的方式对o进行迭代更新:
其中,μτ是随迭代次数改变地步长;
Figure BDA00022325794600001017
表示在集合
Figure BDA00022325794600001018
Figure BDA00022325794600001019
的交集中取1,不在时取0;
Figure BDA00022325794600001020
Figure BDA00022325794600001021
分别表示为:
Figure BDA0002232579460000111
Figure BDA0002232579460000112
对于公式(11)的迭代模型,其初始迭代向量的求取与Wirtinger Flow算法大致相似,不同的地方在于对矩阵Y的范围加了一个阈值:
Figure BDA0002232579460000113
其中,
Figure BDA0002232579460000114
αh,αy为设定的阈值。算法流程如图表2所示,表2为散斑相关成像的Truncated Wirtinger Flow算法流程。
表2
Figure BDA0002232579460000115
根据本发明实施例提出的基于Wirtinger Flow算法的散斑相关成像方法,通过搭建散斑相关成像系统,在光学记忆效应范围内照射目标,获取目标的散斑图像;对获取的散斑图像做相关计算;根据维纳-辛钦定理,目标的功率谱是散斑图像的自相关的傅里叶变换,通过对散斑图像进行自相关运算并进行傅里叶变换,可得到目标的功率谱;利用Wirtinger Flow相关算法建立目标功率谱与目标图像相关的代价函数,并求取函数对目标图像的类梯度并利用类似梯度下降的方法进行优化,最终得到最优解。由此,可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
其次参照附图描述根据本发明实施例提出的基于Wirtinger Flow算法的散斑相关成像装置。
图4为根据本发明一个实施例的基于Wirtinger Flow算法的散斑相关成像装置结构示意图。
如图4所示,该基于Wirtinger Flow算法的散斑相关成像装置包括:变换模块100和成像模块200。
变换模块100,用于获取目标散斑图像,根据维纳-辛钦定理对目标散斑图像进行自相关的傅里叶变换,得到目标功率谱。
成像模块200,用于通过Wirtinger Flow算法建立目标图像和目标功率谱的代价函数,通过优化算法对代价函数进行优化并求解代价函数的最优解,根据最优解重建目标图像。
该装置可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
进一步地,在本发明的一个实施例中,获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取目标散斑图像。
进一步地,在本发明的一个实施例中,还包括:处理模块;
处理模块,用于通过窗口函数对目标功率谱进行处理,其中,窗函数包括矩形窗口函数或塔基窗口函数。
需要说明的是,前述对基于Wirtinger Flow算法的散斑相关成像方法实施例的解释说明也适用于该实施例的装置,此处不再赘述。
根据本发明实施例提出的基于Wirtinger Flow算法的散斑相关成像装置,通过搭建散斑相关成像系统,在光学记忆效应范围内照射目标,获取目标的散斑图像;对获取的散斑图像做相关计算;根据维纳-辛钦定理,目标的功率谱是散斑图像的自相关的傅里叶变换,通过对散斑图像进行自相关运算并进行傅里叶变换,可得到目标的功率谱;利用Wirtinger Flow相关算法建立目标功率谱与目标图像相关的代价函数,并求取函数对目标图像的类梯度并利用类似梯度下降的方法进行优化,最终得到最优解。由此,可以提高目标重建过程对于测量噪声和系统畸变的鲁棒性,不需要目标的任何先验信息,提高重建目标的质量。
此外,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括至少一个该特征。在本发明的描述中,“多个”的含义是至少两个,例如两个,三个等,除非另有明确具体的限定。
在本说明书的描述中,参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外,在不相互矛盾的情况下,本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。
尽管上面已经示出和描述了本发明的实施例,可以理解的是,上述实施例是示例性的,不能理解为对本发明的限制,本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (10)

1.一种基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,包括以下步骤:
S1,获取目标散斑图像,根据维纳-辛钦定理对所述目标散斑图像进行自相关的傅里叶变换,得到目标功率谱;
S2,通过Wirtinger Flow算法建立目标图像和所述目标功率谱的代价函数,通过优化算法对所述代价函数进行优化并求解所述代价函数的最优解,根据所述最优解重建所述目标图像。
2.根据权利要求1所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,所述获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取所述目标散斑图像。
3.根据权利要求1所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,还包括:
通过窗口函数对所述目标功率谱进行处理,其中,所述窗函数包括矩形窗口函数或塔基窗口函数。
4.根据权利要求1所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,所述S2进一步包括:
通过Wirtinger Flow算法建立所述目标图像和所述目标功率谱的最小二乘代价函数,通过梯度优化算法对所述最小二乘代价函数进行迭代优化,求解所述最小二乘代价函数的最优解。
5.根据权利要求4所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,所述S2进一步包括:
S211,建立所述目标图像和所述目标功率谱的所述最小二乘代价函数
Figure FDA0002232579450000011
其中,ai为所述目标功率谱,bi为傅里叶变换向量,是bi的共轭向量,o为所述目标图像;
S212,设定常数λ,其中,
Figure FDA0002232579450000013
n为所述目标图像的维数,根据所述目标功率谱与傅里叶变换向量建立矩阵
Figure FDA0002232579450000014
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S213,根据初始向量o0值和梯度优化算法对所述最小二乘代价函数进行迭代优化使得f(o)取最小值,其中,
Figure FDA0002232579450000021
优化公式为:
Figure FDA0002232579450000022
其中,μτ+1是随迭代次数变化的迭代步长,τ0和μmax为经验值,为所述最小二乘代价函数对所述目标图像的类梯度。
6.根据权利要求1所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,所述S2进一步包括:
通过Truncated Wirtinger Flow算法建立所述目标图像和所述目标功率谱的log-似然代价函数,通过梯度优化算法对所述log-似然代价函数进行迭代优化,求解所述log-似然代价函数的最优解。
7.根据权利要求6所述的基于Wirtinger Flow算法的散斑相关成像方法,其特征在于,所述S2进一步包括:
S221,建立所述目标图像和所述目标功率谱的所述log-似然代价函数
Figure FDA0002232579450000025
其中,ai为所述目标功率谱,bi为傅里叶变换向量,是bi的共轭向量,o为所述目标图像;
S222,设定常数λ,其中,n为所述目标图像的维数,根据所述目标功率谱与傅里叶变换向量建立矩阵
Figure FDA0002232579450000028
设定初始向量o0为矩阵Y的最大特征值对应的特征向量,并使||o0||=λ,求解初始向量o0值;
S223,根据初始向量o0值和梯度优化算法对所述log-似然代价函数进行迭代优化使得l(o;ai)取最小值,其中,
Figure FDA0002232579450000029
优化公式为:
Figure FDA0002232579450000031
其中,
Figure FDA0002232579450000032
为所述log-似然代价函数对所述目标图像的类梯度,μτ是随迭代次数变化的迭代步长,
Figure FDA0002232579450000033
表示在集合
Figure FDA0002232579450000034
Figure FDA0002232579450000035
的交集中取1,不在时取0,
Figure FDA0002232579450000036
Figure FDA0002232579450000037
分别表示为:
Figure FDA0002232579450000038
Figure FDA0002232579450000039
Figure FDA00022325794500000310
αh,αy为设定的阈值。
8.一种基于Wirtinger Flow算法的散斑相关成像装置,其特征在于,包括:
变换模块,用于获取目标散斑图像,根据维纳-辛钦定理对所述目标散斑图像进行自相关的傅里叶变换,得到目标功率谱;
成像模块,用于通过Wirtinger Flow算法建立目标图像和所述目标功率谱的代价函数,通过优化算法对所述代价函数进行优化并求解所述代价函数的最优解,根据所述最优解重建所述目标图像。
9.根据权利要求6所述的基于Wirtinger Flow算法的散斑相关成像装置,其特征在于,所述获取目标散斑图像,包括:
在预先搭建的无透镜的散斑相关成像光学系统中,通过激光器发出光波,经过旋转毛玻璃产生非相干赝热光,非相干赝热光经过孔径光阑后在光学记忆效应范围内照射目标,再透过强散射介质被探测器接收,通过探测器获取所述目标散斑图像。
10.根据权利要求6所述的基于Wirtinger Flow算法的散斑相关成像装置,其特征在于,还包括:处理模块;
所述处理模块,用于通过窗口函数对所述目标功率谱进行处理,其中,所述窗函数包括矩形窗口函数或塔基窗口函数。
CN201910972577.1A 2019-10-14 2019-10-14 基于Wirtinger Flow算法的散斑相关成像方法及装置 Active CN110807822B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910972577.1A CN110807822B (zh) 2019-10-14 2019-10-14 基于Wirtinger Flow算法的散斑相关成像方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910972577.1A CN110807822B (zh) 2019-10-14 2019-10-14 基于Wirtinger Flow算法的散斑相关成像方法及装置

Publications (2)

Publication Number Publication Date
CN110807822A true CN110807822A (zh) 2020-02-18
CN110807822B CN110807822B (zh) 2022-03-22

Family

ID=69488364

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910972577.1A Active CN110807822B (zh) 2019-10-14 2019-10-14 基于Wirtinger Flow算法的散斑相关成像方法及装置

Country Status (1)

Country Link
CN (1) CN110807822B (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111724328A (zh) * 2020-06-30 2020-09-29 苏州兴钊防务研究院有限公司 一种光电协同散射介质成像系统及其方法
CN112161953A (zh) * 2020-08-25 2021-01-01 西安电子科技大学 一种基于散射介质的宽光谱单帧散射成像方法
CN112634380A (zh) * 2020-12-01 2021-04-09 西安电子科技大学 一种单帧超光学记忆效应的多目标彩色散射成像方法
CN112907481A (zh) * 2021-03-13 2021-06-04 南京大学 一种对于噪声鲁棒的高质量无透镜成像方法及系统
CN113781592A (zh) * 2021-07-01 2021-12-10 杭州电子科技大学 一种基于复杂度引导相位恢复的散斑成像重建方法
CN113804101A (zh) * 2020-06-11 2021-12-17 中国科学院上海光学精密机械研究所 透过散射介质的无侵入光学成像和定位的装置和方法
CN114022365A (zh) * 2021-11-25 2022-02-08 中国科学院光电技术研究所 一种梯度下降散斑照明超分辨目标重建方法
CN115508899A (zh) * 2022-10-21 2022-12-23 中铁二院工程集团有限责任公司 基于最优化理论的航空大地电磁功率谱估算方法及设备
CN115696041A (zh) * 2022-10-26 2023-02-03 清华大学 基于波前调制迭代的非侵入式散射介质内部聚焦成像方法

Citations (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120014611A1 (en) * 2010-07-19 2012-01-19 Dean Bruce H System and method for determining phase retrieval sampling from the modulation transfer function
CN103744061A (zh) * 2014-01-15 2014-04-23 西安电子科技大学 基于迭代最小二乘方法的mimo雷达doa估计方法
WO2015063779A1 (en) * 2013-11-04 2015-05-07 Yeda Research And Development Co. Ltd. System and method for phase retrieval in lensless imaging
US20170006218A1 (en) * 2014-03-31 2017-01-05 Fujifilm Corporation Image processing device, imaging device, image processing method, and image processing program
CN108227187A (zh) * 2018-01-24 2018-06-29 深圳大学 一种扩展光学成像景深的方法及系统
CN108550108A (zh) * 2017-09-28 2018-09-18 武汉大学 一种基于相位迭代最小化的傅里叶叠层成像图像重建方法
CN108983579A (zh) * 2018-09-05 2018-12-11 南京大学 无透镜数字全息显微成像相位恢复和重建的方法及其装置
CN109828371A (zh) * 2019-03-28 2019-05-31 清华大学深圳研究生院 一种基于移动散斑光源的大视场散射成像方法
CN110132901A (zh) * 2019-05-21 2019-08-16 北京理工大学 合成孔径穿散射介质成像的系统和方法
CN110132175A (zh) * 2019-05-30 2019-08-16 北京理工大学 基于幅度调制的单像素相位成像方法和装置
CN110274877A (zh) * 2019-05-21 2019-09-24 西安电子科技大学 一种基于散射介质的3d光谱成像系统及方法

Patent Citations (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120014611A1 (en) * 2010-07-19 2012-01-19 Dean Bruce H System and method for determining phase retrieval sampling from the modulation transfer function
WO2015063779A1 (en) * 2013-11-04 2015-05-07 Yeda Research And Development Co. Ltd. System and method for phase retrieval in lensless imaging
CN103744061A (zh) * 2014-01-15 2014-04-23 西安电子科技大学 基于迭代最小二乘方法的mimo雷达doa估计方法
US20170006218A1 (en) * 2014-03-31 2017-01-05 Fujifilm Corporation Image processing device, imaging device, image processing method, and image processing program
CN108550108A (zh) * 2017-09-28 2018-09-18 武汉大学 一种基于相位迭代最小化的傅里叶叠层成像图像重建方法
CN108227187A (zh) * 2018-01-24 2018-06-29 深圳大学 一种扩展光学成像景深的方法及系统
CN108983579A (zh) * 2018-09-05 2018-12-11 南京大学 无透镜数字全息显微成像相位恢复和重建的方法及其装置
CN109828371A (zh) * 2019-03-28 2019-05-31 清华大学深圳研究生院 一种基于移动散斑光源的大视场散射成像方法
CN110132901A (zh) * 2019-05-21 2019-08-16 北京理工大学 合成孔径穿散射介质成像的系统和方法
CN110274877A (zh) * 2019-05-21 2019-09-24 西安电子科技大学 一种基于散射介质的3d光谱成像系统及方法
CN110132175A (zh) * 2019-05-30 2019-08-16 北京理工大学 基于幅度调制的单像素相位成像方法和装置

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
EMMANUEL J. CANDÈS等: "Phase retrieval via Wirtinger flow Theory and algorithms", 《 IEEE TRANSACTIONS ON INFORMATION THEORY》 *
MENG LI等: "Noise-robust coded-illumination imaging with low computational complexity", 《OPTICS EXPRESS》 *
YUXIN CHEN等: "Solving Random Quadratic Systems of Equations Is Nearly as Easy as Solving Linear Systems", 《COMMUNICATIONS ON PURE AND APPLIED MATHEMATICS》 *
刘欢: "基于Wirtinger_Flow算法的相位恢复", 《中国优秀硕士学位论文全文数据库基础学科辑》 *
肖晓等: "基于赝热光照明的单发光学散斑成像", 《物理学报》 *
董昕等: "基于压缩感知和散斑相关法的散射介质成像方法研究", 《激光生物学报》 *
贾辉等: "透过散射介质对直线运动目标的全光成像及追踪技术", 《物理学报》 *

Cited By (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113804101A (zh) * 2020-06-11 2021-12-17 中国科学院上海光学精密机械研究所 透过散射介质的无侵入光学成像和定位的装置和方法
CN111724328B (zh) * 2020-06-30 2024-03-05 苏州兴钊防务研究院有限公司 一种光电协同散射介质成像系统及其方法
CN111724328A (zh) * 2020-06-30 2020-09-29 苏州兴钊防务研究院有限公司 一种光电协同散射介质成像系统及其方法
CN112161953A (zh) * 2020-08-25 2021-01-01 西安电子科技大学 一种基于散射介质的宽光谱单帧散射成像方法
CN112634380B (zh) * 2020-12-01 2023-08-15 西安电子科技大学 一种单帧超光学记忆效应的多目标彩色散射成像方法
CN112634380A (zh) * 2020-12-01 2021-04-09 西安电子科技大学 一种单帧超光学记忆效应的多目标彩色散射成像方法
CN112907481A (zh) * 2021-03-13 2021-06-04 南京大学 一种对于噪声鲁棒的高质量无透镜成像方法及系统
CN112907481B (zh) * 2021-03-13 2023-12-12 南京大学 一种对于噪声鲁棒的高质量无透镜成像方法及系统
CN113781592A (zh) * 2021-07-01 2021-12-10 杭州电子科技大学 一种基于复杂度引导相位恢复的散斑成像重建方法
CN113781592B (zh) * 2021-07-01 2024-04-23 杭州电子科技大学 一种基于复杂度引导相位恢复的散斑成像重建方法
CN114022365A (zh) * 2021-11-25 2022-02-08 中国科学院光电技术研究所 一种梯度下降散斑照明超分辨目标重建方法
CN114022365B (zh) * 2021-11-25 2023-05-26 中国科学院光电技术研究所 一种梯度下降散斑照明超分辨目标重建方法
CN115508899B (zh) * 2022-10-21 2023-08-29 中铁二院工程集团有限责任公司 基于最优化理论的航空大地电磁功率谱估算方法及设备
CN115508899A (zh) * 2022-10-21 2022-12-23 中铁二院工程集团有限责任公司 基于最优化理论的航空大地电磁功率谱估算方法及设备
CN115696041A (zh) * 2022-10-26 2023-02-03 清华大学 基于波前调制迭代的非侵入式散射介质内部聚焦成像方法
CN115696041B (zh) * 2022-10-26 2023-11-14 清华大学 基于波前调制迭代的非侵入式散射介质内部聚焦成像方法

Also Published As

Publication number Publication date
CN110807822B (zh) 2022-03-22

Similar Documents

Publication Publication Date Title
CN110807822B (zh) 基于Wirtinger Flow算法的散斑相关成像方法及装置
EP3506209A1 (en) Image processing method, image processing device and storage medium
Pantin et al. Deconvolution and blind deconvolution in astronomy
Sheng et al. A constrained variable projection reconstruction method for photoacoustic computed tomography without accurate knowledge of transducer responses
Roitner et al. Deblurring algorithms accounting for the finite detector size in photoacoustic tomography
CN112634380A (zh) 一种单帧超光学记忆效应的多目标彩色散射成像方法
González et al. Compressive optical deflectometric tomography: A constrained total-variation minimization approach
Yu et al. A multiplicative Nakagami speckle reduction algorithm for ultrasound images
Nguyen et al. An adaptive filter to approximate the Bayesian strategy for sonographic beamforming
Fedorczak et al. Tomographic reconstruction of tokamak plasma light emission from single image using wavelet-vaguelette decomposition
Churchill et al. Sampling-based spotlight SAR image reconstruction from phase history data for speckle reduction and uncertainty quantification
CN111369638B (zh) 激光反射层析成像欠采样的重建方法、存储介质和系统
Maji et al. A multifractal-based wavefront phase estimation technique for ground-based astronomical observations
Churchill Synthetic aperture radar image formation with uncertainty quantification
Dong et al. A study of multi-static ultrasonic tomography using propagation and back-propagation method
Balodi et al. Comparison of despeckle filters for ultrasound images
Paltauf et al. Artifact removal in photoacoustic section imaging by combining an integrating cylindrical detector with model-based reconstruction
Sahlström et al. Utilizing variational autoencoders in photoacoustic tomography
John et al. Compressive Sensing Based Algorithms for Limited-View PAT Image Reconstruction
Oral et al. Plug-and-play reconstruction with 3D deep prior for complex-valued near-field MIMO imaging
Pal Multi-scale methods for deconvolution from wavefront sensing
US20220262050A1 (en) Systems and methods for opto-acoustic image reconstruction with multiple acquisitions
Gerg et al. Deep adaptive phase learning: Enhancing synthetic aperture sonar imagery through coherent autofocus,”
Chira et al. Blind deconvolution for ultrasound sequences using a noninverse greedy algorithm
Majee High speed imaging via advanced modeling

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant