CN106909724B - 基于调和场的各向异性特定患者巩膜有限元建模方法 - Google Patents
基于调和场的各向异性特定患者巩膜有限元建模方法 Download PDFInfo
- Publication number
- CN106909724B CN106909724B CN201710085497.5A CN201710085497A CN106909724B CN 106909724 B CN106909724 B CN 106909724B CN 201710085497 A CN201710085497 A CN 201710085497A CN 106909724 B CN106909724 B CN 106909724B
- Authority
- CN
- China
- Prior art keywords
- sclera
- field
- harmonic
- anisotropic
- gradient
- 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.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics
- G16H50/50—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics for simulation or modelling of medical disorders
Landscapes
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Physics & Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Public Health (AREA)
- Epidemiology (AREA)
- Pathology (AREA)
- Databases & Information Systems (AREA)
- Data Mining & Analysis (AREA)
- General Health & Medical Sciences (AREA)
- Primary Health Care (AREA)
- Biomedical Technology (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Prostheses (AREA)
Abstract
本发明公开一种基于调和场的各向异性特定患者巩膜有限元建模方法,包括以下步骤:采用高精度激光扫描得到特定患者巩膜立体模型;对所述立体模型指派巩膜权值、设计伞算子,得到调和场;为调和场设置约束条件,分布调和场的梯度场和等值线场;通过梯度场和等值线场平滑地分布规则六面体网格在非规则巩膜上;对立体模型设置基于调和场的各向异性材料参数;对比实验结果。本发明的基于调和场的各向异性特定患者巩膜有限元建模方法可以为特定患者巩膜进行有效建模,加入了各向异性材质方向场,不但保证了建模精确质量高,而且可以对非规则半球巩膜建模进行个性化参数设置。
Description
技术领域
本发明涉及三维有限元建模和医学成像的交叉技术领域,具体地,涉及一种基于调和场的各向异性特定患者巩膜有限元建模方法。
背景技术
最近,研究人员专注于理解巩膜生物力学的特性,因为巩膜和筛板的生物力学特性决定了视神经乳头的生物力学变化,这种特性在处理由于眼内压增大所导致的视网膜神经细胞缺失和视神经损伤上扮演着着重的角色。
眼球后部包括三层:巩膜、脉络膜、视网膜,其中巩膜最厚,视网膜最薄。对巩膜、脉络膜、视网膜施加同样的压力,它们的切线模量在不同数量级上,其中巩膜最高。因此,巩膜在维持眼球形状上扮演着至关重要的角色。巩膜包裹着眼球,它由纤维组织组成,由几乎完全平行且交错的稠密带状原蛋白组成,这些原蛋维持着巩膜的生物力学特征。研究发现,大部分动物巩膜组织的许多相似的特征在结构上都是各向异性的。在后部和视神经乳头周围区域,巩膜纤维几乎是呈环状面的,但不与前部和赤道部区域对齐。环状面巩膜纤维可能充当一个增强环去防止视神经乳头的变形。巩膜胶原纤维的生物力学特性可展示各向异性对外力的反应。
一些研究采用规则的映射六面体网格去生成相同厚度的规则半球型巩膜,这一策略并不能适用于特定几何结构的巩膜。Pandolfi等人开发了一种基于参数的人类角膜网格生成器。它基于一个二维的网格生成算法,这个网络生成器构造出了角膜的结构。这个双锥形函数的输入限定在几个几何参数,这些参数描述了角膜的内、外表面。但是真实的巩膜形状并非规则的,并且形状对巩膜内部区域的厚度变化有重要影响。
发明内容
本发明的目的是提供一种可以为特定患者巩膜建模的方法。
本发明所采用的技术方案为:一种基于调和场的各向异性特定患者巩膜有限元建模方法,包括以下步骤:
S1、采用激光扫描得到特定患者巩膜立体模型;
S2、对所述立体模型指派巩膜权值、设计伞算子,得到调和场;
S3、为所述调和场设置约束条件,分布所述调和场的梯度场和等值线场;
S4、通过所述梯度场和所述等值线场平滑地分布规则六面体网格在所述巩膜上;
S5、对所述立体模型设置基于所述调和场的各向异性材料参数;
S6、通过IA-FEMesh生成器在所述巩膜上生成目标网格,将所述目标网格与所述规则六面体网格进行实验结果对比。
优选的,在步骤S3中,设置约束条件全局最小值赋值给所有约束最小值,全局最大值赋值给所有约束最大值。
优选的,所述步骤S5中,在所述巩膜的立体模型的每一个单元设置一个本地坐标系,其中,将环状方向作为X轴方向,梯度方向作为Y轴方向,环状方向与梯度方向的交点方向作为Z轴方向,所述Z轴方向用于表示所述巩膜的厚度方向,对于所述巩膜的各向异性弹性参数,设Ex=8.6MPa,Ey=6MPa,Ez=2.5MPa。
与相关技术相比,本发明提供的一种基于调和场的各向异性特定患者巩膜有限元建模方法可以在特定患者的非规则巩膜上建立模型,并基于调和场在巩膜上分布规则六面体网格,同时设置各向异性材料参数,从而可以为特定患者巩膜进行有效建模,这种建模方法还加入了各向异性材质方向场,不但保证了建模精确质量高,而且可以对非规则半球巩膜建模进行个性化参数设置。
附图说明
为了更清楚地说明本发明实施例中的技术方案,下面将对实施例描述中所需要使用的附图作简单地介绍,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其它的附图,其中:
图1为本发明的提供调和场的各向异性有限元建模方法的流程图;
图2(a)为巩膜模型的后部示意图;
图2(b)为巩膜模型的赤道部示意图;
图2(c)为巩膜模型的前部示意图;
图2(d)为带等厚线的后部巩膜视图;
图3(a)为巩膜厚度颜色映射示意图;
图3(b)为后巩膜的六面体网格图;
图3(c)为后巩膜和视神经乳头的网格图;
图4(a)和(b)为完美的周围形的巩膜纤维示图;
图4(c)为视神经乳头纤维旋转30°,外围纤维旋转40°后示图;
图5(a)为IA-FEMesh多模块结构的手动编辑网格化;
图5(b)为IA-FEMesh后巩膜六面体网格。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅是本发明的一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例,都属于本发明保护的范围。
本发明提供的基于调和场的各向异性特定患者巩膜有限元建模方法,是在非规则的特定患者巩膜上,平滑地分布一个规则六面体网格。其中,调和场是最有效的平滑分布工具之一,它的梯度向量场和等值线场平滑川流在模型表面,并且它们互相垂直,可以用它们来表示巩膜的主要纤维方向,也能驱动六面体网格化。
另外,用于表面网格的拉普拉斯-贝尔特拉米算子接近于模型的正常平均曲率,以致于调和场分布图像能很好的符合模型形状,本发明采样了梯度向量场和等值线场的流线去生成六面体网格。这些采样的流线精确地川流在模型表面,它可以保留巩膜的最初形状去生成“特定患者”网格,这样可以有效地处理特定患者巩膜的规则形状和显著的巩膜内部区域厚度变化。
请参阅图1,图1为本发明的提供调和场的各向异性有限元建模方法的流程图。本发明提供了一种基于调和场的各向异性特定患者巩膜有限元建模方法,包括以下步骤:
S1、采用激光扫描得到特定患者巩膜立体模型;
S2、对所述立体模型指派巩膜权值、设计伞算子,得到调和场;
S3、为所述调和场设置约束条件,分布所述调和场的梯度场和等值线场;
S4、通过所述梯度场和所述等值线场平滑地分布规则六面体网格在所述巩膜上;
S5、对所述立体模型设置基于所述调和场的各向异性材料参数;
S6、通过IA-FEMesh生成器在所述巩膜上生成目标网格,将所述目标网格与所述规则六面体网格进行实验结果对比。
对于上述六个步骤,具体包括如下内容:
1.采用高精度激光扫描得到特定患者巩膜立体模型。
请一并参阅图2(a)-2(c)。所述巩膜的立体模型包括外表面和内表面,用三角分段线性表面网格来表示。围绕视神经乳头的后巩膜是一个特殊兴趣区域,提取后巩膜半球模型,计算并添加等厚线。最厚的巩膜区域是眼睛的后极部,为1.1毫米;最薄的区域出现在赤道部上,为0.38毫米,具体如图2(d)所示,为带等厚线的后部巩膜视图。
2.对所述立体模型指派巩膜权值、设计伞算子,得到调和场。
通过在视神经乳头和赤道部的最后端环中定位最小值和最大值约束条件(如图3(a)中的箭头所指区域),生成调和场用以驱动网格化策略。对于后巩膜的外表面和内表面网格,构建一个调和场f使得Δf=0,其中,Δ是拉普拉斯算子,受限于狄利克雷边界条件。标准拉普拉斯算子是定义在分段线性表面网格M的伞算子,所述伞算子式子如下:
其中,j∈N(i)是表示与点i相邻的顶点集,wij是指派给边(i,j)的巩膜权值。巩膜权值wij的标准选择为离散的调和权值这里αij和βij表示相对边的角度。将顶点函数值fi组装成n维向量f,拉普拉斯算子可以写成Lf=0,其中,的L由以下式子确定:
消除约束点对应的行和列,把它们移到等式右边得到一个线性系统形式Ax=b,包含一个正定稀疏矩阵A和右边的向量b。预处理共轭梯度这种迭代算法能有效地得到线性系统的解。
3.为所述调和场设置约束条件,分布所述调和场的梯度场和等值线场。
对于后巩膜的外表面和内表面,一组相同的采样种子点放置在赤道部的后圈上,顺着由调和场标量域导出的梯度向量场方向。这些梯度流线平滑地川流和精确的收敛在视神经乳头的最末端圈上。同时,用同样的标量值采样一组调和巩膜场等值线。通过实践尝试,将后巩膜半球的梯度场和等值线场的线条数分别设置为60和15,可以很好的拟合出后巩膜的形状。等值线分布密度是根据巩膜的薄厚程度来的,靠近乳头区域线条分布的比赤道部的要密。这些梯度流线和等值线形成了一个完美的巩膜外表面和内表面四边形网格,且额外地在外表面与内表面之间插入了一个中间层。
4.通过所述梯度场和所述等值线场平滑地分布规则六面体网格在所述巩膜上。
这些梯度流线和等值线形成了一个完美的巩膜外表面和内表面四边形网格,且额外地在外表面与内表面之间插入了一个中间层。四边形网格的这3层有同样的拓扑结构,并且可以自动连接去生成完整的巩膜六面体网格,详见图3(b)所示。最后,这些8个结点的线性六面体元被转化成20结点的非线性六面体元,以解释大型材料形变的原因并提高有限元的精确性,具体如图3(c)所示,为后巩膜和视神经乳头的网格示图。
值得注意的是,尽管这个特定患者的后巩膜形状是非规则的半球型结构,并且视神经乳头也不在半球的中心,但是这个六面体网格结果仍然是非常规则的排列着。此外,这个网络元是自适合地分布,密度大的在视神经乳头区域,密度小的在赤道部区域。
5.对所述立体模型设置基于所述调和场的各向异性材料参数。
根据以前的研究,巩膜具有各向异性的特性,并且巩膜纤维几乎都在视神经乳头圆周面。在巩膜视神经乳头区域,纤维很少与前部和赤道部对齐。换而言之,需要对所述规则六面体网格中的每一个单元定义一个本地坐标系,作为各向异性材料的设置。如果使用的是规则的半球型巩膜模型的话,这个步骤很简单,但是对于非规则的特定患者巩膜,生成点到点的内部组织变化正交轴向量场是比较困难。
为了解决上述难题,本发明使用已有的调和场。很明显,等值线的切线方向与巩膜周围的纤维轨迹是兼容的,梯度向量的方向与经线纤维轨迹是兼容的。具体来说,为所述巩膜的立体模型中的每一个单元设置一个本地坐标系,将环状方向作为X轴方向,梯度方向作为Y轴方向,详见图4(a)和图4(b)。环状方向与梯度方向的交点方向作为Z轴方向,表示巩膜的厚度方向。
对于巩膜的各向异性弹性参数,我们基于之于的研究,设Ex=8.6MPa,Ey=6MPa,Ez=2.5MPa。假设组织是不可压缩的,为了避免非收敛数值的情况,将泊松率设为0.49。
为了确定哪种程度的各向异性特性会影响巩膜的应力和应变分布,创建了一个各向同性的巩膜有限元模型。为了兼容以前的研究,取上界与下界的平均值作为各向同性弹性参数E=3.8MPa。此外,正如之前研究指出的,巩膜纤维几乎都在视神经乳头圆周面,纤维很少与前部和赤道部对齐,我们创建了一组模型来作比较。我们首先在“r1”模型中,将外围区域纤维旋转10°;然后在“r2”模型中,将神经头区域纤维旋转10°,外围区域旋转20°;将模型“r3”的神经头区纤维旋转20°,外围区纤维旋转30°;在模型“r10”中,将神经头区域纤维旋转80°,外围区域旋转90°;图4(c)中模型为“r4”,“r4”模型中视神经乳头纤维旋转30°,外围纤维旋转40°后示图。
6.通过IA-FEMesh生成器在所述巩膜上生成目标网格,将所述目标网格与所述规则六面体网格进行实验结果对比。
对于图3(a)-3(c)中的巩膜模型,我们的基于调和场的六面体巩膜网格生成器,产生5410个六面体元。总共的计算耗时938ms。采用灵敏度分析去衡量计算结果的质量。结果表明,最大位移变化小于0.1%,然而,最大的主要应力少于2.5%,这证明了这种网格的合理性。
为了对我们的模型进行比较,我们使用了IA-FEMesh生成器,它采用了多模块网格化策略来生成目标网格,这种目标网格同样为六面体网格。这里需要非平凡边界,“块结构”技术的运用,使用户可以手动打破域而进入到拓扑块,正如图5(a)所示。这种块生成步骤花费了3分钟时间。计算过程主要包括两个步骤。首先,最近点投影用于将直线非结构化网格表面结点转化为兴趣区的底部表面。当表面节点建立之后,用户可以使用椭圆或超限插值来计算内部结点。使用同样的巩膜模型,总共的映射计算耗时大概2分钟,图5(b)为生成的后巩膜六面体网格。
从完美六面体元的双面角变形统计分析结果表明,我们的方法可以产生更高质量的网格元,如表1所示,表1为单元格双面角扭曲统计分析结果对照表。
表1
Our Method | Mesh Matching Method | |
0°~10° | 41.2% | 19.7% |
10°~20° | 43.6% | 30.2% |
20°~30° | 11.3% | 20.3% |
30°~40° | 2.3% | 14.1% |
40°~50° | 1.6% | 9.3% |
>50° | 0% | 6.4% |
Maximum | 43° | 86° |
Average | 7.8° | 17.5° |
可以观察到,在本发明提供的方法中,双面角变形小于20°的占85%,但在IA-FEMesh生成器中只占50%。这种方法的平均双面角扭曲分别是7.8°和17.5°,并且最大的扭曲角度是43°和86°。换句话说,IA-FEMesh生成器会产生翻转或自交元。
IA-FEMesh生成器中,在一些网格顶点中,最近点映射算法结果偏离了输入模型,尤其在较尖锐的边上。相反,我们的调和场网格生成器要求所有的网格输入模型是精确的,并且这样可以更好的保留特定患者巩膜的原始几何形状。
与相关技术相比,本发明提供的一种基于调和场的各向异性特定患者巩膜有限元建模方法可以在特定患者的非规则巩膜上建立模型,并基于调和场在巩膜上分布规则六面体网格,同时设置各向异性材料参数,从而可以为特定患者巩膜进行有效建模,这种建模方法还加入了各向异性材质方向场,不但保证了建模精确质量高,而且可以对非规则半球巩膜建模进行个性化参数设置。
以上所述仅为本发明的实施例,并非因此限制本发明的专利范围,凡是利用本发明说明书及附图内容所作的等效结构或等效流程变换,或直接或间接运用在其它相关的技术领域,均同理包括在本发明的专利保护范围内。
Claims (3)
1.一种基于调和场的各向异性特定患者巩膜有限元建模方法,其特征在于,包括以下步骤:
S1、采用激光扫描得到特定患者巩膜立体模型;
S2、对所述立体模型指派巩膜权值、设计伞算子,得到调和场;
通过在视神经乳头和赤道部的最后端环中定位最小值和最大值约束条件,生成调和场用以驱动网格化策略;
S3、为所述调和场设置约束条件,分布所述调和场的梯度场和等值线场;
S4、通过所述梯度场和所述等值线场平滑地分布规则六面体网格在所述巩膜上,所述梯度场的梯度流线和所述等值线场的等值线形成了巩膜的外表面和内表面四边形网格,且额外地在外表面与内表面之间插入了一个中间层,四边形网格的这三层有同样的拓扑结构,并且可自动连接去生成完整的巩膜六面体网格;
S5、对所述立体模型设置基于所述调和场的各向异性材料参数,在所述巩膜的立体模型的每一个单元设置一个本地坐标系,其中,将环状方向作为X轴方向,梯度方向作为Y轴方向,环状方向与梯度方向的交点方向作为Z轴方向,所述Z轴方向用于表示所述巩膜的厚度方向,对于所述巩膜的各向异性弹性参数,设Ex=8.6MPa,Ey=6MPa,Ez=2.5MPa;
S6、通过IA-FEMesh生成器在所述巩膜上生成目标网格,将所述目标网格与所述规则六面体网格进行实验结果对比。
3.根据权利要求1所述的基于调和场的各向异性特定患者巩膜有限元建模方法,其特征在于,在步骤S3中,设置约束条件全局最小值赋值给所有约束最小值,全局最大值赋值给所有约束最大值。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710085497.5A CN106909724B (zh) | 2017-02-17 | 2017-02-17 | 基于调和场的各向异性特定患者巩膜有限元建模方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710085497.5A CN106909724B (zh) | 2017-02-17 | 2017-02-17 | 基于调和场的各向异性特定患者巩膜有限元建模方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN106909724A CN106909724A (zh) | 2017-06-30 |
CN106909724B true CN106909724B (zh) | 2020-06-23 |
Family
ID=59208298
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710085497.5A Active CN106909724B (zh) | 2017-02-17 | 2017-02-17 | 基于调和场的各向异性特定患者巩膜有限元建模方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN106909724B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110070096B (zh) * | 2019-05-31 | 2021-01-12 | 中国科学院自动化研究所 | 针对非刚性形状匹配的局部频域描述子生成方法及装置 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102654923A (zh) * | 2012-04-20 | 2012-09-05 | 浙江大学 | 一种由内向外的六面体网格生成方法 |
CN103258343A (zh) * | 2013-05-07 | 2013-08-21 | 中南大学 | 一种基于图像编辑的眼部图像处理方法 |
-
2017
- 2017-02-17 CN CN201710085497.5A patent/CN106909724B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102654923A (zh) * | 2012-04-20 | 2012-09-05 | 浙江大学 | 一种由内向外的六面体网格生成方法 |
CN103258343A (zh) * | 2013-05-07 | 2013-08-21 | 中南大学 | 一种基于图像编辑的眼部图像处理方法 |
Non-Patent Citations (1)
Title |
---|
"医学三维模型重建和体网格生成与变形研究";廖胜辉;《中国博士学位论文全文数据库 信息科技辑 I138-73》;20080815;第21-84页 * |
Also Published As
Publication number | Publication date |
---|---|
CN106909724A (zh) | 2017-06-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US5891131A (en) | Method and apparatus for automated simulation and design of corneal refractive procedures | |
US8812280B2 (en) | Method for modeling biomechanical properties of an eye | |
De Greef et al. | Effects of model definitions and parameter values in finite element modeling of human middle ear mechanics | |
Ariza-Gracia et al. | Automatized patient-specific methodology for numerical determination of biomechanical corneal response | |
US20070162265A1 (en) | Method and apparatus for automated simulation and design of corneal refractive procedures | |
EP3158537B1 (en) | Removing sharp cusps from 3d shapes for additive manufacturing | |
Karimi et al. | Analysis of the effects of finite element type within a 3D biomechanical model of a human optic nerve head and posterior pole | |
JP4654028B2 (ja) | 視力を改善する方法およびシステム | |
Griffiths et al. | Analysis of cornea curvature using radial basis functions–Part I: Methodology | |
KR102630585B1 (ko) | 시각적 작업을 위한 시각적 장비의 성능을 평가하기 위한 장치 및 방법 | |
CN106909724B (zh) | 基于调和场的各向异性特定患者巩膜有限元建模方法 | |
Simonini et al. | Numerical estimation of stress and refractive power maps in healthy and keratoconus eyes | |
CN109222884A (zh) | 一种角膜塑形镜配适状态评估方法及装置 | |
Tang et al. | Mean curvature mapping for detection of corneal shape abnormality | |
Dai et al. | Simulating the effects of elevated intraocular pressure on ocular structures using a global finite element model of the human eye | |
CN108763709B (zh) | 一种基于有限元分析的再造耳支架拓扑结构优化方法 | |
CN113171172B (zh) | 一种角膜术后情况的模拟方法 | |
Doll et al. | Which feature influences on-eye power change of soft toric contact lenses: Design or corneal shape? | |
Jia et al. | Anisotropic Finite Element Modeling Based on a Harmonic Field for Patient‐Specific Sclera | |
Crouch et al. | Finite element model of cornea deformation | |
Dai et al. | Finite element analysis of the mechanical characteristics of glaucoma | |
Kazmitcheff et al. | Registration of a validated mechanical atlas of middle ear for surgical simulation | |
Griffiths et al. | Analysis of cornea curvature using radial basis functions–Part II: Fitting to data-set | |
Ballesta et al. | Variable complexity corneal surfaces characterization by modal geometrical reconstruction methods: comparative study | |
Katani et al. | Numerical solution of nonlinear integral equation arising in optometry |
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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |