Rock core image registration method based on multi-modal image segmentation and feature matching

Abstract

The invention provides a core image registration method based on multi-modal image segmentation and feature matching. The method is characterized in that the imaging difference between a rock core laser confocal image and a CT image is reduced through preprocessing color conversion and a normalization strategy, a marker is segmented by using deep learning, an algorithm target area is optimized in combination with morphological operation, feature points are detected under the guidance of an ROI area by using an ORB algorithm, and a target image is obtained. Matching of feature points is achieved by combining an FLANN feature point matching algorithm and a random K-d tree acceleration algorithm, image alignment is completed through a rigid transformation model, and registration precision is quantified through Dice. According to the method provided by the invention, the registration problem of the multi-modal core image caused by imaging difference is solved, the average Dice reaches 0.9210, and the registration efficiency and accuracy are remarkably improved.

Description

Technical Field

[0001] This disclosure relates to the field of core image processing technology, and more specifically, to a core image processing method based on multimodal image segmentation and feature matching. Background Technology

[0002] High-precision registration of core images is a key technology in geological exploration and hydrocarbon reservoir analysis, supporting pore structure and 3D reconstruction by integrating multimodal data such as laser confocal imaging and CT scans. However, existing core image registration methods are not perfect and currently suffer from the following technical shortcomings that need to be addressed:

[0003] Resolution mismatch: CT images have a large field of view (>10cm) 3 However, the resolution is low (micrometer level), while laser image resolution reaches the nanometer level but the field of view is small (<1mm). 3 This leads to the failure of cross-modal pixel-level mapping, as shown in the literature: Schirra Jr, Randall T., and Peijun Zhang. "Correlative fluorescence and electron microscopy." Current protocols in cytometry 70.1(2014):12-36;

[0004] Artifact interference: CT metal artifacts cause HU value shifts >300 units, and laser images produce brightness tomography due to fluorescence attenuation (gradient change >50%), as shown in the literature: Puvanasunthararajah, Sathyathas et al. “The application of metal artifact reduction methods on computed tomography scans for radiotherapy applications: A literature review.” Journal of applied clinical medical physics vol.22,6(2021):198-223.

[0005] Uneven feature distribution: The coverage of matching areas such as cracks / cement is less than 30% of the entire image, and there is nonlinear distortion in the morphology (curvature change >15°).

[0006] The aforementioned defects lead to the following verifiable technical bottlenecks:

[0007] Gray-level inconsistency: The statistical correlation between the optical reflectance intensity of laser images (range 0-255) and the material density of CT (HU value range -1000 to 3000) is <0.3, as shown in the literature: Zitova, Barbara, and Jan Flusser. "Image registration methods: a survey." Image and vision computing 21.11(2003):977-1000.

[0008] Feature extraction failure: The morphological and textural features of markers (such as artificial markers or mineral particles) in multimodal images exhibit nonlinear distortions. Conventional feature extraction algorithms (such as SIFT and SURF) are sensitive to noise and have low computational efficiency, with a false matching rate exceeding 60%, as shown in the literature: Lowe, David G. "Distinctive image features from scale-invariant keypoints." International journal of computer vision 60 (2004): 91-110.

[0009] Insufficient process coupling: In the existing registration process, there is little overlap in information between the preprocessing and feature extraction stages, and the texture style differences are large, making registration difficult. This leads to problems such as blurred edges and uneven contrast, which affects the robustness of subsequent transformation matrix estimation.

[0010] The defects in the above image registration methods directly lead to the following problems:

[0011] The preprocessing method is too simplistic: traditional grayscale normalization only adjusts the brightness range by linear stretching, without considering the differences in noise types in multimodal images (CT radiation noise variance > 100, laser optical noise variance > 50), resulting in a decrease in the segmentation Dice coefficient of more than 0.2, as shown in the literature: Gu,njia,Yimin Chen, and Weiqin Tong. "KARAN: Mitigating Feature Heterogeneity and Noise for Efficient and Accurate Multimodal Medical Image Segmentation." Electronics 13.23(2024):4594;

[0012] Low feature matching efficiency: Although deep learning-based semantic segmentation can improve the accuracy of marker localization, its output binary mask still needs to rely on morphological post-processing to extract stable edges. However, the traditional dilatation-erosion operation parameters are fixed, making it difficult to handle the structural diversity of different core samples, as shown in the literature: Fang, Chao, **aopeng Wang, and Qingsheng Wang. "Adaptive morphology structural element construction algorithm based on local pixel density and symmetry." Multimedia Tools and Applications 82.1(2023):195-215;

[0013] The matching algorithm lacks adaptability: Existing methods (such as affine transformation based on RANSAC) are easily affected by local optima when dealing with rigid registration, and they do not combine multi-scale feature descriptors (such as binary encoding) to accelerate matching, resulting in a sharp increase in the time consumption of large-scale core data registration.

[0014] Registration evaluation lacks quantitative indicators: Most studies rely on manual visual inspection or mean squared error (MSE) to evaluate registration results, which cannot accurately quantify the overlap of marker regions, thus restricting the automated optimization of algorithm parameters, as shown in the reference: Rohlfing, Torsten. "Image similarity and tissue overlaps as surrogates for image registration accuracy: widely used but unreliable." IEEE transactions on medical imaging 31.2(2011):153-163.

[0015] To address the aforementioned issues, a registration method integrating multimodal image enhancement, adaptive morphological processing, and efficient feature matching is needed to improve the accuracy and efficiency of cross-modal core image registration and provide a reliable data foundation for geological modeling. Summary of the Invention

[0016] To address the aforementioned technical challenges, this disclosure proposes a core image registration method based on multimodal image segmentation and feature matching. This registration method integrates multimodal image enhancement, adaptive morphological processing, and efficient feature matching, thereby improving the accuracy and efficiency of cross-modal core image registration. It provides a reliable data foundation for geological modeling and solves problems such as the reliance on a single preprocessing method, low feature matching efficiency, insufficient algorithm adaptability, and lack of quantitative evaluation metrics.

[0017] The technical solution and implementation steps adopted in this disclosure are as follows:

[0018] First, a core image registration method based on multimodal image segmentation and feature matching is presented, including the following steps:

[0019] S1: Preprocess the laser confocal image and CT slice image, including grayscale normalization, noise removal and contrast enhancement, to unify the grayscale distribution and resolution of the two modes, and obtain laser confocal image A and CT slice image B.

[0020] S2: Train a deep learning semantic segmentation model using a core image dataset containing pixel-level annotations of markers in laser confocal image A and slice image B. The trained deep learning semantic segmentation model is then used to segment the preprocessed images, generating a binary mask file that includes laser confocal image marker segmentation mask C and CT image marker segmentation mask D.

[0021] S3: Morphological processing is performed on the laser confocal image marker segmentation mask C and the CT image marker segmentation mask D generated by S2, including dilation, erosion and difference operations, to obtain the accurate marker edge region;

[0022] S4: Based on the precise marker edge feature region obtained in S3, key points are detected and key feature points are generated using the ORB algorithm;

[0023] S5: Using laser confocal image A as a fixed image and CT slice image B as a floating image, the key feature points obtained in S4 are matched using the FLANN matching algorithm and combined with a random Kd tree to accelerate the search, thus obtaining the corresponding feature point pairs between laser confocal image A and CT slice image B.

[0024] S6: Based on the corresponding feature point pairs obtained in S5, the RANSAC algorithm is applied to remove mismatched point pairs. In addition, the reliable matching point pairs selected by RANSAC are used to calculate the rigid transformation matrix of the floating image B to the fixed image A, and the rigid transformation matrix is ​​obtained.

[0025] S7: Based on the rigid transformation matrix obtained in S6, perform a spatial transformation on the floating image to align it with the fixed image. Simultaneously, transform the segmentation file and perform the same spatial transformation on the CT image corresponding marker segmentation mask file D to obtain the transformed CT slice image marker segmentation mask E, maintaining consistency with the transformed image, and obtain the best CT slice image after registration.

[0026] S8: The best CT slice image obtained after registration based on S7 is spatially aligned and saved with the original laser confocal image A according to the same coordinates. At the same time, the corresponding marker mask files are saved, which contain the marker segmentation mask E of the laser confocal image A and the CT slice image. The parameters of the final transformation matrix used for the slice registration are recorded, thus completing the generation of the core registration result.

[0027] S9: Finally, the registration quality is evaluated by calculating the Dice coefficient between the laser confocal image marker segmentation mask C and the registered CT slice image marker segmentation mask E. The CT slice with the highest Dice coefficient is selected as the final registration result of the core sample.

[0028] Furthermore, the preprocessing in step S1 includes:

[0029] Gray-level normalization: unifies the gray-level range of laser confocal images and CT images to [0, 255];

[0030] Noise removal: Gaussian filtering (σ = 1.0) was used for laser confocal images, and non-local means denoising (NL-Means) was used for CT images;

[0031] Contrast enhancement: Histogram equalization (CLAHE algorithm, 8×8 grid size);

[0032] Resolution unification: The laser confocal image is downsampled to the CT slice resolution using bicubic interpolation;

[0033] Format conversion: Convert CLSM image to 8-bit grayscale image;

[0034] Output: Preprocessed laser confocal image A and CT slice image B.

[0035] Further, step S2 includes:

[0036] Data augmentation operations were performed on the preprocessed laser confocal images and CT images, including random cropping, rotation, flipping, and brightness adjustment, to enhance the model's ability to generalize to complex textures in core images.

[0037] A semantic segmentation model based on selective state space is adopted. The model includes a texture feature enhancement module, a feature selection module, and a visual Mamba layer. The texture feature enhancement module extracts microstructural features such as particles, pores, and cracks from core images through multi-scale convolution operations. The feature selection module selects key features related to markers based on an attention mechanism to suppress the influence of background noise. The visual Mamba layer uses a selective state space mechanism to efficiently model long sequence features and generate high-precision pixel-level segmentation results.

[0038] The model was trained using a labeled core image dataset containing pixel-level annotations of markers in laser confocal images and CT images. The model parameters were optimized using the cross-entropy loss function and the Adam optimizer during training.

[0039] The output of the semantic segmentation model is a pixel-level probability map of the same size as the input image. The value of each pixel in the probability map represents the predicted probability that the pixel belongs to the category of the marker.

[0040] The probability map output by the model is post-processed by setting a threshold to convert the probability map into a binary mask file, which contains the laser confocal image marker segmentation mask C and the CT image marker segmentation mask D.

[0041] Further, step S3 includes:

[0042] The dilation operation enhances the connectivity of a target object in an image by expanding its boundaries. It is used to expand the target region in an image, preparing it for subsequent erosion and edge extraction. The specific mathematical formula is as follows:

[0043]

[0044] Where A is a binary image, B is a structuring element, and B z This indicates that the structuring element B is translated to position Z. This represents the erosion operator. This represents the inflated image; when the area covered by structuring element B intersects with the foreground region of A, the center pixel Z is marked as a foreground pixel;

[0045] Etching reduces noise and stray areas by narrowing the boundary of the target region, as shown in the formula:

[0046]

[0047] in, This means that when all pixels of structuring element B are completely contained within the foreground region of A, the center pixel Z is retained as a foreground pixel. This represents the image after erosion; through the erosion operation, the boundaries of the target object are shrunk, effectively reducing small stray objects.

[0048] The difference operation between dilation and erosion is used to extract the edge region of the target object. It can clearly extract the edge texture features of the marker, thereby more accurately locating the boundary of the ROI region. The formula is:

[0049] Edge Roi =A dilated -A eroded

[0050] Among them, A dilated Indicating the result of the expansion, A eroded Indicates the corrosion result, Edge Roi This indicates the edge area of ​​the sign.

[0051] Furthermore, in step S4, the ORB algorithm is used to detect key points and generate key feature points, following the path below:

[0052] Within the ROI area, the FAST algorithm is used to detect key points. The criteria for determining a point as a corner are based on the continuity of pixel grayscale changes within a certain range around the candidate point. When the grayscale values ​​of multiple adjacent pixels are simultaneously higher or lower than the grayscale value of the center pixel, the point can be considered a corner.

[0053] For the detected key points, the BRIEF algorithm is used to generate feature descriptors. This algorithm selects pixel pairs in the local neighborhood and compares their grayscale relationships to construct descriptive information in binary form. This results in binary feature descriptors generated from the key points, thus realizing the detection of key points and generation of key feature points using the ORB algorithm.

[0054] Furthermore, in step S5, using the laser confocal image A as the fixed image and the CT slice image B as the floating image, the matching is performed using the FLANN feature matching algorithm, including:

[0055] Calculate feature point similarity based on Euclidean distance from feature descriptors;

[0056] Calculate the ratio of the nearest neighbor distance to the second nearest neighbor distance for each feature point;

[0057] When the ratio is less than a set threshold (range 0.5-0.8), the nearest neighbor matching point is accepted as a reliable matching point;

[0058] A random Kd-tree data structure is used to accelerate the nearest neighbor search of feature points;

[0059] Match the corresponding feature point pairs between laser confocal image A and CT slice image B.

[0060] Furthermore, the rigid transformation matrix obtained in step S6 is executed according to the following path:

[0061] Transformation parameter calculation: Based on reliable matching point pairs selected by RANSAC, the optimal transformation parameters are solved using the least squares method.

[0062]

[0063] Where N is the number of interior points, and pCLSM and pCT are the coordinates of the matching points in the fixed image and the floating image, respectively;

[0064] Based on the obtained reliable matching point pairs, a similarity transformation model is adopted, and the specific transformation matrix formula is as follows:

[0065]

[0066] In the transformation matrix, S is the scaling factor, θ is the rotation angle, and t x and t y These are the translation amounts; this transformation matrix defines how to align a floating image to a fixed image through rotation, scaling, and translation operations, with laser confocal image A as the fixed image and CT slice image B as the floating image.

[0067] Furthermore, step S7 also includes:

[0068] Spatial transformation is performed: The floating image is spatially transformed using bilinear interpolation according to the following formula to ensure precise alignment with its fixed image:

[0069]

[0070] Where (x′, y′) are the transformed coordinates, w ij These are bilinear weighting coefficients;

[0071] The identifier segmentation mask D is synchronously transformed using nearest neighbor interpolation:

[0072] F(x′,y′)=D(round(x),round(y))

[0073] Finally, the transformed CT slice image identifier segmentation mask E was obtained, maintaining consistency with the transformed image.

[0074] Furthermore, in step S8, the generation of core registration results also includes:

[0075] Registration data storage:

[0076] Save the registered best CT slice and the original laser confocal image in a unified coordinate system, and store them in the lossless compressed format TIFF or PNG format. At the same time, save the corresponding marker mask file, which includes the marker segmentation mask E of the laser confocal image A and the CT slice image.

[0077] Transformation parameter record:

[0078] The parameters of the rigid transformation matrix are stored in a machine-readable text format.

[0079] Furthermore, the registration quality assessment in step S9 includes:

[0080] Dice coefficient calculation: Dice is an index that measures the similarity between two sets. The calculation formula is as follows:

[0081]

[0082] Where A and B represent the segmentation results of laser confocal map markers and CT image markers, respectively; A∩B represents the overlapping area of ​​the two; and |A| and |B| represent the number of pixels in the mask A and B regions, respectively.

[0083] Optimal slice selection rule: Iterate through the Dice coefficients of all CT slices; select the slice that satisfies Dice ≥ 0.90 and takes the maximum value as the final registration result;

[0084] If all slices have a Dice < 0.90, then initiate iterative optimization: readjust the FLANN matching threshold in step S5 between 0.5 and 0.8 and adjust the RANSAC tolerance in step S6 between 2.0 and 5.0 pixels. Then, repeat steps S5 to S9.

[0085] Compared with the prior art, this disclosure has the following advantages and technical effects:

[0086] This invention provides a core image registration method based on multimodal image segmentation and feature matching. It addresses the problems of large differences in imaging mechanisms, insufficient preprocessing, poor robustness of feature matching, and inaccurate registration evaluation in existing technologies. It proposes a complete set of executable and reproducible image registration procedures and has achieved clear results in experimental verification.

[0087] First, addressing the issue of inconsistent grayscale statistical distributions between laser confocal images and CT images, making direct grayscale registration impossible using traditional mutual information or SIFT algorithms, this invention introduces a deep learning semantic segmentation model to automatically identify common marker regions in both images. This bypasses direct alignment at the grayscale level, instead using region structure alignment as the registration basis, thereby improving the feasibility and stability of cross-modal registration. This method is particularly suitable for situations where the physical mechanisms of the two imaging modalities differ significantly and their grayscale information is uncorrelated.

[0088] Secondly, to address the issues of low image contrast, blurred boundaries, and high noise interference caused by insufficient consideration of modal differences in the image preprocessing stage, this invention introduces a strategy combining grayscale normalization, median filtering, and Gaussian blurring in the preprocessing stage. The combination of Gaussian blurring (σ = 1.5) and median filtering improves the PSNR of CT images by >8dB and the SSIM of laser images by >0.15. Simultaneously, morphological dilation, erosion, and differential operations enhance the edge regions of the segmentation mask, making subsequent feature point extraction more concentrated in clearly defined effective structural regions, thus avoiding matching failures caused by uneven feature distribution across the entire image.

[0089] In terms of feature matching, addressing the issues of high computational cost, noise sensitivity, and low matching accuracy of traditional feature algorithms (such as SIFT and SURF) in core images, this invention employs a lightweight ORB algorithm to extract key points in the ROI region, combines it with BRIEF to generate binary feature descriptors, and then uses a FLANN matcher with ratio filtering and Kd-tree acceleration mechanisms to achieve fast and accurate feature point matching. This strategy improves feature matching efficiency and reduces the risk of mismatches. Experimental results show that this feature point matching scheme can complete the registration and matching stage of a set of images in less than 1.5 seconds on average. The ORB+FLANN combination takes less than 1.5 seconds to match in images, which is 3 times faster than SIFT (the measured average SIFT time is 4.8 seconds).

[0090] In the image transformation and evaluation stage, to address the lack of structural quantitative evaluation standards after registration and the inadequacy of relying solely on manual visual inspection or simple mean square error to reflect regional consistency, this invention employs Dice to quantitatively evaluate the overlap of marker regions in the registration results. Experimental data show that in core sample testing, this method achieves a Dice value of 0.9210 without manual intervention, and the Dice coefficient of 0.921±0.042 is significantly higher than that of traditional methods (mutual information method Dice = 0.712±0.098), making it particularly suitable for scenarios such as pore structure analysis and core digital modeling.

[0091] In summary, this invention constructs an automated multimodal registration system suitable for core images by introducing modal adaptive preprocessing, edge structure enhancement feature extraction, a lightweight matching strategy guided by structural regions, and a quantitative evaluation mechanism for structural overlap rate. This system significantly reduces the need for manual intervention while ensuring registration accuracy, making it suitable for large-scale batch processing of core images and providing a stable data registration foundation for geological modeling, reservoir analysis, and rock property research. Attached Figure Description

[0092] The accompanying drawings, which form part of this disclosure, are used to provide a further understanding of this disclosure. The illustrative embodiments of this disclosure and their descriptions are used to explain this disclosure and do not constitute an undue limitation of this disclosure. In the drawings:

[0093] Figure 1 This is an overall flowchart of a specific embodiment of the present disclosure.

[0094] Figure 2 This is a flowchart illustrating the assisted registration process based on segmentation results in a specific embodiment of this disclosure.

[0095] Figure 3 This is a schematic diagram illustrating the ORB feature point detection principle in a specific embodiment of this disclosure.

[0096] Figure 4 This is a schematic diagram illustrating the principle of laser confocal ORB feature point detection in a specific embodiment of this disclosure.

[0097] Figure 5 This is a schematic diagram illustrating the ORB feature point detection principle in a specific embodiment of this disclosure for CT images. Detailed Implementation

[0098] The implementation process of the technical solution will be described in detail below with reference to the embodiments.

[0099] A core image registration method based on multimodal image segmentation and feature matching specifically includes the following steps:

[0100] S1: Preprocess the laser confocal image and CT slice image, including grayscale normalization, noise removal and contrast enhancement, to unify the grayscale distribution and resolution of the two modes, and obtain laser confocal image A and CT slice image B.

[0101] S2: Train a deep learning semantic segmentation model using a core image dataset containing pixel-level annotations of markers in laser confocal image A and slice image B. The trained deep learning semantic segmentation model is then used to segment the preprocessed images, generating a binary mask file that includes laser confocal image marker segmentation mask C and CT image marker segmentation mask D.

[0102] S3: Morphological processing is performed on the laser confocal image marker segmentation mask C and the CT image marker segmentation mask D generated by S2, including dilation, erosion and difference operations, to obtain the accurate marker edge region;

[0103] S4: Based on the precise marker edge feature region obtained in S3, key points are detected and key feature points are generated using the ORB algorithm;

[0104] S5: Using laser confocal image A as a fixed image and CT slice image B as a floating image, the key feature points extracted by ORB based on S4 are matched by the FLANN matching algorithm and combined with random Kd tree to accelerate the search, so as to obtain the corresponding feature point pairs between laser confocal image A and CT slice image B.

[0105] S6: Based on the feature point matching pairs obtained in S5, the RANSAC algorithm is applied to remove mismatched point pairs. In addition, the reliable matching point pairs selected by RANSAC are used to calculate the rigid transformation matrix of the floating image B to the fixed image A, and the rigid transformation matrix is ​​obtained.

[0106] S7: Based on the rigid transformation matrix obtained in S6, perform spatial transformation on the floating image B to make it precisely aligned with the fixed image A. At the same time, transform the segmentation file and perform the same spatial transformation on the CT image corresponding marker segmentation mask file D to obtain the transformed CT slice image marker segmentation mask E, maintaining consistency with the transformed image.

[0107] S8: Based on S7, the best CT slice after registration is spatially aligned and saved with the original laser confocal image A according to the same coordinates. At the same time, the corresponding marker mask files are saved, which contain the marker segmentation mask E of the laser confocal image A and the CT slice image. The parameters of the final transformation matrix used for the slice registration are recorded, thus completing the generation of the core registration result.

[0108] S9: Finally, the registration quality is evaluated by calculating the Dice coefficient between the laser confocal image marker segmentation mask C and the registered CT slice image marker segmentation mask E. The CT slice with the highest Dice coefficient is selected as the final registration result of the core sample.

[0109] Further specifying, the preprocessing of the laser confocal image and CT slice image is a fundamental operation of this invention. Its purpose is to unify the distribution of grayscale, noise features, and resolution in the multimodal images, providing consistent input for subsequent segmentation and registration. Specifically, the image preprocessing includes the following sub-operations: the grayscale normalization step unifies the pixel value range of the laser confocal image and CT image to [0, 255]; the RGB image of the laser confocal image is converted to a grayscale image, as shown in the following formula.

[0110] I=0.2989·R+0.5870·G+0.1140·B

[0111] Where I represents the pixel intensity value of the output grayscale image, and R, G, B represent the pixel values ​​of each channel in the original RGB image, typically ranging from [0, 255], to remove color interference. In the noise removal process, the CLSM image uses Gaussian filtering (σ = 1.0), while the CT image uses the Non-Local Means (NL-Means) algorithm. The contrast enhancement operation uses CLAHE (8×8 grid) to improve local image contrast. Resolution unification involves scaling the CLSM image to the CT image size using bicubic interpolation and converting the CLSM image to an 8-bit grayscale image through format conversion. After this process, a uniformly formatted laser confocal image A and a CT slice image B are output as input images for subsequent segmentation.

[0112] Further specifying, the step involves introducing a novel semantic segmentation network based on selective state space to perform pixel-level structural recognition on both the laser confocal image A and the CT image B. This network uses a classic U-structure as its backbone, superimposed with a texture feature enhancement module, an attention feature selection module, and a multi-layer visual Mamba sequence modeling structure, enabling effective identification of structures such as pores, matrix, cracks, and cement in core images.

[0113] Specifically: the texture enhancement module performs boundary detection based on first-order and second-order gradient information, and captures small features by combining shallow texture convolution; the feature selection module activates important channels through Tanh nonlinear transformation and introduces residual skip connections to solve the gradient vanishing problem; the visual Mamba module models high-resolution long sequence features through selective state space, taking into account both global perception capability and efficient inference.

[0114] To improve the robustness of the network on multimodal datasets, this invention employs various data augmentation strategies, including rotation, translation, brightness perturbation, and random cropping, to amplify training samples, with some of these augmentations utilizing the MedAugment scheme. The final output consists of two binary masks: a laser image marker mask C and a CT image mask D, used to accurately identify regions of interest and provide ROI support for subsequent keypoint detection.

[0115] Further defining the structural boundaries of the markers extracted from the segmentation masks C and D ensures that subsequent key point extraction focuses on the true edges of the physical structure, improving the accuracy of feature point matching. The initial morphological processing involves a series of operations based on image structural changes, mainly including dilation, erosion, and difference operations. These operations play a crucial role in edge extraction and feature matching. The dilation operation enhances the connectivity of the target object by expanding its boundary in the image; the specific mathematical formula is as follows:

[0116]

[0117] Where A is a binary image, B is a structuring element, and B z This indicates that the structuring element B is translated to position Z. This represents the erosion operator. This represents the expanded image. When the area covered by structuring element B intersects with the foreground region of A, the center pixel Z is marked as a foreground pixel.

[0118] Dilation helps expand the target region in an image, preparing it for subsequent erosion and edge extraction. Erosion, on the other hand, reduces noise and clutter by shrinking the boundaries of the target region; the formula is:

[0119]

[0120] in, This means that when all pixels of structuring element B are completely contained within the foreground region of A, the center pixel Z is retained as a foreground pixel. This represents the image after erosion. Through the erosion operation, the boundaries of the target object are shrunk, effectively reducing small stray objects. The difference operation between dilation and erosion is used to extract the edge regions of the target object; its formula is:

[0121] Edge Roi =A dilated -A eroded

[0122] Among them, A dilated Indicating the result of the expansion, A eroded Indicates the corrosion result, Edge Roi This method represents the edge region of the marker. Through this differential operation, the edge texture features of the marker can be clearly extracted, thereby more accurately locating the boundary of the ROI region. This method significantly reduces the interference of background noise, providing higher accuracy and robustness for subsequent feature extraction and matching. This processing can accurately locate the spatial boundaries of important structures such as cement, pores, and cracks, ensuring that registration keypoints do not fall into the background or unstructured areas, improving the overall matching accuracy and stability, and thus obtaining accurate marker edge regions.

[0123] Further specifying, the method of using the ORB algorithm to detect key points and generate key feature points within the feature area of ​​the precise marker includes two processes: FAST corner detection: determining whether it is a corner point based on the continuous jump of local gray level in the image, and quickly locating areas with drastic structural changes; BRIEF feature description: generating a binary feature vector for each corner point, and compressing and encoding the gray level relationship through pixels to improve the matching speed.

[0124] Further specifying, the FAST corner detection involves using the FAST algorithm to detect key points within the ROI region. The criteria for identifying a point as a corner are based on the continuity of pixel grayscale changes within a certain range around the candidate point. When the grayscale values ​​of multiple adjacent pixels are simultaneously higher or lower than the grayscale value of the center pixel, the point is considered a corner. In the edge region of the cement, the grayscale changes are more pronounced, thus it is identified as a feature point. This method is computationally efficient and can quickly detect key points in the image. In the feature description stage, the BRIEF algorithm is used to generate binary descriptors for the extracted key points, as shown in the following formula:

[0125]

[0126] p: The position of the currently detected image pixel; Ip: The gray value of pixel p; xi: The i-th pixel in the circular neighborhood; Ixi: The gray value of the neighboring pixel xi; t: The gray difference threshold, which determines the corner sensitivity; n: The number of pixels that continuously meet the condition (usually 12 or higher to increase stability).

[0127] To enhance rotation invariance, ORB performs orientation normalization by calculating the principal direction for each corner point. The principal direction is calculated based on the gray-level centroid of the corner point's neighborhood image and is defined as follows:

[0128]

[0129] (xp,yp): Coordinates of corner point p; (xi,yi): Coordinates of pixels in the neighborhood of corner point; I(xi,yi): Gray value of pixel (xi,yi); θ: Main direction of corner point, used for subsequent rotation normalization processing of BRIEF descriptor.

[0130] After completing corner detection and orientation estimation, ORB uses the BRIEF algorithm to generate a compact binary feature vector for each corner. BRIEF constructs features by comparing the grayscale values ​​of a predefined set of pixel pairs (pi, qi) within the corner's neighborhood, with each comparison result being a one-bit binary value. Its basic formula is:

[0131]

[0132] p is the current corner point, pi is the pixel pair selected in the neighborhood of p (the relative position coordinates can be predefined or randomly sampled), I(pi) is the gray value of pixel point pi and qi, and τ is a binary bit of the BRIEF descriptor.

[0133] The results of all τ are concatenated to form an N-dimensional binary feature vector, typically 256 bits in length. ORB performs orientation compensation on the original BRIEF by rotating the pixel pair (pi,qi) according to the principal direction θ of the corner point, making the feature description invariant to rotation. This step can quickly obtain the initial feature point matching pair between the laser confocal image A and the CT image B, providing support for rigid transformation estimation.

[0134] Further specifying, the obtained feature point matching pairs are processed using the RANSAC algorithm to remove mismatched point pairs, first performing FLANN-based feature point matching:

[0135] For each ORB feature point extracted from laser image A and CT image B, construct a 256-dimensional binary feature descriptor;

[0136] Use the FLANN matching algorithm to compare the Hamming distance between two sets of feature descriptors;

[0137] For each laser image feature point, find its nearest and second nearest neighbors in the CT image, and calculate their descriptor distances d1 and d2.

[0138] The ratio test method is used to eliminate unstable matching pairs:

[0139]

[0140] The threshold τ = 0.75. When the ratio is less than the threshold, the match is considered reliable and the pair is retained; otherwise, it is discarded.

[0141] Furthermore, it is optimized to accelerate the search by using a random Kd-tree:

[0142] To improve matching efficiency, the feature vector index is constructed using a randomized Kd-tree structure from FLANN:

[0143] Multiple random Kd-trees are constructed in parallel, with each tree using a different random subset;

[0144] When querying, you can specify the maximum number of backtracks and the maximum number of checkpoints to balance matching speed and accuracy;

[0145] In this embodiment, the number of trees is set to 4, the maximum number of search checks is 64, and the average matching time is approximately 1.2 seconds per sample.

[0146] Finally, RANSAC is used to remove mismatches:

[0147] There may be mismatches in the set of matching point pairs. To remove outlier point pairs, the RANSAC algorithm is used to estimate the initial rigid transformation model, as follows:

[0148] Randomly sample the smallest subset (2 pairs of points) from the set of matching points;

[0149] Calculate the affine rigid transformation model based on this subset;

[0150] Substitute all matching points into the model and calculate its projection error (Euclidean distance):

[0151] ε i =‖P i CT -T(P i CLSM )‖

[0152] Among them, P i CLSM P represents the coordinates of key points in the laser image. i CT Here are the coordinates of the matching points on the CT image, T is the transformation function, and ε is the coordinates of the matching points. iIf the number of pixels is less than 2, it is included in the inner point set. Steps 1 to 4 are repeated for a total of 1000 iterations. The model with the most inner points is selected as the final model. This yields the reliable matching points selected by ANSAC.

[0153] Further specifying, the reliable matching point pairs selected by RANSAC are used to calculate the rigid transformation matrix. A rigid transformation model is employed to calculate the transformation matrix between images. This model includes translation, rotation, and scaling operations, enabling spatial transformation of images while preserving their structural integrity. This transformation effectively aligns identical regions within the image while avoiding any distortion to the image structure. The transformation matrix is ​​calculated based on the coordinates of the selected matching point pairs, using the least squares method to solve for the optimal transformation parameters. The specific transformation matrix formula is as follows:

[0154]

[0155] In the transformation matrix, S is the scaling factor, θ is the rotation angle, and t x and t y These are the translation amounts. This transformation matrix defines how to align a floating image (CT slice) to a fixed image (laser confocal image) through rotation, scaling, and translation operations. By calculating the transformation matrix, the floating image can be spatially transformed to precisely align it with the fixed image, thus obtaining the rigid transformation matrix.

[0156] Further, the rigid transformation matrix is ​​used to spatially transform the floating image B so that it is precisely aligned with the fixed image A, and the CT image mask file is updated synchronously.

[0157] Spatial transformation execution: Bilinear interpolation is used to perform a spatial transformation on the floating image (CT slice image) to ensure precise alignment with its fixed image A.

[0158]

[0159] Where Baligned(x′,y′) is the pixel value of the transformed image at (x′,y′), where (x′,y′) are the transformed coordinates, and B(x′,y′) is the pixel value of the transformed image at (x′,y′). i ,y j ) is the floating image B in (x i ,y j The pixel value of w ij These are the bilinear interpolation weighting coefficients;

[0160] The identifier segmentation mask D is synchronously transformed using nearest neighbor interpolation:

[0161] F(x′,y′)=D(round(x),round(y))

[0162] Where F(x′,y′) represents the transformed mask, D(x,y) represents the pixel value of the segmented mask, and round(x) and round(y) represent rounding the floating-point coordinates to the nearest integer coordinates;

[0163] Finally, the transformed CT slice image identifier segmentation mask E was obtained, maintaining consistency with the transformed image.

[0164] Further specifying, the registered optimal CT slice image and the laser confocal image A are spatially aligned and saved using a unified coordinate system, stored in a lossless compressed format (TIFF or PNG). Simultaneously, two corresponding marker mask files are saved, containing the marker segmentation mask E for the laser confocal image A and the CT slice image. Transformation parameter recording: rigid transformation matrix parameters are stored in a machine-readable text format.

[0165] Further refine the process by evaluating the quality of the registration results and calculating the Dice coefficient between the laser confocal image marker segmentation mask C and the registered CT image marker segmentation mask E.

[0166] Dice coefficient calculation: Dice is an index that measures the similarity between two sets. The calculation formula is as follows:

[0167]

[0168] Here, A and B represent the segmentation results of the laser confocal image markers and the CT image markers, respectively. A∩B represents the overlapping area between the two, and |A| and |B| represent the number of pixels in the mask regions A and B, respectively. The value of Dice ranges from 0 to 1. A value close to 1 indicates high registration accuracy between the two images, meaning the marker regions highly overlap; a value close to 0 indicates poor registration with significant deviation. This metric effectively evaluates the quality of image registration and can be used to optimize registration algorithms to improve registration accuracy and stability.

[0169] Optimal slice selection rule: Iterate through the Dice coefficients of all CT slices; select the slice that satisfies Dice ≥ 0.90 and takes the maximum value as the final registration result;

[0170] If all slices have a Dice < 0.90, then initiate iterative optimization: readjust the FLANN matching threshold in step S5 between 0.5 and 0.8 and adjust the RANSAC tolerance in step S6 between 2.0 and 5.0 pixels. Then, repeat steps S5 to S9.

[0171] Example 1:

[0172] The samples used in this embodiment were of diverse sizes, and 26 representative core samples were selected, including 16 conventional sandstone samples, 10 dense sandstone samples, and 6 high-porosity, high-permeability sandstone samples (porosity greater than 20%). The experiment used a nanoVoxel-3000CT instrument, which employs a microfocus X-ray source and can achieve micrometer-level spatial resolution, with a physical resolution between 4 and 20 micrometers. Laser confocal imaging was performed using a Leica Thunder DMi8 imaging instrument, scanning within a single 20×20 pixel area to acquire high-resolution, large-field-of-view microscopic images of approximately 15000×15000 pixels.

[0173] like Figure 1 As shown, this embodiment provides a core image registration method system based on multimodal image segmentation and feature matching, the process of which includes the following steps:

[0174] Image preprocessing module (S1)

[0175] Segmentation and mask generation module (S2)

[0176] Mask edge extraction module (S3)

[0177] ORB Feature Extraction Module (S4)

[0178] Feature point matching (S5)

[0179] RANSAC Transform Estimation Module (S6)

[0180] Image spatial transformation module (S7)

[0181] Registration result output module (S8)

[0182] Accuracy assessment module (S9)

[0183] The specific images and masks upon which this disclosure relies during execution include:

[0184] Laser confocal image A;

[0185] CT slice image B;

[0186] Laser image segmentation mask C;

[0187] CT image segmentation mask D;

[0188] CT mask space transformation result E.

[0189] 1) Image preprocessing (see...) Figure 2 (Left side process)

[0190] The image preprocessing module (S1) performs the following operations on the input laser confocal image A and CT image B respectively:

[0191] Standardized format: CLSM images are converted from RGB to grayscale, and the RGB images from laser confocal microscopy are converted to grayscale images, as shown in the following formula.

[0192] I=0.2989·R+0.5870·G+0.1140·B

[0193] Where I is the pixel intensity value of the output grayscale image, and R, G, B represent the pixel values ​​of each channel in the original RGB image, typically in the range [0, 255], to remove color interference;

[0194] Gray-level normalization: This unifies the gray-level range of the two images to [0, 255]. Gray-level normalization standardizes the gray-level values ​​of the laser confocal image and the CT image to the range of 0 to 255, eliminating brightness differences. The formula is as follows:

[0195]

[0196] Where I norm (x, y) represents the normalized gray value, and I represents the gray value of a pixel in the original image. max I represents the maximum grayscale value in the image. min This represents the minimum grayscale value in the image;

[0197] Noise suppression: Gaussian filtering (σ = 1.5) was used for laser images, and nonlocal mean denoising was used for CT images;

[0198] Contrast enhancement: The CLAHE algorithm (8×8 grid) is used to improve the visibility of boundary structures;

[0199] Resolution registration: The CLSM image is scaled to the CT image size using bicubic interpolation.

[0200] The output images A and B from this step are used as inputs for subsequent segmentation.

[0201] 2) Deep semantic segmentation and mask generation (see...) Figure 2 (Central)

[0202] A deep learning model was used to perform semantic segmentation on images A and B. The network uses a classic U-structure as the backbone, superimposed with a texture feature enhancement module, an attention feature selection module, and a multi-layer visual Mamba sequence modeling structure, which can effectively identify structures such as pores, matrix, cracks, and cement in core images.

[0203] Specifically: the texture enhancement module performs boundary detection based on first-order and second-order gradient information, and captures small features by combining shallow texture convolution; the feature selection module activates important channels through Tanh nonlinear transformation and introduces residual skip connections to solve the gradient vanishing problem; the visual Mamba module models high-resolution long sequence features through selective state space, taking into account both global perception capability and efficient inference.

[0204] To improve the robustness of the network on multimodal datasets, this invention employs various data augmentation strategies, including rotation, translation, brightness perturbation, and random cropping, to augment training samples, partially utilizing the MedAugment augmentation scheme; the final output consists of two binary masks:

[0205] The laser image marker mask C and the CT image mask D are used to accurately mark the region of interest, providing ROI support for subsequent key point detection.

[0206] This module marks the Region of Interest (ROI) in the image for subsequent feature extraction.

[0207] 3) Mask edge extraction processing (see...) Figure 2 (Right side)

[0208] Perform morphological operations on masks C and D:

[0209] Morphological processing involves a series of operations based on changes in image structure, primarily including dilation, erosion, and difference operations. These operations play a crucial role in edge extraction and feature matching. Dilation enhances the connectivity of objects in an image by expanding their boundaries. The specific mathematical formula is as follows:

[0210]

[0211] Where A is a binary image, B is a structuring element, and B z This indicates that the structuring element B is translated to position Z. This represents the erosion operator. This represents the expanded image. When the area covered by structuring element B intersects with the foreground region of A, the center pixel Z is marked as a foreground pixel.

[0212] Dilation helps expand the target region in an image, preparing it for subsequent erosion and edge extraction. Erosion, on the other hand, reduces noise and clutter by shrinking the boundaries of the target region; the formula is:

[0213]

[0214] in, This means that when all pixels of structuring element B are completely contained within the foreground region of A, the center pixel Z is retained as a foreground pixel. This represents the image after erosion. Through the erosion operation, the boundaries of the target object are shrunk, effectively reducing small stray objects. The difference operation between dilation and erosion is used to extract the edge regions of the target object; its formula is:

[0215] Edge Roi =A dilated -A eroded

[0216] Among them, A dilated Indicating the result of the expansion, A eroded Indicates the corrosion result, Edge Roi This method represents the edge region of a marker. Through this difference operation, the edge texture features of the marker can be clearly extracted, thereby more accurately locating the boundary of the ROI region. This method significantly reduces the interference of background noise, providing higher accuracy and robustness for subsequent feature extraction and matching.

[0217] like Figure 3 As shown, this process generates an edge region E, which is used to precisely define the search range for key points.

[0218] 4) Keypoint detection and feature description (see...) Figure 3 , Figure 4 , Figure 5 )

[0219] The ORB feature extraction module (S4) detects key points in the edge region E, and the steps are as follows:

[0220] Corner detection using the FAST algorithm:

[0221] FAST corner detection uses the FAST algorithm to detect key points within the ROI region. The criteria for identifying a point as a corner are based on the continuity of pixel grayscale changes within a certain range around the candidate point. When the grayscale values ​​of multiple adjacent pixels are simultaneously higher or lower than the grayscale value of the center pixel, the point is considered a corner. In the edge region of the cement, the grayscale changes are more pronounced, thus it is identified as a feature point. This method is computationally efficient and can quickly detect key points in an image. In the feature description stage, the BRIEF algorithm is used to generate binary descriptors for the extracted key points, as shown in the following formula:

[0222]

[0223] p is the position of the currently detected image pixel, Ip is the gray value of pixel ppp, xi is the i-th pixel in the circular neighborhood, Ixi is the gray value of the neighboring pixel xi, t is the gray difference threshold, which determines the corner sensitivity, and n is the number of pixels that continuously meet the condition.

[0224] To enhance rotation invariance, ORB performs orientation normalization by calculating the principal direction for each corner point. The principal direction is calculated based on the gray-level centroid of the corner point's neighborhood image and is defined as follows:

[0225]

[0226] (x i ,y i Let I(x) be the coordinates of the pixel in the neighborhood of the corner point. i ,y i ) is pixels (x i ,y i The grayscale value of ) ∑y i I(x i ,y i ∑x is a weighted sum of the y-coordinates of each pixel in its neighborhood multiplied by its grayscale value. i I(x i ,y i ) is the weighted sum of the x-coordinate of each pixel in the neighborhood multiplied by its gray value, and θ is the main direction of the corner point, which is used for the rotation normalization processing of the BRIEF descriptor in the subsequent process.

[0227] Generate a 256-dimensional binary feature descriptor using BRIEF:

[0228] After completing corner detection and orientation estimation, ORB uses the BRIEF algorithm to generate a compact binary feature vector for each corner. BRIEF constructs features by comparing the grayscale values ​​of a predefined set of pixel pairs (pi, qi) within the corner's neighborhood, with each comparison result being a one-bit binary value. Its basic formula is:

[0229]

[0230] p is the current corner point, pi is the pixel pair selected in the neighborhood of p (the relative position coordinates can be predefined or randomly sampled), I(pi) is the gray value of pixel point pi and qi, and τ is a binary bit of the BRIEF descriptor.

[0231] The results of all τ values ​​are concatenated to form an N-dimensional binary feature vector, typically 256 bits. ORB performs orientation compensation on the original BRIEF, rotating the pixel pair (pi,qi) according to the principal direction θ of the corner points, making the feature description rotation-invariant. This step quickly obtains initial feature point matching pairs between the laser confocal image A and the CT image B, providing support for rigid transformation estimation. Figure 4 This is a schematic diagram showing the distribution of key points extracted from the laser image and CT image respectively (green dots are FAST corner points).

[0232] 5) Feature point matching

[0233] Feature point matching and the RANSAC module (S5) perform the following steps:

[0234] First, FLANN-based feature point matching is performed:

[0235] For each ORB feature point extracted from laser image A and CT image B, construct a 256-dimensional binary feature descriptor;

[0236] Use the FLANN matching algorithm to compare the Hamming distance between two sets of feature descriptors;

[0237] For each laser image feature point, find its nearest and second nearest neighbors in the CT image, and calculate their descriptor distances d1 and d2.

[0238] Perform a ratio test to eliminate unstable matching pairs:

[0239]

[0240] The threshold τ = 0.75. When the ratio is less than the threshold, the match is considered reliable and the pair is retained; otherwise, it is discarded.

[0241] Furthermore, it is optimized to accelerate the search by using a random Kd-tree:

[0242] To improve matching efficiency, the feature vector index is constructed using a randomized Kd-tree structure from FLANN:

[0243] Multiple random Kd-trees are constructed in parallel, with each tree using a different random subset;

[0244] When querying, you can specify the maximum number of backtracks and the maximum number of checkpoints to balance matching speed and accuracy;

[0245] In this embodiment, the number of trees is set to 4, the maximum number of search checks is 64, and the average matching time is approximately 1.2 seconds per sample.

[0246] Finally, RANSAC is used to remove mismatches:

[0247] There may be mismatches in the set of matching point pairs. To remove outlier point pairs, the RANSAC algorithm is used to estimate the initial rigid transformation model, as follows:

[0248] Randomly sample the smallest subset (2 pairs of points) from the set of matching points;

[0249] Calculate the affine rigid transformation model based on this subset;

[0250] Substitute all matching points into the model and calculate its projection error (Euclidean distance):

[0251] ε i =‖P i CT -T(P i CLSM )‖

[0252] Among them, P i CLSM P represents the coordinates of key points in the laser image. i CT Here are the coordinates of the matching points on the CT image, T is the transformation function, and ε is the coordinates of the matching points. i If the number of pixels is less than 2, it is included in the inner point set. Steps 1 to 4 are repeated for a total of 1000 iterations. The model with the most inner points is selected as the final model. This yields the reliable matching points selected by RANSAC.

[0253] 7) Transformation estimation

[0254] Calculation of rigid transformation matrix:

[0255] Based on the reliable matching points selected by RANSAC, a rigid transformation model is used to calculate the transformation matrix between images. The rigid transformation model includes translation, rotation, and scaling operations, enabling spatial transformation of images while preserving their structural integrity. This transformation effectively aligns identical regions within the images while avoiding any distortion of the image structure. The transformation matrix is ​​calculated based on the coordinates of the selected matching point pairs, using the least squares method to find the optimal transformation parameters. The specific transformation matrix formula is as follows:

[0256]

[0257] In the transformation matrix, S is the scaling factor, θ is the rotation angle, and t x and t y These are the translation amounts. This transformation matrix defines how to align a floating image (CT slice) to a fixed image (laser confocal image) through rotation, scaling, and translation operations. By calculating the transformation matrix, the floating image can be spatially transformed to precisely align it with the fixed image, thus obtaining the rigid transformation matrix.

[0258] 8) Image spatial transformation

[0259] Based on the obtained rigid transformation matrix, a spatial transformation is performed on the floating image (CT slice image) using bilinear interpolation to ensure precise alignment with its fixed image A:

[0260]

[0261] Where Baligned(x′,y′) is the pixel value of the transformed image at (x′,y′), where (x′,y′) are the transformed coordinates, and B(x′,y′) is the pixel value of the transformed image at (x′,y′). i ,y j ) is the floating image B in (x i ,y j The pixel value of w ij These are the bilinear interpolation weighting coefficients;

[0262] The identifier segmentation mask D is synchronously transformed using nearest neighbor interpolation:

[0263] F(x′,y′)=D(round(x),round(y))

[0264] Where F(x′,y′) represents the transformed mask, D(x,y) represents the pixel value of the segmented mask, and round(x) and round(y) represent rounding the floating-point coordinates to the nearest integer coordinates;

[0265] Finally, the transformed CT slice image identifier segmentation mask E was obtained, maintaining consistency with the transformed image.

[0266] 9) Registration result output

[0267] The registered optimal CT slice image and the laser confocal image A are spatially aligned and saved using a unified coordinate system in a lossless compressed format (TIFF or PNG). Simultaneously, two corresponding marker mask files are saved, containing the marker segmentation mask E for the laser confocal image A and the CT slice image. Transformation parameter records: rigid transformation matrix parameters are stored in a machine-readable text format.

[0268] 10) Registration quality assessment

[0269] Calculate the Dice coefficient between the laser confocal image marker segmentation mask C and the registered CT image marker segmentation mask E;

[0270] Dice coefficient calculation: Dice is an index that measures the similarity between two sets. The calculation formula is as follows:

[0271]

[0272] Here, A and B represent the segmentation results of the laser confocal image markers and the CT image markers, respectively. A∩B represents the overlapping area between the two, and |A| and |B| represent the number of pixels in the mask regions A and B, respectively. The value of Dice ranges from 0 to 1. A value close to 1 indicates high registration accuracy between the two images, meaning the marker regions highly overlap; a value close to 0 indicates poor registration with significant deviation. This metric effectively evaluates the quality of image registration and can be used to optimize registration algorithms to improve registration accuracy and stability.

[0273] Optimal slice selection rule: Iterate through the Dice coefficients of all CT slices; select the slice that satisfies Dice ≥ 0.90 and takes the maximum value as the final registration result;

[0274] If all slices have a Dice < 0.90, then initiate iterative optimization: readjust the FLANN matching threshold in step S5 between 0.5 and 0.8 and adjust the RANSAC tolerance in step S6 between 2.0 and 5.0 pixels. Then, repeat steps S5 to S9.

[0275] It should be noted that the above embodiments are only used to illustrate the present invention and are not intended to limit the present invention. After understanding the technical solution and core concept of the present invention, those skilled in the art can make various modifications, equivalent substitutions, or simple variations to the embodiments within the scope of protection defined in the claims. These modifications, substitutions, or variations should all be considered to be included within the scope of protection of the present invention.

Claims

1. A core image registration method based on multimodal image segmentation and feature matching, characterized in that, Includes the following steps: S1: Preprocess the laser confocal image and CT slice image, including grayscale normalization, noise removal and contrast enhancement, to unify the grayscale distribution and resolution of the two modes, and obtain laser confocal image A and CT slice image B. S2: Train a deep learning semantic segmentation model using a core image dataset containing pixel-level annotations of markers in laser confocal image A and slice image B. The trained deep learning semantic segmentation model is then used to segment the preprocessed images, generating a binary mask file that includes laser confocal image marker segmentation mask C and CT image marker segmentation mask D. S3: Morphological processing is performed on the laser confocal image marker segmentation mask C and the CT image marker segmentation mask D generated by S2, including dilation, erosion and difference operations, to obtain the accurate marker edge region; S4: Based on the precise marker edge feature region obtained in S3, key points are detected and key feature points are generated using the ORB algorithm; S5: Using laser confocal image A as a fixed image and CT slice image B as a floating image, the key feature points obtained in S4 are matched using the FLANN matching algorithm and combined with a random Kd tree to accelerate the search, thus obtaining the corresponding feature point pairs between laser confocal image A and CT slice image B. S6: Based on the corresponding feature point pairs obtained in S5, the RANSAC algorithm is applied to remove mismatched point pairs. In addition, the reliable matching point pairs selected by RANSAC are used to calculate the rigid transformation matrix of the floating image B to the fixed image A, and the rigid transformation matrix is ​​obtained. S7: Based on the rigid transformation matrix obtained in S6, perform a spatial transformation on the floating image to align it with the fixed image. Simultaneously, transform the segmentation file and perform the same spatial transformation on the CT image corresponding marker segmentation mask file D to obtain the transformed CT slice image marker segmentation mask E, maintaining consistency with the transformed image, and obtain the best CT slice image after registration. S8: The best CT slice image obtained after registration based on S7 is spatially aligned and saved with the original laser confocal image A according to the same coordinates. At the same time, the corresponding marker mask files are saved, which contain the marker segmentation mask E of the laser confocal image A and the CT slice image. The parameters of the final transformation matrix used for the slice registration are recorded, thus completing the generation of the core registration result. S9: Finally, the registration quality is evaluated by calculating the Dice coefficient between the laser confocal image marker segmentation mask C and the registered CT slice image marker segmentation mask E. The CT slice with the highest Dice coefficient is selected as the final registration result of the core sample.

2. The core image registration method based on multimodal image segmentation and feature matching according to claim 1, characterized in that, The preprocessing in step S1 includes: Gray-level normalization: unifies the gray-level range of laser confocal images and CT images to [0, 255]; Noise removal: Gaussian filtering (σ = 1.0) was used for laser confocal images, and non-local means denoising (NL-Means) was used for CT images; Contrast enhancement: Histogram equalization (CLAHE algorithm, 8×8 grid size); Resolution unification: The laser confocal image is downsampled to the CT slice resolution using bicubic interpolation; Format conversion: Convert CLSM image to 8-bit grayscale image; Output: Preprocessed laser confocal image A and CT slice image B.

3. The core image registration method based on multimodal image segmentation and feature matching according to claim 2, characterized in that, Step S2 includes: Data augmentation operations were performed on the preprocessed laser confocal images and CT images, including random cropping, rotation, flipping, and brightness adjustment, to enhance the model's ability to generalize to complex textures in core images. A semantic segmentation model based on selective state space is adopted. The model includes a texture feature enhancement module, a feature selection module, and a visual Mamba layer. The texture feature enhancement module extracts microstructural features such as particles, pores, and cracks from core images through multi-scale convolution operations. The feature selection module selects key features related to markers based on an attention mechanism to suppress the influence of background noise. The visual Mamba layer uses a selective state space mechanism to efficiently model long sequence features and generate high-precision pixel-level segmentation results. The model was trained using a labeled core image dataset containing pixel-level annotations of markers in laser confocal images and CT images. The model parameters were optimized using the cross-entropy loss function and the Adam optimizer during training. The output of the semantic segmentation model is a pixel-level probability map of the same size as the input image. The value of each pixel in the probability map represents the predicted probability that the pixel belongs to the category of the marker. The probability map output by the model is post-processed by setting a threshold to convert the probability map into a binary mask file, which contains the laser confocal image marker segmentation mask C and the CT image marker segmentation mask D.

4. The core image registration method based on multimodal image segmentation and feature matching according to claim 3, characterized in that, Step S3 includes: The dilation operation enhances the connectivity of a target object in an image by expanding its boundaries. It is used to expand the target region in an image, preparing it for subsequent erosion and edge extraction. The specific mathematical formula is as follows: Where A is a binary image, B is a structuring element, and B z This indicates that the structuring element B is translated to position Z. The erosion operator is represented by A⊕B, which represents the dilated image. When the area covered by the structuring element B intersects with the foreground region of A, the center pixel Z is marked as a foreground pixel. Etching reduces noise and stray areas by narrowing the boundary of the target region, as shown in the formula: in, This means that when all pixels of structuring element B are completely contained within the foreground region of A, the center pixel Z is retained as a foreground pixel. This represents the image after erosion; through the erosion operation, the boundaries of the target object are shrunk, effectively reducing small stray objects. The difference operation between dilation and erosion is used to extract the edge region of the target object. It can clearly extract the edge texture features of the marker, thereby more accurately locating the boundary of the ROI region. The formula is: Edge Roi =A dilated -A eroded Among them, A dilated Indicating the result of the expansion, A eroded Indicates the corrosion result, Edge Roi This indicates the edge area of ​​the sign.

5. The core image registration method based on multimodal image segmentation and feature matching according to claim 4, characterized in that, In step S4, the ORB algorithm is used to detect key points and generate key feature points, following the path below: Within the ROI area, the FAST algorithm is used to detect key points. The criteria for determining a point as a corner are based on the continuity of pixel grayscale changes within a certain range around the candidate point. When the grayscale values ​​of multiple adjacent pixels are simultaneously higher or lower than the grayscale value of the center pixel, the point can be considered a corner. For the detected key points, the BRIEF algorithm is used to generate feature descriptors. This algorithm selects pixel pairs in the local neighborhood and compares their grayscale relationships to construct descriptive information in binary form. This results in binary feature descriptors generated from the key points, thus realizing the detection of key points and generation of key feature points using the ORB algorithm.

6. The core image registration method based on multimodal image segmentation and feature matching according to claim 5, characterized in that, In step S5, using laser confocal image A as a fixed image and CT slice image B as a floating image, the matching process using the FLANN feature matching algorithm includes: Calculate feature point similarity based on Euclidean distance from feature descriptors; Calculate the ratio of the nearest neighbor distance to the second nearest neighbor distance for each feature point; When the ratio is less than a set threshold (range 0.5-0.8), the nearest neighbor matching point is accepted as a reliable matching point; A random Kd-tree data structure is used to accelerate the nearest neighbor search of feature points; Match the corresponding feature point pairs between laser confocal image A and CT slice image B.

7. The core image registration method based on multimodal image segmentation and feature matching according to claim 6, characterized in that, The rigid transformation matrix obtained in step S6 is executed according to the following path: Transformation parameter calculation: Based on reliable matching point pairs selected by RANSAC, the optimal transformation parameters are solved using the least squares method. Where N is the number of interior points, and pCLSM and pCT are the coordinates of the matching points in the fixed image and the floating image, respectively; Based on the obtained reliable matching point pairs, a similarity transformation model is adopted, and the specific transformation matrix formula is as follows: In the transformation matrix, S is the scaling factor, θ is the rotation angle, and t x and t y These are the translation amounts; The transformation matrix defines how to align a floating image to a fixed image through rotation, scaling, and translation operations, with laser confocal image A as the fixed image and CT slice image B as the floating image.

8. The core image registration method based on multimodal image segmentation and feature matching according to claim 7, characterized in that, Step S7 further includes: Spatial transformation is performed: The floating image is spatially transformed using bilinear interpolation according to the following formula to ensure precise alignment with its fixed image: Where (x′, y′) are the transformed coordinates, w ij These are bilinear weighting coefficients; The identifier segmentation mask D is synchronously transformed using nearest neighbor interpolation: F(x′,y′)=D(round(x),round(y)) Finally, the transformed CT slice image identifier segmentation mask E was obtained, maintaining consistency with the transformed image.

9. The core image registration method based on multimodal image segmentation and feature matching according to claim 8, characterized in that, In step S8, the generation of core registration results also includes: Registration data storage: Save the registered best CT slice and the original laser confocal image in a unified coordinate system, and store them in the lossless compressed format TIFF or PNG format. At the same time, save the corresponding marker mask file, which includes the marker segmentation mask E of the laser confocal image A and the CT slice image. Transformation parameter record: The parameters of the rigid transformation matrix are stored in a machine-readable text format.

10. The core image registration method based on multimodal image segmentation and feature matching according to claim 9, characterized in that, The registration quality assessment in step S9 includes: Dice coefficient calculation: Dice is an index that measures the similarity between two sets. The calculation formula is as follows: Where A and B represent the segmentation results of laser confocal map markers and CT image markers, respectively; A∩B represents the overlapping area of ​​the two; and |A| and |B| represent the number of pixels in the mask A and B regions, respectively. Optimal slice selection rule: Iterate through the Dice coefficients of all CT slices; select the slice that satisfies Dice ≥ 0.90 and takes the maximum value as the final registration result; If all slices have a Dice < 0.90, then initiate iterative optimization: readjust the FLANN matching threshold in step S5 between 0.5 and 0.8 and adjust the RANSAC tolerance in step S6 between 2.0 and 5.0 pixels. Then, repeat steps S5 to S9.

Images