A geographic information acquisition method and system based on remote sensing images

Abstract

The present application relates to remote sensing image processing technical field, specifically to a kind of geographic information acquisition method and system based on remote sensing image, comprising the following steps: narrowband spectral sensor obtains multiband image data, peak detection extracts reflectivity peak value coordinate and calculates band peak displacement data, Kalman filter analyzes stability and filters stable band set, based on set construction standard spectral feature library, difference calculation spectral similarity value, bayesian classification outputs land use type and generates geographic information distribution map, in the present application, by extracting band peak displacement and carrying out time series stability analysis, stable band is filtered, based on stable reflection feature, spectral expression is constructed, the influence of imaging condition fluctuation on sample is reduced, the type distinguishing degree is strengthened using stable band difference value symbol matching, the accuracy of class judgment based on likelihood inference is improved, so that classification keeps coherence under illumination change scene, the reliability and consistency of spatial distribution result are improved.

Description

Technical Field

[0001] This invention relates to the field of remote sensing image processing technology, and in particular to a method and system for acquiring geographic information based on remote sensing images. Background Technology

[0002] The field of remote sensing image processing technology involves the technical activities of identifying, extracting, classifying, and spatially modeling surface targets using multi-source remote sensing data acquired from aerospace platforms, including multispectral, panchromatic, hyperspectral, radar, and laser point cloud data. Its core aspects include image radiometric correction, geometric correction, mosaicking and fusion, feature extraction, spatial positioning, and attribute inference. Through systematic processing of the spectral, textural, and spatial structural features of remote sensing data, it provides fundamental support for surface information analysis and geospatial analysis. This technical field comprehensively covers the entire process from remote sensing data acquisition, preprocessing, feature construction to feature identification and spatial information representation.

[0003] Traditional geographic information acquisition methods refer to the process of using remote sensing imagery as the main data source, identifying the outlines of ground features and recording their spatial location and attribute information through manual interpretation or rule-based interpretation. This type of method relies on operators to compare the shape, tone, shadow and texture of ground features one by one on the image, and then fill in the classification table and draw the boundary lines in combination with existing ground feature interpretation marks. At the same time, the range of ground features is determined by manual measurement or semi-automatic segmentation based on thresholds. Then, the corresponding point, line and surface data are entered into the geographic information system according to the determined outline coordinates to construct basic spatial elements. When necessary, the interpretation content is corrected through on-site verification. The whole process is completed by means of manual feature recognition, rule comparison and threshold segmentation.

[0004] Existing technologies rely on manual interpretation and rule comparison. The processing is affected by subjective differences in image tone and texture, making it difficult to maintain the stability of spectral performance under illumination fluctuations. Land feature characteristics are prone to drift with changes in imaging conditions. Threshold segmentation is sensitive to temporal noise. Category boundaries shift with changes in local brightness. Spatial representation accuracy is affected by band inconsistencies. Classification results lack coherence across multiple image periods, affecting the reliable representation of land feature distribution. Summary of the Invention

[0005] To address the technical problems of existing technologies that rely on manual interpretation and rule comparison, are affected by subjective differences in image tone and texture, struggle to maintain spectral stability under fluctuating illumination, exhibit drift in ground feature characteristics with changes in imaging conditions, are sensitive to temporal noise in threshold segmentation, have category boundaries shift with changes in local brightness, suffer from interference in spatial representation accuracy due to band inconsistencies, and lack consistency in classification results across multiple image periods, thus affecting the reliable representation of ground feature distribution, this invention provides a geographic information acquisition method based on remote sensing imagery.

[0006] To achieve the above objectives, the present invention employs a geographic information acquisition method based on remote sensing imagery, comprising the following steps:

[0007] S1: Acquire multi-band remote sensing image data of the target area using a narrowband spectral sensor, extract the peak position coordinates of the band reflectance based on the peak detection algorithm, calculate the difference between the peak position coordinates and the reference time node coordinates, and generate the peak displacement data of the band.

[0008] S2: Call the peak displacement data of the band, perform time-series stability analysis using a Kalman filter, calculate the band displacement variance, compare it with the illumination change threshold, filter out stable bands below the threshold, and generate a set of stable bands.

[0009] S3: Based on the stable band set, perform band filtering on remote sensing image data to obtain spectral reflectance data of sample pixels of different land use types, calculate the average reflectance value for multiple types of samples, and construct a standard spectral feature library.

[0010] S4: Call the standard spectral feature library, calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in a band-by-band manner, count the number of difference sign matches, and calculate the spectral similarity value;

[0011] S5: Call the spectral similarity value, perform maximum likelihood estimation based on the Bayesian classifier, select the land use type with the highest similarity value as the pixel classification label, and generate a geographic information distribution map.

[0012] As a further aspect of the present invention, the peak displacement data of the bands includes wavelength displacement, reflectance difference and spectral shift index; the stable band set includes signal-to-noise ratio, radiation stability index and spectral consistency level; the standard spectral feature library includes category mean vector, typical reflectance curve and feature band value; the spectral similarity value includes Euclidean distance, spectral angle and correlation coefficient; and the geographic information distribution map includes classification layer, spatial partition unit and coverage type structure.

[0013] As a further aspect of the present invention, the specific steps of S1 are as follows:

[0014] S101: Based on a narrowband spectral sensor, acquire multi-band remote sensing image data, extract the band reflectance sequence, perform amplitude difference comparison between adjacent points according to the amplitude of the sampling point, and determine the peak inflection position based on the rise and fall changes. Record and aggregate the indexes with amplitude inflection characteristics to generate a band peak index set.

[0015] S102: Based on the band peak index set, call the coordinates of the sampling points corresponding to the index and the coordinates of the reference time nodes to perform item-by-item coordinate difference calculation, and serialize and arrange the difference components according to the index order. Then, perform a one-way consistency check on the direction of the sequence values ​​to obtain the band displacement vector set.

[0016] S103: Based on the band displacement vector set, the difference components in the vector are called and the corresponding band positions are aggregated. The aggregated sequence is then judged according to the numerical segment execution range and encoded in sequence order to generate band peak displacement data.

[0017] As a further aspect of the present invention, the specific steps of S2 are as follows:

[0018] S201: Based on the peak displacement data of the band, perform adjacent difference on the band sequence in the data and smooth the difference sequence through a Kalman filter. Then, matrix-pair the smoothed sequence with the original sequence and extract the variance to generate a band displacement variance matrix.

[0019] S202: Call the band displacement variance matrix, perform amplitude comparison between the variance components corresponding to the bands in the matrix and the illumination change threshold, record the band indices that satisfy the condition that the variance components are lower than the threshold, perform aggregation processing on the indices, and obtain the threshold comparison index set.

[0020] S203: Based on the threshold comparison index set, extract the corresponding index sequence in the peak displacement data of the band, reconstruct the order of the extracted sequence, and simultaneously check the displacement continuity between adjacent samples to generate a stable band set.

[0021] As a further aspect of the present invention, the source of the illumination change threshold is determined by performing discretization processing on continuously collected brightness samples in illumination change scenarios and statistically analyzing the fluctuation range of the discrete samples.

[0022] As a further aspect of the present invention, the specific steps of S3 are as follows:

[0023] S301: Based on the stable band set, the full-band reflectance matrix of the remote sensing image data is filtered. The reflectance sequence of the band position in the full-band reflectance matrix is ​​retrieved according to the stable band index, and the reflectance sequence that does not match the index is removed to obtain the filtered band reflectance matrix.

[0024] S302: Call the filtered band reflectance matrix to perform aggregation processing on the reflectance sequence of land use type sample pixels, classify the reflectance sequence according to the land use type label, and perform numerical accumulation and counting on the same type of sequence at the band position to generate the land type mean reflectance sequence;

[0025] S303: Based on the mean reflectance sequence of the land type, perform key-value aggregation on the mean sequence corresponding to the land use type, using the type label as the key and the mean sequence as the value, construct a mapping structure and perform structured sequence organization to establish a standard spectral feature library.

[0026] As a further aspect of the present invention, the specific steps of S4 are as follows:

[0027] S401: Call the standard spectral feature library, and based on the spectral reflectance value of the pixel in the region to be classified, perform a difference operation on the corresponding band values ​​of the reflectance values ​​and the standard reflectance values ​​of the land use type in the standard spectral feature library, and record the difference sign to generate a band difference sign sequence.

[0028] S402: Call the band difference symbol sequence, perform consistency determination on the two types of symbols at the corresponding band positions for the band symbol reference sequence in the standard reflectance value of land use type, convert the consistency status into binary code and accumulate it into a quantity value, and obtain the symbol matching quantity accumulation value;

[0029] S403: Based on the cumulative value of the symbol matching quantity, call the number of bands of the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type to perform a ratio operation on the quantity value and the number of bands, and use the ratio as the similarity characterization value to obtain the spectral similarity value.

[0030] As a further aspect of the present invention, the specific steps of S5 are as follows:

[0031] S501: Based on the spectral similarity value, perform numerical normalization on the similarity parameter group of land use type, and calculate the ratio of each normalized similarity parameter group to the spectral response benchmark value of the same source pixel. Then, perform weighted aggregation on the ratio results to generate a likelihood probability matrix.

[0032] S502: Call the conditional probability parameter group in the likelihood probability matrix, and perform maximum likelihood determination on the multi-type conditional probability parameter group of each pixel based on the Bayesian classifier. Record the maximum numerical index corresponding to the determination result as the land use type number and perform rasterization encoding to obtain the pixel type index set.

[0033] S503: Based on the land use type sequence number in the pixel type index set, perform spatial mapping on the coordinate group of pixels with the same sequence number, perform coordinate reprojection, and perform regional stitching on the reprojected pixel partitions to obtain a geographic information distribution map.

[0034] As a further aspect of the present invention, the spectral response benchmark value is based on standard sample pixels verified by field testing of land cover categories, the statistical results of stable reflectance of land use types in multispectral bands, and the calibrated spectral response value is obtained by interval averaging of the spectral reflectance of sample pixels of the same category after removing spectral noise interference.

[0035] A geographic information acquisition system based on remote sensing imagery, comprising:

[0036] The spectral positioning module acquires multi-band remote sensing image data of the target area through a narrowband spectral sensor, extracts the position coordinates of the peak reflectance of the band based on the peak detection algorithm, calculates the difference between the peak reflectance of the band and the coordinates of the reference time node, and generates the peak displacement data of the band.

[0037] The band stabilization module calls the peak displacement data of the band, performs time-series stability analysis using a Kalman filter, calculates the band displacement variance, compares it with the illumination change threshold, filters out stable bands below the threshold, and generates a set of stable bands.

[0038] The type feature module performs band filtering on remote sensing image data based on the stable band set, obtains spectral reflectance data of sample pixels of different land use types, calculates the average reflectance value for multiple types of samples, and constructs a standard spectral feature library.

[0039] The similarity calculation module calls the standard spectral feature library to calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in a band-by-band manner, counts the number of difference sign matches, and calculates the spectral similarity value.

[0040] The classification output module calls the spectral similarity value, performs maximum likelihood estimation based on a Bayesian classifier, selects the land use type with the highest similarity value as the pixel classification label, and generates a geographic information distribution map.

[0041] Compared with the prior art, the advantages and positive effects of the present invention are as follows:

[0042] In this invention, by extracting the peak displacement of bands and screening bands with steady-state performance through temporal stability analysis, the subsequent spectral expression is constructed based on relatively stable reflection characteristics, reducing the disturbance of sample statistics to the fluctuation of imaging conditions. By matching the difference sign of stable bands, the distinguishability of different types is enhanced, and the certainty of spectral discrimination is improved in the category judgment based on likelihood inference. This ensures that the classification process remains consistent under the scenario of changing illumination and improves the reliability and consistency of spatial distribution results. Attached Figure Description

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

[0044] Figure 1 This is a schematic diagram of the steps of the present invention;

[0045] Figure 2 This is a detailed schematic diagram of S1 of the present invention;

[0046] Figure 3 This is a detailed schematic diagram of S2 of the present invention;

[0047] Figure 4 This is a detailed schematic diagram of S3 of the present invention;

[0048] Figure 5 This is a detailed schematic diagram of S4 of the present invention;

[0049] Figure 6 This is a detailed schematic diagram of S5 of the present invention;

[0050] Figure 7 This is a system module diagram of the present invention. Detailed Implementation

[0051] The technical solution of the present invention will now be described with reference to the accompanying drawings.

[0052] In embodiments of the present invention, words such as "exemplarily," "for example," etc., are used to indicate that something is an example, illustration, or description. Any embodiment or design described as "exemplary" in the present invention should not be construed as being more preferred or advantageous than other embodiments or designs. Specifically, the use of the word "exemplary" is intended to present the concept in a concrete manner. Furthermore, in embodiments of the present invention, the meaning expressed by "and / or" can be both, or either one.

[0053] In the embodiments of this invention, the terms "image" and "picture" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning. Similarly, the terms "of," "corresponding (relevant)," and "corresponding" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning.

[0054] In this embodiment of the invention, sometimes the subscript such as W1 is written in a non-subscript form such as W1. When the difference is not emphasized, the meaning they express is the same.

[0055] To make the technical problems, technical solutions and advantages of the present invention clearer, a detailed description will be given below in conjunction with the accompanying drawings and specific embodiments.

[0056] Please see Figure 1 This invention provides a method for geographic information acquisition based on remote sensing imagery, comprising the following steps:

[0057] S1: Acquire multi-band remote sensing image data of the target area using a narrowband spectral sensor, extract the peak position coordinates of the band reflectance based on the peak detection algorithm, calculate the difference between the peak position coordinates and the reference time node coordinates, and generate the peak displacement data of the band.

[0058] S2: Call the peak displacement data of the band, use the Kalman filter to perform time-series stability analysis, calculate the variance of the band displacement, compare it with the illumination change threshold, filter out stable bands below the threshold, and generate a set of stable bands.

[0059] S3: Based on a stable band set, band filtering is performed on remote sensing image data to obtain spectral reflectance data of sample pixels of different land use types. The average reflectance value is calculated for multiple types of samples, and a standard spectral feature library is constructed.

[0060] S4: Call the standard spectral feature library, calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in each band, count the number of difference sign matches, and calculate the spectral similarity value;

[0061] S5: Call the spectral similarity value, perform maximum likelihood estimation based on the Bayesian classifier, select the land use type with the highest similarity value as the pixel classification label, and generate a geographic information distribution map.

[0062] The peak displacement data includes wavelength displacement, reflectance difference, and spectral shift index; the stable band set includes signal-to-noise ratio, radiometric stability index, and spectral consistency level; the standard spectral feature library includes category mean vector, typical reflectance curve, and characteristic band values; the spectral similarity values ​​include Euclidean distance, spectral angle, and correlation coefficient; and the geographic information distribution map includes classification layers, spatial partitioning units, and coverage type structure.

[0063] Please see Figure 2 The specific steps of S1 are as follows:

[0064] S101: Based on a narrowband spectral sensor, acquire multi-band remote sensing image data, extract the band reflectance sequence, perform amplitude difference comparison between adjacent points according to the amplitude of the sampling point, and determine the peak inflection position based on the rise and fall changes. Record and aggregate the indexes with amplitude inflection characteristics to generate a band peak index set.

[0065] A pushbroom hyperspectral imager mounted on a low-altitude UAV platform was selected as the narrowband spectral sensor, with a spectral response range of 400nm to 1000nm and a spectral resolution better than 5nm. Pushbroom imaging of the target area was performed under operating conditions at a flight altitude of 200 meters to acquire the raw digital quantized data matrix. Radiometric calibration was then performed on this raw data, using preset absolute radiometric calibration coefficients. and dark current bias Convert digital quantization values ​​to entrance pupil radiance. Combined with the downlink irradiance data acquired simultaneously Using the formula Calculate pixel-by-pixel, and extract the band reflectance sequence corresponding to each pixel. This sequence contains 120 bands. Adjacent difference comparisons based on sampling point amplitudes are performed on this sequence, with a sliding window size of 3. For the [number]th band in the sequence... reflectivity of each band Retrieve the reflectivity of its preceding bands respectively With subsequent band reflectivity Set the noise suppression threshold for determining the effective peak value. The threshold It is not a fixed value, but is dynamically set based on 5% of the current pixel's average reflectance across all wavelengths. For example, the current pixel's average reflectance across all wavelengths... Then calculate Execute the peak-shaped transition judgment logic; if the condition is met... That is, to determine that the position is the turning point of the wave peak, if the condition is met. This indicates a trough inflection point; here, only peak features are extracted from sequence fragments. For example, the difference between the median value of 0.18 and the values ​​on both sides is 0.03 and 0.02 respectively, both of which are greater than the threshold of 0.012. Therefore, it is determined that there is a value inflection feature at the index position. Traversing the entire The sequence identifies the band index values ​​that meet the peak conditions. For example, peaks are detected at bands 15, 42, and 88. The index values ​​are then... Record and aggregate the band peak indexes to construct a one-dimensional array. This completes the structured extraction of spectral features for single pixels.

[0066] S102: Based on the band peak index set, the coordinate difference between the sampling point coordinates corresponding to the index and the reference time node coordinates is calculated item by item, and the difference components are serialized and arranged in the index order. Then, a one-way consistency check is performed on the direction of the sequence values ​​to obtain the band displacement vector set.

[0067] Call the generated band peak index set Set reference time nodes The baseline image acquisition time is [time], and the current acquisition time is [time]. For each index in the index set ,exist Extract the pixel coordinates of the sampling point corresponding to this band from the geometric correction model of the remote sensing image at that time. At the same time Retrieve reference coordinates for corresponding ground features from the reference image database at a given time. Perform a step-by-step coordinate difference calculation, the formula is as follows: The calculation yields results including lateral deviation. Longitudinal deviation To verify the usability of the data, specific implementation data as shown in Table 1 needs to be introduced. Table 1 lists the coordinate difference data extracted based on the band peak index. Based on the data in the table, for index 15, the difference is calculated. , Similarly, calculate the difference at index 42. The difference at index 88 The difference components are serialized and arranged in index order to obtain the initial vector sequence. Then, a one-way consistency check is performed on the numerical direction of the sequence, and the arithmetic mean offset of the vectors in the sequence in the X and Y axes is calculated. , Set a consistency judgment tolerance Item-by-item testing and Is it less than For the vector corresponding to index 42, the X-axis deviation Y-axis deviation To meet consistency requirements, if a vector deviates beyond the tolerance, it is removed from the sequence. Finally, the vectors that pass the test are retained to form the band displacement vector set. .

[0068] Table 1: Band Sampling Coordinate Difference Data Table

[0069]

[0070] As shown in Table 1, this table details the spatial location offset corresponding to the peak indices of three typical bands. By comparing the current coordinates... Reference time coordinates It can be intuitively observed that the transverse and longitudinal difference components of band 15 are... The difference component in band 42 increases to The difference component of band 88 is The specific values ​​reveal that the geometric distortions occurring at different spectral feature locations during imaging are non-uniform, providing accurate data support for subsequent one-way consistency checks and the construction of displacement vector sets.

[0071] S103: Based on the band displacement vector set, the difference components in the vector are called and the corresponding band positions are aggregated. The aggregated sequence is then judged according to the numerical segment range and the sequence order is encoded to generate the band peak displacement data.

[0072] Based on the output band displacement vector set Iterate through each element in the set and extract the difference components from the vector. Corresponding original band position index Perform tuple aggregation to form aggregated items. The aggregated sequence values ​​are used to determine the range of values. Three value segments are set to quantify the degree of displacement. Segment 1 is set as... The corresponding code is "L", and segment 2 is set to... The corresponding code is "M", and segment 3 is set to... The corresponding code is "H". Outliers exceeding 1.0 are marked as "E" for aggregate items. Encode, where It falls into segment 2, coded as "M". Falling into segment 1, the code is "L", and combining them yields the code string "15-ML". Similarly, for... , Belongs to segment 3 (boundary value is assigned to higher segment), coded "H". Belonging to segment 2, coded as "M", we get "42-HM". Encoded as "88-ML", the above independent codes are connected using a sequential encoding method with a specific separator "|" to generate the final band peak displacement data string. The results indicate that at specific spectral feature points in a certain band, ground targets have undergone low to medium degree of geometric distortion or displacement relative to the reference time, and the magnitude of displacement between different bands exhibits a non-uniform distribution characteristic. This string will be written into a specific tag field of the image metadata file as a control parameter for subsequent geometric fine correction.

[0073] Please see Figure 3 The specific steps of S2 are as follows:

[0074] S201: Based on the peak displacement data of the bands, the band sequence in the data is subjected to adjacent difference and smoothed by Kalman filter. Then, the smoothed sequence is matrix-paired with the original sequence and the variance is extracted to generate the band displacement variance matrix.

[0075] Retrieve the generated peak displacement data of the band Extract the displacement amplitude sequence contained therein, assuming the displacement amplitude sequence is... Perform adjacent difference operation on the sequence, the calculation formula is as follows: The difference sequence is obtained. Construct a Kalman filter model and define the state transition matrix. Observation matrix Noise covariance Measure noise covariance Initialize state estimates and error covariance For the first element in the difference sequence Perform the time update step using the formula Calculate prior state estimation Using the formula Calculate the prior error covariance The advantage of the formula lies in its ability to introduce process noise. To prevent the filter from over-converging and thus maintain sensitivity to dynamic changes, a measurement update step is then performed, using the formula... Calculate the Kalman gain and substitute it into the numerical calculation. Using the formula Update state estimate Using the formula Update error covariance Repeat the above iterative process to smooth the entire sequence and obtain a smoothed sequence. The smoothed difference sequence is matrix-paired with the original difference sequence, and the residual sequence between them is calculated. For each band, the residual subset is used with the formula Extracting variance, assuming the residual sequence for a certain band is... Calculate the mean Substitute into the formula to calculate the variance The results indicate that there are minute random jitters in the displacement changes of this band between adjacent time points, and the amplitude of these jitters has been quantified into specific statistical values. The variances calculated for the bands are arranged by band index to generate a band displacement variance matrix. .

[0076] S202: Call the band displacement variance matrix, perform amplitude comparison between the variance components corresponding to the bands in the matrix and the illumination change threshold, record the band indices that satisfy the condition that the variance components are lower than the threshold, perform aggregation processing on the indices to obtain the threshold comparison index set;

[0077] Calling the band displacement variance matrix This includes the displacement stability metric of the waveband, such as... Set the threshold for changes in illumination. The threshold needs to be obtained based on the data collected in actual lighting change scenarios. A set of image brightness mean samples continuously collected under dynamically changing lighting conditions is selected, as shown in Table 2. Brightness sequences are then extracted. Perform fluctuation range statistical analysis on the discrete sample sequence and calculate the sample standard deviation:

[0078] ;

[0079] Set threshold decision coefficient Using the formula Calculate the baseline threshold, substitute the values ​​to obtain... This result indicates that the upper limit of allowable normal background noise fluctuations caused by ambient light has been quantized to 0.00311. Subsequently, an amplitude comparison was performed, comparing the variance components in the matrix one by one with this threshold. For band index 15, the variance... If the condition is met, record index 15; for band index 42, the variance... If the condition is not met, it is considered an unstable band severely affected by light interference; for band index 88, the variance... If the condition is met, record index 88. After traversing all the indices, record the indices that meet the condition. Perform aggregation processing to obtain a threshold comparison index set. .

[0080] Table 2: Sample Data Table of Light Intensity Variation

[0081]

[0082] Table 2 lists sample data of light intensity changes. Brightness sequences are extracted based on the data in the table, and statistical analysis of fluctuation range is performed on the discrete sample sequences.

[0083] S203: Based on the threshold comparison index set, extract the corresponding index sequence from the peak displacement data of the band, reconstruct the order of the extracted sequence, and simultaneously check the displacement continuity between adjacent samples to generate a stable band set.

[0084] Compare index sets based on thresholds From the original peak displacement data of the band Extract the corresponding displacement vector data to obtain the sequence. ,in Corresponding displacement vector , Corresponding displacement vector The extracted sequence is reconstructed in ascending order based on the wavelength values ​​of the bands, and the continuity check between adjacent samples is performed simultaneously, defining a continuity check function. The advantage of the formula lies in its ability to accurately quantify the degree of displacement consistency between two discrete waveband feature points within a spatial geometric domain using the Euclidean distance principle, and to set a continuity judgment threshold. This threshold is set based on the spatial resolution limits of adjacent bands of the sensor. and Substitute into the formula to calculate the spacing Execute the judgment logic, due to the calculation result This indicates that the displacement change between band 15 and band 88 maintains a continuous and smooth transition characteristic in spatial geometry, without any abrupt tortuosity. This result shows that the selected band group has a high degree of inherent consistency in geometric imaging characteristics. Therefore, both bands are determined to be effective stable bands. If the calculated result is greater than the threshold, distance abrupt points need to be removed. Finally, the band indices that pass the verification are confirmed, and a stable band set is generated. .

[0085] Please see Figure 4 The specific steps of S3 are as follows:

[0086] S301: Based on the stable band set, the full-band reflectance matrix of the remote sensing image data is filtered. The reflectance sequence of the band position in the full-band reflectance matrix is ​​retrieved according to the stable band index, and the reflectance sequence that does not match the index is removed to obtain the filtered band reflectance matrix.

[0087] Call the generated stable band set As a key feature index source, the full-band reflectance matrix corresponding to the original hyperspectral remote sensing image of the currently acquired area is loaded. Let the dimension of this matrix be 1. ,in Represents the number of pixels in the image height. Represents the image width in pixels. Represents the total number of original bands of the sensor, and the position of any pixel in the matrix. The corresponding full-band spectral vector is expressed as:

[0088] ;

[0089] Initialize an empty filter band matrix container. Its preset dimension is Execute the band index matching logic, traversing each band index in the original matrix. For the current index Combine it with the stable band set The existence of elements in the data is compared, and the determination logic is as follows: ,like If true, then the band is locked. exist Complete two-dimensional slice data in, if For example, false. or If the data in that band is invalid, it will be marked as invalid and removed; it will not be written into the new matrix. The data will be stored in pixels. For example, its original spectral vector contains 120 data values, of which After indexing and After matching and filtering, the extracted retained vector is Map the vector to Corresponding coordinate position The channel depth direction, for the whole map The above dimensionality reduction mapping is performed pixel by pixel to complete the structured transformation from high-dimensional redundant data to key feature data, ultimately obtaining a reflectance matrix of selected bands containing only stable spectral features. .

[0090] S302: Call the filter band reflectance matrix to perform aggregation processing on the reflectance sequence of land use type sample pixels, classify the reflectance sequence according to the land use type label, and perform numerical accumulation and counting on the band position of the same type of sequence to generate the land type mean reflectance sequence.

[0091] Call the filter band reflectivity matrix Based on the regions of interest (ROIs) delineated using ground-based survey data, a set of sample pixels with clear land use type labels was extracted. Specific sample data is shown in Table 3. Based on the data in the table, two typical land cover labels, "coniferous forest" and "lake / water body," were identified. For the "coniferous forest" type, the spectral sequences of its associated sample pixels were extracted, and an aggregation accumulator was set. and counter Iterate through the data rows in Table 3 with sample IDs 1 and 2. For sample 1, the reflectance of band 15 is 0.12 and the reflectance of band 88 is 0.45. Perform an accumulation operation. The counter increments automatically. For sample 2, the reflectance of band 15 is 0.14, and the reflectance of band 88 is 0.47. This is continued to accumulate. The counter is updated to 2, and the mean calculation formula is executed:

[0092] ;

[0093] in the formula The mean reflectance vector representing a specific land cover type. The total number of samples, For the first The reflectance vector of each sample is used in this formula to eliminate single-point sampling noise and extract general spectral characteristics of ground features through arithmetic averaging. The cumulative value is then substituted into the formula to calculate the mean sequence of "coniferous forest". Similarly, for the "lake water body" type, samples 3 and 4 were traversed, and data from band 15 were accumulated. Accumulated data from 88 bands Calculate the mean sequence The results show that coniferous forests have high reflectance characteristics in band 88 (near-infrared region), while water bodies exhibit low reflectance characteristics in both bands. The numerical results accurately quantify the differences in physical attributes of different land features. The calculated mean vectors of multiple types are arranged in order to generate a land type mean reflectance sequence.

[0094] Table 3: Pixel Reflectance Data for Land Use Types

[0095]

[0096] As shown in Table 3, Table 3 lists the reflectance values ​​of sample pixels of different land use types at the selected bands.

[0097] S303: Based on the mean reflectance sequence of land type, perform key-value aggregation on the mean sequence corresponding to land use type, use type label as key and mean sequence as value, construct mapping structure and perform structured sequence processing to establish a standard spectral feature library;

[0098] Based on the generated land type mean reflectance sequence, a structured aggregation operation in key-value pair form is performed. The data structure of the standard spectral feature library is defined as a hash map, with "land use type label" as the unique key and "mean reflectance sequence" as the corresponding value. For the "coniferous forest" category, the key-value pair Entry1={Key:"Coniferous_Forest", Value:[0.13, 0.46]} is constructed; for the "lake water" category, the key-value pair Entry2={Key:"Lake_Water", Value:[0.025, 0.045]} is constructed. To ensure the robustness of the feature library, a normalization check needs to be performed on the sequence values, and the effective range of reflectance is set to... If any value exceeds this range, it is marked as an anomaly and correction logic is triggered. Since the current calculated values ​​all fall within the valid range, no correction is needed. Next, the sequence is formatted and encoded, retaining the values ​​to four decimal places to standardize precision. Updated to , Updated to The key-value pairs are sorted and organized according to the coding order of land cover types. The organized entries are then encapsulated into a JSON object: Library={"100":[0.1300, 0.4600], "200":[0.0250, 0.0450]}. Finally, this structured object is serialized into a binary stream file, and header metadata including a creation timestamp and a list of applicable band indices is added. Including sensor model identification, the standard spectral feature library was established.

[0099] Please see Figure 5 The specific steps of S4 are as follows:

[0100] S401: Call the standard spectral feature library, and based on the spectral reflectance values ​​of the pixels in the region to be classified, perform a difference operation on the corresponding band values ​​of the reflectance values ​​and the standard reflectance values ​​of land use types in the standard spectral feature library, and record the difference sign to generate a band difference sign sequence.

[0101] The established standard spectral feature library is used to extract the standard mean reflectance sequence corresponding to the "coniferous forest" to be matched. Simultaneously, the spectral reflectance value of a specific pixel in the region to be classified in the corresponding spectral band is obtained. Assuming the pixel is affected by atmospheric scattering, its observed value is Perform a step-by-step difference calculation of the band values ​​using the formula. Calculate the difference; for band 15, calculate... For band 88, calculate The sign of the difference is recorded based on the positive or negative value of the calculation result, and a sign function is defined. ,when The time marker is "+", when The time marker is "-", thereby generating a band difference symbol sequence containing directional information. This sequence visually reflects the deviation direction of the pixel under test relative to the standard spectral curve. "+" indicates that the observation has an overestimation bias, and "-" indicates that the observation has an underestimation bias. This provides basic data support for subsequent matching and fault tolerance analysis. The symbol sequence is temporarily stored in the memory buffer for use in the next step.

[0102] S402: Call the band difference symbol sequence, perform consistency determination on the two types of symbols at the corresponding band positions for the band symbol reference sequence in the standard reflectance value of land use type, convert the consistency status into binary code and accumulate it into a quantity value, and obtain the cumulative value of symbol matching quantity;

[0103] Calling the band difference symbol sequence It also reads the "band symbol reference sequence" for this land cover type from the feature library metadata. This reference sequence defines the allowed deviation patterns, and the ideal reference sequence is set as... Where "0" represents a state of no significant difference within the tolerance range, the consistency determination logic is executed, and a dynamic tolerance threshold needs to be introduced. To bridge the logical relationships between symbols, and set The threshold setting is based on a comprehensive consideration of the sensor signal-to-noise ratio and the variance within land cover types. For example, the variance within coniferous forest types is 0.015, so the threshold is slightly higher than this value to ensure robustness. For band 15, the absolute value of the difference is... ,because The difference in this band is determined to be within the acceptable range, and the consistency status is determined to be "compliant". For band 88, the absolute value of its difference is... ,because Similarly, if the difference exceeds 0.02, it is judged as "compliant". If the difference exceeds 0.02, it is judged as "incompatible". The above consistency status is converted into a binary code, and compliance is marked as... The non-compliance mark is To obtain the matching state sequence .

[0104] Table 4: Spectral Feature Matching Calculation Data Table

[0105]

[0106] As shown in Table 4, Table 4 lists the matching operation data between the pixels to be classified and the standard spectrum. Based on the data in the table, numerical accumulation operation is performed on the binary encoded sequence using the formula... Calculate the cumulative value of symbol matching. This value quantifies the degree of fit between the measured pixel and the standard spectrum in the key feature bands.

[0107] S403: Based on the cumulative value of symbol matching, the number of bands of the spectral reflectance value of the pixels in the area to be classified and the standard reflectance value of the land use type are called to perform a ratio calculation on the number of bands and the number of bands, and the ratio is used as the similarity representation value to obtain the spectral similarity value;

[0108] Cumulative value based on the number of symbol matches in the output Obtain the total number of bands involved in the calculation. Perform ratio calculations to normalize the similarity measure, using the formula... Calculate the spectral similarity value and substitute it into the numerical calculation. The results indicate that, within the selected key band set, the pixel to be classified exhibits 100% characteristic agreement with the standard spectrum of "coniferous forest." To verify the exclusivity of the classification, another set of comparative calculations needs to be introduced, assuming that the pixel is compared with the standard value of "lake water body." The comparison was performed, and the differences were calculated as follows: and All of them far exceed the threshold of 0.02, therefore the corresponding cumulative values ​​are... Similarity By comparing the two sets of similarity values The spectral attributes of the pixels in this area are determined to significantly belong to the coniferous forest category. The final calculated similarity value will then be used to determine the forest type. The classification confidence parameter is output to complete the attribute recognition process of a single pixel.

[0109] Please see Figure 6 The specific steps of S5 are as follows:

[0110] S501: Based on spectral similarity values, perform numerical normalization on the similarity parameter group of land use types, and calculate the ratio of each normalized similarity parameter group to the spectral response benchmark value of the same source pixel. Then, perform weighted aggregation on the ratio results to generate a likelihood probability matrix.

[0111] Based on the output spectral similarity values, the similarity components of the pixels to be classified for "coniferous forest" are extracted. Similarity component for "lake water body" First, numerical normalization is performed to construct the basis of the probability distribution, using the formula... Since the sum is 1 in this example, the normalization result remains unchanged. Then, the spectral response baseline value of the homologous pixels is obtained. The establishment of this value must strictly follow statistical laws, and samples are collected within the "coniferous forest" standard sample area verified through field surveys. Calculate the average reflectance of each pixel of the same type in the characteristic band using its spectral data. and standard deviation Set the elimination interval as After removing five noisy samples that fall outside the specified interval, the interval average is calculated for the remaining samples to obtain a benchmark value characterizing the spectral stability of this type of land cover. Similarly, the baseline value of the water body is obtained. Perform ratio calculations to introduce stability weights, using the formula The normalized similarity is corrected, and the advantage of the formula is that it uses the inherent spectral stability parameters of ground features to weight and correct the observational similarity, reducing the risk of misjudgment caused by the spectral variability of ground features themselves. Substituting the values ​​into the numerical calculation yields... , Weighted aggregation is performed on the corrected ratio results. In practice, this involves normalization again to generate the likelihood probability, using the formula... Calculate the final likelihood and substitute the values. The results indicate that, after incorporating prior knowledge of the spectral stability of ground features, the confidence level of this pixel belonging to coniferous forest is further confirmed to be 100%, generating a likelihood probability matrix. .

[0112] Table 5: Parameters for Likelihood Probability Calculation

[0113]

[0114] As shown in Table 5, Table 5 lists the parameters for calculating the likelihood probability, and details the derivation process from the original similarity to the final likelihood probability.

[0115] S502: Call the conditional probability parameter set in the likelihood probability matrix, and perform maximum likelihood determination on the multi-type conditional probability parameter set of each pixel based on the Bayesian classifier. Record the maximum numerical index corresponding to the determination result as the land use type number and perform rasterization encoding to obtain the pixel type index set.

[0116] Call the likelihood probability matrix Extract the conditional probability parameter set contained therein. Construct a Bayesian classifier decision model and set prior probabilities. This probability value is set based on the area proportion of different land use types in the original land use map. It is assumed that the original proportion of coniferous forest in the region is 60%, and water bodies account for 40%, i.e. Execute the maximum likelihood decision logic and compute the posterior discriminant function. , here This refers to the likelihood probability calculated in the previous steps, which is then substituted into the numerical calculation of the discriminant score for coniferous forests. Water body identification Perform numerical comparison and determination using formulas. Extract the category index corresponding to the maximum value, because The determination result points to the "coniferous forest" category, and its standard classification number is recorded. Perform rasterization encoding to initialize a classification result matrix with the same resolution as the original image. Locate the row and column coordinates of the current cell. , classify the serial number If a case with equal discriminant scores is found when writing data into this matrix cell, it is marked as a category to be determined. It then triggers the neighborhood-assisted determination logic, traverses all image pixels to complete the above calculations and assignments, and finally generates a pixel type index set containing complete land cover category distribution information.

[0117] S503: Based on the land use type sequence number in the pixel type index set, perform spatial mapping on the coordinate group of pixels with the same sequence number, perform coordinate reprojection, and perform regional stitching on the reprojected pixel partitions to obtain a geographic information distribution map.

[0118] Based on cell type index set Identify and extract the tags as The pixel coordinate sequence of (coniferous forest), for example Spatial mapping and coordinate reprojection are performed on the coordinate group of pixels with the same sequence number. The georeferenced transformation parameters (GeoTransform) of the original remote sensing image are read, including the projection coordinates of the upper left corner. and pixel resolution Meters, using the affine transformation formula:

[0119] ;

[0120] Calculate the geospatial coordinates of the center of each pixel in the formula. These are the projected geographic coordinates. These are the matrix row and column indices. Matrix operations aim to accurately map discrete grid points in image space to their physical locations on the Earth's surface. Substitute the example cell... Calculated rice, After reprojecting a single image, region stitching is performed on multiple blocks of data covering a wide area. A weighted average fusion algorithm is used to process overlapping areas, and an overlap bandwidth is set. Pixels, using the inverse distance weighting formula The pixel category values ​​at the smooth seams are transitioned, ultimately generating a land use geographic information distribution map with a unified geographic coordinate system and seamless connections.

[0121] Please see Figure 7A geographic information acquisition system based on remote sensing imagery, comprising:

[0122] The spectral positioning module acquires multi-band remote sensing image data of the target area through a narrowband spectral sensor, extracts the position coordinates of the peak reflectance of the band based on the peak detection algorithm, calculates the difference between the peak reflectance of the band and the coordinates of the reference time node, and generates the peak displacement data of the band.

[0123] The band stabilization module calls the peak displacement data of the band, performs time-series stability analysis using a Kalman filter, calculates the band displacement variance, compares it with the illumination change threshold, filters out stable bands below the threshold, and generates a set of stable bands.

[0124] The type feature module performs band filtering on remote sensing image data based on the stable band set, obtains spectral reflectance data of sample pixels of different land use types, calculates the average reflectance value for multiple types of samples, and constructs a standard spectral feature library.

[0125] The similarity calculation module calls the standard spectral feature library to calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in a band-by-band manner, counts the number of difference sign matches, and calculates the spectral similarity value.

[0126] The classification output module calls the spectral similarity value, performs maximum likelihood estimation based on a Bayesian classifier, selects the land use type with the highest similarity value as the pixel classification label, and generates a geographic information distribution map.

[0127] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A method for acquiring geographic information based on remote sensing imagery, characterized in that, Includes the following steps: S1: Acquire multi-band remote sensing image data of the target area using a narrowband spectral sensor, extract the peak position coordinates of the band reflectance based on the peak detection algorithm, calculate the difference between the peak position coordinates and the reference time node coordinates, and generate the peak displacement data of the band. S2: Call the peak displacement data of the band, perform time-series stability analysis using a Kalman filter, calculate the band displacement variance, compare it with the illumination change threshold, filter out stable bands below the threshold, and generate a set of stable bands. S3: Based on the stable band set, perform band filtering on remote sensing image data to obtain spectral reflectance data of sample pixels of different land use types, calculate the average reflectance value for multiple types of samples, and construct a standard spectral feature library. S4: Call the standard spectral feature library, calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in a band-by-band manner, count the number of difference sign matches, and calculate the spectral similarity value; S5: Call the spectral similarity value, perform maximum likelihood estimation based on the Bayesian classifier, select the land use type with the highest similarity value as the pixel classification label, and generate a geographic information distribution map; The specific steps for S3 are as follows: S301: Based on the stable band set, the full-band reflectance matrix of the remote sensing image data is filtered. The reflectance sequence of the band position in the full-band reflectance matrix is ​​retrieved according to the stable band index, and the reflectance sequence that does not match the index is removed to obtain the filtered band reflectance matrix. S302: Call the filtered band reflectance matrix to perform aggregation processing on the reflectance sequence of land use type sample pixels, classify the reflectance sequence according to the land use type label, and perform numerical accumulation and counting on the same type of sequence at the band position to generate the land type mean reflectance sequence; S303: Based on the mean reflectance sequence of the land type, perform key-value aggregation on the mean sequence corresponding to the land use type, using the type label as the key and the mean sequence as the value, construct a mapping structure and perform structured sequence organization to establish a standard spectral feature library; The specific steps of S4 are as follows: S401: Call the standard spectral feature library, and based on the spectral reflectance value of the pixel in the region to be classified, perform a difference operation on the corresponding band values ​​of the reflectance values ​​and the standard reflectance values ​​of the land use type in the standard spectral feature library, and record the difference sign to generate a band difference sign sequence. S402: Call the band difference symbol sequence, perform consistency determination on the two types of symbols at the corresponding band positions for the band symbol reference sequence in the standard reflectance value of land use type, convert the consistency status into binary code and accumulate it into a quantity value, and obtain the symbol matching quantity accumulation value; S403: Based on the cumulative value of the symbol matching quantity, call the number of bands of the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type to perform a ratio operation on the quantity value and the number of bands, and use the ratio as the similarity characterization value to obtain the spectral similarity value.

2. The geographic information acquisition method based on remote sensing imagery according to claim 1, characterized in that, The peak displacement data of the bands includes wavelength displacement, reflectance difference and spectral shift index; the stable band set includes signal-to-noise ratio, radiation stability index and spectral consistency level; the standard spectral feature library includes category mean vector, typical reflectance curve and feature band value; the spectral similarity value includes Euclidean distance, spectral angle and correlation coefficient; and the geographic information distribution map includes classification layer, spatial partition unit and coverage type structure.

3. The geographic information acquisition method based on remote sensing imagery according to claim 1, characterized in that, The specific steps of S1 are as follows: S101: Based on a narrowband spectral sensor, acquire multi-band remote sensing image data, extract the band reflectance sequence, perform amplitude difference comparison between adjacent points according to the amplitude of the sampling point, and determine the peak inflection position based on the rise and fall changes. Record and aggregate the indexes with amplitude inflection characteristics to generate a band peak index set. S102: Based on the band peak index set, call the coordinates of the sampling points corresponding to the index and the coordinates of the reference time nodes to perform item-by-item coordinate difference calculation, and serialize and arrange the difference components according to the index order. Then, perform a one-way consistency check on the direction of the sequence values ​​to obtain the band displacement vector set. S103: Based on the band displacement vector set, the difference components in the vector are called and the corresponding band positions are aggregated. The aggregated sequence is then judged according to the numerical segment execution range and encoded in sequence order to generate band peak displacement data.

4. The geographic information acquisition method based on remote sensing imagery according to claim 3, characterized in that, The specific steps of S2 are as follows: S201: Based on the peak displacement data of the band, perform adjacent difference on the band sequence in the data and smooth the difference sequence through a Kalman filter. Then, matrix-pair the smoothed sequence with the original sequence and extract the variance to generate a band displacement variance matrix. S202: Call the band displacement variance matrix, perform amplitude comparison between the variance components corresponding to the bands in the matrix and the illumination change threshold, record the band indices that satisfy the condition that the variance components are lower than the threshold, perform aggregation processing on the indices, and obtain the threshold comparison index set. S203: Based on the threshold comparison index set, extract the corresponding index sequence in the peak displacement data of the band, reconstruct the order of the extracted sequence, and simultaneously check the displacement continuity between adjacent samples to generate a stable band set.

5. The geographic information acquisition method based on remote sensing imagery according to claim 4, characterized in that, The threshold for illumination change is determined by discretizing continuously collected brightness samples in illumination change scenarios and statistically analyzing the fluctuation range of the discrete samples.

6. The geographic information acquisition method based on remote sensing imagery according to claim 1, characterized in that, The specific steps of S5 are as follows: S501: Based on the spectral similarity value, perform numerical normalization on the similarity parameter group of land use type, and calculate the ratio of each normalized similarity parameter group to the spectral response benchmark value of the same source pixel. Then, perform weighted aggregation on the ratio results to generate a likelihood probability matrix. S502: Call the conditional probability parameter group in the likelihood probability matrix, and perform maximum likelihood determination on the multi-type conditional probability parameter group of each pixel based on the Bayesian classifier. Record the maximum numerical index corresponding to the determination result as the land use type number and perform rasterization encoding to obtain the pixel type index set. S503: Based on the land use type sequence number in the pixel type index set, perform spatial mapping on the coordinate group of pixels with the same sequence number, perform coordinate reprojection, and perform regional stitching on the reprojected pixel partitions to obtain a geographic information distribution map.

7. The geographic information acquisition method based on remote sensing imagery according to claim 6, characterized in that, The spectral response benchmark value is based on standard sample pixels verified by field testing of land cover categories, the statistical results of stable reflectance of land use types in multispectral bands, and the calibrated spectral response value is obtained by interval averaging of the spectral reflectance of sample pixels of the same category after removing spectral noise interference.

8. A geographic information acquisition system based on remote sensing imagery, characterized in that, The system is used to implement the geographic information acquisition method based on remote sensing imagery as described in any one of claims 1-7, and the system comprises: The spectral positioning module acquires multi-band remote sensing image data of the target area through a narrowband spectral sensor, extracts the position coordinates of the band reflectance peak based on the peak detection algorithm, calculates the difference with the coordinates of the reference time node, and generates band peak displacement data. The band stabilization module calls the peak displacement data of the band, performs time-series stability analysis using a Kalman filter, calculates the variance of the band displacement, compares it with the illumination change threshold, filters out stable bands below the threshold, and generates a set of stable bands. The type feature module performs band filtering on remote sensing image data based on the stable band set, obtains spectral reflectance data of sample pixels of different land use types, calculates the average reflectance value for multiple types of samples, and constructs a standard spectral feature library. The similarity calculation module calls the standard spectral feature library to calculate the difference between the spectral reflectance value of the pixel in the area to be classified and the standard reflectance value of the land use type in a band-by-band manner, counts the number of difference sign matches, and calculates the spectral similarity value. The classification output module calls the spectral similarity value, performs maximum likelihood estimation based on a Bayesian classifier, selects the land use type with the highest similarity value as the pixel classification label, and generates a geographic information distribution map.

Images