CN103914826A - 一种中心矢状面的提取方法 - Google Patents

Abstract

本发明提供一种中心矢状面的提取方法，包括：提供医学图像及矢状面方程，所述矢状面方程包含待定系数，通过梯度下降法迭代系数，以确定中心矢状面方程；其中对任一次所述迭代包括：获取当前系数，通过所述当前系数对应的当前矢状面将所述医学图像分为两部分，获取所述两部分灰度的均方差测度函数及所述均方差测度函数对所述系数的偏导数，以迭代所述系数。本发明采用适用于求解非线性的梯度下降法对代价函数，即均方差测度函数进行迭代，并使用迭代终止时的系数对应的矢状面作为中心矢状面，提高中心矢状面提取的精准度。

Description

一种中心矢状面的提取方法

技术领域

本发明涉及一种医学扫描图像的处理方法，尤其涉及一种中心矢状面的提取方法。

背景技术

中心矢状面(mid-sagittal plane，MSP)是医学图像处理领域经常使用的一个平面，该平面是一个通过身体腹部与背部的垂直平面。通过所述中心矢状面将身体分成对称的左右两部分。提取初次扫描的3D医疗图像的中心矢状面是确定再次扫描方向的前提。

文献1：Babak A.Ardekani，etc，‘Automatic Detection of theMid-Sagittal Plane in 3-D Brain Images’，IEEE transactions on medicalimaging，vol.，16，N).6，December 1997提供的中心矢状面提取方法，采用单纯形下降法优化求解了三维头部图像的中心矢状面方程MSP，代价函数选取为互相关函数。考虑到医学图像的复杂性，并且代价函数一般为MSP中心矢状面方程参数的非线性函数，单纯形法是求解线性规划问题的常用方法，所以用于求MSP中心矢状面方程并不是非常适合的。

文献2：Bhanu Prakash K.N.，‘Locating a Mid-Sagittal Plane’，United States Patent Application Publication，Application number：US2007/0276219 A1提供的中心矢状面提取方法，通过用优化方法寻找相互正交的三个面，其中第一个平面为所求的中心矢状面。代价函数选取三个面与各自相平行的临近数个面对应点的灰度熵值，同样采用了基于单纯形的Nelder-Mead优化方法，并且需要统计三个平面附件的熵值，计算量大。

文献3：Yiwen Wang et al，‘System and method for consistentdetection of mid-sagittal planes for magnetic resonance brain scans’，United States Patent Application Publication，Application number：US2011/7986823B2提供的中心矢状面提取方法，为了保证从多次扫描检测出的中心矢状面的一致性(consistent)，把医生或者设备操作者标定的中心矢状面作为金标准，然后以互信息为代价函数衡量自动检测的中心矢状面同金标准的差异。实现上需要专家检测中心矢状面，如果不同病人的MSP第一次扫描，则均需人工检测时，操作上会非常麻烦。

发明内容

本发明解决的问题是提供一种中心矢状面的提取方法，用以简单快捷地提取中心矢状面。为了解决上述问题，本发明提供了一种中心矢状面的提取方法，包括：

提供医学图像及矢状面方程，所述矢状面方程包括至少一个待定系数，所述待定系数通过梯度下降法迭代确定，以确定中心矢状面方程，实现中心矢状面的提取；

其中对任一次所述迭代包括：获取当前系数，通过所述当前系数对应的当前矢状面将所述医学图像分为两部分；获取所述两部分灰度的均方差测度函数，及计算所述均方差测度函数对所述系数的偏导数，用以迭代系数；迭代后的系数作为下一次迭代的当前系数。

可选的，对于第一次迭代，通过对所述系数进行初始化，以获得所述当前系数；对于第一次后的迭代，将上一次迭代后获得的系数作为当前系数。

可选的，所述梯度下降法包括：提供初始的修正步长，并通过设置松弛因子降低所述修正步长，以迭代所述修正步长。

可选的，还包括：提供步长阈值和最大迭代次数，当所述修正步长不大于所述步长阈值，或者所述迭代次数达到所述最大迭代次数，所述迭代终止，并提取当前迭代后的系数对应的矢状面为中心矢状面。

可选的，所述初始修正步长的范围为0.01～0.03，所述松弛因子的范围为0.6～0.8，所述步长阈值的范围为所述初始修正步长的0.5％～1.5％，所述最大迭代次数的范围为50～150。

可选的，所述初始修正步长为0.02，所述松弛因子为0.6，所述步长阈值为所述初始修正步长的1％，所述最大迭代次数为100次。

可选的，所述均方差测度函数为I(p_i)和I(R(p_i))为分别位于所述当前矢状面两侧互对称的像素点p_i和R(p_i)的灰度值，L为所述当前矢状面的任一部分像素点p_i集合。

可选的，包括：读取所述像素p_i的灰度值，以所述当前矢状面为对称面，计算所述像素p_i的对称点的坐标作为R(p_i)。

可选的，若所述对称点的坐标为整数，则将对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值；若所述对称点的坐标为非整数，则在所述医学图像中进行插值，以获得对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值。

可选的，所述矢状面方程为x＝My+Nz+B，其中，所述M、N、B为所述待定的系数，其中，第i次迭代的当前系数对应为M_i、N_i及B_i，第一次迭代的当前系数为M₀、N₀及B₀。

可选的，还包括获取所述医学图像的质心坐标，并根据迭代后的系数M、N及所述质心坐标，获取对应的系数B。

可选的，若所述质心坐标为(x₀、y₀、z₀)。，所述系数的初始值为M₀＝0，N₀＝0，对应的初始矢状面方程为x＝B₀＝x₀。

可选的，所述偏导数为所述均方差测度函数对系数M、N的偏导数，所述均方差测度函数为C，所述均方差测度函数C对M和N的偏导数分别为 $S=\frac{\∂C}{\∂M},$ $T=\frac{\∂C}{\∂N}.$

可选的，所述第i次迭代后的系数为 其中，i为从1开始的整数，D_i为第i次优化的修正步长。

可选的，包括：i大于1的第i次迭代的当前迭代偏导数为S_i，T_i，上次迭代偏导数为S_i-1，T_i-1，

若所述(S_i-S_i-1)*(T_i-T_i-1)＜0，则对上次迭代的修正步长进行下降后作为当前迭代的修正步长；

若所述(S_i-S_i-1)*(T_i-T_i-1)＞0，则继续采用上次迭代使用的修正步长作为当前迭代的修正步长。

与现有技术相比，本发明具有以下优点：

因为中心矢状面的代价函数一般为非线性函数，采用求解线性规划问题的单纯形法求解所述代价函数并不合适，本发明则采用适用于求解非线性的梯度下降法对代价函数，即均方差测度函数进行迭代，并使用迭代终止时的系数对应的矢状面作为中心矢状面，提高中心矢状面提取的精准度；

进一步地，利用中心矢状面的质心特性作为约束条件，减少需要迭代的系数，降低迭代计算量，加快迭代的速度；

最后，所述均方差测度函数借鉴图像配准中的迭代优化方法，把中心矢状面图像分为左右两部分，通过比较两部分的灰度的均方差实现配准，对所述均方差测度函数进行迭代，获取待定系数及对应的中心矢状面，提高中心矢状面提取的精准度。

附图说明

图1所示为本发明一个实施例的中心矢状面的提取方法的流程示意图。

具体实施方式

在下面的描述中阐述了很多具体细节以便于充分理解本发明。但是本发明能够以很多不同于在此描述的其它方式来实施，本领域技术人员可以在不违背本发明内涵的情况下做类似推广，因此本发明不受下面公开的具体实施的限制。

其次，本发明利用示意图进行详细描述，在详述本发明实施例时，为便于说明，所述示意图只是实例，其在此不应限制本发明保护的范围。

本发明提供了一种中心矢状面的提取方法，包括：提供医学图像及矢状面方程，所述矢状面方程包含若干待定的系数，通过梯度下降法迭代所述系数，以确定系数进而提取中心矢状面；其中对任一次所述迭代包括：获取当前系数，通过所述当前系数对应的当前矢状面将所述医学图像分为两部分，获取所述两部分灰度的均方差测度函数及所述均方差测度函数对所述系数的偏导数，以迭代所述系数。

因为中心矢状面的代价函数一般为非线性函数，采用求解线性规划问题的单纯形法求解所述代价函数并不合适，本发明则采用适用于求解非线性的梯度下降法对代价函数，即均方差测度函数进行迭代，并使用迭代终止时的系数对应的矢状面作为中心矢状面，提高中心矢状面提取的精准度。

其中，所述梯度下降法包括：设置初始修正步长，并设置松弛因子以降低所述修正步长，获得满足迭代结束的系数确定为系数。

具体地，对于第一次迭代，通过对所述系数进行初始化，以获得所述当前系数；对于第一次后的迭代，将上一次迭代后获得的系数作为当前系数。

进一步地，所述两部分为左侧部分和右侧部分，其中所述左侧部分坐标小于所述当前矢状面上的任意一点坐标；所述右侧部分大于所述当前矢状面上的任意一点坐标。以下实施例，以所述两部分为左侧部分和右侧部分为例进行说明。

图1所示为本发明一个实施例的中心矢状面的提取方法的流程示意图，包括：

执行步骤S1，输入图像，计算灰度质心。所述图像为包含有背景和骨骼等信息的脑组织的三维图像。若忽略所述背景信息和骨骼信息，仅利用所述脑组织图像进行计算，可以提高计算速度，但是精度下降较严重。较佳地，在后续处理中，一直保留所述背景信息和骨骼信息。

计算质心坐标x₀，y₀，z₀；通过质心计算，获得所述图像的灰度质心坐标I(x₀、y₀、z₀)。

本发明中的矢状面的方程形式为x＝My+Nz+B，所述M、N、B为所述待定的系数。由于矢状面必然经过医学图像的质心，所以只有两个斜率项是未知的，故本文仅选取M，N两个优化变量，而每次迭代后由两个斜率的更新值计算新的截距B。

执行步骤S2，输入系数的初始值。所述初始值为(M₀、N₀、B₀)。其中，设定被优化系数的初始值M₀＝0，N₀＝0，接着，根据上步骤S1计算获得的质心坐标，及图像必经过质心的条件计算出B₀＝x₀，从而初始的MSP方程为x＝f(x，y，z|M₀，N₀，B₀)＝0y+0z+B₀＝x₀。

作为其他实施例，所述系数的初始值还可以为其它数值。也不限定经过某种先验知识求得的MSP方程作为初始。设定所述系数的初始值后，通过后续的优化对所述系数进行优化，以获得最终的优化系数。

执行步骤S3，更新当前系数，并获取对应的当前矢状面。若为第一次优化，则接收的系数为所述系数初始值(M₀、N₀、B₀)，对应的矢状面的方程为初始的矢状面方程x＝x₀。若为第一次以后的优化，则接收的系数为上次迭代结果作为当前更新的系数，并根据更新的系数获得当前矢状面。

执行步骤S4，读取未处理的像素点。设所述像素点坐标为p_i(x_i、y_i、z_i)，所述读取的内容包括所述像素点的位置坐标及对应的灰度值。其中，后续循环需要读尽所有像素。

执行步骤S5，像素点是否位于当前矢状面左侧，所述当前矢状面为根据步骤S3的系数确定的矢状面。本步骤为第一循环的判断条件，所述第一循环用于读取所述医学图像中位于同一部分的所有像素点。本实施例中，所述同一部分为左侧部分，作为其他实施例，所述同一部分为右侧部分。

其中，判断所述像素点是否位于所述当前矢状面的左侧，包括：计算获得所述像素点到所述矢状面的投射点坐标(x_t、y_t、z_t)；判断所述投射点与所述像素点的相对位置，获得所述像素点与所述矢状面的相对位置。如判断比较像素点和投射点对应的坐标大小，若所述投射点的坐标大于所述像素点坐标，则所述像素点位于所述投射点的左侧，即所述像素点位于所述矢状面的左侧；若所述投射点的坐标小于所述像素点坐标，则所述像素点位于所述投射点的右侧，即所述像素点位于所述矢状面的右侧。

进一步，还可能所述像素点位于所述矢状面上，则所述像素点坐标等于所述投射点的坐标。本发明采用图像配准的迭代的方法，通过矢状面将图像分成左右两部分，并分别将左右部分看作图像配准中的参考图像和浮动图像。通过对所述参考图像和浮动图像进行配准，并且比较两部分像素点的灰度均值差，以找到所述图像的中心矢状面。步骤S5及后续地步骤S6则为寻找位于矢状面两侧对称的像素点。

若所述步骤S5的判断结果为是，则执行步骤S6，以当前矢状面为对称面，计算像素点的右侧对称点。若所述步骤S5的判断结果为否，则执行步骤S4，继续读取未处理的像素点，并放弃上步骤读取的像素点信息。

具体地，包括读取所述像素p_i的灰度值，以所述当前矢状面为对称面，计算所述像素p_i的对称点的坐标。若所述对称点的坐标为整数，则将对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值；若所述对称点的坐标为非整数，则在所述医学图像中进行插值，以获得对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值。

具体地，计算所述像素p_i右侧镜像对称点R(p_i|M，N，B)的坐标为p′_xi，p′_yi，p′_zi，其中， $p_{\mathrm{xi}}^{\′}=p_{\mathrm{xi}}-\frac{2F}{\mathrm{KF}},$ $p_{\mathrm{yi}}^{\′}=p_{\mathrm{yi}}-\frac{2\mathrm{MF}}{\mathrm{KF}},$ $p_{\mathrm{zi}}^{\′}=p_{\mathrm{zi}}-\frac{2\mathrm{NF}}{\mathrm{KF}}$

式中，F＝p_xi-Mp_yi-Np_zi-B，KF＝1+M²+N²；p_xi，p_yi，p_zi为在MSP左侧图像L像素p_i的坐标。

执行步骤S7，累加计算均方差测度函数对系数的偏导数。其中，所述均方差测度函数为C，所述均方差测度函数C对M和N的偏导数分别为 $S=\frac{\∂C}{\∂M},$ $T=\frac{\∂C}{\∂N}.$

进一步地，均方差测度函数表达式为I(p_i)和I(R(p_i))为分别位于所述当前矢状面两侧互对称的像素点p_i和R(p_i)的灰度值，L为所述当前矢状面的任一部分像素点p_i集合。将上述两像素点的灰度值的差值进行平方后累加，获得所述均方差。所述均方差测度用于衡量以所述矢状面为对称面，互为镜像的两侧像素点的灰度的差值水平。

接着，按照上述的均方差测度函数，计算像素p_i对均方差对M和N的偏导数 的贡献，其中，均方差对M和N的偏导数(以下简称偏导数)的计算公式如下：

$\frac{\∂C}{\∂M}=\frac{2}{{\mathrm{KF}}^2}\underset{p_i\∈L}{\mathrm{\Σ}}{\mathrm{diff}}_i(G_{\mathrm{xi}}({\mathrm{dy}}_i\·\mathrm{KF}-2\mathrm{MF})+G_{\mathrm{yi}}((F+{\mathrm{Mdy}}_i)\mathrm{KF}-{2M}^{2}F)+G_{\mathrm{zi}}N\left(({\mathrm{dy}}_i\·\mathrm{KF}-2\mathrm{MF})\right))$

$\frac{\∂C}{\∂N}=\frac{2}{{\mathrm{KF}}^2}\underset{p_i\∈L}{\mathrm{\Σ}}{\mathrm{diff}}_i(G_{\mathrm{xi}}({\mathrm{dz}}_i\·\mathrm{KF}-2\mathrm{NF})+G_{\mathrm{yi}}M({\mathrm{dz}}_i\·\mathrm{KF}-2\mathrm{NF})+G_{\mathrm{zi}}((F+{\mathrm{Ndz}}_i)\mathrm{KF}-{2N}^2\·F))$

式中dy_i＝y₀-p_yi，dz_i＝z₀-p_zi；

G_xi，G_yi，G_zi为像素点p_i的右侧镜像点图像的梯度。进一步地，研究发现采用简单的一阶相邻点的均方差分作为梯度是最合适的。

按照上述步骤，获取当前系数的偏导数值。

具体地，执行完步骤S7后，则对应执行步骤S8的判断，是否有未处理的像素点。若是，则执行步骤S4，继续读取未处理的像素点。直至获取所有的矢状面的左侧的像素点的偏导数的累加值。本步骤S8为第二循环的判断条件，所述第二循环用于对所述医学图像中同一部分的像素点计算所述均方差测度函数的偏导数。

若步骤S8的判断结果为否，则执行步骤S9，根据梯度下降法迭代系数。具体地，当按上述步骤对全部像素遍历完毕，得到均方差对系数的偏导数，并根据梯度下降法更新系数。所述梯度下降法包括：设置初始修正步长，并设置松弛因子以降低所述修正步长，获得满足迭代结束的系数确定为系数。

若当前迭代为第1次迭代，初始系数为(M₀、N₀)，则采用所述梯度下降法计算获得的迭代系数M₁、N₁为 $M_1=M_0-\frac{{\mathrm{SD}}_1}{\sqrt{S^2+T^2}},$ $N_1=N_0-\frac{{\mathrm{TD}}_1}{\sqrt{S^2+T^2}},$ 其中，D₁为初始修正步长。所述初始修正步长的范围为0.01～0.03，优选地，所述初始修正步长为0.02。

若当前迭代为不为1的第i次迭代，上一次迭代后的系数为(M_i-1、N_i-1)，则采用所述梯度下降法计算获得的迭代系数为 其中，D_i为第i次优化的修正步长，所述均方差测度函数为C，所述均方差测度函数C对M和N的偏导数分别为

在首次优化后的循环优化中，所述修正步长按梯度变化向量的内积和与松弛因子计算，呈现几何级数下降。所述修正步长根据当前和上次优化获得的偏导数决定，具体地：若当前获取的偏导数为 上次优化获取的偏导数为 $S_{i-1}=\frac{\∂C}{\∂M},$ $T_{i-1}=\frac{\∂C}{\∂N}.$

若所述(S_i-S_i-1)*(T_i-T_i-1)＜0，则对上次迭代的修正步长进行下降后作为当前迭代的修正步长，具体地对所述修正步长乘以松弛因子，所述松弛因子的范围为(0.6～0.8)，优选地，所述松弛因子为0.6。

若所述(S_i-S_i-1)*(T_i-T_i-1)＞0，则不进行改变，继续采用上次优化使用的修正步长作为当前迭代的修正步长。

通过所述步骤S7计算获得一组优化后的(M_i、N_i、B_i)，执行步骤S10，判断当前迭代后是否满足迭代终止条件。若是，则执行步骤S11，固定当前系数对应的矢状面为中心矢状面。即迭代到此终止，输出当前迭代的系数对应的矢状面。本步骤S10为第三循环的判断条件，所述第三循环为采用不同修正步长对所述系数进行迭代。

所述第三循环的终止条件即步骤S10的迭代终止条件为提供步长阈值及最大迭代次数，当所述修正步长不大于所述步长阈值，或者所述迭代次数达到所述最大迭代次数，所述迭代终止，并提取当前迭代后的系数对应的矢状面为中心矢状面。

进一步地，所述初始修正步长的范围为0.01～0.03，所述松弛因子的范围为0.6～0.8，所述步长阈值的范围为所述初始修正步长的0.5％～1.5％，所述最大迭代次数的范围为50～150。

较佳地，所述初始修正步长为0.02，所述步长阈值为所述初始修正步长的1％，所述最大迭代次数为100次，所述松弛因子为0.6。

本发明还可以通过等距离或者多分辨率采样的方法来使所述迭代的收敛加速，此处不详述。

若步骤S10的判断结果为否，即所述修正步长大于所述步长阈值或者所述迭代次数未达到最大迭代次数，则执行步骤S3，更新当前系数值，并获取当前矢状面。其中，所述当前矢状面为根据所述当前系数设置。

与现有技术相比，本发明具有以下优点：

因为中心矢状面的代价函数一般为非线性函数，采用求解线性规划问题的单纯形法求解所述代价函数并不合适，本发明则采用适用于求解非线性的梯度下降法对代价函数，即均方差测度函数进行迭代，并使用迭代终止时的系数对应的矢状面作为中心矢状面，提高中心矢状面提取的精准度；

进一步地，利用中心矢状面必然通过医学图像的质心的特性作为约束条件，减少需要迭代的系数，降低迭代计算量，加速迭代的速度；

最后，所述均方差测度函数借鉴图像配准中的迭代优化方法，把中心矢状面图像分为左右两部分，并将所述左右两部分视为参考图像和浮动图像，通过比较两部分的灰度的均方差实现配准，对所述均方差测度函数进行迭代，获取待定系数及对应的中心矢状面，提高中心矢状面提取的精准度。

本发明虽然已以较佳实施例公开如上，但其并不是用来限定本发明，任何本领域技术人员在不脱离本发明的精神和范围内，都可以利用上述揭示的方法和技术内容对本发明技术方案做出可能的变动和修改，因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化及修饰，均属于本发明技术方案的保护范围。

Claims (15)

1.一种中心矢状面的提取方法，其特征在于，包括：

提供医学图像及矢状面方程，所述矢状面方程包括至少一个待定系数，所述待定系数通过梯度下降法迭代确定，进而确定中心矢状面方程，实现中心矢状面的提取；

其中对任一次所述迭代包括：获取当前系数，通过所述当前系数对应的当前矢状面将所述医学图像分为两部分；获取所述两部分医学图像灰度的均方差测度函数，计算所述均方差测度函数对所述当前系数的偏导数，利用所述偏导数进行迭代，迭代后的系数作为下一次迭代的当前系数。

2.如权利要求1所述的中心矢状面的提取方法，其特征在于，对于第一次迭代，通过对所述系数进行初始化，以获得所述当前系数；对于第一次后的迭代，将上一次迭代后获得的系数作为当前系数。

3.如权利要求1所述的中心矢状面的提取方法，其特征在于，所述梯度下降法包括如下过程：提供初始的修正步长，并通过设置松弛因子降低所述修正步长，以迭代所述修正步长。

4.如权利要求3所述的中心矢状面的提取方法，其特征在于，还包括：提供步长阈值和最大迭代次数，当所述修正步长不大于所述步长阈值，或者所述迭代次数达到所述最大迭代次数，所述迭代终止，并提取当前迭代后的系数对应的矢状面为中心矢状面。

5.如权利要求4所述的中心矢状面的提取方法，其特征在于，所述初始的修正步长的范围为0.01～0.03，所述松弛因子的范围为0.6～0.8，所述步长阈值的范围为所述初始修正步长的0.5％～1.5％，所述最大迭代次数的范围为50～150次。

6.如权利要求5所述的中心矢状面的提取方法，其特征在于，所述初始修正步长为0.02，所述松弛因子为0.6，所述步长阈值为所述初始修正步长的1％，所述最大迭代次数为100次。

7.如权利要求1所述的中心矢状面的提取方法，其特征在于，所述均方差测度函数为其中I(p_i)和I(R(p_i))为分别位于所述当前矢状面两侧互对称的像素点p_i和R(p_i)的灰度值，L为所述当前矢状面的任一部分的像素点p_i集合。

8.如权利要求7所述的中心矢状面的提取方法，其特征在于，包括：读取所述像素p_i的灰度值，以所述当前矢状面为对称面，计算所述像素p_i的对称点的坐标作为R(p_i)的坐标。

9.如权利要求8所述的中心矢状面的提取方法，其特征在于，若所述对称点的坐标为整数，则将对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值；若所述对称点的坐标为非整数，则在所述医学图像中进行插值，以获得对应坐标的灰度值作为所述像素点R(p_i)对应的灰度值。

10.如权利要求1所述的中心矢状面的提取方法，其特征在于，所述矢状面方程为x＝My+Nz+B，其中，所述M、N及B为所述待定的系数，其中，第i次迭代后的系数对应为M_i、N_i及B_i，第一次迭代的当前系数为M₀、N₀及B₀。

11.如权利要求10所述的中心矢状面的提取方法，其特征在于，还包括获取所述医学图像的质心坐标，并根据确定的系数M、N及所述质心坐标，获取待定系数B。

12.如权利要求11所述的中心矢状面的提取方法，其特征在于，若所述质心坐标为(x₀、y₀、z₀)，，所述系数的初始值为M₀＝0，N₀＝0，对应的初始矢状面方程为x＝B₀＝x₀。

13.如权利要求10所述的中心矢状面的提取方法，其特征在于，所述偏导数为所述均方差测度函数C分别对系数M、N的偏导数，分别记为S、T，表达式为： $S=\frac{\∂C}{\∂M},$ $T=\frac{\∂C}{\∂N}.$

14.如权利要求13所述的中心矢状面的提取方法，其特征在于，包括：第i次迭代后的系数为 $M_i=M_{i-1}-\frac{{\mathrm{SD}}_i}{\sqrt{S^2+T^2}},$ $N_i=N_{i-1}-\frac{{\mathrm{TD}}_i}{\sqrt{S^2+T^2}};$ 其中，i为从1开始的整数，D_i为第i次优化的修正步长。

15.如权利要求13所述的中心矢状面的提取方法，其特征在于，包括：i大于1的第i次迭代的当前迭代偏导数为S_i，T_i，上次迭代偏导数为S_i-1，T_i-1，

若所述(S_i-S_i-1)*(T_i-T_i-1)＜0，则对上次迭代的修正步长进行下降后作为当前迭代的修正步长；

若所述(S_i-S_i-1)*(T_i-T_i-1)＞0，则继续采用上次迭代使用的修正步长作为当前迭代的修正步长。