CN100561378C - 一种像场弯曲的数值校正方法 - Google Patents

Abstract

本发明公开了一种像场弯曲的数值校正方法，对数字全息图进行数值重建，并进行相位去包裹运算；对相位去包裹之后的结果进行一次拟合；根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果。本发明可以有效地消除数字全息显微术中因为显微物镜的存在所引入的像差，仅需要拍摄一幅全息图即可以得到无像差的再现像，适用于动态过程的测量与监测，具有简单，快速，准确，方便等优点，在数字全息显微术等领域具有一定的应用前景。

Description

一种像场弯曲的数值校正方法

技术领域

本发明涉及数字全息技术领域，特别涉及一种数字全息显微术中像场弯曲的数值校正方法。

背景技术

数字全息技术以固体成像器件CCD(Charge Coupled Device，电荷耦合器件)或CMOS(Complementary Metal-Oxide Semiconductor，金属氧化物半导体)替代传统的银盐干板来记录全息图，通过计算机数值模拟数字全息图的衍射再现过程，利用快速傅里叶变换算法及频谱滤波处理等手段获得物光波场的振幅和相位信息，实现三维物光波场重构，避免了传统的湿化学处理过程，具有快速，实时等优点。将该技术与光学显微术相结合，不仅可以对振幅型样品成像，而且无需经过染色处理即可以直接对相位型样品成像，获得其三维相位分布。不过，在显微成像过程中，显微物镜会使再现像产生像场弯曲，影响数值再现结果的相位分布。针对这一问题，Cuche等提出了利用单幅全息图对像场弯曲进行校正的数值方法(CucheE，Marquet P，Depeursinge C.Simultaneous amplitude-contrast and quantitative phase-contrastmicroscopy by numerical reconstruction of Fresnel off-axis holograms.Appl.Opt.，1999，38(34)：6994～7001)，Ferraro等提出一种类似二次曝光法的像场弯曲校正方法(Ferraro P，DeNicola S，Finizio A et al..Compensation of the inherent wave front curvature in digital holographic coherentmicroscopy for quantitative phase-contrast imaging.Appl.Opt.，2003，42(11)：1938～1946)，Miccro等提出利用泽尼克多项式的二维拟合方法以消除数值再现像的像场弯曲(Miccio L，Alfieri D，Grilli S et al..Direct full compensation of the aberrations in quantitative phase microscopy of thinobjects by a single digital hologram.Appl.Phys.Lett.，2007，90(4)：041104)。然而，Cuche的方法需要逐步调整校正参数，难以保证最终实验结果的可靠性，Ferraro的方法主要针对利用反射式光路测量MEMS器件的过程，在后续的相位去包裹运算中会引入误差，Miccro的方法需要经过多次计算，运算过程较为繁琐、复杂。

发明内容

为了克服现有技术不能有效消除数字全息显微术中的像场弯曲的不足，本发明提出数字全息显微术中一种基于曲面拟合的像场弯曲的数值校正方法，主要是在全息图数值再现过程中采用曲面拟合的思想利用计算机通过数值方法校正数字全息显微术中的像场弯曲，获得理想的无像差的数值再现像。

本发明解决其技术问题所采用的技术方案包括以下步骤：

1.对数字全息图进行数值重建，并进行相位去包裹运算；

2.对相位去包裹之后的结果进行一次拟合；

3.根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果。

所述步骤1包含以下步骤：

a.针对数字全息显微术中所拍摄的数字全息图，对其进行数值重建；

b.对数值重建结果进行相位去包裹运算，得到包含有像场弯曲的波前相位分布矩阵P″_i{x″_i，y″_i，z″_i}₁ ⁿ，其中，x″，y″分别为波前相位分布矩阵中垂直于光轴方向的采样点的位置坐标，z″为此矩阵中采样点的相位值坐标。

所述步骤2包含以下步骤：

a.根据所使用CCD的型号，确定波前相位分布矩阵中垂直于光轴方向的采样点的横向位置，即波前相位分布矩阵中位置坐标x″，y″的横向位置x′，y′，其中，x′＝x″×Δx，y′＝y″×Δy，Δx和Δy分别为CCD的横向和纵向像素尺寸；

b.根据记录全息图的激光器波长，确定波前相位分布矩阵中沿光轴方向的采样点的纵向位置，即波前相位分布矩阵中相位值坐标z″的纵向位置z′，其中z′＝z″×λ/2π，λ为激光波长；

c.得到波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ，并对其进行均匀抽样，得到波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.由上述波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ得到超定方程组Ax＝b，其中，

$A=\left[\begin{array}{cccc}{x}_1& y_1& z_1& 1\\ x_2& y_2& z_2& 1\\ M& M& M& M\\ x_n& y_n& z_n& 1\end{array}\right],$ x＝(a，b，c，d)^T， $b=\left[\begin{array}{c}{x}_1^2+y_1^2+z_1^2\\ x_2^2+y_2^2+z_2^2\\ M\\ x_n^2+y_n^2+z_n^2\end{array}\right];$

e.应用奇异值分解法求得上述超定方程组的最小二乘解x，并经过坐标转换公式获得拟合球面的球心坐标(x₀，y₀，z₀)和曲率半径D，其中，

$x_0=-\frac{a}{2},$ $y_0=-\frac{b}{2},$ $z_0=-\frac{c}{2},$ $D=\sqrt{\frac{a^2+b^2+c^2-4d}{2}};$

f.根据拟合球面方程 $\mathrm{\Φ}(x,y)=\exp [-\frac{\mathrm{i\π}}{\mathrm{\λD}}({\left({x-x}_0\right)}^2+{(y-y_0)}^2)],$ 得到拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ；

g.将波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ与拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ作差，获得经过

一次拟合的再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ。

所述步骤3，分为三种不同方法。其中，方法一包含以下步骤：

a.借助软件自动寻找再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中的z_i的最大值z_max，并将最大值z_max乘以系数K作为阈值点，K取值为10％～100％；

b.将再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中z_i大于所选取阈值点的部分矩阵元去除，并确定再现像矩阵中剩余点的位置；

c.根据再现像矩阵中剩余点的位置，对波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

所述步骤3的方法二包含以下步骤：

a.根据再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ得到的强度再现结果，选取再现结果中相邻矩阵元素灰度差值小于整个图像上灰度分布差的10％灰度分布的部分作为相对均匀部分；

b.根据上述灰度分布均匀部分位置对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

c.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合3的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

所述步骤3的方法三包含以下步骤：

a.选择再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ，对其进行灰度直方图统计；

b.根据统计结果，确定灰度直方图统计中灰度分布概率最大点位置，并选择其周围L比例范围内的矩阵元素，L取值为10％～100％；

c.根据所选矩阵元素对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

本发明一种像场弯曲的数值校正方法的步骤3中的三种方法可以根据实际情况选择使用，当样品上有多处较为平坦的区域时，选择方法一或者是方法二的效果比较好。当样品整体不平坦，对透射或反射自它的光波波前的相位调制的区域很多时，选择方法三的效果较好。

此外，当样品自身对透射或反射自它的光波波前的相位起伏影响较小时，经过一次校正获得的再现像的结果即比较理想，可以不进行本发明一种像场弯曲的数值校正方法的步骤3中的二次拟合过程。

本发明的有益效果是：由于采用基于曲面拟合的数值校正方法可以有效地消除数字全息显微术中因为显微物镜的存在所引入的像差。该方法有别于二次曝光法需要拍摄两幅全息图并进行相减操作，仅需要拍摄一幅全息图即可以得到无像差的再现像，适用于动态过程的测量与监测。此外，该方法完全在计算机的数值计算过程结合数学方法，利用曲面拟合的思想对像差进行校正，与Cuche、Ferraro和Miccro等人的方法相比，具有简单，快速，准确，方便等优点，在数字全息显微术等领域具有一定的应用前景。

下面结合附图和实施例对本发明进一步说明。

附图说明

图1是本发明一种像场弯曲的数值校正方法的实施例1中大蒜表皮细胞的光学显微镜示意图；

图2是本发明一种像场弯曲的数值校正方法中未经过像差校正的大蒜表皮细胞的数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图3是本发明一种像场弯曲的数值校正方法中经过一次拟合像差校正的大蒜表皮细胞数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图4是本发明一种像场弯曲的数值校正方法中采用方法一经过二次拟合像差校正的大蒜表皮细胞数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图5是本发明一种像场弯曲的数值校正方法中实施例2中蝉翼的光学显微镜照片；

图6是本发明一种像场弯曲的数值校正方法中未经过像差校正的蝉翼的数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图7是本发明一种像场弯曲的数值校正方法中经过一次拟合像差校正的大蒜表皮细胞数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图8是本发明一种像场弯曲的数值校正方法中采用方法二经过二次拟合像差校正的蝉翼数值再现像的相位分布结果，(a)二维相位分布结果，(b)三维相位分布结果；

图9是本发明一种像场弯曲的数值校正方法中实施例3中针对蝉翼的步骤3a的直方图统计结果。

图10是本发明一种像场弯曲的数值校正方法中采用方法三经过二次拟合像差校正的蝉翼数值再现像的二维相位分布结果。

具体实施方式

方法实施例1：如图1到图4所示，其是本发明一种像场弯曲的数值校正方法的实施效果图，包含步骤如下：

1.对数字全息图进行数值重建，并进行相位去包裹运算；

2.对相位去包裹之后的结果进行一次拟合；

3.根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果。

上述对数字全息图进行数值重建，并进行相位去包裹运算的步骤1，包含以下步骤：

a.针对图1中所示的大蒜表皮细胞利用数字全息显微术拍摄数字全息图，对其进行数值重建；

b.对数值重建结果进行相位去包裹运算，得到包含有像场弯曲的波前相位分布矩阵P″_i{x″_i，y″_i，z″_i}₁ ⁿ，其中，x″，y″分别为波前相位分布矩阵中垂直于光轴方向的采样点的位置坐标，z″为此矩阵中采样点的相位值坐标。图2所示即为本发明一种像场弯曲的数值校正方法中未经过像差校正的大蒜表皮细胞的数值再现像的相位分布结果，其中，(a)为二维相位分布结果，

(b)为三维相位分布结果。

上述对相位去包裹之后的结果进行一次拟合的步骤2，包含以下步骤：

a.根据所使用CCD的型号(Pulnix1400CL)，确定波前相位分布矩阵中垂直于光轴方向的采样点的横向位置，即波前相位分布矩阵中位置坐标x″，y″的横向位置x′，y′，其中，x′＝x″×Δx，y′＝y″×Δy，Δx和Δy分别为CCD的横向和纵向像素尺寸，取Δx＝Δy＝0.45μm；

b.根据记录全息图的激光器波长(He-Ne激光器)，确定波前相位分布矩阵中沿光轴方向的采样点的纵向位置，即波前相位分布矩阵中相位值坐标z″的纵向位置z′，其中z′＝z″×λ/2π，λ为激光波长，为632.8nm；

c.得到波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ，并按照比例对其进行均匀抽样，得到波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.由上述波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ得到超定方程组Ax＝b，其中， $A=\left[\begin{array}{cccc}{x}_1& y_1& z_1& 1\\ x_2& y_2& z_2& 1\\ M& M& M& M\\ x_n& y_n& z_n& 1\end{array}\right],$ x＝(a，b，c，d)^T， $b=\left[\begin{array}{c}{x}_1^2+y_1^2+z_1^2\\ x_2^2+y_2^2+z_2^2\\ M\\ x_n^2+y_n^2+z_n^2\end{array}\right];$

e.应用奇异值分解法求得上述超定方程组的最小二乘解x，并经过坐标转换公式获得拟合球面的球心坐标(x₀，y₀，z₀)和曲率半径D，其中， $x_0=-\frac{a}{2},$ $y_0=-\frac{b}{2},$ $z_0=-\frac{c}{2},$ $D=\sqrt{\frac{a^2+b^2+c^2-4d}{2}};$

f.根据拟合球面方程 $\mathrm{\Φ}(x,y)=\exp [-\frac{\mathrm{i\π}}{\mathrm{\λD}}({\left({x-x}_0\right)}^2+{(y-y_0)}^2)],$ 得到拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ；

g.将波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ与拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ作差，获得经过一次拟合的再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ，图3所示即为经过一次拟合像差校正的大蒜表皮细胞数值再现像的相位分布结果，其中，(a)为二维相位分布结果，(b)为三维相位分布结果。

上述根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果的步骤3，包含以下步骤：

a.借助软件自动寻找再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中的z_i的最大值z_max，并将最大值z_max乘以系数K作为阈值点，K取值为60％；

b.将再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中z_i大于所选取阈值点的部分矩阵元去除，并确定再现像矩阵中剩余点的位置；

c.根据再现像矩阵中剩余点的位置，对波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ，图4所示即为本发明一种像场弯曲的数值校正方法中本实施例中方法经过二次拟合像差校正得到的大蒜表皮细胞数值再现像的相位分布结果，其中，(a)为二维相位分布结果，(b)为三维相位分布结果。

此外，对于实施例1中的大蒜表皮细胞，由于其对光波波前的相位起伏影响较为均匀，一次拟合像差校正获得的再现像的结果即比较理想，可以不进行本发明一种像场弯曲的数值校正方法中步骤3中的二次拟合过程。

方法实施例2：如图5到图8所示，其是本发明一种像场弯曲的数值校正方法的实施效果图，包含步骤如下：

1.对数字全息图进行数值重建，并进行相位去包裹运算；

2.对相位去包裹之后的结果进行一次拟合；

3.根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果。

上述对数字全息图进行数值重建，并进行相位去包裹运算的步骤1，包含以下步骤：

a.针对图5中所示的蝉翼利用数字全息显微术拍摄数字全息图，对其进行数值重建；

b.对数值重建结果进行相位去包裹运算，得到包含有像场弯曲的波前相位分布矩阵P″_i{x″_i，y″_i，z″_i}₁ ⁿ，其中，x″，y″分别为波前相位分布矩阵中垂直于光轴方向的采样点的位置坐标，z″为此矩阵中采样点的相位值坐标。图6所示即为本发明一种像场弯曲的数值校正方法中未经过像差校正的蝉翼的数值再现像的相位分布结果，其中，(a)为二维相位分布结果，(b)为三维相位分布结果。

上述对相位去包裹之后的结果进行一次拟合的步骤2，包含以下步骤：

a.根据所使用CCD的型号(Pulnix1400CL)，确定波前相位分布矩阵中垂直于光轴方向的采样点的横向位置，即波前相位分布矩阵中位置坐标x″，y″的横向位置x′，y′，其中，x′＝x″×Δx，y′＝y″×Δy，Δx和Δy分别为CCD的横向和纵向像素尺寸，取Δx＝Δy＝0.45μm；

b.根据记录全息图的激光器波长(He-Ne激光器)，确定波前相位分布矩阵中沿光轴方向的采样点的纵向位置，即波前相位分布矩阵中相位值坐标z″的纵向位置z′，其中z′＝z″×λ/2π，λ为激光波长，为632.8nm；

c.得到波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ，并按照比例对其进行均匀抽样，得到波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.由上述波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ得到超定方程组Ax＝b，其中， $A=\left[\begin{array}{cccc}{x}_1& y_1& z_1& 1\\ x_2& y_2& z_2& 1\\ M& M& M& M\\ x_n& y_n& z_n& 1\end{array}\right],$ x＝(a，b，c，d)^T， $b=\left[\begin{array}{c}{x}_1^2+y_1^2+z_1^2\\ x_2^2+y_2^2+z_2^2\\ M\\ x_n^2+y_n^2+z_n^2\end{array}\right];$

e.应用奇异值分解法求得上述超定方程组的最小二乘解x，并经过坐标转换公式获得拟合球面的球心坐标(x₀，y₀，z₀)和曲率半径D，其中， $x_0=-\frac{a}{2},$ $y_0=-\frac{b}{2},$ $z_0=-\frac{c}{2},$ $D=\sqrt{\frac{a^2+b^2+c^2-4d}{2}};$

f.根据拟合球面方程 $\mathrm{\Φ}(x,y)=\exp [-\frac{\mathrm{i\π}}{\mathrm{\λD}}({\left({x-x}_0\right)}^2+{(y-y_0)}^2)],$ 得到拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ；

g.将波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ与拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ作差，获得经过一次拟合的再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ，图7所示即为本发明一种像场弯曲的数值校正方法中经过一次拟合像差校正的大蒜表皮细胞数值再现像的相位分布结果，其中，(a)为二维相位分布结果，(b)为三维相位分布结果。

上述根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果的步骤3，包含以下步骤：

a.根据再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ得到的强度再现结果，选取再现结果中相邻矩阵元素灰度差值小于整个图像上灰度分布差的10％灰度分布的部分作为相对均匀部分；

b.根据上述灰度分布均匀部分位置对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

c.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合3的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ，图8所示即为本发明一种像场弯曲的数值校正方法中采用本实施例中方法经过二次拟合像差校正的蝉翼数值再现像的相位分布结果，其中，(a)为二维相位分布结果，(b)为三维相位分布结果。

实施例3：与实施例2类似，所述一种像场弯曲的数值校正方法包含步骤如下：

1.对数字全息图进行数值重建，并进行相位去包裹运算；

2.对相位去包裹之后的结果进行一次拟合；

3.根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果。

不同的是上述根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果的步骤3，包含以下步骤：

a.选择再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ，对其进行灰度直方图统计；

b.根据统计结果，确定灰度直方图统计中灰度分布概率最大点位置，并选择其周围80％比例范围内的矩阵元素；

c.根据所选矩阵元素对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

d.针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤2中的d到g步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

图9是实施例3中针对蝉翼的步骤3a的直方图统计结果。图10是采用方法三经过二次拟合像差校正的蝉翼数值再现像的二维相位分布结果。

Claims (3)

1、一种像场弯曲的数值校正方法，其特征在于包括下述步骤：

(1)对数字全息图进行数值重建，并进行相位去包裹运算，具体包括以下步骤：

(a)针对数字全息显微术中所拍摄的数字全息图，对其进行数值重建；

(b)对数值重建结果进行相位去包裹运算，得到包含有像场弯曲的波前相位分布矩阵P″_i{x″_i，y″_i，z″_i}₁ ⁿ；

(2)对相位去包裹之后的结果进行一次拟合，具体包括以下步骤：

(a)根据所使用CCD的型号，确定波前相位分布矩阵中垂直于光轴方向的采样点的横向位置，即波前相位分布矩阵中位置坐标x″，y″的横向位置x′，y′，其中，x′＝x″×Δx，y′＝y″×Δy，Δx和Δy分别为CCD的横向和纵向像素尺寸；

(b)根据记录全息图的激光器波长，确定波前相位分布矩阵中沿光轴方向的采样点的纵向位置，即波前相位分布矩阵中相位值坐标z″的纵向位置z′，其中z′＝z″×λ/2π，λ为激光波长；

(c)得到波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ，并对其进行均匀抽样，得到波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

(d)由上述波前相位分布的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ得到超定方程组Ax＝b，其中，

$A=\left[\begin{array}{cccc}{x}_1& y_1& z_1& 1\\ x_2& y_2& z_2& 1\\ M& M& M& M\\ x_n& y_n& z_n& 1\end{array}\right],$ x＝(a，b，c，d)^T， $b=\left[\begin{array}{c}{x}_1^2+y_1^2+z_1^2\\ x_2^2+y_2^2+z_2^2\\ M\\ x_n^2+y_n^2+z_n^2\end{array}\right];$

(e)应用奇异值分解法求得上述超定方程组的最小二乘解x，并经过坐标转换公式获得拟合球面的球心坐标(x₀，y₀，z₀)和曲率半径D，其中，

$x_0=-\frac{a}{2},$ $y_0=-\frac{b}{2},$ $z_0=-\frac{c}{2},$ $D=\frac{\sqrt{a^2+b^2+c^2-4d}}{2};$

(f)根据拟合球面方程 $\mathrm{\Φ}(x,y)=\exp [-\frac{\mathrm{i\π}}{\mathrm{\λD}}({(x-x_0)}^2+{(y-y_0)}^2)],$ 得到拟合球面矩阵

M_i{x_i，y_i，Φ_i}₁ ⁿ；

(g)将波前相位分布的位置矩阵P″_i{x′_i，y′_i，z′_i}₁ ⁿ与拟合球面矩阵M_i{x_i，y_i，Φ_i}₁ ⁿ作差，获得经过一次拟合的再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ；

(3)根据一次拟合效果对相位去包裹之后的结果进行二次拟合，得到理想的无像差的再现结果，具体包括以下步骤：

(a)借助软件自动寻找再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中的z_i的最大值z_max，并将最大值z_max乘以系数K作为阈值点，K取值为10％～100％；

(b)将再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ中z_i大于所选取阈值点的部分矩阵元去除，并确定再现像矩阵中剩余点的位置；

(c)根据再现像矩阵中剩余点的位置，对波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

(d)针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤(2)中的(d)到(g)步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

2、根据权利要求1的一种像场弯曲的数值校正方法，其特征在于所述步骤(3)可以替换为下述步骤：

(a)根据再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ得到的强度再现结果，选取再现结果中相邻矩阵元素灰度差值小于整个图像上灰度分布差的10％灰度分布的部分作为相对均匀部分；

(b)根据上述灰度分布均匀部分位置对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

(c)针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤(2)中的(d)到(g)步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

3、根据权利要求1的一种像场弯曲的数值校正方法，其特征在于所述步骤(3)可以替换为下述步骤：

(a)选择再现像矩阵A_i{x_i，y_i，z_i}₁ ⁿ，对其进行灰度直方图统计；

(b)根据统计结果，确定灰度直方图统计中灰度分布概率最大点位置，并选择其周围L比例范围内的矩阵元素，L取值为10％～100％；

(c)根据所选矩阵元素对波前相位分布位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ重新进行抽样，重新得到元素数目较波前相位分布的位置矩阵P′_i{x′_i，y′_i，z′_i}₁ ⁿ更少的拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ；

(d)针对拟合矩阵P_i{x_i，y_i，z_i}₁ ⁿ，重复步骤(2)中的(d)到(g)步，即可以获得经过二次拟合的理想的无像差的再现像矩阵A′_i{x_i，y_i，z_i}₁ ⁿ。

Images