A preoperative virtual reconstruction system for thoracic surgery

Abstract

The present application relates to the technical field of three-dimensional reconstruction and intelligent analysis of medical images, and specifically discloses a preoperative virtual reconstruction system for thoracic surgery, which constructs a spatial graph structure with voxels as nodes and introduces respiratory motion displacement trajectories as dynamic attributes; extracts and fuses static structural features and dynamic motion features to generate context feature vectors rich in anatomical and functional information; in a serialized topological constraint segmentation module, the graph nodes are sorted according to topological geodetic distance and feature saliency, and hierarchical segmentation from trunk to periphery is realized by combining anatomical topological constraints based on a rule base and a sequence prediction model; connectivity test, fracture repair and geometric reasonableness verification are performed to output a high-quality preoperative virtual reconstruction three-dimensional model that is accurate, error-free and consistent with anatomical rules.

Description

Technical Field

[0001] This invention relates to the field of three-dimensional reconstruction and intelligent analysis technology of medical images, specifically to a preoperative virtual reconstruction system for thoracic surgery. Background Technology

[0002] In preoperative planning for thoracic surgery, accurate 3D reconstruction of organs and lesions based on computed tomography (CT) images is crucial. Existing techniques typically employ segmentation algorithms such as graph cut, which achieve tissue separation by mapping images to a graph structure and minimizing an energy function. However, when dealing with complex regions with highly overlapping grayscale values ​​and tightly adhered structures, these methods heavily rely on image grayscale gradients and spatial smoothness constraints, often making it difficult to accurately distinguish the boundaries of different tissues, resulting in limited segmentation accuracy.

[0003] The existing technology has the following shortcomings:

[0004] When dealing with adhesions exhibiting similar grayscale values ​​and complex topologies, energy-minimizing graph cut algorithms may exhibit a mathematical "short-circuiting" effect in their smoothing constraints. This leads them to prioritize segmentation paths that traverse a few high-contrast regions rather than actual tortuous boundaries, resulting in non-physical topological breaks in tubular structures such as blood vessels in the 3D reconstruction model. While this type of error is less common in simple, routine cases, it is particularly prominent in extreme scenarios such as "deadlock" adhesions in the hilar region. The resulting topologically distorted model can severely mislead surgical planning and even lead to incorrect surgical procedures. Traditional methods lack explicit mechanisms for preserving anatomical topological connectivity, making it impossible to effectively mitigate this risk. Summary of the Invention

[0005] The purpose of this invention is to provide a virtual reconstruction system for preoperative thoracic surgery to solve the problems mentioned above.

[0006] The objective of this invention can be achieved through the following technical solutions:

[0007] A virtual reconstruction system for preoperative thoracic surgery includes:

[0008] The multimodal data fusion and dynamic graph construction module is used to simultaneously acquire chest CT image data and respiratory CT image data; an initial graph structure is constructed based on the chest CT image data, voxel units are defined as graph nodes, and spatially adjacent nodes are connected to form edges; the displacement trajectory of each voxel unit in respiratory motion is extracted from the respiratory CT image data and added to the corresponding graph node as a dynamic attribute;

[0009] The cross-modal context-aware feature encoding module extracts the local 3D image block and corresponding displacement trajectory sequence of each image node in the CT image data; processes the image block and trajectory sequence separately to obtain static structural features and dynamic motion features; integrates the static structural features and dynamic motion features to generate a fused context feature vector for each image node.

[0010] The serialization topology constraint segmentation module sorts graph nodes according to the fused context feature vector to form a processing sequence. Based on the processing sequence, a sequence prediction model is used to process each graph node in turn. The sequence prediction model performs label prediction by combining the current graph node features, the state of the processed neighboring graph nodes, and a predefined anatomical topology rule base. The rule base specifies the morphological and connectivity constraints of key anatomical structures.

[0011] The segmentation result verification and model reconstruction module maps the labels obtained from sequence prediction to three-dimensional space to generate an initial segmentation mask; it performs global anatomical rationality verification and correction on the initial segmentation mask, and outputs a preoperative virtual reconstructed three-dimensional model for preoperative clinical use in thoracic surgery.

[0012] As a further aspect of the present invention: the construction of the initial image structure based on chest CT image data specifically includes:

[0013] The spatial coordinates of each voxel unit in the CT image data are mapped to graph nodes. Based on the preset 26 adjacency rules, all voxel units are traversed. If two voxel units are connected by a face, edge or point in three-dimensional space, an undirected edge is established between the graph nodes corresponding to the two voxel units.

[0014] Calculate the gradient magnitude of image grayscale values ​​between adjacent voxel units; assign initial weights to each edge based on the gradient magnitude of image grayscale values;

[0015] Centered on each graph node, image patch features at different scales are extracted; the image patch features are fused with the initial weights to generate enhanced edge attributes, thus completing the optimized construction of the initial graph structure.

[0016] As a further aspect of the present invention: the extraction of the displacement trajectory of each voxel unit in respiratory motion from respiratory CT image data and its addition as a dynamic attribute to the corresponding graph node specifically includes:

[0017] The respiratory CT image data is divided into respiratory phases, and the end of inspiration or end of expiration phase is selected as the reference frame; the dense three-dimensional motion vector field of each phase image relative to the reference frame is calculated.

[0018] For each voxel unit, based on its spatial position in each phase image, a spatiotemporal displacement trajectory spanning a complete respiratory cycle is formed by interpolating its displacement vector in the continuous phase motion vector field.

[0019] Three feature quantities are extracted from the displacement trajectory: trajectory length, average displacement amplitude, and main direction of motion. These three features together form a trajectory feature vector. The trajectory feature vector is then used as a dynamic attribute and attached to the corresponding graph node in the initial graph structure.

[0020] As a further aspect of the present invention: the processing of image blocks and trajectory sequences to obtain static structural features and dynamic motion features specifically includes:

[0021] For local 3D image blocks, a set of directionally separable 3D convolutional kernels are used for feature extraction. The 3D convolutional kernels have different weights along the axial, sagittal, and coronal directions to anisotropically perceive structural information of different anatomical planes and output multi-directional feature maps.

[0022] For the displacement trajectory sequence, its displacement vector is decomposed into three motion primitive components: translation, rotation, and scaling; the amplitude and phase change spectrum of each motion primitive component in each respiratory cycle is calculated to form a dynamic motion feature vector;

[0023] Principal component analysis was performed to reduce the dimensionality of the multi-directional feature maps and dynamic motion feature vectors, and the reduced features were mapped to a unified numerical range to obtain standardized static structural features and dynamic motion features.

[0024] As a further aspect of the present invention: the integration of static structural features and dynamic motion features to generate a fused context feature vector for each graph node specifically includes:

[0025] Construct a fully connected temporary graph using all graph nodes in the current batch as vertices; use the static structural features and dynamic motion features of each node as the two attributes of the corresponding vertex;

[0026] Calculate the attention weight between any two vertices based on their static structural features similarity; use the attention weight to weighted aggregate the dynamic motion features of the target vertex, and generate a dynamic feature enhanced with structural context for each vertex;

[0027] The original static structural features and the dynamic features enhanced by structural context are input into a gated fusion unit. The gated fusion unit learns a dynamic weight to control the fusion ratio of the two types of features and finally outputs a unique fusion context feature vector for each graph node.

[0028] As a further aspect of the present invention: the step of sorting graph nodes based on the fused context feature vector to form a processing sequence specifically includes:

[0029] Based on the fusion context feature vector, graph nodes belonging to the starting point of the main trunk of the key anatomical structure are identified and used as the initial root node set;

[0030] Starting from the root node set, calculate the topological geodesic distance from all other graph nodes to the nearest root node on the initial graph structure; the topological geodesic distance refers to the shortest number of hops along the graph edge connection path;

[0031] All graph nodes are sorted in ascending order of topological geodesic distance. Graph nodes with the same distance are then sorted in a secondary order according to the magnitude of their fused context feature vectors. Finally, a processing sequence radiating from the core of the structure to the edge is generated.

[0032] As a further aspect of the present invention: the sequential prediction model is used to process each graph node in turn, specifically including:

[0033] Before predicting the label of the current graph node, query the label status of the processed neighboring graph nodes to determine whether they have triggered a specific rule;

[0034] If a rule is triggered, the attention weights in the sequence prediction model are adjusted according to the rule type.

[0035] The sequence prediction model integrates the features of the current node, the state context of the regulated neighboring nodes, and finally decides the label of the current node through a gated recurrent unit.

[0036] As a further aspect of the present invention: the global dissection rationality check and correction of the initial segmentation mask specifically includes:

[0037] Extract the topological skeleton network of key anatomical structures from the initial segmentation mask; calculate the number of connected components in the topological skeleton network and compare it with the number of normal connected components specified in the predefined anatomical rule base to identify whether there are abnormal disconnections;

[0038] For the identified abnormal breakage areas, determine the fracture endpoints and calculate the geometric characteristics of the three-dimensional spatial gap at the fracture point;

[0039] Using the fracture endpoints and their neighboring contour points as constraints, a smooth transition repair surface is generated within the fracture gap to restore the correct topological connectivity of the anatomical structure.

[0040] As a further aspect of the present invention: the output of the preoperative virtual reconstruction 3D model specifically includes:

[0041] The total volume and surface area of ​​the repaired 3D segmentation mask are calculated and compared with the range of typical organ geometric parameters recorded in the predefined anatomical topology rule base to verify its geometric rationality.

[0042] Multi-resolution analysis was performed on the initial surface mesh model obtained by the mask conversion after restoration. At the coarse resolution level, the unnatural surface undulations introduced by interpolation restoration were smoothed, while at the fine resolution level, the real anatomical details from the original image data were preserved.

[0043] The optimized curved mesh is subjected to watertightness processing and error-free patch inspection, and finally a preoperative virtual reconstructed 3D model is output. The preoperative virtual reconstructed 3D model can be directly used for 3D visualization, printing or surgical navigation.

[0044] The beneficial effects of this invention are:

[0045] (1) Through multimodal data fusion, dynamic graph construction, and cross-modal context-aware feature encoding, this system can comprehensively consider the static structural information of chest CT images and the dynamic characteristics of respiratory motion. This makes the segmentation of key anatomical structures more accurate, especially when structures that move with respiration are involved. The sequential topological constraint segmentation module, combined with an anatomical rule base, ensures that the segmentation results are not only accurate in local image features but also reasonable under global anatomical topological constraints. This high-precision segmentation and reconstruction helps surgeons to understand the patient's anatomy more accurately, thereby developing safer and more effective surgical plans.

[0046] (2) The preoperative virtual reconstruction 3D model generated by the system has a high-quality mesh structure that meets clinical use standards such as watertightness and error-free surface areas. It can be output in multiple formats, making it easy to integrate into different medical processes, such as 3D visualization software, 3D printing, or surgical navigation systems. This provides surgeons with an intuitive visual aid, supporting more detailed preoperative planning, simulating the surgical process, and providing real-time guidance during surgery. In addition, by performing global anatomical rationality checks and corrections on the segmentation results, the geometric rationality and connectivity of the final model are ensured, further improving its clinical applicability. Attached Figure Description

[0047] The invention will now be further described with reference to the accompanying drawings.

[0048] Figure 1 This is a flowchart of the system of the present invention. Detailed Implementation

[0049] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0050] Please see Figure 1As shown, the present invention is a virtual reconstruction system for preoperative thoracic surgery, comprising:

[0051] The multimodal data fusion and dynamic graph construction module is used to simultaneously acquire chest CT image data and respiratory CT image data; an initial graph structure is constructed based on the chest CT image data, voxel units are defined as graph nodes, and spatially adjacent nodes are connected to form edges; the displacement trajectory of each voxel unit in respiratory motion is extracted from the respiratory CT image data and added to the corresponding graph node as a dynamic attribute;

[0052] The cross-modal context-aware feature encoding module extracts the local 3D image block and corresponding displacement trajectory sequence of each image node in the CT image data; processes the image block and trajectory sequence separately to obtain static structural features and dynamic motion features; integrates the static structural features and dynamic motion features to generate a fused context feature vector for each image node.

[0053] The serialization topology constraint segmentation module sorts graph nodes according to the fused context feature vector to form a processing sequence. Based on the processing sequence, a sequence prediction model is used to process each graph node in turn. The sequence prediction model performs label prediction by combining the current graph node features, the state of the processed neighboring graph nodes, and a predefined anatomical topology rule base. The rule base specifies the morphological and connectivity constraints of key anatomical structures.

[0054] The segmentation result verification and model reconstruction module maps the labels obtained from sequence prediction to three-dimensional space to generate an initial segmentation mask; it performs global anatomical rationality verification and correction on the initial segmentation mask, and outputs a preoperative virtual reconstructed three-dimensional model for preoperative clinical use in thoracic surgery.

[0055] In the multimodal data fusion and dynamic graph construction module, data acquisition and graph structure initialization are performed. The system simultaneously acquires two types of computed tomography (CT) image datasets of the target object: first, conventional chest-enhanced CT image data, which provides high spatial resolution static anatomical information; second, four-dimensional dynamic respiratory CT image data, which captures the motion changes of the thoracic cavity structure throughout the respiratory cycle at continuous time points, resulting in chest CT image data and respiratory CT image data, respectively. The process of constructing the initial graph structure begins with the static enhanced CT image data. The system reads the digital imaging and communication in medicine format file of this data and parses out the three-dimensional volume data matrix. In this matrix, each smallest volume unit, i.e., voxel unit, is uniquely mapped to a graph node according to its index coordinates in the three-dimensional matrix. Subsequently, the system uses a 26-adjacency rule to define the connection relationships between graph nodes. Specifically, the algorithm traverses each voxel unit in three-dimensional space and checks its spatial relationship with its 26 directly adjacent voxel units (including all cases of face adjacency, edge adjacency, and point adjacency). For any two voxel units that satisfy the 26 adjacency condition and are both located within the effective image area, the system will establish an undirected edge between their corresponding two graph nodes, thus initially forming a graph structure skeleton with voxels as nodes and spatial adjacency relationships as edges.

[0056] After establishing the basic connectivity of the graph, each edge needs to be assigned an initial weight value, which should reflect the likelihood of the existence of tissue boundaries belonging to adjacent voxel units. The system calculates the difference in grayscale values ​​between two adjacent voxel units in the acquired computed tomography images, and further calculates the magnitude of this difference, i.e., the magnitude of the image grayscale gradient. This magnitude intuitively represents the degree of intensity change between adjacent voxels. Initial weights are assigned to each edge based on this magnitude, with the principle that when the grayscale values ​​of two adjacent voxel units are very close, i.e., the gradient magnitude is small, they are likely to belong to the same homogeneous tissue, so the edge between them is given a larger weight, indicating a tendency not to cut at this point during segmentation; conversely, if the gradient magnitude is large, it indicates the possible existence of tissue boundaries, and the edge is given a smaller weight. However, relying solely on local grayscale gradients may be affected by noise and cannot capture broader contextual information. Therefore, the system further performs multi-scale feature enhancement. Three-dimensional image patch features at different scales are extracted from each graph node. Smaller-scale image patches capture local details, while larger-scale image patches contain broader anatomical context information. The system processes these image patches through specific feature extraction operations (such as convolutional pooling) to obtain feature vectors representing information at different scales. Finally, these multi-scale feature vectors are fused with initial edge weights calculated based on gray-level gradients (e.g., through a small neural network or weighted concatenation) to generate enhanced edge attributes that simultaneously reflect local boundary information and regional context information, thereby completing the optimized construction of the initial graph structure.

[0057] The system processes four-dimensional dynamic respiratory computed tomography (CT) image data to extract dynamic attributes. First, the acquired image sequence at consecutive time points is divided into respiratory phases. This can be based on the motion curves of the diaphragm or specific anatomical points, dividing the entire sequence into several key phases such as inspiratory and expiratory phases. An image of a respiratory extreme phase (e.g., end of inspiration or end of expiration) is selected as a fixed reference frame. Then, for each phase image outside the reference frame, a non-rigid registration algorithm is used to calculate the three-dimensional motion displacement vector of each voxel position relative to its corresponding position in the reference frame. This process generates a dense three-dimensional motion vector field for each phase, describing the overall deformation from the reference frame to that phase frame. For each voxel unit in the three-dimensional volume data, the system, based on its registered spatial position in each phase image, uses algorithms such as trilinear interpolation to interpolate in the motion vector field of consecutive phases, tracing the voxel's motion path throughout the respiratory cycle, thus forming a spatiotemporal displacement trajectory describing its spatial position changing over time. To improve data reliability, the system performs a smoothing filter on the calculated original displacement trajectory based on respiratory mechanics priors to suppress abnormal jitter caused by acquisition noise and ensure that the trajectory is smooth and conforms to physiological movement laws.

[0058] From the smoothed displacement trajectory, the system quantifies and extracts several features that characterize its motion. These features include at least: trajectory length, i.e., the total length of the voxel's path throughout the entire respiratory cycle; average displacement amplitude, i.e., the average displacement distance relative to the reference frame position; and principal direction of motion, i.e., the main direction vector of the voxel's motion determined by methods such as principal component analysis. These features are combined to form a quantitative trajectory feature vector. This vector encapsulates the dynamic behavior characteristics of the voxel unit during respiration. Finally, the system attaches this trajectory feature vector as a key dynamic attribute to the corresponding graph node in the initial graph structure constructed in the first step.

[0059] In the cross-modal context-aware feature encoding module, for each map node, its local 3D image block and corresponding displacement trajectory sequence in the CT image data are extracted; the image block and trajectory sequence are processed separately to obtain static structural features and dynamic motion features; the static structural features and dynamic motion features are integrated to generate a fused context feature vector for each map node, specifically including:

[0060] Features that distinguish different anatomical structures are extracted from static computed tomography (CT) images. For each graph node (corresponding to a voxel unit) in the graph structure, a local 3D image block of a fixed size (e.g., extending a certain number of voxel distances in each of the three spatial directions from the voxel center) is extracted centered on its coordinates in 3D space. This 3D image block is then convolved using a set of pre-designed 3D convolutional kernels to extract features. The special feature of this set of convolutional kernels is their separable orientation; that is, for the three standard anatomical planes of the human body—axial, sagittal, and coronal—three 3D convolutional kernels with different weight parameters are used respectively. In practice, the system uses three independent sets of 3D convolutional kernels to perform convolution calculations on the input 3D image block. The axial convolutional kernel primarily senses structural changes on the transverse plane, the sagittal convolutional kernel primarily senses the structure of the lateral sections in the left-right direction, and the coronal convolutional kernel primarily senses the structure of the coronal plane in the anterior-posterior direction. Each convolution calculation involves sliding the convolution kernel across the image block, multiplying elements at corresponding locations, summing the results, adding a bias term, and finally applying a non-linear activation function (such as ReLU) to generate the kernel's feature response at that location. This sliding window calculation generates a feature map for each set of directional convolution kernels. Finally, these three feature maps are concatenated along the channel dimension to output a multi-directional feature map that integrates structural information from multiple directions. This process enables the system to anisotropically perceive structural information from different anatomical planes, thus more effectively capturing key anatomical features such as blood vessel cross-sections and the long axis of the bronchi.

[0061] Quantitative features describing tissue motion patterns are extracted from dynamic respiratory computed tomography (CT) data. For the same graph node, its corresponding displacement trajectory sequence is obtained, which consists of three-dimensional displacement vectors at multiple respiratory phase points. The system first decomposes each three-dimensional displacement vector into three basic motion primitive components: translation, rotation, and scaling. The translation component directly represents the displacement of the centroid; the rotation component is characterized by calculating the azimuth and elevation angles of the displacement vector relative to a reference coordinate system (e.g., with the voxel's position in the reference frame as the origin); and the scaling component is characterized by calculating the rate of change of the displacement vector's magnitude relative to the initial state. After completing the motion primitive decomposition at each time point, the system analyzes these primitive components over a complete respiratory cycle. For each motion primitive component (the amplitude of translation), its variation within the respiratory cycle is calculated, and its spectrum is analyzed (using the Fast Fourier Transform algorithm to calculate its spectrum), extracting its main frequency components (i.e., the phase change spectrum) and the corresponding amplitudes of these frequency components (i.e., the amplitude change spectrum). The key parameters (such as dominant frequency, dominant frequency amplitude, and spectral energy) in the amplitude and phase spectra of the three motion primitives (translation, rotation, and scaling) are extracted and sequentially combined into a high-dimensional vector, forming the "dynamic motion feature vector". This vector comprehensively describes the complex motion pattern of the voxel point during respiration.

[0062] Because multi-directional feature maps and dynamic motion feature vectors have high dimensionality and may contain redundant information, directly using them is not only computationally intensive but may also affect the performance of subsequent models. Therefore, the system performs dimensionality reduction on these two types of high-dimensional features. Principal component analysis (PCA) is used to achieve this. For the multi-directional feature map, the eigenvalues ​​of all spatial locations are first expanded into a one-dimensional long vector. Then, the covariance matrix of this vector set (from a batch of training samples) is calculated, and the eigenvalues ​​and eigenvectors of the covariance matrix are calculated. Next, the eigenvectors corresponding to the k largest eigenvalues ​​are selected as principal component directions in descending order of eigenvalues, with the value of k determined based on a preset variance contribution rate (e.g., 95%). Finally, the original high-dimensional feature vectors are projected onto these k principal component directions to obtain the dimensionality-reduced low-dimensional static structure feature vectors. The same operation is performed on the dynamic motion feature vectors to obtain the dimensionality-reduced low-dimensional dynamic motion feature vectors. To eliminate the impact of differences in physical dimensions and numerical ranges of different features, the system further normalizes the two types of feature vectors after dimensionality reduction. For example, the Z-score standardization method is used, which subtracts the mean of the feature over all samples and then divides it by its standard deviation, so that the processed feature values ​​are distributed in a relatively uniform range (e.g., the mean is 0 and the standard deviation is 1), thereby obtaining standardized static structural features and dynamic motion features.

[0063] To effectively integrate static structural information and dynamic motion information, the system constructs a temporary, fully connected undirected graph. The vertex set of this temporary graph is the set of graph nodes currently being processed. Each vertex has two attributes: a normalized static structural feature vector obtained in the previous step and a normalized dynamic motion feature vector. On this temporary graph, the system performs an attention-based feature interaction operation. First, the attention weight between any two vertices (denoted as vertex i and vertex j) is calculated. This weight is calculated based on the similarity of their static structural feature vectors. Specifically, the static structural feature vectors of vertices i and j are linearly transformed (i.e., multiplied by a learnable weight matrix), then the two transformed vectors are dot-producted. The dot-product result is then divided by a scaling factor (the square root of the feature vector dimension), and finally normalized using a Softmax function so that the sum of the attention weights of all edges connected to vertex i is 1. This weight value reflects the importance of the static structural information of vertex j for understanding vertex i. Then, using this attention weight, the dynamic motion feature vectors of all other vertices (j) are weighted and summed (also after a linear transformation) to generate a "structural context-enhanced dynamic feature" for vertex i. This new feature incorporates dynamic information from nodes with similar global static structures.

[0064] After obtaining the dynamic features enhanced by the structural context, they need to be fused with the original static structural features. The system uses a gated fusion unit to accomplish this task. The input to this unit is the original static structural feature vector and the enhanced dynamic feature vector obtained in the previous step. The core of the gated fusion unit is a learnable gating weight, which is a scalar value ranging from 0 to 1. It determines the contribution ratio of static and dynamic features in the final fusion result. The calculation process of the gating weight is as follows: first, the static feature vector and the enhanced dynamic feature vector are concatenated, and then the concatenated vector is input into a small feedforward neural network (containing a hidden layer and an output layer with a sigmoid activation function). The output of this neural network is the gating weight value. The final fusion context feature vector is generated as follows: the enhanced dynamic feature vector is weighted using the gating weight, and the original static structural feature vector is weighted using (1 minus the gating weight). Then, the two weighted vectors are added together to form a unique feature representation for each graph node that integrates static form and dynamic function, i.e., the fusion context feature vector.

[0065] In the serialization topology constraint segmentation module, the system determines which nodes in the graph structure correspond to the starting points of key anatomical structures (such as the trunks of vascular trees and bronchial trees) that require priority processing. This task is accomplished based on the fused context feature vector from each graph node. The fused context feature vector integrates static anatomical morphology and dynamic motion information. The system utilizes a pre-trained classifier, a support vector machine (SVM). The classifier's input is the fused context feature vector of each node, and its output is a binary classification result, determining whether the node belongs to the root node of a key anatomical structure. The classifier's training data consists of root node positions labeled by experts on sample data. The system traverses all nodes in the graph structure, collecting all nodes identified as root nodes by the classifier to form an initial root node set. This set provides a reliable starting point for subsequent sequence sorting.

[0066] After obtaining the initial set of root nodes, the system needs to calculate the proximity of all other nodes in the graph to these root nodes, using topological geodesic distance as the metric. Topological geodesic distance is defined here as the minimum number of edges required to reach any nearest root node from a non-root node along the established graph edge connection path on the initial graph structure; that is, the shortest path hop count. The system uses a breadth-first search algorithm for this calculation. During algorithm initialization, the topological geodesic distance of all root nodes is set to 0, and they are added to a queue for processing. The initial distance values ​​of all other non-root nodes are set to a very large number (such as infinity). Then, the algorithm removes a node from the queue and checks all its adjacent nodes (nodes directly connected by edges). For each adjacent node, if the distance from the current node to that adjacent node (i.e., the distance of the current node plus 1) is less than the currently recorded distance value of that adjacent node, the distance value of the adjacent node is updated to this smaller value, and the adjacent node is added to the queue. This process continues until the queue is empty, at which point each node has obtained the shortest path hop count to the nearest root node, i.e., its topological geodesic distance. This process ensures that even if there are multiple root nodes, each node can correctly calculate the distance to the nearest root node.

[0067] Based on the calculated topological geodesic distances, the system sorts all graph nodes to determine the order in which the sequence prediction model processes them. The primary criterion for sorting is the topological geodesic distance, arranged in ascending order. This means that nodes closer to the root node (i.e., those in the core or main regions of the anatomical structure) are processed first, while nodes farther away (i.e., those on the edges or at the ends) are processed later. This "inside-out" processing order aligns with the natural hierarchy of the anatomical structure, helping to first establish the correct labels for the main structures before gradually inferring the detailed regions. When multiple nodes have the same topological geodesic distance, the system employs a secondary sorting criterion to break the tie. This secondary criterion compares the magnitude (Euclidean length) of the fused context feature vectors of these nodes. The magnitude of the feature vector reflects, to some extent, the salience or information intensity of the feature at that node. The system sorts nodes with the same distance in descending order of magnitude, prioritizing nodes with more salient features. Finally, the system outputs a linear, ordered list of graph nodes, i.e., the processing sequence.

[0068] Before the sequence prediction model begins processing each node in the sequence, it needs to determine whether specific anatomical topology rules need to be activated to constrain the label prediction of the current node based on the states of the already processed neighboring nodes. The system predefines an anatomical topology rule base, in which the morphological and connectivity constraints that key anatomical structures must satisfy are encoded in the form of logical conditions. For example, a rule might state: "If a node has been labeled as a pulmonary artery, and one of its neighboring nodes (in the original graph structure) has also been labeled as a pulmonary artery, then nodes on the path between these two nodes, if their features are ambiguous, should also be labeled as pulmonary arteries to maintain vascular connectivity." When predicting the label of the current node, the system queries the label states of all its adjacent nodes in the original graph structure that have already been processed in the sequence (i.e., processed neighboring nodes). These label states are then matched against the conditions in the rule base. If the neighborhood state of the current node satisfies the triggering condition of one (or more) rules, the rule is determined to be triggered. The output of rule triggering is usually an indication signal, which may contain the type of rule triggered and the label tendency suggested by the rule.

[0069] When a rule is triggered, this prior anatomical knowledge needs to be incorporated into the decision-making process of the sequence prediction model. Sequence prediction models (typically based on gated recurrent units or long short-term memory network architectures) utilize attention mechanisms to weigh the influence of the hidden states of already processed nodes (especially their neighboring nodes) on the current node's decision. Normally, attention weights are automatically learned by the model, based on factors such as node feature similarity. When a rule is triggered, the system dynamically overrides or modifies these automatically learned attention weights according to the type of rule. For example, if the triggered rule requires maintaining vascular connectivity, the system will forcibly increase the attention weights for the current node's already marked vascular neighboring nodes. Specifically, the rule trigger indication signal can be converted into a weight adjustment vector, which is combined with the model's original attention score (e.g., through addition or multiplication) to generate a new set of adjusted attention weights that conform to anatomical constraints. This allows the model to refer more to information from neighboring nodes that are crucial for maintaining topological correctness when making decisions.

[0070] The core of the sequence prediction model is a gated recurrent unit (GRU). The input to this unit mainly includes the following: first, the fused context feature vector of the current node; second, the weighted sum of the hidden states of all processed nodes after adjustment in the previous step (i.e., the regulated neighborhood node state context), which uses adjusted attention weights; and third, the rule-triggered indicator signal can optionally be included as part of the input. Internally, the GRU uses its update and reset gates to determine how much historical information (from the processed sequence) to retain and how much new input information to accept. It synthesizes all this input information, updates its internal hidden states, and finally calculates the probability distribution of the current node belonging to various possible labels (such as background, blood vessel, tumor, etc.) through an output layer (usually a fully connected layer with a Softmax activation function). The label with the highest probability is determined as the final predicted label for that node. Because this decision-making process integrates the features of the current node, the rule-guided neighborhood context, and the sequence history information, it can make robust judgments that conform to both local image features and global anatomical topological constraints. The model then passes the predicted label and updated hidden state of the current node down to process the next node in the sequence.

[0071] In the segmentation result verification and model reconstruction module, the topological correctness of key anatomical structures is quantitatively evaluated from the initial segmentation mask (a three-dimensional binary or label array) output by the serialized topological constraint segmentation module. The system first extracts the topological skeleton network from the corresponding binary segmentation region for anatomical structures requiring guaranteed connectivity (such as the pulmonary artery, pulmonary vein, and bronchi). The extraction process employs a three-dimensional thinning algorithm (also known as skeletonization), which iteratively removes voxels from the object's surface while maintaining its topological structure (such as connectivity and the number of holes) until only a single voxel-wide line structure remains; this structure is the object's central skeleton. Subsequently, the system performs connected component analysis on this skeleton network. Connected component analysis traverses all voxels (nodes) in the skeleton network, marking interconnected nodes as the same connected region. The system calculates the total number of connected components. This number is then compared to the number of normal connected components specified for that specific anatomical structure in a predefined anatomical topology rule base. For example, a normal left pulmonary artery trunk should be a single connected component on the skeleton network. If the number of connected components calculated is greater than 1, it indicates that the structure may have abnormal breaks (i.e., non-physiological breaks) in the segmentation result; if the number is abnormally low, it may indicate that there are false connections in parts that should not be connected. The system will record the location of these topologically abnormal regions.

[0072] For identified areas with abnormal breaks, the system needs to accurately locate the break endpoints and analyze the geometric characteristics of the fracture gap to provide a basis for subsequent repair. The system identifies the break location on the skeleton network and locates the skeleton endpoints at both ends of the break. Since the skeleton is only a single voxel wide, these endpoints can be precisely located. Next, the system returns to the original 3D segmentation mask and analyzes the solid region around the break endpoints. Surface contour points within the local 3D neighborhood of the break endpoint are extracted on the 3D repair mask. The system calculates the geometric features of the 3D spatial gap at the break, primarily including the Euclidean distance (straight-line distance) between the two ends of the break, and the average gap width estimated by calculating the average of the shortest distances from all voxels along the path from one endpoint to the other to the opposite surface. Furthermore, the system analyzes the orientation of local surfaces at the break, for example, by calculating the normal vector direction of the surface point cloud at the openings at both ends of the break through principal component analysis to determine the relative orientation of the fracture surface. These geometric features are crucial for selecting appropriate repair algorithms and parameters.

[0073] After obtaining the precise location of the fracture endpoints and the geometric features of the gap, the system performs a repair operation to restore the correct topological connectivity of the anatomical structure. The goal of the repair is to generate a smooth transition surface within the fracture gap, naturally connecting the two ends of the fracture. This invention employs a constraint-based surface interpolation algorithm, such as the radial basis function interpolation method represented by thin plate splines. This algorithm uses the contour points identified at both ends of the fracture (including the endpoints and their neighboring contour points) as constraints for interpolation. Specifically, the three-dimensional coordinates of these contour points are used as known data points, requiring the generated surface to pass through these points. The thin plate spline interpolation algorithm finds a smooth function (i.e., the interpolation surface) by solving a system of linear equations. This function satisfies the condition of passing through all constraint points while minimizing its bending energy (usually defined by the square integral of the second derivative), thus ensuring that the generated surface is the smoothest. By solving this system of equations, a function defined in the entire three-dimensional space can be obtained. Then, the system extracts the zero isosurface of this function, which is the smooth transition surface connecting the two ends of the fracture. Finally, this repaired surface is voxelized and merged into the original 3D segmentation mask, thereby geometrically filling the gaps and topologically restoring connectivity.

[0074] After topology restoration, the macroscopic geometric parameters of the entire 3D segmentation mask need to be verified to ensure not only correct connectivity but also that its overall size and shape are within an anatomically reasonable range. The system calculates the total volume and total surface area of ​​the restored 3D segmentation mask. Volume is calculated by counting all voxels belonging to the target tissue and then multiplying by the physical volume of each voxel (determined by the image's slice thickness and pixel pitch). Surface area is calculated by summing the areas of the triangular patches generated after surface reconstruction of the segmentation mask. The system then compares these calculated volume and surface area values ​​with the ranges of typical organ geometric parameters recorded in a predefined anatomical topology rule base. This rule base contains normal ranges of organ volume and surface area derived from a large amount of anatomical statistics from normal populations. If the restored model's geometric parameters fall within this typical range, its geometric rationality is considered validated; if they significantly exceed the range (e.g., abnormally small or large volume), it may indicate a large error in the segmentation or restoration process, requiring re-examination or parameter adjustment.

[0075] The initial surface mesh obtained through direct conversion may contain unnatural smooth areas introduced by interpolation restoration, and may also lose some true anatomical details from the original image data. To address this issue, the system performs multi-resolution analysis on the initial surface mesh. First, the mesh is simplified to generate a low-resolution (coarse-resolution) version. At this coarse-resolution level, the Laplacian smoothing algorithm is applied to smooth out unnatural small-scale surface undulations that may be caused by restoration. Then, the smoothed low-resolution mesh is fused with the original high-resolution mesh for detail fusion. Specifically, the geometric details (such as displacement vectors) of the original high-resolution mesh relative to the low-resolution mesh are calculated, and these details are selectively added back to the smoothed low-resolution mesh. This selection process is based on the scale of the details, retaining medium to fine-scale details that represent true anatomical features, while filtering out extremely high-frequency details that may be noise.

[0076] The final step is to ensure that the optimized curved mesh meets the quality requirements for clinical use. First, watertightness is ensured, guaranteeing the mesh is a completely closed manifold structure without boundary edges. This is achieved by checking and filling any holes in the mesh. Then, error-free facet checks are performed to ensure that the normals of all triangular facets are aligned (usually outwards) and that there are no self-intersecting or severely distorted facets. Problematic facets are automatically repaired or require manual intervention. Finally, the system outputs a 3D model with correct topological connectivity, reasonable geometry, and a high-quality mesh structure—the preoperative virtual reconstructed 3D model. This model is stored in standard 3D file formats (such as STL, PLY, or OBJ), which can be directly imported into 3D visualization software for rotation and scaling, printed as a physical model, or integrated into a surgical navigation system for intraoperative registration and guidance, providing accurate anatomical references for preoperative planning in thoracic surgery.

[0077] The working principle of this invention is as follows: A multimodal data fusion and dynamic graph construction module simultaneously acquires chest CT image data and respiratory CT image data, constructs an initial graph structure with voxels as nodes and spatial adjacency relationships as edges, and extracts the displacement trajectory features of each voxel during the respiratory cycle from the dynamic data as the dynamic attributes of the nodes; A cross-modal context-aware feature encoding module extracts the static structural features and dynamic motion features of the displacement trajectory sequence of each graph node, and generates a context feature vector that integrates static morphology and dynamic function through feature interaction and gating fusion based on a temporary graph attention mechanism; A serialized topological constraint segmentation module sorts the nodes from the core of the anatomical structure to the edge according to the feature vector to form a processing sequence, and uses a sequence prediction model combined with a predefined anatomical topology rule base for label prediction to ensure that the segmentation results conform to the morphological and connectivity constraints of key structures; A segmentation result verification and model reconstruction module performs topological skeleton connectivity verification, fracture repair, geometric rationality verification and multi-resolution mesh optimization on the initial segmentation mask, and outputs a high-quality, topologically correct 3D model, which can be directly used for preoperative planning, navigation or printing.

[0078] The foregoing has provided a detailed description of one embodiment of the present invention, but this description is merely a preferred embodiment and should not be construed as limiting the scope of the invention. All equivalent variations and modifications made within the scope of the claims of this invention should still fall within the patent coverage of this invention.

Claims

1. A virtual reconstruction system for preoperative thoracic surgery, characterized in that, include: The multimodal data fusion and dynamic graph construction module is used to simultaneously acquire chest CT image data and respiratory CT image data; an initial graph structure is constructed based on the chest CT image data, voxel units are defined as graph nodes, and spatially adjacent nodes are connected to form edges; the displacement trajectory of each voxel unit in respiratory motion is extracted from the respiratory CT image data and added to the corresponding graph node as a dynamic attribute; The cross-modal context-aware feature encoding module extracts the local 3D image block and corresponding displacement trajectory sequence of each image node in the CT image data; processes the image block and trajectory sequence separately to obtain static structural features and dynamic motion features; integrates the static structural features and dynamic motion features to generate a fused context feature vector for each image node. The serialization topology constraint segmentation module sorts graph nodes according to the fused context feature vector to form a processing sequence. Based on the processing sequence, a sequence prediction model is used to process each graph node in turn. The sequence prediction model performs label prediction by combining the current graph node features, the state of the processed neighboring graph nodes, and a predefined anatomical topology rule base. The rule base specifies the morphological and connectivity constraints of key anatomical structures. The segmentation result verification and model reconstruction module maps the labels obtained from sequence prediction to three-dimensional space to generate an initial segmentation mask; it performs global anatomical rationality verification and correction on the initial segmentation mask, and outputs a preoperative virtual reconstructed three-dimensional model for preoperative clinical use in thoracic surgery.

2. The virtual reconstruction system for preoperative thoracic surgery according to claim 1, characterized in that, The initial image structure constructed based on chest CT image data specifically includes: The spatial coordinates of each voxel unit in the CT image data are mapped to graph nodes. Based on the preset 26 adjacency rules, all voxel units are traversed. If two voxel units are connected by a face, edge or point in three-dimensional space, an undirected edge is established between the graph nodes corresponding to the two voxel units. Calculate the gradient magnitude of image grayscale values ​​between adjacent voxel units; assign initial weights to each edge based on the gradient magnitude of image grayscale values; Centered on each graph node, image patch features at different scales are extracted; the image patch features are fused with the initial weights to generate enhanced edge attributes, thus completing the optimized construction of the initial graph structure.

3. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The step of extracting the displacement trajectory of each voxel unit during respiratory motion from respiratory CT image data and adding it as a dynamic attribute to the corresponding graph node specifically includes: The respiratory CT image data is divided into respiratory phases, and the end of inspiration or end of expiration phase is selected as the reference frame; the dense three-dimensional motion vector field of each phase image relative to the reference frame is calculated. For each voxel unit, based on its spatial position in each phase image, a spatiotemporal displacement trajectory spanning a complete respiratory cycle is formed by interpolating its displacement vector in the continuous phase motion vector field. Three feature quantities are extracted from the displacement trajectory: trajectory length, average displacement amplitude, and main direction of motion. These three features together form a trajectory feature vector. The trajectory feature vector is then used as a dynamic attribute and attached to the corresponding graph node in the initial graph structure.

4. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The processing of image blocks and trajectory sequences yields static structural features and dynamic motion features, specifically including: For local 3D image blocks, a set of directionally separable 3D convolutional kernels are used for feature extraction. The 3D convolutional kernels have different weights along the axial, sagittal, and coronal directions to anisotropically perceive structural information of different anatomical planes and output multi-directional feature maps. For the displacement trajectory sequence, its displacement vector is decomposed into three motion primitive components: translation, rotation, and scaling; the amplitude and phase change spectrum of each motion primitive component in each respiratory cycle is calculated to form a dynamic motion feature vector; Principal component analysis was performed to reduce the dimensionality of the multi-directional feature maps and dynamic motion feature vectors, and the reduced features were mapped to a unified numerical range to obtain standardized static structural features and dynamic motion features.

5. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The process of integrating static structural features and dynamic motion features to generate a fused context feature vector for each graph node specifically includes: Construct a fully connected temporary graph using all graph nodes in the current batch as vertices; use the static structural features and dynamic motion features of each node as the two attributes of the corresponding vertex; Calculate the attention weight between any two vertices based on their static structural features similarity; use the attention weight to weighted aggregate the dynamic motion features of the target vertex, and generate a dynamic feature enhanced with structural context for each vertex; The original static structural features and the dynamic features enhanced by structural context are input into a gated fusion unit. The gated fusion unit learns a dynamic weight to control the fusion ratio of the two types of features and finally outputs a unique fusion context feature vector for each graph node.

6. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The process of sorting graph nodes based on fused context feature vectors to form a processing sequence specifically includes: Based on the fusion context feature vector, graph nodes belonging to the starting point of the main trunk of the key anatomical structure are identified and used as the initial root node set; Starting from the root node set, calculate the topological geodesic distance from all other graph nodes to the nearest root node on the initial graph structure; the topological geodesic distance refers to the shortest number of hops along the graph edge connection path; All graph nodes are sorted in ascending order of topological geodesic distance. Graph nodes with the same distance are then sorted in a secondary order according to the magnitude of their fused context feature vectors. Finally, a processing sequence radiating from the core of the structure to the edge is generated.

7. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The sequential prediction model is used to process each graph node in turn, specifically including: Before predicting the label of the current graph node, query the label status of the processed neighboring graph nodes to determine whether they have triggered the rules in the predefined dissection topology rule base. If a rule is triggered, the attention weights in the sequence prediction model are adjusted according to the rule type. The sequence prediction model integrates the features of the current node, the state context of the regulated neighboring nodes, and finally decides the label of the current node through a gated recurrent unit.

8. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The process of performing a global dissection rationality check and correction on the initial segmentation mask specifically includes: Extract the topological skeleton network of key anatomical structures from the initial segmentation mask; calculate the number of connected components in the topological skeleton network and compare it with the number of normal connected components specified in the predefined anatomical rule base to identify whether there are abnormal disconnections; For the identified abnormal breakage areas, determine the fracture endpoints and calculate the geometric characteristics of the three-dimensional spatial gap at the fracture point; Using the fracture endpoints and their neighboring contour points as constraints, a smooth transition repair surface is generated within the fracture gap to restore the correct topological connectivity of the anatomical structure.

9. The preoperative virtual reconstruction system for thoracic surgery according to claim 1, characterized in that, The output of the preoperative virtual reconstructed 3D model specifically includes: The total volume and surface area of ​​the repaired 3D segmentation mask are calculated and compared with the range of typical organ geometric parameters recorded in the predefined anatomical topology rule base to verify its geometric rationality. Multi-resolution analysis was performed on the initial surface mesh model obtained by the mask conversion after restoration. At the coarse resolution level, the unnatural surface undulations introduced by interpolation restoration were smoothed, while at the fine resolution level, the real anatomical details from the original image data were preserved. The optimized curved mesh is subjected to watertightness processing and error-free patch inspection, and finally a preoperative virtual reconstructed 3D model is output. The preoperative virtual reconstructed 3D model can be directly used for 3D visualization, printing or surgical navigation.

Images