A wireless ultrasound multi-modal image registration and fusion system

Abstract

This invention belongs to the field of medical image processing technology. It discloses a wireless ultrasound multimodal image registration and fusion system. The system uses a comprehensive acquisition module to acquire Doppler blood flow data, inertial navigation attitude information, and CT vascular models in real time. An attitude compensation module eliminates zero-bias drift, ensuring probe attitude measurement stability. A Doppler angle calculation module objectively calculates geometric Doppler angles based on CT data and the attitude matrix, eliminating subjective estimation errors. A vector estimation module constructs a three-dimensional blood flow velocity field combined with vascular anatomical constraints. A spatial registration module achieves accurate fusion of ultrasound and CT data. A consistency verification module automatically identifies attitude-sensitive regions. An iterative correction module performs adaptive optimization on unstable regions. This invention improves the stability and accuracy of wireless ultrasound Doppler measurements, overcomes the dimensional limitations of radial velocity measurement, and achieves seamless registration of multimodal data.

Description

Technical Field

[0001] This invention relates to the field of medical image processing technology, and more specifically, to a wireless ultrasound multimodal image registration and fusion system. Background Technology

[0002] Current wireless Doppler ultrasound blood flow measurement technology faces a physical and engineering contradiction, affecting its application value in clinical vascular assessment. The Doppler effect inherently exhibits an angle-dependent (cosine effect), meaning the accuracy of blood flow velocity measurement depends entirely on the angle between the sound beam and the blood flow direction. In clinical practice, this angle dependence becomes an obstacle when physicians hand-hold wireless probes to scan complex blood vessel routes. Especially at critical anatomical locations such as carotid bifurcations, cerebral vascular bends, or abdominal aortic branches, even slight changes in probe posture can cause significant fluctuations in measurement results, making it impossible for physicians to obtain stable and reliable diagnostic parameters. The problem becomes even more complex when attempting to spatially register these unstable Doppler data with CT three-dimensional vascular models. CT imaging, as an omnidirectional isotropic modality, provides accurate vascular anatomy information, but the inertial navigation unit built into the wireless ultrasound probe faces a cumulative zero-bias drift problem, gradually distorting the mapping relationship between the probe coordinate system and the patient's anatomical coordinate system. Traditional Doppler ultrasound can only measure the radial velocity component along the sound beam direction, and cannot directly obtain lateral velocity information perpendicular to the sound beam, creating an "information blind spot" in three-dimensional spatial registration. This "uncertainty principle" manifests clinically as a drastic fluctuation in quantitative blood flow parameters (such as peak systolic velocity, PSV) with changes in probe holding angle, lacking physical consistency. When physicians rely on these unstable data to assess the severity of carotid artery stenosis or peripheral vascular disease, it may lead to misdiagnosis or unnecessary interventional treatment decisions.

[0003] In view of this, the present invention proposes a wireless ultrasound multimodal image registration and fusion system to solve the above problems. Summary of the Invention

[0004] To overcome the aforementioned deficiencies of the prior art and to achieve the above objectives, the present invention provides the following technical solution: a wireless ultrasound multimodal image registration and fusion system, comprising:

[0005] The integrated acquisition module is used to acquire Doppler blood flow data from the wireless ultrasound probe, attitude quaternion sequences output by the inertial navigation unit, and CT three-dimensional vascular model data. The Doppler blood flow data includes radial velocity measurements along the sound beam axis and corresponding nominal Doppler angle values.

[0006] The attitude compensation module is used to perform zero-bias drift compensation on the attitude quaternion sequence and extract the real-time rotation matrix from the probe coordinate system to the patient anatomical coordinate system after compensation, which is denoted as the probe attitude matrix.

[0007] The Doppler angle calculation module is used to extract the tangential direction vector of the blood vessel centerline based on CT three-dimensional blood vessel model data, and project the sound beam direction vector onto the patient's anatomical coordinate system in combination with the probe posture matrix, and calculate the measured angle between the sound beam direction and the blood vessel tangential direction, which is denoted as the geometric Doppler angle.

[0008] The vector estimation module is used to perform cosine correction on radial velocity measurements based on geometric Doppler angles, obtain the axial velocity component along the tangent direction of the blood vessel, and construct a probabilistic constraint model of the transverse velocity component based on the normal vector of the blood vessel section to generate a three-dimensional vector estimate of blood flow velocity.

[0009] The spatial registration module is used to spatially register the estimated three-dimensional vector value of blood flow velocity with the three-dimensional CT vascular model data in the patient's anatomical coordinate system to construct a Doppler-CT fused blood flow field.

[0010] The consistency verification module is used to extract the peak flow velocity distribution at each blood vessel cross section in the Doppler-CT fusion blood flow field, perform consistency verification with the nominal value of the Doppler angle, and identify abnormal areas where the peak flow velocity changes with the probe posture.

[0011] The iterative correction module is used to introduce attitude sensitivity weights into the probe attitude matrix based on the spatial location of the jump abnormal region, and to iteratively correct the geometric Doppler angle to generate stable blood flow quantitative results.

[0012] The technical effects and advantages of the wireless ultrasound multimodal image registration and fusion system of the present invention are as follows:

[0013] This invention effectively solves the problem of velocity measurement jumps caused by probe posture changes in traditional methods by eliminating angle-dependent factors, enabling clinicians to obtain stable and consistent hemodynamic parameters, especially in anatomically complex areas such as carotid bifurcation, intracranial vessels, and the abdominal aorta. This invention overcomes the inherent limitations of radial velocity measurement, achieving the reconstruction and characterization of a complete three-dimensional blood flow vector field, providing a solid foundation for advanced hemodynamic analyses such as blood flow shear force, eddy current distribution, and collateral circulation. By establishing a precise spatial correspondence between ultrasound data and CT models, this invention achieves seamless integration of functional blood flow information and anatomical structural information, enabling clinicians to comprehensively assess vascular lesions within a unified reference framework. This stable blood flow quantification capability significantly improves the accuracy and consistency of stenosis assessment, providing a reliable basis for precise intervention decisions and reducing the misjudgment rate of borderline cases. Furthermore, this invention provides a risk assessment framework for clinical decision-making by intelligently identifying and marking areas of reduced measurement reliability, enhancing medical safety and optimizing resource utilization efficiency. Attached Figure Description

[0014] Figure 1 This is a schematic diagram of a wireless ultrasound multimodal image registration and fusion system according to the present invention. Detailed Implementation

[0015] 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.

[0016] This invention provides a wireless ultrasound multimodal image registration and fusion system. The system's execution entities include, but are not limited to, ultrasound diagnostic equipment, medical imaging workstations, vascular interventional surgery navigation platforms, and intelligent diagnostic auxiliary systems, which can be considered as general computing nodes in this application.

[0017] Please see Figure 1 In this embodiment of the invention, a wireless ultrasound multimodal image registration and fusion system includes:

[0018] The integrated acquisition module acquires Doppler blood flow data from the wireless ultrasound probe, attitude quaternion sequences output by the inertial navigation unit, and CT 3D vascular model data. The Doppler blood flow data includes radial velocity measurements along the sound beam axis and corresponding nominal Doppler angle values, which are the fundamental inputs for blood flow velocity estimation. The inertial navigation unit integrates a three-axis gyroscope, accelerometer, and magnetometer, and can output the probe's attitude quaternion sequences in real time, describing the probe's orientation in three-dimensional space. The CT 3D vascular model data provides precise structural information about the patient's vascular anatomy, including the vessel centerline, lumen contour, and branch topology. These three types of data are acquired in real time through a multi-channel synchronous acquisition interface, providing comprehensive raw data for subsequent registration and fusion, ensuring the completeness and accuracy of the system analysis.

[0019] The attitude compensation module compensates for zero-bias drift in the attitude quaternion sequence, extracting the real-time rotation matrix from the compensated probe coordinate system to the patient anatomical coordinate system, denoted as the probe attitude matrix. Attitude compensation is a crucial step in ensuring spatial registration accuracy, improving the long-term stability of attitude measurements by eliminating zero-bias drift of the inertial sensor. The compensation process first performs segmented analysis on the attitude quaternion sequence, identifying zero-bias features during the stationary phase; then, a drift prediction model is established to achieve dynamic zero-bias compensation; finally, the compensated quaternion is converted into a rotation matrix to accurately describe the spatial relationship between the probe and the patient anatomical coordinate system. This real-time compensation mechanism effectively addresses electromagnetic interference and temperature variations in the clinical environment, ensuring the consistency and reliability of attitude measurements during long-term scanning.

[0020] The Doppler angle calculation module extracts the tangential direction vector of the vessel centerline from CT 3D vascular model data. Combined with the probe posture matrix, it projects the sound beam direction vector onto the patient's anatomical coordinate system, calculating the measured angle between the sound beam direction and the vessel tangential direction, denoted as the geometric Doppler angle. Doppler angle calculation is a core component of blood flow velocity estimation. By calculating the angle between the sound beam and the blood flow direction through geometric spatial relationships, it overcomes the limitations of angle estimation in traditional Doppler ultrasound. The calculation process first accurately extracts the vessel centerline trajectory from the CT model; then, combined with the probe posture matrix, it projects the ultrasound beam direction onto the same coordinate system; finally, it calculates the spatial angle between the two direction vectors to obtain the true geometric Doppler angle. This angle calculation method based on spatial geometry avoids the reliance on operator subjective judgment in traditional Doppler ultrasound, improving the objectivity and accuracy of blood flow velocity estimation.

[0021] The vector estimation module performs cosine correction on radial velocity measurements based on geometric Doppler angles to obtain the axial velocity component along the vessel tangent. It then constructs a probabilistic constraint model for the lateral velocity component based on the vessel cross-section normal vector, generating a three-dimensional vector estimate of blood flow velocity. Vector estimation is a key technology for extending one-dimensional Doppler measurements to three-dimensional blood flow characterization, achieving complete velocity reconstruction through physical models and probabilistic constraints. The estimation process first uses geometric Doppler angles to perform cosine correction on the radial velocity to recover the axial component; then, based on vascular anatomy, a probabilistic model for the lateral velocity is established; finally, by integrating axial and lateral information, a three-dimensional velocity vector and uncertainty assessment are generated. This vector estimation method, which integrates measurement data and anatomical constraints, overcomes the limitation of traditional Doppler technology, which can only measure radial velocity, providing rich information for a comprehensive understanding of hemodynamic characteristics.

[0022] The spatial registration module is used to spatially register the estimated three-dimensional vector values ​​of blood flow velocity with CT three-dimensional vascular model data in the patient's anatomical coordinate system, constructing a Doppler-CT fused blood flow field. Spatial registration is a core step in achieving multimodal data fusion. Through coordinate transformation and interpolation filling, point-like velocity measurements are mapped to a continuous vascular model. The registration process first projects the velocity vector to the nearest vascular centerline node; then interpolates along the vascular axis to fill unsampled areas; finally, the principle of flow conservation is applied for continuity constraints to ensure the physical rationality of the blood flow field. This anatomical-based spatial registration method achieves seamless fusion of discrete ultrasound measurements and continuous CT models, providing an intuitive and comprehensive visualization of the blood flow field for clinical assessment.

[0023] The consistency verification module extracts the peak velocity distribution at each vessel cross-section in the Doppler-CT fusion blood flow field and performs consistency verification with the nominal Doppler angle value, identifying abnormal regions where the peak velocity jumps with probe posture changes. Consistency verification is an important mechanism for evaluating the reliability of fusion results. By comparing and analyzing, it identifies posture-sensitive areas, providing a basis for subsequent corrections. The verification process first extracts the peak velocity distribution at each cross-section from the fusion blood flow field; then analyzes the time correlation between flow velocity and probe angle; finally, it identifies regions with abnormal velocity jumps and marks them as unstable regions requiring special attention. This consistency verification method based on time-series analysis can automatically discover potential problems in the registration and fusion process, ensuring the clinical reliability of the final results and avoiding measurement bias caused by posture changes.

[0024] The iterative correction module introduces attitude sensitivity weights into the probe attitude matrix based on the spatial location of the aberration region, iteratively correcting the geometric Doppler angle to generate stable quantitative blood flow results. Iterative correction is an innovative mechanism to improve the system's anti-interference capability by reducing measurement fluctuations in attitude-sensitive regions through an adaptive weighting strategy. The correction process first analyzes the angular instability characteristics of the aberration region; then, it establishes an attitude sensitivity weight model to weaken the influence of unstable regions; finally, iterative optimization algorithms are used to gradually correct the geometric Doppler angle until the results converge and stabilize. This data-driven iterative correction method effectively addresses the challenges of complex vascular geometries and probe operation variations, significantly improving the stability and repeatability of quantitative blood flow, and providing more reliable hemodynamic assessments for clinical decision-making.

[0025] The modules are connected via wired and / or wireless means to enable data transmission between them.

[0026] In this embodiment of the invention, the detailed implementation steps of performing zero-bias drift compensation on the attitude quaternion sequence and extracting the real-time rotation matrix from the compensated probe coordinate system to the patient anatomical coordinate system, denoted as the probe attitude matrix, include:

[0027] The attitude quaternion sequence is segmented according to time windows. The mean value of the gyroscope output within each stationary period is calculated and recorded as the current zero-bias estimate. Time window segmentation is the fundamental step in identifying stationary periods, dividing motion and stationary phases by analyzing the rate of attitude change. The segmentation process employs a sliding window technique, calculating the standard deviation of the angular velocity within the window. When the standard deviation is below a preset threshold (typically 0.5-1° / s), it is determined to be a stationary period. The identification of stationary periods leverages the characteristic that probes often remain relatively stable during clinical ultrasound examinations, providing an ideal observation window for zero-bias estimation. Within each identified stationary period, the average value of the three-axis gyroscope output is calculated as the zero-bias estimate for the current time period. This zero-bias estimation method based on stationary periods requires no additional calibration steps, can adaptively acquire zero-bias characteristics during normal examination procedures, and adapts to different equipment and environmental conditions.

[0028] The current zero-bias estimate is substituted into a first-order Markov drift model to predict the zero-bias value at the next time step. This prediction is then weighted and fused with the measured zero-bias estimate to obtain a smoothed zero-bias sequence. The first-order Markov drift model is an effective tool for describing the time-varying nature of sensor zero-bias, improving prediction accuracy by considering time correlation. In the model specification, the time correlation coefficient is a key parameter, determining the smoothness of the zero-bias change. The correlation coefficient is determined by analyzing the autocorrelation function of historical zero-bias sequences, with typical values ​​between 0.9 and 0.99, reflecting the slow nature of zero-bias drift. The prediction equation is expressed as follows:

[0029] ;

[0030] in, for The zero bias of the prediction at time t. for The zero-bias estimate at time t. The initial zero bias value, The time correlation coefficient is used. The prediction variance is calculated using the error propagation rule and serves as the basis for weighted fusion. When a new measured zero-bias estimate is available, the predicted and measured values ​​are fused using the Kalman gain formula to obtain an updated zero-bias estimate. The Kalman gain calculation considers the ratio of the prediction variance to the observation noise variance, automatically adjusting the relative weights of the predicted and measured values ​​to achieve the optimal estimate. This zero-bias estimation method based on a state-space model can smoothly filter out noise interference while rapidly responding to actual zero-bias changes, providing a stable and reliable benchmark for attitude measurement.

[0031] The smoothed zero-bias sequence is subtracted frame-by-frame from the attitude quaternion sequence to obtain the zero-bias-compensated quaternion sequence. Zero-bias compensation is a key operation for improving attitude measurement accuracy, enhancing the long-term stability of attitude angles by systematically eliminating the influence of zero bias. The compensation process first converts the smoothed zero-bias sequence into an equivalent angular increment; then, forward integration correction is performed on the original attitude quaternion sequence to eliminate the influence of accumulated errors; finally, the zero-bias-compensated quaternion sequence is obtained, accurately representing the probe's true attitude change. The compensation algorithm operates directly in quaternion space, avoiding the singularity problem in Euler angle representation and ensuring computational stability under arbitrary attitudes. This frame-by-frame compensation strategy effectively handles the cumulative effect of zero bias during long-term inspections, significantly improving the long-term accuracy of attitude measurement and laying the foundation for subsequent spatial registration.

[0032] The quaternion sequence after zero-bias compensation is converted into a 3×3 rotation matrix, which is then output as the probe attitude matrix. Rotation matrix conversion is a preparatory step for spatial registration, transforming the quaternion representation into a matrix form directly applicable to vector transformations. The conversion process follows the standard mapping relationship from quaternions to rotation matrices, ensuring the mathematical equivalence of the two representations. Quaternions Convert to rotation matrix The calculation formula is:

[0033] ;

[0034] Numerical optimization techniques are employed in the transformation implementation to ensure the orthogonality of the rotation matrix and avoid distortion caused by accumulated errors. The output probe attitude matrix directly describes the spatial transformation relationship from the probe coordinate system to the patient anatomical coordinate system, providing an accurate attitude reference for subsequent beam direction projection and spatial registration. This quaternion-based attitude representation method, combined with a zero-bias compensation strategy, effectively solves the challenges of accuracy and stability in attitude measurement in wireless ultrasound systems, creating conditions for multimodal fusion.

[0035] In this embodiment of the invention, the detailed implementation steps include: extracting the tangential direction vector of the blood vessel centerline based on CT three-dimensional blood vessel model data, projecting the sound beam direction vector onto the patient's anatomical coordinate system in conjunction with the probe posture matrix, and calculating the measured angle between the sound beam direction and the blood vessel tangential direction, denoted as the geometric Doppler angle.

[0036] A centerline extraction algorithm is applied to CT 3D vascular model data. A node is sampled along the centerline at preset arc lengths, and cubic spline fitting is performed on the connections between adjacent nodes to obtain the tangent direction vector at each sampling node. Centerline extraction is a fundamental step in vascular geometry analysis, obtaining the central axis representation of the vascular lumen through a skeletonization algorithm. The extraction process first performs a distance transformation on the CT vascular model, calculating the distance from each voxel to the vascular boundary; then, based on the maximum distance criterion, iteratively extracts the center point sequence; finally, through topological connectivity analysis, a complete centerline network is constructed, including the trunk and branch structures. The extracted raw centerlines often contain local noise and irregularities, requiring further optimization. The sampling and fitting steps involve uniform sampling along the centerline at equal arc length intervals (typically 1-3 mm), followed by cubic spline interpolation applied to the sampling point sequence to generate a smooth and continuous centerline curve. At each sampling node, a normalized tangent direction vector is calculated using the first derivative of the spline curve to accurately represent the local vascular orientation. This centerline representation method based on spline fitting preserves the overall morphological characteristics of blood vessels while eliminating the influence of local noise, providing a reliable reference for blood vessel orientation in Doppler angle calculation.

[0037] The nominal direction vector of the emitted sound beam from the wireless ultrasound probe in the probe coordinate system is multiplied by the probe attitude matrix to obtain the direction vector of the sound beam in the patient's anatomical coordinate system. Sound beam direction projection is a key step in achieving coordinate system unification, mapping the sound beam direction in the probe's local coordinates to the patient's global coordinate system through matrix transformation. The projection process first determines the nominal direction vector of the probe's sound beam in the probe coordinate system, typically the normal direction of the ultrasound transducer array, represented as a unit vector [0,0,1]. Then, using the probe attitude matrix, coordinate transformation is achieved through vector-matrix multiplication; finally, the actual direction vector of the sound beam in the patient's anatomical coordinate system is obtained. To improve accuracy, in practical applications, the internal structural parameters of the ultrasound probe are considered, including transducer arrangement and focusing characteristics, and the nominal direction is finely adjusted to ensure consistency with the actual sound beam propagation direction. This real-time attitude-based sound beam projection method can accurately track the continuous changes in the sound beam direction during ultrasound scanning, laying the foundation for accurate calculation of the geometric Doppler angle.

[0038] The dot product of the sound beam direction vector and the tangent direction vector at the nearest node on the vessel centerline to the probe projection point is calculated, and the inverse cosine is used to obtain the geometric Doppler angle. Doppler angle calculation is a core step in velocity estimation, determining the angle correction coefficient through the spatial relationship between the sound beam and blood flow direction. The calculation process first determines the intersection point of the sound beam and the vessel by projecting the sound beam along its direction vector to calculate the point of minimum distance to the vessel centerline; then, the nearest centerline node is identified, and the tangent direction vector at that point is extracted; finally, the dot product of the two normalized direction vectors is calculated, and the spatial angle is obtained using the inverse cosine function. The calculation formula is:

[0039] ;

[0040] in, The geometric Doppler angle, The normalized beam direction vector. This is the normalized tangent vector of the blood vessel. This represents the vector dot product operation. When the calculated geometric Doppler angle is close to 0° or 180°, due to the flatness of the cosine function near its extreme points, small directional errors can lead to large angle changes. To improve stability, the system employs a weighted averaging strategy in this region, fusing angle estimates from adjacent nodes to reduce noise impact. This Doppler angle calculation method based on geometric spatial relationships overcomes the limitations of traditional ultrasound angle estimation, which relies on operator experience, and provides an objective basis for the accurate quantification of blood flow velocity.

[0041] When the geometric Doppler angle exceeds a preset angle range, an angle validity flag is triggered, invalidating the Doppler data of the corresponding frame. Angle validity verification is a crucial mechanism for ensuring data quality; by eliminating unreliable angle measurements, it improves the overall reliability of the results. The verification process sets a reasonable angle range threshold, typically 20°-160°, avoiding the extreme regions where the Doppler effect is insensitive. When the calculated geometric Doppler angle is within the range of 0°-20° or 160°-180°, the Doppler frequency shift is minimal or near its maximum due to the angles approaching parallel or antiparallel angles, significantly reducing the signal-to-noise ratio and easily introducing large velocity errors. The system automatically identifies these unfavorable angle regions, adds an invalid flag to the Doppler data of the corresponding frame, and handles or removes them in subsequent processing. To avoid data fragmentation, the system provides a smooth transition mechanism, assigning decreasing weights to data approaching the threshold boundary to achieve a gradual transition from valid to invalid. This automatic angle validity verification mechanism greatly improves the robustness and reliability of the system, avoids velocity estimation bias caused by poor acquisition angles, and provides quality assurance for clinical applications.

[0042] In this embodiment of the invention, the detailed implementation steps for generating a three-dimensional vector estimate of blood flow velocity include: cosine correction of radial velocity measurements based on geometric Doppler angles to obtain the axial velocity component along the tangential direction of the blood vessel; and constructing a probabilistic constraint model for the transverse velocity component based on the normal vector of the blood vessel cross section.

[0043] The radial velocity measurement is divided by the cosine of the geometric Doppler angle to obtain the axial velocity component. Cosine correction is a fundamental step in Doppler velocity calculation, converting the measured radial velocity into the true axial velocity through angular geometry. The correction process is based on the fundamental principle of the Doppler effect, namely, that the frequency shift detected by ultrasound is proportional to the cosine component of the blood flow direction. The correction formula is:

[0044] ;

[0045] in, This represents the velocity component along the axial direction of the blood vessel. Radial velocity measured by Doppler ultrasound. The geometric Doppler angle is used. In practical applications, when the angle approaches 90°, the cosine value approaches zero, and the correction coefficient increases significantly, amplifying measurement noise. To improve stability, the system employs an interpolation strategy within the angle range of approximately 90° ± 5°, combining measurements from adjacent time points to avoid instability caused by excessively large correction coefficients. The corrected axial velocity directly reflects the main motion component of blood flow along the vessel centerline, serving as a core indicator for hemodynamic analysis and the foundation for constructing a complete three-dimensional velocity field. This precise correction method based on the geometric Doppler angle significantly improves the accuracy and reliability of axial velocity measurement compared to empirical angle estimation in traditional ultrasound.

[0046] Extract the cross-sectional normal vector at the current vessel centerline node. Using this normal vector as a basis, establish a two-dimensional orthogonal basis vector pair in the cross-sectional plane. The construction of the transverse components is a crucial step in the three-dimensional flow field reconstruction. A velocity distribution model within the vessel cross-section is established through anatomical constraints. The construction process first determines the vessel cross-sectional plane, using the tangent direction at the centerline node as the normal. Then, an orthogonal coordinate system is established in the cross-sectional plane as a reference basis for describing the transverse velocity. The construction of the two-dimensional orthogonal basis uses the Gram-Schmidt orthogonalization method. First, an arbitrary vector that is not collinear with the tangent direction is selected (usually [0,1,0] or [1,0,0]). The first basis vector is generated through projection and normalization operations. Then, the cross product of the tangent vector and the first basis vector is calculated to obtain the second orthogonal basis vector. These two basis vectors together span the vessel cross-sectional plane, providing a natural decomposition basis for the transverse velocity components. The construction of the orthogonal basis ensures the independence and integrity of components in different directions, which is the mathematical foundation for accurately expressing the three-dimensional flow field. This method of establishing local coordinates based on vascular anatomy makes full use of prior information about vascular structures, providing an anatomical basis for lateral velocity modeling.

[0047] Based on the Poiseuille flow prior assumptions, a zero-mean Gaussian distribution is used to probabilistically model the transverse velocity component within the plane spanned by two-dimensional orthogonal basis vectors. The standard deviation is determined by the ratio of the vessel diameter to the axial velocity component. Transverse velocity modeling is an innovative step in flow field reconstruction, estimating velocity components that cannot be directly measured through physical constraints and statistical methods. The model is based on the laminar flow assumption, simplifying the flow within the vessel to a Poiseuille flow model, where the velocity distribution is parabolic, with the highest velocity at the center and approaching zero at the edges. Within the cross-section, the transverse velocity component is typically significantly smaller than the axial component and tends to zero in ideal laminar flow. Considering the complexity of actual blood flow and measurement errors, the model uses a Gaussian distribution to probabilistically represent the transverse component, with a mean of zero, and the standard deviation is determined by vessel characteristic parameters. The formula for calculating the standard deviation is:

[0048] ;

[0049] in, The standard deviation of the lateral velocity component. For the axial velocity component, This is the equivalent diameter of the blood vessel at the current location. The lateral diffusion coefficient, typically ranging from 0.1 to 0.3, reflects the degree to which actual blood flow deviates from ideal laminar flow. This probabilistic representation method based on a physical model considers both the main characteristics of blood flow and retains adaptability to complex flow patterns, enabling reasonable estimation of lateral velocity components that cannot be directly measured by ultrasound.

[0050] The expected values ​​of Gaussian sampling of the axial and lateral velocity components are combined to construct a three-dimensional vector estimate of blood flow velocity, while the corresponding covariance matrix is ​​output as an uncertainty characterization. Three-dimensional vector construction is the final step in velocity estimation, generating a complete three-dimensional velocity representation by combining the axial and lateral components. The construction process first determines the components in each direction in the local coordinate system. The axial component uses a cosine-corrected value, and the lateral component uses the expected value (mean of zero) of the probability model. Then, through coordinate transformation, the velocity components in the local coordinate system are converted into three-dimensional vectors in the patient's anatomical coordinate system. Finally, the covariance matrix is ​​calculated to quantify the uncertainty of the velocity estimate. The covariance matrix contains the variance of each component and its correlation, directly reflecting the reliability of the velocity estimate. The variance of the axial component mainly comes from Doppler measurement noise, while the variance of the lateral component is determined by the standard deviation of the probability model. This three-dimensional velocity estimation with uncertainty characterization not only provides the best estimate of the velocity but also quantifies the reliability of the estimate, providing complete information for subsequent analysis and decision-making. Meanwhile, the covariance matrix also provides a theoretical basis for the weighting strategy in the fusion registration process, realizing the organic combination of data quality and processing methods.

[0051] In this embodiment of the invention, the detailed implementation steps for constructing a Doppler-CT fused blood flow field by spatially registering the estimated three-dimensional vector value of blood flow velocity with the CT three-dimensional vascular model data in the patient's anatomical coordinate system include:

[0052] The estimated three-dimensional vector values ​​of blood flow velocity are projected onto the nearest vessel centerline node in the CT three-dimensional vascular model according to their spatial coordinates. Spatial projection is the first step in the registration process, associating discrete measurement points with the continuous vascular model through spatial proximity relationships. The projection process first calculates the shortest distance from each velocity measurement point to the vessel centerline, identifying the nearest centerline node; then, the velocity vector is assigned to that node, establishing a preliminary point-to-point mapping relationship. To improve mapping accuracy, the system adopts a weighted allocation strategy, distributing the contribution of each measurement point inversely proportionally to multiple neighboring nodes to achieve a smooth transition. When multiple measurement points are mapped to the same centerline node, the system performs a weighted average based on the time recentity of the measurement and the uncertainty of the estimate, ensuring the consistency and reliability of the node velocity values. This spatial relationship-based projection registration method effectively solves the difference between the discreteness of ultrasound sampling and the continuity of the CT model, laying the foundation for subsequent fusion processing. At the same time, the projection process also automatically filters out abnormal measurement points that are significantly deviated from the vascular path, improving the reliability of the registration results.

[0053] The projected three-dimensional blood flow velocity vector is interpolated along the vessel centerline using arc-length weighted interpolation to fill the centerline segment not covered by Doppler sampling. Interpolation filling is a necessary step in constructing a continuous flow field, using a mathematical model to compensate for data gaps in under-sampled areas. The filling process is performed in the arc-length parameterized space along the vessel centerline to ensure continuity in the vessel topology. The system employs radial basis function (RBF) interpolation, using known velocity nodes as control points to construct a continuous interpolation function. The selection of the radial basis function considers hydrodynamic characteristics, using a modified Gaussian kernel to maintain a smooth transition while respecting local features. To reflect the influence of sampling density in different regions, the system introduces an arc-length weighting mechanism, with lower weights and higher uncertainties for regions farther from known points. In bifurcation regions, the system considers blood flow conservation and splitting characteristics to ensure that the interpolation results conform to basic hydrodynamic principles. This physically constrained interpolation strategy avoids non-physical results that may be caused by simple mathematical interpolation, providing reasonable velocity estimates for areas without sampling coverage and achieving complete reconstruction of the blood flow field.

[0054] For each vessel centerline node, the interpolated blood flow velocity is multiplied by the cross-sectional flow conservation coefficient derived from the CT vessel cross-sectional area, and then corrected using the continuity equation constraint. Continuity constraint is a crucial step in ensuring the physical rationality of the flow field, adjusting the velocity distribution through the principle of flow conservation to eliminate inconsistencies. The correction process is based on the continuity equation for incompressible fluids, namely the fundamental principle that inflow equals outflow. The system first extracts the vessel cross-sectional area at each node from the CT model; then calculates the ratio of the cross-sectional areas of adjacent nodes as a baseline coefficient for velocity adjustment; finally, the interpolated velocity is adjusted proportionally to ensure spatial continuity of the flow. The adjustment formula is:

[0055] ;

[0056] in, and These represent the axial velocities of two adjacent nodes. and This represents the corresponding cross-sectional area of ​​the blood vessel. In the bifurcation region, the system considers flow distribution rules, allocating the total flow according to the diameter ratio of the branch vessels to ensure flow conservation between the main branch and the branches. This flow field constraint correction based on physical principles significantly improves the physical rationality of the fusion results, avoids flow discontinuities that may be caused by simple interpolation, and provides a reliable basis for clinical flow field analysis.

[0057] The velocity field, corrected for continuity equation constraints, is stored as a Doppler-CT fused blood flow field. Fields include node coordinates, 3D velocity vectors, and uncertainty covariance matrices. Data organization is fundamental for the efficient use of fusion results; structured storage enables a comprehensive representation of the blood flow field. The storage structure employs a node-edge topology framework, fully preserving the anatomical relationships of the vascular network. For each node, the system stores key information such as its 3D coordinates, vessel diameter, velocity vector, and covariance matrix; for each edge, it records the connected node pairs and the corresponding vessel segment attributes. The velocity vector is represented in the patient's anatomical coordinate system, facilitating cross-modal fusion and clinical interpretation; the covariance matrix quantifies the uncertainty of the velocity estimation, providing an important indicator of result reliability. Data organization considers time-series characteristics, supporting the storage and evolution analysis of multi-frame velocity fields, facilitating the assessment of hemodynamic changes within the cardiac cycle. The storage format design of the fused blood flow field balances computational efficiency and clinical application needs, supporting rapid querying, visualization rendering, and hemodynamic parameter extraction, providing a complete and convenient data foundation for subsequent clinical analysis and decision-making.

[0058] In this embodiment of the invention, the detailed implementation steps for extracting the peak flow velocity distribution at each blood vessel cross-section in the Doppler-CT fused blood flow field, verifying its consistency with the nominal Doppler angle value, and identifying abrupt changes in peak flow velocity with probe orientation include:

[0059] In the Doppler-CT fused blood flow field, cross-sections are extracted along the vessel centerline at preset intervals, and the maximum velocity vector modulus within each cross-section is extracted and recorded as the peak velocity of the cross-section. Cross-sectional analysis is a standard method for flow field assessment, providing velocity characteristics of various parts of the vessel through systematic sampling. The analysis process first involves uniformly setting sampling points along the vessel centerline, with a typical spacing of 5-10 mm to ensure sufficient coverage of key areas; then, a cross-sectional plane perpendicular to the centerline is established at each sampling point; finally, the maximum velocity modulus within the cross-section is searched as the representative value of the peak velocity at that location. For vessels with complex shapes, the system employs an adaptive sampling strategy, increasing the sampling density in areas with large curvature changes to improve the accuracy of characterization. Peak velocity, as an important indicator for hemodynamic assessment, directly reflects the local blood flow state and serves as the basis for subsequent consistency analysis. This systematic cross-sectional peak velocity extraction method provides an overview of the velocity distribution throughout the vessel, facilitating the identification of local abnormal areas and overall flow patterns.

[0060] The peak velocity sequence at each cross-section is differencing along the time dimension to obtain the time rate of change sequence of the peak velocity. Time differencing is a fundamental step in dynamic characteristic analysis, reflecting the evolution trend of parameters by calculating the changes at consecutive time points. The differencing process uses a central difference scheme, calculating the change between adjacent frames of the peak velocity time sequence at each cross-section location, and dividing by the time interval to obtain the rate of change. The formula for calculating the rate of change is:

[0061] ;

[0062] in, Let be the peak flow rate change in the i-th frame. and These are the peak flow rates for the next frame and the previous frame, respectively. The time interval between frames is defined as . Compared to forward or backward differencing, center differencing provides a more accurate derivative approximation and reduces the impact of noise. For the first and last frames of the sequence, forward or backward differencing is used as a supplement. The rate of change of time directly reflects the dynamic characteristics of the peak flow velocity; stable flow regions exhibit a lower rate of change, while regions affected by probe attitude may show significant fluctuations in the rate of change. This dynamic analysis method based on time differencing provides a sensitive indicator for capturing non-physiological fluctuations in flow velocity.

[0063] The Doppler angle nominal value change rate sequence at corresponding time points is extracted synchronously. Angle change analysis is a crucial step in detecting the influence of posture, identifying potential interference sources by tracking the temporal changes in angle parameters. Similar to flow velocity differential analysis, the extraction process involves centrally differencing the Doppler angle nominal value sequence to calculate the angle change rate at each time point. The nominal values ​​are derived from the internal estimates of the wireless ultrasound equipment, reflecting angle changes adjusted by the operator or automatically tracked by the system. The angle change rate calculation employs special processing in the angle domain, considering the periodicity of the angle and avoiding the impact of discontinuities near 0° / 360°. This synchronous extraction strategy ensures the temporal correspondence between angle changes and flow velocity changes, providing paired data for subsequent correlation analysis. The angle change rate sequence reflects the dynamic characteristics of probe operation and is an important basis for identifying the influence of human factors. By comparing the temporal patterns of angle changes and flow velocity changes, the system can effectively distinguish between physiological blood flow changes and measurement fluctuations caused by operational factors.

[0064] Calculate the cross-correlation coefficient between the peak flow velocity time rate of change series and the Doppler angle nominal value rate of change series. Cross-correlation analysis is a statistical method for assessing the correlation between variables, quantifying their dependence by calculating the time correlation between the series. The analysis process uses the standard cross-correlation function to evaluate the degree of correlation between the two series under different time lag conditions. The formula for calculating the cross-correlation coefficient is:

[0065] ;

[0066] in, Lag time The cross-correlation coefficients under the following conditions and The first The rate of change of frame velocity and the rate of change of angle. and This represents the mean of the corresponding sequence. The total sequence length is given. The calculation considers a lag range of ±5 frames to capture potential time delay effects. A high cross-correlation coefficient indicates a significant statistical association between flow velocity and angle changes, potentially meaning that the measurement results are affected by probe attitude changes. This time-correlation-based statistical analysis method can objectively assess the degree of interference from attitude changes on the measurement results, providing a quantitative basis for identifying unstable regions.

[0067] When the cross-correlation coefficient at a certain cross-section exceeds a preset correlation threshold, the cross-section is determined to have an angle-dependent jump and is marked as an abnormal region. Threshold judgment is the final step in anomaly detection, identifying areas requiring special attention by setting objective standards. The judgment process assesses whether the maximum value of the flow velocity-angle cross-correlation coefficient exceeds the preset threshold for each vessel cross-section. The threshold is typically set in the range of 0.7-0.8, determined based on extensive clinical validation data, balancing sensitivity and specificity requirements. When the cross-correlation coefficient exceeds the threshold, the system determines that the flow velocity measurement at that cross-section has a significant angle dependence, which may lead to non-physiological measurement fluctuations. These cross-sections are marked as abnormal regions requiring special attention in subsequent processing. The system also considers cross-correlation significance testing to ensure the statistical reliability of the correlation results and avoid accidental associations in small sample data. This statistical threshold-based anomaly detection method provides objective and consistent evaluation criteria, effectively identifying key areas requiring stabilization and providing precise location for subsequent iterative corrections.

[0068] In this embodiment of the invention, the detailed implementation steps of the iterative correction module, which introduces attitude sensitivity weights into the probe attitude matrix based on the spatial location of the jump abnormal region, iteratively corrects the geometric Doppler angle, and generates stable blood flow quantitative results, include:

[0069] For each abrupt change region, the standard deviation of the geometric Doppler angles over multiple frames within that region is calculated and denoted as angular instability. Angle stability analysis is the first step in formulating a correction strategy, determining the correction requirements for different regions by quantifying the degree of angular fluctuation. The analysis process extracts the geometric Doppler angle data for each identified abrupt change region over time and calculates its standard deviation as an instability index. The formula for calculating the standard deviation is:

[0070] ;

[0071] in, For angular instability, For the first The geometric Doppler angle of the frame, The mean of the angle sequence. The standard deviation represents the sequence length. A higher standard deviation indicates significant fluctuations in angle measurements within this region, requiring stronger stabilization. In addition to the standard deviation, the system also considers auxiliary indicators such as the distribution of extreme angle values ​​and the frequency of change to comprehensively assess the characteristics and degree of angle instability. This statistically based instability quantification method provides an objective basis for subsequent adaptive weight design, ensuring that the correction strength matches actual needs and avoiding problems of over-smoothing or under-correction.

[0072] Angle instability is mapped to attitude sensitivity weights within the interval (0,1], with larger instability resulting in smaller weights. Weight mapping is the core mechanism of adaptive correction, transforming statistical indicators into practical control parameters through nonlinear transformation. The mapping uses a negative exponential function to ensure that the weights decrease smoothly with increasing instability while remaining within the effective interval. The mapping function is:

[0073] ;

[0074] in, For attitude sensitivity weights, For angular instability, The mapping coefficient, typically set to 1.5-3, controls the steepness of the weight descent. The mapping result is limited to the (0,1) interval to ensure that all regions retain a certain level of real-time responsiveness, while imposing stronger smoothing constraints on unstable regions. For extremely unstable regions (standard deviation exceeding a preset upper limit), the system sets a minimum weight threshold to avoid excessive suppression of the impact of real-time data. This adaptive weighting strategy based on exponential mapping achieves a precise match between the correction strength and the degree of instability, providing differentiated stabilization processing for different regions and balancing the needs of measurement stability and real-time responsiveness.

[0075] The corrected geometric Doppler angle is obtained by weighting the current frame's geometric Doppler angle with the exponentially moving average of the historical frame's geometric Doppler angle using attitude sensitivity weights. Temporal fusion is the core step in the stabilization process, smoothing the current measurement with historical information and reducing the impact of instantaneous fluctuations. The fusion process uses an exponentially moving average model, which considers both historical accumulation and the influence of new data. The correction formula is:

[0076] ;

[0077] in, The corrected geometric Doppler angle for the current frame. The original geometric Doppler angle of the current frame. The angle is the exponential moving average of the previous frame. The weights are for attitude sensitivity. The exponential moving average is updated recursively over time, and the result after each fusion serves as the historical reference for the next frame. This recursive calculation avoids the need to store the complete historical sequence, while automatically adjusting the influence range of historical data through weight decay. For systems that have just started or have experienced a long interruption, a warm-up strategy is adopted to gradually establish a reliable historical reference, avoiding instability in the initial stage. This weighted fusion-based correction method effectively suppresses non-physical fluctuations in angle while preserving measurement responsiveness, providing a more reliable angle correction basis for blood flow quantification.

[0078] The axial velocity component is recalculated by substituting the corrected geometric Doppler angle into the cosine correction step, and the three-dimensional velocity vector estimate is updated. Velocity recalculation is the implementation step in the correction process, optimizing the velocity estimation result by updating the angle parameters. The recalculation process is consistent with the initial vector estimation procedure. First, cosine correction is performed on the radial velocity based on the corrected geometric Doppler angle to obtain a stable axial velocity component; then, the probability distribution parameters of the lateral component are adjusted according to the updated axial value; finally, the corrected three-dimensional velocity vector and its uncertainty characterization are generated. The recalculation maintains complete physical model constraints, ensuring that the correction result conforms to fluid dynamics principles, while significantly improving the reliability of the velocity estimate through angle stabilization. Compared to directly smoothing velocity values, this angle-correction-based method is more consistent with the physical nature of Doppler measurements, avoiding potential inconsistencies between different velocity components. The corrected velocity estimate not only retains the basic characteristics of the original measurement but also eliminates non-physical jumps caused by attitude fluctuations, providing more accurate and reliable quantitative blood flow results for clinical assessment.

[0079] The correction process is repeated until the change in peak flow velocity in adjacent iterations is less than a preset convergence threshold. The final result is then output as a stable quantitative blood flow result. Iterative optimization is the mechanism to ensure the effectiveness of the correction, achieving incremental optimization of the results through multiple rounds of adjustments. The iterative process uses the corrected angle and velocity as input for the next round, repeatedly executing weight calculation, angle correction, and velocity update steps until convergence is achieved. Convergence is determined based on the change in peak flow velocity in adjacent iterations. When the change at all monitoring points is less than a preset threshold (usually 1-3%), the correction process is considered to have reached a stable state. To avoid over-iteration, the system sets a maximum iteration limit (usually 5-10 times), terminating the iteration process early in special circumstances. The final stable quantitative blood flow result includes the corrected geometric Doppler angle, three-dimensional velocity vector, and uncertainty assessment, providing high-quality hemodynamic data for clinical applications. This iterative optimization strategy based on convergence criteria ensures the reliability and efficiency of the correction process, improving result stability while avoiding information loss due to over-processing, providing an optimal solution for balancing accuracy and reliability in blood flow assessment.

[0080] In this embodiment of the invention, the detailed implementation steps of substituting the current zero-bias estimate into a first-order Markov drift model to predict the zero-bias value at the next time step, and then weighted and fused with the measured zero-bias estimate to obtain a smoothed zero-bias sequence include:

[0081] The time correlation coefficient of the first-order Markov model is set by fitting the autocorrelation function of the historical zero-biased series. Setting the correlation coefficient is a crucial step in model parameterization, determining the temporal characteristics of the zero-biased drift through a data-driven approach. The process involves first collecting historical zero-biased data for long-term series, covering different usage scenarios and environmental conditions; then calculating the autocorrelation function of the series at different time lags; and finally, fitting the autocorrelation curve using an exponential decay model to extract the time correlation coefficient. The formula for calculating the autocorrelation function is:

[0082] ;

[0083] in, For time lag The autocorrelation coefficient, For a moment The zero bias value, The mean of the zero-biased sequence. Let be the sequence length. The fitting model uses an exponential form. ,in For characteristic correlation time, the time correlation coefficient The relationship is , The sampling time interval is defined as 0.9–0.99. Typical time correlation coefficients range from 0.9 to 0.99, reflecting the slow nature of zero-bias drift. The system supports multiple correlation coefficient settings to accommodate the characteristics of different sensor types and environmental conditions. This parameter determination method based on measured data ensures the match between the model and actual sensor characteristics, providing a reliable statistical basis for drift prediction.

[0084] The predicted zero-bias value and prediction variance for the next time step are calculated using the time correlation coefficient and the current zero-bias estimate. Zero-bias prediction is a core step in the application of the drift model, predicting the sensor state at future times through a recursive relationship. The prediction process is based on the state transition equation of a first-order Markov process, combining the current zero-bias estimate and prior model parameters to calculate the predicted value for the next time step. The prediction equation is:

[0085] ;

[0086] in, For based on Time information The zero bias of the prediction at time t. This is the current zero-bias estimate. The initial zero bias value (usually the long-term mean) is taken. This represents the time correlation coefficient. The prediction variance is calculated using the error propagation rule, taking into account both the current estimation uncertainty and process noise.

[0087] ;

[0088] in, To predict the variance of zero bias, The variance of the current zero-bias estimate, The process noise variance describes the random fluctuations of the system. The prediction variance directly reflects the level of prediction uncertainty and provides a basis for weight allocation in subsequent fusion steps. This prediction method based on a state-space model considers both the long-term trend of zero bias and retains the responsiveness to short-term changes, providing a theoretical framework for stable and sensitive zero-bias tracking.

[0089] The Kalman gain formula is used to weight the predicted zero-bias value and the measured zero-bias estimate according to the ratio of the prediction variance to the observation noise variance, thus obtaining the smoothed zero-bias value at that moment. Kalman fusion is a standard method for optimal estimation, integrating prediction and observation information through mathematical optimality criteria. The fusion process first calculates the Kalman gain as a weighting coefficient between the predicted and measured values; then, a weighted average is performed based on the gain values ​​to obtain the updated zero-bias estimate; finally, the estimation variance is updated to provide a basis for prediction at the next moment. The Kalman gain calculation formula is:

[0090] ;

[0091] in, For a moment Kalman gain, To predict variance, To describe the uncertainty of the measured zero-bias estimate, we use the observation noise variance. The zero-bias update formula is:

[0092] ;

[0093] in, For the updated zero-biased estimate, To predict the zero bias, This is the measured zero-biased estimate. The variance update formula is:

[0094] ;

[0095] Kalman gain automatically adjusts the relative weights of predicted and measured values, increasing the prediction weight when the prediction is accurate and increasing the measured weight when the prediction is uncertain, thus achieving dynamic optimal fusion. This fusion method based on a Bayesian framework demonstrates superior performance in dealing with zero-bias variations and measurement noise, providing the system with a smooth and accurate zero-bias estimation sequence.

[0096] The smoothed zero-bias values ​​are arranged chronologically to form a smoothed zero-bias sequence. Sequence construction is the final step in data organization, providing a complete history of zero-bias evolution through temporal sorting. The construction process sorts the smoothed zero-bias values ​​at each moment according to the acquisition timestamp, forming a continuous time series. For multi-axis sensors, the system processes the zero-bias data for each axis separately, constructing a three-dimensional zero-bias vector sequence. The smoothed sequence not only records the instantaneous values ​​of the zero-bias but also preserves the corresponding uncertainty estimates, providing complete information for subsequent compensation. In sequence construction, the system addresses practical issues such as uneven timestamps and missing data, ensuring the continuity and integrity of the sequence through interpolation and extrapolation methods. For long-term operations, the system implements sliding window management, controlling memory usage while retaining key historical information, supporting long-term stable operation. This method of completely recording zero-bias evolution provides a detailed zero-bias benchmark for subsequent attitude compensation, ensuring the long-term accuracy of the attitude matrix and creating conditions for multimodal fusion.

[0097] In this embodiment of the invention, based on the Poiseuille flow prior assumption, a zero-mean Gaussian distribution is used to probabilistically model the transverse velocity component in the plane spanned by two-dimensional orthogonal basis vector pairs, with the standard deviation determined by the ratio of the vessel diameter to the axial velocity component. Detailed implementation steps include:

[0098] The equivalent circular cross-sectional diameter of the blood vessel at the current node is extracted from CT 3D blood vessel model data. Diameter extraction is a fundamental step in blood vessel geometry analysis, obtaining lumen size parameters through cross-sectional processing. The extraction process first establishes a cross-sectional plane perpendicular to the tangent direction at the blood vessel centerline node; then, the intersection of the CT voxel data and this plane is calculated to obtain the blood vessel contour; finally, the equivalent circular diameter is calculated based on the contour shape. The calculation of the equivalent circular diameter considers the actual area and perimeter of the blood vessel cross-section and handles common cases of non-circular cross-sections. The calculation formula is:

[0099] ;

[0100] in, The equivalent circular cross-section diameter is This represents the actual area of ​​the vessel cross-section. For vessels with complex shapes, the system also calculates ellipticity and irregularity indices, evaluates the applicability of the equivalent circle approximation, and employs a more complex cross-sectional model in highly non-circular regions. Vessel diameter information is a key parameter for flow pattern determination, directly affecting the morphological characteristics of the velocity distribution and the probabilistic model of the lateral component, providing crucial geometric constraints for flow field reconstruction.

[0101] The ratio of the axial velocity component to the diameter of the equivalent circular section is denoted as the velocity-to-diameter ratio. Velocity-to-diameter ratio calculation is a core step in flow characteristic analysis, assessing the fundamental properties of the flow through dimensionless parameters. The calculation uses the estimated axial velocity at the current node and the extracted equivalent circular diameter to form the velocity-to-diameter ratio index. The calculation formula is:

[0102] ;

[0103] in, The velocity-to-diameter ratio, For the axial velocity component, The velocity-to-diameter ratio (VDR) is the equivalent circular cross-section diameter. In fluid mechanics, the VDR is crucial, closely related to the Reynolds number, and directly reflects the relative strength of inertial and viscous forces, serving as a key indicator for determining flow type (laminar or turbulent). In blood flow, the VDR is also closely related to the shear rate, affecting the viscoelastic properties of blood flow. The system dynamically adjusts the flow model parameters based on estimates of the VDR and the local Reynolds number to adapt to the flow characteristics of different vascular regions. This flow characteristic analysis based on physical parameters provides a reasonable theoretical foundation for probabilistic modeling of lateral velocities, ensuring consistency between the model and actual blood flow characteristics.

[0104] Multiplying the velocity-to-diameter ratio by a preset lateral diffusion coefficient yields the standard deviation of the Gaussian distribution of the lateral velocity components. Setting the standard deviation is a crucial step in probabilistic model parameterization, determining the amplitude range of lateral fluctuations through fluid dynamic relationships. The setting process is based on the velocity-to-diameter ratio, combined with the preset lateral diffusion coefficient, to calculate the probability distribution parameters of the lateral velocity components. The calculation formula is:

[0105] ;

[0106] in, The standard deviation of the lateral velocity component. The lateral diffusion coefficient is... Velocity-to-diameter ratio. Lateral diffusion coefficient. This reflects the degree to which actual blood flow deviates from ideal laminar flow, with typical values ​​between 0.1 and 0.3, determined through extensive phantom experiments and computational fluid dynamics simulations. The system supports dynamic adjustment of the diffusion coefficient based on local flow characteristics, increasing the coefficient value in high Reynolds number or post-turbulent regions to better capture complex flow patterns. The standard deviation directly determines the possible range of variation of the transverse velocity component and is a key parameter for three-dimensional flow field reconstruction. This parameter setting method based on a physical model considers both the basic characteristics of blood flow and retains adaptability to complex flow patterns, providing theoretical support for accurately reconstructing the transverse component that cannot be directly measured by ultrasound.

[0107] When the equivalent circular cross-section diameter is smaller than a preset minimum diameter threshold, the standard deviation is forcibly set to zero, degenerating into a pure axial flow assumption to avoid divergence in the lateral component estimation in small vessels. Boundary condition handling is a mechanism to ensure model stability, avoiding ill-conditioned solutions by identifying and handling special cases. The processing sets a minimum diameter threshold, typically 1-2 mm, reflecting the significant changes in flow characteristics in small vessels. When the vessel diameter is smaller than the threshold, the Reynolds number decreases significantly, the flow becomes highly laminar, and the lateral component is almost negligible. Simultaneously, the relative error of velocity estimation increases in small vessels, potentially leading to instability in the lateral component estimation based on proportional relationships. In this case, the system forcibly sets the standard deviation of the lateral velocity component to zero, degenerating the model into a pure axial flow assumption, avoiding unnecessary complexity and potential estimation divergence. For medium-diameter vessels, the system implements a smooth transition strategy, gradually reducing the standard deviation as the diameter decreases to ensure the continuity and stability of the model in vessels of different scales. This adaptive boundary condition handling method effectively balances model complexity and result reliability, providing a unified framework for flow field reconstruction of full-scale vascular networks.

[0108] In this embodiment of the invention, the detailed implementation steps for correcting the continuity equation constraint by multiplying the interpolated blood flow velocity at each vessel centerline node by the cross-sectional flow conservation coefficient derived from the CT vessel cross-sectional area include:

[0109] The cross-sectional areas of blood vessels at adjacent nodes are extracted from CT 3D vascular model data and denoted as the upstream and downstream cross-sectional areas. Cross-sectional area extraction is a fundamental step in flow analysis, providing geometric parameters for flow calculation through precise measurement. The extraction process is similar to diameter extraction: a cross-sectional plane perpendicular to the local tangent is established at adjacent nodes along the vessel centerline; the intersection of the vessel contour and the plane is calculated; and the actual area is calculated based on the closed contour. The calculation employs either a polygon area formula or a pixel counting method, selecting the most suitable algorithm based on the CT data resolution and vessel size. To improve accuracy, the system applies morphological processing to the extracted contour to eliminate noise and artifacts, ensuring the reliability of the area calculation. For low-resolution regions, a model-assisted inference strategy is used, combining prior knowledge of vessel diameter and shape to estimate possible cross-sectional area variations. The upstream and downstream concepts are defined based on blood flow direction, using Doppler ultrasound data or anatomical knowledge to determine the local blood flow direction and clarify the relative positional relationship. This precise cross-sectional area extraction method provides the necessary geometric foundation for flow conservation analysis and is a prerequisite for the application of the continuity equation.

[0110] Calculate the ratio of the upstream cross-sectional area to the downstream cross-sectional area, denoted as the rate of change of cross-section. This ratio calculation is a core step in flow relationship analysis, reflecting the flow rate variation pattern in the pipeline through geometric parameter proportions. The calculation uses a simple division operation to obtain the rate of change of cross-section between adjacent nodes:

[0111] ;

[0112] in, The rate of change of cross section The upstream cross-sectional area, The ratio represents the downstream cross-sectional area. A ratio greater than 1 indicates that the blood vessel narrows from upstream to downstream, while a ratio less than 1 indicates that the blood vessel dilates. The rate of change of cross-section is a direct parameter for applying the principle of flow conservation, reflecting the proportional relationship that the velocity needs to be adjusted. The system handles special cases in ratio calculation, such as unstable division caused by extremely small cross-sectional areas. By setting a minimum area threshold and a ratio limit range, abnormal ratios are avoided. For blood vessel bifurcation regions, the system employs a special processing strategy, considering complex situations such as one inflow and multiple outflows or multiple inflows and one outflow, ensuring the application of the principle of total flow conservation. This ratio analysis method based on geometric changes provides a direct mathematical basis for velocity adjustment and is a core parameter for the application of the continuity equation.

[0113] Multiplying the interpolated axial velocity component at the upstream node by the rate of change of cross section yields the corrected velocity at the downstream node, satisfying the continuity equation for incompressible fluids. Velocity correction is a step in implementing the continuity constraint, applying the principle of flow conservation through mathematical transformations. The correction process is based on the fundamental assumptions of incompressible fluids: fluid density remains constant, and flow rate (volume flow rate) is conserved. According to the continuity equation, the volume flow rate is equal across any cross section:

[0114] ;

[0115] Therefore, the formula for calculating the downstream velocity is derived:

[0116] ;

[0117] in, For the corrected axial velocity of the downstream node, The interpolated axial velocity of the upstream node. The cross-sectional change rate is denoted as . The correction is applied only to the axial velocity component; the lateral component is re-estimated according to the new axial value using a probabilistic model. Velocity correction ensures spatial continuity of blood flow, avoiding physical inconsistencies that may result from simple interpolation. For long vessel segments, the system employs a segment-by-segment cumulative correction strategy, applying continuity constraints node by node from the inlet to ensure overall flow conservation. This physics-based velocity correction method significantly improves the physical plausibility of the fused blood flow field, providing a reliable foundation for subsequent hemodynamic analysis.

[0118] When the rate of change of the cross section exceeds a preset threshold, the system determines that the segment contains vascular stenosis or aneurysmal dilatation, triggering local turbulence marking. The system then amplifies each component of the uncertainty covariance matrix of the corresponding node by a preset factor to mark areas of reduced reliability in the Doppler-CT fused blood flow field. Abnormal region identification is a crucial step in assessing the reliability of the results, identifying potentially complex flow regions through geometric anomalies. The identification process sets a threshold for the rate of change of the cross section, typically ±30%; values ​​exceeding this range are considered significant morphological changes. Rapidly narrowing regions may correspond to vascular stenosis, leading to increased local flow velocity and streamline distortion; rapidly dilating regions may correspond to aneurysms or post-stenotic areas, often exhibiting flow separation and vortex structures. The flow patterns in these regions significantly deviate from laminar flow assumptions, and simple axial-lateral decomposition models may no longer be applicable. The system adds local turbulence markings to the identified abnormal regions and reflects the reduced reliability of the estimation results by increasing the amplification factor of the uncertainty matrix (typically 2-5 times). The amplification of the covariance matrix directly affects the weight of the data in this region in subsequent analyses, automatically reducing the impact of unreliable data. This flow complexity assessment method based on geometric features provides an important reliability indicator for fusion results, helping clinicians to correctly interpret blood flow data and avoid over-reliance on and misjudgment of low-reliability regions.

[0119] This invention achieves precise registration and fusion of wireless ultrasound Doppler blood flow data and CT three-dimensional vascular models by integrating a comprehensive acquisition module, an attitude compensation module, a Doppler angle calculation module, a vector estimation module, a spatial registration module, a consistency verification module, and an iterative correction module. The multimodal registration method of this invention can effectively compensate for the influence of ultrasound probe attitude changes, providing stable and accurate three-dimensional vector estimation of blood flow velocity.

[0120] The above are merely preferred embodiments of the present invention and are not intended to limit the present invention. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

[0121] It should be noted that all formulas in this manual are calculated by removing dimensions and taking their numerical values. The formulas are derived from software simulations based on a large amount of collected data to obtain the most recent real-world results. The preset parameters and thresholds in the formulas are set by those skilled in the art according to the actual situation.

[0122] Although embodiments of the invention have been shown and described, those skilled in the art will understand that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the claims and their equivalents.

Claims

1. A wireless ultrasound multimodal image registration and fusion system, characterized in that, include: The integrated acquisition module is used to acquire Doppler blood flow data from the wireless ultrasound probe, attitude quaternion sequences output by the inertial navigation unit, and CT three-dimensional vascular model data. The Doppler blood flow data includes radial velocity measurements along the sound beam axis and corresponding nominal Doppler angle values. The attitude compensation module is used to perform zero-bias drift compensation on the attitude quaternion sequence and extract the real-time rotation matrix from the compensated probe coordinate system to the patient anatomical coordinate system, which is denoted as the probe attitude matrix. The Doppler angle calculation module is used to extract the tangent direction vector of the blood vessel centerline based on the CT three-dimensional blood vessel model data, and project the sound beam direction vector onto the patient's anatomical coordinate system in combination with the probe posture matrix, and calculate the measured angle between the sound beam direction and the blood vessel tangent direction, which is denoted as the geometric Doppler angle. The vector estimation module is used to perform cosine correction on the radial velocity measurement value based on the geometric Doppler angle, obtain the axial velocity component along the tangent direction of the blood vessel, and construct a probabilistic constraint model of the transverse velocity component based on the normal vector of the blood vessel section to generate a three-dimensional vector estimate of the blood flow velocity. The spatial registration module is used to spatially register the estimated three-dimensional vector value of blood flow velocity with the three-dimensional CT vascular model data in the patient anatomical coordinate system to construct a Doppler-CT fused blood flow field. The consistency verification module is used to extract the peak flow velocity distribution at each blood vessel cross section in the Doppler-CT fused blood flow field, perform consistency verification with the nominal value of the Doppler angle, and identify abnormal regions where the peak flow velocity changes with the probe posture. The iterative correction module is used to introduce attitude sensitivity weights into the probe attitude matrix according to the spatial location of the jump abnormal region, and to iteratively correct the geometric Doppler angle to generate a stable blood flow quantitative result.

2. The system according to claim 1, characterized in that, The zero-bias drift compensation is performed on the attitude quaternion sequence, and the real-time rotation matrix from the compensated probe coordinate system to the patient anatomical coordinate system is extracted and denoted as the probe attitude matrix, including: The attitude quaternion sequence is segmented according to a time window, and the mean value of the gyroscope output is calculated during each static period and recorded as the current zero bias estimate. The current zero-bias estimate is substituted into the first-order Markov drift model to predict the zero-bias value at the next moment, and then weighted and fused with the measured zero-bias estimate to obtain a smooth zero-bias sequence. The smoothed zero-bias sequence is subtracted frame by frame from the attitude quaternion sequence to obtain the zero-bias compensated quaternion sequence; The quaternion sequence after zero bias compensation is converted into a 3×3 rotation matrix and output as the probe attitude matrix.

3. The system according to claim 1, characterized in that, The process involves extracting the tangential direction vector of the vessel centerline based on the CT three-dimensional vascular model data, projecting the sound beam direction vector onto the patient's anatomical coordinate system using the probe posture matrix, and calculating the measured angle between the sound beam direction and the vascular tangential direction, denoted as the geometric Doppler angle, including: The centerline extraction algorithm is executed on the CT three-dimensional blood vessel model data. A node is sampled every preset arc length along the centerline, and cubic spline fitting is performed on the line connecting adjacent nodes to obtain the tangent direction vector at each sampled node. The nominal direction vector of the sound beam emitted from the wireless ultrasound probe in the probe coordinate system is multiplied by the probe attitude matrix to obtain the direction vector of the sound beam in the patient anatomical coordinate system. Calculate the inner product of the sound beam direction vector and the tangent direction vector at the blood vessel centerline node closest to the probe projection point, and take the inverse cosine to obtain the geometric Doppler angle; When the geometric Doppler angle exceeds the preset angle range, an angle validity flag is triggered, and the Doppler data of the corresponding frame is set to an invalid state.

4. The system according to claim 1, characterized in that, The process of cosine-correcting the radial velocity measurement based on the geometric Doppler angle to obtain the axial velocity component along the tangential direction of the blood vessel, and constructing a probabilistic constraint model for the transverse velocity component based on the normal vector of the blood vessel cross section to generate a three-dimensional vector estimate of the blood flow velocity includes: Divide the radial velocity measurement by the cosine of the geometric Doppler angle to obtain the axial velocity component; Extract the cross-sectional normal vector at the current vessel centerline node, and establish a two-dimensional orthogonal basis vector pair in the cross-sectional plane using this normal vector as the basis; In the plane spanned by the two-dimensional orthogonal basis vectors, the transverse velocity component is probabilistically modeled using a zero-mean Gaussian distribution, and the standard deviation is determined by the ratio of the vessel diameter to the axial velocity component. The Gaussian sampling expectation values ​​of the axial velocity component and the lateral velocity component are combined to form the three-dimensional vector estimate of the blood flow velocity, and the corresponding covariance matrix is ​​output.

5. The system according to claim 1, characterized in that, The step of spatially registering the estimated three-dimensional blood flow velocity vector with the CT three-dimensional vascular model data in the patient's anatomical coordinate system to construct a Doppler-CT fused blood flow field includes: The estimated three-dimensional vector value of blood flow velocity is projected onto the nearest vessel centerline node corresponding to the CT three-dimensional vascular model data according to its spatial coordinates. For the projected three-dimensional vector of blood flow velocity, interpolate along the direction of the vessel centerline by arc length weighting to fill the centerline segment not covered by Doppler sampling; For each blood flow velocity at the centerline node of the blood vessel, the interpolated blood flow velocity is multiplied by the cross-sectional flow conservation coefficient derived from the CT blood vessel cross-sectional area, and the continuity equation is constrained and corrected. The velocity field, after being corrected by the continuity equation constraint, is stored as the Doppler-CT fused blood flow field, with fields including node coordinates, three-dimensional velocity vector, and uncertainty covariance matrix.

6. The system according to claim 1, characterized in that, The peak flow velocity distribution at each vessel cross-section in the extracted Doppler-CT fused blood flow field is compared with the nominal Doppler angle value for consistency verification, and abnormal regions of peak flow velocity abrupt changes with probe orientation are identified, including: In the Doppler-CT fused blood flow field, cross sections are cut along the blood vessel centerline at preset intervals, and the maximum value of the velocity vector magnitude in each cross section is extracted and recorded as the peak flow velocity of the cross section. The peak flow velocity sequence of the cross section is differentially divided along the time dimension to obtain the peak flow velocity time change rate sequence; The sequence of the rate of change of the nominal Doppler angle at the corresponding time point is extracted simultaneously; Calculate the cross-correlation coefficient between the peak flow velocity time change rate sequence and the Doppler angle nominal value change rate sequence; When the cross-correlation coefficient at a certain cross section exceeds a preset correlation threshold, it is determined that the cross section has an angle-dependent jump and is marked as the jump abnormal region.

7. The system according to claim 1, characterized in that, The step involves introducing attitude sensitivity weights into the probe attitude matrix based on the spatial location of the abrupt change region, iteratively correcting the geometric Doppler angle, and generating stable blood flow quantification results, including: For each of the aforementioned abrupt change regions, the standard deviation of the geometric Doppler angles over multiple frames within that region is calculated and denoted as the angular instability. The angular instability is mapped to attitude sensitivity weights in the interval (0,1). The corrected geometric Doppler angle is obtained by weighting and fusing the exponentially moving average of the geometric Doppler angles of the current frame and the historical frames with the attitude sensitivity weights. Substitute the corrected geometric Doppler angle into the cosine correction step to recalculate the axial velocity component and update the estimated three-dimensional vector value of the blood flow velocity; Repeat the above correction process until the change in peak flow rate in adjacent iterations is less than the preset convergence threshold, and output the final result as the stabilized blood flow quantification result.

8. The system according to claim 2, characterized in that, The step of substituting the current zero-bias estimate into a first-order Markov drift model to predict the zero-bias value at the next time step, and then weighting and fusing it with the measured zero-bias estimate to obtain a smoothed zero-bias sequence includes: The time correlation coefficient of the first-order Markov model is set, which is obtained by fitting the autocorrelation function of the historical zero-biased sequence; The predicted zero bias and prediction variance for the next time step are calculated using the time correlation coefficient and the current zero bias estimate. Using the Kalman gain formula, the predicted zero bias value and the measured zero bias estimate are weighted according to the ratio of the prediction variance to the observation noise variance to obtain the smoothed zero bias value at that moment. The frame-by-frame smoothed zero-bias values ​​are arranged in chronological order to form the smoothed zero-bias sequence.

9. The system according to claim 4, characterized in that, In the plane spanned by the two-dimensional orthogonal basis vectors, a zero-mean Gaussian distribution is used to probabilistically model the transverse velocity component, with the standard deviation determined by the ratio of the vessel diameter to the axial velocity component, including: Extract the equivalent circular cross-sectional diameter of the blood vessel at the current node from the CT three-dimensional blood vessel model data; Calculate the ratio of the axial velocity component to the diameter of the equivalent circular cross section, and denote it as the velocity-to-diameter ratio; Multiply the velocity-to-diameter ratio by a preset lateral diffusion coefficient to obtain the standard deviation of the Gaussian distribution of the lateral velocity component. When the diameter of the equivalent circular cross section is less than the preset minimum diameter threshold, the standard deviation is forcibly set to zero, degenerating into a pure axial flow assumption, in order to avoid divergence in the estimation of the lateral component in small blood vessels.

10. The system according to claim 5, characterized in that, The interpolated blood flow velocity at each vessel centerline node is multiplied by the cross-sectional flow conservation coefficient derived from the CT vessel cross-sectional area, and the continuity equation constraint correction is performed, including: Extract the cross-sectional area of ​​the blood vessel at two adjacent nodes from the CT three-dimensional blood vessel model data, and record it as the upstream cross-sectional area and the downstream cross-sectional area; Calculate the ratio of the upstream cross-sectional area to the downstream cross-sectional area, and denot it as the rate of change of cross-section; Multiply the interpolated axial velocity component at the upstream node by the cross-sectional change rate to obtain the corrected velocity at the downstream node that satisfies the incompressible fluid continuity equation. When the rate of change of the cross section exceeds a preset rate of change threshold, it is determined that there is vascular stenosis or aneurysm-like dilation in the section, triggering local turbulence marking, and expanding each component of the uncertainty covariance matrix of the corresponding node by a preset magnification factor to mark the area of ​​reduced confidence of quantitative results in the Doppler-CT fused blood flow field.

Images