CN108205800B - 一种快速有效的像素化反卷积技术 - Google Patents

Abstract

fiDrizzle像素化反卷积技术是数字图像处理领域为数不多的包含去像素化过程的图像叠加技术。目的是要解决成像仪器的低采样率导致的像素化问题，也顺便解决姿态稳定性差导致的图像形状模糊等问题。我们在前人工作的基础上，改进了实现方法，省略了一些冗余的计算步骤，并直接瞄准输出网格的分辨率。这样做省时省力，而且保留了更多的图像细节，得到的图像质量更高。本技术快速有效地解决了比较基础性的低采样率导致的像素化问题，能获得更高保真的重建图像。所以在对图像的像素分辨率有极高要求的天文学，微观物理学等领域有广泛用途，甚至在微生物学，医学和数字监控领域也有用武之地。

Description

一种快速有效的像素化反卷积技术

September 30，2019

1技术领域

所有的采集技术都不能包含无限的细节信息，因为图像的几何尺寸和数据存储容量都是有限的。对于在像平面上尺度很小的物体(半高宽FWHM＜10像素)来说，细节往往会因为采样率低(undersampling)而丧失，也就是我们说的放大之后全是马赛克-像素化(pixelization)。如何恢复或者重建因为像素化而丢失的细节(也就是像素化反卷积)是我们要解决的问题。

2背景技术

科学家发现，虽然某些观测已经达到了分辨率极限，不可避免地会出现低采样率问题，但通过记录仪器的震颤(dither)在不同位置拍摄目标物体，就会得到很多低采样率样本，而把这些样本通过合理的技术进行叠加就可以得到突破仪器限制的超高分辨率。也就是从某种程度上实现了像素化的反卷积。在通过叠加低采样率图像得到高采样率图像的技术发展过程中涌现出了诸如Interlace，shift-and-add，Drizzle和iDrizzle等方法。Fruchter等人发展的Drizzle和iDrizzle技术已经成功用在哈勃望远镜宽场图像的叠加拼接上。本文中我们主要用我们新发展的像素化反卷积方法跟Drizzle和iDrizzle做对比。

Drizzle的特点是速度快，但内秉的滤波器(filter)并没有被反卷积出来，所以Drizzle只是图像在高分辨率网格下的单纯叠加，不涉及反卷积，因而也没有很好地去像素化。随后在Drizzle的基础上，Fruchter于2011年又发展出了多次迭代结合傅立叶空间平滑的方法来多次提取残留信号，也即iDrizzle技术。iDrizzle首先要通过Drizzle技术把多次观测到原始图像放到一个过采样的网格上，再多次与原始图像对比，从残差图像中提取大部分信号，并加到结果中，最后通过sinc插值到目标分辨率，比如不高于临界采样的分辨率(critical sampling，符合香农-奈奎斯特采样定理的采样)。iDrizzle能比较好地实现反卷积，给出的图像跟真实图像差不太多。

但随后的研究中我们发现iDrizzle会在傅立叶变换和反变换中消耗大量运算资源，另外一开始就把原观测图像叠加到过采样的高分辨率网格上，最终再差值回来，也会加重计算负载。而且傅立叶变换和反变换之间要进行平滑，最后的sinc差值也相当于一次平滑，因此会最终损失图像在高频部分的信息(细节)。下面我们阐述我们新发展的快速有效像素化反卷积技术(a fast iDrizzle method，简称fiDrizzle)。

3发明内容

3.1fiDrizzle技术

我们的目的是尽量保留和发展iDrizzle的反卷积功能，并简化步骤和合理采样，来增加计算速度。fiDrizzle的步骤如下：

步骤(1)：获取原始观测图像

将对准某区域观测获取的N幅震颤过的原始图像标记为{g¹(m，n)，g²(m，n)...g^N(m，n)}，用g^k(m，n)表示其中的第k幅，N，m，n和k为正整数，1≤k≤N，(m，n)为原始图像像素位置；

步骤(2)：重建图像

将所述原始观测图像{g¹(m，n)，g²(m，n)，g³(m，n)…g^N(m，n)}，用Drizzle方法叠加到一个预定的高分辨率标准网格上，得到重建图像f₁(x，y)，也即升采样操作，f的下标1代表这是第1次重建图像的结果，(x，y)为高分辨率的像素位置，x，y为正整数；

步骤(3)：迭代处理

步骤(3.1)：对第i次图像重建所得的高分辨率图像f_i(x，y)进行N次与步骤(1)相同采样条件的重新采样，也即降采样操作，从而获得N幅图像

用

表示N幅重新采样图像中的第k幅，其中i为自然数，f₁(x，y)也就是步骤(2)获得的重建图像；

所述重新采样图像

与原始观测图像g^k(m，n)一一对应，区别在于前者是对高分辨率重建图像f_i(x，y)的重新采样，后者是真实的观测图像，之后进入下一步；

步骤(3.2)：用N幅原始图像g^k(m，n)分别减去对应的重新采样图像

得到N幅残差图像

其中进入下一步；

步骤(3.3)：将N幅残差图像

用Drizzle方法叠加到上次重建图像f_i(x，y)的标准网格上，从而得到新的重建图像f_i+1(x，y)，若f_i+1(x，y)不满足要求，则进入下一步；

步骤(3.4)：返回步骤(3.1)，进行下一次迭代循环操作，直至输出满足要求的图像。

迭代结束后，我们就可以把输出当成最终结果了。对比iDrizzle的步骤，我们的fiDrizzle做到了如下改进：1.省略了正反傅立叶变换和其中的平滑过程；2.一开始就把输出图像精度当成目标分辨率，而不是用过采样的分辨率。3.因为没有用过采样分辨率，所以省去了最终使用sinc插值的过程。4.因为去掉了平滑和插值步骤，我们期待新的fiDrizzle能保留一些高频信息，也就是重建更多的细节。下面我们对Drizzle，iDrizzle和fiDrizzle三种图像重建技术进行了视觉表观方面，定量方面和时间复杂度方面的全面对比。为了控制对比三种方法，我们首先给出一幅真实图像，然后把真实图像通过重新采样得到一些低分辨率原始图像(dithered images)，然后我们用这三种方法对这些低分辨率原始图像进行图像重建，最后输出结果进行对比。这其中为了模拟真实观测，震颤涉及到的仪器的位移，旋转是随机的，另外还考虑了每次曝光0.1％的CCD随机几何形变。

3.2结果1：视觉表观对比

在图1中，我们使用了图像处理领域广为人知的测试图像Miss Lena(尺寸512×512pixel)作为真实图像，左图。我们把真实图像进行重新采样，将5幅原始观测图像分辨率降低了

其中之一就是右图。我们可以看到低采样率或者说低分辨率图像中的Lena的眼睛和睫毛都模糊不清了，帽子边缘锯齿化严重，帽子上的条纹也消失了，这些都是像素化的结果，所有的细节都马赛克化了，也是我们现实中真实观测图像的表现。之后我们用上述三种方法把这5幅原始观测图像进行重建。如图2，左侧上图是真实图像。右侧上图是Drizzle重建的结果，左侧下图是iDrizzle的结果，右侧下图fiDrizzle的结果。这其中iDrizzle和fiDrizzle都进行了5次迭代。很明显Drizzle重建的图像已经比原始观测图像好多了，也增加了更多的细节。然而与iDrizzle和fiDrizzle对比，Drizzle就差多了。主要原因是Drizzle并没有进行像素化反卷积操作。对比iDrizzle和fiDrizzle我们发现fiDrizzle重建的图像包含更多的细节比如竟然重建了帽子上的条纹。而在真实图像中条纹的特征尺度也就是1-4个像素。说明fiDrizzle可以重建尺度在1-2个原始像素的细节了。另一方面iDrizzle的效果似乎也差不太多，主要区别在细节方面。为了进一步展示iDrizzle和fiDrizzle的区别我们又画了残差图像3。其中左上是真实图像，其余是各方法重建结果与真实图像相减后的残差。残差图像里信息越多说明图像重建效果越差。我们明显可以看到fiDrizzle的残差是最小的，除了噪音几乎看不到什么信号。

3.3结果2：定量对比

除了视觉表观对比之外，我们还进行了定量对比，主要是对比重建图像中光源的形状信息。我们从哈勃望远镜的数据库里HST-jc1g03010-drc.fits摘出一幅481×481分辨率的图像，其中星系的中心位置赤经赤纬为R.A.＝195.01238deg.和Dec＝28.023106deg.。我们将该图像作为真实图像，然后重新采样了10幅比真实图像低

分辨率的原始图像。在图4中，图像的布局跟图2是一样的。可以发现iDrizzle和fiDrizzle明显要比Drizzle好，但iDrizzle和fiDrizzle看起来区别不大。我们沿着X＝Y做了一个流量剖面图，来对比在X＝Y像素上的流量。图5中，所有像素上的流量已经除以真实图像中央像素的流量值(也是最大流量值)，其中上图是流量轮廓图，下图是流量残差图，也就是重建图像流量与真实图像流量相减的差值。可以明显地看到中心区域Drizzle的重建是最差的，iDrizzle次之，fiDrizzle是最好的。在***三种方法区别不大。我们再次看到Drizzle和iDrizzle对高频细节信息重建得不够好。

源在像平面的形状轮廓(或等强度轮廓线，contour shape)是我们辨别被观测源形态的主要依据。在微观物理学，微生物学，天文学等领域，被观测物体的成像尺度往往比较小，而形状信息却尤其重要。比如天体物理中弱引力透镜信号的测量尤其依赖像的形状。这里为了直观简洁，我们用描述二维图形的椭率[e₊，e_×，是计算剪切形变(shearestimator)的重要量]和图形主轴的方位角(position angle，P.A)来定量比较drizzle与fiDrizzle对源形状的重建效果。其中对于任意二维图形(像)的椭率和方位角是这样求得的：

1.首先我们用一个椭圆高斯的亮度轮廓来拟合这个像，也就是最小化下式

其中I(x)是像的亮度轮廓，这样就确定了(A，x₀，M)。我们定义椭圆高斯的权重函数：

2.其次最小化E也会得到像的中心和二阶矩：

3.最后像的椭率则是这样定义的：

e₊＝(M_xx-M_yy)/TM

e_×＝2M_xy/TM

TM＝(M_xx+M_yy)； (4)

值得一提的是椭率e₊主要是描述垂直和水平方向的，而e_×则是描述45度和135度方向的。

方位角：

为了避免发散问题，大家都会在测量形状的时候对测量图像乘以一个圆高斯权重函数权重半径是r_w。我们在图6，7和8中把相应的物理量表示为权重半径r_w的函数，以比较不同重建图像在各个尺度上的差异。我们发现在图6和7中，Drizzle几乎不能在全部区域精确重建图形的椭率。iDrizzle在中心区域表现不好，但在***很好。因为原始图像的像素尺寸跟真实图像的像素尺寸比例是3∶1的，所以图中3×3个像素相当于一个原始像素的尺寸。fiDrizzle在

原始像素的尺度上表现完美，另外即使在

原始像素的中心区域也比iDrizzle重建的精度提高2-3倍。在方位角的比较中也可以发现相似的结论，如图8。

4计算复杂度比较

毫无疑问的是Drizzle的速度是最快的。由于一次重新采样的计算量与一次Drizzle的计算量是相当的。所以fiDrizzle实际上比Drizzle多消耗了2(K-1)倍的计算量(K是迭代次数)。在图2中，真实图像具有512×512尺寸，那么原始图像就具有

像素的尺寸。Drizzle的计算量不仅依赖于5×181×181的原始像素数，还依赖输出图像的分辨率，因此Drizzle消耗了T(n)＝5×181×181×320

52,000,000次的计算量(平均来讲，每一个原始像素会被真实像素网格分为10个区域，每个区域会对应多边形切割操作24次和线积分8次运算)。这样它的复杂度是O(200n²)(令n＝512)。而fiDrizzle要多4次迭代，因此复杂度大概是O(2000n²)。

至于iDrizzle，除了复杂度为O(2000n²)的迭代运算量，还有其他的运算支出：

1.首先正反傅立叶变换要消耗的运算量至少是T(n)＝5×2×log₂(512×512)×512×512

47,000,000，也就是说复杂度为O(20n²log₂n).

2.其次，iDrizzle需要一个过采样的输出网格，也是能大量消耗CPU资源的。比如一个2×2倍高的过采样输出，就会带来一个复杂度为O(8000n²)的运算量，另外还有O(80n²log₂n)的运算复杂度在正反傅立叶变换中。

3.最后，由于要从过采样网格插值到临界采样(或低于临界采样)的网格，这将贡献O(4n²)的计算复杂度.

显而易见，随着输入和输出图像尺寸的增加，iDrizzle机制中的正反傅立叶变换和过采样插值消耗了越来越多的计算资源。如果要使用一个比临界采样率高2×2倍分辨率的过采样网格，那么iDrizzle要比fiDrizzle多花3.9倍的计算时间。

4.1讨论和结论

通过删除iDrizzle的一些冗余步骤，并直接叠加原始观测图像到临界采样率网格，fiDrizzle最终以更快的速度重建了比先前方法具有更高保真度的图像。fiDrizzle对比iDrizzle的进步和优势主要体现在如下3个方面：

1.fiDrizzle更加有效。因为人为的平滑插值等因素被去掉了，原始图像中反映真实图像细节的高频信息得以保留到最终的输出。从视觉表观和定量方面的比较中我们也发现fiDrizzle能够很好地重建小尺度上的细节信息，是高保真的。

2.fiDrizzle节省了大量计算资源。因为省去了频域空间的操作，而且是直接输出到目标网格，fiDrizzle最终用较少的运算量完成了高质量的图像重建工作。

3.fiDrizzle因为节省了内存支出，并且去除了正反傅立叶变换，因此能使程序运行得更加稳定可靠。

在不久的将来很多新的望远镜设备要投入观测中来：比如NASA的Wide FieldInfrared Survey Telescope(WFIRST)，欧空局的Euclid，美国国家自然科学基机会资助的Large Synoptic Survey Telescope(LSST)和中国的天宫望远镜(the Chinese SpaceStation optical Telescope，CSST)。届时会产生海量的数据，如何快速有效地处理这些数据将成为当务之急。实际上这恰好给了fiDrizzle很大的用武之地。

天文领域：

1.对于跟踪姿态保持不够好的望远镜，可以使用多次曝光技术，缩短曝光时间，这样就会大大降低由于姿态不稳引起的图像轮廓模糊，然后再用fiDrizzle重建高分辨率观测图像。

2.无论是空间望远镜还是地面望远镜，都可以选择好的时段对准被观测源进行多次曝光，然后用fiDrizzle重建源的高分辨率图像。这就避免了不良天气条件或者空间望远镜被地球月球等遮挡的因素影响。

3.由于天文观测中，暗弱源往往占大多数，有了fiDrizzle就可以使用这些数目可观的以前只能舍弃的暗弱源来做更深空间的天文研究。

4.只要知道拍摄望远镜的硬件情况，位置等参数我们就可以在fiDrizzle的协助下对同一个源的来自不同望远镜不同历史时期的观测数据进行处理并生成高分辨率图像，从而达到充分利用历史数据的目的。

5.天文观测在有些情况下进行多次曝光是必须的，因为在有多个信号源而且亮度相差很大的情况下，为了增强弱源的信号强度同时又要保证强源不饱和，就必须采用多次曝光，这时候再用fiDrizzle快有有效地重建具有更高保真度的源的图像就很必要了，至少比当前最好的iDrizzle要快要好。

数字图像监控领域，fiDrizzle技术能够利用连续的多帧录像画面重建得到被监控事物的高清图像。

微观物理学，微生物，医学成像等领域可以在分辨率有限的情况下用多次曝光结合fiDrizzle得到病毒，细菌或者有机分子的高分辨率图像。

本技术的局限是：多次曝光时段被观测源不能有很明显的形状和亮度轮廓的改变，否则必须舍弃这次曝光。

最后，该技术为数字成像领域提供了这样一条途径：通过多次曝光和fiDrizzle技术一样可以实现更高的图像分辨率。而不必一味在硬件上追求高的物理像素分辨率，另外即使有高的物理像素分辨率但还要克服长时间曝光中的仪器姿态位置稳定性问题。

5具体实施方式

我们用C语言开发了基于fiDrizzle技术的程序软件，只要知道成像仪器在拍摄图像时候的一些参数(如位置，旋转，CCD形变等)就可以重建高分辨率的源图像，给出更多细节信息。

致谢：本发明受国家973项目(No.2015CB857003，No.2015CB857000，No.2013CB834900)，国家基金委项目(No.11333008，No.11233005，No.11273061)，江苏省杰出青年项目(No.BK20140050)和中科院宇宙学结构先导项目(No.XDB09010000)的支持，发明者在此对这些支持项目表示感谢。

6附图说明

本专利共有8幅附图，用来在视觉表观和定量方面对比Dirzzle，iDrizzle和fiDrizzle三种重建方法的不同效果。附图说明如下：

图1：我们用Lena小姐作为真实图像(左图)，通过重新采样生成了原始观测图像(右图)，其中右图比左图在分辨率上降低了

倍。所以基本上相当于在左图中2×2像素的细节会被右图中一个像素代替。

图2：真实图像与三种不同方法重建结果的对比。左上是真实图像，右上是Drizzle重建结果，左下是iDrizzle的结果，右下是fiDrizzle的重建结果。三种方法都使用了5幅原始观测图像，其中iDrizzle和fiDrizzle都进行了5次迭代。三种重建结果都具有比原始观测图像高2×2的分辨率，也就是说跟真实图像分辨率相同，便于直接比较。

图3：三种重建方法的残差图像。也就是图2中左上的真实图像与其余的重建图像分别相减得到的差。其中左上还是真实图像，其他三幅布局跟图2是一样的。对比阴影中Lena肖像的残差，很明显fiDrizzle跟真实图像的差异是最小的。

图4：另一个图像重建测试图。来自哈勃望远镜的数据库。图像布局与图2是一致的。

图5：沿着X＝Y的像素上的流量轮廓和残差图。上图实线代表真实图像，点划线代表Drizzle重建的结果，划线代表iDrizzle的结果，点线则是fiDrizzle的重建流量。下图代表流量残差。注意：这里三个像素的长度与一个原始像素变长是一样的，因为二者分辨率差3×3倍。

图6：与流量图5类似，只是椭率分量e₊作为权重半径的函数r_ω。误差棒来源于图像噪音的估计。为了简化比较，在残差图中省略了误差棒。

图7：与图6类似，只是另一个椭率分量e_×的函数图像。

图8：与图6类似，只是方位角P.A.的函数图像。

Claims (1)

1.一种快速有效的像素化反卷积技术，其特征在于包括以下步骤：

步骤(1)：获取原始观测图像

将对准某区域观测获取的N幅震颤过的原始图像标记为{g¹(m，n)，g²(m，n)...g^N(m，n)}，用g^k(m，n)表示其中的第k幅，N，m，n和k为正整数，1≤k≤N，(m，n)为原始图像像素位置；

步骤(2)：重建图像

将所述原始观测图像{g¹(m，n)，g²(m，n)，g³(m，n)…g^N(m，n)}，用Drizzle方法叠加到一个预定的高分辨率标准网格上，得到重建图像f₁(x，y)，也即升采样操作，f的下标1代表这是第1次重建图像的结果，(x，y)为高分辨率的像素位置，x，y为正整数；

步骤(3)：迭代处理

步骤(3.1)：对第i次图像重建所得的高分辨率图像f_i(x，y)进行N次与步骤(1)相同采样条件的重新采样，也即降采样操作，从而获得N幅图像

用

表示N幅重新采样图像中的第k幅，其中i为自然数，f₁(x，y)也就是步骤(2)获得的重建图像；

所述重新采样图像与原始观测图像g^k(m，n)一一对应，区别在于前者是对高分辨率重建图像f_i(x，y)的重新采样，后者是真实的观测图像，之后进入下一步；

步骤(3.2)：用N幅原始图像g^k(m，n)分别减去对应的重新采样图像

得到N幅残差图像

其中进入下一步；

步骤(3.3)：将N幅残差图像

用Drizzle方法叠加到上次重建图像f_i(x，y)的标准网格上，从而得到新的重建图像f_i+1(x，y)，若f_i+1(x，y)不满足要求，则进入下一步；

步骤(3.4)：返回步骤(3.1)，进行下一次迭代循环操作，直至输出满足要求的图像。

Images