Land use change detection system based on remote sensing cloud computing

Abstract

The application discloses a land use change detection system based on remote sensing cloud computing. The system comprises a remote sensing data acquisition module, a spectral shift compensation module, a response field construction module, a confidence evaluation module, a multi-scale segmentation optimization module and a sub-pixel boundary optimization module. The application firstly realizes accurate spectral alignment by using invariant ground feature spectral reference, then establishes a spatial confidence evaluation model through change direction consistency and amplitude dispersion analysis to filter core change pixels; and based on boundary stability and internal homogeneity, the optimal segmentation scale is automatically determined, and the change area boundary is accurately corrected through sub-pixel spectral decomposition technology, which significantly improves the reliability and accuracy of land use change detection, greatly reduces the false change detection rate, effectively enhances the boundary positioning accuracy, adapts to complex change scenes and improves the processing efficiency.

Description

Technical Field

[0001] This invention relates to the fields of remote sensing technology and land resource monitoring technology, and more specifically, to a land use change detection system based on remote sensing cloud computing. Background Technology

[0002] Existing land use change detection technologies face numerous challenges when processing complex and variable remote sensing images. First, systematic spectral shifts are prevalent across multi-temporal images. Factors such as atmospheric variations, sensor attenuation, seasonal differences, and varying solar angles lead to significant inconsistencies in the spectral responses of the same land cover across different time phases. This is particularly evident in arid-humid seasonal transition zones, such as the vegetation cover change detection in tropical monsoon regions, which often results in some false change areas. Second, traditional methods are overly sensitive to image noise and local anomalies, generating numerous "salt-pepper noise"-like false change points in cloud transition zones, shadowed areas, and sensor-faulted stripes, severely interfering with decision analysis. Insufficient accuracy of change boundaries is particularly prominent in urban-rural transition zones and agroforestry areas. Jagged boundaries generated by pixel-level detection can deviate from actual land cover boundaries by 1-3 pixels, leading to significant errors in change area estimation. Furthermore, fixed-parameter segmentation methods cannot simultaneously adapt to large-scale changes in urban expansion areas and small-scale changes in fragmented farmland, resulting in poor adaptability to change scales in land monitoring practice. Current technologies lack effective spatial confidence assessment mechanisms, making it difficult to distinguish between high-confidence explicit changes and low-confidence potential changes, thus hindering the provision of reliable evidence in land law enforcement and monitoring applications. In the analysis of low-to-medium resolution images where mixed pixels are prevalent, sub-pixel-level change analysis capabilities are limited, particularly in woodland-grassland transition zones and land-water boundary areas, where mixed pixels can account for 25%-40% of the total changed area, leading to severe distortion in change area estimation. Existing methods over-rely on single-pixel spectral features, neglecting the spatial structure and contextual relationships of changes, making it difficult to identify continuous change patterns in complex topographic areas and highly heterogeneous landscapes.

[0003] In view of this, the present invention proposes a land use change detection system based on remote sensing cloud computing to solve the above problems. Summary of the Invention

[0004] To overcome the aforementioned deficiencies of the prior art and to achieve the above objectives, the present invention provides the following technical solution: a land use change detection system based on remote sensing cloud computing, comprising:

[0005] The remote sensing data acquisition module is used to acquire multi-temporal remote sensing image data of the same geographic area and extract spectral reference samples of invariant ground features in each temporal image.

[0006] The spectral drift compensation module is used to calculate the spectral drift compensation vector for each band based on the spectral response differences of the spectral reference sample at different time phases, and to perform pixel-by-pixel spectral alignment processing on multi-temporal remote sensing image data based on the spectral drift compensation vector.

[0007] The response field construction module is used to calculate the pixel-by-pixel spectral difference of the spectrally aligned image, construct the initial change response field, and extract the set of change candidate pixels based on the spectral difference distribution characteristics of the pixels in the initial change response field.

[0008] The confidence assessment module is used to obtain the neighborhood window of each candidate pixel in the set of candidate pixels with changes, calculate the consistency index of the change direction and the dispersion of the change amplitude of the pixels within the neighborhood window, and construct a spatial confidence evaluation model for candidate pixels based on the consistency index of the change direction and the dispersion of the change amplitude, and select candidate pixels with spatial confidence greater than the preset confidence threshold as core changed pixels.

[0009] The multi-scale segmentation optimization module is used to calculate the boundary stability coefficient and internal homogeneity coefficient of the changed region at different scales by using the core changed pixel as the seed point and combining the multi-scale segmentation results. Based on the weighted fusion result of the boundary stability coefficient and internal homogeneity coefficient, the optimal segmentation scale is determined, and a change detection patch is generated at the optimal segmentation scale.

[0010] The sub-pixel boundary optimization module is used to perform sub-pixel-level spectral decomposition on the boundary pixels of change detection patches, obtain the land cover component ratio of the boundary pixels, correct the boundary position of the change detection patches according to the temporal change rate of the land cover component ratio, and generate a land use change detection result map based on the corrected change detection patches.

[0011] The modules are connected via wired and / or wireless means to enable data transmission between them.

[0012] Compared with the prior art, the present invention has at least the following advantages:

[0013] (1) By eliminating systematic errors in traditional methods, this invention reduces the false change detection rate, enabling decision-making departments to obtain real and reliable change information and avoid resource waste and regulatory loopholes.

[0014] (2) This invention can accurately depict the dynamic evolution process of urban and rural boundaries, providing a scientific basis for urban planning; in the monitoring of ecologically fragile areas, it can accurately track the vegetation degradation and restoration process, supporting the evaluation of ecological restoration projects.

[0015] (3) The high adaptability of this invention enables it to simultaneously cope with various changing scenarios, from large-scale urbanization to small-scale farmland structure adjustment, meeting the monitoring needs at different scales. The ability to accurately locate boundaries significantly improves the accuracy of changing area statistics, providing precise data for land resource accounting and carbon sink estimation.

[0016] (4) The efficient processing capability of the present invention makes large-scale, high-frequency monitoring possible, and can support the dynamic monitoring of national land changes on a quarterly or even monthly basis, providing technical support for timely detection and handling of illegal land use. Attached Figure Description

[0017] Figure 1 This is a schematic diagram of the land use change detection system based on remote sensing cloud computing of the present invention. Detailed Implementation

[0018] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0019] This application provides a land use change detection system based on remote sensing cloud computing. The execution entities of the system include, but are not limited to, remote sensing data processing platforms, land resource monitoring centers, land spatial planning systems, and ecological environment monitoring platforms, which can be regarded as general computing nodes of this application. The change detection system includes, but is not limited to, at least one of the following: cloud-based remote sensing processing engine, distributed image analysis system, and intelligent change recognizer.

[0020] Please see Figure 1 In this embodiment of the invention, the land use change detection system based on remote sensing cloud computing includes:

[0021] The remote sensing data acquisition module is used to acquire multi-temporal remote sensing image data of the same geographic area and extract spectral reference samples of invariant features from each temporal image. This module accesses multi-source satellite image resources through a remote sensing cloud platform to acquire historical and current image data of the study area. Spectral reference samples of invariant features refer to ground targets that remain relatively stable during the study period, such as large buildings, stable water bodies, and mature forests. These samples serve as reference benchmarks for spectral changes, providing a basis for subsequent spectral alignment. The extraction process employs multi-temporal image overlay analysis technology, combined with spectral variation coefficient threshold screening, to automatically identify areas with high spectral stability in the time series as reference samples. These samples are spatially distributed throughout the study area, ensuring the regional representativeness and accuracy of spectral compensation.

[0022] The spectral drift compensation module calculates spectral drift compensation vectors for each band based on the differences in spectral response of spectral reference samples across different time phases, and performs spectral alignment processing. This module analyzes the changes in the spectral response of invariant ground features across images from different time phases, quantifying systematic spectral shifts caused by factors such as atmospheric conditions, sensor parameters, and seasonal variations. The compensation vector is calculated based on the statistical difference in spectral reflectance between the reference and target time phases, taking into account weighted corrections and normalization of spatial distribution density. Spectral alignment processing eliminates spectral differences caused by non-changing factors by applying the corresponding compensation vector to each pixel of the target time phase, laying the foundation for the detection of real-world ground feature changes.

[0023] The response field construction module is used to calculate pixel-by-pixel spectral dissimilarity in spectrally aligned images, construct an initial change response field, and extract a set of candidate change pixels. This module calculates pixel-level spectral differences between aligned images to form a continuous response field characterizing the degree of surface change. The spectral dissimilarity calculation employs a hybrid metric of improved spectral angle mapping and Euclidean distance to enhance sensitivity to different types of changes. The initial change response field identifies anomalous regions using adaptive threshold segmentation technology. Combining spectral change direction angle analysis and neighborhood consistency support assessment, pixels with high comprehensive change indication values ​​are selected as candidate change pixels, providing key areas for subsequent spatial analysis.

[0024] The confidence assessment module is used to acquire the neighborhood window of candidate pixels with changes, calculate the consistency index of change direction and the dispersion of change amplitude within the neighborhood, and construct a spatial confidence evaluation model. This module evaluates the spatial consistency and reliability of the change signal by analyzing the relationship between the change characteristics of candidate pixels and their surrounding pixels. The consistency index of change direction reflects the similarity of the change vector directions of neighboring pixels, while the dispersion of change amplitude quantifies the spatial uniformity of change intensity. The spatial confidence model combines these two indicators and considers the local extremum characteristics of candidate pixels in the change response field to comprehensively evaluate the spatial credibility of the change signal, select high-confidence core change pixels, and provide seed points for precise localization of the change region.

[0025] The multi-scale segmentation optimization module uses core change pixels as seed points and combines multi-scale segmentation results to calculate the boundary stability coefficient and internal homogeneity coefficient of change regions at different scales, thus determining the optimal segmentation scale. This module segments the image at multiple scale levels, evaluating the geometric characteristics and internal consistency of segmented objects containing core change pixels at different scales. The boundary stability coefficient is calculated using the boundary overlap rate and shape index between adjacent scales, while the internal homogeneity coefficient is based on the proportion of core change pixels within the segmented object and internal spectral fluctuation analysis. Through a comprehensive evaluation of boundary stability and internal homogeneity, the module identifies the segmentation scale that most accurately represents the change region, generating change detection patches with optimized boundaries.

[0026] The sub-pixel boundary optimization module performs sub-pixel-level spectral decomposition on the boundary pixels of change detection patches to obtain the proportion of land cover components in the boundary pixels and corrects the patch boundary positions based on the component change rate. Addressing the issue of blurred boundaries in change areas, this module employs a spectral mixing model to analyze the land cover composition within the boundary pixels. By extracting clean pixels inside and outside the change patches as endmembers, a spectral library is established. Fully constrained least-squares spectral unmixing is then performed on the boundary pixels to obtain the component proportions of changed and unchanged land cover. Based on temporal change analysis of component proportions and spectral similarity assessment, the attribution of boundary pixels is accurately defined, achieving sub-pixel-scale boundary precision and generating high-precision land use change detection result maps.

[0027] The modules are connected via wired and / or wireless means to enable data transmission between them.

[0028] In this embodiment of the invention, the detailed implementation steps for calculating the spectral drift compensation vector for each band include:

[0029] Statistical analysis was performed on the spectral reflectance values ​​of the spectral reference samples in images from various time phases to extract the median spectral reflectance for each band. Statistical analysis is a fundamental step in quantifying spectral changes, assessing systematic drift by analyzing the differences in the spectral response of the same ground feature across different time phases. The analysis process first extracts the spectral reflectance values ​​of each reference sample at its corresponding location in all time phase images, forming a multidimensional spectral dataset; then, the data is organized by band, and statistical analysis is performed separately for each band; finally, the median reflectance for each band is calculated as a representative statistic. The median, rather than the mean, was chosen to reduce the influence of outliers and improve statistical robustness. The median spectral reflectance intuitively reflects the central trend of each band, providing a reliable benchmark for drift calculation.

[0030] The initial drift is denoted as the shift in median spectral reflectance for each band between the reference and target time phases. Calculating the initial drift is a crucial step in determining the direction and magnitude of the spectral drift, directly reflecting the systematic spectral changes between different time phases. The calculation process selects one time phase as the reference (usually the phase with the best imaging quality or optimal atmospheric conditions), and the other time phases are used as the target time phases. The median difference for each band is calculated separately. The formula for calculating the drift is:

[0031] ;

[0032] in, For band The initial drift amount, For reference phase band The median spectral reflectance, For target time phase band The median spectral reflectance is used. The sign of the offset indicates the direction of drift, and the absolute value reflects the magnitude of drift. This band-independent calculation method takes into account the possible differential drift characteristics in different bands, providing a basis for accurate compensation.

[0033] Based on the spatial distribution density of the spectral reference samples, the initial drift is spatially weighted and corrected to obtain the spatially weighted drift. Spatial weighting correction is a key step in improving the regional adaptability of the drift, adjusting the global drift by considering the spatial distribution characteristics of the reference samples. The correction process first constructs a spatial density map of the reference samples, using kernel density estimation to quantify the concentration of sample distribution within the study area; then, a spatial weighting function is designed to give higher contribution weights to areas with high reference sample density; finally, the initial drift is locally adjusted to generate a spatially varying drift field. The formula for calculating the spatially weighted drift is:

[0034] ;

[0035] in, For position band Spatial weighted drift, This represents the sample density value at that location. For average density, For maximum density, This is the adjustment coefficient (usually taken as 0.2-0.5). This spatially adaptive drift calculation overcomes the limitations of traditional global compensation methods and can effectively cope with heterogeneous spectral variations over a large study area.

[0036] The ratio of the spatially weighted drift to the standard deviation of the spectral reflectance for each band is used as the normalized drift coefficient. Normalization ensures the relative rationality of drift compensation across different bands, taking into account the influence of band-specific variability on drift assessment. The process first calculates the standard deviation of the spectral reflectance of the reference sample for each band, quantifying the band's natural variability level; then, it calculates the ratio of the drift to the standard deviation, forming a dimensionless normalized coefficient. Normalization ensures that drift compensation for bands with high natural variability is relatively conservative, while compensation for bands with low natural variability is more sensitive, conforming to the objective laws of spectral variation and improving the rationality and accuracy of compensation.

[0037] Spectral drift compensation vectors for each band are constructed based on normalized drift coefficients. The direction of these vectors is determined by the sign of the initial drift, while the amplitude is determined by the absolute value of the normalized drift coefficients. Vector construction is the final step in spectral drift compensation, transforming the previously calculated normalized drift coefficients into directly applicable correction values. The construction process first determines the compensation direction, aligning it with the initial drift direction; then, it sets the compensation amplitude based on the absolute value of the normalized drift coefficients, while considering a global adjustment factor to control the overall compensation intensity; finally, it generates a complete compensation vector field for pixel-by-pixel spectral adjustment. The compensation vectors are applied to each pixel of the target phase, achieving precise spectral alignment, eliminating spectral differences caused by non-changing factors, and providing a reliable basis for subsequent change detection.

[0038] In this embodiment of the invention, the detailed implementation steps for extracting the set of candidate pixels with changes include:

[0039] An adaptive threshold segmentation method is used to segment the spectral dissimilarity of the initial change response field, marking pixels with spectral dissimilarity greater than the adaptive threshold as preliminary candidate pixels. Adaptive threshold segmentation is a key step in identifying significantly changed regions from continuous change responses, automatically adapting to the characteristics of different images through a data-driven threshold determination method. The segmentation process first analyzes the statistical distribution characteristics of the change response field, obtaining the distribution pattern through histogram analysis or kernel density estimation; then, the OTSU method or an improved valley search algorithm is used to determine the optimal segmentation threshold, which maximizes the variance between changed and non-changed categories; finally, pixels with spectral dissimilarity greater than the threshold are marked as preliminary candidate pixels as potential change regions. Adaptive thresholding is more flexible than fixed thresholding, allowing adjustment of the segmentation criteria according to specific data characteristics, improving the adaptability and robustness of change detection.

[0040] The angle between the difference vectors of preliminary candidate pixels across multispectral bands is calculated to obtain the spectral change direction angle of the preliminary candidate pixels. The spectral change direction angle is an important characteristic characterizing the nature of land use change; different types of land use change typically exhibit different spectral change directions. The calculation process first constructs the spectral vector of each candidate pixel in the preceding and following time phases, represented as two points in n-dimensional space (n being the number of bands); then, the angle between these two spectral vectors is calculated to quantify the directional characteristics of the spectral change. The formula for calculating the spectral change direction angle is:

[0041] ;

[0042] in, The spectral change direction angle, and These are the spectral vectors for the preceding and following time phases, respectively. Let S denote the vector inner product, and ||S|| denote the vector norm. The directional angle of spectral variation ranges from [0, π]. Different angles correspond to different types of land surface changes, such as vegetation increase, water body expansion, or urbanization, providing directional indicators for the analysis of the nature of the changes.

[0043] The number of pixels with similar spectral change direction angles within the neighborhood of a preliminary candidate pixel is counted and denoted as the directional consistency support. Directional consistency support is a spatial indicator for evaluating the reliability of change detection; it verifies the authenticity of the change signal by analyzing the consistency of change direction within a local region. The statistical process first determines an appropriate neighborhood window size (typically 3×3 or 5×5 pixels); then calculates the difference in spectral change direction angle between the center pixel and each neighboring pixel; finally, it counts the number of pixels whose direction angle difference is less than a threshold (typically π / 8 or π / 6), which is taken as the directional consistency support. Genuine change areas typically exhibit high local directional consistency, while noise or spurious changes show random directional distribution and low support. This indicator effectively improves the noise resistance of change detection and reduces the impact of isolated spurious change points.

[0044] The product of directional consistency support and spectral dissimilarity of the initial candidate pixels is used as the comprehensive change indicator. This comprehensive change indicator is a composite metric that fuses intensity and directional consistency, comprehensively evaluating the significance and reliability of changes. The calculation process multiplies spectral dissimilarity (reflecting change intensity) with directional consistency support (reflecting change reliability) to obtain the comprehensive metric, considering both the magnitude and spatial consistency of the change. The advantage of this indicator is its ability to balance the intensity and spatial structure features of changes, identifying both high-intensity but isolated change points and moderate-intensity but spatially consistent change regions, thus improving the overall performance of change detection.

[0045] Initial candidate pixels are sorted in descending order of their comprehensive change indicator values. A predetermined percentage of initial candidate pixels with the highest comprehensive change indicator values ​​are selected to form a set of change candidate pixels. Sorting and filtering is a crucial step in determining the final candidate pixel set, ensuring the selection of the most significant and reliable change regions through priority sorting. The filtering process first sorts all initial candidate pixels from highest to lowest comprehensive change indicator value; then, pixels ranking at the top are selected according to a predetermined percentage (usually the top 10%-20%) or a fixed threshold; finally, a set of change candidate pixels is formed, serving as input for subsequent spatial confidence assessment. This sorting-based filtering method is more flexible than fixed-threshold filtering, adaptively controlling the number of candidate pixels, balancing the completeness and accuracy of detection, and providing a good candidate foundation for high-quality change detection.

[0046] In this embodiment of the invention, the detailed implementation steps for calculating the consistency index of the change direction and the dispersion of the change amplitude of pixels within the neighborhood window include:

[0047] For each pixel within a neighborhood window, the difference in its temporal spectral vector between preceding and following phases in multispectral space is calculated and denoted as the pixel change vector. The pixel change vector is a fundamental characteristic representation of surface changes, directly reflecting the direction and magnitude of change in multispectral space. The calculation process first extracts the spectral vector (n-dimensional vector, where n is the number of bands) of each pixel between preceding and following phases; then, the vector difference is calculated to obtain the vector representing the change. The formula for calculating the pixel change vector is:

[0048] ;

[0049] in, For position The pixel change vector at that location. and These are the spectral vectors for the preceding and subsequent time phases, respectively. The direction of the change vector indicates the nature of the spectral change (such as increased vegetation or decreased water volume), while the magnitude of the vector reflects the intensity of the change. This vector representation preserves all the multidimensional information of the change, providing fundamental data for subsequent directional and magnitude analysis.

[0050] Using the pixel change vector of a candidate pixel as a reference vector, the cosine similarity between the pixel change vectors of other pixels within the neighborhood window and the reference vector is calculated. Cosine similarity calculation is the core step in evaluating the consistency of change directions, quantifying the degree of coordination of change directions within the neighborhood through vector angular similarity. The calculation process first sets the change vector of the central candidate pixel as the reference; then, it calculates the cosine similarity between the change vector of each pixel within the neighborhood and the reference vector. The cosine similarity ranges from [-1, 1], with values ​​closer to 1 indicating more consistent directions, values ​​closer to -1 indicating more opposite directions, and values ​​close to 0 indicating unrelated directions. This index accurately characterizes the spatial consistency pattern of change directions within the neighborhood.

[0051] The cosine similarity of all pixels within a neighborhood window is averaged and denoted as the direction consistency index. The direction consistency index is a statistical indicator that comprehensively reflects the coordination of change directions in a local area, providing a regional-level assessment by integrating the direction similarity of all pixels within the neighborhood. The statistical process calculates the arithmetic mean of all cosine similarity values ​​within the neighborhood window to obtain a single consistency measure. Typically, the direction consistency index is high (close to 1) within a region of genuine change, indicating that pixel changes within the region have similar directions; while the consistency index is low in noisy or spurious change regions, exhibiting randomness and inconsistency in direction. This index effectively distinguishes between overall change and scattered noise, providing a reliable basis for spatial confidence assessment based on direction consistency.

[0052] Calculate the magnitude of the pixel change vectors for all pixels within the neighborhood window, obtaining a magnitude sequence. Magnitude calculation is a fundamental step in assessing the distribution of change magnitude, quantifying the intensity level of surface change through the magnitude of the change vectors. The calculation process calculates the Euclidean norm for the change vector of each pixel within the neighborhood window, forming a scalar sequence representing the intensity of change. The magnitude calculation formula is:

[0053] ;

[0054] in, For position The magnitude of the vector at the point of change. The first change vector One portion, The vector dimension (number of bands) is used. The modulus length sequence reflects the spatial distribution characteristics of intensity variation within a local region, providing basic data for amplitude consistency analysis.

[0055] The standard deviation of the modulus sequence is calculated, and the ratio of the standard deviation to the mean of the modulus sequence is used as the amplitude dispersion. Amplitude dispersion is a statistical indicator for evaluating the spatial uniformity of variation intensity, quantifying the relative dispersion of variation intensity within a neighborhood using the coefficient of variation. The calculation process first calculates the standard deviation and mean of the modulus sequence; then, the ratio of the two is calculated to obtain the coefficient of variation as a measure of dispersion. The formula for calculating amplitude dispersion is:

[0056] ;

[0057] in, For the dispersion of the change range, The standard deviation of the modulus-length sequence within the neighborhood window. This represents the mean of the modulus length sequence. A smaller dispersion value indicates a more uniform distribution of change intensity, typically corresponding to actual areas of land cover change; a larger dispersion value indicates a more uneven distribution, potentially corresponding to noise or boundary transition zones. This index complements the directional consistency index, together providing a comprehensive assessment of the spatial characteristics of change.

[0058] In this embodiment of the invention, the detailed implementation steps for constructing the spatial confidence evaluation model for candidate pixels include:

[0059] The directional confidence component is obtained by multiplying the directional consistency index by a preset directional consistency weight. This directional confidence component is the first component of spatial confidence assessment, reflecting the contribution of spatial consistency in directional change to the overall confidence score. The calculation process multiplies the directional consistency index by a predetermined weight coefficient, adjusting the relative importance of directional factors in the overall assessment. These weight coefficients are typically preset based on the application scenario and change type, with higher weights assigned to change types sensitive to directional consistency (such as seasonal vegetation change). The directional confidence component directly affects the calculation of spatial confidence; areas with high directional consistency receive higher confidence scores, improving the assessment model's ability to identify real changes.

[0060] The inverse of the dispersion of the amplitude is normalized to obtain the normalized inverse dispersion. This normalized inverse dispersion is then multiplied by a preset amplitude stability weight to obtain the amplitude confidence component. The amplitude confidence component is the second component of spatial confidence assessment, reflecting the contribution of spatial uniformity of change intensity to the overall confidence. The calculation process first takes the inverse of the dispersion, ensuring that uniform regions receive high values ​​and non-uniform regions receive low values; then, it normalizes the values ​​to the [0,1] interval using maximum and minimum values; finally, it multiplies by a preset weight coefficient to adjust the relative importance of the amplitude factor. The amplitude stability weight is typically set considering the characteristics of the change type, assigning higher weights to change types sensitive to uniformity (such as large-area land cover transformation). The amplitude confidence component and the directional confidence component complement each other, jointly improving the comprehensiveness and accuracy of spatial confidence assessment.

[0061] The local maximum indicator factor (ROMA) for candidate pixels in the initial change response field is calculated. The MMAA is determined by whether the spectral dissimilarity of the candidate pixel is the maximum value within its neighborhood. The MMAA is a key indicator for evaluating the significance of a pixel in the change response, identifying local response peaks, which typically correspond to the core region of change. The calculation compares the spectral dissimilarity of a candidate pixel with the dissimilarity of all other pixels within its neighborhood window. If the candidate pixel's value is the maximum, the indicator factor is 1; otherwise, it is 0 or a ranking-based attenuation value. This evaluation method based on local extrema effectively highlights the central change region, reduces the weight of edge and transition regions, and improves the accuracy of identifying the core change region. Simultaneously, MMA analysis also helps suppress the influence of spatially correlated noise, enhancing the spatial consistency of change detection.

[0062] A spatial confidence evaluation model for candidate pixels is constructed by weighting and summing the direction confidence component, amplitude confidence component, and local maximum indicator factor. The spatial confidence evaluation model is a comprehensive evaluation framework integrating multi-dimensional spatial features, consolidating the contributions of various factors through reasonable weight allocation. The construction process first determines the weight coefficients of the three components; direction confidence, amplitude confidence, and local maximum indicator are typically allocated in a 4:3:3 ratio, but can be dynamically adjusted according to specific application scenarios. Then, the weighted sum is calculated to obtain the final spatial confidence score. The spatial confidence calculation formula is:

[0063] ;

[0064] in, For spatial confidence, For the direction confidence component, For the amplitude confidence component, As a local maximum indicator, , and The corresponding weight coefficients and satisfying The spatial confidence score ranges from [0,1], with higher values ​​indicating higher reliability of the change. The spatial confidence evaluation model comprehensively considers the consistency of the direction of change, the uniformity of the amplitude, and the local significance, providing a multi-dimensional evaluation standard for the selection of core change pixels and significantly improving the reliability and accuracy of change detection.

[0065] Candidate pixels with spatial confidence scores greater than a preset confidence threshold are selected as core change pixels. Core change pixel selection is a crucial step in identifying high-confidence change regions, using a confidence threshold to identify the most reliable change signals. The selection process first sets an appropriate confidence threshold (typically 0.7-0.8); then, all candidate pixels with spatial confidence scores exceeding the threshold are marked as core change pixels, serving as seed points for the change regions. These core change pixels exhibit high directional consistency, uniform change amplitude, and local saliency, representing regions most likely to experience real changes. Core change pixels provide reliable starting points for subsequent multi-scale segmentation and region growing, effectively guiding the accurate extraction and boundary determination of change regions.

[0066] In this embodiment of the invention, the detailed implementation steps for calculating the boundary stability coefficient and internal homogeneity coefficient of the changing region at different scales include:

[0067] Within a preset scale range, the image is segmented at multiple scales with a preset scale step size to obtain a set of segmented objects at each scale. Multi-scale segmentation is a key step in generating objects of different granularities, creating spatial units at multiple scales by adjusting segmentation parameters. The segmentation process first determines the scale range (typically ranging from fine-grained to coarse-grained to cover the expected changes in object size) and the scale step size (controlling the degree of difference between adjacent scales); then, region growing or multi-resolution segmentation algorithms are used to segment the image at each scale parameter; finally, a series of segmented object sets with increasing scales are obtained, forming a multi-scale segmentation hierarchy. Segmented objects at different scales capture the spatial structure of the image at different scales; fine-scale segments retain details but objects are fragmented, while coarse-scale segments have smooth boundaries but may mix different ground features. Multi-scale segmentation provides a candidate object set for subsequent optimal scale selection and is the foundation for target-level change detection.

[0068] For each scale, the overlap rate of boundary pixels of segmented objects containing core changed pixels is calculated between adjacent scales and denoted as the boundary overlap rate. The boundary overlap rate is a core indicator for assessing the stability of segmentation boundaries, quantifying the reliability of the boundary by the consistency of boundary positions between adjacent scales. The statistical process first identifies segmented objects containing core changed pixels at each scale; then, it extracts the boundary pixel set of these objects; next, it calculates the ratio of the intersection of the boundary sets at the current scale and those at adjacent scales (usually the previous and next scales) to the current boundary set, obtaining the boundary overlap rate. The boundary overlap rate is the ratio of the number of overlapping boundary pixels to the number of pixels in the current boundary set. A high overlap rate indicates that the boundary remains stable with scale changes, usually corresponding to a true ground feature boundary; a low overlap rate indicates that the boundary is prone to change with scale, possibly indicating an artificially segmented boundary. This indicator provides a direct basis for boundary stability assessment.

[0069] The boundary stability coefficient is obtained by multiplying the boundary coincidence rate by the reciprocal of the perimeter-to-area ratio of the segmented object. The boundary stability coefficient is a composite index that comprehensively considers boundary stability and shape regularity, evaluating boundary quality through a combination of boundary coincidence rate and shape index. The calculation process first calculates the perimeter-to-area ratio of the segmented object, an index reflecting shape complexity, with regular shapes having smaller values ​​and complex shapes having larger values; then, its reciprocal is taken, resulting in a higher value for regular shapes; finally, it is multiplied by the boundary coincidence rate to obtain the final stability coefficient. This evaluation method, combining boundary stability and shape regularity, considers both the temporal continuity of the boundary and the rationality of the spatial morphology, providing a reliable indicator for selecting segmentation scales with stable and natural boundaries.

[0070] The standard deviation of the spectral differences among all pixels within a segmented object is calculated and denoted as the internal spectral fluctuation value. The internal spectral fluctuation value is a fundamental indicator for assessing the homogeneity within a segmented object, quantifying internal consistency through the degree of variation in spectral differences. The calculation process statistically analyzes the standard deviation of the spectral differences among all pixels within the segmented object, reflecting the dispersion of the differences. A smaller standard deviation indicates a more consistent degree of variation within the object, typically corresponding to homogeneous variation regions; a larger standard deviation indicates inconsistent internal variation, potentially mixing varied and non-variable regions. This indicator directly reflects the internal homogeneity of the segmented object as a unit of variation, providing a spectral dimension basis for evaluating segmentation quality.

[0071] The internal homogeneity coefficient is obtained by multiplying the proportion of core variable pixels within a segmented object by the inverse of its internal spectral fluctuation value. The internal homogeneity coefficient is a comprehensive index for evaluating the internal consistency of a segmented object, quantifying its internal characteristics through a combination of core pixel density and spectral uniformity. The calculation process first calculates the ratio of the number of core variable pixels to the total number of pixels in the segmented object, obtaining the core pixel proportion; then, the inverse of the internal spectral fluctuation value is taken, resulting in higher values ​​for objects with lower fluctuations; finally, the two are multiplied to obtain the homogeneity coefficient. This evaluation method, combining spatial distribution and spectral characteristics, considers both the representativeness of the change signal and the consistency of change characteristics, providing a comprehensive index for selecting a segmentation scale with internal homogeneity. An ideal segmented object should have a high proportion of core pixels and low internal fluctuations, corresponding to natural units of actual ground feature changes.

[0072] In this embodiment of the invention, the detailed implementation steps for determining the optimal segmentation scale include:

[0073] Scale sequence normalization is performed on the boundary stability coefficients to obtain normalized boundary stability coefficients. Scale sequence normalization is a key step in ensuring the comparability of indicators across different scales, eliminating the influence of inherent scale differences by unifying the dimensional range. The normalization process uses the maximum-minimum normalization method to map the boundary stability coefficients in the entire scale sequence to the [0,1] interval. The normalization formula is:

[0074] ;

[0075] in, For scale The normalized boundary stability coefficient, These are the original coefficients. and These are the minimum and maximum values ​​in the scale sequence, respectively. Normalization ensures that evaluation indicators at different scales have the same numerical range, facilitating subsequent indicator fusion and comparison, and providing a standardized measurement basis for optimal scale selection.

[0076] The internal homogeneity coefficients are scale-series normalized to obtain normalized internal homogeneity coefficients. Similar to the boundary stability coefficients, the internal homogeneity coefficients also need scale-series normalization to ensure fair comparisons across different scales. Normalization uses the same minimax method to map the original homogeneity coefficients to a uniform interval. This standardization process eliminates systemic differences between different scales, making the assessment of internal homogeneity scale-independent and objectively reflecting the internal characteristic quality of the segmented object, unaffected by differences in the scale parameters themselves.

[0077] The weighted average of the normalized boundary stability coefficient and the normalized internal homogeneity coefficient yields the comprehensive evaluation index for each scale. This comprehensive evaluation index is an integrated indicator that combines boundary and internal dimensions, balancing the importance of boundary quality and internal consistency through reasonable weight allocation. The calculation process first determines the weight coefficients, typically allocating boundary stability and internal homogeneity in a 1:1 ratio, but this can be adjusted according to application requirements; for example, the weight of boundary stability can be increased for scenarios with high boundary accuracy requirements. Then, the weighted average is calculated to obtain the comprehensive evaluation score for each scale. The comprehensive evaluation index comprehensively reflects multiple aspects of segmentation quality; high scores correspond to ideal segmentation objects with stable boundaries and internal homogeneity, providing a comprehensive measure for optimal scale selection.

[0078] The gradient rate of change of the comprehensive evaluation index across the scale sequence is calculated to identify the inflection point scale where the gradient rate of change turns from positive to negative. Gradient rate of change analysis is a key method for discovering inflection points in the evaluation curve, identifying the turning point in the growth pattern of the evaluation index through changes in the derivative. The analysis process first calculates the difference in the comprehensive evaluation index between adjacent scales to obtain the first-order gradient; then, it calculates the rate of change of the gradient, i.e., the second-order gradient; finally, it identifies the inflection point where the gradient turns from positive to negative, i.e., the position where the growth rate of the evaluation index begins to decline. This type of inflection point scale usually corresponds to the optimal segmentation scale, at which the evaluation index has reached a high level and further increases in scale parameters yield diminishing returns. The inflection point analysis method can adaptively identify key turning points in the evaluation curve, avoiding the problem of excessive bias towards large scales that may occur with simply selecting the maximum value.

[0079] If an inflection point scale exists, it is used as the optimal segmentation scale; otherwise, the scale corresponding to the maximum value of the comprehensive evaluation index is used as the optimal segmentation scale. Determining the optimal scale is the final decision-making step in multi-scale analysis, flexibly selecting the most suitable segmentation level through a strategy of prioritizing inflection points and selecting maximum values ​​as alternatives. The decision-making process first checks whether a clear gradient inflection point exists; if so, the scale corresponding to the inflection point is selected as the optimal scale, which is usually the balance point between boundary quality and internal homogeneity; if no clear inflection point exists (e.g., the evaluation curve shows monotonous growth or large fluctuations), it degenerates into a maximum value strategy, selecting the scale with the highest comprehensive evaluation index. This dual-standard decision-making mechanism improves the flexibility and robustness of scale selection, adapting to the characteristics of different types of change regions, ensuring the generation of change detection patches at the optimal scale, and providing a high-quality foundation for subsequent boundary optimization.

[0080] In this embodiment of the invention, the detailed implementation steps for obtaining the proportion of land cover components in boundary pixels include:

[0081] Pure pixels within and outside the detected patch are extracted and designated as changed and unchanged endmembers, respectively. Endmember extraction is a crucial prerequisite for spectral mixture analysis, establishing a spectral reference basis by identifying representative pure pixels. The extraction process first defines the criteria for pure pixels, typically based on the statistical distribution of spectral dissimilarity. The set of pixels with the highest dissimilarity within the patch (e.g., the quartiles above) is selected as changed endmembers, and the set of pixels with the lowest dissimilarity in the external neighborhood (e.g., the quartiles below) is selected as unchanged endmembers. Then, the selected endmembers are evaluated for spectral purity, eliminating potentially mixed pixels to ensure representativeness. This statistically based endmember extraction method automatically adapts to different types of change characteristics, providing a reliable endmember basis for subsequent spectral decomposition. Endmember quality directly affects the accuracy of decomposition; therefore, the identification of pure pixels is a critical step in sub-pixel analysis.

[0082] An endmember spectral library is constructed, comprising both variable and unvariable endmembers. This construction is a systematic step in organizing the spectral reference set, providing a standardized set of endmembers for spectral demixing. The process begins by collecting the spectral vectors of all identified variable and unvariable endmembers; then, statistical analysis is performed, clustering each endmember class where possible to reduce redundancy and extract representative spectra; finally, a structured endmember spectral library is formed, with each endmember containing a complete spectral reflectance curve and class identifier. The endmember spectral library typically contains 2-5 variable endmembers and 2-5 unvariable endmembers, dynamically adjusted according to the complexity of the variation. A high-quality endmember spectral library lays the foundation for accurate spectral decomposition, particularly important in the analysis of mixed pixels at the boundaries of variable regions.

[0083] For the spectral vectors of boundary pixels, endmember decomposition is performed using the fully constrained least squares method to obtain the abundance values ​​of each endmember in the boundary pixels. Endmember decomposition is the core algorithm of sub-pixel analysis, which decomposes mixed pixels into linear combinations of pure endmembers through a mathematical model. The decomposition process first assumes that the spectrum of the boundary pixels is a linear combination of the endmember spectra; then, an optimization problem is constructed, aiming to minimize the squared error between the actual spectrum and the endmember combined spectrum; finally, the fully constrained least squares problem is solved (the constraint conditions are that the abundance is non-negative and the sum is 1) to obtain the contribution ratio of each endmember. The mathematical model of the fully constrained least squares method is a linear mixture model, which ensures the physical interpretability of the decomposition results, and the abundance values ​​directly correspond to the area proportion of land cover. Through this method, the proportion of each land cover component in the boundary pixels is accurately quantified, providing sub-pixel level information for fine-tuning the boundary.

[0084] The sum of abundance values ​​of endmembers in the changed category is used as the proportion of changed land cover components in the boundary pixel, while the sum of abundance values ​​of endmembers in the unchanged category is used as the proportion of unchanged land cover components in the boundary pixel. Component proportion calculation is a semantic interpretation step in the endmember decomposition results, obtaining the change information needed for practical applications through endmember category abundance aggregation. The calculation process adds the abundance values ​​of all endmembers in the same category to obtain the overall proportion of that category. This abundance aggregation method simplifies the complex endmember structure, transforming the decomposition results into an intuitive binary classification (changed and unchanged), facilitating subsequent boundary decisions. The proportion of changed land cover components directly reflects the actual area of ​​change within a pixel, providing a quantitative basis for sub-pixel-level boundary precision and serving as a key indicator for determining the attribution of mixed pixels.

[0085] In this embodiment of the invention, the detailed implementation steps for correcting the boundary position of the change detection patch include:

[0086] The difference in the proportion of changed land cover components in a boundary pixel between different time phases is calculated and denoted as the component change. The component change is a direct indicator for assessing the degree of change in boundary pixels, quantifying the magnitude of sub-pixel-scale changes through temporal comparison. The calculation process analyzes the changing proportion of land cover components at the same boundary pixel location between different time phases and calculates the difference. The larger the absolute value of the component change, the more significant the change in the land cover composition within the pixel; a change close to 0 indicates that the land cover composition within the pixel remains essentially unchanged. This indicator directly reflects the degree of change experienced by the boundary pixel, providing a primary basis for boundary determination.

[0087] A threshold for component variation is set. Boundary pixels with component variation greater than the threshold are marked as defined boundary pixels, while those with component variation less than the threshold are marked as undefined boundary pixels. Boundary pixel classification is a crucial first step in boundary optimization, distinguishing boundary pixels into defined and undefined categories based on the threshold. The classification process first sets an appropriate component variation threshold (typically 0.3-0.5); then, pixels are marked based on the comparison between the component variation and the threshold. Defined boundary pixels represent clearly defined boundary areas and are directly included in the final variation region; undefined boundary pixels require further analysis to determine their classification, representing the fuzzy boundary of the variation region. This preliminary classification based on the degree of variation provides a framework for boundary refinement, particularly beneficial for determining the boundaries of gradual transition zones.

[0088] For uncertain boundary pixels, the spectral similarity between them and the core region inside the change detection patch, as well as their spectral similarity with the outer neighborhood of the patch, is calculated. Spectral similarity calculation is a key basis for boundary attribution determination, identifying the attribute of an uncertain pixel by comparing its spectral relationship with the inner and outer regions. The calculation process first defines the core region inside the patch (usually the inner region excluding the boundary) and the outer neighborhood (usually a certain range outside the boundary); then, it extracts the average or representative spectra of these two regions; finally, it calculates the spectral similarity between the uncertain boundary pixel and both the core region and the outer neighborhood, typically using spectral angle measurement or Euclidean distance as the similarity metric. This two-way similarity comparison fully utilizes spectral information, providing a spectral dimensional criterion for the final attribution of boundary pixels, and is particularly suitable for change types with obvious spectral features.

[0089] Based on the comparison results of spectral similarity, uncertain boundary pixels with higher spectral similarity to the internal core region are retained within the change detection patch, while uncertain boundary pixels with higher spectral similarity to the external neighborhood are removed from the change detection patch. Boundary assignment decision is the final step in boundary optimization, determining the final assignment of boundary pixels through similarity comparison. The decision process compares the spectral similarity between uncertain boundary pixels and the internal and external regions, using the "proximity principle" to determine their assignment. This spectral similarity-based decision method effectively utilizes spectral information to address uncertainties in component analysis, ensuring that boundary determination considers both the degree of sub-pixel change and the overall similarity of spectral features, thus improving the accuracy and reliability of boundary determination.

[0090] Based on the determined change boundary pixels and the retained uncertain boundary pixels, the boundary locations of change detection patches are reconstructed. Boundary reconstruction is the final output step of boundary optimization, forming a corrected change region by integrating all retained pixels. The reconstruction process combines all determined change boundary pixels and the uncertain boundary pixels to be retained, generating new change detection patch boundaries. Compared to the original segmentation boundaries, the reconstructed boundaries more accurately reflect the actual range of the change region, especially in the gradual transition region, where sub-pixel-level analysis significantly improves the accuracy of boundary location. This boundary optimization method based on component analysis and spectral similarity achieves sub-pixel precision in determining change boundaries, overcoming the jagged effect and mixed pixel misjudgment problems of traditional pixel-level boundaries, providing a reliable guarantee for high-precision land use change detection.

[0091] In this embodiment of the invention, the detailed implementation steps for calculating the positional overlap rate of boundary pixels of a segmented object containing core change pixels between adjacent scales include:

[0092] The system extracts the set of boundary pixels for the segmented object at the current scale, denoted as the current boundary set. Boundary pixel extraction is a fundamental step in boundary stability analysis, determining the object boundary by identifying the surrounding pixels of the segmented object. The extraction process first identifies the segmented object containing core changing pixels at the current scale; then, a morphological boundary extraction algorithm is used to determine the edge positions of the object. Specifically, eight-neighbor connectivity analysis is used to mark pixels within the object that are directly adjacent to external pixels of the segmented object as boundary pixels. To improve the accuracy of boundary determination, the system employs an improved Laplacian edge detection method, considering gray-level gradient changes inside and outside the object to ensure the accuracy of the boundary position. The boundary pixel set is typically stored using a spatial index structure, such as an R-tree or quadtree, to support efficient spatial queries and set operations. The current boundary set completely records the boundary information of changing objects at a specific scale, providing a reference benchmark for comparisons at adjacent scales.

[0093] The set of boundary pixels for segmented objects containing core changed pixels at adjacent scales is extracted and denoted as the adjacent boundary set. The process of extracting boundaries at adjacent scales is similar to that at the current scale, but focuses on other scale levels adjacent to the current scale. The extraction process first identifies adjacent scales, typically including the next scale with a smaller step size than the current scale and the previous scale with a larger step size. Then, within each adjacent scale, segmented objects containing the same core changed pixels are identified. Finally, the same boundary extraction method as at the current scale is applied to determine the set of boundary pixels for these objects. In multi-scale analysis, segmented objects at different scales may have inclusion or intersection relationships; therefore, boundary extraction needs to ensure the correct correspondence of objects. The adjacent boundary set reflects the boundary characteristics of the changed region at different scale parameters, providing a comparative object for boundary stability assessment.

[0094] The system calculates the intersection of the current boundary set and the adjacent boundary sets to obtain the set of overlapping boundary cells. This intersection calculation is a crucial step in quantifying boundary position consistency, assessing boundary stability by identifying boundary cells that maintain the same position at different scales. The calculation process employs efficient set operation algorithms, such as hash table lookup or spatial index matching, to identify cells that appear simultaneously in both the current and adjacent boundary sets. To handle minor cell position shifts, the system implements a tolerance matching mechanism, allowing positional correspondence determination within a preset small range (e.g., a distance of 1 pixel), enhancing the robustness of boundary matching. The set of overlapping boundary cells represents the boundary portion that remains stable during scale changes, typically corresponding to stable features of actual land cover boundaries, such as clear land cover transition areas or distinct terrain boundaries.

[0095] The boundary overlap rate is the ratio of the number of pixels in the overlapping boundary pixel set to the number of pixels in the current boundary set. Calculating the overlap rate is the final step in boundary stability quantification, visually representing the degree of boundary consistency through a proportional relationship. The calculation process first counts the number of pixels in the overlapping boundary pixel set and the number of pixels in the current boundary set; then, it calculates the ratio between the two, obtaining an overlap rate within the range [0,1]. A overlap rate closer to 1 indicates a more stable boundary across different scales, typically corresponding to a real ground feature boundary; a overlap rate closer to 0 indicates a significant change in the boundary with scale, possibly corresponding to an artificial boundary generated by a segmentation algorithm. As a crucial component of the boundary stability coefficient, the boundary overlap rate directly influences the selection of the optimal scale, providing an objective basis for generating change detection patches with natural boundaries.

[0096] In practical applications, the system typically calculates the average overlap rate with adjacent scales for each segmented object, comprehensively evaluating the overall stability of the boundary across the scale sequence. When the overlap rate at a certain scale is observed to be significantly higher than that at other scales, that scale is usually considered a candidate scale for capturing real-world feature boundaries. Boundary overlap rate analysis effectively distinguishes between stable natural boundaries and unstable artificially segmented boundaries, providing a reliable quantitative basis for subsequent boundary stability coefficient calculations and optimal scale determination.

[0097] This invention achieves accurate detection of land use change in multi-temporal remote sensing images through spectral benchmark sample extraction, spectral drift compensation, change response field construction, spatial confidence assessment, multi-scale segmentation optimization, and sub-pixel boundary correction. The multi-level optimization detection method of this invention can accurately identify change areas and provide sub-pixel precision boundaries, effectively addressing the challenges of detecting complex surface changes and providing a systematic technical solution for land resource monitoring and national spatial planning.

[0098] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

[0099] It should be noted that all formulas in this manual are calculated by removing dimensions and taking their numerical values. The formulas are derived from software simulations based on a large amount of collected data to obtain the most recent real-world results. The preset parameters and thresholds in the formulas are set by those skilled in the art according to the actual situation.

[0100] Although embodiments of the invention have been shown and described, those skilled in the art will understand that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the claims and their equivalents.

Claims

1. A land use change detection system based on remote sensing cloud computing, characterized in that, include: The remote sensing data acquisition module is used to acquire multi-temporal remote sensing image data of the same geographic area and extract spectral reference samples of invariant ground features in each temporal image. The spectral drift compensation module is used to calculate the spectral drift compensation vector for each band based on the spectral response difference of the spectral reference sample in different time phases, and to perform pixel-by-pixel spectral alignment processing on the multi-temporal remote sensing image data based on the spectral drift compensation vector. The response field construction module is used to calculate the pixel-by-pixel spectral difference of the image after spectral alignment, construct an initial change response field, and extract a set of change candidate pixels based on the spectral difference distribution characteristics of the pixels in the initial change response field. The confidence assessment module is used to obtain the neighborhood window of each candidate pixel in the set of candidate pixels with change, calculate the consistency index of change direction and the dispersion of change amplitude of the pixels in the neighborhood window, and construct a spatial confidence evaluation model of candidate pixels based on the consistency index of change direction and the dispersion of change amplitude, and select candidate pixels with spatial confidence greater than a preset confidence threshold as core change pixels. The multi-scale segmentation optimization module is used to calculate the boundary stability coefficient and internal homogeneity coefficient of the change region at different scales by using the core change pixel as the seed point and combining the multi-scale segmentation results, and to determine the optimal segmentation scale based on the weighted fusion result of the boundary stability coefficient and the internal homogeneity coefficient, and to generate change detection patches at the optimal segmentation scale. The sub-pixel boundary optimization module is used to perform sub-pixel-level spectral decomposition on the boundary pixels of the change detection patch, obtain the land cover component ratio of the boundary pixels, correct the boundary position of the change detection patch according to the temporal change rate of the land cover component ratio, and generate a land use change detection result map based on the corrected change detection patch. The step of determining the optimal segmentation scale based on the weighted fusion result of the boundary stability coefficient and the internal homogeneity coefficient includes: The boundary stability coefficients are normalized by scale sequence to obtain normalized boundary stability coefficients; The internal homogeneity coefficients are normalized by scale sequence to obtain normalized internal homogeneity coefficients; The normalized boundary stability coefficient and the normalized internal homogeneity coefficient are weighted and averaged to obtain the comprehensive evaluation index for each scale. Calculate the gradient rate of change of the comprehensive evaluation index on the scale sequence, and identify the inflection point scale where the gradient rate of change turns from positive to negative; If the inflection point scale exists, then the inflection point scale is used as the optimal segmentation scale; if the inflection point scale does not exist, then the scale corresponding to the maximum value of the comprehensive evaluation index is used as the optimal segmentation scale. The step of correcting the boundary position of the change detection patch based on the temporal rate of change of the proportion of the land feature components includes: The difference in the proportion of the land feature components in the boundary pixel between the preceding and following time phases is calculated and denoted as the component change. Set a component change threshold, mark boundary pixels whose component change is greater than the component change threshold as defined change boundary pixels, and mark boundary pixels whose component change is less than the component change threshold as undefined boundary pixels; For the uncertain boundary pixel, calculate its spectral similarity with the core region inside the change detection patch and its spectral similarity with the outer neighborhood of the patch; Based on the comparison results of the spectral similarity, the uncertain boundary pixels with higher spectral similarity to the internal core region are retained in the change detection patch, while the uncertain boundary pixels with higher spectral similarity to the external neighborhood are removed from the change detection patch; Based on the determined boundary pixels and the retained uncertain boundary pixels, the boundary position of the change detection patch is reconstructed.

2. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The step of calculating the spectral drift compensation vector for each band based on the spectral response differences of the spectral reference sample at different time phases includes: Statistical analysis was performed on the spectral reflectance values ​​of the spectral reference sample in each temporal image to extract the median spectral reflectance of each band; Calculate the offset of the median spectral reflectance of each band between the reference time phase and the target time phase, and denote it as the initial drift. Based on the spatial distribution density of the spectral reference sample, the initial drift amount is spatially weighted and corrected to obtain the spatially weighted drift amount; The ratio of the spatially weighted drift to the standard deviation of the spectral reflectance of each band is used as the normalized drift coefficient. The spectral drift compensation vector for each band is constructed based on the normalized drift coefficient.

3. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The step of extracting a set of candidate pixels based on the spectral difference distribution characteristics of pixels in the initial change response field includes: An adaptive threshold segmentation is performed on the spectral difference value of the initial change response field, and pixels with spectral difference greater than the adaptive threshold are marked as preliminary candidate pixels; Calculate the angle between the difference vectors of the preliminary candidate pixels in the multispectral bands to obtain the spectral change direction angle of the preliminary candidate pixels; The number of pixels in the neighborhood of the preliminary candidate pixel that have similar spectral change direction angles is counted and denoted as the directional consistency support. The product of the directional consistency support and the spectral difference of the preliminary candidate pixels is used as the comprehensive change indicator value; The preliminary candidate pixels are arranged in descending order according to the comprehensive change indicator value, and the preliminary candidate pixels with the highest percentage of comprehensive change indicator values ​​are selected to form the change candidate pixel set.

4. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The calculation of the consistency index of the change direction and the dispersion of the change magnitude of pixels within the neighborhood window includes: For each pixel within the neighborhood window, calculate the difference between its preceding and following temporal spectral vectors in the multispectral space, denoted as the pixel change vector; Using the pixel change vector of the candidate pixel as a reference vector, calculate the cosine similarity between the pixel change vector of other pixels in the neighborhood window and the reference vector. The mean of the cosine similarity of all pixels within the neighborhood window is calculated and denoted as the consistency index of the direction of change. Calculate the magnitude of the pixel change vector of all pixels within the neighborhood window, and obtain the magnitude sequence; The standard deviation of the modulus sequence is calculated, and the ratio of the standard deviation to the mean of the modulus sequence is used as the dispersion of the variation range.

5. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The step of constructing a spatial confidence evaluation model for candidate pixels based on the consistency index of the direction of change and the dispersion of the magnitude of change includes: Multiply the change direction consistency index by a preset direction consistency weight to obtain the direction confidence component; The reciprocal of the variation amplitude dispersion is normalized to obtain the normalized reciprocal of dispersion, and the normalized reciprocal of dispersion is multiplied by a preset amplitude stability weight to obtain the amplitude confidence component. Calculate the local maximum indicator factor of the candidate pixel in the initial change response field; The spatial confidence evaluation model for the candidate pixel is constructed by weighting and summing the direction confidence component, the amplitude confidence component, and the local maximum indicator factor. The output value of the spatial confidence evaluation model is the spatial confidence of the candidate pixel.

6. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The calculation of the boundary stability coefficient and internal homogeneity coefficient of the changing region at different scales includes: Within a preset scale range, the image is segmented at multiple scales with a preset scale step size to obtain a set of segmented objects at each scale. For each scale, the overlap rate of the boundary pixels of the segmented object containing the core changed pixels between adjacent scales is calculated and denoted as the boundary overlap rate. The boundary stability coefficient is obtained by multiplying the reciprocal of the perimeter-area ratio of the segmented object by the boundary overlap rate. Calculate the standard deviation of the spectral differences of all pixels within the segmented object, and denot it as the internal spectral fluctuation value; The internal homogeneity coefficient is obtained by multiplying the reciprocal of the internal spectral fluctuation value by the proportion of the core change pixels within the segmented object.

7. The land use change detection system based on remote sensing cloud computing according to claim 6, characterized in that, The statistics include the overlap rate of boundary pixels of the segmented object containing the core changed pixels between adjacent scales, including: Extract the set of boundary pixels of the segmented object at the current scale, and denote it as the current boundary set; Extract the set of boundary pixels of the segmented object containing the core changed pixels at adjacent scales, and denote it as the adjacent boundary set; Calculate the intersection of the current boundary set and the adjacent boundary set to obtain the set of overlapping boundary cells; The ratio of the number of pixels in the overlapping boundary pixel set to the number of pixels in the current boundary set is used as the boundary overlap rate.

8. The land use change detection system based on remote sensing cloud computing according to claim 1, characterized in that, The step of performing sub-pixel-level spectral decomposition on the boundary pixels of the detected change patches to obtain the proportion of ground cover components in the boundary pixels includes: The clean pixels inside the change detection patch and the clean pixels in the outer neighborhood of the patch are extracted and used as the change-class endmembers and the unchange-class endmembers, respectively. Construct an endmember spectral library containing the changed endmembers and the unchanged endmembers; The spectral vector of the boundary pixel is decomposed into endmembers using the fully constrained least squares method to obtain the abundance value of each endmember in the boundary pixel. The sum of the abundance values ​​of the changed endmembers is taken as the proportion of the changed land cover components in the boundary pixel, and the sum of the abundance values ​​of the unchanged endmembers is taken as the proportion of the unchanged land cover components in the boundary pixel.

Images