Automatic rock stratum identification method and system based on multispectral remote sensing

Abstract

The invention relates to the technical field of artificial intelligence and data processing, in particular to an automatic rock stratum identification method and system based on multispectral remote sensing, and the method specifically comprises the following steps: employing a continuum removal-based adaptive normalization method for a collected multispectral remote sensing data sample; decomposing the normalized sample based on a mixed pixel spectrum of a physical mechanism of a rock mass, and calculating an abundance coefficient of each end member; calculating a first derivative value of each normalized sample in each wave band, and calculating a fusion feature value of each dimension by combining an abundance coefficient; constructing a rock stratum identification model based on deep learning, and inputting fusion feature values of all dimensions to carry out rock stratum category probability prediction; performing optimization training on the rock stratum identification model; newly collected data are processed and then input into the trained model, and a final recognition result is obtained. The accuracy, the stability and the geological rationality of the automatic rock stratum identification result can be improved.

Description

Technical Field

[0001] This invention relates to the fields of artificial intelligence and data processing technology, and in particular to an automatic rock strata identification method and system based on multispectral remote sensing. Background Technology

[0002] In geological exploration and mineral resource surveys, accurate identification of rock strata distribution is crucial for guiding field investigations, resource development, and geological modeling. Multispectral remote sensing technology, due to its wide coverage and non-contact acquisition capabilities, has become an important tool for rock strata identification. However, existing multispectral remote sensing-based rock strata identification methods still face numerous challenges in practical applications. First, different rock strata exhibit significant differences in background reflectance, and traditional min-max normalization methods struggle to preserve local absorption characteristics, making it difficult to distinguish between dark rocks and low-reflectance minerals. Second, remote sensing pixels often exhibit mixed reflectance phenomena of various lithologies, and existing linear mixing models, assuming fixed endmembers, cannot characterize the nonlinear mixing effects caused by variations in mineral composition. Third, conventional feature construction methods simply concatenate the original spectral and abundance information without considering their physical origins and expression differences, resulting in insufficient model recognition capabilities for key lithological features. Furthermore, most deep learning models rely solely on classification accuracy during training, neglecting the mineral mixing physical mechanisms that rock strata identification should satisfy, leading to a lack of model interpretability. Therefore, there is an urgent need for a novel method that integrates physical constraints and adaptive expression of spectral features to improve the accuracy, stability, and geological rationality of rock strata identification.

[0003] Chinese invention patent CN117830203A discloses a method for extracting the attitude of rock strata based on DEM and high-resolution optical remote sensing data in the field of basic geological integrated remote sensing research. The main contents include: automatically extracting contour lines of appropriate precision from DEM data of the study area using a computer; automatically extracting rock strata boundaries from high-resolution optical remote sensing data of the study area and performing manual correction; extracting the intersection data of rock strata boundaries and contour lines; automatically extracting the intersection data of a certain rock stratum with two adjacent contour lines having a number of intersections greater than or equal to 3; calculating the dip and dip angle parameters of different rock strata using the three-point coplanarity method; marking the attitude of each point along the boundary line of each rock stratum in the study area and outputting it in graphical or tabular form; efficiently and quickly acquiring rock strata attitude data over a large area; and outputting it in graphical or tabular form for further analysis and utilization.

[0004] Existing technologies have the following objective drawbacks: conventional min-max normalization ignores local absorption features, which can easily lead to the compression of dark rock features and a decrease in classification accuracy; conventional mixed pixel processing relies on fixed endmembers and cannot cope with the fluctuations in actual mineral composition and nonlinear mixing effects; conventional feature stitching methods do not distinguish the physical differences between derivatives and abundance, making it difficult for the model to focus on key lithological features; and using only cross-entropy loss ignores physical constraints, resulting in prediction results that lack physical consistency with the rock layer mixing mechanism.

[0005] Therefore, this invention proposes an automatic rock strata identification method and system based on multispectral remote sensing to solve the above problems. Summary of the Invention

[0006] This invention addresses the shortcomings of existing technologies by developing an automatic rock strata identification method based on multispectral remote sensing. This invention improves the accuracy, stability, and geological plausibility of the automatic rock strata identification results.

[0007] On the one hand, the technical solution of this invention to solve the technical problem is an automatic rock strata identification method based on multispectral remote sensing, as follows: An adaptive normalization method based on continuum removal was applied to the collected multispectral remote sensing data samples to obtain normalized samples. The collected multispectral remote sensing data samples were from open-source datasets or self-collected datasets, and all collected multispectral remote sensing data samples were labeled with lithological categories. The specific adaptive normalization method based on continuum removal is as follows: For the original reflectance of each sample in each band after normalization, the effects of atmospheric scattering and illumination changes are first eliminated by fitting the spectral envelope. Specifically, a cubic polynomial fitting function is used. The center wavelength of the band and the original reflectance value are input to obtain the fitted value. Then, the maximum value of the original reflectance value and the fitted value is taken as the continuum function value of each sample in each band, which is used to characterize the spectral envelope. Next, the local absorption depth is normalized based on the difference between the continuum function value and the minimum reflectance, resulting in the normalized reflectance value of each sample in each band, i.e., the normalized sample. Based on the physical mechanism of the rock mass, the mixed pixel spectrum is used to decompose the normalized sample. Each endmember of the decomposed sample is initialized, and the abundance coefficient of each endmember is calculated. The calculation process of the abundance coefficient of each endmember is as follows: For each decomposed endmember, a clustering method is used to generate a sample cluster. Based on the normalized samples, an indicator function is used to determine whether the current band matches the known mineral spectral library. The initial reflectance value of each endmember in each band is calculated. Then, with the initial reflectance values ​​of fixed endmembers, the abundance coefficient of each sample in each endmember is calculated by minimizing the L2 norm squared of the linear combination of the normalized reflectance value and each endmember, combined with the coefficient regularization term. Here, endmembers represent different rock strata types.

[0008] The first derivative of each normalized sample in each band is calculated, and then combined with the abundance coefficient and attention-weighted feature fusion to obtain multi-dimensional fused feature values. The calculation process of the multi-dimensional fused feature values ​​is as follows: For each band of each sample, the difference in normalized reflectance values ​​between adjacent bands is calculated, and the first derivative value is obtained by combining the difference in center wavelength between adjacent bands with the derivative scaling factor. Then, for each sample, the attention weight of the first derivative value and the attention weight of the abundance coefficient in each dimension are calculated. Then, the weighted first derivative feature value and the weighted abundance coefficient are summed to obtain the fusion feature value in each dimension. The dimensions include the total number of bands and the total number of endmembers.

[0009] A deep learning-based rock strata identification model is constructed. Based on the model's spectral adaptive activation function and multi-scale residual block architecture, the fused feature values ​​of each dimension are used to predict the probability of rock strata categories. The specific operations in the deep learning-based rock strata identification model are as follows: The fused feature values ​​of each dimension are input into the model. These fused feature values ​​are then input into the linear layer. Combined with the projection weight matrix and bias vector of the input layer, the initial activation values ​​of each sample at each hidden layer node are obtained. These initial activation values ​​are then input into the multi-scale residual block architecture. Specifically, the activation function value of the previous layer residual block is processed by the spectral adaptive activation function and combined with the activation function values ​​of the previous two layers residual blocks multiplied by the cross-layer skip connection weight matrix to obtain the activation function value of the current layer residual block. Through layer-by-layer calculation, the activation function value of the final layer residual block is obtained. This activation function value is then combined with the output layer weight matrix and the output layer bias vector, and the predicted probability of the rock stratum category is calculated using the Softmax function.

[0010] By incorporating the physical constraints of spectral decomposition into the loss function, the rock strata identification model is optimized, and the model parameters are updated using an adaptive learning rate mechanism of parameter grouping, resulting in a well-trained rock strata identification model. The physical constraints of spectral decomposition are incorporated into the loss function, and the total loss is calculated by combining the cross-entropy classification loss term and the physical regularization term. Specifically, the cross-entropy loss term is calculated by measuring the deviation between the predicted probability of the rock stratum category obtained by the quantification model and the actual category. For each sample, the L2 norm square of the linear combination of the normalized reflectance value, the predicted abundance coefficient, and the endmember reflectance value is calculated, and the physical regularization term is obtained by combining the absolute value of the abundance coefficient. A parameter grouping adaptive learning rate mechanism is adopted, which groups the trainable parameters according to their physical role in the model. A dynamically adjusted cosine annealing learning rate is applied to the weight parameters that control feature transformation to stabilize the convergence of deep networks. At the same time, a fixed small learning rate is applied to the bias parameters that control the offset to suppress overfitting. The training process of the model is as follows: The collected data is forward-propagated, processed, and then input into the rock strata identification model to obtain the predicted probability of the rock strata category. The loss function is then calculated, and the model parameters are updated through backpropagation. The model is iteratively optimized by the Adam optimizer until the iteration termination condition is met, at which point the optimization ends and a trained rock strata identification model is obtained.

[0011] The newly acquired multispectral remote sensing data is processed and then input into the trained model to obtain the final rock stratum category probability prediction result. The lithology category corresponding to the maximum probability is taken as the identification result.

[0012] On the other hand, the present invention also provides an automatic rock strata identification system based on multispectral remote sensing, including modules for executing processing instructions for each step in the automatic rock strata identification method based on multispectral remote sensing, including the following modules: Data collection module: Collects multispectral remote sensing data samples; Data preprocessing module: performs normalization, mixed pixel decomposition, and special fusion operations on the collected samples; Rock strata identification module: Input the preprocessed sample into the rock strata identification module to predict the lithology of the sample; Optimize the training module: Optimize and update the parameters of the rock strata identification module through the model's total loss function and adaptive learning rate mechanism to obtain a well-trained rock strata identification module; Recognition result output module: After the sample to be identified is preprocessed, it is input into the trained rock stratum recognition module and the final recognition result is output.

[0013] The effects described in the invention are merely those of the embodiments, and not all the effects of the invention. The above technical solutions have the following advantages or beneficial effects: This invention takes into account that multispectral remote sensing data is presented in the form of one-dimensional spectral curves, which are composed of reflectance values ​​of multiple bands. It has the characteristics of high dimensionality, strong noise interference, and scale inconsistency caused by illumination and atmospheric conditions. The invention adopts the continuum removal normalization method. Through the calculation of continuum function, the difference in background reflectance can be eliminated, that is, the influence of atmospheric scattering and illumination changes can be eliminated, the integrity of absorption characteristics is preserved, the key absorption characteristics of lithology are highlighted, and the ability to identify dark rocks and low reflectance minerals is improved. A physically constrained hybrid pixel decomposition model is constructed to process the mixed reflection signals of multiple rock layers contained in a single pixel. The decomposed endmembers are initialized and combined with key absorption bands of the mineral spectral library to enhance the geological interpretability of the decomposition results. This can avoid the problem that the endmember spectra will change dynamically due to the variation of mineral composition, resulting in insufficient decomposition accuracy and failure to capture nonlinear mixing effects. This invention employs a dual-path attention fusion mechanism of first derivative and abundance coefficient. By calculating the first derivative of the spectrum, the lithologically specific absorption edge features can be highlighted. Then, the weighted derivative features and the weighted abundance coefficients are summed to adaptively enhance highly discriminative features and avoid feature dimension redundancy and physical semantic confusion. By jointly optimizing the classification loss term and the physical regularization term, spectral reconstruction error and abundance sparsity are incorporated into the training, forcing the model's prediction results to conform to the physical mechanism of rock and mineral mixing. This improves the model's physical consistency and classification robustness, and avoids the problem of model prediction results being disconnected from lithological physical mechanisms, resulting in decomposition and prediction fragmentation.

[0014] In summary, this invention can improve the accuracy, stability, and geological rationality of automatic rock strata identification results. Attached Figure Description

[0015] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used together with the embodiments of the invention to explain the invention and do not constitute a limitation thereof.

[0016] Figure 1 This is a schematic diagram of the method flow of the present invention.

[0017] Figure 2 This is a schematic diagram of the reflectance of the original spectral lines and the continuous envelope.

[0018] Figure 3 This is a schematic diagram of the spectral reflectance after continuous normalization.

[0019] Figure 4 This is a schematic diagram illustrating the impact of the normalization method of this invention on the accuracy of rock strata identification compared to conventional normalization methods.

[0020] Figure 5 This is a schematic diagram of the spectral lines of the physically constrained endmember.

[0021] Figure 6 This is a schematic diagram of the abundance distribution of mixed pixels.

[0022] Figure 7 This is a schematic diagram of the attention weights of the spectral derivative.

[0023] Figure 8 This is a diagram illustrating the attention weights for abundance coefficients.

[0024] Figure 9 This is a schematic diagram of the distribution of rock strata in the feature space. Detailed Implementation

[0025] To clearly illustrate the technical features of this solution, the invention will be described in detail below through specific implementation methods and in conjunction with the accompanying drawings.

[0026] Example 1 like Figure 1 As shown, an automatic rock strata identification method based on multispectral remote sensing includes the following steps: S1. Collected multispectral remote sensing data samples, and label the collected multispectral remote sensing data samples with lithological category labels; S2. An adaptive normalization method based on continuum removal is used to process the collected multispectral remote sensing data samples to obtain normalized samples. S3. Based on the physical mechanism of the rock mass, the mixed pixel spectrum is decomposed into normalized samples, each endmember of the decomposed sample is initialized, and the abundance coefficient of each endmember is calculated. S4. Calculate the first derivative value of each normalized sample in each band, and combine it with the abundance coefficient memory attention weighted feature fusion to obtain multi-dimensional fusion feature values. S5. Construct a rock strata identification model based on deep learning, and predict the probability of rock strata category based on the spectral adaptive activation function and multi-scale residual block architecture of the model and the fused feature values ​​of each dimension. S6. The rock strata identification model is optimized by incorporating the physical constraints of spectral decomposition into the loss function, and the model parameters are updated by an adaptive learning rate mechanism of parameter grouping to obtain a trained rock strata identification model. S7. After processing, the newly acquired multispectral remote sensing data is input into the trained model to obtain the final rock stratum category probability prediction result. The lithology category corresponding to the maximum probability is taken as the identification result.

[0027] In a specific implementation, S1 is as follows: The collected data can be from open-source standard datasets or self-collected datasets; The open-source standard dataset can be the Landsat-8-OLI open-source dataset, which contains seven visible to shortwave infrared bands, with a wavelength range of 430–2300 nanometers and a spatial resolution of 30 meters, meeting the requirements of rock strata identification for spectral resolution and spatial coverage. This open-source standard dataset already has labels, including labels for igneous rocks such as granites and basalts, sedimentary rocks such as carbonate rocks and clastic rocks, and metamorphic rocks such as gneiss and marble. If the dataset is self-collected, it must be obtained from a standardized multispectral remote sensing image dataset covering the target area via satellite or airborne platform. The data source must meet the following core requirements: 1) Band integrity: It must include key bands that are sensitive to lithology (such as Band 1–Band 7), with particular emphasis on shortwave infrared bands, as they are sensitive to the absorption characteristics of silicate and carbonate minerals; 2) Radiometric calibration and atmospheric correction: The raw data needs to be converted into apparent reflectance through radiometric calibration to eliminate atmospheric scattering effects and ensure the accuracy of the physical meaning of reflectance; 3) Regional representativeness: The data should cover typical exposed rock areas, such as mountains and mining areas, and avoid invalid pixels with vegetation coverage greater than 20% or cloud pollution greater than 5%. Self-collected datasets need to be labeled. The labeling method is manual labeling, specifically based on field-measured geological maps and mineral sample spectral libraries. Pure rock strata areas are delineated on the images as training samples. The labeling categories cover the main lithologies, such as granite, basalt, and limestone, with no less than 500 pixels labeled for each category.

[0028] In a specific implementation, S2 is as follows: Since multispectral remote sensing data is presented in the form of one-dimensional spectral curves, consisting of reflectance values ​​of multiple bands, it has the characteristics of high dimensionality, strong noise interference, and scale inconsistency caused by illumination and atmospheric conditions. Conventional methods often use global min-max normalization to process such data, such as scaling the data to the 0 to 1 range. This method ignores the local absorption characteristics of the spectral curve and its physical meaning, and cannot eliminate the systematic bias caused by atmospheric scattering. Therefore, this invention employs an adaptive normalization method based on continuum removal. First, it eliminates the influence of atmospheric scattering and illumination variations by fitting the spectral envelope. Then, it normalizes the local absorption depth based on the difference between the envelope and the minimum reflectance, thereby suppressing noise interference and scale differences while preserving key spectral features crucial for lithological identification. The specific steps are as follows: S2.1 Calculation of Continuum Functions: For each sample and each band's original reflectance value, a cubic polynomial fitting function is used. The center wavelength of the band and the original reflectance value are input to obtain the fitted value. Then, the maximum value between the original reflectance value and the fitted value is taken as the continuum function value, used to characterize the spectral envelope. This eliminates the influence of atmospheric scattering and illumination variations, and ensures that the envelope is always higher than the original spectrum to avoid underestimating the absorption characteristics. By fitting a cubic polynomial envelope higher than the original spectrum, atmospheric scattering and illumination variations are eliminated, preserving the integrity of the absorption characteristics. This is expressed as: , In the formula, To indicate the first The sample at the th The continuum function value for each band, i.e. the calculated value of the spectral envelope in that band, is used to eliminate the effects of atmospheric scattering and changes in illumination. The sample index identifies different remote sensing pixels or observation points; For band indexing, identify different bands in multispectral data; Indicates the first The sample at the th The original reflectance values ​​for each band are multispectral input data obtained through direct observation; Indicates the first The center wavelength of each band is directly read from the remote sensing data metadata and determined by the sensor parameters; This represents a polynomial fitting function, specifically using a cubic polynomial for a given wavelength. and the corresponding original reflectivity value Perform fitting; This indicates the order of the polynomial fit, i.e., the highest degree of the polynomial is 3; This represents the function that takes the maximum value. It should be noted that spectral absorption characteristics manifest as reflectance dips, and the envelope must cover the upper boundary of these dips to quantify the absorption depth. If the fitted curve is lower than the original data, it will lead to errors in the calculation of the absorption depth. The term ensures that the obtained continuous function value is fitted. Always at or above the original reflectance value To ensure the envelope lies above the original spectrum and avoid underestimating absorption characteristics, polynomial fitting is achieved using the least squares method, aiming to minimize the sum of squared residuals between the observed values ​​and the fitted curve. Specifically, the input wavelength vector... and corresponding reflectivity value Construct the design matrix For example, suppose nm, The design matrix obtained by fitting a cubic polynomial is: Then solve for the sparse vectors. , sparse vector The first element, sparse vector The second element, sparse vector The third element, sparse vector The fourth element, then, outputs the fitted value. For example, setting ; S2.2 Absorption depth normalization: Using the continuum function value and the minimum reflectance value of the sample across all bands, the difference between the original reflectance value and the continuum function value is calculated. Based on this difference, a normalized reflectance value is obtained, highlighting local absorption characteristics and eliminating global scale differences, expressed as: , In the formula, The normalized first The sample at the th The reflectance values ​​of each band represent the spectral characteristics after eliminating global scale differences and highlighting local absorption features; To indicate the first Each sample in all bands The minimum reflectance value in the spectrum is used as the reference point for normalization, which usually corresponds to the deepest absorption feature in the spectrum of that sample. This is a band index, distinct from 'b', used to traverse all bands; To represent the continuum normalization scaling factor, the default value is set to 1. Its function is to control the range of the final normalized value, avoid the absolute value of the normalized result being too large due to the denominator being too small, and thus prevent the dynamic range of the feature from being over-compressed. It should be noted that, unlike conventional normalization which uses global range or statistics, this method relies on local absorption characteristics rather than absolute reflectance intensity. The term represents the difference between the continuum value and the minimum reflectance, based on a sample-level minimum reflectance benchmark. Physically equivalent to mapping absorption depth to the scale of the current sample's potential maximum absorption capacity, it enables comparability of absorption intensity across samples. For example, the strong absorption of iron ore in the 700-900 nm range and the weak absorption of gypsum in the 1500 nm range are compressed to similar numerical ranges after normalization, preventing high-reflectance minerals from masking the key absorption of low-reflectance minerals. More importantly, it addresses the differences in background reflectance between different lithologies, such as dark basalt and light sandstone. Conventional fixed-range normalization would cause all band values ​​of dark rocks to approach 0, losing their separability. As a dynamic benchmark, it can automatically adapt to the differences in background reflection of different lithologies.

[0029] In a specific implementation, S3 is as follows: The normalized spectral curve still faces the problem of mixed pixels, that is, a single pixel contains mixed reflection signals of multiple rock layers. Conventional methods usually use a linear mixing model to process it. This model assumes that the end-member spectrum is fixed. However, in practical applications, the end-member spectrum will change dynamically due to the variation of mineral composition, resulting in insufficient decomposition accuracy, inability to capture nonlinear mixing effects, and the decomposition results are inconsistent with the absorption characteristics of rocks and minerals. Therefore, this invention combines lithological physical mechanisms and adopts a decomposition framework with endmember adaptive updating and abundance sparsity constraints. Through an iterative optimization process, it simultaneously learns endmember spectra and abundance coefficients to ensure that the decomposition results conform to the physical characteristics of rock and mineral absorption. The specific steps are as follows: S3.1 Initialization of physical constraint endpoints: Based on the normalized reflectance values, sample clusters are generated using clustering methods. For each endmember, its initial reflectance value is calculated as the average of the normalized reflectance values ​​across the key absorption band set. This average is used to characterize the spectral characteristics of pure rock layers, ensuring that the endmembers conform to the true mineral absorption characteristics and avoiding random initialization bias. This is expressed as: , In the formula, Indicates the first The terminal in the first The initial reflectance values ​​of each band characterize the spectral features of the pure rock layer; Endmember indexes are used to identify different rock strata types; Indicates the first The sample clusters corresponding to each endmember are all clustered using k-means clustering. of Genesis, where the cluster center corresponds to a physically pure rock layer; This indicates an indicator function that outputs 1 when the condition is met and 0 otherwise. Indicates the first The key absorption band set of rock-like strata is predefined through a priori mineral spectral library to characterize physical constraints; To represent sample clusters The number of samples in the sample; It should be noted that endmembers represent the spectral characteristics of pure rock types, indicating the reflectance curve of a single lithology under ideal conditions; key absorption band sets... The mineral spectral library is determined by matching characteristic absorption bands in a known mineral spectral library. This library can be the open-source ASTER (Advanced-Spaceborne-Thermal-Emission-and-Reflection-Radiometer). For example, if the... The rock strata are composed of hematite, and mineral spectral libraries show that they contain hematite. nm has strong absorption, then Includes band index satisfy nm, through physical constraints, ensures that the endmember initialization conforms to the true absorption characteristics of minerals, avoiding randomness; conventional clustering initialization applies equal weight to all bands, but the spectral characteristics of minerals are concentrated only in specific bands, such as epidote in the 2300 nm iron absorption band, through an indicator function. Limit the calculation to a predefined set of key absorption bands. A physical-driven feature filter is constructed to achieve adaptive suppression of noise and redundant bands. In the near-infrared band where there are no feature changes, the algorithm automatically ignores reflectance fluctuations in this region and extracts only spectral segments strongly correlated with mineral composition. Furthermore, because... By encoding domain knowledge into mathematical constraints from a prior knowledge base, endmember initialization is transformed from a purely data-driven approach to a physical-guided hybrid paradigm. When two minerals have similar reflections in non-characteristic bands, conventional methods may generate confusing endmembers, while band-selective averaging can force the separation of mineral-specific spectral morphology. S3.2, Sparse Constraint Abundance Estimation With fixed endmember reflectance values, the abundance coefficient is calculated for each sample. The abundance coefficient is obtained by minimizing the L2 norm squared of the linear combination of the normalized reflectance value and the endmembers, combined with a sparse regularization term. This reflects the spatial discontinuity of the rock strata and suppresses excessive mixing effects, and is expressed as: , In the formula, Indicates the first In the nth sample The abundance coefficient of each endmember, representing the proportion of each rock layer in the mixed pixel, must satisfy the following conditions: To avoid negative values ​​leading to non-physical interpretations, and simultaneously satisfying that the sum of the abundance coefficients of all endmembers is 1, i.e. ; This indicates that the function to be minimized is the function whose optimization objective is to minimize the subsequent expression. The total number of endmembers; Indicates the first The terminal in the first Reflectance values ​​for each band; Represents the L2 norm; This represents the sparsity regularization coefficient, used to enhance the sparsity of abundance coefficients, such as... =0.01; It should be noted that, since pixels are usually dominated by a few rock layers, The term represents the abundance coefficient. Apply sparsity constraints to force The abundance coefficients of most endmembers are approximately 0 to reflect the spatial discontinuity of rock strata distribution, avoid over-mixing, and improve decomposition accuracy; initial reflectance value It is the initial endmember generated by clustering, and the reflectance value. This is the updated value after iterative optimization; because the abundance needs to satisfy... and Minimizing abundance forces non-dominant endmember coefficients to approach zero, while dominant endmembers are automatically assigned higher weights. The term, as a regularization term, represents the sparsity pressure applied to the abundance and sum. It works synergistically with the non-negative abundance constraint to effectively address the spatial sparsity of rock layer distribution and suppress excessive mixing effects. Furthermore, the objective function is solved through optimization algorithms. To obtain the optimal abundance coefficient and reflectivity value The optimization algorithm can be gradient descent or non-negative least squares.

[0030] In a specific implementation, S4 is as follows: While the abundance coefficients obtained from mixed pixel decomposition reduce the impact of the mixing effect, they lose detailed information from the original spectrum. Conventional methods often simply concatenate the original spectrum and abundance coefficients as features, which can easily lead to feature dimension explosion and redundancy, making it difficult for the model to focus on key bands that are crucial for lithological discrimination. Therefore, this invention adopts a multi-level feature enhancement mechanism. First, the spectral derivative is calculated to highlight the lithologically specific absorption edge features. Then, an attention mechanism is used to weightedly fuse the derivative features with the abundance coefficients, thereby adaptively enhancing lithological features with high discriminative power. The specific steps are as follows: S4.1 Calculation of the first derivative of the spectrum: For each sample and each band, the difference in normalized reflectance values ​​between adjacent bands is calculated. This difference is then combined with the difference in center wavelength between adjacent bands and the derivative scaling factor to obtain the first derivative value. This first derivative value is used to characterize the local slope of the spectral curve, amplify the rate of change in reflectance between bands, highlight lithologically specific absorption edge characteristics, and suppress smoothing noise. This is expressed as: , In the formula, For the first The sample at the th The first derivative values ​​of each band characterize the local slope of the spectral curve and are used to identify lithology-specific absorption valleys. The normalized first The sample at the th The reflectance values ​​of each band represent the spectral characteristics after eliminating global scale differences and highlighting local absorption features; For the first The center wavelength of each band; For the first The center wavelength of each band is directly read from the remote sensing data metadata and determined by the sensor parameters. The unit is nanometers; The derivative scaling factor is used to amplify weak absorbed signals while suppressing noise interference. =1.2; It should be noted that conventional derivative calculations directly use the reflectivity difference between adjacent bands, ignoring the slope distortion caused by uneven band spacing. The term represents the wavelength difference normalization operation, and the derivative scaling factor. Forming a coupling relationship, wavelength difference As the denominator, it physically achieves comparable first derivatives across different bands. The scaling factor can amplify the discriminative power of weak absorption edges. Essentially, The model's sensitivity to spectral variation rates is parameterized and can be dynamically adjusted according to mineral type. For example, in the gradually varying absorption valley of chlorite at 2250 nm, the conventional derivative may be overwhelmed by noise, while setting... The hyperparameters are made significant by linear amplification. At the same time, because the derivative itself is small, it still remains within a reasonable range after amplification. It can form obvious peaks and valleys in the derivative space, thus widening the difference with the surrounding rock. S4.2 Attention-weighted feature fusion: For each sample, the attention weights of the derivative features and the abundance coefficients are calculated. Then, the weighted derivative features and the weighted abundance coefficients are summed to obtain the fused feature value, which is used to adaptively enhance the high-discriminative lithological features, as shown below: , In the formula, For the first The sample at the th The fused feature values ​​of each feature dimension The total dimension of the fused features is calculated as follows: ; Total number of bands; For the first The sample at the th The derivative feature attention weights for each band are calculated as follows: ; It is a natural exponential function; The attention score function is a derivative feature, specifically a learnable linear transformation; For the first The sample at the th The first derivative value of each band; To distinguish it from the band index of b, it is used to traverse all bands; For the first The sample at the th The abundance coefficients and attention weights for each band are calculated as follows: ; The attention score function is the abundance coefficient, specifically a learnable linear transformation; Indicates the first In the nth sample Abundance coefficients of each endmember; To distinguish it from the end-member index of k, it is used to traverse all end-members; Both the attention score function for first-order derivative features and the attention score function for abundance coefficients are learnable linear transformation functions. For example, the implementation method is represented as follows: ,in, The weights of the attention scoring function, The biases of the attention scoring function are all trainable parameters, calculated as follows: Indicates attention score; It should be noted that the attention weights of derivative features Abundance coefficients and attention weights The calculation methods are all based on Softmax normalization, which, while ensuring a weight sum of 1, focuses on high-discriminative bands and strengthens the representation of important strata types. Conventional splicing or weighted averaging cannot distinguish the physical differences between the derivative characteristics of local band variations and the abundance coefficients of global lithology proportions. This is addressed by integrating the characteristic values. During the computation process, dual-path attention achieves a synergistic technical effect on heterogeneous features, through an independent attention mechanism. and The calculation constructs an adaptive weight allocator for feature types, achieving complementary enhancement of discriminative information. For example, when the derivative of a certain band caused by the 700nm absorption edge of hematite is significant, The weight is automatically increased, and when the abundance of a certain end-member is prominent due to the high proportion of granite, This enhances its contribution. The dual-path design avoids feature crosstalk and focuses on the derivative path of the absorption edge of altered minerals and the abundance path of the host rock proportion.

[0031] In a specific implementation, S5 is as follows: Because the fused feature values ​​are high-dimensional and highly nonlinear, conventional fully connected networks use standard activation functions such as ReLU for processing. However, such functions are difficult to adapt to the continuous and smooth characteristics of the spectrum, easily leading to gradient vanishing and insufficient feature representation, thus making it difficult to capture the subtle spectral variation patterns required for lithology identification. Therefore, this invention adopts a spectral adaptive activation function and constructs a multi-scale residual block architecture, which can effectively preserve the continuous information of the spectrum in deep networks, while enhancing the feature representation capability of the model. The specific steps are as follows: S5.1 Input Layer and Feature Projection: The fused feature values ​​are input into the linear layer, and combined with the input projection weight matrix and the input layer bias vector to obtain the initial activation values. This achieves linear compression of the fused feature dimension, retains key lithological information, and reduces subsequent computational complexity, as shown below: , In the formula, For the first The sample at the th The initial activation values ​​of the hidden layer nodes; The index of the hidden layer node identifies the different feature dimensions after dimensionality reduction. , The reduced-dimensionality feature representation provides input for subsequent residual blocks; The input projection weight matrix consists of trainable parameters that are learned through training to compress the feature mapping. The input layer bias vector is a trainable parameter used to adjust the feature offset. For the hidden layer dimension, for example, set to To balance expressive power and computational efficiency; S5.2 Multi-scale residual block connection: For each residual block, the activation value of the previous layer is applied with a spectral adaptive activation function, and then multiplied by the activation values ​​of the previous two layers by the cross-layer skip connection weight matrix to obtain the activation output value of the current layer. This maintains spectral continuity and alleviates gradient vanishing in deep networks, and is expressed as: , In the formula, For the first The sample at the th Layer residual block The activation output value of each hidden layer node; For the layer index of the residual block, ; This represents the total number of residual blocks; For the first The weight matrix of the residual block is a trainable parameter used to perform feature transformation; For the first The sample at the th Layer residual block The output activation values ​​of each hidden layer node; For the first The sample at the th Layer residual block The output activation values ​​of each hidden layer node; For the first The bias vector of the layer residual block is a trainable parameter; For the cross-layer skip connection weight matrix, the first layer... The characteristics of the layer residual block are directly transferred to Layer residual blocks are trainable parameters; Let be the spectral adaptive activation function, and let its input be... The calculation method is expressed as follows: ; The input values ​​for the spectral adaptive activation function are: item; For hyperbolic tangent function, the output range is... It provides fundamental nonlinearity; This is a scaling factor used to adjust the strength of the nonlinear term. ; This is the attenuation coefficient, used to control the rate of Gaussian decay. (Setting...) ; It should be noted that the conventional ReLU activation function faces a trade-off between preserving continuity and detecting abrupt changes in spectral data processing. The hard truncation of the ReLU activation function disrupts continuity, while the Softplus activation function provides excessive smoothing. The term provides basic nonlinearity, with an output range of [value missing]. , The term represents the Gaussian decay term, which is derived from... Control the decay rate, The term makes the function approximately linear in the smooth region of the spectrum and enhances nonlinearity in the region of abrupt changes, which fits the continuous variation characteristics of the spectrum, such as when Hours characterize the spectral smooth region. The spectral adaptive activation function is approximately linear, preserving spectral continuity when Large time-characteristic abrupt change region, , The term plays a dominant role, enhancing nonlinearity to capture absorption edges. Through adaptive smoothing-nonlinear switching, it improves sensitivity to subtle spectral changes. For example, in the quartz-feldspar mixed spectrum, the gradient of the weak absorption valley of feldspar may vanish under the ReLU activation function, while the spectral adaptive activation function maintains near-linear propagation with small inputs, enabling effective backpropagation of the gradient. Simultaneously, the parameters... and By controlling the nonlinear intensity and decay rate respectively, an adjustable spectral response surface can be formed; S5.3, Probability Prediction of Rock Strata Types: The activation values ​​of the final residual block are combined with the output layer weight matrix and the output layer bias vector, and the softmax function is applied to calculate the predicted probability of the rock strata category to support the classification decision, as shown below: , In the formula, For the first The sample belongs to the first The predicted probabilities of each category, ; This represents the total number of rock strata categories. The softmax function converts the input vector into a probability distribution. is the output layer weight matrix, which consists of trainable parameters; The output layer bias vector is a trainable parameter.

[0032] In a specific implementation, S6 is as follows: Conventional methods only use the cross-entropy loss function, focusing solely on classification accuracy and neglecting the physical consistency constraints that the spectral decomposition process should satisfy, such as the sparsity of abundance coefficients and the minimization of endmember reconstruction errors. This can easily lead to a disconnect between the model's prediction results and the physical mechanism of lithology, resulting in a separation between decomposition and prediction, and reducing the physical interpretability of rock strata identification. Therefore, this invention incorporates the physical constraints of spectral decomposition into the loss function, and forces the model's prediction results to conform to the physical mechanism of rock and mineral mixing by jointly optimizing the classification loss term and the physical regularization term. S6.1.1 Constructing a framework for the joint loss function: The total loss function is obtained by combining the cross-entropy classification loss term and the physical regularization term. The classification loss and physical regularization term are jointly optimized to balance model accuracy and the physical consistency of spectral decomposition, expressed as: , In the formula, The total loss function guides the end-to-end training of the model; The cross-entropy classification loss term quantifies the deviation between the predicted lithology category and the true label. As a physical regularization term, it constrains the network output to conform to the physical mechanism of spectral decomposition; The physical regularization weight coefficient controls the strength of physical constraints. =0.05; S6.1.2, Calculation of Cross-Entropy Classification Loss: For all training samples and all categories, the average of the logarithmic products of the true labels and predicted probabilities is calculated to obtain the cross-entropy loss term. This term quantifies the deviation between the predicted probabilities and the true labels, driving the model to accurately identify rock strata categories, and is expressed as: , In the formula, The total number of training samples; This represents the total number of rock strata categories. For the first The sample belongs to the first One-hot encoding of each category, when the sample Belongs to the The value is 1 if there is a specific category, otherwise it is 0; S6.1.3 Construction of Physical Regularization Terms: For each sample, the L2 norm squared of the linear combination of the normalized reflectance value, the predicted abundance coefficient, and the endmember reflectance value is calculated. Combined with the absolute value of the abundance coefficient, a physical regularization term is obtained. Through reconstruction error and abundance sparsity constraints, the network prediction is forced to conform to the mineral mixing physical mechanism, expressed as: , In the formula, The network predicts the first [number]th ... In the nth sample Abundance coefficients of each endmember; For absolute value functions, sparsity constraints are applied to the abundance coefficients; The auxiliary layer can be a fully connected layer, with the input fused feature values. Output the predicted abundance coefficient During training, the physical regularization term Supervision abundance coefficient The generation process is forced to approximate the true abundance distribution, and a Softmax function is applied after the auxiliary layer to satisfy... ; It should be noted that, through reflectivity values... With abundance coefficient The hidden layer encoding of the forced classification network yields physically interpretable spectral decomposition results. The term characterizes the reconstruction error, requiring the network to predict abundance energy that can accurately reconstruct the input spectrum. Ensuring that the predicted abundance energy can accurately reconstruct the original spectrum forces physical consistency decomposition, which is equivalent to adding a self-supervised signal. The term represents the abundance sparsity. Since a few rock layers dominate in a pixel, the abundance sparsity can reflect the sparsity of rock layer distribution. During model training, the reflectance value... As a learnable parameter, the mineral mixing law is backpropagated through the physical regularization term, so that the hidden features of the network carry lithological physical interpretation. For example, when the classifier confuses gabbro and diorite, it can punish the reconstruction deviation of the endmember spectra of the two, driving the network to learn more discriminative mineral absorption feature expressions. When training models with high-dimensional spectral data, it is easy to get stuck in local optima. Conventional Adam optimizers use a fixed learning rate update strategy for all parameters, which cannot distinguish the different roles of weight parameters and bias parameters in the model and their update requirements. This leads to unstable convergence of endmember-related parameters and reduces the model's ability to capture subtle differences in lithology. Therefore, this invention adopts a parameter grouping adaptive learning rate mechanism, which groups parameters according to their physical roles in the model. A dynamically adjusted cosine annealing learning rate is applied to the weight parameters that control feature transformation to stabilize the convergence of deep networks. At the same time, a fixed small learning rate is applied to the bias parameters that control the offset to suppress overfitting. S6.2.1 Define parameter update rules: Using the Adam optimizer, the trainable parameters are updated based on the bias correction terms from the first-order moment estimates, the bias correction terms from the second-order moment estimates, and the group-adaptive learning rate, as follows: , In the formula, For the first The trainable parameters for the next iteration are in set form; For the first The trainable parameters for the next iteration are in set form; For the first The bias correction term for the first-order moment estimate in the next iteration is calculated as follows: ; For the first The first moment of the gradient in the next iteration; The exponential decay rate of the first moment of the gradient Power; For the first The bias correction term for the second-order moment estimate in the next iteration is calculated as follows: ; For the first The second moment of the gradient in the next iteration; The exponential decay rate of the second moment of the gradient Power; To set a numerical stability constant and prevent the denominator from being zero, ; Adaptive learning rate for grouping; It should be noted that the exponential decay rate of the first moment of the gradient... and the exponential decay rate of the second moment of the gradient These are hyperparameters of the Adam optimizer, used to eliminate biases in the initial moment estimation during optimization. They are typically set to... , ; S6.2.2, Group Learning Rate Setting: Based on the physical effects of the parameters, the parameters are divided into two groups. For the weight parameter group, the base learning rate is multiplied by the cosine function value, and cosine annealing is applied to accelerate convergence. For the bias parameter group, the base learning rate is multiplied by a fixed value, and a small fixed learning rate is used to suppress overfitting. This can be expressed as: , In the formula, It is a trainable set of weight parameters that controls feature transformation and information transmission; A set of trainable bias parameters controls the offset of the decision boundary; Set the base learning rate as the update benchmark. =0.001; This represents the current training iteration number. This represents the maximum number of training iterations. It is a cosine function. The cosine annealing strategy is implemented to achieve the weight parameters. These are trainable parameters, in set form; Pi is used for calculation. ; It should be noted that in the initial stage of training, , This accelerates weight convergence; in the later stages of training, , Fine-tuning the parameters stabilizes the model; while the bias parameters are fixed with a small learning rate. To avoid prediction bias caused by overfitting; S6.3 The model training process is as follows: The collected data is forward propagated, processed, and then input into the rock strata identification model to obtain the predicted probability of rock strata categories. Then the loss function is calculated, and the model parameters are updated through back propagation. The model is iteratively optimized through the Adam optimizer until the iteration termination condition is met, at which point the optimization ends and a trained rock strata identification model is obtained. The iteration termination condition is as follows: After each training round, evaluate the accuracy of the validation set. If the validation loss does not decrease or the accuracy fluctuation is less than 0.5% for 10 consecutive rounds, then terminate the training and save the optimal model parameters to avoid overfitting. The maximum number of training iterations is set to 500 rounds, and the base learning rate is fixed at 0.001.

[0033] In a specific implementation, S7 is as follows: A new multispectral remote sensing image is input, and continuous normalization is performed pixel by pixel to eliminate illumination / atmospheric bias. Then, physical constraint hybrid pixel decomposition is performed, and the endmember abundance coefficients are output. Then, attention-weighted feature fusion is performed to fuse derivative features and abundance coefficients. The fused features are input into a trained multiscale residual network, and after linear projection and residual block nonlinear transformation, the lithology category probability distribution of each pixel is output. The lithology category corresponding to the maximum probability is taken as the recognition result of that pixel.

[0034] Example 2 An automatic rock strata identification system based on multispectral remote sensing includes modules for executing processing instructions for each step of the automatic rock strata identification method based on multispectral remote sensing, including the following modules: Data collection module: Collects multispectral remote sensing data samples; Data preprocessing module: performs normalization, mixed pixel decomposition, and special fusion operations on the collected samples; Rock strata identification module: Input the preprocessed sample into the rock strata identification module to predict the lithology of the sample; Optimize the training module: Optimize and update the parameters of the rock strata identification module through the model's total loss function and adaptive learning rate mechanism to obtain a well-trained rock strata identification module; Recognition result output module: After the sample to be identified is preprocessed, it is input into the trained rock stratum recognition module and the final recognition result is output.

[0035] Example 3 like Figures 2 to 3 The figures shown are schematic diagrams of the reflectance of the original spectral curve and the envelope of the continuum, and the reflectance of the spectrum after the continuum is normalized. Figure 2 The original reflectance characteristics of different rock layers (granite, basalt, limestone, gneiss, and sandstone) in Landsat-8-OLI multispectral bands (seven bands from Band 1 to Band 7) are shown. The dashed lines represent the continuum envelope fitted by a cubic polynomial. This envelope is always above the original spectrum to accurately capture spectral absorption characteristics. This processing can effectively eliminate the influence of atmospheric scattering and illumination changes, and provide a physically accurate benchmark for subsequent normalization. Figure 3 The spectral curves are presented after adaptive normalization, which is based on the difference between the continuum envelope and the minimum reflectance. This significantly highlights the key absorption characteristics of different rock layers, such as the strong absorption valley of basalt in Band 5. At the same time, it eliminates the scale inconsistency caused by the difference in background reflectance of different lithologies, making the characteristics of rock layers with large reflectance differences, such as dark basalt and light sandstone, comparable.

[0036] like Figure 4As shown, a comparative analysis of normalization methods was conducted to verify the effect of continuum normalization on improving the accuracy of rock strata identification. The horizontal axis represents five typical rock strata: granite, basalt, limestone, gneiss, and sandstone, while the vertical axis represents the classification accuracy (dimensionless, numerical range 0-1). Figure 4 The medium gray bars represent the conventional min-max normalization method, while the blue bars represent the continuum normalization method of this invention. Experimental data show that the conventional min-max normalization method, by ignoring local absorption characteristics, results in significantly lower accuracy in identifying basalt (dark rocks) and limestone (low-reflectance minerals), proving its inability to adapt to the differences in background reflectance among different lithologies. The continuum normalization method of this invention achieves higher accuracy than the conventional method across all five types of rock formations. Because it employs a dynamic benchmark and spectral envelope fitting, it solves the problem of dark rock characteristics approaching zero, while preserving the key absorption characteristics of silicate and carbonate minerals. Figure 4 The uniform and stable height of the blue columns in the middle verifies that the continuum normalization achieves cross-lithological comparability through wavelength difference standardization.

[0037] Example 4 like Figures 5 to 6 As shown, Figure 5 Three pure rock layer endmember spectra were demonstrated by clustering and key absorption band selection. These endmembers were calculated based on the initial reflectance of the predefined key absorption bands in the mineral spectral library (such as the silicate absorption characteristics of granite in Band 4-5), avoiding random initialization bias and ensuring that the endmembers conform to the true mineral absorption characteristics, thus providing physical mechanism constraints for the decomposition of mixed pixels. Figure 6 The bar chart visually displays the composition ratio of different rock layers in a single pixel. The abundance coefficient is obtained through sparsity constraint optimization, reflecting the spatial discontinuity of rock layer distribution. The dominant rock layer (such as granite) has a high abundance value, while the abundance of non-dominant rock layers is close to zero, which can effectively suppress the excessive mixing effect.

[0038] Example 5 like Figures 7 to 8 As shown, Figure 7 The study revealed the degree of attention the model pays to the derivative features of different bands. Higher weights are concentrated in key bands that are sensitive to lithology, such as Band4-5 (640-880nm). The adaptive attention mechanism can enhance discriminative features such as the 700nm absorption edge of hematite and suppress noise interference in non-characteristic bands. Figure 8 The model demonstrates its differentiated focus on the abundance characteristics of different rock strata, with limestone receiving the highest weight as the dominant stratum. This attention allocation allows the model to focus on spatially significant rock strata types and enhances the spectral contribution of the dominant strata during feature fusion.

[0039] like Figure 9As shown, the feature enhancement mechanism of this invention is used to analyze the separation ability of different rock layers through two-dimensional feature space visualization. Figure 9 In the scatter plot shown, the horizontal axis represents "shortwave infrared characteristics (band 5-7 derivatives)" (dimensionless), and the vertical axis represents "abundance characteristic weighted values" (dimensionless). Figure 9 The five types of rock strata show obvious clustering distribution. Granite (red dot group) is concentrated in the upper right area, basalt (blue dot group) is clustered in the lower left area, limestone (green dot group) is distributed in the middle, gneiss (purple dot group) is located in the upper middle, and sandstone (cyan dot group) occupies the lower right. The experimental results show that: 1) Similar rock strata clusters are closely clustered, proving that the feature representation has high consistency; 2) Dissimilar rock strata have clear boundaries, and the decision boundary can smoothly separate the boundaries of different categories; 3) The feature spatial structure conforms to geological laws, such as gneiss (metamorphic rock) being located in the transition zone between granite (igneous rock) and limestone (sedimentary rock). The experimental results also revealed noteworthy anomalies: approximately 10% of the granite samples intruded into the gneiss region, consistent with migmatization observed in field measurements; another 5% of the limestone samples overlapped with sandstone, reflecting the similarity in carbonate mineral composition between the two. Experimental results validate the effectiveness of attention-weighted feature fusion, which enhances the discriminative power of the original spectrum by strengthening derivative and abundance features.

[0040] Although the specific embodiments of the invention have been described above in conjunction with the accompanying drawings, this is not intended to limit the scope of protection of the invention. Based on the technical solutions of the invention, various modifications or variations that can be made by those skilled in the art without creative effort are still within the scope of protection of the invention.

Claims

1. An automatic rock strata identification method based on multispectral remote sensing, characterized in that, include: An adaptive normalization method based on continuum removal was used to process the collected multispectral remote sensing data samples to obtain normalized samples. The mixed pixel spectrum based on the physical mechanism of the rock mass is used to decompose the normalized sample, initialize each endmember after decomposition, and calculate the abundance coefficient of each endmember. The first derivative of each normalized sample in each band is calculated, and the multi-dimensional fused feature value is obtained by combining the abundance coefficient memory attention weighted feature fusion. A rock strata identification model based on deep learning is constructed, and the rock strata category probability is predicted based on the spectral adaptive activation function and multi-scale residual block architecture in the model and the fused feature value of each dimension. The rock strata identification model is optimized by incorporating the physical constraints of spectral decomposition into the loss function. An adaptive learning rate mechanism with parameter grouping is used to update the model parameters, resulting in a well-trained rock strata identification model. Newly acquired multispectral remote sensing data is processed and input into the trained model to obtain the final rock strata category probability prediction result. The lithology category corresponding to the maximum probability is taken as the identification result.

2. The method for automatic rock strata identification based on multispectral remote sensing according to claim 1, characterized in that, The adaptive normalization method based on continuum removal is as follows: For the original reflectance of each sample in each band after normalization, the effects of atmospheric scattering and illumination changes are first eliminated by fitting the spectral envelope. Specifically, a cubic polynomial fitting function is used. The center wavelength of the band and the original reflectance value are input to obtain the fitted value. Then, the maximum value of the original reflectance value and the fitted value is taken as the continuum function value of each sample in each band, which is used to characterize the spectral envelope. Next, the local absorption depth is normalized based on the difference between the continuum function value and the minimum reflectance, resulting in the normalized reflectance value of each sample in each band, i.e., the normalized sample.

3. The method for automatic rock strata identification based on multispectral remote sensing according to claim 2, characterized in that, The abundance coefficient of each endmember is calculated as follows: For each decomposed endmember, a clustering method is used to generate a sample cluster. Based on the normalized samples, an indicator function is used to determine whether the current band matches the known mineral spectral library. The initial reflectance value of each endmember in each band is calculated. Then, with the initial reflectance values ​​of fixed endmembers, the abundance coefficient of each sample in each endmember is calculated by minimizing the L2 norm squared of the linear combination of the normalized reflectance value and each endmember, combined with the coefficient regularization term. Here, endmembers represent different rock strata types.

4. The method for automatic rock strata identification based on multispectral remote sensing according to claim 3, characterized in that, The calculation process for multi-dimensional fusion feature values ​​is as follows: For each band of each sample, the difference in normalized reflectance values ​​between adjacent bands is calculated, and the first derivative value is obtained by combining the difference in center wavelength between adjacent bands with the derivative scaling factor. Then, for each sample, the attention weight of the first derivative value and the attention weight of the abundance coefficient in each dimension are calculated. Then, the weighted first derivative feature value and the weighted abundance coefficient are summed to obtain the fusion feature value in each dimension. The dimensions include the total number of bands and the total number of endmembers.

5. The method for automatic rock strata identification based on multispectral remote sensing according to claim 4, characterized in that, The specific operations in the deep learning-based rock strata identification model are as follows: The fused feature values ​​of each dimension are input into the model. These fused feature values ​​are then input into the linear layer. Combined with the projection weight matrix and bias vector of the input layer, the initial activation values ​​of each sample at each hidden layer node are obtained. These initial activation values ​​are then input into the multi-scale residual block architecture. Specifically, the activation function value of the previous layer residual block is processed by the spectral adaptive activation function and combined with the activation function values ​​of the previous two layers residual blocks multiplied by the cross-layer skip connection weight matrix to obtain the activation function value of the current layer residual block. Through layer-by-layer calculation, the activation function value of the final layer residual block is obtained. This activation function value is then combined with the output layer weight matrix and the output layer bias vector, and the predicted probability of the rock stratum category is calculated using the Softmax function.

6. The method for automatic rock strata identification based on multispectral remote sensing according to claim 5, characterized in that, The physical constraints of spectral decomposition are incorporated into the loss function, and the total loss is calculated by combining the cross-entropy classification loss term and the physical regularization term. Specifically, the cross-entropy loss term is calculated by measuring the deviation between the predicted probability of the rock stratum category obtained by the quantification model and the actual category. For each sample, the L2 norm square of the linear combination of the normalized reflectance value, the predicted abundance coefficient, and the endmember reflectance value is calculated, and the physical regularization term is obtained by combining the absolute value of the abundance coefficient.

7. The method for automatic rock strata identification based on multispectral remote sensing according to claim 6, characterized in that, A parameter grouping adaptive learning rate mechanism is adopted, which groups the trainable parameters according to their physical role in the model. A dynamically adjusted cosine annealing learning rate is applied to the weight parameters that control feature transformation to stabilize the convergence of deep networks, while a fixed small learning rate is applied to the bias parameters that control the offset to suppress overfitting.

8. The method for automatic identification of rock strata based on multispectral remote sensing according to claim 7, characterized in that, The training process of the model is as follows: The collected data is forward-propagated, processed, and then input into the rock strata identification model to obtain the predicted probability of the rock strata category. The loss function is then calculated, and the model parameters are updated through backpropagation. The model is iteratively optimized by the Adam optimizer until the iteration termination condition is met, at which point the optimization ends and a trained rock strata identification model is obtained.

9. The automatic rock strata identification method based on multispectral remote sensing according to claim 1, characterized in that: The collected multispectral remote sensing data samples are sourced from open-source datasets or self-collected datasets; All collected multispectral remote sensing data samples are labeled with lithological categories.

10. An automatic rock strata identification system based on multispectral remote sensing, characterized in that: The system includes modules for executing processing instructions for each step in the automatic rock strata identification method based on multispectral remote sensing as described in any one of claims 1-9, including the following modules: Data collection module: Collects multispectral remote sensing data samples; Data preprocessing module: performs normalization, mixed pixel decomposition, and special fusion operations on the collected samples; Rock strata identification module: Input the preprocessed sample into the rock strata identification module to predict the lithology of the sample; Optimize the training module: Optimize and update the parameters of the rock strata identification module through the model's total loss function and adaptive learning rate mechanism to obtain a well-trained rock strata identification module; Recognition result output module: After the sample to be identified is preprocessed, it is input into the trained rock stratum recognition module and the final recognition result is output.

Images