CN113945544A - 面向空间频域测量的动态单像素成像及退化图像复原方法 - Google Patents
面向空间频域测量的动态单像素成像及退化图像复原方法 Download PDFInfo
- Publication number
- CN113945544A CN113945544A CN202111193686.7A CN202111193686A CN113945544A CN 113945544 A CN113945544 A CN 113945544A CN 202111193686 A CN202111193686 A CN 202111193686A CN 113945544 A CN113945544 A CN 113945544A
- Authority
- CN
- China
- Prior art keywords
- measurement
- image
- sampling
- imaging
- pixel
- 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
Links
- 238000005259 measurement Methods 0.000 title claims abstract description 108
- 238000003384 imaging method Methods 0.000 title claims abstract description 95
- 238000000034 method Methods 0.000 title claims abstract description 50
- 238000005070 sampling Methods 0.000 claims abstract description 92
- 238000001228 spectrum Methods 0.000 claims abstract description 33
- 230000008569 process Effects 0.000 claims abstract description 22
- 230000010363 phase shift Effects 0.000 claims abstract description 17
- 238000001914 filtration Methods 0.000 claims abstract description 12
- 230000033001 locomotion Effects 0.000 claims description 40
- 230000015556 catabolic process Effects 0.000 claims description 26
- 238000006731 degradation reaction Methods 0.000 claims description 26
- 238000006073 displacement reaction Methods 0.000 claims description 20
- 230000008859 change Effects 0.000 claims description 4
- 238000011084 recovery Methods 0.000 claims description 4
- 230000003595 spectral effect Effects 0.000 claims description 3
- 230000021615 conjugation Effects 0.000 abstract 1
- 230000003287 optical effect Effects 0.000 description 8
- 238000010586 diagram Methods 0.000 description 5
- 230000009466 transformation Effects 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 4
- 238000003745 diagnosis Methods 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 102000001554 Hemoglobins Human genes 0.000 description 1
- 108010054147 Hemoglobins Proteins 0.000 description 1
- 208000026062 Tissue disease Diseases 0.000 description 1
- QVGXLLKOCUKJST-UHFFFAOYSA-N atomic oxygen Chemical compound [O] QVGXLLKOCUKJST-UHFFFAOYSA-N 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000000295 complement effect Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000005284 excitation Effects 0.000 description 1
- 238000005286 illumination Methods 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 229910044991 metal oxide Inorganic materials 0.000 description 1
- 150000004706 metal oxides Chemical class 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 229910052760 oxygen Inorganic materials 0.000 description 1
- 239000001301 oxygen Substances 0.000 description 1
- 239000004065 semiconductor Substances 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 230000003068 static effect Effects 0.000 description 1
- 230000002194 synthesizing effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
- G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
- G01N21/47—Scattering, i.e. diffuse reflection
- G01N21/4738—Diffuse reflection, e.g. also for testing fluids, fibrous materials
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/90—Dynamic range modification of images or parts thereof
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/20—Analysis of motion
- G06T7/262—Analysis of motion using transform domain methods, e.g. Fourier domain methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20056—Discrete and fast Fourier transform, [DFT, FFT]
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Immunology (AREA)
- Biochemistry (AREA)
- General Health & Medical Sciences (AREA)
- Analytical Chemistry (AREA)
- Chemical & Material Sciences (AREA)
- Pathology (AREA)
- Mathematical Physics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Multimedia (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Health & Medical Sciences (AREA)
- Image Analysis (AREA)
Abstract
本发明涉及一种面向空间频域测量的动态单像素成像及退化图像复原方法，包括如下步骤：采用三步相移傅里叶变换基函数作为单像素采样模板，分别以零频和共轭对称的两个调制频率为中心对傅里叶频谱进行采样，确定待测傅里叶系数，并生成对应的单像素采样模板；生成的单像素采样模板分别对应不同的采样空间频率，计算每个单像素采样模板与三个中心频率的最小欧式距离，根据最小欧式距离确定测量优先级，将采样模板上传到数字微镜器电路，在动态成像过程中按顺序循环加载采样模板进行连续测量；第一次测量完采样模板后基于三步相移公式计算傅里叶系数，计算逆滤波复原频谱，重建第一帧复原图像；实现动态成像。
Description
技术领域
本发明属于空间频域成像(Spatial Frequency Domain imaging,SFDI)技术领域，具体涉及到一种面向空间频域测量的动态单像素成像及退化图像复原方法。
背景技术
空间频域成像采用空间正弦调制光照射生物组织，并测量表面漫反射图像。由于漫反射光信号与组织体内部分子发生了复杂的相互作用，携带了组织体在特定光源激励下丰富的生理学信息，因而结合光学传输模型可以重构出生物组织的光学参数，进而可以得到与之关联的血红蛋白浓度、氧饱和度等生理参数。作为一种快速非接触式宽场成像技术，SFDI非常适合于皮肤等薄层组织疾病的诊断及治疗监测，近年来受到了广泛关注[1,2]。
传统SFDI系统一般采用基于电荷耦合器件或互补金属氧化物半导体的面阵传感器采集组织体表面的漫反射图像。面阵探测器在可见波长光谱中的应用中具有廉价、高性能的技术优势，但在非可见光波段，如光学诊断的近红外窗口，面阵探测器价格昂贵或难以实现。相比于面阵探测器，单像素探测器在非可见波段的灵敏度以及高精度的时间分辨率方面具有显著的性能优势，结合锁相光子计数技术，单像素探测器还可以实现面阵探测器难以实现的多波长并行测量，大幅度提升采集效率[3]。因此，基于单像素探测器测量的单像素成像(Single Pixel Imaging,SPI)技术逐渐被应用于空间频域成像领域[4]。
SPI使用一系列具有空间分辨的采样模板对待测物体的空间信息进行编码与压缩，采用单像素探测器记录压缩后的强度测量值，随后将测量结果与采样模板相结合，重建待测物体图像[5]。由于重建每帧图像需要多次测量，因此成像速度是SPI目前面临的主要问题之一。原则上，重建像素总数为N的图像需要M＝N个不同的采样模板。基于压缩感知(Compressive sensing,CS)原理，即大多数自然图像在适当的基表示中是稀疏的(有许多接近或等于零的系数)，可以用更少的测量实现图像重建。传统CS图像重建算法比较耗时，难以实现实时成像。为了提高图像重建速度，逐渐发展了基于正交变换的单像素成像方法，如hadamard变换、离散傅里叶变换和离散余弦变换等。针对空间频域漫反射图像为正弦谐波构成的特点，图像在余弦基中具有非常稀疏的表示，因此在SFDI中应用傅里叶单像素成像方法可以有效降低采样率[3]。
然而，目前单像素成像获取一张较高分辨率的图像，仍然需要数百次测量。同时，单个采样模板测量速度还受限于硬件系统(如数字微镜器件)的刷新频率，特别是傅里叶单像素成像成像的采样模板采用灰度图像，难以实现高速测量。此外，成像过程中目标体运动导致的图像退化会极大地降低成像质量。因此需发展有效的动态单像素空间频域成像方法，以及针对运动退化发展有效的图像复原方法，在保证成像质量的基础上进一步提高成像速度。
[参考文献]
[1]D.J.Cuccia,F.Bevilacqua,A.J.Durkin,F.Ayers,and B.J.Tromberg,“Quantitation and mapping of tissue optical properties using modulatedimaging,”J.Biomed.Opt.14(2),024012(2009).
[2]G.T.Kennedy,R.Stone II,A.C.Kowalczewski,R.Rowland,J.H.Chen,M.L.Baldado,A.Ponticorvo,N.Bernal,R.J.Christy,and A.J.Durkin,“Spatialfrequency domain imaging:aquantitative,noninvasive tool for in vivomonitoring of burn wound and skin graft healing,”J.Biomed.Opt.24(7),071615(2019).
[3]T.X.Li,Z.P.Qin,X.Hou,M.Dan,J.LI,L.M.Zhang,Z.X.Zhou and F.Gao,“Multi-wavelength spatial frequency domain diffuse optical tomography usingsingle-pixel imaging based on lock-in photon counting,”Opt.Express 27(16),23138(2019).
[4]M.Torabzadeh,I.Y.Park,R.A.Bartels,A.J.Durkin and B.J.Tromberg,“Compressed single pixel imaging in the spatial frequency domain,”J.Biomed.Opt.2017,22(3):30501
[5]J.Romberg,“Imaging via compressive sampling,”IEEE SignalProcess.Mag.25(2),14–20(2008).
发明内容
为了克服现有单像素成像速度的不足以及物体运动造成的图像质量退化，本发明提出一种面向空间频域测量的动态单像素成像及退化图像复原方法，采用交叠成像窗口更新测量值实现动态单像素成像，在无需改变测量系统硬件条件下，可以有效提高动态单像素成像帧速；基于傅里叶频移特性对单像素测量值进行逆滤波，复原物体运动造成的图像质量退化，消除运动伪影。技术方案如下：
一种面向空间频域测量的动态单像素成像及退化图像复原方法，包括如下步骤：
第一步，采用三步相移傅里叶变换基函数作为单像素采样模板，每个单像素采样模板下的测量值即为漫反射图像的傅里叶系数；基于漫反射图像频谱分布特点，分别以零频和共轭对称的两个调制频率为中心对傅里叶频谱进行采样，确定M个待测傅里叶系数，并生成对应的单像素采样模板；
第二步，生成的单像素采样模板分别对应不同的采样空间频率，中心频率为直流及两个共轭对称的交流调制频率，计算每个单像素采样模板与三个中心频率的最小欧式距离，根据最小欧式距离确定测量优先级，其中距离越小对应的优先级越高，每一个测量优先级中都同时包含表征直流和交流调制图像信息的采样模板，且越高测量优先级的采样模板表征的图像信息越多，基于此测量优先级按从高到低的顺序将采样模板上传到数字微镜器电路，在动态成像过程中按顺序循环加载采样模板进行连续测量，并实时传输测量数据；
第三步，第一次测量完M个采样模板后基于三步相移公式计算傅里叶系数，计算逆滤波复原频谱，重建第一帧复原图像，方法如下：
一帧图像测量时间内目标体离散为M个运动状态，以第一次测量时的运动状态为参考，设测量第i个傅里叶系数时目标体空间位移为(xi,yi)，则此时漫反射图像表示为
其中xi和yi分别为第i次测量时目标体沿x、y方向的运动位移，IDC(x-xi,y-yi)和IAC(x-xi,y-yi)分别为第i次测量时的直流调制图像和交流调制图像，fx为交流调制频率，
为相位,在动态单像素测量中，同一采样模板的第一帧测量值Yi与新一帧测量值Yi+M满足相移关系：
其中(x0,y0)为一帧图像测量过程中目标体沿x、y方向的运动位移，u、v分别表示测量值所对应的x、y方向的空间频率；利用相邻两个傅里叶系数前后两帧测量值Yi、Yi+M和Yi+1、Yi+1+M，根据公式(2)建立两个方程求解出物体位移参数(x0,y0)；在一帧成像时间内目标体设为匀速运动，线性计算出一帧图像测量时间内每个运动状态对应的位移(xi,yi)；
每一次单像素测量时对应直流调制图像的退化函数为：
每一次单像素测量时对应交流调制图像的退化函数为：
基于退化函数，分别对单像素测量值进行逆滤波获得复原图像频谱
其中Hi(u,v)表示Hi,DC(u,v)或Hi,AC(u,v)；逆滤波后的傅里叶频谱分别复原了直流成分和交流成分，进行傅里叶逆变换得到直流复原漫反射图像和交流复原漫反射图像，采用单次快照解调方法分别提取复原的直流幅度图像和交流幅度图像；
第四步，根据成像帧速需要确定相邻两帧成像窗口测量值交叠个数，按照测量顺序用新测量值替换上一帧中对应采样模板的测量值，得到新一帧图像的傅里叶频谱，在此过程中，处于同一测量优先级的测量值始终同步更新；利用新一帧成像窗口中最新测量的傅里叶系数更新目标体运动的位移参数，重复第三步实现图像复原；在连续测量下，重复第三步和第四步实现动态成像。
进一步地，通过改变相邻两帧成像窗口交叠长度可以调整成像帧速。
进一步地，每新测量一个傅里叶系数都替换前一帧对应的测量值，将成像的时间分辨率提高到单次测量所需的时间，获得连续变化过程中更多细节信息。
本发明的有益效果为：1.本发明通过交叠成像窗口更新采样模板测量值，可以获得传统单像素成像方法重建的相邻两帧图像的中间过程图像，有效提高了单像素成像动态测量的图像帧速；2.本发明在无需改变系统硬件测量速度下，通过改变成像窗口交叠长度可以调整成像帧速；3.本发明直接利用单像素测量值估算目标体运动参数，在连续动态测量中能够有效追踪目标体运动状态，并能够有效消除每帧重建图像的运动伪影；4.本发明可以直接应用于基于锁相光子技术的多波长单像素成像系统，实现多波长空间频域漫反射图像并行动态测量。
附图说明
图1是本发明提出的针对空间频域单像素成像的傅里叶频谱分布及测量优先级示意图
图2是本发明提出的动态单像素交叠成像窗口示意图
图3是本发明使用的基于锁相光子计数的多波长单像素空间频域成像系统结构示意图
具体实现方式
下面结合附图和具体实施例对本发明技术方案作进一步详细描述，所描述的具体实施例仅对本发明进行解释说明，并不用以限制本发明。具体实施方式中，使用附图3所示的基于锁相光子计数的多波长单像素空间频域成像系统实现测量。
通过设置合适的成像窗口交叠长度，可以实现对空间频域漫反射图像的动态测量与重建。由于漫反射图像频谱中同时包含直流和交流成分，针对此频谱特点，本发明采用傅里叶单像素成像方法，并基于一种椭圆形采样方式对直流及交流频率中心特定半径范围内的傅里叶系数进行测量，实现最大化采样效率。由于基于交叠成像窗口的动态单像素成像方法重建新一帧图像只需更新部分最新测量值，为保证每一帧重建的漫反射图像都同时反映漫反射图像中的直流和交流调制图像的变化信息(即直流和交流成分始终被同步采样)，本发明为所有单像素采样模板设计了测量优先级以保证待测漫反射图像中的直流和交流成分在动态成像过程中始终被同步采样。同时，针对动态成像过程中物体运动导致的图像退化，发展了一种基于傅里叶频移特性的图像复原算法，可以有效解决运动伪影。
本发明采用的技术方案如下。在动态单像素成像过程中，M个采样模板以设定的优先级顺序循环测量，所有测量值组成一个时间序列。其中任意M次连续测量都构成一个成像窗口，可以重建出一帧图像。当测量完第一个完整帧时，交叠移动成像窗口，更新采样模板测量值，可以快速重建出新一帧的图像。通过改变相邻两帧成像窗口的交叠长度，可以动态调整成像帧速。同时，利用每一帧成像窗口中新测量的傅里叶系数与被替换的傅里叶系数间的相位关系，计算目标体的运动位移。随后利用傅里叶频移特性计算测量的每个傅里叶系数对应的运动退化函数。结合退化函数采用逆滤波消除直流和交流调制图像运动伪影，进行傅里叶逆变换重建漫反射图像并解调得到直流和交流调制图像。包括如下步骤：
1.采用三步相移傅里叶变换基函数作为单像素采样模板，每个单像素采样模板下的测量值即为漫反射图像的傅里叶系数。为实现采样效率最大化，即使用最少的单像素采样模板获得足够多的图像信息，本发明针对漫反射图像傅里叶频谱的分布特点，以椭圆形采样方式确定待测傅里叶系数，并生成对应的单像素采样模板。本发明对采样方式进行了改进，常规方法假设待测图像只有零频一个中心频率，形成圆形采样，本发明针对待测图像的频谱特点(有三个中心频率)，形成椭圆形采样。具体说明如下。
由于空间频域成像采用正弦结构宽场照明，光在组织传播后，表面漫反射图像表示为：
其中IDC(x,y)和IAC(x,y)分别为直流调制图像和交流调制图像，fx为交流调制频率。因此空间频域漫反射图像同时包含直流(零频)和频率为fx的交流调制成分，其傅里叶频谱能量集中在零频(位于频谱中心)及调制频率fx(关于频谱中心左右对称)，如附图1左图所示。基于此频谱分布特点，在成像时分别以零频及调制频率fx为中心，对一定半径范围内的傅里叶系数进行测量，形成附图1右图所示的椭圆形采样方式。此椭圆形采样方式可以最大化对空间频域漫反射图像的采样效率。其中，椭圆形采样半径根据具体应用场景确定，结构越复杂的生理组织成像时所需的采样半径越大。随后，根据所选择的M个待测傅里叶系数，基于三步相移方法生成对应的单像素采样模板。
2.生成的单像素采样模板分别对应不同的采样空间频率。其中，对应低频的采样模板主要表征直流调制图像傅里叶系数，对应高频的采样模板主要表征交流调制图像的傅里叶系数。由于动态单像素过程中只更新部分测量值来重建新一帧的漫反射图像，为保证每一帧重建的漫反射图像中直流调制图像和交流调制图像的变换信息都被同步更新，需要为单像素采样模板设计合适的测量优先级。采用如下方法：在步骤1所述的椭圆形采样方式下，计算每个采样模板与三个中心频率(直流以及左右对称的两个交流调制频率)的最小欧式距离，根据最小欧式距离确定测量优先级，其中距离越小对应的优先级越高。在此设置下，每一个优先级中都同时包含表征直流和交流调制图像信息的采样模板，且越高优先级的采样模板表征的图像信息越多。基于此测量优先级按从高到低的顺序将采样模板上传到数字微镜器电路，在动态成像过程中按顺序循环加载采样模板进行连续测量，并实时传输测量数据。
3.第一次测量完M个采样模板后基于三步相移公式计算傅里叶系数，重建第一帧图像。但当测量过程中目标体运动时，单像素测量的是不同运动状态下目标图像的傅里叶系数，直接重建会存在运动伪影，因此需要先进行图像复原。图像复原的关键是计算测量过程中物体运动导致的退化函数。由于单像素成像过程中物体运动导致的图像退化机制与传统基于阵列相机的方式有很大不同，因此传统图像复原方法在单像素成像中无法应用。本发明针对单像素图像的运动退化机制，提出一种基于傅里叶频谱相移特性的图像复原方法，用以解决本发明中单像素图像的运动退化问题。具体实现如下：
由于单个傅里叶系数测量时间很短(<0.1s)，此时段内物体可视为静止状态，因此一帧图像测量时间内目标体可离散为M个运动状态。以第一次测量时的运动状态为参考，假设测量第i个傅里叶系数时目标体空间位移为(xi,yi)，则此时漫反射图像表示为
在动态单像素测量中，同一采样模板的第一帧测量值Yi与新一帧测量值Yi+M满足相移关系：
其中(x0,y0)为一帧图像测量过程中目标体运动位移。利用相邻两个傅里叶系数前后两帧测量值Y1、Y2和Y1+M、Y2+M，根据公式(3)可以建立方程组求解两个未知数(x0,y0)。在一帧成像时间内目标体假设为匀速运动，可线性计算出一帧图像测量时间内每个运动状态对应的位移(xi,yi)。
针对公式(2)中的直流调制图像，直接利用傅里叶频谱的相移特性，第i个运动状态直流调制图像的傅里叶系数与初始状态直流调制图像的傅里叶系数满足：
因此，每一次单像素测量时对应直流调制图像的退化函数为：
针对公式(2)中的交流调制图像，利用傅里叶频谱的相移特性和余弦调制性质，第i个运动状态交流调制图像的傅里叶系数与初始状态交流调制图像的傅里叶系数满足：
其中u0对应交流调制频率fx。由于交流调制图像频谱成分F(u-u0,v)和F(u+u0,v)分别以(-u0,0)和(u0,0)为中心，关于原点左右对称。对于大部分信息集中在低频的生物组织成像，F(u-u0,v)≈0,u＞0(F(u+u0,v)≈0,u＜0)。因此，每一次单像素测量时对应交流调制图像的退化函数可处理为：
基于公式(5)和(7)推导的退化函数，分别对单像素测量值进行逆滤波获得复原图像频谱：
其中Hi(u,v)表示Hi,DC(u,v)或Hi,AC(u,v)。逆滤波后的傅里叶频谱分别复原了直流成分和交流成分，进行傅里叶逆变换得到两幅复原图像。采用单次快照解调方法分别从两幅图像中提取出直流复原图像和交流复原图像。通过查表法可进一步重建目标体的光学特性图像。
4.在随后的成像中，根据成像帧速需要确定相邻两帧成像窗口测量值交叠个数，按照测量顺序用新测量值替换上一帧中对应采样模板的测量值，得到新一帧图像的傅里叶频谱。在此过程中，处于同一测量优先级的测量值始终同步更新。利用新一帧成像窗口中最新测量的傅里叶系数更新目标体运动的位移参数，重复步骤3实现图像复原。在连续测量下，重复步骤3、4实现动态成像。通过改变相邻两帧成像窗口交叠长度可以调制成像帧速。最快情况下，每新测量一个傅里叶系数都替换前一帧对应的测量值，将成像的时间分辨率提高到单次测量所需的时间，获得连续变化过程中更多细节信息。
下面结合实施例对本发明进行说明。
第一步，三波长方波调制光源通过DMD1向待测目标体投照空间频率为fx＝0.2mm-1的正弦空间调制光。确定成像区域大小为40×40mm2，单像素成像分辨率为128×128，在此条件下傅里叶频谱图像的中心频率分别为0(65,65)，-fx(57,65)和fx(73,65)，括号内为中心频率在频谱图像中的像素位置，fx对应的。基于此频谱分布特点，采用椭圆形采样方式确定待测傅里叶系数，采样半径为9.5，共确定M＝485个待测傅里叶系数。随后，基于三步相移方法生成对应的单像素采样模板Φi，i＝1，2…，M。
第二步，在上述椭圆形采样方式下，计算每个采样模板与三个中心频率的最小欧式距离，根据最小欧式距离确定测量优先级，其中距离越小对应的优先级越高。随后将生成的采样模板以测量优先级从高到低的顺序依次上载到DMD2控制电路。对于同一测量优先级的采样模板，以直流、交流正、交流负的顺序依次上载。控制DMD2依次循环加载采样模板，每个采样模板的测量时间为16ms(对应三步相移的总测量时间为48ms)。锁相光子计数模块同步工作，检测并分离各波长下每个采样模板的单像素测量数据，并实时传输各波长的测量数据。
第三步，第一次完整测量M个采样模板后，需要对测量数据进行逆滤波以消除物体运动导致的图像退化。首先利用相邻两个傅里叶系数前后两帧测量值，根据公式(3)建立方程组求解一帧图像测量过程中目标体运动位移(x0,y0)。对应第一帧图像使用的测量值是Y2、Y3和Y2+M、Y3+M(由于第一个采样模板对应的原点，代入公式(3)求解结果将为0，因此跳过第一个采样模板的测量值)，建立方程组求解。
此步骤中，综合使用三个波长的数据分别计算(x0,y0)并取均值，以降低测量噪声引起的误差。基于匀速运动假设线性计第一帧图像测量过程中第i个傅里叶系数Yi对应的位移(xi,yi)。根据公式(5)和(7)分别计算直流调制图像和交流调制图像的退化函数。其中直流调制图像对应的退化函数为：
交流调制图像对应的退化函数为：
其中，调制频率fx＝0.2mm-1对应的傅里叶离散采样频率u0＝0.625。
对于各波长测量数据，分别上述退化函数进行逆滤波。逆滤波后的傅里叶频谱分别复原了直流成分和交流成分，进行傅里叶逆变换可得到两幅复原图像。采用单次快照解调方法分别从两幅图像中提取出直流复原图像和交流复原图像，通过查表法可进一步重建目标体的光学特性图像。
第四步，在随后的成像中，以附图2中所示的成像窗口交叠更新顺序，即同一测量优先级测量值同时更新，替换前一帧中对应采样模板的测量值，得到新一帧图像的傅里叶频谱。同时利用新一帧成像窗口中最新测量的两个傅里叶系数更新位移参数(x0,y0)，计算新一帧图像每个傅里叶系数对应的直流和交流退化函数，对各波长的新一帧测量数据进行逆滤波，傅里叶逆变换并解调后得到直流和交流复原图像。采用查表法从复原的直流和交流漫反射图像中重建物体的光学特性图像。重复步骤3、4实现连续动态成像。在16ms的单个采样模板测量时间下，可并行实现每个波长约10帧/s的动态成像帧速。
尽管结合以上测量实例对本发明进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨的情况下，还可以做出很多变形，这些均属于本发明的保护之内。
Claims (3)
1.一种面向空间频域测量的动态单像素成像及退化图像复原方法，包括如下步骤：
第一步，采用三步相移傅里叶变换基函数作为单像素采样模板，每个单像素采样模板下的测量值即为漫反射图像的傅里叶系数；基于漫反射图像频谱分布特点，分别以零频和共轭对称的两个调制频率为中心对傅里叶频谱进行采样，确定M个待测傅里叶系数，并生成对应的单像素采样模板；
第二步，生成的单像素采样模板分别对应不同的采样空间频率，中心频率为直流及两个共轭对称的交流调制频率，计算每个单像素采样模板与三个中心频率的最小欧式距离，根据最小欧式距离确定测量优先级，其中距离越小对应的优先级越高，每一个测量优先级中都同时包含表征直流和交流调制图像信息的采样模板，且越高测量优先级的采样模板表征的图像信息越多，基于此测量优先级按从高到低的顺序将采样模板上传到数字微镜器电路，在动态成像过程中按顺序循环加载采样模板进行连续测量，并实时传输测量数据；
第三步，第一次测量完M个采样模板后基于三步相移公式计算傅里叶系数，计算逆滤波复原频谱，重建第一帧复原图像，方法如下：
一帧图像测量时间内目标体离散为M个运动状态，以第一次测量时的运动状态为参考，设测量第i个傅里叶系数时目标体空间位移为(xi,yi)，则此时漫反射图像表示为
其中xi和yi分别为第i次测量时目标体沿x、y方向的运动位移，IDC(x-xi,y-yi)和IAC(x-xi,y-yi)分别为第i次测量时的直流调制图像和交流调制图像，fx为交流调制频率，
为相位,在动态单像素测量中，同一采样模板的第一帧测量值Yi与新一帧测量值Yi+M满足相移关系：
其中(x0,y0)为一帧图像测量过程中目标体沿x、y方向的运动位移，u、v分别表示测量值所对应的x、y方向的空间频率；利用相邻两个傅里叶系数前后两帧测量值Yi、Yi+M和Yi+1、Yi+1+M，根据公式(2)建立两个方程求解出物体位移参数(x0,y0)；在一帧成像时间内目标体设为匀速运动，线性计算出一帧图像测量时间内每个运动状态对应的位移(xi,yi)；
每一次单像素测量时对应直流调制图像的退化函数为：
每一次单像素测量时对应交流调制图像的退化函数为：
基于退化函数，分别对单像素测量值进行逆滤波获得复原图像频谱
其中Hi(u,v)表示Hi,DC(u,v)或Hi,AC(u,v)；逆滤波后的傅里叶频谱分别复原了直流成分和交流成分，进行傅里叶逆变换得到直流复原漫反射图像和交流复原漫反射图像，采用单次快照解调方法分别提取复原的直流幅度图像和交流幅度图像；
第四步，根据成像帧速需要确定相邻两帧成像窗口测量值交叠个数，按照测量顺序用新测量值替换上一帧中对应采样模板的测量值，得到新一帧图像的傅里叶频谱，在此过程中，处于同一测量优先级的测量值始终同步更新；利用新一帧成像窗口中最新测量的傅里叶系数更新目标体运动的位移参数，重复第三步实现图像复原；在连续测量下，重复第三步和第四步实现动态成像。
2.根据权利要求1所述的动态单像素成像及退化图像复原方法，其特征在于，通过改变相邻两帧成像窗口交叠长度可以调整成像帧速。
3.根据权利要求2所述的动态单像素成像及退化图像复原方法，其特征在于，每新测量一个傅里叶系数都替换前一帧对应的测量值，将成像的时间分辨率提高到单次测量所需的时间，获得连续变化过程中更多细节信息。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111193686.7A CN113945544B (zh) | 2021-10-13 | 2021-10-13 | 面向空间频域测量的动态单像素成像及退化图像复原方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111193686.7A CN113945544B (zh) | 2021-10-13 | 2021-10-13 | 面向空间频域测量的动态单像素成像及退化图像复原方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113945544A true CN113945544A (zh) | 2022-01-18 |
| CN113945544B CN113945544B (zh) | 2023-10-17 |
Family
ID=79329573
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111193686.7A Active CN113945544B (zh) | 2021-10-13 | 2021-10-13 | 面向空间频域测量的动态单像素成像及退化图像复原方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113945544B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114859377A (zh) * | 2022-04-22 | 2022-08-05 | 中国科学院合肥物质科学研究院 | 一种运动目标实时捕获单像素成像方法及设备 |
Citations (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1326403A (zh) * | 1998-11-09 | 2001-12-12 | 西尔弗布鲁克研究有限公司 | 打印机及其操作方法 |
| EP2063260A1 (en) * | 2007-11-19 | 2009-05-27 | Lambda-X | Fourier transform deflectometry system and method |
| CN103440644A (zh) * | 2013-08-08 | 2013-12-11 | 中山大学 | 一种基于最小描述长度的多尺度图像弱边缘检测方法 |
| CN103617617A (zh) * | 2013-12-05 | 2014-03-05 | 淮海工学院 | 基于功率谱描述的水下图像质量评价测量方法 |
| WO2016015516A1 (zh) * | 2014-07-29 | 2016-02-04 | 暨南大学 | 一种使用单像素探测器的光学成像方法 |
| CN105841925A (zh) * | 2016-03-22 | 2016-08-10 | 中国科学院国家空间科学中心 | 一种基于探测器像素响应傅里叶谱获取的图像重建方法 |
| US10009587B1 (en) * | 2017-08-14 | 2018-06-26 | Christie Digital Systems Usa, Inc. | Real-time spatial-based resolution enhancement using shifted superposition |
| CN112465922A (zh) * | 2020-10-09 | 2021-03-09 | 浙江大学 | 一种结合色差模糊成像与图像重建技术的高光谱成像系统 |
| CN113114882A (zh) * | 2021-03-26 | 2021-07-13 | 暨南大学 | 一种高采样效率的傅里叶单像素成像方法 |
| CN113395501A (zh) * | 2021-06-15 | 2021-09-14 | 北京理工大学 | 一种基于互补性质的傅里叶单像素成像方法 |
-
2021
- 2021-10-13 CN CN202111193686.7A patent/CN113945544B/zh active Active
Patent Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1326403A (zh) * | 1998-11-09 | 2001-12-12 | 西尔弗布鲁克研究有限公司 | 打印机及其操作方法 |
| CN1508014A (zh) * | 1998-11-09 | 2004-06-30 | ��������³���о�����˾ | 数字打印系统 |
| EP2063260A1 (en) * | 2007-11-19 | 2009-05-27 | Lambda-X | Fourier transform deflectometry system and method |
| CN103440644A (zh) * | 2013-08-08 | 2013-12-11 | 中山大学 | 一种基于最小描述长度的多尺度图像弱边缘检测方法 |
| CN103617617A (zh) * | 2013-12-05 | 2014-03-05 | 淮海工学院 | 基于功率谱描述的水下图像质量评价测量方法 |
| WO2016015516A1 (zh) * | 2014-07-29 | 2016-02-04 | 暨南大学 | 一种使用单像素探测器的光学成像方法 |
| CN105841925A (zh) * | 2016-03-22 | 2016-08-10 | 中国科学院国家空间科学中心 | 一种基于探测器像素响应傅里叶谱获取的图像重建方法 |
| US10009587B1 (en) * | 2017-08-14 | 2018-06-26 | Christie Digital Systems Usa, Inc. | Real-time spatial-based resolution enhancement using shifted superposition |
| CN112465922A (zh) * | 2020-10-09 | 2021-03-09 | 浙江大学 | 一种结合色差模糊成像与图像重建技术的高光谱成像系统 |
| CN113114882A (zh) * | 2021-03-26 | 2021-07-13 | 暨南大学 | 一种高采样效率的傅里叶单像素成像方法 |
| CN113395501A (zh) * | 2021-06-15 | 2021-09-14 | 北京理工大学 | 一种基于互补性质的傅里叶单像素成像方法 |
Non-Patent Citations (3)
| Title |
|---|
| 张玉叶;李明珠;王春歆;李开端;: "运动模糊图像复原技术研究", 中国设备工程, no. 06 * |
| 毕祥丽;: "反射式压缩编码孔径超分辨率成像技术研究", 光电技术应用, no. 03 * |
| 赵宽心;李同心;侯茜;但迈;高峰;: "基于锁相光子计数测量的单像素空间频率域成像系统", 光学学报, no. 04 * |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114859377A (zh) * | 2022-04-22 | 2022-08-05 | 中国科学院合肥物质科学研究院 | 一种运动目标实时捕获单像素成像方法及设备 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113945544B (zh) | 2023-10-17 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110231292B (zh) | 一种单像素太赫兹波成像方法和系统 | |
| CN109124615B (zh) | 一种可选区高动态激光散斑血流成像装置及方法 | |
| US11368608B2 (en) | Compressed sensing based object imaging system and imaging method therefor | |
| CN106485764B (zh) | Mri图像的快速精确重建方法 | |
| CN109828371B (zh) | 一种基于移动散斑光源的大视场散射成像方法 | |
| US10438346B2 (en) | Method for extraction of spatial frequency information for quantitative tissue imaging | |
| CN113945544A (zh) | 面向空间频域测量的动态单像素成像及退化图像复原方法 | |
| CN111986118B (zh) | 一种加权核范数最小化的水下计算鬼成像图像去噪方法及系统 | |
| CN116051423B (zh) | 基于空间频域滤波的激光散斑衬比血流成像方法及系统 | |
| CN113709325B (zh) | 基于Hadamard频域变换矩阵阈值滤波的单像素成像方法 | |
| CN104207752A (zh) | 一种高速扫频光学相干断层成像系统 | |
| CN110533734A (zh) | 基于传统单能ct的多能谱分段稀疏扫描迭代重建方法 | |
| CN115690254A (zh) | 用于光相干层析成像的色散校正处理方法及电子设备 | |
| Wang et al. | Fast demodulation algorithm for multi-wavelength LED frequency-division modulation transmission hyperspectral imaging | |
| CN110044474B (zh) | 一种基于插值算法的关联成像优化方法 | |
| CN108230267B (zh) | 一种激光散斑血流成像图像处理方法 | |
| CN113395501A (zh) | 一种基于互补性质的傅里叶单像素成像方法 | |
| CN117968562A (zh) | 一种空间频域成像方法及系统 | |
| CN116309195B (zh) | 基于三维傅里叶变换的全场激光散斑血流成像方法及系统 | |
| CN115639198B (zh) | 一种全场光学时空相干编码动态体成像装置及方法 | |
| CN116659414B (zh) | 一种基于改进HiIbert变换的结构光解调方法 | |
| Ye et al. | A ghost imaging method based on multi-frequency fusion | |
| Hualong | A Visible and Infrared Image Fusion Method Based on Ghost Imaging | |
| Dan et al. | A lock-in photon-counting based single pixel imaging toward real-time multi-wavelength SFD-DOT | |
| CN114710615B (zh) | 高效单像素成像方法及系统 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |