一种频谱修复与自适应寻优的显微照明误差校正方法

本发明涉及计算成像与显微成像，尤其涉及一种频谱修复与自适应寻优的显微照明误差校正方法。

背景技术：

1、傅里叶叠层显微成像（fourier ptychographic microscopy，fpm）是一种融合合成孔径与相位恢复的计算成像技术。该技术通常采用可编程光源阵列从不同空间位置依次照明样品，使样品频谱中处于物镜通带之外的高频信息在倾斜照明作用下移动至物镜通带内。相机采集不同照明位置对应的低分辨率图像后，通过频谱拼接和迭代相位恢复，重建样品的高分辨率振幅信息和相位信息。

2、在傅里叶叠层显微成像过程中，各光源单元的空间位置用于确定相应的照明波矢以及物体频谱子区域的位置。受光源阵列安装偏移、阵列旋转、器件制造误差和装配误差等因素影响，光源单元的实际位置与标称位置之间容易产生几何失调。若仍按照标称位置执行频谱拼接，相邻频谱子区域之间将出现位置偏差，进而使迭代重建结果产生褶皱、重影和相位失真。

3、现有照明误差校正方式包括根据成像装置的机械尺寸进行几何标定，以及根据低分辨率图像或其频域特征反演光源位置参数。机械标定方式依赖测量设备和测量精度，难以充分补偿光源阵列的整体平移误差、旋转误差及装配误差。基于频域特征的校正方式能够利用倾斜照明形成的圆形光瞳特征确定照明波矢，但低分辨率图像相对边界之间的不连续性会在离散傅里叶变换结果中形成沿两个正交频率轴分布的十字型伪影，从而干扰圆形光瞳特征的提取以及光源位置参数的计算。

4、常规周期性平滑分解方法虽然能够减弱由图像边界不连续引起的十字型伪影，但在提取周期性分量的过程中可能造成非直流频谱能量衰减。与此同时，采用固定搜索步长求解光源阵列的平移误差和旋转误差时，较大的搜索步长虽然能够扩大搜索范围，但难以获得较高的局部定位精度；较小的搜索步长虽然有利于精细定位，但容易受初始参数影响，难以兼顾全局搜索范围和局部定位精度。

5、因此，如何在保留有效频谱信息的基础上修复正交轴区域的异常频谱，补偿频谱修复造成的非直流频谱能量变化，并根据寻优状态分别调整平移参数和旋转参数的搜索步长，是显微照明误差校正过程中需要解决的技术问题。

技术实现思路

1、针对上述问题，本发明提供一种频谱修复与自适应寻优的显微照明误差校正方法，以解决正交轴频谱伪影干扰圆形频谱特征提取、频谱修复造成非直流频谱能量变化，以及固定搜索步长难以兼顾搜索范围和定位精度的技术问题。

2、为解决上述技术问题，本发明提供一种频谱修复与自适应寻优的显微照明误差校正方法，包括以下步骤：

3、s1、频谱修复及能量补偿：控制光源阵列中的不同光源单元依次照明样品，获取各照明位置对应的低分辨率图像；对所述低分辨率图像进行二维傅里叶变换，获得原始频谱；对所述低分辨率图像进行周期性平滑分解，获得周期性分量，并对所述周期性分量进行二维傅里叶变换，获得周期性频谱；对所述周期性频谱中位于正交轴区域的异常频谱像素进行邻域插值，根据所述原始频谱的非直流频谱能量与插值频谱的非直流频谱能量确定能量补偿系数，并在保持直流分量不变的条件下对所述插值频谱的非直流分量进行能量补偿，获得恢复频谱特征图像。

4、s2、圆形频谱特征提取：将所述恢复频谱特征图像输入预训练特征提取网络，获得表征光瞳位置的圆形频谱特征。

5、s3、照明几何参数自适应寻优：分别设置平移搜索步长和旋转搜索步长，基于当前平移搜索步长和当前旋转搜索步长生成候选几何失调参数；根据各候选几何失调参数确定对应的仿真光瞳位置，生成与步骤s2所得圆形频谱特征具有相同数据形式的仿真圆形频谱特征；计算所述仿真圆形频谱特征与步骤s2所得圆形频谱特征之间的特征误差，根据所述特征误差更新当前最优几何失调参数，并在连续预设次数未获得更小特征误差时分别衰减所述平移搜索步长和所述旋转搜索步长；

6、s4、照明误差校正及傅里叶叠层显微重建：根据所得几何失调参数修正各光源单元的空间位置及对应照明波矢，并根据修正后的照明波矢进行傅里叶叠层显微重建。

7、进一步地，在所述步骤s1中，所述周期性平滑分解满足：

8、；

9、；

10、其中，为所述低分辨率图像，为周期性分量，为平滑分量，为根据所述低分辨率图像的边界不连续性构建的边界项，为二维离散拉普拉斯算子。

11、进一步地，所述二维离散拉普拉斯算子满足：

12、；

13、其中，分别表示当前像素左侧、右侧、上侧、下侧的平滑分量值。

14、进一步地，将所述边界项初始化为零，并按照下式对其边界位置进行赋值：

15、；

16、；

17、；

18、；

19、其中，为所述低分辨率图像的行数，为所述低分辨率图像的列数；和用于表征相对行边界之间的不连续性，和用于表征相对列边界之间的不连续性，符号“”表示在当前边界项数值的基础上进行更新。

20、进一步地，对离散泊松方程两侧进行二维傅里叶变换，得到：

21、；

22、根据二维离散拉普拉斯算子的频域性质，平滑分量的频域表示为：

23、；

24、其中，，为垂直频率方向的离散索引，为水平频率方向的离散索引，为平滑分量的频域表示，表示二维傅里叶变换。

25、所述周期性频谱满足：

26、；

27、其中，为周期性分量对应的周期性频谱。

28、进一步地，对所述周期性频谱中位于正交轴区域的异常频谱像素进行邻域插值，满足：

29、；

30、其中，为插值频谱中坐标为的频谱值，为频谱中心的坐标；当频谱像素位于经过频谱中心的垂直频率轴且不位于频谱中心时，采用其沿行方向相邻的两个频谱像素进行平均插值；当频谱像素位于经过频谱中心的水平频率轴且不位于频谱中心时，采用其沿列方向相邻的两个频谱像素进行平均插值；当频谱像素位于频谱中心时，保持所述频谱中心的直流分量不变。

31、进一步地，分别计算原始频谱与插值频谱的非直流频谱能量：

32、；

33、；

34、其中，为低分辨率图像的原始频谱，为扣除频谱中心直流能量后的原始非直流频谱能量，为扣除频谱中心直流能量后的插值后非直流频谱能量。

35、根据所述原始频谱非直流能量与所述插值频谱非直流能量确定能量补偿系数：

36、；

37、其中，为能量补偿系数。

38、进一步地，保持插值频谱的直流分量不变，并利用所述能量补偿系数对插值频谱的非直流分量进行补偿：

39、；

40、其中，为恢复频谱特征图像中坐标为的频谱值。

41、进一步地，所述特征提取网络为u-net网络。所述u-net网络包括编码器、瓶颈层和解码器；编码器的卷积层通道数依次为8、16、32和64，瓶颈层具有128个通道，解码器的卷积层通道数依次为64、32、16和8；所述u-net网络采用3×3卷积核，并通过2×2最大池化操作进行下采样。

42、进一步地，所述u-net网络的输入为尺寸为256×256的单通道恢复频谱特征图像；训练样本包括300张恢复频谱特征图像，并按照7：3的比例划分为训练集和验证集；训练过程采用二元交叉熵作为损失函数，采用平均绝对误差作为评价指标，并利用反向传播更新网络参数。

43、进一步地，将平移搜索步长初始化为200个位置偏移坐标单位，将旋转搜索步长初始化为0.5度，将当前最优第一平移量、当前最优第二平移量和当前最优旋转量均初始化为零，将局部未优化计数器和累计未优化计数器均初始化为零，并将历史最小特征误差初始化为无穷大。

44、进一步地，基于当前平移搜索步长和当前旋转搜索步长，利用均匀分布随机数生成多组候选第一平移偏移量、候选第二平移偏移量和候选旋转偏移量；将所述候选第一平移偏移量、所述候选第二平移偏移量和所述候选旋转偏移量，分别叠加至当前最优第一平移量、当前最优第二平移量和当前最优旋转量，获得多组候选几何失调参数。

45、进一步地，基于各组候选几何失调参数更新照明波矢坐标，根据更新后的照明波矢坐标确定对应的仿真光瞳位置，并根据所述仿真光瞳位置生成与步骤s2所得圆形频谱特征具有相同图像尺寸和数据形式的仿真圆形频谱特征；计算所述仿真圆形频谱特征与步骤s2所得圆形频谱特征之间的平均绝对误差。

46、当候选几何失调参数对应的平均绝对误差小于历史最小特征误差时，将该候选几何失调参数更新为当前最优几何失调参数，将对应的平均绝对误差更新为历史最小特征误差，并将局部未优化计数器和累计未优化计数器均清零。

47、当候选几何失调参数对应的平均绝对误差不小于所述历史最小特征误差时，保持当前最优几何失调参数不变，并分别增加所述局部未优化计数器和所述累计未优化计数器的计数值。

48、当所述局部未优化计数器的计数值达到5时，将平移搜索步长衰减至衰减前的40%，将旋转搜索步长衰减至衰减前的55%，并将所述局部未优化计数器清零；所述累计未优化计数器在搜索步长衰减时保持当前计数值不变。随后，基于衰减后的平移搜索步长和旋转搜索步长继续生成候选几何失调参数。

49、进一步地，当所述历史最小特征误差不大于预设误差阈值，或者所述累计未优化计数器的计数值达到50时，终止寻优，并将当前最优几何失调参数作为光源阵列几何失调参数输出。

50、进一步地，将所述光源阵列几何失调参数补偿至各光源单元的标称空间位置，获得各光源单元的校正空间位置；根据所述校正空间位置确定对应的照明波矢，并采用自适应相位恢复算法对不同照明位置对应的低分辨率图像进行迭代重建，获得样品的振幅信息和相位信息。

51、进一步地，在所述步骤s3中，所述几何失调参数包括光源阵列沿阵列平面内相互正交的第一方向和第二方向产生的第一平移量和第二平移量，以及所述光源阵列绕光轴产生的旋转量；以当前最优几何失调参数为中心，根据当前平移搜索步长分别沿所述第一方向和所述第二方向的正向及反向生成候选平移参数，并根据当前旋转搜索步长沿绕光轴的正向及反向生成候选旋转参数，由所述候选平移参数和所述候选旋转参数组成多组候选几何失调参数。

52、进一步地，在所述步骤s4中，根据所述第一平移量、所述第二平移量和所述旋转量，对各光源单元的标称空间位置进行平移变换和旋转变换，获得各光源单元的校正空间位置；根据所述校正空间位置重新计算各光源单元对应的照明波矢，并根据重新计算的照明波矢确定各低分辨率图像对应的频谱子区域，对各频谱子区域依次进行相位恢复迭代，获得样品的高分辨率振幅信息和相位信息。

53、与现有技术相比，本发明至少具有以下有益效果：

54、（1）本发明通过周期性平滑分解减弱图像相对边界不连续性造成的频谱伪影，并对正交轴区域的异常频谱像素进行邻域插值，能够在抑制十字型伪影的同时保留正交轴附近的频谱连续性，降低伪影对光瞳圆形特征提取的干扰。

55、（2）本发明根据原始频谱与插值频谱的非直流频谱能量确定能量补偿系数，在保持插值频谱的直流分量不变的条件下，对所述插值频谱的非直流分量进行能量补偿，能够减少异常频谱像素邻域插值造成的有效频谱能量衰减，使恢复频谱特征图像更准确地表征光瞳位置。

56、（3）本发明分别设置平移搜索步长和旋转搜索步长，并在连续预设次数未获得更小特征误差时按照不同衰减比例缩小两类搜索步长，使参数搜索过程能够由较大范围的粗搜索逐步过渡至较小范围的精细搜索。

57、（4）本发明将频谱修复、光瞳特征提取、光源阵列几何失调参数寻优和傅里叶叠层显微重建相结合，能够对光源阵列沿阵列平面内两个相互正交方向的平移误差以及绕光轴的旋转误差进行校正，减少光源位置偏差对重建振幅信息和相位信息的影响。