Dual-modality tomographic image registration method and system based on multi-directional structural features

Abstract

The application discloses a kind of based on multi-direction structural feature dual-mode tomographic image registration method and system, the method includes: according to each direction angle corresponding structure response diagram selects effective direction angle, and according to each effective direction angle in first modality image determines closed path;With closed path mapping to second modality image obtains initial mapping closed path and carries out position correction, obtains target mapping closed path;Along closed path, extract first response value sequence, along target mapping closed path, extract second response value sequence;With second response value sequence after transformation, calculate and first response value sequence's degree of coincidence, with all effective direction angle degree of coincidence's weighted sum as global registration degree of coincidence;According to global registration degree of coincidence adjustment space transformation parameter, output registration result.Dual-mode image gray response opposite physical characteristics are used to establish cross-modal structure corresponding relationship, can realize the stable automatic registration of dual-mode tomographic image.

Description

Technical Field

[0001] This invention belongs to the field of multimodal image registration technology, and particularly relates to a dual-modal tomographic image registration method and system based on multi-directional structural features. Background Technology

[0002] Existing bimodal image registration methods mainly fall into two categories: registration methods based on local features and registration methods based on gray-level statistics. The former typically uses local features such as Scale Invariant Feature Transform (SIFT), Speed-Up Robust Features (SURF), and Oriented Fast Rotation Binary Features (ORB) to establish cross-modal correspondences; the latter typically uses similarity measures such as Mutual Information (MI), Normalized Cross-Correlation (NCC), Enhanced Correlation Coefficient (ECC), Modal Independent Neighborhood Descriptor (MIND), and Normalized Gradient Field (NGF), and combines them with affine models, rigid body models, or multi-resolution pyramid strategies for optimization.

[0003] For general multimodal images, the methods described above can address translation, rotation, scale changes, and some structural mismatches to a certain extent. However, due to significant differences in imaging mechanisms, neutron / X-ray dual-modal tomographic images often exhibit complex nonlinear relationships in the grayscale response of the same structure in both modes, with local areas potentially showing contrast reversal. Furthermore, neutron and X-ray images differ significantly in edge width, texture representation, and noise characteristics. Therefore, traditional feature-point-based methods are susceptible to weak features, repetitive textures, and mismatches, while traditional methods based on full-image grayscale statistics are prone to getting trapped in local optima, making it difficult to simultaneously achieve global geometric consistency and accurate alignment of internal structures.

[0004] Especially for samples such as cylindrical lithium-ion batteries, the metal shell boundary in X-ray images has a strong response, while neutron images are more sensitive to the internal electrolyte and rolled layered structure. This makes existing full-image grayscale registration methods prone to being dominated by the outer contour, thus weakening the registration accuracy of the internal structural regions. Traditional local feature methods are easily affected by near-cylindrical symmetry, repetitive rolled textures, and insufficient local recognition, resulting in mismatches, reliance on manual point selection, or insufficient automation.

[0005] Currently, most publicly available studies on neutron / X-ray dual-mode tomography focus on joint imaging and subsequent structural characterization. There is relatively little systematic research specifically on the registration method for neutron / X-ray dual-mode tomography images, and a relatively complete and stable methodological system has not yet been formed. In particular, there is a lack of a registration method that can both enhance the comparability of internal structures across modes and take into account the correspondence between slices of the entire tomographic sequence and the transfer of continuous parameters. Summary of the Invention

[0006] The technical problem to be solved by this invention is: to address the problems of large differences in grayscale response, obvious local contrast reversal, and repetitive and directional internal structures between neutron / X-ray dual-mode tomographic images, a registration method is proposed that can both establish a stable structural correspondence by utilizing the complementary response characteristics of dual modes and achieve automatic global affine registration.

[0007] In a first aspect, the present invention provides a dual-modal tomographic image registration method based on multi-directional structural features, comprising: Acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; Anisotropic Gaussian derivative kernels are constructed using multiple preset orientation angles, and the first modal image and the second modal image are filtered according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image. The multi-directional structural response tensor contains the structural response map corresponding to each orientation angle. Effective orientation angles are selected from the preset multiple orientation angles based on the structural response diagrams corresponding to each orientation angle, and closed paths are determined in the first modal image based on each effective orientation angle, wherein each effective orientation angle corresponds to a closed path; The closed path is mapped to the second modal image after applying the current spatial transformation parameters to obtain the initial mapped closed path. The position of the initial mapped closed path is then corrected according to the second multi-directional structural response tensor to obtain the target mapped closed path. Extract a first response value sequence in the first modal image along the closed path, and extract a second response value sequence in the second modal image along the target mapping closed path; The second response value sequence is transformed to obtain the inverted response value sequence. The sequence matching degree between the first response value sequence and the inverted response value sequence is calculated, and the weighted sum of the sequence matching degrees corresponding to all valid direction angles is used as the global registration matching degree. The current spatial transformation parameters are adjusted according to the global registration fit, and the adjusted current spatial transformation parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

[0008] Secondly, the present invention provides a dual-modal tomographic image registration system based on multi-directional structural features, comprising: The acquisition module is configured to acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; The filtering module is configured to construct anisotropic Gaussian derivative kernels at multiple preset direction angles, and filter the first modal image and the second modal image according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image, wherein the multi-directional structural response tensor contains the structural response map corresponding to each direction angle; The determination module is configured to select effective direction angles from the preset multiple direction angles based on the structural response map corresponding to each direction angle, and determine a closed path in the first modal image based on each effective direction angle, wherein each effective direction angle corresponds to a closed path; The mapping module is configured to apply the current spatial transformation parameters to the closed path and map it to the second modal image to obtain an initial mapped closed path, and to perform position correction on the initial mapped closed path according to the second multi-directional structural response tensor to obtain a target mapped closed path; The extraction module is configured to extract a first response value sequence in the first modal image along the closed path, and to extract a second response value sequence in the second modal image along the target mapping closed path; The calculation module is configured to perform a sequence transformation on the second response value sequence to obtain an inverted response value sequence, calculate the sequence matching degree between the first response value sequence and the inverted response value sequence, and use the weighted sum of the sequence matching degrees corresponding to all effective direction angles as the global registration matching degree. The output module is configured to adjust the current spatial transformation parameters according to the global registration fit, apply the adjusted current spatial transformation parameters to the second modal image, align the second modal image with the first modal image, and output the registration result.

[0009] Thirdly, an electronic device is provided, comprising: at least one processor, and a memory communicatively connected to the at least one processor, wherein the memory stores instructions executable by the at least one processor, the instructions being executed by the at least one processor to enable the at least one processor to perform the steps of the dual-modal tomographic image registration method based on multi-directional structural features according to any embodiment of the present invention.

[0010] Fourthly, the present invention also provides a computer-readable storage medium having a computer program stored thereon, wherein when the program instructions are executed by a processor, the processor performs the steps of the dual-modal tomographic image registration method based on multi-directional structural features according to any embodiment of the present invention.

[0011] The dual-modal tomographic image registration method and system based on multi-directional structural features proposed in this application have the following beneficial effects: By utilizing the physical characteristic that dual-modal images have opposite gray-level responses to the same physical structure, the registration basis is changed from finding positions with similar responses in two images to finding positions with opposite responses in two images, fundamentally solving the problem of traditional gray-level statistical methods being affected by nonlinear gray-level relationships; By extracting the structural orientation response of the image through a multi-directional anisotropic Gaussian derivative kernel, the image is mapped from the grayscale domain to the structural response domain, so that the registration basis is expanded from absolute grayscale values ​​to edge, layered texture and directional structural information, effectively highlighting the internal structural information with cross-modal consistency. Effective orientation angles are selected based on cross-modal response consistency metrics, so that registration depends only on the orientation channel with the strongest structural response, avoiding interference from weak structural orientations and noise responses on the optimization process; By combining the area preservation of the region enclosed by the closed path with the statistical analysis of response polarity, the mapping path determination method can accurately find the structural boundary position corresponding to the closed path of the first mode in the second modal image, thus achieving a stable correspondence between cross-modal structural boundaries.

[0012] By inverting the response polarity and calculating the sequence fit, the dual-modal complementary response characteristics are quantified into a numerical registration quality metric, providing a clear direction for optimizing spatial transformation parameters. Attached Figure Description

[0013] To more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings used in the following description of the embodiments will be briefly introduced. Obviously, the drawings described below are some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0014] Figure 1 A flowchart illustrating a dual-modal tomographic image registration method based on multi-directional structural features, provided in an embodiment of the present invention; Figure 2 A structural block diagram of a dual-modal tomographic image registration system based on multi-directional structural features is provided in an embodiment of the present invention. Figure 3 This is a schematic diagram of the structure of an electronic device provided in an embodiment of the present invention. Detailed Implementation

[0015] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, 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, 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] Please see Figure 1 The diagram shows a flowchart of a dual-modal tomographic image registration method based on multi-directional structural features according to this application.

[0017] like Figure 1 As shown, the dual-modal tomographic image registration method based on multi-directional structural features specifically includes the following steps: Step S101: Obtain a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure.

[0018] In this step, step S1011, dual-mode tomographic images of the same sample are acquired. The sample is scanned using a neutron imaging system to obtain a neutron tomographic image of the sample, which serves as the first modal image. The scanning parameters of the neutron imaging system are set according to the sample type and detection requirements, including the neutron source energy, detector exposure time, projection angle, etc. Neutrons have high sensitivity to light elements (such as hydrogen, lithium, carbon, etc.) and can clearly display the structure of light elements in the sample.

[0019] The same sample is scanned using an X-ray imaging system to obtain an X-ray tomographic image of the sample, which serves as the second modal image. The scanning parameters of the X-ray imaging system are set according to the sample type and detection requirements, including the X-ray source voltage, X-ray source current, detector exposure time, projection angle, etc. X-rays have high sensitivity to heavy elements (such as metal elements) and can clearly display the structure of heavy elements in the sample.

[0020] When neutron imaging and X-ray imaging systems scan the same sample, the scanning geometry parameters of the two systems (such as the position of the rotation center, the size of the detector pixel, the distance from the X-ray source to the detector, etc.) are set independently according to the imaging requirements of their respective systems, and do not need to be completely consistent.

[0021] Step S1012: Verify whether the grayscale responses of the two modal images to the same physical structure are opposite. For the same physical structure boundary in the first modal image and the second modal image, the gray value change direction of the first modal image and the gray value change direction of the second modal image are extracted along the normal direction of the boundary, respectively.

[0022] The grayscale value change direction is determined as follows: at the boundary of the physical structure, along the direction from the outside of the physical structure to the inside, the grayscale value sequence of the first modality image is obtained. If the grayscale value sequence increases with the position, the grayscale transition direction of the first modality image at the boundary of the physical structure is recorded as positive; if the grayscale value sequence decreases with the position, the grayscale transition direction is recorded as negative. The grayscale transition direction of the second modality image at the boundary of the physical structure is determined in the same way.

[0023] When the gray-level transition direction of the first modal image at the physical structure boundary is opposite to that of the second modal image at the physical structure boundary, it is determined that the gray-level response of the first modal image and the second modal image to the physical structure is opposite, and the current image pair is confirmed to meet the subsequent processing conditions.

[0024] Step S1013: Standardize the grayscale response direction. If the grayscale responses of the first modal image and the second modal image to the same physical structure are not completely opposite, but rather there are cases where the response directions of some structures are the same and the response directions of some structures are opposite, then the following unification process is performed: In the first modal image and the second modal image, select multiple representative physical structure boundaries respectively, and count the proportion of the gray-level transition directions of the first modal image at each boundary that are opposite to those of the gray-level transition directions of the second modal image at each boundary. If the proportion is greater than a preset proportion threshold, such as 80%, then the gray-level values ​​of the second modal image are reversed as a whole, so that the gray-level responses of the two images at all physical structure boundaries satisfy the opposite relationship.

[0025] The method for inverting the overall grayscale value is as follows: Subtract the sum of the maximum and minimum grayscale values ​​of the image from the grayscale value of each pixel in the second modality image, so that the range of grayscale values ​​after inversion is the same as the range of grayscale values ​​before inversion. That is, for the grayscale value I(x,y) of each pixel in the image, the inverted grayscale value I'(x,y) = I max +I min -I(x,y), where I max I is the maximum grayscale value of the image. min This is the minimum grayscale value of the image.

[0026] If the ratio is not greater than the preset ratio threshold, recheck the calibration status of the imaging system or adjust the position of the sample in the imaging system, reacquire the image, and verify again.

[0027] Step S1014: Establish an image spatial coordinate system and record the spatial correspondence. After verifying and unifying the inverse relationship of grayscale response, spatial coordinate systems are established for the first and second modal images respectively. The spatial coordinate system takes the upper left corner pixel position of the image as the origin, the horizontal direction as the first coordinate axis direction, the vertical direction as the second coordinate axis direction, and the vertical direction as the third coordinate axis direction (i.e., the axial direction of the tomographic sequence).

[0028] Record the initial spatial correspondence between the first modal image and the second modal image. The initial spatial correspondence includes the position offset of the two images in their respective spatial coordinate systems, the pixel size ratio, and the angle deviation. The initial spatial correspondence is used as the initial value of the current spatial transformation parameters.

[0029] Step S102: Construct anisotropic Gaussian derivative kernels for multiple preset direction angles, and filter the first modal image and the second modal image according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image. The multi-directional structural response tensor contains the structural response map corresponding to each direction angle.

[0030] In this step, the first direction angle of each pixel position in the first modal image along the direction of the maximum change in gray value in the neighborhood is obtained, and the first direction angles of all pixel positions are statistically analyzed according to the angle value to obtain a first direction angle histogram. The second direction angle of each pixel position in the second modal image along the direction of the maximum change in gray value in the neighborhood is obtained, and the second direction angles of all pixel positions are statistically analyzed according to the angle value to obtain a second direction angle histogram. Each local peak angle in the first orientation angle histogram is marked as a first candidate angle, and it is determined whether the distribution frequency at the same angle as a certain first candidate angle in the second orientation angle histogram is greater than a preset frequency threshold. If the distribution frequency is greater than the preset frequency threshold, then a certain first candidate angle is retained as a second candidate angle, and all first candidate angles are traversed to obtain a set of candidate direction angles composed of all second candidate angles. Perform a clustering operation on each second candidate angle in the candidate orientation angle set, and take the angle corresponding to each cluster center as a preset orientation angle to obtain a preset orientation angle set; Using each preset direction angle in the preset direction angle set as a rotation angle, rotate the two-dimensional Gaussian kernel around an axis perpendicular to the plane where the two-dimensional Gaussian kernel is located to the direction corresponding to the preset direction angle. Calculate the first-order partial derivative of the rotated two-dimensional Gaussian kernel along the direction corresponding to the preset direction angle to obtain the anisotropic Gaussian derivative kernel corresponding to each preset direction angle. The anisotropic Gaussian derivative kernels corresponding to each preset direction angle in the preset direction angle set are convolved with the first modal image to obtain the first structural response map under each preset direction angle. All the first structural response maps under the preset direction angles are arranged in the order of the preset direction angles to form the first multi-directional structural response tensor. The anisotropic Gaussian derivative kernels corresponding to each preset orientation angle in the preset orientation angle set are convolved with the second modality image to obtain the second structural response map under each preset orientation angle. All the second structural response maps under the preset orientation angles are arranged in the order of the preset orientation angles to form the second multi-directional structural response tensor.

[0031] In one specific embodiment, the optimal set of orientation angles is determined by utilizing the consistency of the orientation distribution of the physical structure in the dual-modal image. Then, an anisotropic Gaussian derivative kernel is constructed using this set of orientation angles, and the dual-modal images are filtered to obtain a multi-directional structural response tensor, specifically: Step S1021: Obtain the orientation angle histogram of the first modality image. Obtain the first direction angle of each pixel position in the first modality image along the direction of the maximum change in gray value in the neighborhood. The neighborhood is a window of a preset size centered on the pixel position, such as a 5×5 pixel window. The direction of the maximum change in gray value is determined by calculating the image gradient vector at the pixel position in the neighborhood. The direction of the gradient vector is the direction of the maximum change in gray value. The gradient vector can be calculated using the Sobel operator or the structure tensor method.

[0032] The first orientation angle of all pixel positions is statistically analyzed according to the angle value to obtain a first orientation angle histogram. The horizontal axis of the histogram is the angle value, ranging from 0° to 180°; the vertical axis is the number of pixels falling into each angle interval. The width of the angle interval division of the histogram can be set according to actual needs, for example, it can be divided into 180 intervals, with each interval width being 1°.

[0033] In the same way, the second direction angle of each pixel position in the second modality image along the direction of the largest change in gray value in the neighborhood is obtained, and the second direction angles of all pixel positions are statistically analyzed according to the angle value to obtain the second direction angle histogram.

[0034] Step S1022: Determine the candidate orientation angle set. Each local peak angle in the first direction angle histogram is marked as a first candidate angle. The method for determining the local peak angle is as follows: In the first direction angle histogram, for a certain angle position, if the frequency of the angle position is greater than the frequency of the adjacent angle positions on both sides, and the frequency of the angle position exceeds a preset multiple of the average frequency of all angle positions in the first direction angle histogram, such as 2 times, then the angle position is determined as a local peak angle.

[0035] For each first candidate angle, query the distribution frequency of the same angle as the first candidate angle in the second direction angle histogram, and determine whether the distribution frequency is greater than a preset frequency threshold. The preset frequency threshold is a preset multiple of the average frequency of all angle positions in the second direction angle histogram, such as 1.5 times.

[0036] If the distribution frequency is greater than the preset frequency threshold, the first candidate angle is retained as the second candidate angle. This judgment step utilizes the complementary characteristics of bimodal images: the same physical structure will produce directional responses in both modes. Therefore, if a certain angle is the peak of the directional distribution in the first modal image, and the distribution frequency of that angle is also high in the second modal image, then that angle corresponds to the physical structure direction common to both modes and has registration value.

[0037] Iterate through all first candidate angles to obtain a set of candidate direction angles consisting of all second candidate angles.

[0038] Step S1023: Cluster the candidate orientation angles to obtain a preset orientation angle set. Perform a clustering operation on each second candidate angle in the candidate orientation angle set. The clustering operation can be either mean-shift clustering or K-means clustering.

[0039] The specific process of clustering is as follows: take each angle value in the candidate orientation angle set as a data point, perform mean drift clustering on the angle axis, and cluster multiple candidate angles with a distance less than the preset clustering radius (e.g., 10°) into one class, and take the class center angle of each class as the cluster center.

[0040] The angle corresponding to each cluster center after clustering is taken as a preset direction angle, and a preset direction angle set is obtained. The number of direction angles in the preset direction angle set is the number of direction channels for subsequent filtering. This number can be set to 8 to 12 according to actual needs.

[0041] Step S1024: Construct an anisotropic Gaussian derivative kernel Using each preset direction angle in the preset direction angle set as a rotation angle, construct the anisotropic Gaussian derivative kernel corresponding to that preset direction angle.

[0042] First, define a standard two-dimensional Gaussian kernel function G(x,y)=(1 / (2πσ)). x σ y ))·exp(-(x² / (2σ x ²)+y² / (2σ y ²))), where σ x σ is the standard deviation in the horizontal direction. yThe standard deviation is defined as follows: the horizontal and vertical directions are the two coordinate axes of the plane containing the 2D Gaussian kernel. The kernel is rotated about an axis perpendicular to the plane of the kernel, with a preset rotation angle, so that the principal axis of the kernel aligns with the preset rotation angle. The standard deviation of the rotated kernel along the preset rotation angle is the first standard deviation σ. parallel The standard deviation in the direction perpendicular to the preset direction angle is the second standard deviation σ. perpendicular The first standard deviation is greater than the second standard deviation, for example, σ. parallel =2.0, σ perpendicular =0.8, which makes the anisotropic Gaussian kernel stretched along the preset direction angle and compressed in the direction perpendicular to the preset direction angle.

[0043] The first-order partial derivative of the rotated and stretched two-dimensional Gaussian kernel along the preset direction angle is obtained to obtain the anisotropic Gaussian derivative kernel corresponding to the preset direction angle. The first-order partial derivative is the derivative along the preset direction angle. Its function is to enable the kernel to produce a strong response to the structural edge extending along the direction, while the response to the structure not extending along the direction is weak. For the preset direction angle θ, the anisotropic Gaussian derivative kernel is calculated by taking the partial derivative of the rotated anisotropic Gaussian kernel along the θ direction.

[0044] Step S1025: Filter the first modality image to obtain the first multi-directional structural response tensor. The anisotropic Gaussian derivative kernel corresponding to each preset orientation angle in the preset orientation angle set is convolved with the first modality image to obtain the first structural response map under each preset orientation angle. The convolution operation is performed pixel by pixel in the image spatial domain. The convolution result at that position is obtained by multiplying the pixel value in the area covered by the convolution kernel with the coefficient at the corresponding position of the convolution kernel and summing the results.

[0045] For the i-th preset direction angle θ i The corresponding first structural response diagram R 1,i (x,y)=(Kθ i * I1)(x,y), where Kθ i For the i-th preset direction angle θ i The corresponding anisotropic Gaussian derivative kernel, I1 is the first modality image, * is the convolution operator, and (x,y) are the spatial coordinates in the plane where the two-dimensional Gaussian kernel is located.

[0046] Arrange the first structural response maps under all preset orientation angles in the order of the preset orientation angles (e.g., from smallest to largest) to form the first multi-directional structural response tensor. The first multi-directional structural response tensor is a three-dimensional data volume with dimensions of height × width × number of directional channels, where the first two dimensions are spatial dimensions and the third dimension is the directional dimension.

[0047] Step S1026: Filter the second modality image to obtain the second multi-directional structural response tensor. The anisotropic Gaussian derivative kernel corresponding to each preset orientation angle in the preset orientation angle set is convolved with the second modality image to obtain the second structural response map under each preset orientation angle.

[0048] For the i-th preset direction angle θ i The corresponding second structural response diagram R 2,i (x,y)=(Kθ i *I2)(x,y), where I2 is the second modality image.

[0049] Arrange all the second structural response maps under the preset orientation angles in the same orientation angle order as in step S1025 to form the second multi-directional structural response tensor.

[0050] The method for constructing and filtering multi-directional structural response tensors provided in this step has the following technical advantages compared to existing technologies: The preset set of orientation angles is determined by extracting local peaks from the orientation angle histogram and performing bimodal consistency verification. This allows the selection of orientation angles to be driven by the structural features of the image itself, rather than by using a fixed and uniform angle division. The local peaks of the orientation angle histogram represent the concentrated distribution direction of structural edges in the image. The bimodal consistency verification ensures that the selected orientation angles have significant structural responses in both modal images, so that subsequent filtering operations can focus on the orientation channels with the most registration value in the image.

[0051] By performing clustering operations on candidate direction angles, adjacent candidate angles are merged into cluster centers, avoiding redundancy between preset direction angles. The clustering operation can adaptively group spatially similar candidate angles into one category, so that the final set of preset direction angles is evenly distributed on the angle axis and does not overlap with each other, ensuring the independence and representativeness of each direction channel.

[0052] The constructed anisotropic Gaussian derivative kernel exhibits direction selectivity, enhancing the structural edge response along a preset directional angle while suppressing structures and noise extending outside this direction. The anisotropic Gaussian kernel has a larger standard deviation (stretching) along the directional angle and a smaller standard deviation (compression) in the vertical direction, allowing the kernel's coverage to extend along the structural direction and compress along the vertical direction. This results in a more accurate response to structural boundaries extending along this direction. Calculating the first-order partial derivative of the rotated and stretched Gaussian kernel along the directional angle reveals positive and negative polarity responses from the anisotropic Gaussian derivative kernel to structural edges along this direction. These polarity responses directly correspond to the grayscale transition directions of the physical structural boundaries, providing a foundation for establishing cross-modal structural correspondences using complementary response relationships.

[0053] The first and second modal images are filtered separately. The two modal images use the same set of anisotropic Gaussian derivative kernels to ensure the consistent correspondence between the channels in each direction in the structural response tensor between the two modes. Since the grayscale responses of the two modal images to the same physical structure are opposite, the response values ​​on the same direction channel after filtering show opposite polarities in the two modes. This polarity difference directly reflects the complementary response characteristics of the dual modal images.

[0054] The obtained multi-directional structural response tensor maps the image from the grayscale domain to the structural orientation domain, transforming the image representation from absolute grayscale values ​​to edge structural responses along each orientation. Since the nonlinear grayscale relationship between the neutron image and the X-ray image does not affect the distribution of structural orientations, this mapping effectively reduces the interference of cross-modal grayscale differences on registration, providing a structural feature representation with cross-modal consistency for subsequent steps.

[0055] Step S103: Select effective direction angles from the preset multiple direction angles according to the structural response diagrams corresponding to each direction angle, and determine closed paths in the first modal image according to each effective direction angle, wherein each effective direction angle corresponds to a closed path.

[0056] In this step, the response values ​​of a first structural response map in the first multi-directional structural response tensor and a second structural response map in the second multi-directional structural response tensor at the same pixel position are multiplied, and the product results of all pixel positions are summed to obtain a cross-modal response consistency metric value of a certain preset orientation angle. The first structural response map and the second structural response map are both corresponding to the preset orientation angle. All preset orientation angles are sorted from largest to smallest according to the cross-modal response consistency metric value, and the orientation angles that are in the top preset number of sorts are selected as valid orientation angles; In the first structural response map corresponding to a certain effective direction angle, the pixel positions where the response value along the direction of the certain effective direction angle presents a local maximum and the distance between two adjacent local maxima is greater than a preset distance threshold are extracted. The extracted pixel positions are connected in sequence to form a structural center line under the certain effective direction angle. Taking each pixel position on the center line of a certain structure as the center, extend to both sides in a direction perpendicular to the center line of the certain structure, and find the position with the largest change rate of response value in the first structural response map corresponding to a certain effective direction angle as the boundary point to obtain the left boundary point set and the right boundary point set; The left boundary point set is connected in sequence to form the left boundary line, the right boundary point set is connected in sequence to form the right boundary line, and the left boundary line and the right boundary line are connected at both ends to form a closed path, wherein the closed path surrounds the physical structure boundary corresponding to the effective direction angle.

[0057] In one specific embodiment, the accumulated value of the product of response values ​​in the same directional channel of the dual-modal structural response is used to evaluate the structural corresponding quality of each directional angle, and the directional angle with the highest registration value is selected; then, the structural centerline is extracted from the first structural response map corresponding to each effective directional angle and a closed path is formed, specifically as follows: Step S1031: Calculate the cross-modal response consistency metric. The response values ​​at the same pixel location are multiplied by the response values ​​at the same pixel location in the first multi-directional structural response tensor corresponding to a preset orientation angle in the first multi-directional structural response tensor and the second structural response tensor corresponding to the same preset orientation angle in the second multi-directional structural response tensor. With the i-th preset direction angle θ i For example, let the response value at pixel position (x, y) in the first structural response map be R. 1,i (x,y), the response value at the same pixel position (x,y) in the second structural response map is R. 2,i If (x, y), then the product of the responses at pixel position (x, y) is R. 1,i (x,y)×R 2,i (x,y), Summing the product values ​​of the responses at all pixel locations yields the cross-modal response consistency metric C(θ_i) for the preset orientation angle θ_i. i The physical meaning of this cross-modal response consistency metric is as follows: when the first and second structural response maps exhibit opposite polarities at the same physical structural boundary, i.e., the grayscale responses of the two modal images to the physical structure are opposite, the product value of the responses corresponding to this boundary region is negative; when the two structural response maps exhibit responses of the same polarity, the product value is positive. Because the grayscale responses of the two modal images to the physical structure are opposite, at the correct structural boundary location, R... 1,i (x,y) and R 2,i The signs of (x,y) are opposite, and the product is negative. Therefore, C(θ) i The larger the absolute value of ) and the more negative it is, the higher the consistency of the structural response of the two modes under that orientation angle, and the more registration value that orientation angle has.

[0058] Step S1032: Select the effective direction angle All preset orientation angles are sorted in ascending order of their cross-modal response consistency metric values. Since a negative cross-modal response consistency metric value indicates that the two modal structural responses are opposite, and a smaller value indicates a higher quality structure, the orientation angles ranked higher in the ascending order correspond to the optimal orientation angles. Select the first preset number of direction angles in the sorting as valid direction angles. The preset number can be set according to actual needs, such as 4 to 6. Alternatively, a consistency threshold can be set, and direction angles with cross-modal response consistency metric values ​​less than the consistency threshold can be selected as valid direction angles.

[0059] Step S1033, extract the structural centerline For a certain effective direction angle θ eff Obtain the first structural response map R corresponding to the effective orientation angle. 1,eff In the first structural response diagram, along a certain effective direction angle θ eff The direction of the scan extracts the pixel positions where the response value presents a local maximum and the distance between two adjacent local maxima is greater than a preset distance threshold.

[0060] Specifically, along a certain effective direction angle θ eff The direction, for each scan line, when θ eff Scanning along the line when approaching the horizontal direction, when θ eff When the direction is close to vertical, the scanning is performed along the column; otherwise, the scanning is performed at an angle. The pixel position on the scan line whose response value is greater than the response value of its adjacent pixels is identified as the candidate maximum position. For two adjacent candidate maximum positions, the distance between them along the scan line is calculated. If the distance is greater than a preset distance threshold, such as 5 pixels, both candidate maximum positions are retained. If the distance is not greater than the preset distance threshold, only the one with the larger response value is retained.

[0061] All the candidate maximum positions are connected sequentially according to their order on the scan line to obtain one or more line segments. If multiple line segments are extracted under a certain effective direction angle, the longest line segment is selected as the structural center line under the effective direction angle.

[0062] Each pixel position on the center line of the structure represents the center position of the physical structure boundary in the effective angular direction.

[0063] Step S1034: Determine the left boundary point and the right boundary point. Taking each pixel position on the structural center line as the center, extend to both sides in a direction perpendicular to the structural center line. The extension method is as follows: for each pixel position on the structural center line, in the first structural response map, move from that pixel position to the first side and the second side pixel by pixel in a direction perpendicular to the structural center line, and record the response value of each moving position until the moving distance reaches the preset extension distance, such as 20 pixels.

[0064] During the extension process, the location with the largest rate of change of response value is identified as the boundary point. The rate of change of response value is the absolute value of the difference between the response values ​​of adjacent pixel locations. The location with the largest rate of change of response value corresponds to the transition position between the physical structure boundary and the surrounding area, i.e., the boundary point.

[0065] For each pixel position on the center line of the structure, find the position with the largest rate of change of response value on the path extending to the first side and record it as the first side boundary point; find the position with the largest rate of change of response value on the path extending to the second side and record it as the second side boundary point. All first side boundary points constitute the left boundary point set, and all second side boundary points constitute the right boundary point set.

[0066] Step S1035: Form a closed path Connect the left boundary points in order of their corresponding positions on the structure's center line to form the left boundary line, and connect the right boundary points in order of their corresponding positions on the structure's center line to form the right boundary line.

[0067] Connect the first endpoint of the left boundary line to the first endpoint of the right boundary line, and connect the second endpoint of the left boundary line to the second endpoint of the right boundary line to form a closed path around the physical structure boundary.

[0068] The effective orientation angle selection and closed path determination method provided in this step has the following technical advantages compared to existing technologies: Valid orientation angles are selected by using a cross-modal response consistency metric. The cross-modal response consistency metric is the sum of the products of the response values ​​of the first and second structural response maps in the same directional channel. Since the grayscale responses of the two-modal images to the same physical structure are opposite, the correctly corresponding structural boundary exhibits opposite response polarities in the structural response maps of the two modalities, resulting in a negative product. This makes the cross-modal response consistency metric of the orientation angle a large negative value. By selecting the orientation angles that rank in the top preset number as valid orientation angles, it is ensured that the orientation angles used for subsequent registration have the highest cross-modal structural correspondence quality.

[0069] When selecting effective orientation angles, the complementary characteristics of structural orientation distribution in bimodal images are fully utilized. The local peaks in the orientation angle histogram correspond to the concentrated distribution direction of structural edges in the image. Bimodal consistency verification is performed, i.e. the frequency of the same angle distribution in the second orientation angle histogram is greater than the preset frequency threshold, to ensure that the selected orientation angle has a significant structural response in both modal images. Through clustering operations, multiple local peaks are merged into cluster centers, so that the effective orientation angles are evenly distributed on the angle axis and do not overlap, thus avoiding redundancy of orientation angles.

[0070] The structural centerline is formed by extracting the local maxima locations along the effective direction angle from the first structural response map. The local maxima condition ensures that the extracted pixel locations correspond to the centerline of the physical structural boundary. The condition that the distance between adjacent maxima is greater than a preset distance threshold prevents over-extraction caused by noise or small fluctuations, thus ensuring the continuity and stability of the structural centerline.

[0071] Boundary points are determined by searching for the location with the largest rate of change of response value along the direction perpendicular to the center line of the structure. The location with the largest rate of change of response value corresponds to the transition position between the physical structure boundary and the surrounding background, i.e. the actual outline of the boundary. By searching to both sides and connecting them to form the left and right boundary lines respectively, and then connecting them at both ends to form a closed path, a complete closed enclosure of the physical structure boundary is achieved, providing a clear spatial reference for finding the corresponding mapping closed path in the second modal image.

[0072] For each valid orientation angle, a corresponding closed path is determined, so that each orientation angle corresponds one-to-one with the closed path. This correspondence is maintained in subsequent steps, so that the registration contribution of each orientation angle can be calculated and weighted independently, providing a structural basis for the weighted summation of the global registration fit.

[0073] Step S104: Apply the current spatial transformation parameters to the closed path and map it to the second modal image to obtain an initial mapped closed path. Then, perform position correction on the initial mapped closed path according to the second multi-directional structural response tensor to obtain the target mapped closed path.

[0074] In this step, the current spatial transformation parameters are obtained. The current spatial transformation parameters include horizontal translation, vertical translation, rotation angle, uniform scaling factor, and shearing parameter. The initial value of the current spatial transformation parameters is determined by the global coarse registration result between the first modal image and the second modal image. During the fine registration iteration process, the current spatial transformation parameters are updated with each iteration. Obtain a first internal region enclosed by a certain closed path, and calculate the first area value of the first internal region in the first modal image; After applying the current spatial transformation parameters to a certain closed path, the initial position of the mapping path is obtained, and a second internal region is formed in the second modality image based on the initial position of the mapping path; Centered on the centroid of the second internal region, and with the first area value as the target area, the region is expanded radially outward layer by layer in the second modal image until the area of ​​the expanded region is equal to the first area value. The boundary of the expanded region is then used as the first candidate mapping closed path. Determine the response polarity of each pixel position on the first candidate mapping closed path in the direction of the effective direction angle corresponding to a certain closed path in the second multi-directional structure response tensor, and count the number of first pixel positions with negative response polarity and the number of second pixel positions with positive response polarity, wherein the response polarity is the sign of the response value; If the absolute value of the difference between the number of first pixel positions and the number of second pixel positions is less than a preset difference threshold, then the first candidate mapping closed path is used as the initial mapping closed path. If the absolute value of the difference between the number of first pixel positions and the number of second pixel positions is not less than the preset difference threshold, then it is determined whether the number of first pixel positions is greater than the number of second pixel positions. If yes, then the centroid is moved by a preset step along the distribution centroid direction of the pixel positions with negative response polarity on the boundary of the extended region to obtain the updated centroid position; if no, then the centroid is moved by a preset step along the distribution centroid direction of the pixel positions with positive response polarity on the boundary of the extended region to obtain the updated centroid position. Using the updated centroid position as the center and the first area value as the target area, the layer-by-layer expansion is performed again to obtain the second candidate mapping closed path. The response polarity is counted again. If the condition is still not met, the centroid is moved and expanded again until the condition is met or the preset maximum number of moves is reached. The candidate mapping closed path corresponding to the preset maximum number of moves is taken as the initial mapping closed path.

[0075] Further, obtain a certain effective direction angle corresponding to a certain initial mapping closed path, and the first response value of each pixel position on the closed path corresponding to the certain initial mapping closed path in the first multi-directional structure response tensor along the direction of the certain effective direction angle, and obtain the first response polarity of the first response value, wherein the first response polarity is the positive or negative sign of the first response value; Invert the first response polarity of each pixel position on the initial mapping closed path to obtain the expected response polarity of each pixel position. Centered on a pixel position on a certain initial mapping closed path, a first direction and a second direction search are performed along the direction perpendicular to the direction indicated by a certain effective direction angle. Within the search range, the second response value of each search position in the second multi-directional structural response tensor along the direction of the certain effective direction angle is obtained, and the second response polarity of the second response value is obtained. The search range is a first preset distance extended from a pixel position in the first direction and a second preset distance extended in the second direction, wherein the first direction and the second direction are opposite to each other. All search positions within the search range whose second response polarity matches the expected response polarity are formed into a candidate position set. The search position with the largest absolute value of the second response value is selected from the candidate position set as the optimized candidate position for a certain pixel position. Connect all the candidate positions of the optimized pixel positions in sequence according to the original order of the pixel positions on the initial mapping closed path to obtain a target mapping closed path corresponding to the initial mapping closed path.

[0076] In one specific embodiment, the closed path in the first modal image is first mapped to the second modal image using the current spatial transformation parameters to obtain an initial mapped closed path; then, the mapped path is corrected pixel-by-pixel according to the second multi-directional structural response tensor, so that the corrected path is aligned with the actual physical structure boundary in the second modal image, specifically as follows: Step S1041: Obtain the current spatial transformation parameters Obtain the current spatial transformation parameters, which are affine transformation parameters, including horizontal translation, vertical translation, rotation angle, uniform scaling factor, and shearing parameter. The current spatial transformation parameters are used to describe the spatial mapping relationship from the second modality image to the first modality image.

[0077] The initial values ​​of the current spatial transformation parameters are determined by the global coarse registration results between the first modality image and the second modality image. The global coarse registration is a multi-resolution pyramid optimization performed in the grayscale domain using mutual information as a similarity metric. Its output is the affine transformation parameters after coarse registration. During the fine registration iteration process, the current spatial transformation parameters are updated with each iteration. That is, after each iteration, the current spatial transformation parameters are updated to the parameter values ​​optimized in this iteration, which serve as the initial values ​​for the next iteration.

[0078] Step S1042: Calculate the area of ​​the region enclosed by the closed path. For a certain closed path C, obtain the first internal region enclosed by the closed path. The first internal region is the set of pixels enclosed by the closed path. Calculate the first area value A1 of the first internal region in the first modal image. The area value is the number of pixels contained in the first internal region multiplied by the area of ​​a single pixel.

[0079] Step S1043: Generate the initial mapping closed path The initial position of the mapped path is obtained by applying the current spatial transformation parameters to a closed path. Specifically, for each pixel position on a closed path, the horizontal coordinate x and vertical coordinate y of the pixel position are used as the first coordinate value and the second coordinate value, and the mapped horizontal coordinate x′ and mapped vertical coordinate y′ of the pixel position are calculated according to the following formula: x′=s·(x·cosθ-y·sinθ)+s h ·(x·sinθ+y·cosθ)+t x , y′=s·(x·sinθ+y·cosθ)+t y , In the formula, θ is the rotation angle in the current spatial transformation parameters, and t x t represents the horizontal translation amount in the current spatial transformation parameters. y s is the vertical translation in the current spatial transformation parameters, and s is the uniform scaling factor in the current spatial transformation parameters. h This refers to the shearing parameter in the current spatial transformation parameters.

[0080] Where x·cosθ-y·sinθ represents the first coordinate value after rotating the pixel position around the origin by an angle θ, and x·sinθ+y·cosθ represents the second coordinate value after rotation; s·(x·cosθ-y·sinθ) and s·(x·sinθ+y·cosθ) represent the results of uniformly scaling the rotated coordinate values, respectively; s h ·(x·sinθ+y·cosθ) represents the effect of the shearing operation on the x-axis, with the shearing amount applied along the rotated second coordinate direction; plus t x and t y This means that after rotation, scaling, and shearing are performed, the results will be translated horizontally and vertically, respectively.

[0081] The mapped coordinates of all pixel positions on a closed path are connected sequentially according to the original order of the pixel positions on the closed path to obtain the mapped path, which encloses the second inner region in the second modality image.

[0082] Step S1044: Generation of candidate mapping closed paths based on area preservation Centered on the centroid of the second inner region, and with the first area value as the target area, the region expands radially outward layer by layer in the second modal image. The centroid is the average value of the position coordinates of all pixels in the second inner region. The expansion method is as follows: based on the current boundary, expand one layer of pixels radially outward to form a new expansion region; calculate the area of ​​the new expansion region. If the area is still less than the first area value, continue to expand one layer outward until the area of ​​the expansion region is equal to the first area value.

[0083] The boundary of the extended region is used as the first candidate mapping closure path, the area of ​​which is the first area value, but the position and shape may not yet be fully aligned with the actual structural boundary in the second modality image.

[0084] Step S1045: Verification of candidate mapping closed paths based on response polarity statistics Determine the response polarity of each pixel position on the first candidate mapping closed path in the effective direction angle direction corresponding to a certain closed path in the second multi-directional structure response tensor. The response polarity is the sign of the response value: the response polarity is positive when the response value is positive and negative when the response value is negative.

[0085] The number N of the first pixel positions with negative response polarity. neg The number N of second pixel positions with positive response polarity pos .

[0086] If N neg With N pos The absolute value of the difference is less than the preset difference threshold T diff If the distribution of positive and negative response polarities on the first candidate mapping closed path is relatively balanced, and the path is basically aligned with the structural boundary in the second modality image, then the first candidate mapping closed path is taken as the initial mapping closed path.

[0087] If the absolute value is not less than the preset difference threshold, it indicates that the distribution of positive and negative response polarities is unbalanced, and the path is not yet aligned with the structural boundary in the second modality image, so the centroid position needs to be adjusted.

[0088] Step S1046, Centroid Position Adjustment and Iterative Expansion When N neg With N pos If the absolute value of the difference is not less than a preset difference threshold, then determine N. neg Is it greater than N? pos .

[0089] If N neg >N posIf the first candidate mapping closed path has a majority of negative polarity responses, the centroid is located on one side of the actual structure boundary, that is, the side with negative response polarity. The centroid is moved by a preset step along the centroid direction of the distribution of pixels with negative response polarity on the boundary of the extended region to obtain the updated centroid position.

[0090] If N neg ≤N pos If the positive polarity response is the majority on the first candidate mapping closed path, the centroid is moved by a preset step along the centroid direction of the distribution of pixels with positive response polarity on the boundary of the extended region to obtain the updated centroid position.

[0091] Using the updated centroid position as the center and the first area value as the target area, perform layer-by-layer expansion again to obtain the second candidate mapping closed path, and then recalculate the response polarity. If the condition is still not met, i.e., N... neg With N pos If the absolute value of the difference is still not less than the preset difference threshold, then continue to move the centroid and expand again until the condition is met or the preset maximum number of moves is reached.

[0092] The candidate mapping closed path corresponding to the condition is used as the initial mapping closed path. If the condition is still not met after reaching the preset maximum number of moves, the candidate mapping closed path corresponding to the preset maximum number of moves is used as the initial mapping closed path.

[0093] Step S1047: Pixel-by-pixel correction of the mapping path based on response polarity inversion Obtain the closed path C of a given initial mapping init The corresponding effective direction angle θ eff Obtain the first response value of each pixel position on the closed path corresponding to the initial mapping closed path in the first multi-directional structural response tensor along the effective directional angle direction, and obtain the first response polarity of the first response value, which is the positive or negative sign of the first response value.

[0094] The first response polarity of each pixel position on the initial mapping closed path is inverted to obtain the expected response polarity of each pixel position. The physical basis of this inversion operation is that the grayscale response of the two modal images to the same physical structure is opposite. Therefore, the response polarity at the boundary of a certain physical structure in the first modal image should be opposite in the second modal image.

[0095] Centered on a pixel position P on the initial mapping closed path, a first direction and a second direction search are performed along the direction perpendicular to the direction indicated by the effective direction angle. The first direction and the second direction are opposite to each other and are both perpendicular to the direction indicated by the effective direction angle.

[0096] Within the search range, the second response value of each search position in the second multi-directional structural response tensor along the effective directional angle direction is obtained, and the second response polarity of the second response value is obtained. The search range is a first preset distance extended from a certain pixel position along the first direction and a second preset distance extended along the second direction.

[0097] All search positions within the search range whose second response polarity matches the expected response polarity are formed into a candidate position set. The search position with the largest absolute value of the second response value is selected from the candidate position set as the optimized candidate position for a certain pixel position.

[0098] If there is no search position within the search range that matches the expected response polarity, then the pixel position itself is used as the optimization candidate position. This processing method ensures that each pixel position has a corresponding optimization candidate position, and the path will not be broken due to local search failure.

[0099] Step S1048: Form a closed path for target mapping. Connect all the candidate positions for optimization of pixels in the original order of the pixel positions on the initial mapping closed path to obtain the target mapping closed path corresponding to the initial mapping closed path.

[0100] The closed path mapping and correction method provided in this step has the following technical advantages compared to existing technologies: The strategy of "area preservation + centroid adjustment + layer-by-layer expansion" is adopted to determine the initial mapping closed path in the second modal image. Area preservation ensures that the area enclosed by the mapping path is consistent with the area enclosed by the corresponding closed path in the first modal image, so that the mapping path matches the original closed path in scale. Centroid adjustment is driven by the statistical distribution of response polarity, so that the mapping path gradually moves closer to the actual structural boundary. Layer-by-layer expansion ensures the continuity of the candidate mapping closed path.

[0101] The introduced centroid adjustment mechanism based on response polarity determines the path offset direction by statistically analyzing the distribution of positive and negative response polarities on the candidate mapping closed path. When negative polarity responses predominate, the centroid moves along the direction of the centroid of the negative polarity distribution; when positive polarity responses predominate, the centroid moves along the direction of the centroid of the positive polarity distribution. This adjustment direction is determined entirely by structural response data, without relying on any external reference points or human intervention, thus ensuring the objectivity and accuracy of the adjustment direction.

[0102] In the pixel-by-pixel correction stage, the first response polarity is inverted and taken as the expected response polarity. Based on the physical characteristic of opposite gray-scale responses in bimodal images, a judgment criterion for cross-modal structure correspondence is established. When the response polarity of a certain search position in the second modal image is consistent with the expected response polarity, the position is included in the candidate position set. Furthermore, the position with the largest absolute value of the response value is selected from the candidate positions as the optimization candidate position, that is, the position with the strongest structural response is selected as the final mapping position. This criterion makes the path correction process driven by physical constraints and has interpretability.

[0103] By setting a fallback mechanism during the correction phase, which states that "if there is no search position within the search range that matches the expected response polarity, then the current pixel position itself will be used as the optimization candidate position," the integrity of the path is ensured. This mechanism prevents the correction process from being interrupted due to local search failures, guarantees the continuity of the corrected path, and avoids the problem of subsequent steps being unable to be executed due to path breakage.

[0104] The determination of the initial mapping closure path and the pixel-by-pixel correction are divided into two stages. The first stage determines the initial mapping closure path based on area preservation and centroid adjustment. The second stage performs pixel-by-pixel position correction based on response polarity inversion. The two stages work together to ensure that the final target mapping closure path can be accurately aligned with the actual boundary position of the physical structure in the second modality image, thereby providing a spatially accurate path for the extraction of the subsequent response sequence.

[0105] Step S105: Extract a first response value sequence from the first modal image along the closed path, and extract a second response value sequence from the second modal image along the target mapping closed path.

[0106] In this step, for a certain closed path, the effective direction angle corresponding to the certain closed path is obtained; The total number of pixels on a closed path is determined. Using the total number of pixels as the number of traversals, the Euclidean distances between adjacent pixel positions are calculated sequentially and accumulated to obtain the total length of the closed path. Starting from any point on the closed path, sampling is performed in a clockwise direction with equal step sizes to obtain multiple sampling points arranged in the path order. The step size is the total length divided by the preset number of sampling points. For each sampling point, the response value of the sampling point along the effective direction angle in the first multi-directional structural response tensor is obtained, and the response values ​​are arranged according to the path order of the sampling points to obtain the first response value sequence corresponding to a certain closed path.

[0107] In one specific embodiment, sampling is performed at equal steps along a closed path. At each sampling point, the response values ​​along the effective direction angles in the multi-directional structural response tensor are obtained, forming an ordered sequence of response values. Specifically: Step S1051: Obtain the effective direction angle of the closed path. For a certain closed path, obtain the effective direction angle corresponding to the closed path. The effective direction angle is the direction angle used to determine the closed path in step S103, and the direction angle corresponds one-to-one with the closed path.

[0108] Step S1052, determine sampling parameters The total number of pixels on a closed path is determined. The total number of pixels is the number of consecutive pixel positions that make up a closed path. A closed path is a closed curve whose starting point and ending point coincide. Therefore, the starting point is not counted repeatedly in the total number of pixels.

[0109] Using the total number of pixels as the number of traversals, calculate the Euclidean distance between adjacent pixel positions and sum them up to obtain the total length L of a closed path.

[0110] Get the preset number of sampling points M. The preset number of sampling points M is a pre-set positive integer used to control the sampling density along the closed path. The value of M can be set according to actual needs, for example, set to 100.

[0111] Calculate the step size Δs = L / M, where the step size is the path length between every two adjacent sampling points along a closed path.

[0112] Step S1053, sampling with equal step size Starting from any point on a closed path, for example, taking the first pixel position P1 in the pixel position sequence as the starting point, proceed clockwise along a closed path.

[0113] Starting from the starting point, a sampling point is recorded for each step Δs. Specifically, the path length is accumulated in a clockwise direction from the starting point. When the accumulated length reaches Δs, the current position is recorded as the first sampling point. Continue moving forward, and when the accumulated length reaches 2×Δs, the current position is recorded as the second sampling point. This process continues until the entire closed path is completed.

[0114] If the theoretical landing point of a sampling point falls between two adjacent pixel positions, the precise coordinates of the sampling point are calculated using linear interpolation. The linear interpolation method is as follows: Let the two adjacent pixel positions be A and B, and the path length of the sampling point from A be d. A The path length from B is d. B Then the coordinates of the sampling point are (d B ·x A +d A ·x B ) / (d A +d B ), (d B ·yA +d A ·y B ) / (d A +d B (x) A y A (x) represents the coordinates of pixel position A in the coordinate plane of the first modality image. B y B ) represents the coordinates of pixel position B in the coordinate plane of the first modality image.

[0115] After sampling is completed, M sampling points are obtained in the order of the path.

[0116] Step S1054: Obtain response values ​​and form a sequence. For each sampling point, obtain the effective direction angle θ of the sampling point in the first multi-directional structural response tensor. eff The response value in the direction, specifically, with respect to sampling point S. k Using the spatial location as an index, locate the position on the spatial dimension of the first multi-directional structural response tensor; at this position, extract the effective direction angle θ. eff The response value R on the corresponding directional channel 1,eff(Sk) , Arrange the response values ​​corresponding to all sampling points in the order of the sampling points to obtain the first response value sequence corresponding to a certain closed path. The length of the first response value sequence is M, which is consistent with the preset number of sampling points M. Step S1055: Extract the second response value sequence in the same manner. The target mapping closed path is processed in the same way, and the second response value sequence is extracted along the target mapping closed path.

[0117] Step S106: Perform sequence transformation on the second response value sequence to obtain the inverted response value sequence, calculate the sequence matching degree between the first response value sequence and the inverted response value sequence, and use the weighted sum of the sequence matching degrees corresponding to all valid direction angles as the global registration matching degree.

[0118] In this step, a certain second response value sequence corresponds to a certain closed path, and a certain first response value sequence corresponds to a certain closed path. The certain closed path is obtained by mapping a certain closed path.

[0119] Multiply each second response value in the given second response value sequence by -1 to obtain a reverse response value sequence.

[0120] The physical basis of this inversion operation is that the grayscale responses of the first modal image and the second modal image to the same physical structure are opposite. Therefore, the response values ​​at corresponding positions of the first response value sequence and the second response value sequence extracted along the boundary of the same physical structure have opposite signs. After multiplying the second response value sequence by -1, its sign is consistent with the sign of the first response value sequence, so that the two sequences show the same direction of change at the boundary of the same physical structure, thereby allowing direct calculation of the positive correlation between the two.

[0121] The first response value sequence and the reverse response value sequence are normalized respectively.

[0122] Calculate a Pearson correlation coefficient between a normalized first response value sequence and a normalized inverse response value sequence, and use the Pearson correlation coefficient as the directional fit of the effective direction angle corresponding to the closed path.

[0123] The global registration match is obtained by weighted summation of the directional match rates of all effective directional angles. The weight of each effective directional angle in the weighted summation is the proportion of the energy value of the structure response map corresponding to the effective directional angle to the total energy value of the structure response maps corresponding to all effective directional angles. The energy value is the sum of the squares of the response values ​​of each pixel position in the structure response map.

[0124] The physical meaning of this weighted summation method is as follows: the higher the structural response energy of the directional angle, the richer the structural information contained in its structural response map, and the higher its importance in the registration optimization. Therefore, it occupies a larger weight in the global registration fit. The lower the structural response energy of the directional angle, the less structural information contained in its structural response map, and the smaller its weight in the global registration fit. Through this weighting method, the registration optimization process automatically focuses on the directional channel with the richest structural information.

[0125] Step S107: Adjust the current spatial transformation parameters according to the global registration matching degree, apply the adjusted current spatial transformation parameters to the second modal image, align the second modal image with the first modal image, and output the registration result.

[0126] In this step, with the goal of maximizing the global registration fit, the Powell algorithm is used to iteratively update the current spatial transformation parameters. In each iteration step, after updating the current spatial transformation parameters, the updated global registration fit is redefined; If the updated global registration match is greater than the original global registration match, then the updated current spatial transformation parameters are accepted. If the updated global registration match is not greater than the original global registration match, then the updated current spatial transformation parameters are rejected, and the original current spatial transformation parameters remain unchanged. When the difference in global registration fit between two adjacent iterations is less than a preset convergence threshold, the iteration is terminated, and the current spatial transformation parameters at the time of termination are used as the final registration parameters. The final registration parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

[0127] In one specific embodiment, with the goal of maximizing global registration fit, the Powell algorithm is used to iteratively update the current spatial transformation parameters. After the convergence condition is met, the final registration parameters are applied to the second modality image to achieve image alignment. Step S1071, set the initial optimization parameters The current spatial transformation parameters are obtained as the initial values ​​for the optimization iteration. These initial values ​​are determined by the global coarse registration results between the first and second modal images. The global coarse registration is a multi-resolution pyramid optimization performed in the grayscale domain using mutual information as a similarity metric. Its output is the affine transformation parameters after coarse registration, including initial estimates of horizontal translation, vertical translation, rotation angle, uniform scaling factor, and shearing parameters. Set the maximum number of iterations K for the Powell algorithm. max For example, K max =100, and a preset convergence threshold ε, for example, ε=10 -4 The maximum number of iterations is used to prevent the optimization process from looping indefinitely, and the preset convergence threshold is used to determine whether the optimization has reached the vicinity of the optimal solution.

[0128] Step S1072: Perform Powell's algorithm iterative update. With the goal of maximizing global registration fit, the Powell algorithm is used to iteratively update the current spatial transformation parameters. The Powell algorithm is a derivative-free optimization method that does not require explicit gradient calculation, and its basic process is as follows: At the beginning of the k-th iteration, a set of linearly independent search directions D={d1,d2,d3,d4,d5} is determined. The five search directions correspond to the five parameters of the affine transformation: horizontal translation, vertical translation, rotation angle, uniform scaling factor, and shearing parameter. The search directions are unit vectors in five-dimensional space, and each search direction indicates the direction of change of a parameter.

[0129] For each search direction d i From the current parameter point P kStarting from this point, perform a one-dimensional search along this direction to find the step size α that maximizes the global registration match. i One-dimensional search can be performed using the golden section search or parabolic interpolation method. The search interval is gradually narrowed within the preset search interval until the optimal step size is found. The global registration fit is the global registration fit value obtained after substituting the affine parameters corresponding to the current parameter point and executing steps S102 to S106.

[0130] Along the search direction d i Movement step size α i Then, the new parameter point P is obtained. k' =P k +α i ·d i And calculate the global registration fit corresponding to the new parameter point.

[0131] The above one-dimensional search is performed sequentially in the five search directions to complete one complete iteration. After one iteration, the parameter change ΔP = P before and after this iteration is calculated. k' -P k The parameter change is added as a new search direction to the search direction set, replacing the original search direction that causes the largest change in global registration fit, in order to accelerate the convergence of subsequent iterations.

[0132] Step S1073: Determine whether to accept parameter updates. If the updated global registration match is greater than the original global registration match, then accept the updated current spatial transformation parameters and set P... k' As the starting point for the next iteration, If the updated global registration match is greater than the original global registration match, then the updated spatial transformation parameters are rejected, and the original spatial transformation parameters are kept unchanged. This ensures that the global registration match remains monotonically constant during the optimization process, guaranteeing that the algorithm eventually converges to a local optimum. Step S1074: Determine the convergence condition When the difference in global registration accuracy between two consecutive iterations is less than a preset convergence threshold ε, the optimization process is considered to have converged, and the iteration is terminated. Alternatively, when the number of iterations k reaches a preset maximum number of iterations K... max When the iteration ends, the iteration is terminated.

[0133] The current spatial transformation parameters at the time of termination are used as the final registration parameters P*.

[0134] Step S1075: Apply the final registration parameters and output the registration result. The final registration parameter P* is applied to the second modality image. Specifically, for each pixel position (x, y) in the second modality image, its mapping coordinates (x', y') are calculated according to the affine transformation formula described in step S104, and the gray value of the pixel position is assigned to the mapping position (x', y') in the coordinate plane of the first modality image.

[0135] If the mapped position (x', y') falls on a non-integer coordinate position, then bilinear interpolation or nearest neighbor interpolation is used to determine the gray value at that position. The bilinear interpolation method is to take a weighted average of the gray values ​​of the four integer coordinate pixels around the mapped position, and the weights are determined by the distance from the mapped position to each integer coordinate pixel.

[0136] After mapping all pixel positions, a second modal image aligned with the first modal image is obtained. The registration result is output, which includes the aligned second modal image, the final registration parameter P*, and the final global registration fit value, for subsequent 3D visualization, dual-modal fusion, or quantitative analysis.

[0137] In summary, the method of this application acquires a first modal image and a second modal image of the same sample, with opposite grayscale responses to the same physical structure. Anisotropic Gaussian derivative kernels are constructed using multiple preset orientation angles, and the two modal images are filtered to obtain a first multi-directional structural response tensor and a second multi-directional structural response tensor. Effective orientation angles are selected based on the structural response maps corresponding to each orientation angle, and closed paths are determined in the first modal image based on each effective orientation angle. The closed paths are mapped to the second modal image to obtain an initial mapped closed path, and position correction is performed to obtain a target mapped closed path. A first response value sequence is extracted along the closed path, and a second response value sequence is extracted along the target mapped closed path. The second response value sequence undergoes a sequence transformation, and its consistency with the first response value sequence is calculated. The weighted sum of the consistency of all effective orientation angles is used as the global registration consistency. Spatial transformation parameters are adjusted based on the global registration consistency to align the second modal image with the first modal image, and the registration result is output. By utilizing the physical property that the grayscale responses of bimodal images are opposite, a cross-modal structural correspondence can be established, enabling stable automatic registration of bimodal tomographic images.

[0138] Please see Figure 2 The diagram shows a structural block diagram of a dual-modal tomographic image registration system based on multi-directional structural features according to this application.

[0139] like Figure 2 As shown, the dual-modal tomographic image registration system 200 includes an acquisition module 210, a filtering module 220, a determination module 230, a mapping module 240, an extraction module 250, a calculation module 260, and an output module 270.

[0140] The acquisition module 210 is configured to acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; the filtering module 220 is configured to construct anisotropic Gaussian derivative kernels with preset multiple orientation angles, and filter the first modal image and the second modal image according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image, wherein the multi-directional structural response tensor contains a structural response map corresponding to each orientation angle; the determination module 230 is configured to select effective orientation angles from the preset multiple orientation angles according to the structural response maps corresponding to each orientation angle, and determine a closed path in the first modal image according to each effective orientation angle, wherein each effective orientation angle corresponds to a closed path; the mapping module 240 is configured to apply the closed path to the current The spatial transformation parameters are mapped to the second modal image to obtain an initial mapping closed path. The position of the initial mapping closed path is corrected according to the second multi-directional structural response tensor to obtain a target mapping closed path. The extraction module 250 is configured to extract a first response value sequence in the first modal image along the closed path and extract a second response value sequence in the second modal image along the target mapping closed path. The calculation module 260 is configured to perform a sequence transformation on the second response value sequence to obtain an inverted response value sequence, calculate the sequence matching degree between the first response value sequence and the inverted response value sequence, and use the weighted sum of the sequence matching degrees corresponding to all effective direction angles as the global registration matching degree. The output module 270 is configured to adjust the current spatial transformation parameters according to the global registration matching degree, apply the adjusted current spatial transformation parameters to the second modal image to align the second modal image with the first modal image, and output the registration result.

[0141] It should be understood that Figure 2 The modules and references described in the document Figure 1 The steps described in the text correspond to those in the method described above. Therefore, the operations, features, and corresponding technical effects described above also apply to the method described in the text. Figure 2 The various modules in the document will not be described in detail here.

[0142] In other embodiments, the present invention also provides a computer-readable storage medium having a computer program stored thereon, wherein when the program instructions are executed by a processor, the processor performs the dual-modal tomographic image registration method based on multi-directional structural features in any of the above method embodiments. In one embodiment, the computer-readable storage medium of the present invention stores computer-executable instructions, which are configured as follows: Acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; Anisotropic Gaussian derivative kernels are constructed using multiple preset orientation angles, and the first modal image and the second modal image are filtered according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image. The multi-directional structural response tensor contains the structural response map corresponding to each orientation angle. Effective orientation angles are selected from the preset multiple orientation angles based on the structural response diagrams corresponding to each orientation angle, and closed paths are determined in the first modal image based on each effective orientation angle, wherein each effective orientation angle corresponds to a closed path; The closed path is mapped to the second modal image after applying the current spatial transformation parameters to obtain the initial mapped closed path. The position of the initial mapped closed path is then corrected according to the second multi-directional structural response tensor to obtain the target mapped closed path. Extract a first response value sequence in the first modal image along the closed path, and extract a second response value sequence in the second modal image along the target mapping closed path; The second response value sequence is transformed to obtain the inverted response value sequence. The sequence matching degree between the first response value sequence and the inverted response value sequence is calculated, and the weighted sum of the sequence matching degrees corresponding to all valid direction angles is used as the global registration matching degree. The current spatial transformation parameters are adjusted according to the global registration fit, and the adjusted current spatial transformation parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

[0143] Computer-readable storage media may include a program storage area and a data storage area, wherein the program storage area may store an operating system and an application program required for at least one function; the data storage area may store data created based on the use of a dual-modal tomographic image registration system based on multi-directional structural features. Furthermore, the computer-readable storage medium may include high-speed random access memory, and may also include memory, such as at least one disk storage device, flash memory device, or other non-volatile solid-state storage device. In some embodiments, the computer-readable storage medium may optionally include memory remotely disposed relative to a processor, which can be connected to the dual-modal tomographic image registration system based on multi-directional structural features via a network. Examples of such networks include, but are not limited to, the Internet, corporate intranets, local area networks, mobile communication networks, and combinations thereof.

[0144] Figure 3 This is a schematic diagram of the structure of the electronic device provided in the embodiment of the present invention, such as... Figure 3 As shown, the device includes a processor 310 and a memory 320. The electronic device may also include an input device 330 and an output device 340. The processor 310, memory 320, input device 330, and output device 340 can be connected via a bus or other means. Figure 3 Taking a bus connection as an example, the memory 320 is the computer-readable storage medium described above. The processor 310 executes various server functions and data processing by running non-volatile software programs, instructions, and modules stored in the memory 320, thereby implementing the dual-modal tomographic image registration method based on multi-directional structural features as described in the above embodiment. The input device 330 can receive input digital or character information and generate key signal inputs related to user settings and function control of the dual-modal tomographic image registration system based on multi-directional structural features. The output device 340 may include a display screen or other display device.

[0145] The aforementioned electronic device can execute the method provided in the embodiments of the present invention, and has the corresponding functional modules and beneficial effects for executing the method. Technical details not described in detail in this embodiment can be found in the method provided in the embodiments of the present invention.

[0146] In one implementation, the above-described electronic device is applied in a dual-modal tomographic image registration system based on multi-directional structural features, serving as a client, and includes: at least one processor; and a memory communicatively connected to the at least one processor; wherein the memory stores instructions executable by the at least one processor, the instructions being executed by the at least one processor to enable the at least one processor to: Acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; Anisotropic Gaussian derivative kernels are constructed using multiple preset orientation angles, and the first modal image and the second modal image are filtered according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image. The multi-directional structural response tensor contains the structural response map corresponding to each orientation angle. Effective orientation angles are selected from the preset multiple orientation angles based on the structural response diagrams corresponding to each orientation angle, and closed paths are determined in the first modal image based on each effective orientation angle, wherein each effective orientation angle corresponds to a closed path; The closed path is mapped to the second modal image after applying the current spatial transformation parameters to obtain the initial mapped closed path. The position of the initial mapped closed path is then corrected according to the second multi-directional structural response tensor to obtain the target mapped closed path. Extract a first response value sequence in the first modal image along the closed path, and extract a second response value sequence in the second modal image along the target mapping closed path; The second response value sequence is transformed to obtain the inverted response value sequence. The sequence matching degree between the first response value sequence and the inverted response value sequence is calculated, and the weighted sum of the sequence matching degrees corresponding to all valid direction angles is used as the global registration matching degree. The current spatial transformation parameters are adjusted according to the global registration fit, and the adjusted current spatial transformation parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

[0147] Through the above description of the embodiments, those skilled in the art can clearly understand that each embodiment can be implemented by means of software plus necessary general-purpose hardware platforms, and of course, it can also be implemented by hardware. Based on this understanding, the above technical solutions, in essence or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product can be stored in a computer-readable storage medium, such as ROM / RAM, magnetic disk, optical disk, etc., including several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute the methods of various embodiments or some parts of embodiments.

[0148] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A dual-modal tomographic image registration method based on multi-directional structural features, characterized in that, include: Acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; Anisotropic Gaussian derivative kernels are constructed using multiple preset orientation angles, and the first modal image and the second modal image are filtered according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image. The multi-directional structural response tensor contains the structural response map corresponding to each orientation angle. Effective orientation angles are selected from the preset multiple orientation angles based on the structural response diagrams corresponding to each orientation angle, and closed paths are determined in the first modal image based on each effective orientation angle, wherein each effective orientation angle corresponds to a closed path; The closed path is mapped to the second modal image after applying the current spatial transformation parameters to obtain the initial mapped closed path. The position of the initial mapped closed path is then corrected according to the second multi-directional structural response tensor to obtain the target mapped closed path. Extract a first response value sequence in the first modal image along the closed path, and extract a second response value sequence in the second modal image along the target mapping closed path; The second response value sequence is transformed to obtain the inverted response value sequence. The sequence matching degree between the first response value sequence and the inverted response value sequence is calculated, and the weighted sum of the sequence matching degrees corresponding to all valid direction angles is used as the global registration matching degree. The current spatial transformation parameters are adjusted according to the global registration fit, and the adjusted current spatial transformation parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

2. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, The step of constructing anisotropic Gaussian derivative kernels with multiple preset orientation angles, and filtering the first modal image and the second modal image according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image includes: The first direction angle of each pixel position in the first modal image along the direction of the maximum change in gray value in the neighborhood is obtained, and the first direction angles of all pixel positions are statistically analyzed according to the angle value to obtain a first direction angle histogram. The second direction angle of each pixel position in the second modal image along the direction of the maximum change in gray value in the neighborhood is obtained, and the second direction angles of all pixel positions are statistically analyzed according to the angle value to obtain a second direction angle histogram. Each local peak angle in the first orientation angle histogram is marked as a first candidate angle, and it is determined whether the distribution frequency at the same angle as a certain first candidate angle in the second orientation angle histogram is greater than a preset frequency threshold. If the distribution frequency is greater than the preset frequency threshold, then a certain first candidate angle is retained as a second candidate angle, and all first candidate angles are traversed to obtain a set of candidate direction angles composed of all second candidate angles. Perform a clustering operation on each second candidate angle in the candidate orientation angle set, and take the angle corresponding to each cluster center as a preset orientation angle to obtain a preset orientation angle set; Using each preset direction angle in the preset direction angle set as a rotation angle, rotate the two-dimensional Gaussian kernel around an axis perpendicular to the plane where the two-dimensional Gaussian kernel is located to the direction corresponding to the preset direction angle. Calculate the first-order partial derivative of the rotated two-dimensional Gaussian kernel along the direction corresponding to the preset direction angle to obtain the anisotropic Gaussian derivative kernel corresponding to each preset direction angle. The anisotropic Gaussian derivative kernels corresponding to each preset direction angle in the preset direction angle set are convolved with the first modal image to obtain the first structural response map under each preset direction angle. All the first structural response maps under the preset direction angles are arranged in the order of the preset direction angles to form the first multi-directional structural response tensor. The anisotropic Gaussian derivative kernels corresponding to each preset orientation angle in the preset orientation angle set are convolved with the second modality image to obtain the second structural response map under each preset orientation angle. All the second structural response maps under the preset orientation angles are arranged in the order of the preset orientation angles to form the second multi-directional structural response tensor.

3. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, The step of selecting effective direction angles from the preset plurality of direction angles based on the structural response map corresponding to each direction angle, and determining the closed path in the first modal image based on each effective direction angle includes: Multiply the response values ​​at the same pixel position of a first structural response map in the first multi-directional structural response tensor and a second structural response map in the second multi-directional structural response tensor, and sum the product results over all pixel positions to obtain a cross-modal response consistency metric value for a certain preset directional angle, wherein the first structural response map and the second structural response map both correspond to the preset directional angle. All preset orientation angles are sorted from largest to smallest according to the cross-modal response consistency metric value, and the orientation angles that are in the top preset number of sorts are selected as valid orientation angles; In the first structural response map corresponding to a certain effective direction angle, the pixel positions where the response value along the direction of the certain effective direction angle presents a local maximum and the distance between two adjacent local maxima is greater than a preset distance threshold are extracted. The extracted pixel positions are connected in sequence to form a structural center line under the certain effective direction angle. Taking each pixel position on the center line of a certain structure as the center, extend to both sides in a direction perpendicular to the center line of the certain structure, and find the position with the largest change rate of response value in the first structural response map corresponding to a certain effective direction angle as the boundary point to obtain the left boundary point set and the right boundary point set; The left boundary point set is connected in sequence to form the left boundary line, the right boundary point set is connected in sequence to form the right boundary line, and the left boundary line and the right boundary line are connected at both ends to form a closed path, wherein the closed path surrounds the physical structure boundary corresponding to the effective direction angle.

4. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, The step of applying the current spatial transformation parameters to the closed path and mapping it to the second modality image to obtain the initial mapped closed path includes: Obtain the current spatial transformation parameters, which include horizontal translation, vertical translation, rotation angle, uniform scaling factor, and shearing parameter. The initial value of the current spatial transformation parameters is determined by the global coarse registration result between the first modal image and the second modal image. During the fine registration iteration process, the current spatial transformation parameters are updated with each iteration. Obtain a first internal region enclosed by a certain closed path, and calculate the first area value of the first internal region in the first modal image; After applying the current spatial transformation parameters to a certain closed path, the initial position of the mapping path is obtained, and a second internal region is formed in the second modality image based on the initial position of the mapping path; Centered on the centroid of the second internal region, and with the first area value as the target area, the region is expanded radially outward layer by layer in the second modal image until the area of ​​the expanded region is equal to the first area value. The boundary of the expanded region is then used as the first candidate mapping closed path. Determine the response polarity of each pixel position on the first candidate mapping closed path in the direction of the effective direction angle corresponding to a certain closed path in the second multi-directional structure response tensor, and count the number of first pixel positions with negative response polarity and the number of second pixel positions with positive response polarity, wherein the response polarity is the sign of the response value; If the absolute value of the difference between the number of first pixel positions and the number of second pixel positions is less than a preset difference threshold, then the first candidate mapping closed path is used as the initial mapping closed path. If the absolute value of the difference between the number of first pixel positions and the number of second pixel positions is not less than the preset difference threshold, then it is determined whether the number of first pixel positions is greater than the number of second pixel positions. If yes, then the centroid is moved by a preset step along the distribution centroid direction of the pixel positions with negative response polarity on the boundary of the extended region to obtain the updated centroid position; if no, then the centroid is moved by a preset step along the distribution centroid direction of the pixel positions with positive response polarity on the boundary of the extended region to obtain the updated centroid position. Using the updated centroid position as the center and the first area value as the target area, the layer-by-layer expansion is performed again to obtain the second candidate mapping closed path. The response polarity is counted again. If the condition is still not met, the centroid is moved and expanded again until the condition is met or the preset maximum number of moves is reached. The candidate mapping closed path corresponding to the preset maximum number of moves is taken as the initial mapping closed path.

5. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, The step of correcting the position of the initial mapping closed path based on the second multi-directional structural response tensor to obtain the target mapping closed path includes: Obtain a valid direction angle corresponding to a certain initial mapping closed path, and a first response value of each pixel position on the closed path corresponding to the certain initial mapping closed path in the first multi-directional structure response tensor along the direction of the certain valid direction angle, and obtain the first response polarity of the first response value, wherein the first response polarity is the positive or negative sign of the first response value; Invert the first response polarity of each pixel position on the initial mapping closed path to obtain the expected response polarity of each pixel position. Centered on a pixel position on a certain initial mapping closed path, a first direction and a second direction search are performed along the direction perpendicular to the direction indicated by a certain effective direction angle. Within the search range, the second response value of each search position in the second multi-directional structural response tensor along the direction of the certain effective direction angle is obtained, and the second response polarity of the second response value is obtained. The search range is a first preset distance extended from a pixel position in the first direction and a second preset distance extended in the second direction, wherein the first direction and the second direction are opposite to each other. All search positions within the search range whose second response polarity matches the expected response polarity are formed into a candidate position set. The search position with the largest absolute value of the second response value is selected from the candidate position set as the optimized candidate position for a certain pixel position. Connect all the candidate positions of the optimized pixel positions in sequence according to the original order of the pixel positions on the initial mapping closed path to obtain a target mapping closed path corresponding to the initial mapping closed path.

6. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, Extracting the first response value sequence from the first modality image along the closed path includes: For a given closed path, obtain the effective direction angle corresponding to that closed path; The total number of pixels on a closed path is determined. Using the total number of pixels as the number of traversals, the Euclidean distances between adjacent pixel positions are calculated sequentially and accumulated to obtain the total length of the closed path. Starting from any point on the closed path, sampling is performed in a clockwise direction with equal step sizes to obtain multiple sampling points arranged in the path order. The step size is the total length divided by the preset number of sampling points. For each sampling point, the response value of the sampling point along the effective direction angle in the first multi-directional structural response tensor is obtained, and the response values ​​are arranged according to the path order of the sampling points to obtain the first response value sequence corresponding to a certain closed path.

7. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 1, characterized in that, in, A certain second response value sequence corresponds to a certain closed mapping path, and a certain first response value sequence corresponds to a certain closed mapping path, wherein the certain closed mapping path is obtained by mapping a certain closed path; The step of performing a sequence transformation on the second response value sequence to obtain an inverted response value sequence, calculating the sequence matching degree between the first response value sequence and the inverted response value sequence, and using the weighted sum of the sequence matching degrees corresponding to all valid direction angles as the global registration matching degree includes: Multiply each second response value in the given second response value sequence by -1 to obtain a reverse response value sequence; The first response value sequence and the inverted response value sequence are respectively normalized; Calculate a Pearson correlation coefficient between a normalized first response value sequence and a normalized inverse response value sequence, and use the Pearson correlation coefficient as the directional fit degree of the effective direction angle corresponding to the closed path. The global registration match is obtained by weighted summation of the directional match rates of all effective directional angles. The weight of each effective directional angle in the weighted summation is the proportion of the energy value of the structure response map corresponding to the effective directional angle to the total energy value of the structure response maps corresponding to all effective directional angles. The energy value is the sum of the squares of the response values ​​of each pixel position in the structure response map.

8. The dual-modal tomographic image registration method based on multi-directional structural features according to claim 4, characterized in that, The step of adjusting the current spatial transformation parameters based on the global registration fit, applying the adjusted current spatial transformation parameters to the second modal image to align the second modal image with the first modal image, and outputting the registration result includes: With the goal of maximizing the global registration matching accuracy, the Powell algorithm is used to iteratively update the current spatial transformation parameters. In each iteration step, after updating the current spatial transformation parameters, the updated global registration fit is redefined; If the updated global registration match is greater than the original global registration match, then the updated current spatial transformation parameters are accepted. If the updated global registration match is not greater than the original global registration match, then the updated current spatial transformation parameters are rejected, and the original current spatial transformation parameters remain unchanged. When the difference in global registration fit between two adjacent iterations is less than a preset convergence threshold, the iteration is terminated, and the current spatial transformation parameters at the time of termination are used as the final registration parameters. The final registration parameters are applied to the second modal image to align the second modal image with the first modal image, and the registration result is output.

9. A dual-modal tomographic image registration system based on multi-directional structural features, characterized in that, include: The acquisition module is configured to acquire a first modal image and a second modal image of the same sample, wherein the first modal image and the second modal image have opposite grayscale responses to the same physical structure; The filtering module is configured to construct anisotropic Gaussian derivative kernels at multiple preset direction angles, and filter the first modal image and the second modal image according to the anisotropic Gaussian derivative kernels to obtain a first multi-directional structural response tensor corresponding to the first modal image and a second multi-directional structural response tensor corresponding to the second modal image, wherein the multi-directional structural response tensor contains the structural response map corresponding to each direction angle; The determination module is configured to select effective direction angles from the preset multiple direction angles based on the structural response map corresponding to each direction angle, and determine a closed path in the first modal image based on each effective direction angle, wherein each effective direction angle corresponds to a closed path; The mapping module is configured to apply the current spatial transformation parameters to the closed path and map it to the second modal image to obtain an initial mapped closed path, and to perform position correction on the initial mapped closed path according to the second multi-directional structural response tensor to obtain a target mapped closed path; The extraction module is configured to extract a first response value sequence in the first modal image along the closed path, and to extract a second response value sequence in the second modal image along the target mapping closed path; The calculation module is configured to perform a sequence transformation on the second response value sequence to obtain an inverted response value sequence, calculate the sequence matching degree between the first response value sequence and the inverted response value sequence, and use the weighted sum of the sequence matching degrees corresponding to all effective direction angles as the global registration matching degree. The output module is configured to adjust the current spatial transformation parameters according to the global registration fit, apply the adjusted current spatial transformation parameters to the second modal image, align the second modal image with the first modal image, and output the registration result.

Images