CN109147061B - Method for carrying out human-computer interaction editing on segmentation result of volume data - Google Patents

Abstract

The invention relates to a method for carrying out human-computer interaction editing on a segmentation result of volume data, and belongs to the technical field of digital medical treatment. The method comprises the following steps: performing triangular mesh model generation and surface drawing on the volume data segmentation result to be edited; performing volume rendering on the volume data; mixing the volume rendering and surface rendering results to obtain a uniform view; recording user mouse interaction to obtain a pen touch line, acquiring intersection of an extension plane and volume data by using a ray projection method to obtain a cross section grid, and displaying the volume data gray level; intersecting the section grids with the segmentation result grids to obtain intersecting lines; carrying out planarization on the cross-section grids and the intersecting lines, and carrying out interactive editing on the intersecting lines on the obtained two-dimensional view; and establishing a triangular mesh optimization model based on volume data gradient, and decomposing and iteratively optimizing the objective function according to the cross line editing result to obtain the optimized triangular mesh vertex coordinates, namely the interactive editing result. The method improves the segmentation precision and the interactive editing precision of the volume data.

Description

Method for carrying out human-computer interaction editing on segmentation result of volume data

Technical Field

The invention relates to a method for carrying out human-computer interaction editing on a segmentation result of volume data, and belongs to the technical field of digital medical treatment.

Background

In the process of digital medical treatment, the volume data of a patient is acquired through diagnosis and treatment means such as CT and MRI, and a target region needs to be acquired through a volume data segmentation technology, which is helpful for further analysis and diagnosis. Due to the complexity of medical volume data, the existing automatic volume data segmentation technology has poor universality, the segmentation result usually has problems, and the segmentation result needs to be further edited by using a human-computer interaction technology. However, it is difficult to perform digitized human-computer interaction editing on three-dimensional volume data, and the existing technology cannot well solve the problems of high degree of freedom, poor interaction precision and the like faced by three-dimensional human-computer interaction. In order to efficiently edit the volume data division result, a well-designed interaction method is required. In the interactive volume data segmentation method adopted by software such as Seg3D, the user marks a segmentation target on a two-dimensional slice of volume data using a brush. In order to complete the segmentation task, a large number of two-dimensional slices need to be marked, and the interaction efficiency is low. Or by means of an automatic segmentation method provided by software, the obtained segmentation result may have certain errors. The existing interactive body data segmentation method based on grid editing expresses the segmentation result in the form of a triangular grid model, and the segmentation purpose is completed by constructing the triangular grid model through user interaction or edits the grid shape through interaction to edit the segmentation result. The method does not fully consider the volume data information in the grid shape calculation process, has strong dependence on user interaction, and may need additional interaction times to ensure that a better segmentation result is obtained. The existing volume data segmentation method based on contour editing adopted by TurtleSeg and other software selects a planar slice for volume data, and a user marks a contour line of a segmentation target on the slice as a constraint condition for grid generation. And performing curved surface reconstruction on the result obtained by marking the plurality of contour lines to obtain a mesh model of the segmentation result. The method has the defects of accuracy and efficiency, and more contours need to be marked to ensure the accuracy of the details of the segmentation result. Related methods can be found in Ijiri T, Yokota H.content-based interface for refining volume segment [ J ]. Computer Graphics Forum,2010,29(7), 2153-.

The existing volume data interactive segmentation method based on two-dimensional slice interaction needs an operator to mark a large number of slices to obtain a good result, the interaction efficiency is low, and the time required for segmenting one piece of data can be as long as several hours. In the grid-based segmentation method, the algorithm has a low utilization rate of volume data information and depends heavily on user interaction, and due to the lack of volume data as constraint, only the grid shape is considered, so that the segmentation precision is low. The curved surface reconstruction method based on the two-dimensional contour line is difficult to accurately segment the fine structure and needs more user interaction times.

Disclosure of Invention

The invention aims to provide a method for carrying out human-computer interaction editing on a segmentation result of volume data, so as to overcome the technical problems that the advantages of three-dimensional interaction and two-dimensional interaction cannot be combined and the volume data information cannot be fully utilized in the prior art, and further improve the interaction efficiency and the segmentation precision.

The invention provides a method for carrying out human-computer interaction editing on a segmentation result of volume data, which comprises the following steps:

(1) reading the segmentation result of the volume data stored in the form of binary data from the computer, and performing isosurface extraction on the segmentation result by using a cube stepping method to obtain a triangular mesh model M of the segmentation result₀Triangular mesh model M₀Representing the shape and volume characteristics of the segmented object, the target value of the extracted isosurface in the cube stepping method is 0.5, and the triangle is usedMesh model M₀Storing in a mode of vertex coordinates and triangular patch indexes; adopting a graphical processing unit rasterization method to carry out the triangular mesh model M₀Performing surface drawing, and recording the volume data coordinates of the positions of all the surface patches in the triangular mesh to obtain a surface drawing two-dimensional image;

(2) reading volume data corresponding to the segmentation result of the step (1) from a computer, respectively performing volume rendering on the volume data from a spatial position corresponding to each computer screen pixel by using a ray projection method according to the position of an observer, namely, emitting rays to the volume data along the sight line direction to obtain projected rays, judging the rays projected onto the volume data according to the volume data coordinates of the position of each patch in the triangular grid of the step (1), if the rays are projected onto the triangular grid patch in the step (1), ending the ray projection, if the rays are not projected onto the patch, intersecting the projected rays and the volume data to obtain an intersecting line segment, obtaining a plurality of sampling points at equal intervals on the intersecting line segment, reading the volume data gray value at the sampling points, and accumulating the volume data gray value along the ray projection direction, obtaining a brightness value accumulated in the transmission process of the projection light in the volume data, and obtaining a color value of each screen pixel according to the brightness value to obtain a two-dimensional image drawn by the volume;

(3) in the same view, mixing the surface rendering two-dimensional image in the step (1) and the volume rendering two-dimensional image in the step (2) by adopting a transparency mixing method to obtain a unified display image;

(4) recording a pen touch line for mouse interaction of a user in the unified display image in the step (3), storing the pen touch line in a coordinate point sequence form of screen pixels, sequentially connecting every two coordinate points in the pen touch line to form a broken line segment, extending the screen pixels corresponding to all points on the broken line segment in the direction of sight of an observer by using the light projection method in the step (2) to obtain a plurality of parallel rays, intersecting the plurality of parallel rays with the volume data in the step (2) to obtain a plurality of parallel line segments, forming a quadrangle by two adjacent parallel line segments, and enabling the quadrangle to be located according to the pen touch lineThe punctuation sequences are connected in sequence to obtain an extension surface, and the extension surface is stored in a triangular mesh form and is recorded as a section mesh M_SMaking the section grid M_SIntersecting the volume data and meshing M according to the section_SThe grid vertex coordinates obtain corresponding volume data coordinates, further obtain the gray level of the volume data, and grid M on the section_SDisplaying a gray scale of the volume data;

(5) the section grid M obtained in the step (4)_SAnd (2) performing the step (1) on the triangular mesh model M₀One triangular patch is taken from each triangular patch, two triangles corresponding to the two triangular patches are subjected to space intersection test to obtain three-dimensional space coordinates and barycentric coordinates of two end points of two intersected triangular patch cross-section in space, all triangular patch selection modes in two triangular meshes are traversed to obtain all cross-section lines of the triangular patches in the two triangular meshes, all cross-section lines are connected into continuous broken line segments according to space adjacent relations and are recorded as cross-section meshes M_SAnd triangular mesh model M₀Cross line C of₀；

(6) The section grid M of the step (4) is processed_SA plurality of quadrangles connected in sequence are arranged on a two-dimensional plane according to the stroke line coordinate point sequence in the step (4) to obtain a section grid P on the two-dimensional plane_SRasterizing the two-dimensional cross-section grid P by adopting a graphic processing unit_SSurface drawing is carried out to obtain volume data gray scale on a section grid P_SThe display result image of (1); according to the above cross-section mesh P_SSetting the intersecting line C of the step (5)₀In the cross section of the grid P_SThe barycentric coordinates on the middle triangular patch are not changed to obtain an intersection line C on a two-dimensional plane_P(ii) a Making the above two-dimensional section mesh P_SIntersecting the two dimensions C_PThe views are displayed together in the same view and are recorded as a two-dimensional view;

(7) the visual boundary of the volume data gray scale on the two-dimensional section grid in the two-dimensional view of the step (6) and the two-dimensional intersection line C in the two-dimensional view_PComparing the positions of the two lines, and if the intersection line deviates from the above steps(1) When the boundary of the object is divided, the point on the intersecting line is moved by using the mouse, and the point on the intersecting line is moved to the boundary of the divided object to obtain the intersecting line C_NAnd obtaining an intersecting line C_NThe barycentric coordinates of the points above;

(8) the intersection line C according to the above step (7)_NBarycentric coordinates of points thereon and the cross-sectional mesh M in the above-described step (4)_SObtaining the intersecting line C by the grid vertex coordinates_NThe three-dimensional space coordinate C of the upper point;

(9) establishing a triangular mesh optimization model based on volume data gradient, wherein an objective function E of the triangular mesh optimization model is as follows:

E＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(V′)

wherein V' is an independent variable of the objective function, namely the optimized triangular mesh model M₀The vertex coordinates are solved to obtain V' by minimizing an objective function E,

‖L(V′)-T(V′)δ₀‖²is a laplace shape constraint, where L is the laplace coordinate operator,

l (V') is an optimized triangular mesh model M₀Laplace coordinate of (d)₀Is an initial grid M₀T (V') is a local transformation of vertices represented by a linear combination of vertex coordinates;

ω‖MV′-C‖²for user interaction constraint, ω is a weight parameter, and C is the intersection line C obtained in the above step (8)_NThe three-dimensional space coordinates C, M of the upper point are parameter matrixes, namely the intersection line C_NThe three-dimensional space coordinate C of the upper point is expressed as a triangular mesh model M₀Linear combination of the vertex coordinates is carried out, and the linear combination is expressed into a matrix form, so that a parameter matrix M can be obtained;

- κ G (V ') is a volume data gradient constraint, κ is a weight parameter of the term, G is a volume data gradient function, G (V ') represents a volume data gradient magnitude at the vertex V ', the gradient function G is numerically calculated using a difference;

(10) rewriting the optimization function E of the above step (9) to E':

E′＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(P)+θ||V′-P||²

wherein P is an auxiliary variable and represents a triangular mesh model M in the optimization process₀Theta is a relaxation variable representing | V' -P |²The weight coefficient of (2) is gradually increased in the optimization process;

(11) extracting the related terms of E ' in the step (10) and taking V ' as a variable to obtain an optimization sub-problem S1, wherein the optimization sub-problem S1 is a known auxiliary variable P, and solving a vertex coordinate V ', and an optimization function is as follows:

the optimization subproblem S1 is a linear equation, which is solved by using a linear least square method to obtain V ', and V' is updated;

(12) extracting the relevant items of E 'in the step (10) by taking P as a variable to obtain an optimization sub-problem S2, wherein the optimization sub-problem S2 is a known vertex coordinate V', and the auxiliary variable P is solved, and the optimization function is as follows:

the optimization sub-problem S2 is a gradient optimization problem, the objective function being independent for each vertex, for the triangular mesh M₀Solving the gradient optimization problem by using a neighborhood searching method to traverse the triangular mesh M₀Obtaining P and updating P at the top point of the graph;

(13) setting the maximum value of the relaxation variable theta to theta_maxSetting the initial value of the relaxation variable theta to theta at initialization₀Taking V' ═ V₀In which V is₀A triangular mesh M of the step (1)₀Then the iteration starts:

(13-1) let P be V', traverse the triangular mesh M₀Performing neighborhood search on the vertex of the P, solving the optimization subproblem S2 to obtain and update P;

(13-2) solving the optimization sub-problem S1 by using a linear least square method according to the P in the step (13-1), and obtaining and updating V';

increasing the value of theta, repeating the steps (13-1) and (13-2), judging the relaxation variable theta, and if theta is larger than or equal to theta_maxStopping iteration and carrying out the step (14) if theta is less than theta_maxAnd then returning to (13-1);

(14) and (5) displaying a triangular mesh corresponding to the vertex coordinates V' in the step (13), wherein the triangular mesh is a result obtained by carrying out one-time human-computer interaction editing on the volume data segmentation result.

The method for carrying out human-computer interaction editing on the segmentation result of the volume data, provided by the invention, has the advantages that:

according to the method, the construction body drawing and the surface drawing are combined and are uniformly displayed in an interactive environment, so that the segmentation target is convenient to identify; the segmentation method based on stroke drawing is used for solving the intersecting lines of the segmentation results, so that the segmentation method can adapt to the segmentation targets in different shapes, and the method is favorable for more intuitively confirming the problems in the edited segmentation results; the three-dimensional slice is flattened to generate a two-dimensional slice, so that interaction can be performed more accurately, and the interactive editing precision is improved; the two-dimensional editing result is transmitted to the three-dimensional by the grid optimization algorithm based on the volume data gradient, so that the number of user interaction is saved, and meanwhile, the partition grid boundary is close to the volume data boundary, so that the partition precision is improved.

Detailed Description

The invention provides a method for carrying out human-computer interaction editing on a segmentation result of volume data, which comprises the following steps:

(1) reading the segmentation result of the volume data stored in the form of binary data from the computer, and performing isosurface extraction on the segmentation result by using a cube stepping (MarchingCubes) method to obtainTriangular mesh model M to segmentation result₀Triangular mesh model M₀The method comprises the steps of representing the shape and the volume characteristics of a segmentation target, wherein the target value of an extracted isosurface in a cube stepping method is 0.5, in the process, an input volume data segmentation result is stored in a binary data form, 0 or 1 is used for representing whether the position belongs to the segmentation target or not at each position, if 0 represents that the position does not belong to the segmentation target, if 1 represents that the position belongs to the segmentation target, all three-dimensional objects composed of 1 positions are the segmentation results represented by the segmentation data, and therefore the target value of the extracted isosurface in the cube stepping method is set to be 0.5. Using the triangular mesh model M₀Storing in a mode of vertex coordinates and triangular patch indexes; rasterizing the triangular mesh model M by a Graphic Processing Unit (GPU)₀And (4) performing surface drawing, wherein in the surface drawing process, the depth cache of the GPU equipment is used for shielding and removing, and a correct surface drawing result is obtained. In addition, in the surface drawing process, volume data coordinates of the positions of all the surface patches in the triangular mesh are recorded; and obtaining a surface drawing two-dimensional image which can be displayed on a screen.

(2) Reading volume data corresponding to the segmentation result of the step (1) from a computer, respectively performing volume rendering on the volume data from a spatial position corresponding to each computer screen pixel by using a ray projection method according to the position of an observer, namely, emitting rays to the volume data along the sight line direction to obtain projected rays, judging the rays projected onto the volume data according to the volume data coordinates of the position of each patch in the triangular grid of the step (1), if the rays are projected onto the triangular grid patch in the step (1), ending the ray projection, if the rays are not projected onto the patch, intersecting the projected rays and the volume data to obtain an intersecting line segment, obtaining a plurality of sampling points at equal intervals on the intersecting line segment, reading the volume data gray value at the sampling points, and accumulating the volume data gray value along the ray projection direction, obtaining a brightness value accumulated in the transmission process of the projection light in the volume data, and obtaining a color value of each screen pixel according to the brightness value to obtain a two-dimensional image drawn by the volume; the two-dimensional image may be displayed on a screen.

(3) In the same view, mixing the surface rendering two-dimensional image in the step (1) and the volume rendering two-dimensional image in the step (2) by adopting a transparency mixing method to obtain a unified display image;

(4) and (4) recording a pen touch line for mouse interaction of a user in the unified display image in the step (3), storing the pen touch line in a coordinate point sequence form of screen pixels, sequentially connecting every two coordinate points in the pen touch line to form a broken line segment, wherein the broken line segment represents a part needing to be edited in a segmentation result by observing the unified display image, and intercepting the section shape of the three-dimensional volume data according to the shape of the part to be edited. Using the light projection method in the step (2), extending the screen pixels corresponding to all points on the broken line segment along the sight line direction of the observer to the inner direction of the screen to obtain a plurality of parallel rays, intersecting the plurality of parallel rays with the volume data in the step (2) to obtain a plurality of parallel line segments, forming a quadrangle by two adjacent parallel line segments, sequentially connecting the quadrangle according to the coordinate point sequence of the stroke line to obtain an extension surface, storing the extension surface in the form of triangular grids and marking as a cross-section grid M_SMaking the section grid M_SIntersecting the volume data and meshing M according to the section_SThe grid vertex coordinates obtain corresponding volume data coordinates, further obtain the gray level of the volume data, and grid M on the section_SDisplaying a gray scale of the volume data; the gray scale condition of the disassembled data can be observed through the observation of the cross section grids, and the segmentation target corresponding to the segmentation result can be distinguished through the gray scale.

(5) The section grid M obtained in the step (4)_SAnd (2) performing the step (1) on the triangular mesh model M₀Respectively taking one triangular patch, carrying out space intersection test on two triangles corresponding to the two triangular patches to obtain three-dimensional space coordinates and barycentric coordinates of two end points of two intersected triangular patch sections in space, traversing all triangular patch selection modes in the two triangular meshes to obtain all intersected triangular patch sections in the two triangular meshes, and carrying out space intersection test on all the intersected patch sections according to the intersection line sectionsThe adjacent space relations are connected into a continuous broken line segment which is marked as a section grid M_SAnd triangular mesh model M₀Cross line C of₀；

(6) The section grid M of the step (4) is processed_SA plurality of quadrangles connected in sequence are arranged on a two-dimensional plane according to the stroke line coordinate point sequence in the step (4) to obtain a section grid P on the two-dimensional plane_SRasterizing the two-dimensional cross-section grid P by adopting a graphic processing unit_SSurface drawing is carried out to obtain volume data gray scale on a section grid P_SThe display result image of (1); according to the above cross-section mesh P_SSetting the intersecting line C of the step (5)₀In the cross section of the grid P_SThe barycentric coordinates on the middle triangular patch are not changed to obtain an intersection line C on a two-dimensional plane_P(ii) a Making the above two-dimensional section mesh P_SIntersecting the two dimensions C_PThe views are displayed together in the same view and are recorded as a two-dimensional view;

(7) the visual boundary of the volume data gray scale on the two-dimensional section grid in the two-dimensional view of the step (6) and the two-dimensional intersection line C in the two-dimensional view_PComparing the positions, when the intersecting line deviates from the boundary of the segmentation target in the step (1), moving the point on the intersecting line by using a mouse, and moving the point on the intersecting line to the boundary of the segmentation target to obtain an intersecting line C_NAnd obtaining an intersecting line C_NThe barycentric coordinates of the points above;

(8) the two-dimensional intersection line C of the step (7) is formed_PTo the new intersection line C_NIs converted into the position transformation of the three-dimensional mesh vertex. The method is according to the intersecting line C of the step (7) above_NBarycentric coordinates of points thereon and the cross-sectional mesh M in the above-described step (4)_SObtaining the intersecting line C by the grid vertex coordinates_NThe three-dimensional space coordinate C of the upper point;

(9) establishing a triangular mesh optimization model based on volume data gradient for transmitting the two-dimensional interactive editing result in the step (7) to three-dimensional, so that the two-dimensional interactive editing operation can drive the three-dimensional mesh M₀And deforming, thereby completing the editing of the volume data segmentation result. TheThe objective function E of the triangular mesh optimization model is as follows:

E＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(V′)

wherein V' is an independent variable of the objective function, namely the optimized triangular mesh model M₀The vertex coordinates are solved to obtain V' by minimizing an objective function E,

‖L(V′)-T(V′)δ₀‖²is a laplace shape constraint, where L is the laplace coordinate operator,

l (V') is an optimized triangular mesh model M₀Laplace coordinate of (d)₀Is an initial grid M₀T (V') is a local transformation of vertices represented by a linear combination of vertex coordinates; to solve the problem that laplace coordinates are sensitive to rotation scaling transformations.

ω‖MV′-C‖²For user interaction constraint, ω is a weight parameter, and the value thereof is adjusted empirically, which in one embodiment of the present invention is 0.5, and C is the intersection line C obtained in the above step (8)_NThe three-dimensional space coordinates C, M of the upper point are parameter matrixes, namely the intersection line C_NThe three-dimensional space coordinate C of the upper point is expressed as a triangular mesh model M₀Linear combination of the vertex coordinates is carried out, and the linear combination is expressed into a matrix form, so that a parameter matrix M can be obtained;

- κ G (V ') is a volume data gradient constraint, κ is a weight parameter of the term, the value of which is empirically adjusted, and in one embodiment of the invention is 1.0, G is a volume data gradient function, G (V ') represents the volume data gradient magnitude at the vertex V ', and the gradient function G is numerically calculated using a difference;

(10) rewriting the optimization function E of the above step (9) to E':

E′＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(P)+θ||V′-P||²

wherein P is an auxiliary variable and represents a triangular mesh model M in the optimization process₀Theta is a relaxation variable representing | V' -P |²The weight coefficient of (2) is gradually increased in the optimization process;

the use of the auxiliary variable P in the above equation replaces V ' in the non-convex function G (V '), with the addition of | V ' -P |²The terms are constrained so that the auxiliary variable P is as close as possible to V' being replaced. When θ → + ∞, there is E' → E;

(11) extracting the related terms of E ' in the step (10) and taking V ' as a variable to obtain an optimization sub-problem S1, wherein the optimization sub-problem S1 is a known auxiliary variable P, and solving a vertex coordinate V ', and an optimization function is as follows:

the optimization subproblem S1 is a linear equation, which is solved by using a linear least square method to obtain V ', and V' is updated;

(12) extracting the relevant items of E 'in the step (10) by taking P as a variable to obtain an optimization sub-problem S2, wherein the optimization sub-problem S2 is a known vertex coordinate V', and the auxiliary variable P is solved, and the optimization function is as follows:

the optimization sub-problem S2 is a gradient optimization problem, the objective function being independent for each vertex, for the triangular mesh M₀Solving the gradient optimization problem by using a neighborhood searching method to traverse the triangular mesh M₀Obtaining P and updating P at the top point of the graph; selecting some points in a certain size of area near each vertex, and calculating the objective function of the optimization subproblem S2 to obtain a vertex coordinate with the minimum objective function value;

(13) the optimization sub-problem S1 in step (11) is iterated over the optimization sub-problem S2 in step (12).

Setting the maximum value of the relaxation variable theta to theta_maxSetting the initial value of the relaxation variable theta to theta at initialization₀Taking V' ═ V₀In which V is₀A triangular mesh M of the step (1)₀Then the iteration starts:

(13-1) let P be V', traverse the triangular mesh M₀Performing neighborhood search on the vertex of the P, solving the optimization subproblem S2 to obtain and update P;

(13-2) solving the optimization sub-problem S1 by using a linear least square method according to the P in the step (13-1), and obtaining and updating V';

increasing the value of theta, repeating the steps (13-1) and (13-2), judging the relaxation variable theta, and if theta is larger than or equal to theta_maxStopping iteration and carrying out the step (14) if theta is less than theta_maxAnd then returning to (13-1);

(14) and (5) displaying a triangular mesh corresponding to the vertex coordinates V' in the step (13), wherein the triangular mesh is a result obtained by carrying out one-time human-computer interaction editing on the volume data segmentation result.

Claims (1)

1. A method for carrying out human-computer interaction editing on a segmentation result of volume data is characterized by comprising the following steps:

(1) reading the segmentation result of the volume data stored in the form of binary data from the computer, and performing isosurface extraction on the segmentation result by using a cube stepping method to obtain a triangular mesh model M of the segmentation result₀Triangular mesh model M₀Representing the shape and volume characteristics of the segmentation target, extracting the isosurface with a target value of 0.5 in the cube stepping method, and modeling the triangular mesh model M₀Storing in a mode of vertex coordinates and triangular patch indexes; adopting a graphical processing unit rasterization method to carry out the triangular mesh model M₀Performing surface drawing, and recording the volume data coordinates of the positions of all the surface patches in the triangular mesh to obtain a surface drawing two-dimensional image;

(2) reading volume data corresponding to the segmentation result of the step (1) from a computer, respectively performing volume rendering on the volume data from a spatial position corresponding to each computer screen pixel by using a ray projection method according to the position of an observer, namely, emitting rays to the volume data along the sight line direction to obtain projected rays, judging the rays projected onto the volume data according to the volume data coordinates of the position of each patch in the triangular grid of the step (1), if the rays are projected onto the triangular grid patch in the step (1), ending the ray projection, if the rays are not projected onto the patch, intersecting the projected rays and the volume data to obtain an intersecting line segment, obtaining a plurality of sampling points at equal intervals on the intersecting line segment, reading the volume data gray value at the sampling points, and accumulating the volume data gray value along the ray projection direction, obtaining a brightness value accumulated in the transmission process of the projection light in the volume data, and obtaining a color value of each screen pixel according to the brightness value to obtain a two-dimensional image drawn by the volume;

(3) in the same view, mixing the surface rendering two-dimensional image in the step (1) and the volume rendering two-dimensional image in the step (2) by adopting a transparency mixing method to obtain a unified display image;

(4) recording a pen touch line for mouse interaction of a user in the unified display image in the step (3), storing the pen touch line in a coordinate point sequence form of screen pixels, sequentially connecting every two coordinate points in the pen touch line to form a broken line segment, using the light projection method in the step (2), extending the screen pixels corresponding to all points on the broken line segment along the sight direction of an observer to the inner direction of the screen to obtain a plurality of parallel rays, intersecting the plurality of parallel rays with the volume data in the step (2) to obtain a plurality of parallel line segments, forming a quadrangle by two adjacent parallel line segments, sequentially connecting the quadrangle according to the coordinate point sequence of the pen touch line to obtain an extension surface, storing the extension surface in a triangular grid form, and recording the extension surface as a cross-section grid M_SMaking the section grid M_SIntersecting the volume data and meshing M according to the section_SThe grid vertex coordinates obtain corresponding volume data coordinates, further obtain the gray level of the volume data, and grid M on the section_SUpper displayThe gray scale of the volume data;

(5) the section grid M obtained in the step (4)_SAnd (2) performing the step (1) on the triangular mesh model M₀One triangular patch is taken from each triangular patch, two triangles corresponding to the two triangular patches are subjected to space intersection test to obtain three-dimensional space coordinates and barycentric coordinates of two end points of two intersected triangular patch cross-section in space, all triangular patch selection modes in two triangular meshes are traversed to obtain all cross-section lines of the triangular patches in the two triangular meshes, all cross-section lines are connected into continuous broken line segments according to space adjacent relations and are recorded as cross-section meshes M_SAnd triangular mesh model M₀Cross line C of₀；

(6) The section grid M of the step (4) is processed_SA plurality of quadrangles connected in sequence are arranged on a two-dimensional plane according to the stroke line coordinate point sequence in the step (4) to obtain a section grid P on the two-dimensional plane_SRasterizing the two-dimensional cross-section grid P by adopting a graphic processing unit_SSurface drawing is carried out to obtain volume data gray scale on a section grid P_SThe display result image of (1); according to the above cross-section mesh P_SSetting the intersecting line C of the step (5)₀In the cross section of the grid P_SThe barycentric coordinates on the middle triangular patch are not changed to obtain an intersection line C on a two-dimensional plane_P(ii) a Making the above two-dimensional section mesh P_SIntersecting the two dimensions C_PThe views are displayed together in the same view and are recorded as a two-dimensional view;

(7) the visual boundary of the volume data gray scale on the two-dimensional section grid in the two-dimensional view of the step (6) and the two-dimensional intersection line C in the two-dimensional view_PComparing the positions, when the intersecting line deviates from the boundary of the segmentation target in the step (1), moving the point on the intersecting line by using a mouse, and moving the point on the intersecting line to the boundary of the segmentation target to obtain an intersecting line C_NAnd obtaining an intersecting line C_NThe barycentric coordinates of the points above;

(8) the intersection line C according to the above step (7)_NBarycentric coordinates of points thereon and the cross-sectional net in the above step (4)Grid M_SObtaining the intersecting line C by the grid vertex coordinates_NThe three-dimensional space coordinate C of the upper point;

(9) establishing a triangular mesh optimization model based on volume data gradient, wherein an objective function E of the triangular mesh optimization model is as follows:

E＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(V′)

wherein V' is an independent variable of the objective function, namely the optimized triangular mesh model M₀The vertex coordinates are solved to obtain V' by minimizing an objective function E,

||L(V′)-T(V′)δ₀||²is a laplace shape constraint, where L is the laplace coordinate operator,

l (V') is an optimized triangular mesh model M₀Laplace coordinate of (d)₀Is an initial grid M₀T (V') is a local transformation of vertices represented by a linear combination of vertex coordinates;

ω||MV′-C||²for user interaction constraint, ω is a weight parameter, and C is the intersection line C obtained in the above step (8)_NThe three-dimensional space coordinates C, M of the upper point are parameter matrixes, namely the intersection line C_NThe three-dimensional space coordinate C of the upper point is expressed as a triangular mesh model M₀Linear combination of the vertex coordinates is carried out, and the linear combination is expressed into a matrix form, so that a parameter matrix M can be obtained;

- κ G (V ') is a volume data gradient constraint, κ is a weight parameter of the term, G is a volume data gradient function, G (V ') represents a volume data gradient magnitude at the vertex V ', the gradient function G is numerically calculated using a difference;

(10) rewriting the optimization function E of the above step (9) to E':

E′＝||L(V′)-T(V′)δ₀||²+ω||MV′-C||²-κG(P)+θ||V′-P||²

wherein P is an auxiliary variable and represents a triangular mesh model M in the optimization process₀The vertex coordinate of (a) is a relaxation variable and represents | V' -P | Y²The weight coefficient of (2) is gradually increased in the optimization process;

(11) extracting the related terms of E ' in the step (10) and taking V ' as a variable to obtain an optimization sub-problem S1, wherein the optimization sub-problem S1 is a known auxiliary variable P, and solving a vertex coordinate V ', and an optimization function is as follows:

the optimization subproblem S1 is a linear equation, which is solved by using a linear least square method to obtain V ', and V' is updated;

(12) extracting the relevant items of E 'in the step (10) by taking P as a variable to obtain an optimization sub-problem S2, wherein the optimization sub-problem S2 is a known vertex coordinate V', and the auxiliary variable P is solved, and the optimization function is as follows:

the optimization sub-problem S2 is a gradient optimization problem, the objective function being independent for each vertex, for the triangular mesh M₀Solving the gradient optimization problem by using a neighborhood searching method to traverse the triangular mesh M₀Obtaining P and updating P at the top point of the graph;

(13) setting the maximum value of the relaxation variable theta to theta_maxSetting the initial value of the relaxation variable theta to theta at initialization₀Taking V' ═ V₀In which V is₀A triangular mesh M of the step (1)₀Then the iteration starts:

(13-1) let P be V', traverse the triangular mesh M₀Performing neighborhood search on the vertex of the P, solving the optimization subproblem S2 to obtain and update P;

(13-2) solving the optimization sub-problem S1 by using a linear least square method according to the P in the step (13-1), and obtaining and updating V';

increasing the value of theta, repeating the steps (13-1) and (13-2), judging the relaxation variable theta, and if theta is larger than or equal to theta_maxStopping iteration and carrying out the step (14) if theta is less than theta_maxAnd then returning to (13-1);

(14) and (5) displaying a triangular mesh corresponding to the vertex coordinates V' in the step (13), wherein the triangular mesh is a result obtained by carrying out one-time human-computer interaction editing on the volume data segmentation result.