CN100386052C - 一种获得声门下压力值及计算发声效率的方法 - Google Patents

Abstract

本发明公开了一种基于声带的三维形状重建和发声的空气动力学建模技术的，获得声门下压力值与计算发声效率及其它空气动力学参数值的方法。具体包括：提出了一种定义声带三维几何形状特征值的方法；在此基础上给出了声带表面三维形状重建的方法；同时提出了一种基于三维有限元算法的发声的空气动力学建模方法；并进而给出了获得声门下压力值、计算发声效率及其它空气动力学参数值的方法。本发明的方法对于单侧或双侧声带麻痹、息肉、杓状肌固定、厚边、组织水肿、上皮损伤、喉癌等临床喉病病症的辅助诊断及治疗均有广泛的应用前景，同时为发声功能评价、发声生理和基础研究，以及艺术嗓音、语音声学、语音信号处理和语言学等学科都有重要意义。

Description

一种获得声门下压力值及计算发声效率的方法

技术领域

本发明属于正常及病变语音的空气动力学参量及发声特性参量技术领域，特别涉及一种获得声门下压力值及发声效率的方法。

背景技术

语音是人类进行交流最基本、最有效、最重要的手段。在语音的产生过程中，声带是最重要的发音器官，但它同时也是最易受损伤和发生病变的部位。声带在发声过程中呈现一种复杂的三维振动，这种三维振动模式通过语音的空气动力学参量和发声声学特性参量表现出来，在根本上决定了语音产生的效率和质量。因此，对于发声过程中的空气动力学参量和发声声学特性参量进行研究分析，不仅对于临床上喉病检查、喉功能评价、发声生理和基础研究有重要意义，而且对艺术嗓音、语音声学、语音信号处理和语言学等学科领域都有着深远的影响。

在语音产生的空气动力学参量中，声门下压力值由于直接控制和影响声门阻抗、声门腔内气流流量、气流状态、声门上压力、气流分离点等重要参量，而成为语音产生过程中的一个非常重要的空气动力学参量。另一个语音产生过程中的重要参数就是发声效率，它是衡量人类发声过程中空气动力学能量向语音声学能量转化能力的重要参数。在所有有关声带的基础研究和实际使用中，以及对病理语音信号进行分析处理和临床辅助诊断中，这两个参数都是十分基础而且重要的。

目前，国内外通常采用两类方法来获得声门下压力值和估算发声效率。一类是侵入性的，甚至是有创的方法来获取声门下压。如采用跨声门压导管和食管气球的方法进行侵入性测量声门下压(VandenBerg，1956)，或利用脊椎麻醉针穿过气管测量声门下压(Isshiki，1964)，另一类是非侵入性的，如利用呼吸气流面罩和压力传感器通过测量口腔内压来估算发音时的声门下压(Ronthenberg，1973)，应用反滤波或者无反射管测量出气流的交直流分量，用它们的比率来估算发声效率和声门下压(Isshiki，1981)，用“波流指数”法测试声门下压和发声效率(Kakita，1986)。前一类方法由于是侵入性和有创的，给测试对象带来许多不便，而后一类方法则是估算值，其精确性受到测量操作及计算方法的局限。

发明内容

针对现有技术存在的缺陷或不足，本发明的目的在于，基于申请人对语音的产生及其空气动力学建模研究已取得的成果，提出一种非侵入的，利用X线断层图、口腔语音气流压力信号、图像边缘提取、声带表面形状三维建模和三维有限元计算方法等多种技术获得较精确声门下压力值及发声效率的方法。

本发明主要面向喉病检查、喉功能评价和发声生理研究，同时兼顾艺术嗓音、语音声学、语音信号处理和语言学等学科领域的需要。

为了实现上述任务，本发明采取如下的技术解决方案：

一种获得声门下压力值及计算发声效率的方法，其特征在于，包括下列步骤：

1)首先对声带表面三维形状特征用一组形状参数来定义声带运动过程中处于各种状态的形状和参数，图2给出了这些形状参数的具体意义示意图；

2)利用图像的边缘检测和提取技术对人体喉部的X线的横向断层及纵向断层图进行边缘检测和提取，将提取出来的声带边缘曲线与声带运动过程中的参数进行拟合，得到最小声门直径、声门入口曲线半径、声门出口曲线半径、声门上表面角、声门下表面角、声门倾角、声门长度、声门上腔长度、声门下腔长度和声门厚度参数；

3)将步骤2)所获得的各项声门形状参数，按照步骤1)进行声带表面三维形状重建，由发声对象发/a:/音时测得的语音压力信号和气流流量值作为边界条件，采用补偿有限元算法进行空气动力学建模，即可获得声门下压力值及其它空气动力学参数值；

假设气流均质、无重力、且其流动为定常、不可压缩，则三维Navier-Stokes方程表达为：

$\frac{\∂v_x}{\∂t}+v_x\frac{\∂v_x}{\∂x}+v_y\frac{\∂v_x}{\∂y}+v_z\frac{\text{\∂}{v}_x}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂x}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_x}{\∂x^2}+\frac{{\∂}^2{v}_x}{\∂y^2}+\frac{{\∂}^2{v}_x}{\∂z^2}]=0---\left(1\right)$

$\frac{\∂v_y}{\∂t}+v_x\frac{\∂v_y}{\∂x}+v_y\frac{\∂v_y}{\∂y}+v_z\frac{\∂v_y}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂y}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_y}{\∂x^2}+\frac{{\∂}^2{v}_y}{\∂y^2}+\frac{{\∂}^2{v}_y}{\∂z^2}]=0---\left(2\right)$

$\frac{\∂v_z}{\∂t}+v_x\frac{\∂v_z}{\∂x}+v_y\frac{\∂v_z}{\∂y}+v_z\frac{\∂v_z}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂z}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_z}{\∂x^2}+\frac{{\∂}^2{v}_z}{\∂y^2}+\frac{{\∂}^2{v}_z}{\∂z^2}]=0---\left(3\right)$

以及

$\frac{\∂v_x}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=0---\left(4\right)$

其中v_x，v_y，v_z分别为x，y和z方向的速度分量，ρ是空气密度，p是空气压力，μ为空气的分子粘度系数；

在补偿有限体积方法中，连续性方程可表达为：

$\frac{{\∂v}_x}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=-\frac{1}{\mathrm{\λ}}p---\left(5\right)$

其中λ是补偿数，当补偿数趋于无穷大时，满足质量守恒定律；

有限体积方程是由Galerkin方法构造的，速度元使用四阶的形状函数进行内插值计算，压力元则采用四面体单元上定义的线性形状函数来内插替换计算，这种压力内插函数能产生高精度的快速收敛结果；在整个气流域中，采用十节点的速度元和四节点的压力元；

整体有限体积方程采用直接迭代法进行求解，每次迭代的解都用亚松驰方法校正以确保迭代的稳定性，直到标准误差比收敛标准小时，整体有限体积方程的迭代就完成，标准误差定义为：

$e_i=\underset{j=1,N}{\max }\left[\right|\frac{a_{i,j}^n-a_{i,j}^{n-1}}{a_{i,\max }^n}\left|\right]---\left(6\right)$

其中下标i(i＝v_x，v_y，v_z，或者p)表示流体变量的每个分量，n是迭代级数，a_i，j表示在第j个节点上的第i个气流变量的节点值，N表示节点的总数；

4)利用步骤3)所获得的声门下压力和通过呼吸面罩获得的口腔气流流量和口腔语音压力的直流分量值来计算获得发声效率：

发声效率通过下式计算：

$\mathrm{VE}=\frac{P_r}{P_a}=\frac{P_s{A}_m}{P_{\mathrm{sb}}U}---\left(7\right)$

其中：P_r为唇端辐射的声学信号功率；

P_a为声门下空气动力功率；

P_s为测量呼吸面罩内麦克风得到的声压信号；

A_m为测量呼吸面罩的内表面面积；

P_sb为声门下平均压力；

U为声门气流体积速度波平均值。

本发明首次提出了一种非侵入的，利用人体喉部的横向与纵向X线断层图和从人体口腔面罩所获得的语音气流压力信号，结合图像的边缘检测和提取、声带的三维形状特征值定义方法、声带的三维形状重建方案和发声的空气动力学建模方法，来获得声门下压力值、发声效率值及其它空气动力学参数值的方法。经过对临床喉病患者的检测和计算机编程计算，不同喉病病人的声门下压力值、发声效率值及其它空气动力学参数值差异较大，经与临床治疗效果对比，本发明所提出的方法对于单侧或双侧声带麻痹、息肉、杓状肌固定、厚边、组织水肿、上皮损伤、喉癌等临床喉病病症的辅助诊断及治疗均有广泛的应用前景，同时为发声功能评价提供了有效的方法。

附图说明

图1为本发明的整体结构框图。

图2为声带表面形状特征参数定义方法。

图3为一个10°声门倾角条件下声带表面三维形状重建的示例。

图4为声门角分别取-40°，-20°，-10°，-5°，0°，5°，10°，20°，40°时，9个重建的声门形状及其它参数。

以下结合附图对本发明人给出的实施方式作进一步的详细说明。

具体实施方式

采用图1所示流程图，用本方法获得声门下压力值、发声效率值及其它空气动力学参数值的步骤如下：

1)定义声带表面三维形状特征值

用一组形状参数来定义声带运动过程中处于各种状态的形状，为声带形状的三维重建打下基础。

图2给出了这些形状参数的具体意义示意图，这些形状参数包括：

①最小声门直径。指声门腔内声门最窄处的间距值，即图2中D_min所示的长度。

②声门入口曲线半径。指声门入口处弧线的曲率半径，即图2中R_in所示的半径。

③声门出口曲线半径。指声门出口处弧线的曲率半径，即图2中R_out所示的半径。

④声门上表面角。指声门上表面切线斜率与水平线之间的夹角，即图2中θ_sup所示的角度。

⑤声门下表面角。指声门下表面切线斜率与水平线之间的夹角，即图2中θ_sub所示的角度。

⑥声门倾角。指声带两个表面某一瞬时时刻的夹角，即图2中T所示的长度。

⑦声门长度。指声带入口曲线与出口曲线之间的直线的长度，即图2中T所示的长度。

⑧声门厚度。指声带上沿与下沿之间的长度。

⑨声门上腔长度。指上声道的长度，即图2中Q₁所示的长度。

⑩声门下腔长度。指下声道的长度，即图2中Q₂所示的长度。

2)声带表面三维形状重建

利用图像的边缘检测和提取技术对人体喉部的X线横向及纵向断层图进行边缘检测和提取，将提取出来的声带边缘曲线与步骤1)所述的形状参数进行拟合，从而重建声带表面的三维形状。

重建声带表面三维形状采用如下步骤即可：

(1)获取研究对象喉部横向与纵向的X线断层图，同时同步采集研究对象的口腔气流流量信号和口腔语音压力信号。

在医生的协助下，在研究对象以胸声区发声方式发元音/a:/时，获取研究对象喉部横向与纵向的X线断层图，一般应分别采集至少十张断层图以获得较为准确的声门形状参数，同时利用口腔面罩同步采集研究对象的口腔气流流量信号和口腔语音压力信号。

(2)对获得的一帧帧喉部横向断层图图像，利用蛇形算法(Karim，1998)检测声带边缘曲线。

(3)将采用步骤(1)所获得的声带横向表面边缘曲线与声门形状参数进行拟合，得到如下参数：①最小声门直径；②声门入口曲线半径；③声门出口曲线半径；④声门上表面角；⑤声门下表面角；⑥声门倾角；⑦声门长度；⑧声门上腔长度；⑨声门下腔长度；

(4)对获得的一帧帧喉部纵向断层图图像，利用蛇形算法(Karim，1998)检测声带边缘曲线。

(5)将获得的声带纵向表面边缘曲线与声门形状参数进行拟合，得到声带的厚度参数。

(6)利用获得的声门形状参数和声带的三维形状特征值定义方法，对声带表面三维形状进行重建。

3)将所重建的声带表面三维形状中的声门角参数进行变化，分别取-40°，-20°，-10°，-5°，0°，5°，10°，20°，40°，其它参数不变，得到9个重建的声门形状。图4给出了这9个重建声带的形状及其它参数，在研究过程中假设声门形状完全对称，因此只给出了一半的形状参数。

4)以口腔气流流量和口腔语音压力信号的直流分量为边界条件，利用发声的三维有限元建模方法，对于重建的9个声带表面三维形状分别计算声门下压力值，及其它空气动力学参数值。将计算所得的各参数值取平均后得到的就是希望获得的参数值。

声门下压力值的计算基于三维有限元算法，其方法是：

由于人类发声过程中Reynolds数较低(一般不大于2000)，因此可认为气流均质、无重力、且其流动为定常、不可压缩，则三维Navier-Stokes方程可表达为：

$\frac{\∂v_x}{\∂t}+v_x\frac{\∂v_x}{\∂x}+v_y\frac{\∂v_x}{\∂y}+v_z\frac{\text{\∂}{v}_x}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂x}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_x}{\∂x^2}+\frac{{\∂}^2{v}_x}{\∂y^2}+\frac{{\∂}^2{v}_x}{\∂z^2}]=0---\left(1\right)$

$\frac{\∂v_y}{\∂t}+v_x\frac{\∂v_y}{\∂x}+v_y\frac{\∂v_y}{\∂y}+v_z\frac{\∂v_y}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂y}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_y}{\∂x^2}+\frac{{\∂}^2{v}_y}{\∂y^2}+\frac{{\∂}^2{v}_y}{\∂z^2}]=0---\left(2\right)$

$\frac{\∂v_z}{\∂t}+v_x\frac{\∂v_z}{\∂x}+v_y\frac{\∂v_z}{\∂y}+v_z\frac{\∂v_z}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂z}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_z}{\∂x^2}+\frac{{\∂}^2{v}_z}{\∂y^2}+\frac{{\∂}^2{v}_z}{\∂z^2}]=0---\left(3\right)$

以及

$\frac{\∂v_x}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=0---\left(4\right)$

其中v_x，v_y，v_z分别为x，y和z方向的速度分量，，ρ是空气密度，p是空气压力，μ为空气的分子粘度系数。

在补偿有限体积方法中，连续性方程(4)可表达为：

$\frac{{\∂v}_x}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=-\frac{1}{\mathrm{\λ}}p---\left(5\right)$

其中λ是补偿数.当补偿数趋于无穷大时，满足质量守恒定律。

有限体积方程是由Galerkin方法构造的。速度元使用四阶的形状函数进行内插值计算，压力元则采用四面体单元上定义的线性形状函数来内插替换计算，这种压力内插函数能产生高精度的快速收敛结果。在整个气流域中，我们采用了十节点的速度元和四节点的压力元。

整体有限体积方程采用直接迭代法进行求解。每次迭代的解都用亚松驰方法校正以确保迭代的稳定性。直到标准误差比收敛标准小时，整体有限体积方程的迭代就完成。标准误差定义为：

$e_i=\underset{j=1,N}{\max }\left[\right|\frac{a_{i,j}^n-a_{i,j}^{n-1}}{a_{i,\max }^n}\left|\right]---\left(6\right)$

其中下标i(i＝v_x，v_y，v_z，或者p)表示流体变量的每个分量，n是迭代级数，a_i，j表示在第j个节点上的第i个气流变量的节点值，N表示节点的总数。

5)利用获得的口腔气流流量和口腔语音压力信号的直流分量，和获得的声门下压力值，通过下式计算发声效率：

$\mathrm{VE}=\frac{P_r}{P_a}=\frac{P_s{A}_m}{P_{\mathrm{sb}}U}---\left(7\right)$

其中：P_r为唇端辐射的声学信号功率；

P_a为声门下空气动力功率；

P_s为测量呼吸面罩内麦克风得到的声压信号；

A_m为测量呼吸面罩的内表面面积；

P_sb为声门下平均压力；

U为声门气流体积速度波平均值；

6)利用获得的其它空气动力学参数值进行发声声学特性分析。

采用本方法，除了能计算出声门下压力值以外，还能计算出以下其它空气动力学参数值，主要有：

①声门阻抗，由跨声门压力与相应的气流流量值相比得到；

②跨声门压力，由声门下压力值与声门上压力值相减得到；

③压力极小值，由声门腔内压力分布场值得到；

④速度极大值，由声门腔内速度分布场值得到；

⑤气流分离点，由声门腔内压力分布场得到；

⑥上声带表面压力与下声带表面压力比率，由声门腔内压力分布场值得到；

⑦压力极低值与跨声门压的比率，由压力极小值与跨声门压值相比得到；

⑧声门上压与声门下压的比率；由声门上压值与声门下压值相比得到。

发声的声学特性主要包括声强、基频和声调等参量，发声的空气动力学参数与发声参量间存在着密切的、非常重要的关系。利用本发明的研究成果还可以定性分析发声的声学特性，主要包括：

①跨声门压的增加会使发声强度增加；

②跨声门压的增加会使发声基频同时增加；

③跨声门压的增加同时会使发声时的声调变高，这是因为声带需要更高的跨声门压来维持同样幅度的振动；

④声门阻抗的增加会降低发声的声强；

另外对于本发明所采用的方法、技术手段和结果分析作以下说明：

1.本发明叙述获得研究对象喉部的横向与纵向断层图的手段时主要提出的是利用X射线断层成像技术(CT)，但本发明的使用范围并不局限于此，本发明对于磁共振成像(MRI)，高能聚焦超声(HIFU)成像，正电子发射断层成像(PET)等其它断层成像方法也有较高的适应性。

2.本发明叙述图像的边缘检测和提取技术时提出了蛇形算法(snake)(Karim，1998)，主要是因为这种方法具有动态灵活、计算量小的特点，但对于图像的边缘检测和提取并不局限于本方法，其它众多有效的图像边缘检测和提取方法也都适用于本发明。

3.本发明所述方法的结果也并不局限于获取声门下压力值、发声效率值及其它空气动力学参数值，正如上述步骤6)所述，发声空气动力学参数与发声声学特性有密切关系，对它们之间关系的深入研究将揭示出更多的研究对象发声特性方面的成果。

本发明不仅对于临床喉病检查、喉功能评价有现实指导意义，而且对发声生理和基础研究，以及艺术嗓音、语音声学、语音信号处理和语言学等学科领域都具有重要意义。更重要的是，它为进一步对于喉部病理状态下语音产生机理的研究、分析、建模，以及病理信号的采集、分析、分类、处理，尤其是由于声门几何形状变化所引起的语音畸变研究打下了良好基础，它同时对嗓音康复、语音合成、语音处理、语音识别的发展，都有着深远的影响。

以下是发明人给出的实施例，但并不局限于这些实施例。

实施例1：男性正常轻声发声时声门下压力值及发声效率获得方法

实施例的研究对象为一名25岁的青年男性，采用如下步骤得到声门下压力、发声效率值及其它空气动力学参数值：

1)在医生的协助下，在研究对象以胸声区发声方式轻声发元音/a:/时，分别采集十张研究对象喉部横向与纵向的X线断层图，同时利用口腔面罩同步采集研究对象的口腔气流流量信号和口腔语音压力信号，所采集到的口腔气流流量为189.4cm³/s，口腔语音压力为0.51Pa.

2)对获得的喉部横向与纵向断层图图像，利用蛇形算法(Karim，1998)检测声带边缘曲线。

3)将获得的声带横向与纵向表面边缘曲线与声门形状参数进行拟合，综合十个断层图像得到以下参数值：①平均最小声门直径＝0.06cm；②声门入口曲线半径＝0.15cm；③声门出口曲线半径＝0.0987cm；④声门上表面角＝0°；⑤声门下表面角＝40°；⑥声门倾角＝10°；⑦声门长度＝0.3cm；⑧声门上腔长度＝0.6cm；⑨声门下腔长度＝1.5cm；⑩声门厚度＝1.2cm。

4)利用获得的声门形状参数和声带的三维形状特征值定义方法，对声带表面三维形状进行重建。图3示出了声门角为10°时，声带的三维形状重建图。

5)将所重建的声带表面三维形状中的声门角参数进行变化，分别取-40°，-20°，-10°，-5°，0°，5°，10°，20°，40°，其它参数不变，得到9个重建的声门形状。图4示出了这9个重建声带的形状及其它参数，假设声门形状完全对称，因此只给出了一半的形状参数。

6)以步骤1)获得的口腔气流流量(189.4cm³/s)和和声门重建形状为边界条件，设模型出口压力为零，利用权利要求1第3)步所述发声的三维有限元建模方法，对于步骤5)重建的9个声带表面三维形状分别计算声门下压力值，得到的9个声门下压力值依次为：7.86，7.87，8.21，8.53，9.81，7.49，6.93，7.55，7.68cmH₂O，其均值7.99cmH₂O即为所求的声门下压力值。

7)利用步骤1)获得的口腔气流流量和口腔语音压力信号，和步骤6)获得的声门下压力值，利用权利要求1中的公式(7)计算发声效率。其中，呼吸面罩内麦克风得到的声压信号为0.51Pa，气流量为189.4cm³/s，呼吸面罩的内表面面积为100cm³，声门下平均压力由步骤6)得到为7.99cmH₂O，由公式可得该男性轻声发声时的发声效率为3.37*10^-4。

8)利用得到的声门下压力值和气流压力速度场，结合权利要求2所述的方法，计算得其它参量值。本例中：①声门阻抗值为4.22gcm^-4s^-1；②跨声门压力为7.99cmH₂O；③压力极小值为-866.08Pa；④速度极大值为7.097m/s；⑤气流分离点在最小声门直径下游0.08cm处；⑥上声带表面压力与下声带表面压力比率为1(由于声带及其压力分布完全对称)；⑦压力极低值与跨声门压的比率为1.08；⑧声门上压与声门下压的比率在出口压力设为大气压的情况下为零。

实施例2：男性正常大声发声时声门下压力值及发声效率获得方法

采用与实施例1相同的研究对象和研究步骤，在该研究对象以胸声区发声方式大声发元音/a:/时，所采集到的口腔气流流量为502.8cm³/s，口腔语音压力为2.6Pa。又对其发声时的X线断层图进行分析，提取的参数值与轻声发声时相同，故声带的三维形状图相同，将声门角参数进行与实施例1相同的变化后计算得到9个声门下压力值依次为：56.62，48.22，45.55，44.94，43.9，43.01，42.55，43.38，44.18cmH₂O，其均值45.87cmH₂O，即为所求的声门下压力值。采用同样方法可计算得到该男性大声发声时的发声效率为1.127*10^-4。再利用得到的声门下压力值和气流压力速度场，得到其它参量值如下：①声门阻抗值为9.12gcm^-4s^-1；②跨声门压力为45.87cmH₂O；③压力极小值为-52.12cmH₂O；④速度极大值为15.3m/s；⑤气流分离点在最小声门直径下游0.10cm处；⑥上声带表面压力与下声带表面压力比率为1(由于声带及其压力分布完全对称)；⑦压力极低值与跨声门压的比率为1.14；⑧声门上压与声门下压的比率在出口压力设为大气压的情况下为零。

Claims (2)

1.一种获得声门下压力值及计算发声效率的方法，其特征在于，包括下列步骤：

1)首先对声带表面三维形状特征用一组形状参数来定义声带运动过程中处于各种状态的形状和参数；

2)利用图像的边缘检测和提取技术对人体喉部的X线的横向断层及纵向断层图进行边缘检测和提取，将提取出来的声带边缘曲线与声带运动过程中的参数进行拟合，得到最小声门直径、声门入口曲线半径、声门出口曲线半径、声门上表面角、声门下表面角、声门倾角、声门长度、声门上腔长度、声门下腔长度和声门厚度参数；

3)将步骤2)所获得的各项声门形状参数，按照步骤1)进行声带表面三维形状重建，由发声对象发/a:/音时测得的语音压力信号和气流流量值作为边界条件，采用补偿有限元算法进行空气动力学建模，即可获得出声门下压力值及其它空气动力学参数值；

假设气流均质、无重力、且其流动为定常、不可压缩，则三维Navier-Stokes方程表达为：

$\frac{\∂v_x}{\∂t}+v_x\frac{\∂v_x}{\∂x}+v_y\frac{\∂v_x}{\∂y}+v_z\frac{\∂v_x}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂x}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_x}{\∂x^2}+\frac{{\∂}^2{v}_x}{\∂y^2}+\frac{{\∂}^2{v}_x}{\∂z^2}]=0---\left(1\right)$

$\frac{\∂v_y}{\∂t}+v_x\frac{\∂v_y}{\∂x}+v_y\frac{\∂v_y}{\∂y}+v_z\frac{\∂v_y}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂y}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_y}{\∂x^2}+\frac{{\∂}^2{v}_y}{\∂y^2}+\frac{{\∂}^2{v}_y}{\∂z^2}]=0---\left(2\right)$

$\frac{\∂v_z}{\∂t}+v_x\frac{\∂v_z}{\∂x}+v_y\frac{\∂v_z}{\∂y}+v_z\frac{\∂v_z}{\∂z}+\frac{1}{\mathrm{\ρ}}\frac{\∂p}{\∂z}-\frac{\mathrm{\μ}}{\mathrm{\ρ}}[\frac{{\∂}^2{v}_z}{\∂x^2}+\frac{{\∂}^2{v}_z}{\∂y^2}+\frac{{\∂}^2{v}_z}{\∂z^2}]=0---\left(3\right)$

以及

$\frac{\∂v_z}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=0---\left(4\right)$

其中v_x，v_y，v_z分别为x，y和z方向的速度分量，ρ是空气密度，p是空气压力，μ为空气的分子粘度系数；

在补偿有限体积方法中，连续性方程可表达为：

$\frac{\∂v_x}{\∂x}+\frac{\∂v_y}{\∂y}+\frac{\∂v_z}{\∂z}=-\frac{1}{\mathrm{\λ}}p---\left(5\right)$

其中λ是补偿数，当补偿数趋于无穷大时，满足质量守恒定律；

有限体积方程是由Galerkin方法构造的，速度元使用四阶的形状函数进行内插值计算，压力元则采用四面体单元上定义的线性形状函数来内插替换计算，这种压力内插函数能产生高精度的快速收敛结果；在整个气流域中，采用十节点的速度元和四节点的压力元；

整体有限体积方程采用直接迭代法进行求解，每次迭代的解都用亚松驰方法校正以确保迭代的稳定性，直到标准误差比收敛标准小时，整体有限体积方程的迭代就完成，标准误差定义为：

$e_i=\begin{array}{c}\mathrm{Max}\\ j=1,N\end{array}\left[\right|\frac{a_{i,j}^n-a_{i,j}^{n-1}}{a_{i,\max }^n}\left|\right]---\left(6\right)$

其中下标i＝v_x，v_y，v_z或者p，表示流体变量的每个分量，n是迭代级数，α_i，j表示在第j个节点上的第i个气流变量的节点值，N表示节点的总数；

4)利用步骤3)所获得的声门下压力和通过呼吸面罩获得的口腔气流流量和口腔语音压力的直流分量值来计算获得发声效率：

发声效率通过下式计算：

$\mathrm{VE}=\frac{P_r}{P_a}=\frac{P_s{A}_m}{P_{\mathrm{sb}}U}---\left(7\right)$

其中：P_r为唇端辐射的声学信号功率；

P_a为声门下空气动力功率；

P_s为测量呼吸面罩内麦克风得到的声压信号；

A_m为测量呼吸面罩的内表面面积；

P_sb为声门下平均压力；

U为声门气流体积速度波平均值。

2.如权利要求1所述的方法，其特征在于，所述的其它空气动力学参数值，主要有：

①声门阻抗，由跨声门压力与相应的气流流量值相比得到；

②跨声门压力，由声门下压力值与声门上压力值相减得到；

③压力极小值，由声门腔内压力分布场值得到；

④速度极大值，由声门腔内速度分布场值得到；

⑤气流分离点，由声门腔内压力分布场得到；

⑥上声带表面压力与下声带表面压力比率，由声门腔内压力分布场值得到；

⑦压力极低值与跨声门压的比率，由压力极小值与跨声门压值相比得到；

⑧声门上压与声门下压的比率；由声门上压值与声门下压值相比得到。

Images