CN1270672C

CN1270672C - Method for correcting brain tissue deformation in navigation system of neurosurgery

Info

Publication number: CN1270672C
Application number: CNB200410024847XA
Authority: CN
Inventors: 宋志坚; 刘翌勋; 李文生; 王满宁; 谢震中; 杜文健
Original assignee: Fudan University
Current assignee: Fudan University
Priority date: 2004-06-01
Filing date: 2004-06-01
Publication date: 2006-08-23
Anticipated expiration: 2024-06-01
Also published as: CN1582863A

Abstract

The present invention relates to a method for correcting brain tissue deformation in a navigation system of a neurosurgery operation, which belongs to the field of medical image treatment and application. Firstly, the present invention uses a three-dimensional automatic partitioning algorithm based on MRI to obtain target tissue (brain tissue), and then, lattices are formed on the partitioned brain tissue; subsequently, a physical model of the brain tissue is established by endowing each lattice unit with a corresponding biomechanics property on the basis of a line elasticity theory; the movement of a bare brain cortex is tracked by three-dimensional laser scanning equipment in a tracking algorithm and used as a boundary condition for finite element calculation by combining the physical model, and thus, the deformation of any position of the whole brain tissue is obtained; finally, a three-dimensional data field is updated by an insertion and return algorithm for guiding operation before the operation. The method has the advantages of simple implementation, reliable precision and convenient clinical application, and can be integrated in the existing navigational systems, and thus, the present invention greatly improves the precision of the navigation systems.

Description

A method for correcting brain tissue deformation in a neurosurgery navigation system

技术领域technical field

本发明属医学图像处理及应用领域，涉及一种外科手术导航系统精度校正方法，具体涉及一种神经外科手术导航系统中脑组织变形校正的方法。The invention belongs to the field of medical image processing and application, and relates to a precision correction method for a surgical navigation system, in particular to a method for correcting brain tissue deformation in a neurosurgery navigation system.

背景技术Background technique

临床外科手术中脑组织变形是影响神经外科手术导航系统精度的重要因素。目前的解决方法包括术中影像校正如术中MRI、术中US等和物理模型校正两种方法。其中术中影像校正是精度最高的方法，但缺点是费用昂贵，且容易造成术中感染，故当前研究主要集中在基于物理模型的校正方法上。基于物理模型的方法能够通过脑组织的生物力学属性约束脑组织的运动，计算量小、精度可靠、实施简单，便于临床应用，因此，它是当前国际上的研究热点。Yale大学的ImageProcessing and Analysis课题组在2001年11月的IEEE Workshop on MathematicalMethods in Biomedical Image Analysis上提出一种线弹性物理模型，他们将表面位移作为边界条件，并采用双立体相机获取边界条件，这种获取边界条件的方法需要改变现有的导航设备来固定两个相机，同时两个相机的校正、同步使得术中操作更为复杂，因此不能方便地应用于临床。Vanderbilt大学的Biomedical Modeling实验室在《IEEE Transactions on medical imaging》(Vol.18，No.10：866-874，1999.)发表的“Model-Updated Image Guidance：Initial Clinical Experiencewith Gravity-Induced Brain Deformation”一文中提出一种固结理论模型，但由于该模型需要将体力作为边界条件，即需要在术中获得脑脊液的位置，然而脑脊液在临床上是难以测量的，使得该模型目前也无法实现临床应用。Brain tissue deformation in clinical surgery is an important factor affecting the accuracy of neurosurgery navigation system. Current solutions include intraoperative image correction such as intraoperative MRI, intraoperative US, etc., and physical model correction. Among them, intraoperative image correction is the method with the highest accuracy, but the disadvantage is that it is expensive and easy to cause intraoperative infection. Therefore, current research mainly focuses on the correction method based on physical models. The method based on the physical model can constrain the movement of brain tissue through the biomechanical properties of the brain tissue. It has the advantages of small calculation, reliable accuracy, simple implementation, and convenient clinical application. Therefore, it is a current research hotspot in the world. The Image Processing and Analysis research group of Yale University proposed a linear elastic physical model on the IEEE Workshop on Mathematical Methods in Biomedical Image Analysis in November 2001. They used surface displacement as the boundary condition and used dual stereo cameras to obtain the boundary condition. The method of obtaining boundary conditions needs to change the existing navigation equipment to fix the two cameras. At the same time, the correction and synchronization of the two cameras make the intraoperative operation more complicated, so it cannot be conveniently applied clinically. "Model-Updated Image Guidance: Initial Clinical Experience with Gravity-Induced Brain Deformation" published by the Biomedical Modeling Laboratory of Vanderbilt University in "IEEE Transactions on medical imaging" (Vol.18, No.10:866-874, 1999.) A consolidation theoretical model is proposed in this paper, but because the model requires physical strength as a boundary condition, that is, the position of the cerebrospinal fluid needs to be obtained during the operation. However, the cerebrospinal fluid is difficult to measure clinically, making this model currently unavailable for clinical application.

发明内容Contents of the invention

本发明的目的是为临床应用提供一种实施简单，操作灵活，无需更改现有导航设备的外科手术导航系统精度校正方法，具体涉及一种神经外科手术导航系统中脑组织变形校正的方法。The purpose of the present invention is to provide a method for correcting the accuracy of a surgical navigation system that is simple to implement and flexible in operation without changing the existing navigation equipment for clinical applications, and specifically relates to a method for correcting brain tissue deformation in a neurosurgery navigation system.

本发明采用基于线弹性理论的物理模型，同时为了方便有效的获取边界条件，采用三维激光成像设备LRS(Laser Range Scanner)实时获取术中裸露的大脑皮层，通过跟踪算法跟踪脑皮层运动，从而获得用于驱动模型的边界条件。将该边界条件与线弹性物理模型相结合，可以有效地对术中的脑组织变形进行校正。本发明方法实施简单，操作灵活，无需更改现有的导航设备，便于临床应用。The present invention adopts a physical model based on the theory of linear elasticity, and at the same time, in order to obtain boundary conditions conveniently and effectively, a three-dimensional laser imaging device LRS (Laser Range Scanner) is used to obtain real-time exposure of the cerebral cortex during the operation, and the movement of the cerebral cortex is tracked by a tracking algorithm, thereby obtaining Boundary conditions used to drive the model. Combining this boundary condition with a linear elastic physical model can effectively correct for intraoperative brain tissue deformation. The method of the invention is simple to implement, flexible to operate, does not need to change the existing navigation equipment, and is convenient for clinical application.

本发明的技术方案是：首先采用基于MRI的三维自动分割算法，获得目标区域如脑组织，随后将分割出的脑组织网格化，在线弹性理论的基础上通过对每一网格单元赋予相应的生物力学属性，建立脑组织的物理模型。借助三维激光扫描设备，通过跟踪算法跟踪裸露脑皮层的运动，将其作为边界条件并结合物理模型进行有限元计算，获得整个脑组织任意位置的变形，最后采用一种插回算法更新术前三维数据场用于指导手术。The technical solution of the present invention is: firstly adopt the three-dimensional automatic segmentation algorithm based on MRI to obtain the target area such as brain tissue, then grid the segmented brain tissue, and assign corresponding biomechanical properties of brain tissue to establish a physical model. With the help of 3D laser scanning equipment, the movement of the exposed cerebral cortex is tracked through a tracking algorithm, which is used as a boundary condition and combined with a physical model for finite element calculations to obtain the deformation of any position of the entire brain tissue, and finally an interpolation algorithm is used to update the preoperative 3D The data field is used to guide the operation.

本发明方法通过下述步骤实现，The inventive method realizes by following steps,

1.三维自动分割算法1. 3D automatic segmentation algorithm

解决脑组织变形的第一步是分割出目标区域如脑组织，本发明采用一种基于灰度直方图原理结合形态学特征的自动分割算法。所述算法首先采用高斯曲线拟合直方图曲线，根据拟合高斯曲线自动判定门限值，门限化后的结果可将脑组织与大部分的颅骨相分离。所采用的门限化判定公式为：The first step to solve the deformation of brain tissue is to segment the target area such as brain tissue. The present invention adopts an automatic segmentation algorithm based on the gray histogram principle combined with morphological features. The algorithm first adopts a Gaussian curve to fit the histogram curve, and automatically determines the threshold value according to the fitted Gaussian curve, and the thresholded result can separate the brain tissue from most of the skull. The thresholding judgment formula adopted is:

t_s＝t_BF+4/5(u_GM-t_BF)t _s ＝t _BF +4/5(u _GM -t _BF )

其中t_BF是背景与前景之间的门限，u_GM为脑白质灰度均值，分别对应高斯曲线的第一和第二个峰值。在门限化后，脑组织和颅骨之间还会存在细微的连接，此时可采用形态学中的腐蚀算法，去除所存在的细微连接；通过设定种子点寻求最大连通域可将脑组织从皮肤、骨骼中分离出来；最后针对分割出的脑组织采用形态学中的膨胀算法，恢复腐蚀操作中所损失的脑组织。Among them, t _BF is the threshold between the background and the foreground, and u _GM is the mean gray value of the brain white matter, corresponding to the first and second peaks of the Gaussian curve, respectively. After thresholding, there will still be subtle connections between the brain tissue and the skull. At this time, the erosion algorithm in morphology can be used to remove the existing subtle connections; by setting the seed point to seek the largest connected domain, the brain tissue can be removed from the The skin and bones are separated; finally, the morphological expansion algorithm is used for the segmented brain tissue to restore the brain tissue lost in the corrosion operation.

2.网格化2. Gridding

针对分割出的脑组织，采用四面体单元进行离散化。线弹性模型的特点在于每一单元都是常应变单元，本发明方法针对在脑组织边界、脑室边界处相对内部组织变化较大的现象，采用一种具有多分变率的网格划分算法，边界处单元密度大，而内部单元密度小，能精确表达表面，较好的模拟出边界处的变形情况，而对于内部变化相对较小的区域，单元则划分的较大些，减少了计算量。For the segmented brain tissue, the tetrahedron unit is used for discretization. The characteristic of the linear elastic model is that each unit is a constant strain unit. The method of the present invention aims at the phenomenon that the brain tissue boundary and the ventricle boundary have relatively large changes in the internal tissue, and adopts a grid division algorithm with multi-point variability. The unit density is high at the location, while the internal unit density is small, which can accurately express the surface and better simulate the deformation at the boundary. For the area with relatively small internal changes, the units are divided into larger units, which reduces the amount of calculation.

所述网格划分算法结合类八叉树和MT(Marching Tetrahedron)算法，采用case表提高算法效率。算法组织如下述，The grid division algorithm combines the octree-like and MT (Marching Tetrahedron) algorithms, and uses the case table to improve the efficiency of the algorithm. The algorithm is organized as follows,

1)采用类八叉树算法将三维图像空间划分为由均匀的六面体单元所组成的网格；1) Using the octree-like algorithm to divide the three-dimensional image space into grids composed of uniform hexahedral units;

2)针对每一六面体单元将其划分为五个四面体单元；2) divide it into five tetrahedral units for each hexahedral unit;

3)对2)所获得均匀四面体网格，在边界处加以细化，获得多分变率网格；3) The uniform tetrahedral grid obtained in 2) is refined at the boundary to obtain a multi-fraction variable rate grid;

4)采用类MT算法对网格进行切割，去除背景，获得最终的仅包含脑组织的多分变4) Use the MT-like algorithm to cut the grid, remove the background, and obtain the final multivariate model that only contains brain tissue.

率网格。rate grid.

3.物理建模3. Physical modeling

针对网格化后的脑组织，通过对每一单元赋予相应的生物力学属性建立脑组织的物理模型。由于脑组织变形过程缓慢而且形变小，应变与应力成线性关系，本发明将其模拟为基于线弹性理论的弹性体。模型的数学表达如下：For the meshed brain tissue, a physical model of the brain tissue is established by assigning corresponding biomechanical properties to each unit. Since the deformation process of the brain tissue is slow and the deformation is small, and the strain and the stress have a linear relationship, the present invention simulates it as an elastic body based on the theory of linear elasticity. The mathematical expression of the model is as follows:

${&dtri; &dtri;}^{22} {u u}_{x x} + + \frac{11}{11 - - 22 v v} \frac{&PartialD; &PartialD;}{&PartialD; &PartialD; x x} ((\frac{{&PartialD; &PartialD; u u}_{x x}}{&PartialD; &PartialD; x x} + + \frac{{&PartialD; &PartialD; u u}_{y the y}}{&PartialD; &PartialD; y the y} + + \frac{{&PartialD; &PartialD; u u}_{z z}}{&PartialD; &PartialD; z z})) + + \frac{{F f}_{x x}}{μ μ} = = 00$

${&dtri; &dtri;}^{22} {u u}_{y the y} + + \frac{11}{11 - - 22 v v} \frac{&PartialD; &PartialD;}{&PartialD; &PartialD; y the y} ((\frac{{&PartialD; &PartialD; u u}_{x x}}{&PartialD; &PartialD; x x} + + \frac{&PartialD; &PartialD; {u u}_{y the y}}{&PartialD; &PartialD; y the y} + + \frac{&PartialD; &PartialD; {u u}_{z z}}{&PartialD; &PartialD; z z})) + + \frac{{F f}_{y the y}}{μ μ} = = 00$

${&dtri; &dtri;}^{22} {u u}_{z z} + + \frac{11}{11 - - 22 v v} \frac{&PartialD; &PartialD;}{&PartialD; &PartialD; z z} ((\frac{&PartialD; &PartialD; {u u}_{x x}}{&PartialD; &PartialD; x x} + + \frac{&PartialD; &PartialD; {u u}_{y the y}}{&PartialD; &PartialD; y the y} + + \frac{{&PartialD; &PartialD; u u}_{z z}}{&PartialD; &PartialD; z z})) + + \frac{{F f}_{z z}}{μ μ} = = 00$

其中u为位移矢量，F为力矢量。v为泊松比，μ＝E/2(1+v)，E为弹性模量。Where u is the displacement vector and F is the force vector. v is Poisson's ratio, μ=E/2(1+v), and E is the modulus of elasticity.

4.边界条件4. Boundary conditions

物理模型完成后，采用合适的边界条件驱动模型。术中最易观察到的是裸露的大脑皮层，所以可通过获得脑皮层变形前和变形后的位置计算出表面位移。将该表面位移作为驱动模型的边界条件。本方法只需引入三维成像设备LRS，在开颅后，手术前进行扫描后即可移开，操作简单，无需更改现有导航设备且不占用术中空间，所述设备能够提供丰富的几何和纹理信息用于跟踪皮层的运动。脑皮层变形前的初始位置可通过刚体配准算法获得，After the physical model is complete, the model is driven with appropriate boundary conditions. The exposed cerebral cortex is the easiest to observe during the operation, so the surface displacement can be calculated by obtaining the pre-deformation and post-deformation positions of the cerebral cortex. Use this surface displacement as a boundary condition to drive the model. This method only needs to introduce the three-dimensional imaging device LRS, which can be removed after craniotomy and preoperative scanning. The operation is simple, no need to change the existing navigation equipment and does not occupy intraoperative space. Texture information is used to track the movement of the cortex. The initial position of the cortex before deformation can be obtained by the rigid body registration algorithm,

本发明采用一种简单的坐标转换来实现上述刚体配准算法。首先借助于术前的点配准(PBR)算法实现图像空间到参考架空间的变换，然后借助于Polaris实现参考架空间到Polaris空间以及Polaris空间到被跟踪器空间的变换，最后借助于校准器，实现被跟踪器空间到LRS空间的变换。经过上述一系列变换，最终可以实现图像空间到LRS空间的变换，从而可以获得脑皮层在LRS空间的初始位置。The present invention adopts a simple coordinate transformation to realize the above-mentioned rigid body registration algorithm. Firstly, the transformation from the image space to the reference frame space is realized by means of the preoperative point registration (PBR) algorithm, and then the transformation from the reference frame space to the Polaris space and from the Polaris space to the tracked space is realized by means of Polaris, and finally, the calibrator is used , realizing the transformation from the tracked space to the LRS space. After the above series of transformations, the transformation from the image space to the LRS space can be realized finally, so that the initial position of the cerebral cortex in the LRS space can be obtained.

变形后脑皮层位置可通过基于活动轮廓模型的非刚体面配准技术实现。所述模型所对应的能量函数定义如下述，The post-deformation cortical position can be achieved by a non-rigid body surface registration technique based on an active contour model. The energy function corresponding to the model is defined as follows,

$E E. ((v v)) = = {&Integral; &Integral; &Integral; &Integral; w w}_{1010} ((r r,, s the s)) {| | \frac{&PartialD; &PartialD; v v ((r r,, s the s))}{&PartialD; &PartialD; r r} | |}^{22} + + {w w}_{0101} ((r r,, s the s)) {| | \frac{&PartialD; &PartialD; v v ((r r,, s the s))}{&PartialD; &PartialD; s the s} | |}^{22} + + P P ((v v ((r r,, s the s)))) drds drds$

其中v(r，s)＝(v₁(r，s)，v₂(r，s)，v₃(r，s))为双参数曲面方程，w₁₀(r，s)为弹性系数函数，w₀₁(r，s)为刚性系数函数，P(v(r，s))为势能函数。Among them, v(r, s)=(v ₁ (r, s), v ₂ (r, s), v ₃ (r, s)) is a two-parameter surface equation, and w ₁₀ (r, s) is an elastic coefficient function , w ₀₁ (r, s) is the rigid coefficient function, and P(v(r, s)) is the potential energy function.

引入时间参数t，该能量函数的极值满足欧拉Euler-Lagrange方程：The time parameter t is introduced, and the extreme value of the energy function satisfies the Euler-Lagrange equation:

$\frac{&PartialD; &PartialD; v v}{&PartialD; &PartialD; t t} \{\begin{matrix} - - \frac{&PartialD; &PartialD;}{&PartialD; &PartialD; s the s} (({w w}_{1010} \cdot &Center Dot; \frac{&PartialD; &PartialD; v v}{&PartialD; &PartialD; s the s})) - - \frac{&PartialD; &PartialD;}{&PartialD; &PartialD; r r} (({w w}_{0101} \cdot \cdot \frac{&PartialD; &PartialD; v v}{&PartialD; &PartialD; r r})) = = - - F f ((v v)) \\ v v ((00,, s the s,, r r)) = = {v v}_{00} ((s the s,, r r)) \end{matrix}$

其中v₀(s，r)为脑皮层初始位置，可由刚体配准获得。F(v)为图像力，定义如下：Among them, v ₀ (s, r) is the initial position of the cerebral cortex, which can be obtained by rigid body registration. F(v) is the image force, defined as follows:

$F f ((v v)) = = &dtri; &dtri; ((P P ((v v)))) = = &dtri; &dtri; ((\frac{11}{11 + + {((&dtri; &dtri; I I))}^{22}}))$

其中I为三维激光扫描设备对裸露脑皮层扫描所成的三维表面。Wherein I is the three-dimensional surface formed by scanning the exposed cerebral cortex by the three-dimensional laser scanning equipment.

Euler-Lagrange方程可采用FEM(Finite Element Method)求解，从而获得各个节点的变形后位置，结合初始位置，可以获得脑皮层表面节点位移，将该表面位移作为边界条件驱动物理模型。The Euler-Lagrange equation can be solved by FEM (Finite Element Method), so as to obtain the deformed position of each node. Combined with the initial position, the displacement of the surface node of the cerebral cortex can be obtained, and the surface displacement can be used as a boundary condition to drive the physical model.

5.有限元计算5. Finite element calculation

根据获得的脑皮层位移，结合物理模型，采用有限元计算出任意节点处的位移，再结合形函数就可获得脑组织任意位置的变形。According to the obtained cerebral cortex displacement, combined with the physical model, the displacement at any node is calculated by using finite element, and then combined with the shape function, the deformation of any position of the brain tissue can be obtained.

6.更新原始三维数据场6. Update the original 3D data field

采用插回算法更新术前三维数据场。从变形后的网格单元出发，寻找出单元内的显示坐标点，利用形函数获得该点在未变形前的位置，再利用三线性插值获得该点的灰度值。针对所有变形后的单元进行上述处理，就可以用变形后的三维数据场更新术前三维数据场。The preoperative three-dimensional data field is updated by inserting back algorithm. Starting from the deformed grid unit, find the display coordinate point in the unit, use the shape function to obtain the position of the point before deformation, and then use trilinear interpolation to obtain the gray value of the point. By performing the above processing on all deformed units, the preoperative three-dimensional data field can be updated with the deformed three-dimensional data field.

本发明引入三维激光扫描设备获取边界条件，并结合线弹性物理模型预测整个脑组织变形情况，既保证了模型预测精度，又解决了边界条件难以测量这一难题，从而可以在临床上实施，大幅度提高手术导航系统的精度。The present invention introduces three-dimensional laser scanning equipment to obtain boundary conditions, and combines the linear elastic physical model to predict the deformation of the entire brain tissue, which not only ensures the prediction accuracy of the model, but also solves the problem that the boundary conditions are difficult to measure, so that it can be implemented clinically. Significantly improve the accuracy of surgical navigation systems.

附图说明Description of drawings

图1为解决脑组织变形的流程图。Figure 1 is a flowchart for addressing brain tissue deformation.

图2为基于256×256×48MRI数据场的三维自动分割算法结果。Figure 2 shows the results of the 3D automatic segmentation algorithm based on the 256×256×48 MRI data field.

图3为针对256×256×48MRI数据场多分辨率网格化的结果，Figure 3 is the result of multi-resolution gridding for 256×256×48 MRI data field,

其中，共划分18485个四面体单元，5410个节点；左图为两维断层图像，中间为该断层的网格化结果，右图为该三维数据场的网格化结果。Among them, a total of 18,485 tetrahedral units and 5,410 nodes are divided; the left image is a two-dimensional fault image, the middle image is the gridded result of the fault, and the right image is the gridded result of the three-dimensional data field.

图4为刚体配准算法，采用坐标变换实现。Figure 4 shows the rigid body registration algorithm, which is realized by coordinate transformation.

图5为边界条件，其中，左图为裸露脑皮层，用于跟踪术中脑皮层运动；中间为不加约束的脑皮层；右图为脑皮层底部，根据临床经验可认为固定不动。Figure 5 shows the boundary conditions. The left picture shows the exposed cortex, which is used to track the movement of the cerebral cortex during the operation; the middle picture shows the unconstrained cerebral cortex; the right picture shows the bottom of the cerebral cortex, which can be considered fixed according to clinical experience.

图6为三维变形结果，其中，左图是未变形前的脑组织，右图是重力所引起的变形后的脑组织。Figure 6 shows the results of three-dimensional deformation, in which the left picture is the brain tissue before deformation, and the right picture is the brain tissue after deformation caused by gravity.

图7为两维变形结果，其中，左图是未变形前的脑组织断层，右图是变形后的相应断层。Figure 7 shows the results of two-dimensional deformation, in which the left image is the brain tissue section before deformation, and the right image is the corresponding section after deformation.

图8为位移场的三维可视化结果，其中，箭头表示位移方向，颜色表示位移大小，颜色由红到篮表示位移由小到大。Figure 8 is the three-dimensional visualization result of the displacement field, where the arrow indicates the displacement direction, the color indicates the displacement size, and the color from red to blue indicates the displacement from small to large.

图9为3D Digital公司的三维激光扫描设备，型号为Model 200，Figure 9 shows the 3D laser scanning equipment of 3D Digital Company, the model is Model 200,

其中200型有一个彩色摄像机，可以将彩色照片和数据点云叠加，有效测量范围是200mm至750mm。仪器在距离被测物体450mm处的测量精度是±125um。Among them, the 200 type has a color camera, which can superimpose color photos and data point clouds, and the effective measurement range is 200mm to 750mm. The measurement accuracy of the instrument at a distance of 450mm from the measured object is ±125um.

图10为脑皮层的扫描结果Figure 10 shows the scan results of the cerebral cortex

具体实施方式Detailed ways

实施例1Example 1

1.针对256×256×48的三维MRI数据场采用三维自动分割算法，获得脑组织。门限值分别对应拟合高斯曲线的第一和第二个峰值，腐蚀元素采用半径为5个象素的球状元素，膨胀元素采用半径为6个象素的球状元素。1. For the 256×256×48 3D MRI data field, use the 3D automatic segmentation algorithm to obtain brain tissue. The threshold values correspond to the first and second peaks of the fitted Gaussian curve respectively, the corrosion element adopts a spherical element with a radius of 5 pixels, and the expansion element adopts a spherical element with a radius of 6 pixels.

2.采用多分辨率网格化算法，将所分割出的脑组织离散为18485个四面体，节点个数为5410。边界处最大四面体为7.5×7.5×7.5mm³(用四面体外接六面体大小衡量)，内部最大四面体为15×15×15mm³；2. Using a multi-resolution grid algorithm, the segmented brain tissue is discretized into 18485 tetrahedrons with 5410 nodes. The largest tetrahedron at the boundary is 7.5×7.5×7.5mm ³ (measured by the size of the hexahedron circumscribed by the tetrahedron), and the largest tetrahedron inside is 15×15×15mm ³ ;

3.对每一单元设置脑组织生物力学属性参数。杨氏模量＝3Kpa，泊松比＝0.45；3. Set the biomechanical property parameters of brain tissue for each unit. Young's modulus=3Kpa, Poisson's ratio=0.45;

4.采用坐标变换实现刚体配准，获得脑皮层在LRS空间的初始位置。跟踪器固定在LRS上，Polaris通过跟踪器跟踪LRS，实现LRS在Polaris监控范围内任意移动。4. Use coordinate transformation to realize rigid body registration, and obtain the initial position of the cerebral cortex in the LRS space. The tracker is fixed on the LRS, and Polaris tracks the LRS through the tracker, so that the LRS can move freely within the monitoring range of Polaris.

5.在开颅后，手术前利用LRS进行一次扫描，LRS距离病人开颅位置的距离控制在450mm左右，扫描时间5～7秒，分辨率512(水平)×500(垂直)。根据设备提供的变形后的三维表面信息以及由步骤4获得的变形前三维表面信息，结合非刚体面配准技术，跟踪脑皮层运动，获得表面节点位移。5. After the craniotomy, use LRS to scan once before the operation. The distance between the LRS and the patient's craniotomy position is controlled at about 450mm, the scanning time is 5-7 seconds, and the resolution is 512 (horizontal) × 500 (vertical). According to the deformed 3D surface information provided by the device and the pre-deformed 3D surface information obtained in step 4, combined with non-rigid body surface registration technology, the movement of the cerebral cortex is tracked to obtain the surface node displacement.

6.采用有限元法将物理模型方程化为矩阵形式：Au＝b(A为总刚矩阵，u为待求位移矢量，b为节点力矢量)，采用共扼梯度法求解含有5410(节点数)×3(维数)＝16230个方程的线性方程组。引入由步骤5所获得的边界条件，消除矩阵A的奇异性，求出位移矢量u，获得脑组织内部任意节点处的位移，再结合有限元中的形函数就可以插值出脑组织内部任意位置的变形。6. Use the finite element method to convert the physical model equation into a matrix form: Au=b (A is the total stiffness matrix, u is the displacement vector to be obtained, and b is the nodal force vector), and the conjugate gradient method is used to solve the equation containing 5410 (number of nodes) )×3 (dimension)=16230 linear equations of equations. Introduce the boundary conditions obtained in step 5, eliminate the singularity of the matrix A, obtain the displacement vector u, and obtain the displacement at any node inside the brain tissue, and then combine the shape function in the finite element to interpolate any position inside the brain tissue deformation.

7.由变形后的数据场，根据形函数计算出所有对显示有贡献的坐标点(整数坐标点)变形前的位置，利用三线性插值计算出该点的灰度值。最后对所有这些整数坐标点组成的三维数据场采用经典的Raycasting算法加以可视化，用于指导手术。7. From the deformed data field, calculate the positions of all coordinate points (integer coordinate points) that contribute to the display before deformation according to the shape function, and use trilinear interpolation to calculate the gray value of the point. Finally, the classic Raycasting algorithm is used to visualize the three-dimensional data field composed of all these integer coordinate points, which is used to guide the operation.

在奔IV2.6G、1G内存、Win2000操作系统下运行，需要2分钟左右。It takes about 2 minutes to run under Ben IV2.6G, 1G memory, and Win2000 operating system.

Claims

1. the gauged method of neurosurgery navigation system midbrain metaplasia is characterized in that comprising the steps: the three-dimensional automatically partitioning algorithm of (1) employing based on MRI, obtains cerebral tissue;

(2) adopt the multi-resolution grid algorithm that the cerebral tissue that is partitioned into is carried out gridding;

(3) give corresponding biomechanics attribute to each grid cell, set up the physical model of cerebral tissue;

(4) obtain the distortion front surface by the rigid body registration Algorithm, open cranium after, adopt 3 D laser scanning equipment, scan strained cortex and obtain the distortion rear surface, by the cortex track algorithm, obtain boundary condition;

(5) carry out FEM calculation in conjunction with boundary condition and physical model, obtain the distortion of all cell node positions, calculate the distortion of cerebral tissue optional position in conjunction with shape function;

(6) employing turns back to the algorithm preceding 3 d data field of distortion renewal art that calculates.

2, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1, it is characterized in that, wherein the three-dimensional partitioning algorithm automatically of step (1) is, adopt the gaussian curve approximation histogram curve, according to the automatic decision gate limit value of fitted Gaussian curve, described thresholding judgement formula is

t _s＝t _BF+4/5(u _GM-t _BF)

T wherein _BFBe the thresholding between background and the prospect, u _GMBe the alba gray average, respectively first and second peak values of corresponding Gaussian curve.

3, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1, it is characterized in that, wherein the multi-resolution grid algorithm of step (2) is, after with class Octree algorithm the 3-D view spatial division being the grid of hexahedral element, each hexahedral element is divided into five tetrahedron elements, to the tetrahedral grid boundary refinement that is obtained, obtain multi-resolution grid, adopt class MarchingTetrahedron algorithm cutting grid, remove background, obtain final multi-resolution grid.

4, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1 is characterized in that wherein the mathematical expression of the physical model of step (3) is,

{&dtri;}^{2} u_{x} + \frac{1}{1 - 2 v} \frac{&PartialD;}{&PartialD; x} (\frac{{&PartialD; u}_{x}}{&PartialD; x} + \frac{{&PartialD; u}_{y}}{&PartialD; y} \frac{{&PartialD; u}_{z}}{&PartialD; z}) + \frac{F_{x}}{μ} = 0

{&dtri;}^{2} u_{y} + \frac{1}{1 - 2 v} \frac{&PartialD;}{&PartialD; y} (\frac{{&PartialD; u}_{x}}{&PartialD; x} + \frac{{&PartialD; u}_{y}}{&PartialD; y} \frac{{&PartialD; u}_{z}}{&PartialD; z}) + \frac{F_{y}}{μ} = 0

{&dtri;}^{2} u_{z} + \frac{1}{1 - 2 v} \frac{&PartialD;}{&PartialD; z} (\frac{{&PartialD; u}_{x}}{&PartialD; x} + \frac{{&PartialD; u}_{y}}{&PartialD; y} \frac{{&PartialD; u}_{z}}{&PartialD; z}) + \frac{F_{z}}{μ} = 0

Wherein u is a displacement vector, and F is a force vector.μ=E/2 (1+v), v are Poisson's ratio, and E is an elastic modelling quantity.

5, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1 is characterized in that, wherein the described cortex track algorithm of step (4) is,

At first adopt the rigid body registration Algorithm to obtain the cortex initial position, algorithm is realized by means of coordinate transform, adopts the non-rigid body surface adjustment Technical Follow-Up cortex motion based on movable contour model then, and the energy function of described movable contour model correspondence defines and is,

E (v) = {&Integral; &Integral; w}_{10} (r, s) {| \frac{&PartialD; v (r, s)}{&PartialD; r} |}^{2} + w_{01} (r, s) {| \frac{&PartialD; v (r, s)}{&PartialD; s} |}^{2} + P (v (r, s)) drds

Wherein v (r, s)=(v ₁(r, s), v ₂(r, s), v ₃(r, s)) is two-parameter surface equation, w ₁₀(r s) is the coefficient of elasticity function, w ₀₁(r s) is the stiffness coefficient function, and (v (r, s)) is a potential-energy function to P;

Introduce time parameter t, the extreme value of described energy function satisfies Euler Euler-Lagrange equation,

\{\begin{matrix} \frac{&PartialD; v}{&PartialD; t} - \frac{&PartialD;}{&PartialD; s} (w_{10} \cdot \frac{&PartialD; v}{&PartialD; s}) - \frac{&PartialD;}{&PartialD; r} (w_{01} \cdot \frac{&PartialD; v}{&PartialD; r}) = - F (v) \\ v (0, s, r) = v_{0} (s, r) \end{matrix}

V wherein ₀(s r) is the cortex initial position, is obtained by the rigid body registration Algorithm;

Wherein F (v) be image force, definition is,

F (v) = &dtri; (P (v)) = &dtri; (\frac{1}{1 + {(&dtri; I)}^{2}})

Wherein I is the imaging of three-dimensional laser device scan.

6, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1, it is characterized in that, wherein step (5) is a boundary condition with the cortex motion, the model equation that adopts Finite Element Method to find the solution the partial differential equation form, obtain the distortion of all cell node positions, calculate the distortion of cerebral tissue optional position again in conjunction with shape function.

7, the gauged method of neurosurgery navigation system midbrain metaplasia according to claim 1, it is characterized in that, wherein the described algorithm that turns back to of step (6) is, from the distortion after grid cell, seek out the rounded coordinate point in the unit, utilize shape function to obtain this position before not being out of shape, utilize Tri linear interpolation to obtain the gray value of this point again, to the distortion after the unit carry out above-mentioned processing after, the 3 d data field after the available distortion upgrades original 3 d data field.