Method and system for reconstructing a three-dimensional image based on two-dimensional x-ray images

Abstract

The application discloses a method and system for reconstructing three-dimensional images based on two-dimensional X-ray images, relates to the technical field of medical images, and comprises a micro X-ray tube array arranged according to a known geometric rule, each tube being independently controllable for exposure; a high-voltage control module providing each tube with an independently stable and switchable high-voltage power supply, supporting sequential pulse exposure; an image acquisition module receiving and converting X-ray signals into multiple frames of original two-dimensional digital images; and an image processing module registering multiple frames of images according to the geometric positions of the tubes, fusing to generate a single composite two-dimensional image, and finally reconstructing a three-dimensional image based on the composite image and multi-view geometric information of the tube array. The technology replaces mechanical movement with a static tube array to obtain multi-view data, and simplifies data processing and improves imaging efficiency and quality through a new process of fusing first and reconstructing later.

Description

Technical Field

[0001] This invention belongs to the field of medical imaging technology, specifically a method and system for reconstructing three-dimensional images based on two-dimensional X-ray images. Background Technology

[0002] In the field of X-ray 3D imaging, conventional techniques primarily rely on a single X-ray tube. This tube must move around the object being inspected using a complex mechanical rotation mechanism to acquire projection data from different angles. This method's data acquisition speed is limited by the mechanical motion, and the system suffers from problems such as large size, high cost, and reliability affected by moving parts. Mechanical motion can also introduce vibrations, leading to errors in the geometric relationships between projected images. Another technical solution is to use multiple fixed X-ray sources to attempt to eliminate mechanical motion. However, in existing solutions, multiple sources often operate simultaneously, and the generated X-rays interfere with each other, resulting in severe scattering noise and reducing the contrast and quality of the projected images. Attempting to have multiple sources operate sequentially in a time-sharing manner places extremely high demands on the high-voltage power supply system, requiring it to achieve extremely fast and stable switching between different sources.

[0003] When reconstructing a 3D image from a series of acquired 2D projected images, the conventional method is to directly input these raw images into a reconstruction algorithm, such as filtered backprojection or iterative reconstruction. These algorithms directly process large amounts of raw data, making the computation complex and time-consuming. Algorithm performance is highly dependent on the precise geometric correspondence between the projected images and the quality of the images themselves. In existing technologies, geometric position errors caused by mechanical motion or image quality degradation due to simultaneous exposure from multiple sources are directly incorporated into the reconstruction process. These errors are amplified by the reconstruction algorithm, ultimately forming various artifacts in the 3D image, impairing the accuracy and clarity of the imaging, and limiting the reliable application of 3D imaging. Summary of the Invention

[0004] This invention aims to solve at least one of the technical problems existing in the prior art; Therefore, this invention proposes a system for reconstructing three-dimensional images based on two-dimensional X-ray images, comprising: The X-ray tube array module contains multiple miniature X-ray tubes arranged in a known geometric pattern. The focal points of all the miniature X-ray tubes together cover the target imaging area, and each miniature X-ray tube can independently control the exposure. A high-voltage control module is used to provide an independent and stable high-voltage power supply for each miniature X-ray tube, supporting the sequential pulsed exposure of the tube array module; The image acquisition module is used to receive X-ray signals after penetrating the object being detected, and to convert the X-ray signals into multiple frames of raw two-dimensional digital image signals. An image processing module is used to receive the multiple frames of original two-dimensional digital image signals and process the multiple frames of original two-dimensional digital image signals. The image processing module processes the multi-frame original two-dimensional digital image signals, specifically including: image registration of the images corresponding to the multi-frame original two-dimensional digital image signals according to the known geometric position of the X-ray tube array module, and image fusion of the registered multi-view images to generate a single synthetic two-dimensional image; the system reconstructs a three-dimensional image based on the single synthetic two-dimensional image and the multi-view geometric information of the X-ray tube array module.

[0005] Preferably, the image processing module performs image registration on the images corresponding to the multiple frames of original two-dimensional digital image signals, specifically including: Establish a projection geometric model between the spatial coordinates of the focal point of each miniature X-ray tube in the X-ray tube array module and the imaging plane coordinates of the image acquisition module; Based on the projection geometry model, calculate the perspective transformation matrix between the multiple frames of original two-dimensional digital image signals corresponding to each miniature X-ray tube; The perspective transformation matrix is ​​used to transform the images corresponding to the multiple frames of original two-dimensional digital image signals to a unified reference coordinate system, thereby completing image registration.

[0006] Preferably, the image processing module performs image fusion on the registered multi-view images, specifically including: Each registered frame image is decomposed into multiple scales to obtain the sub-band coefficients at multiple scale levels corresponding to each frame image. For the same scale layer, the weight map of each frame image at the scale layer is calculated based on the local energy features of the subband coefficients corresponding to each frame image. Based on the weight map, the sub-band coefficients corresponding to all images at the same scale layer are weighted and fused to obtain the fused sub-band coefficients at the scale layer. A multi-scale inverse transform is performed on the fused subband coefficients of all scale layers to generate the single synthesized two-dimensional image.

[0007] Preferably, the system reconstructs a three-dimensional image based on the single synthesized two-dimensional image and the multi-view geometric information of the X-ray tube array module, specifically including: The contour information of the detected object is extracted from the registered multi-view images, and the spatial constraints of the detected object are constructed by combining the spatial ray equation between the focal point of each micro X-ray tube in the X-ray tube array module and the image acquisition module. Using the grayscale gradient information of pixels in the single synthesized two-dimensional image, a density projection relationship in three-dimensional voxel space is established; The system of equations consisting of the spatial constraints and the density projection relationship is solved iteratively to gradually update the density value of each voxel in the three-dimensional voxel space, thereby obtaining the initial three-dimensional volume data.

[0008] Preferably, iteratively solving the system of equations formed by the spatial constraints and the density projection relationship specifically includes: Using the single synthesized two-dimensional image as the projection target, the reconstructed three-dimensional volume data under the multi-view geometric information is compared with the projection target in the forward simulated projection, and the projection difference is calculated. Based on the projection difference and the spatial constraints, calculate the density correction for each voxel in the three-dimensional voxel space; The density value of each voxel in the three-dimensional voxel space is updated according to the density correction amount. The steps of forward simulation projection and calculation of correction amount are repeated until the projection difference is less than a preset threshold or the preset number of iterations is reached, and the initial three-dimensional volume data is obtained.

[0009] Preferably, based on the projection difference and the spatial constraints, the density correction amount for each voxel in the three-dimensional voxel space is calculated, specifically including: For the current voxel to be calculated in the three-dimensional voxel space, determine all spatial rays that pass through the current voxel and extend from the focal point of any miniature X-ray tube in the X-ray tube array module to the imaging plane of the image acquisition module; For each of the spatial rays, the projection difference corresponding to the spatial ray is obtained. The projection difference is the difference between the gray value of the corresponding pixel in the projection target and the forward simulated projection value of the current three-dimensional volume data along the direction of the spatial ray. Calculate the projection contribution weight of the current voxel to each space ray passing through it, the contribution weight being determined based on the shortest distance from the center of the current voxel to the space ray and the current density value of the current voxel; Multiply the projection difference corresponding to each space ray by the contribution weight of the current voxel to the space ray to obtain the preliminary correction component of the space ray to the current voxel. The weighted sum of all the preliminary correction components of the current voxel is then adjusted according to a preset relaxation factor to obtain the first intermediate correction amount. Based on the adjacent voxel smoothness condition in the spatial constraint, calculate the density gradient between the current voxel and its neighboring voxels, and calculate a spatial smoothness constraint factor based on the density gradient. Multiply the first intermediate correction by the spatial smoothing constraint factor to obtain the final density correction of the current voxel.

[0010] Preferably, the system further includes: The 3D image optimization module is used to process the initial 3D volume data to generate the final 3D image.

[0011] Preferably, the 3D image optimization module processes the initial 3D volume data, specifically including: For each voxel in the initial three-dimensional volume data, calculate the magnitude of the density gradient of its neighboring voxels; An anisotropic diffusion model is constructed based on the density gradient magnitude, and the diffusion coefficient in the anisotropic diffusion model is negatively correlated with the density gradient magnitude. The anisotropic diffusion model is applied to smooth the initial three-dimensional volume data, which suppresses noise while preserving the contour edges, resulting in optimized three-dimensional volume data.

[0012] Preferably, the 3D image optimization module is further used for: A density threshold is set, and voxels with density values ​​higher than the density threshold in the optimized three-dimensional volume data are segmented to form a three-dimensional surface point cloud of the detected object. The three-dimensional surface point cloud is processed into a triangular mesh to construct a three-dimensional surface model of the detected object, which is the final three-dimensional image.

[0013] Preferably, the present invention also includes a method for reconstructing three-dimensional images based on two-dimensional X-ray images, the method comprising all the modules and method flow of the system for reconstructing three-dimensional images based on two-dimensional X-ray images as described above.

[0014] Compared with the prior art, the beneficial effects of the present invention are: By employing a fixed geometric arrangement of miniature X-ray tube arrays for independent sequential exposure, the time delay and positional errors caused by mechanical movement are directly eliminated, achieving true static multi-angle data acquisition and preparation. The high-voltage control module supports independent and stable switching of power supply to each X-ray tube, enabling the array to perform sequential pulse exposures with extremely high temporal resolution. This avoids scattering interference caused by simultaneous exposure from multiple sources, ensuring that each frame of projected image has clear contrast and low noise levels. This design allows the system to instantly capture the projection information of the detected object from multiple fixed viewing angles, providing a high-quality data source with well-defined geometric relationships for subsequent processing.

[0015] Multiple original 2D images from different fixed viewpoints are first registered based on the known precise spatial position of the X-ray tube, and then fused to generate a single synthetic 2D image. This step essentially integrates and enhances multi-angle information within a 2D plane. The synthetic image combines projection features from various directions, improving the signal-to-noise ratio and highlighting the key structural information of the detected object. This preprocessing process differs from traditional workflows, generating an intermediate image with higher information density and more concentrated features.

[0016] Based on the aforementioned single synthesized 2D image, and combined with the inherent and known multi-view geometric constraints of the X-ray tube array, a new reconstruction path is constructed for 3D image reconstruction. This method transforms the complex 3D reconstruction problem, which traditionally directly processes massive multi-view projection data, into an optimization problem centered on the fused and enhanced image and constrained by fixed geometry. This simplifies the computational process, reduces the stringent dependence on strict consistency between the original multi-frame images, and thus reduces reconstruction artifacts caused by data errors at the algorithm level, enabling a more stable and efficient restoration of the 3D internal structure of the detected object. Attached Figure Description

[0017] Figure 1 This is a timing diagram of the three-dimensional image reconstruction system based on two-dimensional X-ray images described in this invention; Figure 2 Flowchart for image registration; Figure 3 A flowchart for multi-scale image fusion; Figure 4 This is a graph showing the relationship between the relaxation factor and the iteration convergence rate. Figure 5 This is a graph showing the relationship between point cloud density and triangular meshing quality. Detailed Implementation

[0018] The technical solution of the present invention will be clearly and completely described below with reference to the embodiments. 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.

[0019] See Figure 1The system for reconstructing 3D images based on 2D X-ray images includes an X-ray tube array module, a high-voltage control module, an image acquisition module, and an image processing module. The X-ray tube array module contains multiple miniature X-ray tubes arranged according to known geometric rules. The focal points of all miniature X-ray tubes jointly cover the target imaging area. Each miniature X-ray tube can independently control its exposure. The high-voltage control module provides an independent and stable high-voltage power supply for each miniature X-ray tube, supporting sequential pulsed exposure of the X-ray tube array module. The image acquisition module receives the X-ray signals after penetrating the object being detected and converts the X-ray signals into multiple frames of raw 2D digital image signals. The image processing module receives and processes the multiple frames of raw 2D digital image signals. When processing the multiple frames of raw 2D digital image signals, the image processing module performs image registration on the images corresponding to the multiple frames of raw 2D digital image signals according to the known geometric positions of the X-ray tube array module. The registered multi-view images are then fused to generate a single synthetic 2D image. The system reconstructs 3D images based on the single synthetic 2D image and the multi-view geometric information of the X-ray tube array module.

[0020] In one embodiment of the present invention, see [reference] Figure 2 In practical implementation, the image processing module establishes a projection geometric model between the spatial coordinates of the focal point of each miniature X-ray tube in the X-ray tube array module and the imaging plane coordinates of the image acquisition module. The projection geometric model is described using a pinhole camera model. The process of projecting a three-dimensional point in the world coordinate system onto the image plane is defined by the following formula:

[0021] Where: column vector A column vector representing the homogeneous coordinates of a 3D point in the world coordinate system. This represents the homogeneous pixel coordinates of the point projected onto the imaging plane of the image acquisition module. It is a non-zero scale factor, matrix It is the intrinsic parameter matrix of the image acquisition module. It is a 3x3 rotation matrix and 3x1 translation vector The combined extrinsic parameter matrix describes the rigid body transformation from the world coordinate system to the image acquisition module coordinate system. In some embodiments, the world coordinates of all the focal points of the miniature X-ray tubes in the X-ray tube array module... It is the known quantity pre-calibrated and stored through precision mechanical assembly drawings or optical tracking systems, the intrinsic parameter matrix of the image acquisition module. Similarly, it is obtained and stored through a pre-calibrated camera process. Optionally, for cameras numbered... The miniature X-ray tube, its corresponding extrinsic parameter matrix It is uniquely determined that it associates the world coordinates of the focal point of the miniature X-ray tube with the imaging plane of the image acquisition module.

[0022] In practice, the perspective transformation matrix between multiple frames of original two-dimensional digital image signals corresponding to each miniature X-ray tube is calculated based on the aforementioned projection geometry model. and miniature X-ray tube The image processing module needs to calculate the values ​​from the two acquired images. pixel plane to image homography matrix of pixel plane This is understandable; homography matrix It is a 3x3 matrix that satisfies the following relation for objects on the same plane in space: ,in: and These are the two-dimensional objects in the image. and images Homogeneous pixel coordinates on the image. The image processing module utilizes a known miniature X-ray tube. and miniature X-ray tube extrinsic matrix , and the public intrinsic parameter matrix Based on the assumption that the imaged object is mainly located in an approximately planar region, the homography matrix is ​​calculated through mathematical derivation. In some embodiments, the image processing module selects a specific micro-X-ray tube viewpoint as a unified reference coordinate system, for example, selecting the viewpoint of the micro-X-ray tube located at the geometric center of the tube array as the reference viewpoint. The image processing module calculates the perspective transformation matrix from all other micro-X-ray tube viewpoint images to this reference viewpoint image. ,in An index representing a non-referenced perspective.

[0023] In practice, a perspective transformation matrix is ​​used to transform the images corresponding to multiple frames of original two-dimensional digital image signals to a unified reference coordinate system. The image processing module applies the corresponding perspective transformation matrix to each image from a non-reference viewpoint. Perform a geometric transformation. Optionally, the geometric transformation process resamples the original image, and for each integer pixel position in the transformed image, it is performed using a perspective transformation matrix. The inverse transform finds the corresponding floating-point coordinates of the pixel in the original image. A bilinear interpolation algorithm is then used to obtain the grayscale value of that location from the original image and assign it to the pixel in the transformed image. It can be understood that, after the above geometric transformation and resampling operations, all images acquired by different miniature X-ray tubes are corrected to the same reference imaging plane coordinate system, eliminating the viewing angle differences and image geometric distortions caused by the different focal points of the miniature X-ray tubes, thus completing image registration.

[0024] In one embodiment of the present invention, see [reference] Figure 3 In specific implementation, the image processing module performs image fusion on the registered multi-view images and performs multi-scale decomposition on each registered frame to obtain sub-band coefficients at multiple scale levels corresponding to each frame. Multi-scale decomposition can be performed using the Laplacian pyramid decomposition method or the discrete wavelet transform method. In the Laplacian pyramid decomposition, each registered image is decomposed into a set of sub-band images with different spatial frequencies. In the discrete wavelet transform, each registered image is decomposed into a set of low-frequency approximation coefficients and high-frequency detail coefficients in multiple directions. In some embodiments, the image processing module uses discrete wavelet transform for multi-scale decomposition, for each registered input image... ,in Indicates the image sequence number. Representing pixel coordinates, after a first-order discrete wavelet transform, four sub-bands are obtained: low-frequency approximate sub-band. Horizontal high-frequency detail sub-band Vertical high-frequency detail sub-band and diagonal high-frequency detail subband For low-frequency approximate subbands Continue with the next level of discrete wavelet transform to construct the image. The multi-scale pyramid representation, where each level of decomposition produces a set of sub-band coefficients corresponding to different scale layers and directions.

[0025] In practical implementation, for the same scale layer, the weight map of each frame image at that scale layer is calculated based on the local energy features of the subband coefficients corresponding to each frame image. The local energy features reflect the activity level or information richness of the image subband coefficients within the local window. It can be understood that for the scale layer... and direction ,image In position Local energy at the point It can be calculated by the sum of squares of the sub-band coefficient amplitudes within a local window centered at that location, and the calculation method is defined by the formula:

[0026] in: Representing an image At the scale layer ,direction Above Sub-band coefficient of location, Indicates a A local neighborhood window centered on the region, such as a 3x3 or 5x5 rectangular area. These are relative coordinates within the window. Based on the calculated local energy. ,image At the scale layer ,direction ,Location Fusion weights at the point All can be accessed at this location The local energy of the image is obtained by normalization, for example, using the formula... The calculations are performed to generate a weight value for each position, and the set of weight values ​​for all positions constitutes the image. The weight map at this scale and in this direction. Optionally, to improve the stability of the fusion result, the calculated initial weight map can be subjected to guided filtering or Gaussian smoothing to obtain a spatially smoother final weight map.

[0027] In practice, the sub-band coefficients of all images at the same scale are weighted and fused according to the weight map. The image processing module reads all images at the same scale and in the same direction. Subband coefficients of an image and the corresponding weighted graph For each position in the coefficient matrix Subband coefficients after fusion Calculated by weighted summation: This operation can be understood as being performed independently at each scale level and in each direction, resulting in the fused sub-band coefficient matrix for each scale level and direction. In some embodiments, for the lowest frequency approximate subband, the fusion strategy may employ an averaging method or select the image coefficients with the highest sharpness, rather than a weighted fusion based on local energy.

[0028] In practice, a multi-scale inverse transform is performed on the fused subband coefficients across all scale layers to generate a single synthetic 2D image. The image processing module performs multi-scale reconstruction on the fused subband coefficients at each scale layer and in each direction, following the order from the coarsest scale to the finest scale. Optionally, when using discrete wavelet transform for decomposition, the image processing module performs an inverse discrete wavelet transform, progressively upsampling and synthesizing the fused low-frequency approximate subband coefficients and high-frequency detail subband coefficients at each level, ultimately reconstructing a single image in the spatial domain. This reconstruction process is the inverse of the aforementioned multi-scale decomposition; it integrates the fused subband coefficients representing different frequency information into a complete image, which is the single synthetic 2D image that incorporates multi-view information.

[0029] In practice, the Laplacian pyramid decomposition method is used to decompose each registered frame of the image into multiple scales during the image fusion process. The construction of the Laplacian pyramid begins with the Gaussian pyramid of the input image. For each registered image to be fused, the image processing module first applies a low-pass filter to the image for convolution processing to smooth the image details. Then, the smoothed image is downsampled at a fixed ratio, usually by reducing the number of rows and columns of the image by half. The image obtained after one downsampling becomes the next layer of the Gaussian pyramid. This new image is used as input to repeat the low-pass filtering and downsampling operations. After several iterations, a set of images with progressively lower resolution is obtained, which is the Gaussian pyramid of the input image. Next, a Laplacian pyramid is generated based on the Gaussian pyramid. The image processing module starts from the bottom layer of the Gaussian pyramid and upsamples the image of a certain layer to restore its size to the size of the previous layer. Upsampling usually uses interpolation algorithms such as bilinear interpolation to estimate the gray value of the newly added pixels. Then, the upsampled image is subtracted pixel by pixel from the image of the original resolution of the previous layer in the Gaussian pyramid. This difference image captures the high-frequency detail information lost between two adjacent scales. This difference image is the sub-band image of the Laplacian pyramid at that level. This upsampling and subtraction process is repeated for each layer of the Gaussian pyramid to obtain a set of sub-band images of different scales, i.e., the Laplacian pyramid. Each sub-band image represents the information of the original image in a specific spatial frequency band.

[0030] In the context of image fusion, each registered input image independently generates a Laplacian pyramid. Each layer of the pyramid corresponds to a scale layer, and the data of that layer is the sub-band coefficient of the image at that scale. These coefficients clearly separate the different frequency components of the image from coarse to fine, providing a frequency-division representation basis for subsequent fusion based on local energy features.

[0031] In practical implementation, the core basis for calculating the weight map of each frame at a specific scale level is the local energy features of the subband coefficients. Obtaining these local energy features first requires defining a sliding neighborhood window on the subband coefficient matrix. This window is typically a rectangular region centered on the current pixel, with an odd number of rows and columns, such as a 3x3 or 5x5 pixel block. For a given scale level and direction, the image processing module iterates through each pixel position in the matrix. For the current position, the algorithm extracts the values ​​of all subband coefficients within the neighborhood window and then calculates the sum of squares of these values. This sum of squares is defined as the local energy value of that pixel position at that scale level and direction.

[0032] The calculation of local energy values ​​reflects the activity level or information intensity of the signal in the neighborhood of a pixel. Regions with higher energy values ​​usually correspond to edges, textures, and other structurally rich parts of the image. After calculating the local energy of all pixels in a single subband of a single image, a local energy map with the same size as the coefficient matrix of that subband is obtained. To fuse multiple images, the image processing module needs to align the local energy maps of all input images at the same scale and in the same direction. For the same pixel position in the fusion weight map, the system reads the local energy values ​​of all input images at that position, sums these energy values ​​to obtain the total energy, and then divides the energy value of each image at that position by the total energy. The resulting ratio is the initial fusion weight of that image at that pixel position.

[0033] This calculation ensures that for each pixel, the sum of the weights of all input images at that point is one, and images with higher energy receive greater weights at that point. To obtain a more spatially coherent weight allocation that avoids noise interference, the image processing module also performs post-processing on the initial weight map. Post-processing typically uses an additional smoothing filter, such as a Gaussian filter, to convolve the initial weight map. The filtering operation can eliminate small abrupt changes in the weight map caused by noise, making the spatial transition of weights more natural and smooth. The smoothed weight map is ultimately used to guide the weighted fusion of subband coefficients.

[0034] In one embodiment of the present invention, the system reconstructs a three-dimensional image based on a single synthesized two-dimensional image and the multi-view geometric information of the X-ray tube array module. Contour information of the detected object is extracted from the registered multi-view images. Contour extraction is accomplished by applying an edge detection algorithm to each registered image. The edge detection algorithm identifies pixel positions with drastic grayscale changes in the image and connects these positions to form a two-dimensional contour line of the detected object from the corresponding viewpoint. In some embodiments, for each registered image, the image processing module performs a binarization operation, marking the area within the contour line as the foreground and the area outside the contour line as the background, thereby obtaining a binary mask image of the object from that viewpoint. The binary mask images from all viewpoints collectively describe the two-dimensional shape information of the detected object projected from different directions. The spatial constraints of the detected object are constructed by combining the spatial ray equation between the focal point of each micro-X-ray tube in the X-ray tube array module and the image acquisition module. The spatial ray equation is jointly determined by the three-dimensional coordinates of the focal point of the micro-X-ray tube and the three-dimensional coordinates of the corresponding contour pixels on the imaging plane of the image acquisition module.

[0035] In specific implementation, for the first in the X-ray tube array module A miniature X-ray tube, the coordinates of its focal point in the world coordinate system are: A contour pixel in the registration image corresponding to the X-ray tube. The back-projection ray can be calculated using imaging geometry, and any point on this ray... Satisfy the equation:

[0036] Where: vector From the focus Pointing to image pixels The direction vector corresponding to the direction in space. Spatial constraints mean that the detected object must lie within the common intersection region of all such back-projected rays originating from different foci and passing through its corresponding contour pixels. A density projection relationship in three-dimensional voxel space is established using the gray-level gradient information of pixels in a single synthesized two-dimensional image. The gray-level gradient of the single synthesized two-dimensional image reflects the cumulative effect of material density changes along the X-ray penetration path. It can be understood that the density projection relationship links the voxel density value along a projection path in three-dimensional voxel space with the gray-level value of the corresponding pixel in the synthesized two-dimensional image. For each pixel in the synthesized image... Its grayscale value This can be modeled as starting from a virtual projection source, traversing three-dimensional voxel space, and reaching a pixel. The linear weighted sum of all voxel density values ​​along the path.

[0037] The initial 3D volume data is obtained by iteratively solving the equations formed by spatial constraints and density projection relationships to gradually update the density value of each voxel in the 3D voxel space. The image processing module defines a 3D voxel grid containing the detected object, with each voxel having an initial density value to be solved. Optionally, the initial density value can be set to zero or a uniform background value. In some embodiments, the iterative algorithm adopts an algebraic reconstruction technique framework. In each iteration, the system calculates a simulated synthetic 2D image in the forward direction according to the density projection relationship based on the current 3D voxel density value. The simulated image is compared with a real single synthetic 2D image to obtain the error. Spatial constraints play a role by limiting the range of density updates. For example, for voxels that are determined to be impossible to exist according to spatial constraints, their density value is forced to be zero or a very small value during the iteration process. The image processing module calculates the density correction amount for each voxel based on the error and spatial constraints, and updates the voxel density value with the correction amount. This process is repeated until a preset convergence condition is met, such as the error being less than a threshold or the maximum number of iterations being reached. The resulting 3D voxel density distribution is the initial 3D volume data. It can be understood that the initial three-dimensional volume data is a three-dimensional matrix, and the value of each element in the matrix represents the material density or linear decay coefficient within the corresponding small volume element in three-dimensional space.

[0038] In one embodiment of the present invention, the iterative solution of the equations relating spatial constraints and density projection is achieved through the following steps: a single synthesized two-dimensional image is used as the projection target; the forward simulated projection of the reconstructed three-dimensional volume data under multi-view geometric information is compared with the projection target to calculate the projection difference. The forward simulated projection process is based on the known position of the focal point of each micro-X-ray tube in the X-ray tube array module and the imaging plane geometry of the image acquisition module. For the current three-dimensional volume data, the system calculates the cumulative voxel density values ​​along the ray path from each micro-X-ray tube focal point, through the three-dimensional voxel space, to each pixel on the imaging plane; this cumulative value is the simulated projection value. The projection difference is the difference between the gray value of the corresponding pixel in the projection target and the forward simulated projection value of the current three-dimensional volume data along the spatial ray direction. In some embodiments, for each pixel, the projection difference can be expressed as the difference between the gray value of the projection target and the current forward simulated projection value.

[0039] In practice, based on the projection difference and spatial constraints, the density correction for each voxel in the three-dimensional voxel space is calculated. For the current voxel to be calculated in the three-dimensional voxel space, all spatial rays passing through the current voxel and originating from the focal point of any miniature X-ray tube in the X-ray tube array module to the imaging plane of the image acquisition module are identified. For each spatial ray, the corresponding projection difference is obtained. The projection contribution weight of the current voxel to each spatial ray passing through it is calculated. The contribution weight is determined based on the shortest distance from the center of the current voxel to the spatial ray and the current density value of the current voxel.

[0040] Multiply the projection difference corresponding to each spatial ray by the contribution weight of the current voxel to the spatial ray to obtain the initial correction component of the spatial ray to the current voxel. Sum all the initial correction components of the current voxel with weights and adjust them according to a preset relaxation factor to obtain the first intermediate correction amount. Calculate the density gradient between the current voxel and its neighboring voxels based on the smoothness condition of adjacent voxels in the spatial constraints, and calculate a spatial smoothing constraint factor based on the density gradient. Multiply the first intermediate correction amount by the spatial smoothing constraint factor to obtain the final density correction amount of the current voxel. The calculation of the final density correction amount is defined by the formula:

[0041] in: It is a voxel The final density correction amount, It is a preset relaxation factor. Based on voxels The spatial smoothing constraint factor for calculating the density gradient magnitude. It is all that passes through voxels A collection of spatial rays It is a ray Additional weights, It is a ray The corresponding projection difference multiplied by the voxel The initial correction component is obtained after weighting the contribution of this ray.

[0042] In practice, the density value of each voxel in the 3D voxel space is updated according to the density correction amount; for example, the new density value equals the old density value plus the final density correction amount. The steps of forward simulation projection and correction calculation are repeated until the projection difference is less than a preset threshold or a preset number of iterations is reached, thus obtaining the initial 3D volume data. Optionally, during the iteration process, the threshold for the projection difference can be set as a certain proportion of the average absolute difference of all pixels. This can be understood as a spatial smoothing constraint factor. The calculation depends on the density gradient magnitude of the current voxel, which can be obtained by calculating the square root of the sum of the squares of the density differences between the current voxel and its six neighboring voxels. In some embodiments, the calculation parameters for the contribution weight can be adjusted according to the characteristics of the imaging system; see Table 1, which shows one such contribution weight calculation parameter.

[0043]

[0044] See Figure 4 This is a graph showing the relationship between the relaxation factor and the iteration convergence speed, a professional analysis chart belonging to the iterative optimization stage of the X-ray 3D reconstruction algorithm. The curve exhibits a single-peak shape, with the convergence speed reaching its peak (close to 98 iterations / second) when the relaxation factor is approximately 0.6, indicating that this is the optimal relaxation factor value for the algorithm. In practical X-ray 3D reconstruction systems, setting the relaxation factor to around 0.6 can maximize the convergence speed while ensuring iterative stability, thereby reducing the reconstruction time of the 3D image. The convergence speed directly determines the reconstruction time of the 3D image. In scenarios with high real-time requirements, such as medical and industrial inspection, using the optimal relaxation factor can significantly improve the system's processing efficiency and throughput.

[0045] In one embodiment of the present invention, in a specific implementation, the three-dimensional image optimization module processes the initial three-dimensional volume data to generate a final three-dimensional image. For each voxel in the initial three-dimensional volume data, the density gradient magnitude of its neighboring voxels is calculated. The density gradient magnitude reflects the rate of density change in a local region of the three-dimensional volume data. For voxels located in three-dimensional grid coordinates... The density gradient magnitude of the voxel at that location The calculation method is defined by the formula:

[0046] in: This indicates the initial 3D volume data at the grid position. The density value at that location is used. An anisotropic diffusion model is constructed based on the density gradient magnitude. The diffusion coefficient in the anisotropic diffusion model is negatively correlated with the density gradient magnitude, and the diffusion coefficient is a function that decreases as the gradient magnitude increases. This design allows for weaker diffusion in regions with large density gradient magnitudes, such as the edges of object contours, thus preserving details, while stronger diffusion in regions with small density gradient magnitudes, such as within homogeneous tissue, thus smoothing noise.

[0047] In practice, an anisotropic diffusion model is applied to smooth the initial 3D volume data. The anisotropic diffusion process is achieved by iteratively updating the density value of each voxel. In each iteration, the new density value of a voxel is determined by a weighted average of its current density value and the density values ​​of its neighboring voxels, with the weights controlled by the diffusion coefficient. The 3D image optimization module performs multiple iterations until a preset number of iterations is reached or the density change is less than a set threshold, thereby suppressing noise while preserving the contour edges and obtaining optimized 3D volume data. Optionally, the diffusion coefficient function can be a mathematical expression negatively correlated with the magnitude of the density gradient. The 3D image optimization module also sets a density threshold to segment voxels in the optimized 3D volume data with density values ​​higher than the density threshold, forming a 3D surface point cloud of the detected object.

[0048] In practical implementation, the density threshold can be determined based on the grayscale histogram analysis of the optimized 3D volume data, for example, by automatically calculating a segmentation threshold using the Otsu's method. The 3D image optimization module traverses the optimized 3D volume data, extracting the 3D coordinates of the voxel center points with density values ​​greater than the threshold, forming a discrete set of points, i.e., a 3D surface point cloud. The 3D surface point cloud is then triangulated to construct a 3D surface model of the detected object. In some embodiments, triangulation can employ the moving cube algorithm or the Delaunay triangulation algorithm. The moving cube algorithm generates isosurface triangular patches within the voxel by examining the relationship between the density values ​​of the eight corner points of the voxel and the threshold; the Delaunay triangulation algorithm directly triangulates the discrete 3D point cloud, forming a closed triangular mesh surface. It can be understood that through triangulation, the discrete 3D surface point cloud is connected into a continuous surface model composed of numerous small triangular patches, which is the final 3D image. Optionally, the 3D image optimization module can also perform smoothing or simplification post-processing on the generated triangular mesh.

[0049] See Figure 5This is a graph showing the relationship between point cloud density and triangulation quality. It's a professional analysis chart used in the X-ray 3D image optimization stage to evaluate the impact of point cloud density on the final 3D surface model. As point cloud density increases, mesh accuracy rises rapidly, but the increase slows significantly after reaching approximately 300 points / cm², eventually stabilizing at around 0.95. This indicates that once the point cloud density exceeds a certain threshold, further increasing sampling points has very limited marginal benefits for improving model accuracy. The overall trend is upward, but with significant fluctuations. The generation time peaks (approximately 5.5 seconds) in the point cloud density range of 400-450 points / cm², then slightly decreases. This shows that the computational complexity of the algorithm increases with point cloud density, and computational efficiency may fluctuate in certain density ranges due to data distribution characteristics. In practical applications, a point cloud density of approximately 300 points / cm² is recommended. This represents a balance between accuracy and computational efficiency, ensuring mesh accuracy above 0.93 while keeping the generation time within 4 seconds.

[0050] The above embodiments are only used to illustrate the technical methods of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical methods of the present invention without departing from the spirit and scope of the technical methods of the present invention.

Claims

1. A system for reconstructing three-dimensional images based on two-dimensional X-ray images, characterized in that, The system includes: The X-ray tube array module contains multiple miniature X-ray tubes arranged in a known geometric pattern. The focal points of all the miniature X-ray tubes together cover the target imaging area, and each miniature X-ray tube can independently control the exposure. A high-voltage control module is used to provide an independent and stable high-voltage power supply for each miniature X-ray tube, supporting the sequential pulsed exposure of the tube array module; The image acquisition module is used to receive X-ray signals after penetrating the object being detected, and to convert the X-ray signals into multiple frames of raw two-dimensional digital image signals. An image processing module is used to receive the multiple frames of original two-dimensional digital image signals and process the multiple frames of original two-dimensional digital image signals. The image processing module processes the multi-frame original two-dimensional digital image signals, specifically including: image registration of the images corresponding to the multi-frame original two-dimensional digital image signals according to the known geometric position of the X-ray tube array module, and image fusion of the registered multi-view images to generate a single synthetic two-dimensional image; the system reconstructs a three-dimensional image based on the single synthetic two-dimensional image and the multi-view geometric information of the X-ray tube array module.

2. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 1, characterized in that, The image processing module performs image registration on the images corresponding to the multiple frames of original two-dimensional digital image signals, specifically including: Establish a projection geometric model between the spatial coordinates of the focal point of each miniature X-ray tube in the X-ray tube array module and the imaging plane coordinates of the image acquisition module; Based on the projection geometry model, calculate the perspective transformation matrix between the multiple frames of original two-dimensional digital image signals corresponding to each miniature X-ray tube; The perspective transformation matrix is ​​used to transform the images corresponding to the multiple frames of original two-dimensional digital image signals to a unified reference coordinate system, thereby completing image registration.

3. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 2, characterized in that, The image processing module performs image fusion on the registered multi-view images, specifically including: Each registered frame image is decomposed into multiple scales to obtain the sub-band coefficients at multiple scale levels corresponding to each frame image. For the same scale layer, the weight map of each frame image at the scale layer is calculated based on the local energy features of the subband coefficients corresponding to each frame image. Based on the weight map, the sub-band coefficients corresponding to all images at the same scale layer are weighted and fused to obtain the fused sub-band coefficients at the scale layer. A multi-scale inverse transform is performed on the fused subband coefficients of all scale layers to generate the single synthesized two-dimensional image.

4. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 3, characterized in that, The system reconstructs a three-dimensional image based on the single synthesized two-dimensional image and the multi-view geometric information of the X-ray tube array module, specifically including: The contour information of the detected object is extracted from the registered multi-view images, and the spatial constraints of the detected object are constructed by combining the spatial ray equation between the focal point of each micro X-ray tube in the X-ray tube array module and the image acquisition module. Using the grayscale gradient information of pixels in the single synthesized two-dimensional image, a density projection relationship in three-dimensional voxel space is established; The system of equations consisting of the spatial constraints and the density projection relationship is solved iteratively to gradually update the density value of each voxel in the three-dimensional voxel space, thereby obtaining the initial three-dimensional volume data.

5. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 4, characterized in that, Iteratively solving the system of equations formed by the spatial constraints and the density projection relationship specifically includes: Using the single synthesized two-dimensional image as the projection target, the reconstructed three-dimensional volume data under the multi-view geometric information is compared with the projection target in the forward simulated projection, and the projection difference is calculated. Based on the projection difference and the spatial constraints, calculate the density correction for each voxel in the three-dimensional voxel space; The density value of each voxel in the three-dimensional voxel space is updated according to the density correction amount. The steps of forward simulation projection and calculation of correction amount are repeated until the projection difference is less than a preset threshold or the preset number of iterations is reached, and the initial three-dimensional volume data is obtained.

6. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 5, characterized in that, Based on the projection difference and the spatial constraints, the density correction for each voxel in the three-dimensional voxel space is calculated, specifically including: For the current voxel to be calculated in the three-dimensional voxel space, determine all spatial rays that pass through the current voxel and extend from the focal point of any miniature X-ray tube in the X-ray tube array module to the imaging plane of the image acquisition module; For each of the spatial rays, the projection difference corresponding to the spatial ray is obtained. The projection difference is the difference between the gray value of the corresponding pixel in the projection target and the forward simulated projection value of the current three-dimensional volume data along the direction of the spatial ray. Calculate the projection contribution weight of the current voxel to each space ray passing through it, the contribution weight being determined based on the shortest distance from the center of the current voxel to the space ray and the current density value of the current voxel; Multiply the projection difference corresponding to each space ray by the contribution weight of the current voxel to the space ray to obtain the preliminary correction component of the space ray to the current voxel. The weighted sum of all the preliminary correction components of the current voxel is then adjusted according to a preset relaxation factor to obtain the first intermediate correction amount. Based on the adjacent voxel smoothness condition in the spatial constraint, calculate the density gradient between the current voxel and its neighboring voxels, and calculate a spatial smoothness constraint factor based on the density gradient. Multiply the first intermediate correction by the spatial smoothing constraint factor to obtain the final density correction of the current voxel.

7. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 6, characterized in that, The system also includes: The 3D image optimization module is used to process the initial 3D volume data to generate the final 3D image.

8. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 7, characterized in that, The 3D image optimization module processes the initial 3D volume data, specifically including: For each voxel in the initial three-dimensional volume data, calculate the magnitude of the density gradient of its neighboring voxels; An anisotropic diffusion model is constructed based on the density gradient magnitude, and the diffusion coefficient in the anisotropic diffusion model is negatively correlated with the density gradient magnitude. The anisotropic diffusion model is applied to smooth the initial three-dimensional volume data, which suppresses noise while preserving the contour edges, resulting in optimized three-dimensional volume data.

9. The system for reconstructing three-dimensional images based on two-dimensional X-ray images according to claim 8, characterized in that, The 3D image optimization module is also used for: A density threshold is set, and voxels with density values ​​higher than the density threshold in the optimized three-dimensional volume data are segmented to form a three-dimensional surface point cloud of the detected object. The three-dimensional surface point cloud is processed into a triangular mesh to construct a three-dimensional surface model of the detected object, which is the final three-dimensional image.

10. A method for reconstructing three-dimensional images based on two-dimensional X-ray images, characterized in that, The system includes all modules and method flows of the system for reconstructing three-dimensional images based on two-dimensional X-ray images as described in any one of claims 1 to 9.

Images