CN111795949A - 抗散射成像方法与装置 - Google Patents

Abstract

本申请提出一种抗散射成像方法与装置，包括以下步骤：获取散斑图；计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位；其中，N为正整数；根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位；根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。由此，从单张散斑图即可重建原始目标信息，且幅值和相位恢复过程相互独立，相位恢复精度高，抗噪性强。

Description

抗散射成像方法与装置

技术领域

本申请涉及计算摄像学技术领域，尤其涉及一种抗散射成像方法与装置。

背景技术

当自由空间存在阻碍的散射介质时，如大雾、生物组织等，由于介质的不均匀性使得光子与其互相作用后改变传播方向。这种相互作用往往是多次且不确定的，最终会导致光路变得随机化，即散射。从视觉角度来讲，光通过强散射介质后会失去所携带的目标信息。设备捕获散射光后所成的图像由无数杂乱无章的斑点组成，称为散斑图。

目前已经发展出了许多透过散射介质成像的方法。但大多数方法都需要侵入原始目标平面内或者放置“引导星”，且要用到多张的散斑图。

2014年O.Katz等人提出了一种基于单张散斑自相关的非侵入式实时散射成像方案，发现物体散斑图的自相关本质上和物体强度图的自相关是相等的。通过计算散斑图的自相关函数，获得目标傅里叶域幅值，之后使用相位恢复算法从傅里叶域幅值重建出原始目标信息。该方案极大简化了散斑数据的采集以及成像流程，在一定程度上能实现实时成像。一经提出引起了许多研究者的关注，后续对其改进与提高成像质量的工作成果层出不穷。

O.Katz等人基于散斑自相关获取原始目标傅里叶域幅值信息后，使用相位恢复算法恢复傅里叶域相位信息。这个过程中存在一些缺陷。首先，相位恢复算法依赖于幅值信息，若幅值不准确则相位信息也不会准确；其次，大部分相位恢复算法恢复结果并不理想，只能重建结构简单的原始目标，且无法知道目标在背景中的具体位置；最后，相位恢复算法的抗噪性能不理想，由于拍摄硬件限制导致散斑图噪声较大时，会极大影响重建质量。

发明内容

本申请旨在至少在一定程度上解决相关技术中的技术问题之一。

为此，本申请的一个目的在于提出一种抗散射成像方法和装置，基于双谱分析和高斯牛顿优化对相位进行恢复，能够在一定程度上解决上述问题，双谱分析将目标幅值与相位的估计过程互相独立，且具有良好的抗噪性能，高斯牛顿法进一步优化使相位精度更高。

本申请的另一个目的在于提出一种抗散射成像装置。

为达到上述目的，本申请实施例提出了一种抗散射成像方法，包括以下步骤：

获取散斑图；

计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；

对所述散斑图进行分块处理得到N张子散斑图，并利用所述N张子散斑图计算双谱得到双谱相位；其中，N为正整数；

根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位；

根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像。

本申请实施例的抗散射成像方法，获取散斑图；计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位；其中，N为正整数；根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位；根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。由此，从单张散斑图即可重建原始目标信息，且幅值和相位恢复过程相互独立，相位恢复精度高，抗噪性强。

另外，根据本申请上述实施例的抗散射成像方法，还可以具有如下附加的技术特征：

进一步地，在本申请的一个实施例中，所述获取散斑图，包括：获取具有原始目标信息光线经过散射介质形成的多重散射光束；使用探测器采集所述多重散射光束生成所述散斑图。

进一步地，在所述计算所述散斑图的自相关信息之前，还包括：将所述散斑图除以其自身低频版本；或，对所述散斑图进行多项式拟合；采用高斯核卷积算法对所述散斑图进行平滑滤波处理。

进一步地，在本申请的一个实施例中，所述计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息，包括：

通过以下公式(1)计算所述散斑图的自相关信息：

[I★I]＝[(O*S)★(O*S)]＝[(O★O)*(S★S)]；

其中，I为探测器平面捕获的散斑图，O为原始目标图像，S是散射系统的点扩散函数；I＝O*S，“*”是卷积算子，“★”是自相关算子；

根据维纳辛钦定理，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息为|F{(O★O)}|＝|F{O}|²。

进一步地，在本申请的一个实施例中，对所述公式(1)进行简化处理得到公式(2)：

[I★I]＝[(O★O)]+C；

其中，C是在计算散斑图的自相关过程中产生的常数项。

进一步地，在本申请的一个实施例中，在所述利用所述N张子散斑图计算双谱得到双谱相位之前，还包括：采用高斯窗函数对每块子散斑进行空间滤波；其中，所述高斯窗函数为G＝exp[-(x²+y²/ω²)，x和y表示空间坐标，ω用来控制窗函数窗口大小。

进一步地，在本申请的一个实施例中，利用所述N张子散斑图计算双谱得到双谱相位，对于子散斑图o(x)，相对移位(x₁，x₂)的三阶相关为公式(3)：

对所述公式(3)式进行傅里叶变换并结合卷积的性质可以得到所述双谱相位为公式(4)：

O⁽³⁾(u，v)＝O(u)O(v)O^*(u+v)；

其中，O⁽³⁾(u，v)是所述双谱，O(u)是目标图像的傅里叶变换，u，

是空间频率，所述公式(4)表示所述双谱中含有傅里叶信息。

进一步地，在本申请的一个实施例中，所述根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位，包括：

对所述公式(4)取相位形式为公式(5)：

β(u，v)＝φ(u)+φ(v)-φ(u+v)；

其中，φ(u)是O(u)的相位，β(u，v)是所述双谱相位，通过所述公式(5)与递归算法从散斑图的双谱相位中恢复出原始目标的傅里叶相位信息φ(u)，使用φ(u)作为高斯牛顿优化的初始值φ₀(u)；

所述公式(5)表示为线性结构：β＝Aφ，其中，

表示需要恢复的目标傅里叶相位，

表示双谱相位，

是每行只有三个非零元素的稀疏矩阵，它的三个非零元素包含两个1和一个-1，对应于(5)式等号右边三个相位元素的符号；

构造最小二乘目标函数来优化相位，从而降低相位与散斑数据双谱之间的误差，目标函数为公式(6)：

其中，

是一个对角加权矩阵，对所述公式(6)进行转换得到公式(7)：

使用高斯牛顿迭代法对(6)或(7)式求解得到所述高精度傅里叶域相位。

进一步地，在本申请的一个实施例中，所述根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像，包括：

对所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位进行逆傅里叶变换或在高斯牛顿优化傅里叶域相位时将所述高精度傅里叶域相位更改为所述原始目标图像。

为达到上述目的，本申请第二方面实施例提出了一种抗散射成像装置，包括：

获取模块，用于获取散斑图；

计算变换模块，用于计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；

分块计算模块，用于对所述散斑图进行分块处理得到N张子散斑图，并利用所述N张子散斑图计算双谱得到双谱相位；其中，N为正整数；

计算优化模块，用于根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位；

重建模块，用于根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像。

本申请实施例的抗散射成像装置，获取散斑图；计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位；其中，N为正整数；根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位；根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。由此，从单张散斑图即可重建原始目标信息，且幅值和相位恢复过程相互独立，相位恢复精度高，抗噪性强。

本申请附加的方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本申请的实践了解到。

附图说明

本申请上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解，其中：

图1为本申请实施例的一种抗散射成像方法的流程示例图；

图2散射成像整体流程图；

图3基于散斑自相关获取傅里叶幅值流程图；

图4对单张散斑图分块示意图；

图5基于双谱相位恢复和高斯牛顿优化获取傅里叶域相位流程图；

图6重建原始目标图像流程图；

图7为本申请实施例提供的一种抗散射成像装置的结构示意图。

具体实施方式

下面详细描述本申请的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，旨在用于解释本申请，而不能理解为对本申请的限制。

下面参照附图描述根据本申请实施例提出的抗散射成像方法与装置。

图1为本申请实施例的一种抗散射成像方法的流程示例图。如图1所示，该抗散射成像方法，包括：

步骤101，获取散斑图。

具体地，使用探测器采集目标场景经过散射介质后的单张散斑图。

在本申请的一个实施例中，获取散斑图，包括：获取具有原始目标信息光线经过散射介质形成的多重散射光束；使用探测器采集多重散射光束生成所述散斑图。

步骤102，计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息。

在本申请的一个实施例中，在计算散斑图的自相关信息之前，还包括：将散斑图除以其自身低频版本；或，对散斑图进行多项式拟合；采用高斯核卷积算法对散斑图进行平滑滤波处理。

需要说明的是，对散斑图进行预处理，包括但不限于消除渐变光晕和平滑去噪，预处理会使原始目标重建质量更高。

具体地，获取原始目标傅里叶域幅值，对预处理后的散斑图计算自相关函数，由于原始目标的自相关近似等于散斑图的自相关，能够从散斑图自相关中获取原始目标的信息，若探测器采集的散斑图分辨率高，其自相关函数仅在中心小范围内存在有效信息，可以使用窗函数截取中心的小范围区域，之后，根据维纳辛钦定理，对自相关函数进行傅里叶变换后取模值，可获得原始目标的功率谱密度即傅里叶域幅值信息。

在本申请的一个实施例中，计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息，包括：

通过以下公式(1)计算散斑图的自相关信息：

[I★I]＝[(O*S)★(O*S)]＝[(O★O)*(S★S)]；

其中，I为探测器平面捕获的散斑图，O为原始目标图像，S是散射系统的点扩散函数；I＝O*S，“*”是卷积算子，“★”是自相关算子；

根据维纳辛钦定理，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息为|F{(O★O)}|＝|F{O}|²。

在本申请的一个实施例中，对公式(1)进行简化处理得到公式(2)：

[I★I]＝[(O★O)]+C；

其中，C是在计算散斑图的自相关过程中产生的常数项。

步骤103，对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位，其中，N为正整数。

步骤104，根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位。

步骤105，根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。

具体地，双谱分析恢复相位需要用到多张散斑图，本方法使用空间代替时间的思想，对单张散斑图进行分块，获得若干张子散斑图。每张子散斑图都包含原始目标的全部信息。

可以理解的是，计算子散斑图的双谱，获取双谱相位，双谱和三阶相关性是一对傅里叶变换对，计算双谱的方法包括但不限于计算三阶相关性后进行傅里叶变换或基于双谱定义式：O⁽³⁾(u，v)＝O(u)O(v)O^*(u+v)，从目标的不同空间频率u，v，u+v中直接计算双谱，之后计算若干张子散斑图的平均双谱相位。

另外，由于散斑图的双谱近似等于原始目标的双谱，能够从散斑图双谱中获取原始目标的双谱信息，利用双谱相位β与目标傅里叶域相位φ的关系：β(u，v)＝φ(u)+φ(v)-φ(u+v)，通过递归的方式获取目标傅里叶域相位的初始估计值。

具体地，构造从双谱相位β恢复目标相位φ的线性结构β＝Aφ。其中，A是每行有三个非零元素的稀疏矩阵，三个非零元素包含两个1和一个-1，对应于递归式三个相位元素φ(u)、φ(v)、φ(u+v)的正负符号，之后构造优化目标函数，形式包括但不限于最小二乘。使用高斯牛顿法求解得出目标傅里叶域相位的优化值。

在本申请的一个实施例中，在利用N张子散斑图计算双谱得到双谱相位之前，还包括：采用高斯窗函数对每块子散斑进行空间滤波；其中，高斯窗函数为G＝exp[-(x²+y²/ω²)，x和y表示空间坐标，ω用来控制窗函数窗口大小。

在本申请的一个实施例中，利用N张子散斑图计算双谱得到双谱相位，

对于子散斑图o(x)，相对移位(x₁，x₂)的三阶相关为公式(3)：

对公式(3)式进行傅里叶变换并结合卷积的性质可以得到双谱相位为公式(4)：

O(3⁾(u，v)＝O(u)O(v)O^*(u+v)；

其中，O⁽³⁾(u，v)是双谱，O(u)是目标图像的傅里叶变换，

是空间频率，公式(4)表示双谱中含有傅里叶信息。

在本申请的一个实施例中，根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位，包括：

对公式(4)取相位形式为公式(5)：

β(u，v)＝φ(u)+φ(v)-φ(u+v)；

其中，φ(u)是O(u)的相位，β(u，v)是双谱相位，通过公式(5)与递归算法从散斑图的双谱相位中恢复出原始目标的傅里叶相位信息φ(u)，使用φ(u)作为高斯牛顿优化的初始值φ₀(u)。

公式(5)表示为线性结构：β＝Aφ，其中，

表示需要恢复的目标傅里叶相位，

表示双谱相位，

是每行只有三个非零元素的稀疏矩阵，它的三个非零元素包含两个1和一个-1，对应于(5)式等号右边三个相位元素的符号。

构造最小二乘目标函数来优化相位，从而降低相位与散斑数据双谱之间的误差，目标函数为公式(6)：

其中，

是一个对角加权矩阵，对公式(6)进行转换得到公式(7)：

使用高斯牛顿迭代法对(6)或(7)式求解得到高精度傅里叶域相位。通过高斯牛顿法优化傅里叶相位具有抗噪性强和精度高的优点，相比基于梯度的一阶优化方法具有更快的收敛速度，会缩短求解时间并降低计算成本

在本申请的一个实施例中，根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像，包括：对原始目标傅里叶域幅值信息和高精度傅里叶域相位进行逆傅里叶变换或在高斯牛顿优化傅里叶域相位时将高精度傅里叶域相位更改为原始目标图像。

具体地，原始目标的重建，重建可以通过逆傅里叶变换完成，或者在使用高斯牛顿优化相位信息时将优化目标从相位φ更改为原始图像O，直接在优化中重建原始目标。目标函数的形式包括但不限于：

其中，C是目标O元素的边界约束集。R(o)是正则化算子，α是正则化参数。

需要说明的是，可以利用已获得的原始目标傅里叶域幅值与相位直接进行逆傅里叶变换获得原始目标信息或在高斯牛顿优化傅里叶域相位时将优化目标从相位φ更改为原始图像O，即直接在优化的过程重建原始目标。

本申请实施例的抗散射成像方法，获取散斑图；计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位；其中，N为正整数；根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位；根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。由此，从单张散斑图即可重建原始目标信息，且幅值和相位恢复过程相互独立，相位恢复精度高，抗噪性强。

为了更加清楚描述上述过程，结合图2-图6进行详细描述，图2散射成像整体流程图；图3基于散斑自相关获取傅里叶幅值流程图；图4对单张散斑图分块示意图；图5基于双谱相位恢复和高斯牛顿优化获取傅里叶域相位流程图；图6重建原始目标图像流程图。

具体地，如图1所示，携带有原始目标信息的光线经过散射介质后形成多重散射光束，使用探测器捕获散射光得到散斑图，无法直观的从散斑图中获得原始目标的信息。

利用谱分析技术能够恢复隐藏在其中的信息。具体来讲，原始目标的幅值信息蕴含在散斑图的二阶相关函数中；原始目标的相位信息蕴含在散斑图的三阶相关函数中。他们分别通过功率谱和双谱承载傅里叶域幅值和相位信息。

具体地，获取散斑图后进行预处理，包括消除渐变光晕和平滑滤波。消除渐变光晕有效快速的方法是用散斑图除以其自身低频版本或使用Zernike多项式拟合。平滑滤波采用高斯核卷积的方法。

从预处理后的散斑图获取目标傅里叶域幅值的流程如图2所示。在基于记忆效应的散射成像模型中，探测器捕获的散斑图是目标图像与散射系统的点扩散函数卷积的结果，即：I＝0*S。其中，I为探测器平面捕获的散斑图，O为原始目标图像，S是散射系统的点扩散函数，“*”是卷积算子。对捕获的散斑图计算自相关，得到：

[I★I]＝[(O*S)★(O*S)]＝[(O★O)*(S★S)]； (1)

其中，“★”是自相关算子。上式说明散斑图计算自相关的结果等于原始目标图像自相关与散射系统点扩散函数自相关的卷积运算。由于散斑图的自相关(S★S)是一个6函数(本质上是宽带噪声的自相关)，因此(1)式可以化简成：

[I★I]＝[(O★O)]+C； (2)

其中，C是在计算散斑图的自相关过程中产生的常数项。(2)式说明散斑图的自相关近似等于原始目标图像的自相关。通过散斑图即可获得原始目标图像的自相关信息，再根据维纳辛钦定理，对原始目标图像的自相关进行傅里叶变换后取模值，可以得到目标的功率谱密度(傅里叶域幅值)。即|F{(O★O)}|＝|F{O}|²。

对散斑图分块，获取若干张自散斑图。分块示意图如图3所示。从双谱中恢复原始目标傅里叶相位需要多张散斑图，本方法用空间代替时间的思想，将采集的单张散斑图进行空间上的分块，从而获得多张子散斑图，且每张子散斑都包含有原始目标的全部信息。

在分块的过程中，为充分利用散斑中的散斑颗粒，让相邻子散斑具有一定的重叠区域。为了消除子散斑不连续的边缘在离散二维快速傅里叶变换运算中的影响，还需要采用高斯窗函数对每块子散斑进行空间滤波，滤波函数的形式为：G＝exp[-(x²+y²/ω²)]。其中x和y表示空间坐标，ω用来控制窗函数窗口大小。

用分块后的若干张子散斑图恢复原始目标傅里叶域相位信息，该过程如图4所示。根据散斑图的双谱近似等于原始目标双谱这一性质，对散斑图进行双谱分析也就相当于对目标图像进行双谱分析。双谱的本质是三阶相关性。对于二维图像o(x)，其相对移位(x₁，x₂)的三阶相关表述为：

对(3)式进行傅里叶变换并结合卷积的性质可以得到其双谱：

O⁽³⁾(u，v)＝O(u)O(v)O*(u+v) (4)

其中，O⁽³⁾(u，v)是目标的双谱，O(u)是目标图像的傅里叶变换，

是空间频率。(4)式表示目标双谱中含有傅里叶信息。对(4)式取相位形式：

β(u，v)＝φ(u)+φ(v)-φ(u+v) (5)

其中，φ(u)是O(u)的相位即需要恢复出的傅里叶相位，β(u，v)是与空间三重频率(u，v，u+v)有关的双谱相位。(5)式提供了一种确定性的关系，通过该式与递归算法可以从散斑图的双谱中恢复出原始目标的傅里叶相位信息φ(u)。此处使用该结果作为高斯牛顿优化的初始值φ₀(u)

(5)式可以表示为线性结构：β＝Aφ。其中，

表示需要恢复的目标傅里叶相位，

表示双谱相位，

是每行只有三个非零元素的稀疏矩阵，它的三个非零元素包含两个1和一个-1，对应于(5)式等号右边三个相位元素的符号。可以构造最小二乘目标函数来优化相位，从而降低相位与散斑数据双谱之间的误差。目标函数的形式为：

其中，

是一个对角加权矩阵，其正对角元素由双谱相位β元素的信噪比决定。(6)式存在不可微的问题，利用欧拉公式对其进行修改来避免这一问题，即：

使用高斯牛顿迭代法对(6)或(7)式求解。初始值使用递归方法获得的初始相位φ₀(u)。每次迭代时需要计算迭代方向和步长。迭代的方向通过求解一个线性系统来确定。对于当前迭代y处的目标函数E(y)，求解迭代方向的方程为：

其中

是梯度。步长的选取基于回溯迭代的方法，起初先使用全步长η＝1，之后迭代回溯选择步长0＜η≤1，直到找到保证目标函数绝对减小的步长。目标函数的绝对减小由

来判定。c∈(0，1)是一个约束常数，默认值为c＝1×10^-4。经过高斯牛顿优化方案可以获得原始目标傅里叶域高精度相位信息。

重建原始目标图像的两种方法如图5所示。既可以逆傅里叶变换完成重建，还可以利用优化步骤。对于(6)和(7)式，将优化目标从相位φ更改为原始图像O，在优化过程中完成目标图像的重建。即：

min_o∈C{E₁(o)+αR(o)} (8)

min_o∈C{E₂(o)+αR(o)} (9)

其中，C表示恢复目标O元素的边界约束集。在实际中，通常为非负约束。R(o)是正则化算子，α＞0是正则化参数。

对于(8)和(9)的存在约束条件的优化问题，普通的高斯牛顿法已经无法适用。可以使用投影高斯牛顿法求解。投影高斯牛顿法在每次迭代中将优化变量分成两组，活动集合A和非活动集合I。活动集合A包含约束条件可行域边界上的量，非活动集合I包含了约束条件可行域内的量。对于非活动集合，通过投影版本来确定迭代方向，即：

其中，I_I是投影运算符。对于活动集合，采用投影梯度下降法确定迭代方向，即：

其中，I_A＝I-I_I，

是投影梯度

迭代步长采用采用与普通高斯牛顿法相同的策略，但对于可行域投影Q，判定目标函数绝对减小的条件修正为：

为了实现上述实施例，本申请还提出一种抗散射成像装置。

图7为本申请实施例提供的一种抗散射成像装置的结构示意图。

如图7所示，该装置包括：获取模块701、计算变换模块702、分块计算模块703、计算优化模块704和重建模块705。

获取模块701，用于获取散斑图。

计算变换模块702，用于计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息。

分块计算模块703，用于对所述散斑图进行分块处理得到N张子散斑图，并利用所述N张子散斑图计算双谱得到双谱相位；其中，N为正整数。

计算优化模块704，用于根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位。

重建模块705，用于根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像。

需要说明的是，前述对方法实施例的解释说明也适用于该实施例的装置，此处不再赘述。

本申请实施例的抗散射成像装置，获取散斑图；计算散斑图的自相关信息，对自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；对散斑图进行分块处理得到N张子散斑图，并利用N张子散斑图计算双谱得到双谱相位；其中，N为正整数；根据双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对初始估计值进行迭代优化得到高精度傅里叶域相位；根据原始目标傅里叶域幅值信息和高精度傅里叶域相位重建原始目标图像。由此，从单张散斑图即可重建原始目标信息，且幅值和相位恢复过程相互独立，相位恢复精度高，抗噪性强。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本申请的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。

此外，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者隐含地包括至少一个该特征。在本申请的描述中，“多个”的含义是至少两个，例如两个，三个等，除非另有明确具体的限定。

流程图中或在此以其他方式描述的任何过程或方法描述可以被理解为，表示包括一个或更多个用于实现定制逻辑功能或过程的步骤的可执行指令的代码的模块、片段或部分，并且本申请的优选实施方式的范围包括另外的实现，其中可以不按所示出或讨论的顺序，包括根据所涉及的功能按基本同时的方式或按相反的顺序，来执行功能，这应被本申请的实施例所属技术领域的技术人员所理解。

在流程图中表示或在此以其他方式描述的逻辑和/或步骤，例如，可以被认为是用于实现逻辑功能的可执行指令的定序列表，可以具体实现在任何计算机可读介质中，以供指令执行系统、装置或设备(如基于计算机的系统、包括处理器的系统或其他可以从指令执行系统、装置或设备取指令并执行指令的系统)使用，或结合这些指令执行系统、装置或设备而使用。就本说明书而言，″计算机可读介质″可以是任何可以包含、存储、通信、传播或传输程序以供指令执行系统、装置或设备或结合这些指令执行系统、装置或设备而使用的装置。计算机可读介质的更具体的示例(非穷尽性列表)包括以下：具有一个或多个布线的电连接部(电子装置)，便携式计算机盘盒(磁装置)，随机存取存储器(RAM)，只读存储器(ROM)，可擦除可编辑只读存储器(EPROM或闪速存储器)，光纤装置，以及便携式光盘只读存储器(CDROM)。另外，计算机可读介质甚至可以是可在其上打印所述程序的纸或其他合适的介质，因为可以例如通过对纸或其他介质进行光学扫描，接着进行编辑、解译或必要时以其他合适方式进行处理来以电子方式获得所述程序，然后将其存储在计算机存储器中。

应当理解，本申请的各部分可以用硬件、软件、固件或它们的组合来实现。在上述实施方式中，多个步骤或方法可以用存储在存储器中且由合适的指令执行系统执行的软件或固件来实现。如，如果用硬件来实现和在另一实施方式中一样，可用本领域公知的下列技术中的任一项或他们的组合来实现：具有用于对数据信号实现逻辑功能的逻辑门电路的离散逻辑电路，具有合适的组合逻辑门电路的专用集成电路，可编程门阵列(PGA)，现场可编程门阵列(FPGA)等。

本技术领域的普通技术人员可以理解实现上述实施例方法携带的全部或部分步骤是可以通过程序来指令相关的硬件完成，所述的程序可以存储于一种计算机可读存储介质中，该程序在执行时，包括方法实施例的步骤之一或其组合。

此外，在本申请各个实施例中的各功能单元可以集成在一个处理模块中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个模块中。上述集成的模块既可以采用硬件的形式实现，也可以采用软件功能模块的形式实现。所述集成的模块如果以软件功能模块的形式实现并作为独立的产品销售或使用时，也可以存储在一个计算机可读取存储介质中。

上述提到的存储介质可以是只读存储器，磁盘或光盘等。尽管上面已经示出和描述了本申请的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本申请的限制，本领域的普通技术人员在本申请的范围内可以对上述实施例进行变化、修改、替换和变型。

Claims (10)

1.一种抗散射成像方法，其特征在于，包括以下步骤：

获取散斑图；

计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；

对所述散斑图进行分块处理得到N张子散斑图，并利用所述N张子散斑图计算双谱得到双谱相位；其中，N为正整数；

根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位；

根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像。

2.如权利要求1所述的抗散射成像方法，其特征在于，所述获取散斑图，包括：

获取具有原始目标信息光线经过散射介质形成的多重散射光束；

使用探测器采集所述多重散射光束生成所述散斑图。

3.如权利要求1所述的抗散射成像方法，其特征在于，在所述计算所述散斑图的自相关信息之前，还包括：

将所述散斑图除以其自身低频版本；或，对所述散斑图进行多项式拟合；

采用高斯核卷积算法对所述散斑图进行平滑滤波处理。

4.如权利要求1所述的抗散射成像方法，其特征在于，所述计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息，包括：

通过以下公式(1)计算所述散斑图的自相关信息：

[I★I]＝[(O*S)★(O*S)]＝[(O★O)*(S★S)]；

其中，I为探测器平面捕获的散斑图，O为原始目标图像，S是散射系统的点扩散函数；I＝O*S，“*”是卷积算子，“★”是自相关算子；

根据维纳辛钦定理，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息为|F{(O★O)}|＝|F{O}|²。

5.如权利要求4所述的抗散射成像方法，其特征在于，还包括：

对所述公式(1)进行简化处理得到公式(2)：

[I★I]＝[(O★O)]+C；

其中，C是在计算散斑图的自相关过程中产生的常数项。

6.如权利要求1所述的抗散射成像方法，其特征在于，在所述利用所述N张子散斑图计算双谱得到双谱相位之前，还包括：

采用高斯窗函数对每块子散斑进行空间滤波；其中，所述高斯窗函数为G＝exp[-(x²+y²/ω²)，x和y表示空间坐标，ω用来控制窗函数窗口大小。

7.如权利要求1所述的抗散射成像方法，其特征在于，利用所述N张子散斑图计算双谱得到双谱相位，

对于子散斑图o(x)，相对移位(x₁,x₂)的三阶相关为公式(3)：

对所述公式(3)式进行傅里叶变换并结合卷积的性质可以得到所述双谱为公式(4)：

O⁽³⁾(u,v)＝O(u)O(v)O^*(u+v)；

其中，O⁽³⁾(u,v)是所述双谱，O(u)是目标图像的傅里叶变换，

是空间频率，所述公式(4)表示所述双谱中含有傅里叶信息。

8.如权利要求7所述的抗散射成像方法，其特征在于，所述根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位，包括：

对所述公式(4)取相位形式为公式(5)：

β(u,v)＝φ(u)+φ(v)-φ(u+v)；

其中，φ(u)是O(u)的相位，β(u,v)是所述双谱相位，通过所述公式(5)与递归算法从散斑图的双谱相位中恢复出原始目标的傅里叶相位信息φ(u)，使用φ(u)作为高斯牛顿优化的初始值φ₀(u)；

所述公式(5)表示为线性结构：β＝Aφ，其中，

表示需要恢复的目标傅里叶相位，

表示双谱相位，

是每行只有三个非零元素的稀疏矩阵，它的三个非零元素包含两个1和一个-1，对应于(5)式等号右边三个相位元素的符号；

构造最小二乘目标函数来优化相位，从而降低相位与散斑数据双谱之间的误差，目标函数为公式(6)：

其中，

是一个对角加权矩阵，对所述公式(6)进行转换得到公式(7)：

使用高斯牛顿迭代法对(6)或(7)式求解得到所述高精度傅里叶域相位。

9.如权利要求8所述的抗散射成像方法，其特征在于，所述根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像，包括：

对所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位进行逆傅里叶变换或在高斯牛顿优化傅里叶域相位时将所述高精度傅里叶域相位更改为所述原始目标图像。

10.一种抗散射成像装置，其特征在于，包括：

获取模块，用于获取散斑图；

计算变换模块，用于计算所述散斑图的自相关信息，对所述自相关信息进行傅里叶变换后取模值获取原始目标傅里叶域幅值信息；

分块计算模块，用于对所述散斑图进行分块处理得到N张子散斑图，并利用所述N张子散斑图计算双谱得到双谱相位；其中，N为正整数；

计算优化模块，用于根据所述双谱相位计算原始目标傅里叶域相位的初始估计值，并利用高斯牛顿优化方法对所述初始估计值进行迭代优化得到高精度傅里叶域相位；

重建模块，用于根据所述原始目标傅里叶域幅值信息和所述高精度傅里叶域相位重建原始目标图像。

Images