CN110097634A - 一种自适应多尺度的三维鬼成像方法 - Google Patents

Abstract

本发明公开了一种自适应多尺度的三维鬼成像方法，首先通过鬼成像系统获得目标不同方位的低分辨率鬼像；将低采样率下得到的图像结果进行预处理，滤除噪声；然后利用预处理过后的二维鬼像计算目标表面法向量；以表面法向量的稀疏性作为划分高低分辨率区域的标准划分区域；根据划分出的高分辨率区域生成高分辨率成像的参考光场，并获得高分辨率鬼像；最后将高低分辨率区域的图像拼接，并进行三维重构，从而获得目标三维图像。本发明减少鬼成像系统的采样次数，提高了系统的运算速度；减少了三维成像的测量数据，降低数据冗余，有利于数据的传输和存储。

Description

一种自适应多尺度的三维鬼成像方法

技术领域

本发明属于三维成像技术领域，特别是一种自适应多尺度的三维鬼成像方法。

背景技术

鬼成像(ghost imaging，GI)又称关联成像(correlated imaging)或双光子成像(two-photon imaging)，是一种备受关注的区别于传统成像概念的成像技术。在鬼成像系统中存在两条光路，信号光路以及参考光路。信号光路中，利用一个桶探测器探测通过目标物体的总光强。而在参考光路中，则利用具有空间分辨率的探测器探测光源的空间强度分布。通过两个探测器探测到的光强作交叉关联来重建图像。计算鬼成像是一种改进的鬼成像技术，该技术利用光调制器件如空间光调制器(SLM)和数字微镜阵列(DMD)与鬼成像相结合，实现了仅凭单臂光路即可成像。这种方法简化了鬼成像光路，使得鬼成像从实验室的研究走向实际应用成为了可能。鬼成像相比于传统成像有很多独特的优势和特点，如抗散射和抗大气湍流等，该技术已被应用于许多领域，例如雷达探测、光学相干断层扫描(OCT)和荧光成像。当然，鬼成像技术也存在一些局限。最明显的是该技术的成像时间较长，主要包括前期的数据采集时间和后期的算法处理时间。为了缓解这一问题，基于多尺度思想的图像研究方法也被引入鬼成像。例如Marc Abmann提出了一种基于小波树的自适应压缩计算式鬼成像方案(Compressive adaptive computational ghost imaging[J].Scientificreports,2013,3:1545)，受到了学界的广泛关注。

三维重建技术是当今计算机视觉领域的前沿研究问题，它可以根据自然界中物体的二维信息重建出具有详细几何信息的三维模型。这些模型不仅能为考古、文物保护、城市规划等领域提供可视化场景，更可以为工业测量、游戏、医疗诊断等领域提供数据的测量、分析和存档。

近年来，将鬼成像与其他技术相结合的研究层出不穷，也出现了基于鬼成像的三维成像技术(Improving the noise immunity of3D computational ghost imaging[J].Optics Express,2019,27(3):2344-2353)。但是结合鬼成像的三维重构技术同样受到高采样数以及大量后期运算的影响。同时，获得图像的分辨率越高，系统所需的数据量就越大。因此，三维鬼成像系统的分辨率与采样率是一对矛盾。

发明内容

本发明的目的在于提供一种自适应多尺度的三维鬼成像方法，解决了鬼成像需要大量采样数据以重构目标的高分辨率三维图像这一技术难题，能够缓解分辨率与采样率之间的矛盾。

实现本发明目的的技术解决方案为：一种自适应多尺度的三维鬼成像方法，步骤如下：

步骤1，通过鬼成像系统获得目标不同方位的低分辨率鬼像I_up、I_left、I_right；

步骤2，对低采样率下的鬼成像结果进行预处理，即对图像进行傅立叶变换，然后进行高斯低通滤波，滤除图像高频部分，降低噪声；

步骤3，利用预处理过后的二维鬼像计算目标表面法向量；

步骤4，以表面法向量的稀疏性作为划分高低分辨率区域的标准进行区域划分，即法向量稀疏区域低分辨率成像，法向量不稀疏区域高分辨率成像；

步骤5，根据划分出的高分辨率区域，生成相应的基于Hadamard测量的参考光场；

步骤6，通过鬼成像系统获得目标不同方位的高分辨率鬼像I¹ _up、I¹ _left、I¹ _right；

步骤7，获得分辨率自适应的二维鬼像I² _up、I² _left、I² _right，用I¹ _up、I¹ _left、I¹ _right中有效的高分辨率图像信息替换低采样图像中的对应区域；

步骤8，进行三维重构，获得目标三维图像，根据步骤7中获得的分辨率自适应的二维鬼像I² _i，采用光度立体法恢复目标的三维形态。

本发明与现有的技术相比，其显著优点为：(1)减少三维鬼成像系统的采样次数，提高了系统的速度。相对传统的三维鬼成像方法，自适应多尺度的三维鬼成像技术在保留目标细节的同时，减少了平坦区域的采样数。(2)减少了三维成像的测量数据，在采样的过程中实现了数据压缩，降低数据冗余，有利于后续数据的传输和存储。(3)继承了鬼成像的优势，与其他成像方式相比，具有抗散射能力。(4)将三维成像与鬼成像技术相结合，使鬼成像技术得以实际应用。

下面结合附图对本发明作进一步详细说明。

附图说明

图1为三维鬼像成像系统图。

图2为本发明自适应多尺度的三维鬼成像方法的流程图。

图3为自适应区域划分图像：(a)高低分辨率区域划分二值图像m(x,y)；(b)基于Hadamard的高分辨率参考光场示意图。

图4为成像系统结果图像：(a)低分辨率采样下探测器位于目标上、左、右三个位置时目标的低分辨率二维鬼像I_up，I_left，I_right；(b)高分辨率采样探测器位于目标上、左、右三个位置时目标高分辨率的二维鬼像I¹ _up，I¹ _left，I¹ _right；(c)高低分辨率图像拼接结果I² _up，I² _left，I² _right。

图5为本发明恢复的目标三维结果图。

具体实施方式

本发明结合鬼成像和机器视觉中三维成像技术，将鬼成像系统低采样率下获得的目标不同方位的图像进行预处理，然后计算目标表面法向量，并根据表面法向量划分高低分辨率区域，对划分出的高分辨率区域生成相应的参考光场，从而获得相应区域的高分辨率鬼像，最后拼接图像并进行三维重构，得到目标三维图像。

结合图2，本发明自适应多尺度的三维鬼成像方法的具体实施步骤如下：

步骤1，通过鬼成像系统获得目标不同方位的低分辨率鬼像I_up，I_left，I_right。三维鬼像成像系统图如图1所示，系统主要包括成像目标6，分别位于该成像目标6上方的第一桶探测器1、左方的第二桶探测器2以及右方的第三桶探测器3，还有一个连接计算机8的投影仪4用于投射参考光场5至成像目标6处，采集卡7将各桶探测器输出的模拟信号转化为数字信号后传送给计算机8。三维鬼像成像系统工作过程如下：计算机8预先生产一系列基于Hadamard的参考光场，使用投影仪4将其逐帧投射至成像目标6处，然后用不同方位的第一桶探测器1、第二桶探测器2、第三桶探测器3分别接收成像目标6反射的光强并传送给计算机8，该计算机8根据鬼像公式计算出成像目标6的不同方位的二维阴影图像。图像计算公式如式(1)：

I_i(x,y)＝g⁽²⁾(A(x,y),s_i)＝＜A(x,y)·s_i>-<A(x,y)><s_i> (1)

其中x,y为平面坐标，i为探测器标号(i＝up,left,right)，I(x,y)为目标鬼像，A(x,y)为参考光场，s为桶探测器测得的信号，<…>表示统计平均值。图4(a)为第一桶探测器1、第二桶探测器2、第三桶探测器3测得结果恢复的低分辨率二维鬼像。

基于Hadamard的参考光场生成方法如下：以获取分辨率为2^N×2^N(N＝1,2,3,…)的目标鬼像为例，首先获取大小为2^2N×2^2N的Hadamard矩阵H(i,j)(1≤i≤2^2N,1≤j≤2^2N)。Hadamard矩阵由+1和-1组成，实际使用投影仪4时，通常用亮代表1，暗代表0，因此将Hadamard分为正负图案两部分。将H(i,j)中的-1替换成0，得到Hadamard的正图案H⁺(i,j)；将-H(i,j)中的-1替换成0得到Hadamard的负图案H^-(i,j)。则

H(i,j)＝H⁺(i,j)-H^-(i,j) (2)

正图案中的第i行H⁺(i,:)重排成2^N×2^N的矩阵得到H⁺ _i(:,:)，负图案中的第i行H-(i,:)重排成2^N×2^N的矩阵得到H^- _i(:,:)，则参考光场A(x,y)中的第i帧可表示为H⁺ _i(:,:)-H^- _i(:,:)。正图案投射到参考目标6上后第一桶探测器1、第二桶探测器2、第三桶探测器3测得信号为s⁺，负图案同理可得s^-，鬼像计算所需信号s＝s⁺-s^-。

步骤2，对低采样率下的鬼成像结果进行预处理。低采样率下鬼成像图像信噪比低，不能直接用于目标三维信息的提取，所以需要对图像进行预处理。考虑到噪声在图像中多表现为尖锐的灰度变化，而该部分对其傅里叶变换的高频内容有贡献，因此降低噪声，可先获取图像的傅立叶变换，然后进行高斯低通滤波，衰减高频，再将滤波后的图像傅里叶信息做傅里叶逆变换，恢复空间域的图像，以达到除去噪声的效果，式(3)给出了滤波器的二维形式。

式(3)中u、v为图像在频率域中的坐标，D₀为滤波器的截至频率。图像在频域中进行滤波，需要进行填充，为了使信号强度损失最小，选取D₀的值为填充图像大小的10％。

步骤3，利用预处理过后的二维鬼像计算目标表面法向量。在这里，采用光度立体法与从阴影恢复形状的原理相结合的思想来计算目标表面法向量。光度立体法利用相同视点但由不同光照拍摄的阴影图像实现物体的三维形状，是三维成像领域的一种经典方法。根据从阴影恢复形状的思想，假设物体表现出均匀的朗伯反射率并且光源位于无穷远处。亮度的变化(辐照度方程)变成单纯的局部表面方向的函数：

I(x,y)＝R(p(x,y),q(x,y)) (4)

其中R(p,q)为反射图，(p,q)＝(z_x,z_y)是表面梯度，即深度图的导数。

根据光学原理，探测器(这里的探测器指光学领域的所有探测器，包含本发明使用的桶探测器)获得的图像灰度主要由以下因素决定：目标形状(即目标表面法向量)、目标相对于光源和探测器的位置以及探测器与光源的相对位置。形成的阴影图像中像素的强度可表示为：

其中I_s表示光源的光强，α表示目标表面反照率，对于朗伯体来说是一个常数，代表物体表面指向第k个探测器的单位法向量，则代表物体表面的法向量， 可以表示为：

其中D是第一桶探测器1、第二桶探测器2、第三桶探测器3单位法矢量组成的矩阵，I是对应的图像强度的数组，I＝[I_up,I_left,I_right]^T。

表面梯度(p,q)可以通过获得：

则目标表面法向量用z分量归一化的形式可表示为

步骤4，以表面法向量的稀疏性作为划分高低分辨率区域的标准进行区域划分，即法向量稀疏(平坦)区域低分辨率成像，法向量不稀疏(非平坦)区域高分辨率成像。划分区域的主要依据是局部区域是否平坦，平坦区域可以采用低分辨率成像，非平坦区域则采用高分辨率恢复，以确保保留非平坦区域细节的同时相对地减少采样数，降低数据冗余程度。目标的表面法向量包含了物体表面的平坦度信息，因此目标的平坦程度可以表示为

G_pq＝|p_x|+|p_y|+|q_x|+|q_y| (9)

其中p_x，p_y，q_x和q_y分别代表和用一个3x3(像素)的窗口扫描G_pq的结果图以估计目标局部区域的平坦程度：

ν(x,y)的值越小代表(x,y)附近局部区域相对越平坦。根据目标非平坦区域的相对大小设置阈值T，使得ν(x,y)＞T划分出的区域像素个数为原图像的(k＝1,2,3,…)，可将ν(x,y)的值按从大到小排序之后以第2^2N-k+1个的值作为T。划分所得到的结果以二值矩阵m(x,y)表示：

m(x,y)中值为1的区域为高分辨率区域，值为0的区域则为低分辨率区域。m(x,y)结果如图3(a)所示，其中白色区域值为1，黑色区域值为0。

步骤5，根据划分出的高分辨率区域，生成相应的基于Hadamard的测量(参考)光场。首先，实现参考光场分辨率的提升，即将低分辨率参考光场中的每个像素划分为2^a×2^a(a＝1,2,3,…)个像素，分辨率提升2^2a倍，此时参考光场的分辨率为2^L×2^L(L＝N+a，N为步骤1中所提N)。具体做法为求m(x,y)与维度为2^a的单位矩阵的克罗内克(Kronecker)积，得到高分辨率下划分结果的二值矩阵M(X,Y)，此过程可表示为：

M(X,Y)中值为1的部分为高分辨率区域，即需要进行高分辨率成像的区域；值为0的区域则不需要再进行高分辨率成像。

将矩阵M(X,Y)向量化后得到M(n)，点(X,Y)与n的关系为

然后生成一个大小为2^T×2^T(T＝2L-k，k为步骤4中提到的划分出的区域像素个数为原图像的中的k)Hadamard矩阵H_T，则高分辨率的测量感知矩阵H_L(m,n)为

矩阵H_L的行号m(1≤m≤T)表示参考光场的序号，n对应参考光场上的坐标点(X,Y)，n₁满足：

获得高分辨率的测量感知矩阵H_L(m,n)，仿照步骤1中基于Hadamard矩阵H(i,j)生成参考光场的方法，分正负图案，按行取出重排，得到T帧大小为2^L×2^L(像素)的高分辨率参考光场。图3(b)展示了依据M(X,Y)产生的基于Hadamard的高分辨参考光场的其中三帧，可以看出只有需要进行高分辨率成像的区域(M(X,Y)中值为1的区域)的参考光场有明暗变化；而低分辨率区域(M(X,Y)中值为0的区域)为全黑，即不再对该区域进行采样。

步骤6，通过鬼成像系统获得目标不同方位的高分辨率鬼像I¹ _up，I¹ _left，I¹ _right。成像方法与步骤1描述方法相同，区别在于投影图案以及鬼像计算时的参考图案采用步骤5中所得的更高分辨率的基于Hadamard的参考光场。上、左、右三个方向的成像结果如图4(b)所示。

步骤7，获得分辨率自适应的二维鬼像I² _up，I² _left，I² _right，用I¹ _up，I¹ _left，I¹ _right中有效的高分辨率图像信息替换低采样图像中的对应区域，结果如图4(c)所示。根据区域划分信息M(x,y)找到重叠部分，对高分辨率图像做灰度调整后将其替换至相应位置。首先根据步骤5中的方法调整步骤1中获得的低分辨率图像I_up，I_left，I_right的大小，得到对应的尺寸为2^L×2^L像素的图像I'_up，I'_left，I'_right。然后调整图像灰度，灰度调整系数κ_i可表示为：

其中*表示数组元素对应相乘，则分辨率自适应的二维鬼像I² _i的表达式为：

I² _i＝(1-M(x,y))*I'_i(x,y)+κ_i·M(x,y)*I¹ _i(x,y) (17)

步骤8，进行三维重构，获得目标三维图像，根据步骤7中获得的分辨率自适应的二维鬼像I² _i，采用光度立体法恢复目标的三维形态，即根据步骤3中计算方法得到每个像素的表面梯度(p,q)，即深度图的导数(z_x,z_y)，然后利用得到的梯度图恢复目标表面各点的高度，选取目标区域中的一个像素作为起始点，设定其高度为0，然后向外拓展恢复与其相邻的四个像素的高度。由于每个像素的梯度是已知的，因此可以基于其邻近点的已知高度和像素本身的梯度来估计像素的高度。重复上述过程恢复下一个邻近像素的高度，直到获得所有像素的高度，由此获得目标的三维图像。最终由分辨率自适应的二维鬼像恢复的目标三维形态如图5所示。

Claims (8)

1.一种自适应多尺度的三维鬼成像方法，其特征在于步骤如下：

步骤1，通过鬼成像系统获得目标不同方位的低分辨率鬼像I_up、I_left、I_right；

步骤2，对低采样率下的鬼成像结果进行预处理，即对图像进行傅立叶变换，然后进行高斯低通滤波，滤除图像高频部分，降低噪声；

步骤3，利用预处理过后的二维鬼像计算目标表面法向量；

步骤4，以表面法向量的稀疏性作为划分高低分辨率区域的标准进行区域划分，即法向量稀疏区域低分辨率成像，法向量不稀疏区域高分辨率成像；

步骤5，根据划分出的高分辨率区域，生成相应的基于Hadamard测量的参考光场；

步骤6，通过鬼成像系统获得目标不同方位的高分辨率鬼像I¹ _up、I¹ _left、I¹ _right；

步骤7，获得分辨率自适应的二维鬼像I² _up、I² _left、I² _right，用I¹ _up、I¹ _left、I¹ _right中有效的高分辨率图像信息替换低采样图像中的对应区域；

步骤8，进行三维重构，获得目标三维图像，根据步骤7中获得的分辨率自适应的二维鬼像I² _i，采用光度立体法恢复目标的三维形态。

2.根据权利要求1所述的方法，其特征在于步骤1中，运用鬼成像系统对目标进行成像，计算机预先生产一系列基于Hadamard的参考光场，使用投影仪将其逐帧投射至目标处；使用三个不同方位的桶探测器分别接收目标反射的光强；根据鬼像公式计算出目标的不同方位的二维阴影图像，图像计算公式如式(1)：

I_i(x,y)＝g⁽²⁾(A(x,y),s_i)＝<A(x,y)·s_i>-<A(x,y)><s_i> (1)

其中x,y为平面坐标，i为探测器标号(i＝up,left,right)，I(x,y)为目标鬼像，A(x,y)为参考光场，s为桶探测器测得的信号，<…>表示平均值。

3.根据权利要求2所述的方法，其特征在于基于Hadamard的参考光场生成方法如下：以获取分辨率为2^N×2^N的目标鬼像为例，N＝1,2,3，首先获取大小为2^2N×2^2N的Hadamard矩阵H(i,j)，1≤i≤2^2N,1≤j≤2^2N，Hadamard矩阵由+1和-1组成，使用投影仪时，用亮代表1，暗代表0，因此将Hadamard分为正负图案两部分，将H(i,j)中的-1替换成0，得到Hadamard的正图案H⁺(i,j)；将-H(i,j)中的-1替换成0得到Hadamard的负图案H^-(i,j)，则

H(i,j)＝H⁺(i,j)-H^-(i,j) (2)

正图案中的第i行H⁺(i,:)重排成2^N×2^N的矩阵得到H⁺ _i(:,:)，负图案中的第i行H^-(i,:)重排成2^N×2^N的矩阵得到H^- _i(:,:)，则参考光场A(x,y)中的第i帧可表示为H⁺ _i(:,:)-H^- _i(:,:)，正图案投射到目标上后三个桶探测器测得信号为s⁺，负图案同理可得s^-，鬼像计算所需信号s＝s⁺-s^-。

4.根据权利要求1所述的方法，其特征在于步骤3中，利用预处理过后的二维鬼像计算目标表面法向量，即采用光度立体法与从阴影恢复形状的原理相结合的思想来计算目标表面法向量，假设物体表现出均匀的朗伯反射率并且光源位于无穷远处，亮度的变化变成单纯的局部表面方向的函数：

I(x,y)＝R(p(x,y),q(x,y)) (4)

其中R(p,q)为反射图，(p,q)＝(z_x,z_y)是表面梯度，即深度图的导数；

形成的阴影图像中像素的强度表示为：

其中I_s表示光源的光强，α表示目标表面反照率，对于朗伯体来说是一个常数，代表物体表面指向第k个探测器的单位法向量，则代表物体表面的法向量， 表示为：

其中D是三个方位的桶探测器单位法矢量组成的矩阵，I是对应的图像强度的数组，I＝[I_up,I_left,I_right]^T；

表面梯度(p,q)通过获得：

则目标表面法向量用z分量归一化的形式可表示为

5.根据权利要求1所述的方法，其特征在于步骤4中，根据目标表面法向量，划分高低分辨率区域，平坦区域采用低分辨率成像，非平坦区域则采用高分辨率恢复，目标的表面法向量包含了物体表面的平坦信息，因此目标的平坦程度表示为

G_pq＝|p_x|+|p_y|+|q_x|+|q_y| (9)

其中p_x，p_y，q_x和q_y分别代表和用一个3x3的窗口扫描G_pq的结果图以估计目标局部区域的平坦程度：

根据目标非平坦区域的相对大小设置阈值T，使得ν(x,y)＞T划分出的区域像素个数为原图像的k＝1,2,3，将ν(x,y)的值按从大到小排序之后以第2^2N-k+1个的值作为T，划分所得到的结果以二值矩阵m(x,y)表示：

m(x,y)中值为1的区域为高分辨率区域，值为0的区域则为低分辨率区域。

6.根据权利要求1所述的方法，其特征在于步骤5中，根据划分出的高分辨率区域，生成相应的基于Hadamard测量的参考光场：首先，实现参考光场分辨率的提升，即将低分辨率参考光场中的每个像素划分为2^a×2^a个像素，a＝1,2,3,…，分辨率提升2^2a倍，此时参考光场的分辨率为2^L×2^L，L＝N+a，N为步骤1中所提N，具体做法为求m(x,y)与维度为2^a的单位矩阵的克罗内克积，得到高分辨率下划分结果的二值矩阵M(X,Y)表示：

M(X,Y)中值为1的部分为高分辨率区域，即需要进行高分辨率成像的区域；值为0的区域则不需要再进行高分辨率成像；

将矩阵M(X,Y)向量化后得到M(n)，点(X,Y)与n的关系为

然后生成一个大小为2^T×2^T，Hadamard矩阵H_T，T＝2L-k，则高分辨率的测量感知矩阵H_L(m,n)为

矩阵H_L的行号m表示参考光场的序号，1≤m≤T，n对应参考光场上的坐标点(X,Y)，n₁满足：

获得高分辨率的测量感知矩阵H_L(m,n)，仿照步骤1中基于Hadamard矩阵H(i,j)生成参考光场的方法，分正负图案，按行取出重排，得到T帧大小为2^L×2^L的高分辨率参考光场。

7.根据权利要求1所述的方法，其特征在于步骤7中，用I¹ _up，I¹ _left，I¹ _right中有效的高分辨率图像信息替换低采样图像中的对应区域，获得分辨率自适应的二维鬼像I² _up，I² _left，I² _right，根据区域划分信息找到重叠部分，对高分辨率图像做灰度调整后将其替换至相应位置，即首先根据步骤5中的方法调整低分辨率图像I_up，I_left，I_right的大小，得到对应的尺寸为2^L×2^L像素的图像I'_up，I'_left，I'_right，图像灰度调整系数κ_i表示为：

其中*表示数组元素对应相乘，则分辨率自适应的二维鬼像I² _i的表达式为：

I² _i＝(1-M(x,y))*I'_i(x,y)+κ_i·M(x,y)*I¹ _i(x,y)。 (17)

8.根据权利要求1所述的方法，其特征在于步骤8中，根据步骤7中获得的分辨率自适应的二维鬼像I² _i，采用光度立体法恢复目标的三维形态，即根据步骤3中计算方法得到每个像素的表面梯度(p,q)，即深度图的导数(z_x,z_y)，然后利用得到的梯度图恢复目标表面各点的高度，选取目标区域中的一个像素作为起始点，设定其高度为0，然后向外拓展恢复与其相邻的四个像素的高度；重复上述过程恢复下一个邻近像素的高度，直到获得所有像素的高度。