Rock component identification method and system based on energy spectrum CT

Abstract

The invention discloses a rock component identification method and system based on energy spectrum CT, and belongs to the technical field of rock component identification. The method comprises the following steps: selecting a calibration motif of a to-be-detected sample, and respectively performing CT scanning and reconstruction on the calibration motif and the to-be-detected sample in a high energy region and a low energy region to obtain high-energy region projection data L1 and low-energy region projection data L2 of the calibration motif and high-energy region projection data L '1 and low-energy region projection data L' 2 of the to-be-detected sample; performing substance decomposition based on the high-energy region projection data L1 and the low-energy region projection data L2 of the calibration motif and the high-energy region projection data L '1 and the low-energy region projection data L' 2 of the sample to be detected to obtain a decomposition coefficient ai (r) of each base substance i; and obtaining the electron density rho ei and the effective atomic number Zi of the base substance i, and combining with the decomposition coefficient ai (r) of the base substance i to obtain a rock component identification result of the to-be-detected sample. According to the invention, the internal three-dimensional structure information and material composition information of the rock can be obtained losslessly and efficiently.

Description

Technical Field

[0001] The invention belongs to the technical field of rock component identification and proposes a rock component identification method and system based on energy spectrum CT. Background Art

[0002] The mineral composition and content of rocks are of great significance for lithology classification, sedimentary environment judgment and source direction analysis. At present, the methods for determining rock composition can be mainly divided into qualitative identification and quantitative analysis. Qualitative identification mainly includes specimen identification and thin section observation. This traditional lithology component identification method relies too much on experience and geological professional knowledge accumulation. It is not only time-consuming and highly professional, but also easily affected by subjective factors. Quantitative analysis is to identify rock components by testing the types and contents of different elements or minerals in rocks. There are mainly infrared spectroscopy, X-ray diffraction (XRD), X-ray fluorescence analysis (XRF), scanning electron microscopy, etc. XRF can perform qualitative and quantitative analysis on the elements on the surface of sample slices, and XRD can perform qualitative and quantitative analysis on the minerals of samples. The results are accurate and have high reference value, but the test requires pretreatment such as slicing and grinding, and cannot be non-destructive and rapid detection.

[0003] With the continuous development of computed tomography (CT) imaging technology, it has become possible to describe the material distribution and spatial structure inside rocks with the help of CT's three-dimensional imaging technology and digital image processing technology, which can achieve non-destructive detection and three-dimensional reconstruction of the inside of rocks. Different material components in rocks have different attenuation capabilities for X-rays, which is manifested in different grayscale values ​​of each material in the CT reconstructed image, that is, different linear attenuation coefficients. The difference in linear attenuation coefficients can be used to identify some rock mineral components.

[0004] Conventional CT imaging technology only obtains the linear attenuation information of materials through the intensity and position of X-rays. However, the linear attenuation coefficient is a physical quantity related to the ray energy, the density of the detected object and the atomic number of the detected object. When using the same ray energy for CT imaging, the difference in linear attenuation coefficients may be caused by different densities and atomic numbers of the detected objects, or by differences in both physical quantities. In addition, different detected objects show the same grayscale value due to their different atomic numbers and material densities. Therefore, in conventional CT imaging, different detected objects may show the same grayscale value, and the material composition cannot be accurately identified. Spectral CT collects two sets of ray imaging data of different energies in one scanning imaging, providing energy information for data reconstruction. How to use the attenuation data of materials at different energies to obtain more accurate material composition information for identification and differentiation of rock components is an urgent problem to be solved. Summary of the invention

[0005] In order to make up for the fact that traditional quantitative analysis methods cannot quickly and non-destructively obtain the composition and structure information inside the rock, and that conventional CT's rock composition identification is not accurate enough, the present invention provides a rock composition identification method and system based on energy spectrum CT, which can non-destructively and efficiently obtain the internal three-dimensional structure information and material composition information of the rock.

[0006] To achieve the above-mentioned object of the invention, the technical solution of the present invention includes the following contents.

[0007] A rock component identification method based on energy spectrum CT, the method comprising:

[0008] A calibration phantom of the sample to be detected is selected, and CT scanning and reconstruction of the calibration phantom and the sample to be detected in high and low energy zones are performed respectively to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, as well as high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be detected; wherein the calibration phantom comprises a plurality of calibration samples i composed of a base material i, and the settings of the high energy zone and the low energy zone are obtained based on the attenuation capability of the sample to be detected to X-rays;

[0009] Based on the high-energy projection data L1 and low-energy projection data L2 of the calibration phantom, and the high-energy projection data L′1 and low-energy projection data L′2 of the sample to be tested, the material is decomposed to obtain the decomposition coefficient a of each base material i i (r);

[0010] Get the electron density ρ of base material i ei and the effective atomic number Z i , and combined with the decomposition coefficient a of base substance i i (r), and obtain the electron density ρ of the sample to be tested e and effective atomic number Z eff , and then obtain the rock component identification result of the sample to be tested by looking up the table.

[0011] Furthermore, the step of selecting a calibration phantom of a sample to be tested comprises:

[0012] Perform projection imaging on the sample to be tested and record the maximum attenuation value of the X-ray passing through the sample to be tested;

[0013] Measuring the volume of the sample to be tested to obtain an alternative calibration sample;

[0014] Based on density weighing method and qualitative identification method, the electron density and atomic number of the sample to be tested are roughly estimated;

[0015] A calibration sample 1 composed of a base material 1 is selected from the candidate calibration samples; wherein the difference between the maximum attenuation value of X-rays passing through the calibration sample 1 and the maximum attenuation value of X-rays passing through the sample to be detected is less than a second set threshold, and the electron density and atomic number of the calibration sample 1 are both greater than a rough estimate of the sample to be detected;

[0016] Selecting a calibration sample 2 composed of a base substance 2 from the candidate calibration samples; wherein the electron density and atomic number of the calibration sample 2 are both smaller than the rough estimate of the sample to be tested;

[0017] A calibration sample i composed of base material i is selected from the alternative calibration samples; wherein the maximum attenuation value of X-rays passing through the calibration sample i is less than the maximum attenuation value of X-rays passing through the calibration sample 1, and is greater than the maximum attenuation value of X-rays passing through the calibration sample 2.

[0018] Furthermore, the high energy region and the low energy region are obtained based on the attenuation capability of the sample to be detected to X-rays, including:

[0019] In the high energy region and the low energy region, the photon count of the X-ray passing through the sample to be detected is 25%-40% of the photon count in the empty field;

[0020] and,

[0021] The empty-field photon count of X-rays in the low-energy region is 40%-60% of the empty-field electron count of the X-rays in the high-energy region.

[0022] Furthermore, the high-energy region projection data L1 and low-energy region projection data L2 of the calibration phantom, and the high-energy region projection data L′1 and low-energy region projection data L′2 of the sample to be tested are subjected to material decomposition to obtain the decomposition coefficient a of each base material i. i (r) including:

[0023] Based on the high-energy region projection data L1 and the low-energy region projection data L2 of the calibration phantom, the n basis images f of the calibration phantom are calculated. n (r), the n base images f n The vector of (r) is represented as f(r); where r represents the voxel spatial position;

[0024] By reconstructing the high-energy region projection data L1, the CT reconstruction data of the high-energy region is obtained, and the CT reconstruction data of the high-energy region is subjected to image segmentation processing to obtain the template image t of the base material i i (r);

[0025] Based on the template image t i (r) and the vector representation f(r) are respectively subjected to the least square fitting calculation for each basis material i, and the polynomial coefficient vector c corresponding to the basis material i is obtained.i ;

[0026] Based on the high-energy region projection data L′1 and low-energy region projection data L′2 of the sample to be tested, the n basis images f of the calibration phantom are calculated. n ′(r), the n basis images f n The vector of ′(r) is denoted as f′(r);

[0027] Calculate the decomposition coefficient a of base substance i i (r) = c i ·f′(r).

[0028] Furthermore, the n basis images f of the calibration phantom n (r) = R -1 b n (L1,L2); where R -1 represents the inverse Radon transform, b n (L1, L2) represents a reconstructed image obtained based on the high-energy region projection data L1 and the low-energy region projection data L2.

[0029] Furthermore, the template image of the base material i

[0030] Furthermore, the polynomial coefficient vector corresponding to the base material i

[0031] Furthermore, the electron density ρ of the base substance i is obtained ei , and combined with the decomposition coefficient c of base substance i i , obtaining the rock component identification result of the sample to be tested, including:

[0032] Get the electron density ρ of base material i ei ;

[0033] Based on the decomposition coefficient a i (r) The electron density ρ of the base substance i ei Weighted calculation to obtain the electron density ρ of the sample to be detected e ;

[0034] Calculate the electron density percentage of base substance i

[0035] According to the electron density percentage f of each base substance i i and the effective atomic number Z of the base substance i i , get the effective atomic number Z of the sample to be tested eff ;

[0036] Based on the electron density ρ e and the effective atomic number Z eff A table lookup is performed to obtain the rock component identification result of the sample to be tested.

[0037] Furthermore, after performing CT scanning and reconstruction of the calibration phantom and the sample to be tested in high and low energy zones respectively to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, and high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be tested, the method further includes:

[0038] Median filtering is performed on the high energy area projection data L1, the low energy area projection data L2, the high energy area projection data L′1 and the low energy area projection data L′2

[0039] and / or,

[0040] Performing mean filtering on the high energy region projection data L1, the low energy region projection data L2, the high energy region projection data L′1, and the low energy region projection data L′2;

[0041] and / or,

[0042] Perform pixel merging and projection frame merging on the high energy area projection data L1, the low energy area projection data L2, the high energy area projection data L′1 and the low energy area projection data L′2;

[0043] and / or,

[0044] The CT data reconstruction process is performed on the high energy region projection data L1, the low energy region projection data L2, the high energy region projection data L′1 and the low energy region projection data L′2 respectively through the energy spectrum CT reconstruction algorithm.

[0045] A rock component identification system based on energy spectrum CT, the system comprising:

[0046] A data acquisition module, used for selecting a calibration phantom of a sample to be detected, and performing CT scanning and reconstruction of the calibration phantom and the sample to be detected in high and low energy zones, respectively, to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, and high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be detected; wherein the calibration phantom comprises a plurality of calibration samples i composed of a base material i, and the settings of the high energy zone and the low energy zone are obtained based on the attenuation capability of the sample to be detected to X-rays;

[0047] The material decomposition module is used to perform material decomposition based on the high-energy region projection data L1 and low-energy region projection data L2 of the calibration phantom, and the high-energy region projection data L′1 and low-energy region projection data L′2 of the sample to be tested, and obtain the decomposition coefficient a of each base material i. i (r);

[0048] Component identification module, used to obtain the electron density ρ of the base substance i ei and the effective atomic number Z i , and combined with the decomposition coefficient a of base substance i i (r), and obtain the electron density ρ of the sample to be tested e and effective atomic number Z eff , and then obtain the rock component identification result of the sample to be tested by looking up the table.

[0049] Compared with the existing technology, conventional CT imaging technology only obtains the linear attenuation information of the material through the intensity and position of the X-ray. However, the linear attenuation coefficient is a physical quantity related to the ray energy, the electron density of the detected object and the effective atomic number of the detected object. Different materials may show the same linear attenuation coefficient value under the same ray energy due to the different electron density and effective atomic number of their components. For such materials, it is difficult to distinguish and identify the material components through the gray value analysis method of traditional CT imaging. Therefore, it is not accurate enough to identify the rock components only through the gray value of the linear attenuation coefficient. Obtaining the electron density and effective atomic number information of the sample to be tested through material decomposition can more accurately distinguish and identify different material components. BRIEF DESCRIPTION OF THE DRAWINGS

[0050] Figure 1 Execution block diagram of the rock composition identification method based on spectral CT.

[0051] Figure 2 Implementation flow chart of material decomposition.

[0052] Figure 3 Actual image of the test sample.

[0053] Figure 4 Grayscale image of the test sample CT imaging achieved based on existing technology.

[0054] Figure 5 The effective atomic number image and electron density image of the test sample are realized based on the present invention. DETAILED DESCRIPTION

[0055] The present invention is further described in detail below in conjunction with specific embodiments. The given examples are only for illustrating the present invention, but not for limiting the scope of the present invention.

[0056] The present invention provides a lithology identification method based on spectral CT, which mainly includes spectral CT imaging of calibration phantom and rock sample, CT image preprocessing, material decomposition and rock component identification. The execution block diagram is as follows: Figure 1 shown.

[0057] Step 1: Calibrate spectral CT imaging of the phantom and rock sample.

[0058] In order to identify rock composition by utilizing the attenuation information of radiation at different energies, it is necessary to estimate the effective atomic number and electron density of the object being tested by material decomposition, thereby achieving more accurate rock composition identification.

[0059] For different rock samples, it is necessary to select a suitable calibration phantom for the correction processing of material decomposition. First, it is necessary to perform projection imaging on the rock sample to be tested, record the maximum attenuation value of the X-ray passing through the sample, the volume of the standard sample needs to be similar to the sample to be tested, and the attenuation value of the X-ray passing through it is similar to the maximum attenuation value of the X-ray passing through the sample. Select a standard part that meets the above two conditions as calibration sample 1, and select a standard part with a volume similar to the object to be tested and a smaller density and atomic number as calibration sample 2. If more calibration samples need to be selected, then the maximum attenuation value of the X-ray passing through these calibration samples is less than the maximum attenuation value of the X-ray passing through the calibration sample 1, and greater than the maximum attenuation value of the X-ray passing through the calibration sample 2. Finally, all the calibration samples obtained are combined as the calibration phantom required for material decomposition. Use the same CT imaging conditions to perform CT scanning and reconstruction on the calibration phantom and the sample to be tested respectively. Set the energy range of CT imaging according to the attenuation ability of the sample to be tested to the ray. Specifically, the energy range of CT imaging must meet two conditions at the same time: 1) In the high-energy zone and the low-energy zone, the photon count of X-rays passing through the sample to be detected is 25%-40% of the empty-field photon count, and the preferred value is one-third; 2) The empty-field photon count of X-rays in the low-energy zone is 40%-60% of the empty-field count of X-rays in the high-energy zone, and the preferred value is 50%. After that, two sets of projection data L1 and L2 of the high and low energy zones of the calibration phantom and two sets of projection data L′1 and L′2 of the high and low energy zones of the sample to be detected can be collected respectively.

[0060] Step 2: Image preprocessing.

[0061] There are interference factors such as noise, scattering artifacts, and ring artifacts in CT data that are not conducive to material decomposition and identification. Therefore, it is necessary to preprocess the projection data of the high and low energy areas of the calibration model and the sample to be tested to improve the accuracy of material decomposition. First, median filtering is performed on each projection data to reduce bad pixels in the image, thereby reducing ring artifacts in the reconstructed image. Secondly, mean filtering is performed on each projection data to reduce noise in the image and improve the signal-to-noise ratio of the image. In addition, pixel merging and projection frame merging can also improve the signal-to-noise ratio of the image. The CT data reconstruction is performed on the processed projection image data of the high and low energy areas respectively through the energy spectrum CT reconstruction software to obtain two sets of CT reconstruction data of the high and low energy areas.

[0062] Step 3: Material Decomposition

[0063] The linear attenuation coefficient is a physical quantity related to the ray energy, the atomic number and density of the material, and can be expanded into a linear combination of basis functions:

[0064]

[0065] Among them, μ(E,r) represents the linear attenuation coefficient of the material, and the basis function ψ i (E) is the energy E of the X-ray, and the decomposition coefficient a i (r) is a quantity related to the voxel spatial position r. When the basis function represents the energy-related quantity of the photoelectric effect and Compton scattering, the decomposition coefficient is a quantity related to the atomic number and electron density of the material. When the basis function is expressed as the linear attenuation coefficient of the basis material, the decomposition coefficient is the volume fraction or mass fraction of each basis material. In addition, the linear attenuation coefficient of the material can be expanded into a linear combination of the mass attenuation coefficients of each basis material. At this time, the decomposition coefficient represents the material density of each basis material. The key to material decomposition is to solve the decomposition coefficient a. i (r). In this method, the calibration sample i contained in the calibration model represents the i-th basis material, and the basis function ψ i (E) represents the linear attenuation coefficient of the ith calibration sample, and the decomposition coefficient a i (r) represents the volume fraction of the i-th calibration sample.

[0066] The attenuation of X-rays passing through matter can be expressed as a weighted average of the attenuations at different energies:

[0067]

[0068] Where S(F) represents the energy spectrum of X-rays, D(E) represents the energy response of the detector, and I0 represents the empty field attenuation when the ray does not pass through the object. According to formula (1), the attenuation of the ray passing through the material is:

[0069]

[0070] Among them A i =∫a i (r)dr,A i represents the line integral of the basis function coefficients. A i The reconstruction of A i Radon inverse transform R -1 It is the decomposition coefficient a that is ultimately required to be solved in the material decomposition process. i (r).

[0071] The material decomposition process first requires the basis function coefficient line integral A iThe solution is the numerical inversion process of formula (3). The projection data L1 and L2 of the two energy zones of high and low energy are established using high-order polynomials and the line integral A of the basis function coefficients. i Assume that:

[0072] A i =D i (L1,L2) (4)

[0073] Use higher-order polynomial expansion:

[0074]

[0075] Among them, c i,n b n (L1,L2) represents the nth term of the polynomial, c i,n represents the coefficient of the nth polynomial. The polynomial has a total of N = (J + 1) (K + 1) expansion terms, where J and K represent the order of the polynomial.

[0076] Definition b n The reconstructed image of (L1,L2), i.e. b n Inverse Radon transform R of (L1,L2) -1 For the base image:

[0077] f n (r) = R- 1 b n (L1,L2) (6)

[0078] The reconstruction of the decomposition coefficient line integral can be expressed as a linear combination of the basis images:

[0079]

[0080] c i It is a coefficient vector consisting of a set of polynomial coefficients, which can be solved by least squares estimation:

[0081]

[0082] Template Image i (r) As a priori knowledge, it can represent the density of each base material, or it can represent the volume fraction of the i-th base material in each voxel. In order to estimate the effective atomic number, the volume fraction of the base material is used as the template image:

[0083]

[0084] Taking two base substances as an example, the substance decomposition process firstly processes the data of the calibration model and calculates the coefficient vector c i: For the projection data L1 and L2 of the high and low energy zones of the calibration phantom after the image preprocessing in step 2, the n basis images f of the calibration phantom are obtained based on formula (6): n (r); by reconstructing R from the projection data L1 -1 L1 obtains the CT reconstruction data of the high-energy region, and then obtains two template images t1(r) and t2(r) through image segmentation based on formula (9), corresponding to the two base materials respectively; for the template images t1(r), t2(r) and the base image f n (r), based on formula (8), the least squares fitting calculation is performed to obtain two sets of polynomial coefficient vectors c1 and c2 corresponding to the two base materials.

[0085] For the two sets of projection data L′1 and L′2 of the high and low energy zones of the sample to be detected, the n basis images f′ of the sample to be detected are obtained based on formula (6) n (r); Polynomial coefficient vectors c1 and c2 and basis image f′ n Substituting (r) into formula (7) can obtain the decomposition coefficients a1(r) and a2(r). The implementation process of material decomposition is as follows: Figure 2 shown.

[0086] Step 4: Rock composition identification

[0087] The electron density and atomic number that determine the linear attenuation coefficient of a substance are related to the volume fraction of the base substance. i (r) The electron density ρ for a known basis substance ei The weighting can realize the electron density ρ of the material with unknown element composition e Calculation:

[0088] ρ e =∑a i (r)·ρ ei (10)

[0089] Similarly, the effective atomic number of a substance of unknown composition can be expressed as the electron density percentage f i Weighted calculation of the effective atomic number of the base substance:

[0090]

[0091] Where Z eff is the effective atomic number of the sample to be tested, Z i is the effective atomic number of the known base substance i, Represents the effective atomic number Z i to the 3.4th power.

[0092] Through the above calculation, the volume fraction information of each base substance obtained by decomposing the substance can be converted into the effective atomic number Z eff and the electron density ρ e Finally, based on the effective atomic number and electron density ρ e A query is performed in the pre-constructed rock component identification table to obtain the rock component identification result of the sample to be tested.

[0093] In summary, the present invention estimates the effective atomic number and electron density of a substance through spectral CT material decomposition, thereby achieving more accurate component identification and multiple material identification. The present invention can be applied to various dual-energy CT and spectral CT equipment.

[0094] The following is a specific experiment to illustrate the rock component identification method based on energy spectrum CT provided by the present invention.

[0095] The experiment used 10mmφ*100m cylindrical polymethyl methacrylate (PMMA) and 8mmφ*60m cylindrical titanium (Ti) standard phantoms as calibration objects, and used rock samples with different compositions as test objects, including quartz sandstone, gray mudstone, red mudstone and marl. The highest voltage of CT imaging is 120KVp. By setting the two thresholds of the photon counting detector to 45keV and 75keV respectively, CT imaging data in the two energy intervals of [45keV, 75keV] and [75keV, 120keV] are obtained for rock material decomposition and material composition identification.

[0096] Figure 3 These are the actual images of the test samples, where (a) is quartz sandstone, (b) is gray mudstone, (c) is red mudstone, and (d) is marlstone. Figure 4 The CT imaging grayscale images of different test samples are displayed, which indicates the linear attenuation coefficient information of the test samples. The numerical value represents the different attenuation capabilities of different components to radiation. Figure 4 The grayscale values ​​of different samples are not much different. The grayscale values ​​of (a) quartz sandstone and (c) red mudstone are similar, and the grayscale values ​​of (b) gray mudstone and (d) marlstone are similar. It is difficult to identify and distinguish rocks of different compositions through CT imaging grayscale images.

[0097] Figure 5The effective atomic number images and electron density images of different test samples are shown, among which (a) is the effective atomic number image of quartz sandstone, (b) is the effective atomic number image of gray mudstone, (c) is the effective atomic number image of red mudstone, (d) is the effective atomic number image of marl, (e) is the electron density image of quartz sandstone, (f) is the electron density image of gray mudstone, (g) is the electron density image of red mudstone, and (h) is the electron density image of marl. It can be seen that the effective atomic number of (a) quartz sandstone is significantly smaller than that of (c) red mudstone. The oxides and hydroxides of chemical elements such as iron and aluminum accumulate in the red mudstone, resulting in a relatively high effective atomic number of the red mudstone; while the electron density value of (a) quartz sandstone is significantly greater than that of (c) red mudstone. Similarly, the grayscale values ​​of (b) gray mudstone and (d) marl are similar and cannot be distinguished by grayscale value images, but the effective atomic number of marl is significantly higher than that of gray mudstone and the electron density is significantly lower than that of gray mudstone. Comparing (b) gray mudstone and (c) red mudstone, the effective atomic numbers of the two are similar, but the electron density is significantly different. The rock components can be identified and distinguished by the effective atomic number and electron density values.

[0098] The above embodiments are only for illustrating the technical concept and features of the present invention, and their purpose is to enable people familiar with the technology in this field to understand the content of the present invention and implement it, and they cannot be used to limit the protection scope of the present invention. Any equivalent changes or modifications made according to the spirit of the present invention should be included in the protection scope of the present invention.

Claims

1. A rock component identification method based on spectral CT, characterized in that: The method comprises: A calibration phantom of the sample to be detected is selected, and CT scanning and reconstruction of the calibration phantom and the sample to be detected in high and low energy zones are performed respectively to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, as well as high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be detected; wherein the calibration phantom comprises a plurality of calibration samples i composed of a base material i, and the settings of the high energy zone and the low energy zone are obtained based on the attenuation capability of the sample to be detected to X-rays; Based on the high-energy projection data L1 and low-energy projection data L2 of the calibration phantom, and the high-energy projection data L′1 and low-energy projection data L′2 of the sample to be tested, the material is decomposed to obtain the decomposition coefficient a of each base material i i (r); Get the electron density ρ of base material i ei and the effective atomic number Z i , and combined with the decomposition coefficient a of base substance i i (r), and obtain the electron density ρ of the sample to be tested e and effective atomic number Z eff , and then obtain the rock component identification result of the sample to be tested by looking up the table.

2. The method according to claim 1, characterized in that: The step of selecting a calibration phantom for a sample to be tested comprises: Perform projection imaging on the sample to be tested and record the maximum attenuation value of the X-ray passing through the sample to be tested; Measuring the volume of the sample to be tested to obtain an alternative calibration sample; Based on density weighing method and qualitative identification method, the electron density and atomic number of the sample to be tested are roughly estimated; A calibration sample 1 composed of a base material 1 is selected from the candidate calibration samples; wherein the difference between the maximum attenuation value of X-rays passing through the calibration sample 1 and the maximum attenuation value of X-rays passing through the sample to be detected is less than a second set threshold, and the electron density and atomic number of the calibration sample 1 are both greater than a rough estimate of the sample to be detected; Selecting a calibration sample 2 composed of a base substance 2 from the candidate calibration samples; wherein the electron density and atomic number of the calibration sample 2 are both smaller than the rough estimate of the sample to be tested; A calibration sample i composed of base material i is selected from the alternative calibration samples; wherein the maximum attenuation value of X-rays passing through the calibration sample i is less than the maximum attenuation value of X-rays passing through the calibration sample 1, and is greater than the maximum attenuation value of X-rays passing through the calibration sample 2.

3. The method according to claim 1, characterized in that The high energy region and the low energy region are obtained based on the attenuation capability of the sample to be detected to X-rays, including: In the high energy region and the low energy region, the photon count of the X-ray passing through the sample to be detected is 25%-40% of the photon count in the empty field; and, The empty-field photon count of X-rays in the low-energy region is 40%-60% of the empty-field electron count of the X-rays in the high-energy region.

4. The method according to claim 1, characterized in that: The high-energy region projection data L1 and the low-energy region projection data L2 of the calibration phantom, as well as the high-energy region projection data L′1 and the low-energy region projection data L′2 of the sample to be tested, are used to perform material decomposition to obtain the decomposition coefficient a of each base material i. i (r) including: Based on the high-energy region projection data L1 and the low-energy region projection data L2 of the calibration phantom, the n basis images f of the calibration phantom are calculated. n (r), the n base images f n The vector of (r) is represented as f(r); where r represents the voxel spatial position; By reconstructing the high-energy region projection data L1, the CT reconstruction data of the high-energy region is obtained, and the CT reconstruction data of the high-energy region is subjected to image segmentation processing to obtain the template image t of the base material i i (r); Based on the template image t i (r) and the vector representation f(r) are respectively subjected to the least square fitting calculation for each basis material i, and the polynomial coefficient vector c corresponding to the basis material i is obtained. i ; Based on the high-energy region projection data L′1 and low-energy region projection data L′2 of the sample to be tested, the n basis images f′ of the calibration phantom are calculated. n (r), the n basis images f′ n The vector of (r) is denoted as f′(r); Calculate the decomposition coefficient a of base substance i i (r) = c i ·f′(r).

5. The method according to claim 4, characterized in that The n basis images f of the calibration phantom n (r) = R -1 b n (L1,L2); where R -1 represents the inverse Radon transform, b n (L1, L2) represents a reconstructed image obtained based on the high-energy region projection data L1 and the low-energy region projection data L2.

6. The method according to claim 4, characterized in that The template image of the base material i 7. The method according to claim 4, characterized in that The polynomial coefficient vector corresponding to the base material i 8. The method according to claim 1, characterized in that The electron density ρ of the base substance i is obtained ei , and combined with the decomposition coefficient c of base substance i i , obtaining the rock component identification result of the sample to be tested, including: Get the electron density ρ of base material i ei ; Based on the decomposition coefficient a i (r) The electron density ρ of the base substance i ei Weighted calculation to obtain the electron density ρ of the sample to be detected e ; Calculate the electron density percentage of base substance i According to the electron density percentage f of each base substance i i and the effective atomic number Z of the base substance i i , get the effective atomic number Z of the sample to be tested eff ; Based on the electron density ρ e and the effective atomic number Z eff A table lookup is performed to obtain the rock component identification result of the sample to be tested.

9. The method according to any one of claims 1 to 8, characterized in that: After the calibration phantom and the sample to be tested are subjected to CT scanning and reconstruction in high and low energy zones respectively to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, and high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be tested, the method further includes: Median filtering is performed on the high energy area projection data L1, the low energy area projection data L2, the high energy area projection data L′1 and the low energy area projection data L′2 and / or, Performing mean filtering on the high energy region projection data L1, the low energy region projection data L2, the high energy region projection data L′1, and the low energy region projection data L′2; and / or, Perform pixel merging and projection frame merging on the high energy area projection data L1, the low energy area projection data L2, the high energy area projection data L′1 and the low energy area projection data L′2; and / or, The CT data reconstruction process is performed on the high energy region projection data L1, the low energy region projection data L2, the high energy region projection data L′1 and the low energy region projection data L′2 respectively through the energy spectrum CT reconstruction algorithm.

10. A rock component identification system based on spectral CT, characterized in that: The system comprises: A data acquisition module, used for selecting a calibration phantom of a sample to be detected, and performing CT scanning and reconstruction of the calibration phantom and the sample to be detected in high and low energy zones, respectively, to obtain high energy zone projection data L1 and low energy zone projection data L2 of the calibration phantom, and high energy zone projection data L′1 and low energy zone projection data L′2 of the sample to be detected; wherein the calibration phantom comprises a plurality of calibration samples i composed of a base material i, and the settings of the high energy zone and the low energy zone are obtained based on the attenuation capability of the sample to be detected to X-rays; The material decomposition module is used to perform material decomposition based on the high-energy region projection data L1 and low-energy region projection data L2 of the calibration phantom, and the high-energy region projection data L′1 and low-energy region projection data L′2 of the sample to be tested, and obtain the decomposition coefficient a of each base material i. i (r); Component identification module, used to obtain the electron density ρ of the base substance i ei and the effective atomic number Z i , and combined with the decomposition coefficient a of base substance i i (r), and obtain the electron density ρ of the sample to be tested e and effective atomic number Z eff , and then obtain the rock component identification result of the sample to be tested by looking up the table.

Images