CN114998161A - 一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 - Google Patents
一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 Download PDFInfo
- Publication number
- CN114998161A CN114998161A CN202210624880.4A CN202210624880A CN114998161A CN 114998161 A CN114998161 A CN 114998161A CN 202210624880 A CN202210624880 A CN 202210624880A CN 114998161 A CN114998161 A CN 114998161A
- Authority
- CN
- China
- Prior art keywords
- image
- precision
- pft
- fourier transform
- reconstruction
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 45
- 238000000386 microscopy Methods 0.000 title claims abstract description 25
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 26
- 238000004364 calculation method Methods 0.000 claims description 16
- 238000009499 grossing Methods 0.000 claims description 12
- 210000001747 pupil Anatomy 0.000 claims description 11
- 230000004075 alteration Effects 0.000 claims description 5
- 239000000654 additive Substances 0.000 claims description 3
- 230000000996 additive effect Effects 0.000 claims description 3
- 230000008569 process Effects 0.000 claims description 2
- 238000003325 tomography Methods 0.000 claims 3
- 238000003384 imaging method Methods 0.000 abstract description 22
- 238000005259 measurement Methods 0.000 abstract description 7
- 230000006870 function Effects 0.000 description 21
- 230000000694 effects Effects 0.000 description 18
- 238000001228 spectrum Methods 0.000 description 16
- 238000005516 engineering process Methods 0.000 description 9
- 230000000737 periodic effect Effects 0.000 description 6
- 238000011084 recovery Methods 0.000 description 6
- 238000002474 experimental method Methods 0.000 description 5
- 239000000523 sample Substances 0.000 description 5
- 238000004088 simulation Methods 0.000 description 5
- 238000012937 correction Methods 0.000 description 3
- 238000010586 diagram Methods 0.000 description 3
- 238000012545 processing Methods 0.000 description 3
- 230000008859 change Effects 0.000 description 2
- 238000009795 derivation Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 230000006978 adaptation Effects 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 239000011324 bead Substances 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 239000012472 biological sample Substances 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000000903 blocking effect Effects 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000004134 energy conservation Methods 0.000 description 1
- 230000004907 flux Effects 0.000 description 1
- 238000005286 illumination Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000002906 microbiologic effect Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/10—Image enhancement or restoration using non-spatial domain filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/50—Image enhancement or restoration using two or more images, e.g. averaging or subtraction
-
- 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]
-
- 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/20172—Image enhancement details
- G06T2207/20192—Edge enhancement; Edge preservation
-
- 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/20212—Image combination
- G06T2207/20224—Image subtraction
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Analysis (AREA)
Abstract
本发明公开了一种基于完美傅里叶变换算法(PFT)的傅里叶叠层显微术高精度重构图像方法，该方法包括如下步骤：1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像；2)计算得到平滑的误差图像；3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像。本发明的方法在提高重构图像精度的同时，还保证了与快速傅里叶变换基本一致的效率，还提出了改进的广义傅里叶叠层显微成像模型。实测中，该方法将相位精度的标准差从0.08弧度提高到0.02弧度，对应400分之一波长，精度提高了约4倍。
Description
技术领域
本发明属于光学信息获取与处理技术领域，具体而言涉及傅里叶叠层显微成像技术的高精度相位成像，提供了一种基于PFT的傅里叶叠层显微术高精度重构图像方法，同时也为傅里叶叠层显微成像技术提供了一种新的广义模型，这将进一步为相关定量生物学和测量等应用提供参考。
背景技术
傅里叶叠层显微术(Fourier ptychographic microscopy,FPM)是一种极具前景的计算成像技术，具有高分辨率、大视场和定量相位恢复等特点。其通过获取多角度照明信息(低分辨率图像集)，然后在重构算法中利用这些低分辨率图像集恢复出高分辨率图像。当前，FPM技术已成功应用于定量相位成像[Sun J,Chen Q,Zhang J,et al.Single-shotquantitative phase microscopy based on color-multiplexed Fourier ptychography[J].Optics Letters,2018,43(14):3365-3368.]，且已产生了很多高精度的系统标定方法，如像差去除[Ou X,Zheng G,Yang C.Embedded pupil function recovery forFourier ptychographic microscopy.Opt Express 2014；22(5):4960-4972.]、LED强度波动校正[Bian Z,Dong S,Zheng G.Adaptive system correction for robust Fourierptychographic imaging.Opt Express 2013；21(26):32400–32410.]、LED位置校正[PanA,Zhang Y,Zhao T,Wang Z,Dan D,Lei M,Yao B.System calibration method forFourier ptychographic microscopy.J Biomed Opt 2017；22(9):096005.]、渐晕效应去除[Pan A,Zuo C,Xie Y,Lei M,Yao B.Vignetting effect in Fourier ptychographicmicroscopy.Opt Laser Eng 2019；120:40-48.]、噪声抑制[Zhang Y,Pan A,Lei M,YaoB.Data preprocessing methods for robust Fourier ptychographic microscopy.OptEng 2017；56(12):123107.Yeh L-H,Dong J,Zhong J,Tian L,Chen M,Tang G,Soltanolkotabi M,Waller L.Experimental robustness of Fourier ptychographyphase retrieval algorithms.Opt Express 2015；23(26):33214–33240.]等等。在FPM重构算法中由于快速傅里叶变换(fast Fourier transform,FFT)需要对图像进行周期性扩展[James W.Cooley,John W.Tukey.An algorithm for the machine calculation ofcomplex Fourier series[J].Mathematics of Computation,1965,19(90)：297-301.]，但是非周期性的图像进行FFT时会因为边缘效应产生伪影现象，这种伪影现象对于纯相位样本而言尤为明显，极有可能降低成像精度。该伪影在傅里叶空间中呈现为一种从高频到低频的交叉形状，其中高频伪影对应于空域强度图像视场边缘处较为明显的波动，类似于吉布斯效应或者振铃效应，而低频伪影在重构相位图像上呈现为一种“波纹”状。此外，由于不同的分块操作产生了不同的边缘条件，同一子区域的成像精度取决于分块的大小，这些问题都限制了FPM定量测量的精度提升。
对于边缘效应，直观想，补零插值可能是一种可采取的简单方法，即通过在空域对图像的边缘补零，使图像满足四个边界的连续性，但是采用零填充插值后的图像仍然是非周期的，这无法从根本上解决非周期性图像进行FFT时产生的边缘效应。在强度传输方程等技术中也涉及到边缘效应问题，其通过离散余弦变换(discrete cosine transform,DCT)来解决[Zuo C,Chen Q,Asundi A.Boundary-artifact-free phase retrieval with thetransport of intensity equation:fast solution with use of discrete cosinetransform.Opt Express 2014；22(8):9220-44.；Huang L,Zuo C,Idir M,Qu W,AsundiA.Phase retrieval with the transport-of-intensity equation in anarbitrarilyshaped aperture by iterative discrete cosine transforms.Opt Lett 2015；40(9):1976-1979.]。DCT主要通过对原始图像“添加信息”，使原始图像具有周期性，即先对原始图像进行四倍放大，再通过裁剪其中四分之一的方法来去除FFT后的重复信息。DCT方法可以有效去除伪影，但是也面临着计算时间和内存增加的问题。因此，FPM定量应用方面需要高精度、高效率的边缘效应去除方法。
名词解释
1.傅里叶叠层显微成像技术，英文准确名称为“Fourierptychographicmicroscopy”，简称FPM技术，由美国加州理工学院Yang等人发明于2013年[Zheng G,Horstmeyer R,Yang C.Wide-field,high-resolution Fourierptychographicmicroscopy[J].Nature Photonics,2013,7(9):739-745.]，中文翻译存在多种名称如傅立叶重叠关联成像、傅里叶叠层显微成像技术、傅里叶叠层成像术、傅里叶叠层技术等，本文统一采用傅里叶叠层成像技术这一专业术语。
2.快速傅里叶变换(fast Fourier transform,FFT)[James W.Cooley,JohnW.Tukey.An algorithm for the machine calculation of complex Fourier series[J].Mathematics of Computation,1965,19(90):297-301.]是数字信号处理邻域中一种利用计算机计算离散傅里叶变换的高效、快速计算方法的统称，简称FFT。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少。
3.离散余弦变换(discrete cosine transform,DCT)，是与傅里叶变换相关的一种变换，它类似于离散傅里叶变换，但是只使用实数。
发明内容
本发明的目的是提供一种基于PFT的傅里叶叠层显微术高精度重构图像方法，用于解决现有技术中去除边缘效应算法上对精度和效率的需求。
为了实现上述技术任务，本发明采用如下技术方案予以实现：
一种基于PFT的傅里叶叠层显微术高精度重构图像方法，该方法包括如下步骤：
1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像；
2)计算得到平滑的误差图像；
3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像。
进一步地，所述的1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像，即：
f(x,y)＝g(x,y)+e(x,y)
其中，x∈[0,M-1]and y∈[0,N-1]，f(x,y)表示FPM原始重构结果，其作为初始图像进行后续计算；g(x,y)表示最终所需的高精度分量；e(x,y)表示平滑误差分量。
所述的2)计算得到平滑的误差图像，具体过程如下：
步骤2.1：定义像素集合；
步骤2.2：建立关于高精度分量g(x,y)和平滑误差分量e(x,y)的目标函数E；
步骤2.3：通过最小化步骤2.2中的目标函数E求得平滑误差分量e(x,y)。
其中步骤2.1：定义像素集合，具体为：
定义G(a)为像素a在原始图像所在空间Q中的相邻像素集，符号b表示G(a)的元素，定义Ψ为原始图像FFT变换后的扩展空间，H(a)为Ψ\Q中a的相邻像素集，z表示H(a)的元素，z’定义为像素z映射到Q空间的位置，即：
G(a)＝{b∈Q,|a-b|＝1}
步骤2.2：建立关于高精度分量g(x,y)和平滑误差分量e(x,y)的目标函数E，即：
上述目标函数E约束条件为：
f(x,y)＝g(x,y)+e(x,y)
其中，g(a)-g(z’)为高精度分量g(x,y)边界处的像素计算；e(a)-e(b)为平滑误差分量e(x,y)空间域上的相邻像素差值。
所述的步骤2.3，通过最小化步骤2.2中的目标函数E求得平滑误差分量e(x,y)，具体包括如下步骤：
步骤2.3.1，根据步骤2.1中定义的像素集合，由集合的绝对值计算对应集合的基数；
步骤2.3.2，将最小化步骤2.2中的目标函数E转化为函数K和误差函数e的卷积运算，即：
L(f)＝K*e
其中，*表示卷积运算符，K表示Ψ空间中某个像素到像素a的距离c；L(f)是一个线性运算符，定义为：
步骤2.3.3，对L(f)＝K*e两边进行FFT计算，得到平滑误差分量e(x,y)的FFT形式为：
其中，u(x,y)相当于计算初始图像f(x,y)的上下边缘差值u1和左右边缘差值u2之和，即u＝u1+u2，且u1和u2如下式所示，其中
p∈[0,M-1]，q∈[0,N-1]。
所述的3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像，具体包括：
步骤3.1：利用PFT提取重构图像的傅里叶变换，得到高精度分量g(x,y)的FFT形式，即：
其中x∈[0,M-1]and y∈[0,N-1]，f(x,y)为初始图像，g(x,y)为预期的高精度图像，e(x,y)为平滑误差图像。
步骤3.2：对高精度分量g(x,y)的FFT形式进行逆傅里叶变换，并计算强度，即可得到高精度重构图像。
本发明同时还给出，所述的3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像，其高精度重构图像的原始图像计算模型为：
其中，r＝(x,y)为空间域坐标，k＝(kx,ky)为傅里叶域坐标；P(k)为系统光瞳函数，O(k-ki)为光瞳面的出射波，i为图像的序号，Ii(r)即为原始图像集；权重因子wi对应LED强度波动误差，相位项ejΦ对应光瞳函数中的像差，In,i对应加性噪声，△ki对应位置偏差；F为PFT算法。
优选地，步骤2.3.1的集合的绝对值计算对应集合的基数，即|G(a)|+|H(a)|＝4。
步骤2.3.2中K表示Ψ空间中某个像素到像素a的距离c，则
本发明与现有技术相比，具有以下有益效果：
1)提出了基于PFT的傅里叶叠层显微术高精度重构图像方法，将标准差从0.08弧度提高到0.02弧度，对应400分之一波长，精度提高了4倍。
2)该方法在去除伪影提高成像精度的同时，能够保证效率与原始快速傅里叶变换相当。
3)为FPM高精度测量应用提供一个包含FPM现有系统误差(强度波动、像差、噪声、位置偏差等)的广义模型。
附图说明
图1是边缘效应说明图。其中，(a1-c1)为强度、相位和模拟中真实的傅里叶频谱(无十字伪影)；(a2-c2)为不采用常规带通滤波器重构的频谱；(a3-c3)为带有常规带通滤波器重构的频谱；(a4-c4)为(a2-c2)和(a3-c3)的差分图；(a5-c5)为对(a1)和(b1)中框选区域的重构。
图2是不同初始猜测下的重构结果图。其中，(a1-a4)为不同初始猜测的原始图像；(b1-b4)为初始猜测的傅里叶谱；(c1-c4)为重构的相位；(d1-d4)为重构的傅里叶谱。
图3是FFT、DCT、PFT的原理图。其中，(a)为FFT的隐式周期化假设；(b)为DCT水平、垂直方向均为对称运行；(c)为将恢复后的图像在PFT中分解为高精度图像和平滑误差图像。
图4是PFT推导中变量和集合的说明图。
图5是分别采用FFT、DCT和PFT进行FPM重构的结果图。其中，(a1-c1)真实标记；(a2-c2)为FFT重构结果；(a3-c3)为DCT重构结果；(a4-c4)为PFT重构结果。
图6是美国空军分辨率板的重构结果。其中，(a,a1)为全视场及其特写；(b)为采用常规FFT重构结果；(c)为PFT重构结果；(d)为FFT和PFT方法的差异；(e)为相位精度线廓线。
图7是采用FFT和PFT对Hela细胞的重构实验。其中，(a,a1)为全视场及其特写；(b)为采用常规FFT重构；(c)为PFT重构；(d)为FFT和PFT重构结果的差异。
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
本发明提出了基于PFT的傅里叶叠层显微术高精度重构图像方法。其采用“减法操作”的思路，即在FFT重构完后的图像中减去伪影图像。具体是将FPM重构后非周期的原始图像分成两部分，一部分是满足FFT周期要求的重构图像，也即无伪影的结果，另一部分是一幅平滑的误差图像，也即伪影图像。因为PFT是建立在FFT操作后的一步操作，在去除伪影的同时，能够保证效率基本与原始FFT相当。
FPM中典型的低分辨率图像(low resolution,LR)图像计算模型为
其中r＝(x,y)为样本空间域的坐标，k＝(kx,ky)为傅里叶域坐标；P(k)为系统光瞳函数，O(k-ki)为光瞳面的出射波，i为LR图像的序号；{·}表示二维傅里叶变换，可以定义为
其中g(·)是虚函数，j是虚数单位。
FPM重构算法中使用FFT可以快速计算图像的离散傅里叶变换，但是在FFT中需要对图像进行隐式周期化假设，对于非周期图像，该假设可能会导致傅里叶频谱出现交叉形状的伪影，这种现象称为边缘效应。如图1所示，我们首先对边缘效应的产生原因进行分析，图1(a1-c1)展示了真实的强度、相位及其傅里叶谱，注意其中傅里叶谱中不存在交叉伪影。我们将正入射LR图像的上采样作为初始猜测进行FPM重构后，傅里叶谱(图1(c2))出现了明显的交叉伪影，且对应的强度(图1(a2))和相位(图1(b2))均明显存在精度误差(圆圈和箭头强调区域的对比)。一般，带通滤波器可用来滤除高频误差，但是其使用后的结果(图1(a3-c3))与真实图像相比，仍然存在精度误差，且差分图的强度差(图1(a4))中的差异非常小，视觉上变化不大。对比，图1(c2)与图1(c3)可以发现在傅里叶频谱中消除了高频伪影，但低频伪影始终存在(图1(c3)中的箭头)，而正是低频伪影导致了在强度和相位分布上出现了波纹状的伪影(图1(a3,b3))。如果采用均方根误差(root-mean-square error,RMSE)来客观评价FPM的重构效果，即
其中f(x,y)，g(x,y)是两幅图像，M和N表示图像的大小。计算得到的使用带通滤波器前后恢复的强度RMSE分别为4.36％和4.34％，相位RMSE分别为2.47％和2.47％，即带通滤波器对强度图像的精度有轻微提升效果，但对相位精度几乎没有提升，这说明了相位图像对边缘效应导致的伪影更加敏感。此外，通过对比图1(a3,b3)和选取图1(a1,b2)中子区域恢复的结果图(a5,b5)，可以发现图1(b5)中存在很明显的由高频伪影(图1(c5)中的箭头)引起的边缘波动，虽然在图1(c3)傅里叶谱中也存在类似的伪影，但不是那么明显，因为这些伪影会被傅里叶空间中样本频谱所覆盖，因而与孔径合成的大小有关。同时，成像精度也取决于不同的分块处理方式，因为同一子区域的成像精度与不同边缘条件下的分块大小有关。
比较不同初始猜测和对应重构的傅里叶谱，如图2所示，发现双线性猜测(图2(d1)中箭头)和双三次猜测(图2(d2)中箭头)由于上采样出现了明显的高频伪影，而在所有重构的傅里叶谱(图2(d1-d4))中都存在低频伪影。其中，对双线性猜测、双三次猜测、全一猜测、随机猜测的重构强度图像计算得的RMSE分别为4.36％、4.35％、4.55％、30.75％，重构后的相位图像RMSE分别为2.47％、2.47％、2.40％、4.64％，显然随机猜测与其他三种初始猜测方法相比，更难收敛(图2(c4))。因此，结合图1，我们发现精度误差(图2(c1-c3))主要是由低频伪影造成，且与初始猜测和上采样采用何种方式无关。由于除了随机猜测外，不同的初始猜测会得到相似的结果，我们在接下来的仿真和实验中，使用全一猜测来缓解这些高频伪影。
与强度传输方程技术中使用DCT去除边缘效应不同，我们提出一种基于基于PFT的傅里叶叠层显微术高精度重构图像方法，与FFT、DCT原理对比如图3所示，首先在实施FFT时(图3(a))，数字图像需要周期性地放大，而非周期图像被不精确地视为周期图像，从而由于边缘效应导致了交叉形状的伪影。而DCT方法可以看作是沿着水平轴和垂直轴依次进行对称运算，得到一幅大小为原图像4倍的新图像(图3(b))，然后对新得到的具有空间周期性的对称图像进行FFT快速计算，从而消除边缘效应。但是最终的图像需要对称图像进行四分之一的切割操作，因此DCT需要更多计算时间和内存。而本发明提出的基于PFT的傅里叶叠层显微术高精度重构图像方法(图3(c))是直接去除FPM重构中由于FFT造成的误差，即首先将M×N空间的重构图像分解为两部分(图3(c))，其一为目标高精度结果，其二为交叉伪影对应的平滑，然后在重构时忽略平滑误差图像，即可得到高精度重构结果。
具体地，执行以下操作：
1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像；
2)计算得到平滑的误差图像；
3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像。
所述的步骤1)1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像，即：
f(x,y)＝g(x,y)+e(x,y) (4)
其中，x∈[0,M-1]and y∈[0,N-1]，f(x,y)表示FPM原始重构结果，其作为初始图像进行后续计算；g(x,y)表示最终所需的高精度分量；e(x,y)表示平滑误差分量。
为了计算平滑误差图像e(x,y)：步骤2.1，定义像素集合，Q设置为初始离散图像f(x,y)的空间域，a表示在Q空间的某一个像素，如图4中所示，Ψ空间是M×N的Q空间FFT变换时的扩张，Ψ\Q表示Ψ空间除了Q的部分，定义G(a)为像素a在Q空间中的相邻像素集，符号b表示G(a)的元素:
G(a)＝{b∈Q,|a-b|＝1} (4)
将H(a)定义为Ψ\Q中a的相邻像素集，z表示H(a)的元素，由:
z’定义为像素z映射到Q空间的位置。
将初始图像分解为g(x,y)表示最终所需的高精度分量和e(x,y)表示平滑误差分量两个分量，应满足四个基本规则：1)高精度图像分量g(x,y)的边界应尽可能平滑，以满足隐式周期化假设；2)平滑误差分量e(x,y)对应于低频伪影，应为平滑图像；3)平滑误差分量e(x,y)的平均强度应为零，高精度图像分量g(x,y)的强度应接近f(x,y)，以满足能量守恒。
其次，步骤2.2，建立关于g(x,y)和e(x,y)这两个分量的目标函数E：
上述目标函数E约束条件为：
f(x,y)＝g(x,y)+e(x,y)
步骤2.3，通过最小化目标函数E求得平滑误差分量e(x,y)。
步骤2.3.1，将步骤2.2中约束条件(8)带入目标函数(9)，得到：
使E对e(a)求导，并令导数应为零，推导得：
进一步推导出
步骤2.3.2，根据步骤2.1中定义的像素集合，由集合的绝对值计算对应集合的基数：
|G(a)|+|H(a)|＝4 (12)
步骤2.3.3，将最小化步骤2.2.1中的目标函数E的推导式(11)转化为函数K和误差函数e的卷积运算，即：
L(f)＝K*e (13)
其中，*表示卷积运算符，K表示Ψ空间中某个像素到像素a的距离c，
L(f)是一个线性运算符，定义为：
步骤2.3.4，对L(f)＝K*e两边进行FFT计算，得到平滑误差分量e(x,y)的FFT形式为：
其中，u(x,y)相当于计算初始图像f(x,y)的上下边缘差值u1和左右边缘差值u2之和，即u＝u1+u2，且u1和u2如下式所示，其中
p∈[0,M-1]，q∈[0,N-1]。
参见图4，u(x,y)为初始图像的四个边缘。除了这四条边，所有u(x,y)的值都是0。u(x,y)中每个边缘的值等于它在f(x,y)中的对边减去它在f(x,y)中的对应边。上下边和左右边的值分别定义为u1和u2，因此有u＝u1+u2。
进一步地，所述的步骤3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像，具体包括：
步骤3.1：利用PFT提取重构图像的傅里叶变换，得到高精度分量g(x,y)的FFT形式，即：
步骤3.2：对高精度分量g(x,y)的FFT形式进行逆傅里叶变换，并计算强度，即可得到高精度重构图像。
与此同时，利用本发明的方法得到高精度重构图像，所述的高精度重构图像的原始图像计算模型为：
其中，r＝(x,y)为空间域坐标，k＝(kx,ky)为傅里叶域坐标；P(k)为系统光瞳函数，O(k-ki)为光瞳面的出射波，i为LR图像的序号，Ii(r)即为原始图像集；权重因子wi对应LED强度波动误差，相位项ejΦ对应光瞳函数中的像差，In,i对应加性噪声，△ki对应位置偏差；F为PFT算法。
该原始图像计算模型，我们称之为广义傅里叶叠层显微成像模型，它结合FPM中其他已知的系统误差，可以作为FPM高精度测量的广义模型，进一步为相关定量生物学和测量等应用提供参考。
具体实施例
遵循本发明的整体技术方案，为验证上述模型的可行性与正确性，本发明选用基于高斯牛顿算法的FPM成像系统作为实施例，并对FFT、DCT和PFT分别在仿真组和实验组进行验证。
仿真方面，采用模拟条件：11×11LED阵列(4毫米间距，中心波长630nm)、LED阵列与样品间距为76mm、4×/0.1NA物镜、6.5μm像素尺寸的相机、采集128×128像素大小的图像块；所有结果运行30次迭代以确保收敛。重构的结果如图5所示，图中精度误差用圆圈和箭头突出显示。计算得到的FFT、DCT和PFT重构强度图像的RMSE分别为4.55％、4.49％和0.64％，FFT重构相位图像的RMSE分别为2.4％、1.1％和1.1％。图5(c3,c4)中可以看到，DCT和PFT都可以去除低频伪影，提高相位精度。但是FFT和DCT强度结果相比真实强度图5(a1)都有小的下降，根据RMSE，DCT相比于FFT的精度并没有提高(图5(a3)中的圆圈)；且在DCT产生的傅里叶谱中有一个不明显的黑色实框(图5(c3))。然而PFT的强度和相位图像都非常接近真实值(图5(a4-c4))。同时，PFT的时间代价为5.52s，与FFT的5.44s相当，而DCT的时间代价为12.45s，是PFT和FFT的两倍多。注意，通常计算时间与像素的数量不是线性的。因此，虽然DCT可以在一定程度上降低交叉伪影造成的相位图像和傅里叶谱的精度误差，但它重构的强度图像有一个较小的精度误差，且时间成本更高，而PFT相对而言是一个实现无边缘效应傅里叶变换的快速方法。
实验方面，采用实验条件：32×32的可编程阵列LED阵列(4毫米间距，中心波长标定为631nm，由Arduino控制点亮)，成像使用LED阵列中间的15×15的灯珠、一个4×/0.1NA消色差物镜、一个8位CCD相机(型号为DMK23G445，由德国Imaging Source公司制造，像素为1280×960，像元尺寸为3.75μm)；作为样品的分辨率板放置在距离LED阵列86mm的高度；FPM重构算法采用经典高斯方法，重构算法是在CPU为I7-10700和Matlab2016a的计算机平台上进行。FFT和PFT得到的重构结果如图6所示。分辨率板的整个视场如图6(a)所示，其近距离视图如图6(a1)所示。由FFT与PFT结果计算得到的差分图如图6(d)所示，从相位差可以明显看出本发明提出的PFT取得了很好的效果。如图6(e)所示，我们还提供了定量测量的相位精度曲线，在0附近的虚线表示真值，可以看到常规FFT方法的相位精度标准差为0.08弧度，而PFT算法的相位精度标准差为0.02弧度，即相比较FFT，PFT算法将相位精度标准差从0.08弧度提高到0.02弧度，提高了4倍。此外，在计算时间上，PFT的时间代价为12.33s，与FFT的11.92s相当。
为了进一步验证基于PFT的傅里叶叠层显微术高精度重构图像方法即PFT算法的有效性，我们继续使用未染色的Hela细胞切片作为生物样本，其可以视为纯相位样本。采用20×/0.4NA摄像机进行拍摄，以测试不同目标下的PFT性能；考虑到这个物镜的光通量较小，我们将高度改为55mm以获得足够的亮度。实验结果如图7所示，整个视场如图7(a)所示，其覆盖了单个Hela细胞的特写图如图7(a1)所示。图7(d)所示的FFT和PFT结果的差分图证明了对于纯相位样品，本发明提出的PFT方法效果更加明显，对相位更具优势。同样地，在计算时间上，PFT的时间代价为8.51s，与FFT的7.97s相当。
本方法针对原始的高斯牛顿算法进行FPM重构，成功利用PFT解决了边缘效应引起的伪影，有效提升FPM的成像精度，尤其是相位恢复精度。同理，本方法中的提出的LR的广义计算模型适用于所有FPM定量应用技术，提升成像精度。
Claims (10)
1.一种基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于，该方法包括如下步骤：
1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像；
2)计算得到平滑的误差图像；
3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像。
2.如权利要求1所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的1)将通过FPM重构得到的FFT结果分解为高精度重构图像和平滑的误差图像，即：
f(x,y)＝g(x,y)+e(x,y)
其中，x∈[0,M-1]and y∈[0,N-1]，f(x,y)表示FPM原始重构结果，其作为初始图像进行后续计算；g(x,y)表示最终所需的高精度分量；e(x,y)表示平滑误差分量。
3.如权利要求1所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的2)计算得到平滑的误差图像，具体过程如下：
步骤2.1：定义像素集合；
步骤2.2：建立关于高精度分量g(x,y)和平滑误差分量e(x,y)的目标函数E；
步骤2.3：通过最小化步骤2.2中的目标函数E求得平滑误差分量e(x,y)。
4.如权利要求2所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的步骤2.1：定义像素集合，具体为：
定义G(a)为像素a在Q空间中的相邻像素集，符号b表示G(a)的元素，定义H(a)为Ψ\Q中a的相邻像素集，z表示H(a)的元素，z’定义为像素z映射到Q空间的位置，即：
G(a)＝{b∈Q,|a-b|＝1}
5.如权利要求2所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的步骤2.2：建立关于高精度分量g(x,y)和平滑误差分量e(x,y)的目标函数E，即：
上述目标函数E约束条件为：
f(x,y)＝g(x,y)+e(x,y)
其中，g(a)-g(z’)为高精度分量g(x,y)边界处的像素；e(a)-e(b)为平滑误差分量e(x,y)空间域上的相邻像素差值。
6.如权利要求2或3所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的步骤2.3，通过最小化步骤2.2中的目标函数E求得平滑误差分量e(x,y)，具体包括如下步骤：
步骤2.3.1，根据步骤2.1中定义的像素集合，由集合的绝对值计算对应集合的基数；
步骤2.3.2，将最小化步骤2.2中的目标函数E转化为函数K和误差函数e的卷积运算，即：
L(f)＝K*e
其中，*表示卷积运算符，K表示Ψ空间中某个像素到像素a的距离c；L(f)是一个线性运算符，定义为：
步骤2.3.3，对L(f)＝K*e两边进行FFT计算，得到平滑误差分量e(x,y)的FFT形式为：
其中，u(x,y)相当于计算初始图像f(x,y)的上下边缘差值u1和左右边缘差值u2之和，即u＝u1+u2，且u1和u2如下式所示，其中
p∈[0,M-1]，q∈[0,N-1]。
7.如权利要求1所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像，具体包括：
步骤3.1：利用PFT提取重构图像的傅里叶变换，得到高精度分量g(x,y)的FFT形式，即：
其中x∈[0,M-1]and y∈[0,N-1]，f(x,y)为初始图像，g(x,y)为预期的高精度图像，e(x,y)为平滑误差图像。
步骤3.2：对高精度分量g(x,y)的FFT形式进行逆傅里叶变换，并计算强度，即可得到高精度重构图像。
8.如权利要求1或6所述基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的3)利用基于PFT的FPM重构算法去除平滑的误差图像，得到高精度重构图像，所述的高精度重构图像的原始图像计算模型为：
其中，r＝(x,y)为空间域坐标，k＝(kx,ky)为傅里叶域坐标；P(k)为系统光瞳函数，O(k-ki)为光瞳面的出射波，i为LR图像的序号，Ii(r)即为原始图像集；权重因子wi对应LED强度波动误差，相位项ejΦ对应光瞳函数中的像差，In,i对应加性噪声，△ki对应位置偏差；F为PFT算法。
9.如权利要6所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述步骤2.3.1的集合的绝对值计算对应集合的基数，即|G(a)|+|H(a)|＝4。
10.如权利要6所述的基于PFT的傅里叶叠层显微术高精度重构图像方法，其特征在于：所述的步骤2.3.2中K表示Ψ空间中某个像素到像素a的距离c，即：
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210624880.4A CN114998161A (zh) | 2022-06-02 | 2022-06-02 | 一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210624880.4A CN114998161A (zh) | 2022-06-02 | 2022-06-02 | 一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN114998161A true CN114998161A (zh) | 2022-09-02 |
Family
ID=83030839
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210624880.4A Pending CN114998161A (zh) | 2022-06-02 | 2022-06-02 | 一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114998161A (zh) |
Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20040257530A1 (en) * | 2003-06-20 | 2004-12-23 | Visx, Inc. | Wavefront reconstruction using fourier transformation and direct integration |
| CN103236049A (zh) * | 2013-05-05 | 2013-08-07 | 西安电子科技大学 | 基于序列相似性插值的部分k空间图像重构方法 |
| CN105930002A (zh) * | 2010-05-03 | 2016-09-07 | 平蛙实验室股份公司 | 通过层析成像重构进行触摸确定 |
| CN107219527A (zh) * | 2017-05-27 | 2017-09-29 | 吉林大学 | 一种周期型双阵列通道式成像系统的单快拍快速成像方法 |
| US20180204313A1 (en) * | 2015-07-20 | 2018-07-19 | Okinawa Institute Of Science And Technology School Corporation | 2d discrete fourier transform with simultaneous edge artifact removal for real-time applications |
| CN111062889A (zh) * | 2019-12-17 | 2020-04-24 | 北京理工大学 | 一种用于傅里叶叠层显微成像技术的光强校正方法 |
| CN111158131A (zh) * | 2019-12-31 | 2020-05-15 | 杭州电子科技大学 | 一种基于傅里叶叠层成像的led矩阵校正方法 |
| US20200378745A1 (en) * | 2019-05-28 | 2020-12-03 | University Of Central Florida Research Foundation, Inc. | Iterative Optical Diffraction Tomography (iODT) Method and Applications |
| CN112465703A (zh) * | 2020-12-01 | 2021-03-09 | 大连海事大学 | 一种含频谱偏移量估计的傅里叶叠层成像深度学习图像重构方法 |
| US20210118110A1 (en) * | 2019-10-21 | 2021-04-22 | Illumina, Inc. | Increased Calculation Efficiency for Structured Illumination Microscopy |
| CN113671682A (zh) * | 2021-08-23 | 2021-11-19 | 北京理工大学重庆创新中心 | 一种基于傅里叶叠层显微成像的频域光源位置精确校正方法 |
-
2022
- 2022-06-02 CN CN202210624880.4A patent/CN114998161A/zh active Pending
Patent Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20040257530A1 (en) * | 2003-06-20 | 2004-12-23 | Visx, Inc. | Wavefront reconstruction using fourier transformation and direct integration |
| CN105930002A (zh) * | 2010-05-03 | 2016-09-07 | 平蛙实验室股份公司 | 通过层析成像重构进行触摸确定 |
| CN103236049A (zh) * | 2013-05-05 | 2013-08-07 | 西安电子科技大学 | 基于序列相似性插值的部分k空间图像重构方法 |
| US20180204313A1 (en) * | 2015-07-20 | 2018-07-19 | Okinawa Institute Of Science And Technology School Corporation | 2d discrete fourier transform with simultaneous edge artifact removal for real-time applications |
| CN107219527A (zh) * | 2017-05-27 | 2017-09-29 | 吉林大学 | 一种周期型双阵列通道式成像系统的单快拍快速成像方法 |
| US20200378745A1 (en) * | 2019-05-28 | 2020-12-03 | University Of Central Florida Research Foundation, Inc. | Iterative Optical Diffraction Tomography (iODT) Method and Applications |
| US20210118110A1 (en) * | 2019-10-21 | 2021-04-22 | Illumina, Inc. | Increased Calculation Efficiency for Structured Illumination Microscopy |
| CN111062889A (zh) * | 2019-12-17 | 2020-04-24 | 北京理工大学 | 一种用于傅里叶叠层显微成像技术的光强校正方法 |
| CN111158131A (zh) * | 2019-12-31 | 2020-05-15 | 杭州电子科技大学 | 一种基于傅里叶叠层成像的led矩阵校正方法 |
| CN112465703A (zh) * | 2020-12-01 | 2021-03-09 | 大连海事大学 | 一种含频谱偏移量估计的傅里叶叠层成像深度学习图像重构方法 |
| CN113671682A (zh) * | 2021-08-23 | 2021-11-19 | 北京理工大学重庆创新中心 | 一种基于傅里叶叠层显微成像的频域光源位置精确校正方法 |
Non-Patent Citations (8)
| Title |
|---|
| AN PAN, ET.AL: "Edge effect removal in Fourier ptychographic microscopy via perfect Fourier transformation(PFT)", 《ARXIV:2009.03138V2》 * |
| AN PAN, ET.AL: "Edge effect removal in Fourier ptychographic microscopy via perfect Fourier transformation(PFT)", 《ARXIV》 * |
| LIHENG BIAN, ET.AL: "Fourier ptychographic reconstruction using poisson maximum likelihood and truncated wirtinger gradient", 《OPEN ACCESS》 * |
| LIHENG BIAN, ET.AL: "Fourier ptychographic reconstruction using wirtinger flow optimization", 《OPEN ACCESS》 * |
| MENG XIANG, ET.AL: "Coherent synthetic aperture imaging for visible remote sensing via reflective fourier ptychography", 《ARXIV》 * |
| 刘厚鹏等: "基于PFT算法的输送带钢丝绳芯X线图像的缺陷检测", 《矿山机械》 * |
| 张记龙等: "连续傅里叶变换的高精度PFT算法", 《太远机械学院学报》 * |
| 潘安: "高分辨率大视场快速傅里叶叠层显微成像方法研究", 《中国博士学位论文全文数据库 (信息科技辑)》 * |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107180411B (zh) | 一种图像重构方法及系统 | |
| CN111650738A (zh) | 一种基于深度学习的傅里叶叠层显微图像重构方法及装置 | |
| CN110415201B (zh) | 基于结构光和深度学习的单曝光超分辨率成像方法及装置 | |
| CN110717936B (zh) | 一种基于相机姿态估计的图像拼接方法 | |
| CN108537735B (zh) | 一种焦平面太赫兹成像的图像拼接方法 | |
| CA3206206A1 (en) | Device and method for correspondence analysis in images | |
| WO2019068201A1 (en) | METHODS AND SYSTEMS FOR ESTIMATING 3D STRUCTURE USING NON-UNIFORM REFINEMENT | |
| Pradeep et al. | Improving shape from focus using defocus cue | |
| CN116958419A (zh) | 一种基于波前编码的双目立体视觉三维重建系统和方法 | |
| CN114241072B (zh) | 一种叠层成像重建方法及系统 | |
| Zhang et al. | Group-based sparse representation for Fourier ptychography microscopy | |
| Wu et al. | Adaptive correction method of hybrid aberrations in Fourier ptychographic microscopy | |
| CN111694016A (zh) | 一种非干涉合成孔径超分辨成像重构方法 | |
| CN114998161A (zh) | 一种基于完美傅里叶变换的傅里叶叠层显微术高精度重构图像方法 | |
| Choudhary et al. | Mathematical modeling and simulation of multi-focus image fusion techniques using the effect of image enhancement criteria: a systematic review and performance evaluation | |
| Hui et al. | Image restoration for optical synthetic aperture system via patched maximum–minimum intensity prior and unsupervised DenoiseNet | |
| CN116955900A (zh) | 一种相位解包裹方法 | |
| CN114881907B (zh) | 光学显微图像多景深焦点合成方法、系统及图像处理方法 | |
| CN113554555B (zh) | 相位畸变补偿方法、装置、计算机设备和存储介质 | |
| CN113658317B (zh) | 电子显微镜连拍图像处理方法和装置 | |
| CN115752760A (zh) | 一种适用于微振动环境的相位恢复算法 | |
| CN115601569A (zh) | 一种基于改进piifd的异源图像优化匹配方法及系统 | |
| Qiao et al. | Stainless steel cylindrical pot outer surface defect detection method based on cascade neural network | |
| CN114926357A (zh) | 计算显微成像系统的led阵列光源位姿自校正方法 | |
| CN116485658B (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 | ||
| WD01 | Invention patent application deemed withdrawn after publication | ||
| WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20220902 |