Method and system for improving geometric precision correction of historical geoscience map

Abstract

The invention discloses a geometric precision correction improving method and system for a historical geoscience map. The method comprises the following steps: acquiring a high-resolution panchromatic image and a middle-resolution multispectral image corresponding to a historical geoscience map, performing geometric correction on the high-resolution panchromatic image by taking a standard remote sensing image as a reference, and performing geometric correction on the middle-resolution multispectral image by taking the corrected high-resolution panchromatic image as a reference image; and finally, performing histogram matching and fusion on the single-band intensity images of the high-resolution panchromatic image and the middle-resolution multispectral image to generate a historical earth surface reference color image or further performing three-dimensional analysis for performing geometric correction and pattern spot content correction on a historical geoscience map. According to the method, by reconstructing the high-precision historical reference scene, the problem that high-precision geometric correction and pattern spot content checking and correction are difficult to carry out on the historical geoscience map due to long years and large earth surface change is effectively solved, the high-quality historical geoscience map is obtained, and spatial-temporal change analysis and application can be carried out on the historical geoscience map and the current geoscience map.

Description

Technical Field

[0001] The present invention belongs to the technical fields of cartography and remote sensing, and particularly relates to a method and system for geometric precision correction and improvement of historical geoscience maps. Background Art

[0002] Historical geoscience maps refer to some map documents dating back a long time (such as from the 1970s to the 1980s of the 20th century), and their types include soil maps, land use maps, vegetation maps, geological maps, geomorphological maps, etc. Limited by the data acquisition technology and mapping methods and techniques at that time, there are often a series of inevitable errors in the spatial positions and content information of these geoscience maps. In the fields of agriculture and forestry, environment, ecology, land and resources, water conservancy, etc., there is a wide range of technical needs for the correction, modification, improvement, and perfection of historical geoscience maps.

[0003] Taking the soil map as an example, large numbers of countries around the world have carried out soil surveys, and the most important result is precisely the soil map. The soil map reflects the soil types and spatial distributions in various places, and is the most direct manifestation of the quantity and quality of soil resources. Taking the second national soil survey carried out in China from 1979 to 1986 as an example, most townships completed 1:10,000 soil maps, and counties (districts) completed 1:50,000 soil maps; and the third national soil survey was carried out in 2021. The soil map required by the third soil survey is to be improved and updated on the basis of the soil map obtained from the second soil survey (referred to as the second general survey soil map), and at the same time, analyze and master the changes in the quantity and quality of soil resources in the past nearly 40 years. Therefore, to achieve such soil change comparison, high-quality second general survey soil maps must be used as the basis. However, the general procedure for making the second general survey soil maps in each county, city, and district across the country is as follows: using topographic maps (using aerial photos in places with better conditions) as the base map for on-site surveys, excavating soil profiles, combining laboratory analysis and test data, determining soil types, and demarcating soil type boundaries. If aerial photos are used as the base map, stereoscopic observation of adjacent image pairs of aerial photos can be carried out under a stereoscope, which helps to better judge and determine soil boundaries. The method using aerial photos as the base map has a vivid image of the observation study area and greatly speeds up the soil survey speed and quality, and has obvious improvement compared with the method using topographic maps as the base map. However, the method using aerial photos as the base map also brings some new problems, mainly including:

[0004] 1) Aerial photo distortion: Aerial photos are central projections, with projection errors, tilt errors, large edge distortions, and different scales in different parts. In areas with large topographic undulations, the elevation differences lead to scale changes, further exacerbating the distortion. In some areas with large elevation differences, the scales of different elevation areas are inconsistent, resulting in distortion;

[0005] 2) Drafting error: The area covered by a single aerial photograph is limited. In general counties in China, tens of thousands of aerial photographs are required for complete coverage. The distortions of different photographs vary, and the distortions between adjacent photographs are also different;

[0006] 3) Map sheet edge-matching problem: The distortions between different aerial photographs and different map sheets, as well as the expansion and contraction of the map medium (such as paper), result in misalignment and discontinuity of map patches at the map sheet joints.

[0007] Therefore, due to the limitations of the then aerial photograph soil survey technical mapping and the overall soil survey mapping technical conditions related thereto, there are obviously multiple spatial position errors in the accuracy of the positions of soil boundaries on the finally completed soil map.

[0008] Over time, especially in recent years, with the drastic land use / cover changes brought about by the rapid economic development, it has become extremely difficult to directly use modern high-precision images or topographic maps to register and correct the second national soil survey map, mainly reflected in:

[0009] 1) It is extremely difficult or even infeasible to select corresponding points: This is the fundamental bottleneck in applying the current geometric correction method to the second national soil survey map. The second national soil survey map is a thematic map, mainly presenting soil categories, supplemented by a small amount of information such as rivers, roads, villages, and urban boundaries. It has been nearly 40 years since the completion of the second national soil survey map. During this period, urban expansion has been drastic, road and water systems have changed significantly, farmland has been greatly regularized, and the surface has changed tremendously. It has become extremely difficult to directly find a sufficient number of high-precision, spatially corresponding, and evenly distributed stable corresponding points (i.e., feature points that exist in both periods and whose morphological positions have not changed) between the original second national soil survey map and modern high-precision images / topographic maps, and it cannot be achieved in most areas. Without reliable corresponding points, the traditional point-based geometric registration and correction methods cannot be effectively implemented.

[0010] 2) The accuracy of the existing correction methods for the second national soil survey map is seriously insufficient: Even if attempts are made to correct using the existing technical specifications, it is difficult to improve the accuracy due to the lack of high-quality corresponding points. By randomly checking the second national soil survey maps of some counties and cities, it is found that the average geometric error of soil boundaries reaches 50 meters, and in some places, the error is as high as 200 meters. Such accuracy far from meets the numerous current and future requirements.

[0011] 3) Limitations in the application of other data sources: If only contemporary Landsat / SPOT images are used, although these images are synchronous with the second national soil survey period (around 1979 - 1986), their spatial resolution (20 - 80 meters) is insufficient to support precise positioning and calibration at the meter level, and it is difficult to accurately identify fine feature points for calibration. If only modern high-resolution images are used, there are problems such as the spatio-temporal mismatch of features as mentioned above, and reliable homologous points cannot be found for the precise registration of historical maps; if only original aerial photographs are used, since it is very difficult to obtain original aerial photographs now, and even if the aerial photographs used in that year can be obtained, the processing process is cumbersome, involving a series of steps such as scanning, distortion correction, mosaicking, and transfer mapping. The cost is high, the distortion is complex, and errors will still be introduced when finally transferred to topographic maps. Especially, the coverage area of a single aerial photograph is small, and it takes a huge amount of time and extremely high cost to process a large area.

[0012] Of course, in addition to the second national soil survey maps, historical geoscience maps such as land use maps, vegetation maps, geological maps, and geomorphological maps in the same period (1970s - 1980s) also have similar technical problems that are difficult to perform high-precision geometric calibration and improve the content of map patches due to the long time, huge changes in the surface, and lack of homologous points. Summary of the Invention

[0013] The purpose of the present invention is to solve the technical problems that historical geoscience maps in the prior art are difficult to perform high-precision geometric calibration and improve the content of map patches, and to provide a geometric precise calibration and improvement method and system for historical geoscience maps. The present invention fuses multi-source remote sensing data in historical periods, constructs a high-precision historical surface reference color image as an intermediate transition image information layer, and based on this, establishes the association between the current high-precision standard remote sensing image and the early historical geoscience map, so as to realize the geometric precise calibration and improvement of the content of map patches of the second national soil survey maps.

[0014] The specific technical solutions adopted by the present invention are as follows:

[0015] In the first aspect, the present invention provides a geometric precise calibration and improvement method for historical geoscience maps, which includes:

[0016] S1. According to the spatio-temporal range corresponding to the historical geoscience map to be calibrated, obtain the corresponding high-resolution panchromatic image and medium-resolution multispectral image, and perform preprocessing on them respectively;

[0017] S2. Using the standard geospatial information as the first reference image, select control points on the preprocessed high-resolution panchromatic image and the first reference image, and perform geometric calibration on the high-resolution panchromatic image;

[0018] S3. Use the geometrically corrected high - resolution panchromatic image as the second reference image, select control points on the pre - processed medium - resolution multispectral image and the second reference image, and perform geometric correction on the medium - resolution multispectral image;

[0019] S4. According to the spectral range of the high - resolution panchromatic image, extract the multi - band image within the spectral range from the geometrically corrected medium - resolution multispectral image and convert it into a single - band intensity image. Then, perform histogram matching between the single - band intensity image and the geometrically corrected high - resolution panchromatic image, and generate a historical surface reference color image through image fusion;

[0020] S5. Use the historical surface reference color image as a two - dimensional reference base map to perform geometric correction and patch content correction on the historical geoscience map.

[0021] As a preference of the above - mentioned first aspect, after obtaining the historical surface reference color image, use the historical surface reference color image as the surface color source, use the digital elevation model data within the time - space range as the elevation source, generate a three - dimensional reference base map through three - dimensional visualization technology, and perform geometric correction and patch content correction on the historical geoscience map.

[0022] As a preference of the above - mentioned first aspect, the historical geoscience map is a soil map, a land use map, a vegetation map, a geological map or a geomorphological map.

[0023] As a preference of the above - mentioned first aspect, the historical geoscience map is a second national soil survey soil map, the high - resolution panchromatic image is a Conora satellite image, and the medium - resolution multispectral image is a Landsat MSS image, a Landsat TM image or a SPOT HRV image.

[0024] As a preference of the above - mentioned first aspect, the georeference data uses a topographic map with a scale of not less than 1:10,000 or an orthoimage with a spatial resolution of not less than 1 meter.

[0025] As a preference of the above - mentioned first aspect, the single - band intensity image is obtained by weighted averaging of the multi - band images within the spectral range, or uses the first principal component of the multi - band images within the spectral range, or uses the first component of the Gram - Schmidt transform of the multi - band images within the spectral range.

[0026] As a preference of the above - mentioned first aspect, when performing histogram matching, it is necessary to upsample the single - band intensity image to the same resolution as the geometrically corrected high - resolution panchromatic image, then perform histogram matching between the two, make the cumulative histogram of the high - resolution panchromatic image approximate the upsampled single - band intensity image, and finally perform image fusion on the two types of images after histogram matching through the Gram - Schmidt algorithm to generate a historical surface reference color image.

[0027] Preferably, in the first aspect above, when performing geometric correction and patch content correction on the historical geological map based on the two-dimensional reference base map and the three-dimensional reference base map, it is necessary to select the same-named ground feature points in the same period on the reference base map and the historical geological map as control points, perform geometric correction on the historical geological map, and then compare with the reference base map to correct and improve the soil patch boundary and attribute content on the geometrically corrected historical geological map. Finally, after topological relationship inspection and correction, the final map is output.

[0028] In a second aspect, the present invention provides a computer-readable storage medium, on which a computer program is stored. When the computer program is executed by a processor, it implements the method for geometric precision correction and improvement of historical geological maps as described in any one of the above first aspect solutions.

[0029] In a third aspect, the present invention provides a computer electronic device, which includes a memory and a processor;

[0030] The memory is used to store a computer program;

[0031] The processor is used to implement the method for geometric precision correction and improvement of historical geological maps as described in any one of the above first aspect solutions when executing the computer program.

[0032] The present invention has the following beneficial effects compared with the prior art:

[0033] 1) Solve the problem of historical map correction: The present invention effectively solves the long-existing key technical bottleneck that it is difficult to perform high-precision geometric correction on historical geological maps such as the second general soil map due to the lack of stable same-named points across periods caused by the long time and large surface changes by fusing high and low-resolution remote sensing images of a specific historical period to reconstruct a high-precision historical reference scene.

[0034] 2) Significantly improve geometric accuracy: Compared with the existing method whose correction accuracy is only dozens of meters, or even worse (such as 10 - 200 meters), the method of the present invention can improve the geometric correction accuracy of historical geological maps such as the second general soil map to 2 - 5 meters. This order-of-magnitude improvement in accuracy lays a solid foundation for subsequent refined applications.

[0035] 3) Greatly improve the accuracy and detail of patch content interpretation: The present invention not only solves the positioning problem, but also can significantly improve the content quality of historical geological maps such as soil maps. Using high-precision and information-rich historical reference color images (and three-dimensional views), the original patch boundaries can be accurately adjusted, errors corrected, and omissions supplemented.

[0036] 4) High-efficiency processing: Compared with processing original aerial photos, this method uses satellite images with a large single-scene coverage area (such as a single Corona image covering approximately 700 square kilometers) for processing. Combining modern geographic information system methods and cloud computing technologies (such as the GEE platform), it has high processing efficiency and is very suitable for large-area applications.

[0037] 5) Provide unique historical reference data: The high-resolution (meter-level and sub-meter-level), multi-spectral historical surface reference color images generated by this invention are a highly valuable data product. It restores the true color landscape of the surface corresponding to the historical period of historical geological maps, filling the gap in high-precision color surface data for that period. This specific data combination of Corona + Landsat / SPOT is the key to achieving this goal: Corona provides the necessary meter-level and sub-meter-level geometric accuracy framework and micro-topographic details, and Landsat / SPOT provides the key color and multi-spectral information to enhance the interpretability of ground objects (especially soil-related ground objects). The human eye can generally distinguish differences in 10 gray levels, while through this invention, it is possible to distinguish tens of thousands of different color-level ground objects. Therefore, this combination is not a simple superposition but an organic integration, forming a high-quality, information-rich historical ground object reference.

[0038] 6) Three-dimensional stereoscopic display: By combining color images with meter-level and sub-meter-level geometric accuracy frameworks and micro-topographic details with DEM (Digital Terrain Elevation Model data), this invention can perform three-dimensional stereoscopic display, which is very helpful for the interpretation and judgment of specific soil types and plays an important role in the verification, correction, and adjustment of the map patch types of the second national soil survey.

[0039] 7) Wide application range: This method is not only applicable to the second national soil survey maps but also applicable to correcting, repairing, and improving other historical geological maps of the same period, such as land use maps, vegetation maps, geological maps, geomorphological maps, etc., or directly used for making high-precision thematic maps of that period. BRIEF DESCRIPTION OF THE DRAWINGS

[0040] Figure 1 Schematic diagram of the steps for the geometric precision correction improvement method of historical geological maps;

[0041] Figure 2 Schematic diagram of the structure of a computer electronic device;

[0042] Figure 3 Two-dimensional reference base map and three-dimensional reference base map generated based on the original Conora satellite image and Landsat TM multi-spectral image in the embodiment of the present invention;

[0043] Figure 4The effect of overlaying the second soil survey map on the two-dimensional reference base map in mountainous and plain areas in the embodiments of the present invention;

[0044] Figure 5 Comparison of the second soil survey map before and after correction in the exemplary area in the embodiments of the present invention. Detailed implementation manners

[0045] To make the above objects, features, and advantages of the present invention more obvious and understandable, the following detailed description of the specific implementation manners of the present invention will be given with reference to the accompanying drawings. Many specific details are set forth in the following description to facilitate a thorough understanding of the present invention. However, the present invention can be implemented in many other ways different from those described herein, and those skilled in the art can make similar improvements without departing from the connotation of the present invention. Therefore, the present invention is not limited by the specific embodiments disclosed below. The technical features in various embodiments of the present invention can be combined correspondingly without conflict.

[0046] In the description of the present invention, it should be understood that the terms "first" and "second" are only used for the purpose of distinguishing descriptions, and cannot be understood as indicating or implying relative importance or implicitly specifying the quantity of the indicated technical features. Thus, the features defined with "first" and "second" may explicitly or implicitly include at least one of such features.

[0047] The present invention provides a method for geometric precision correction and improvement of historical geological maps. The method is based on creating a historical surface reference color image by fusing multi-source remote sensing data of a specific historical period as an intermediate transition image information layer, so as to perform precise correction of historical geological maps and high-quality correction and improvement of the content of map patches.

[0048] It should be noted that the historical geological maps applicable to the method of the present invention can be early soil maps, land use maps, vegetation maps, geological maps or geomorphic maps, especially the geological maps collected around the 1980s and earlier. For example, the soil map obtained from the second soil census in China (i.e., the second soil survey map). Due to the limitations of the technical conditions at that time, these geological maps often have relatively large errors.

[0049] In a preferred embodiment of the present invention, as Figure 1 shown, the specific steps of the above method for geometric precision correction and improvement of historical geological maps include:

[0050] S1. According to the spatio-temporal range corresponding to the historical geological map to be corrected, obtain the corresponding high-resolution panchromatic image and medium-resolution multi-spectral image, and perform preprocessing on them respectively.

[0051] It should be noted that the "high resolution" and "medium resolution" in the present invention refer to the relative levels of spatial resolution, and do not limit the absolute spatial resolution value. In theory, the higher the spatial resolution of the high-resolution panchromatic image, the better, and the higher the spatial resolution of the medium-resolution multispectral image, the better. However, limited by the remote sensing images actually existing in the historical period of the historical geoscience maps, it is impossible to obtain panchromatic images and multispectral images with any spatial resolution. Therefore, the highest-resolution remote sensing images can be selected as much as possible within the selectable range.

[0052] It should be noted that for the preprocessing of high-resolution panchromatic images and medium-resolution multispectral images here, the data characteristics of the remote sensing images themselves can be considered and implemented with reference to the conventional preprocessing steps of remote sensing images. General preprocessing operations include contrast stretching, histogram equalization, radiometric calibration, atmospheric correction, cloud removal, stripe noise removal, and bad line processing. Different remote sensing images can be combined and used according to different preprocessing operations as needed.

[0053] It should be noted that for the high-resolution panchromatic image and the medium-resolution multispectral image in the present invention, their corresponding spatio-temporal ranges should be as close as possible to the spatio-temporal range corresponding to the historical geoscience map to be corrected. In theory, it is optimal if the spatio-temporal ranges of the three are exactly the same. However, due to the limitation of remote sensing data resources, it is not possible to fully ensure that the corresponding remote sensing images can be obtained at any time, and there may also be missing image data at some locations spatially. However, since the distribution of surface features (soil distribution, land use type distribution, geomorphic distribution, etc.) does not change in real time, even if there is a certain deviation in the time range, it will not cause a great impact. Therefore, the time range of the high-resolution panchromatic image and the medium-resolution multispectral image can be as close as possible to the time range corresponding to the historical geoscience map to be corrected, and it does not need to be exactly the same. However, the spatial ranges of the high-resolution panchromatic image and the medium-resolution multispectral image should be as consistent as possible with the spatial range corresponding to the historical geoscience map to be corrected. Only when the remote sensing images are incomplete can local images be considered for correction.

[0054] S2. Taking the standard geospatial information as the first reference image, selecting control points on the preprocessed high-resolution panchromatic image and the first reference image, and performing geometric correction on the high-resolution panchromatic image.

[0055] It should be noted that the above-mentioned georeference data is used as a reference image for geometric correction of high-resolution panchromatic images. Therefore, its own accuracy needs to meet the requirements of correction. Moreover, it should be particularly noted that for the aforementioned high-resolution panchromatic images and medium-resolution multispectral images, their corresponding spatio-temporal ranges should be as close as possible to the spatio-temporal range corresponding to the historical geoscience maps to be corrected. However, the above-mentioned georeference data is not restricted by this. The latest accurate data can be used, and theoretically, the higher the accuracy, the better. In the embodiments of the present invention, the above-mentioned standard geospatial information can be a topographic map with a scale of not less than 1:10,000, or a remotely sensed orthoimage with a spatial resolution of not less than 1 meter, such as the ESRI global image base map with a resolution of not less than 1 meter without offset.

[0056] It should be noted that ground control points (GCPs) refer to the homologous feature points (also known as homologous points) collected in pairs on two images. In this step, homologous feature points need to be selected on the preprocessed medium-resolution multispectral image and the second reference image. The selected control points should be evenly distributed within the entire image range, and preference should be given to selecting feature points with clear outlines, stable positions, and easy to accurately identify, such as road intersections, bridges, corner points of large independent buildings, river turning points, characteristic points of reservoir dams, etc., which have existed since the era of the historical geoscience maps and have not changed much so far (relative to the period of the reference image). Different densities of control points need to be selected according to the terrain complexity. For example, in the gentle slope plain area, the recommended density of control points is not less than 0.5 - 1 per square kilometer; while in the hilly mountainous area with large terrain undulations (for example, areas with a height difference exceeding 50 meters and a slope greater than 15 degrees), the density should be increased to 1.5 - 3 per square kilometer. Therefore, the specific control points selected can be determined according to the actual image situation, and no restrictions are imposed on this.

[0057] In addition, in the present invention, geometric correction of the image based on the selected control points can be achieved by using existing geometric correction algorithms, such as the cubic spline interpolation function (Spline) algorithm, the rubber sheet transformation (Rubbersheet) algorithm, etc. However, different images to be corrected have different geometric offset characteristics. Therefore, it is necessary to select a suitable geometric correction algorithm according to the data situation of the image itself. If the high-resolution panchromatic image uses Corona image, it is preferably geometrically corrected using the cubic spline interpolation function (Spline) algorithm and registered to the first reference image (such as the aforementioned ESRI global image base map).

[0058] S3. Using the geometrically corrected high-resolution panchromatic image as the second reference image, select control points on the preprocessed medium-resolution multispectral image and the second reference image, and perform geometric correction on the medium-resolution multispectral image.

[0059] Similarly, in step S3, control points are also required to collect homologous ground feature points by pairing on the preprocessed medium-resolution multispectral image and the second reference image. The specific selection principle can also refer to the principle described above. The geometric correction algorithm used in this step can also be selected according to the data conditions of the medium-resolution multispectral image itself. If the medium-resolution multispectral image uses Landsat image or SPOT image, it is also recommended to use the cubic spline interpolation function (Spline) algorithm for geometric correction and register it to the geometrically corrected high-resolution panchromatic image.

[0060] S4. According to the spectral range of the high-resolution panchromatic image, extract the multi-band image within the spectral range from the geometrically corrected medium-resolution multispectral image and convert it into a single-band intensity image. Then, perform histogram matching on the single-band intensity image and the geometrically corrected high-resolution panchromatic image, and generate a historical surface reference color image through image fusion.

[0061] In the embodiment of the present invention, the above single-band intensity image is obtained by weighted average of the multi-band images within the spectral range of the high-resolution panchromatic image (the weights can also be set to be the same, that is, arithmetic average is performed), or the first principal component (PC1) of the above multi-band images is used, or the first component (GS1) of the Gram-Schmidt (GS) transform of the above multi-band images is used.

[0062] In the embodiment of the present invention, when performing histogram matching on the single-band intensity image and the geometrically corrected high-resolution panchromatic image, the specific method can be implemented as follows: first, upsample the single-band intensity image (bilinear interpolation, cubic convolution interpolation, etc. can be used) to the same resolution as the geometrically corrected high-resolution panchromatic image, then perform histogram matching on the two to make the cumulative histogram of the high-resolution panchromatic image approximate the upsampled single-band intensity image, and finally perform image fusion on the histogram-matched single-band intensity image and the high-resolution panchromatic image through the Gram-Schmidt (GS) algorithm to generate a historical surface reference color image.

[0063] There are various image fusions in the present invention. The above Gram-Schmidt (GS) algorithm is particularly suitable for historical soil maps, such as the second general soil survey map.

[0064] S5. Use the historical surface reference color image as a two-dimensional reference base map to perform geometric correction and patch content correction on the above historical geoscience maps.

[0065] It should be noted that the geometric correction and patch content correction of the above historical geological maps are achieved by interpreting the information shown on the reference base map. Since the historical surface reference color image restores the true color landscape of the surface during that period using the high-resolution panchromatic image and medium-resolution multispectral image of the same period, it can fill the gap in high-precision color surface data for that period. Therefore, through visual interpretation or automatic interpretation using algorithms or neural network techniques, landforms and ground objects can be accurately extracted from the historical surface reference color image. This historical surface reference color image not only ensures a geometric positioning accuracy at the meter level but also provides rich ground object information (color, texture). In actual operation, the historical geological map can be overlaid on this base map of the historical surface reference color image, and then the information in the historical geological map can be visually verified and corrected with reference to the historical surface reference color image, which greatly facilitates and improves the efficiency and quality of accurately adjusting the patch boundary, supplementing missing patches, and accurately verifying and correcting patch attributes.

[0066] In addition, the historical surface reference color image obtained in the above step S4 is two-dimensional. It can be used as a two-dimensional reference base map for geometric correction and patch content correction of historical geological maps, but it is not very intuitive. Therefore, in other embodiments of the present invention, after obtaining the historical surface reference color image, the historical surface reference color image can be used as the source of surface color, and the digital elevation model (DEM) data within the spatio-temporal range corresponding to the historical geological map can be used as the elevation source. Through three-dimensional visualization technology, a three-dimensional reference base map can be generated for geometric correction and patch content correction of historical geological maps. Theoretically, the higher the accuracy of the DEM data, the better. Insufficient DEM accuracy may affect the recognition of micro-topographic features and the accuracy of boundary delineation based on terrain. In the embodiments of the present invention, it is recommended to use DEM data with a spatial resolution better than 30 meters and relatively high elevation accuracy, such as ALOS PALSAR 12.5m DEM, which can be obtained from public channels. Of course, if there is DEM data with higher resolution, the effect will be better.

[0067] Thus, the present invention combines the historical surface reference color image with meter-level and sub-meter-level geometric accuracy framework and micro-topographic details with DEM data, enabling three-dimensional stereoscopic display of the surface color landform. Overlaying the historical geological map on this three-dimensional reference base map is more conducive to the interpretation and reading of specific ground objects in most scenarios compared to the two-dimensional reference base map, and plays an important role in the verification, correction, and adjustment of the patch types of historical geological maps.

[0068] It should be noted that the above two-dimensional reference base map and three-dimensional reference base map in the present invention can be selected alternatively for geometric correction and patch content correction, or can be used alternately in different correction steps, and the specific usage is not limited. When performing geometric correction and patch content correction on the historical geological map based on the two-dimensional reference base map and three-dimensional reference base map, it is necessary to select the same-named ground feature points in the same period on the reference base map and the historical geological map as control points, perform geometric correction on the historical geological map, and then compare with the reference base map to correct and improve the soil patch boundary and attribute content on the geometrically corrected historical geological map (which can be done manually or through an automatic algorithm or a neural network model). Finally, after topological relationship inspection and correction, the final map is output. The content of the topological relationship inspection includes at least: no overlap (Overlap), no gap (Gap), closed boundary between patches, and integrity of patch attributes. If errors are found in the topological relationship inspection, corresponding repairs are required.

[0069] To better understand the specific implementation process and principle of the present invention, the soil map obtained from the second national soil census (i.e., the second general survey soil map) is used below as an exemplary historical geological map to demonstrate the specific practices and technical effects of the above method for improving the geometric precision correction of historical geological maps. Corresponding to the second general survey soil map, the high-resolution panchromatic image used is the Conora satellite image, and the medium-resolution multispectral images are the Landsat MSS image, Landsat TM image or SPOT HRV image. This method includes steps such as data acquisition and preprocessing, geometric correction of high-resolution panchromatic images, geometric calibration of multispectral images, generation of fused historical surface reference color images, three-dimensional visualization by fusing DEMs, geometric correction of the second general survey soil map, and correction and improvement of the content of the second general survey soil map. The specific practices of each step are described below respectively.

[0070] 1) Data acquisition and preprocessing: Obtain panchromatic images with high spatial resolution (0.6 - 2 meters) during the Second National Soil Census (around 1979 - 1986). In this embodiment, declassified Conora satellite images are used; obtain satellite images with medium spatial resolution (20 - 80 meters) and multi - spectral (3 - 6 bands) during the same period. In this embodiment, Landsat MSS / TM images (30 - 80 meters spatial resolution, 4 - 6 spectral bands) or SPOT HRV images (20 meters spatial resolution, 3 bands) are used. Preprocess the acquired various image data. For Corona images, quality inspection of the film - scanned images is required, and contrast enhancement and other processing are carried out if necessary; for Landsat and SPOT images, to reduce the influence of atmosphere and illumination, radiometric calibration and atmospheric correction are required, and they are converted into apparent reflectance or surface reflectance (standard algorithms provided by, for example, FLAASH, QUAC, 6S models or corresponding data processing systems). At the same time, identify and mask clouds, cloud shadows, strip noises or bad lines in the images to avoid affecting subsequent processing.

[0071] 2) Geometric correction of high - resolution panchromatic images: Use existing high - precision georeference data (in this embodiment, topographic maps with a scale larger than 1:10,000 or high - precision orthoimages with a spatial resolution better than 1 meter) as the first reference image to perform precise geometric correction on the acquired high - resolution Corona images.

[0072] 2.1) Control point selection: Control points should be evenly distributed within the entire image range, and preferably select feature points with clear outlines, stable positions, and easy to accurately identify, such as road intersections, bridges, corner points of large independent buildings, river turning points, characteristic points of reservoir dams, etc., which existed in the early 1980s and have changed little so far (relative to the reference data period). Different densities of control points need to be selected according to the terrain complexity. In gently sloping plain areas, the recommended density of control points is not less than 0.5 - 1 per square kilometer; while in hilly mountainous areas with large terrain undulations (for example, areas with a height difference exceeding 50 meters and a slope greater than 15 degrees), the density should be increased to 1.5 - 3 per square kilometer. For a single Corona image (covering approximately 700 square kilometers), it is recommended to select at least 200 - 300 high - quality control points.

[0073] 2.2) Image correction: Considering the possible local non - linear distortion in Corona images, in this embodiment, the cubic spline interpolation function (Spline) algorithm capable of handling local deformation is used as the correction algorithm to perform geometric correction on the pre - processed Corona images, and register them onto the first reference image. The correction target set in this embodiment is to control the root mean square error (RMSE) within 2 - 5 pixels (corresponding to about 2 - 5 meters in the Corona image). After geometric correction, a quasi - ortho Corona image with high geometric accuracy is obtained.

[0074] 3) Geometric calibration of multi - spectral images: Using the quasi - ortho Corona image corrected in the previous step as the second reference image, perform geometric calibration (registration) on the medium - resolution multi - spectral images (Landsat images, SPOT images) of the same period. Referring to the high - resolution panchromatic image geometric correction steps in the previous step, select the homologous ground control points clearly visible on both the Corona image and the medium - resolution multi - spectral images, and use the same correction algorithm (Spline algorithm) for geometric calibration, so that the medium - resolution multi - spectral images are registered onto the second reference image, that is, the Landsat images, SPOT images and the quasi - ortho Corona image are precisely aligned in space.

[0075] 4) Generation of historical surface reference color image fusion: Before fusing the high - resolution panchromatic image (quasi - ortho Corona image) geometrically corrected in step 2) and the medium - resolution multi - spectral images (Landsat / SPOT) geometrically corrected in step 3), it is necessary to first perform radiometric normalization / matching pre - processing to make the radiometric characteristics such as brightness and contrast of the PAN image more consistent with the intensity information of the MS image, so as to better retain spectral information and reduce color distortion in subsequent fusion. The radiometric normalization method adopted in this embodiment is histogram matching, and its specific steps may include:

[0076] (a) Generate the target intensity image: Based on the medium - resolution multi - spectral images (Landsat / SPOT) geometrically corrected in step 3), generate a single - band target intensity image (Target Intensity Image), that is, a single - band intensity image. This single - band intensity image can be used to represent the overall brightness of the image. The methods for generating a single - band intensity image from multi - spectral images include any one of the following a1) - a3):

[0077] a1) According to the spectral range of the high-resolution panchromatic image itself, extract the multispectral bands within this spectral range from the geometrically corrected medium-resolution multispectral image (taking Landsat TM image as an example, the blue, green, red, and near-infrared bands need to be extracted), and perform weighted average or direct arithmetic average on the multispectral bands to calculate a synthetic luminance image as the single-band intensity image;

[0078] a2) Extract the first principal component (PC1) of the principal component analysis from the multispectral bands obtained in a1) above as the single-band intensity image;

[0079] a3) Extract the first component (GS1) of the Gram-Schmidt transform from the multispectral bands obtained in a1) above as the single-band intensity image.

[0080] It should be noted that for different medium-resolution multispectral images, such as Landsat and SPOT, the single-band intensity images need to be extracted respectively. The single-band intensity image extracted corresponding to any medium-resolution multispectral image can be used for the subsequent histogram matching and image fusion of the present invention.

[0081] (b) Perform histogram matching:

[0082] The resolution of the single-band intensity image generated in step (a) is lower than that of the high-resolution panchromatic image (quasi-orthorectified Corona image). Therefore, it is necessary to upsample it to the same spatial resolution as the high-resolution panchromatic image through resampling methods (such as bilinear interpolation, cubic convolution interpolation).

[0083] Then, taking the geometrically corrected high-resolution panchromatic image as the source image and the upsampled single-band intensity image as the target (reference) image, perform the histogram matching algorithm to adjust the pixel values of the high-resolution panchromatic image so that its cumulative histogram approximates the cumulative histogram of the target image, and obtain the radiometrically matched panchromatic image (PAN_matched).

[0084] After completing the radiometric normalization / matching preprocessing, fuse the panchromatic image PAN_matched (high resolution) and the geometrically corrected medium-resolution multispectral image (multispectral bands, maintaining its original resolution) through the image fusion algorithm to generate a high-resolution historical surface reference color image with both high spatial resolution and multispectral information, reflecting the surface conditions during the second national land survey period.

[0085] The image fusion algorithm in the present invention preferably adopts the Gram - Schmidt (GS) transform (or its improved version) fusion algorithm. The reason for choosing the GS algorithm is that while effectively injecting the high - spatial - resolution details of the Corona image, it can better maintain the original spectral characteristics of the Landsat / SPOT image and reduce color distortion during the fusion process. This is crucial for the subsequent interpretation of soil - related features (such as vegetation, land use types) that rely on spectral information. In contrast, some other fusion methods (such as HSV) may sometimes cause more obvious color aberrations, so they are not the preferred methods of the present invention.

[0086] The method for determining the fusion parameters of the image fusion algorithm is as follows: The specific parameters of the Gram - Schmidt fusion (such as the weight coefficient) should be determined according to the specific implementation of the software used (such as ArcGIS Pro, ENVI, or the Google Earth Engine platform). Usually, an optimal selection can be made by combining qualitative visual interpretation of the fusion results (checking spatial clarity, color naturalness, and the presence of obvious artifacts) and quantitative evaluation indicators (see the next step of quality control).

[0087] After completing the image fusion of the single - band intensity image and the geometric - corrected high - resolution panchromatic image through histogram matching and the Gram - Schmidt algorithm, it is necessary to perform image quality control on the fused historical surface reference color image. The specific approach is as follows: Evaluate the quality of the fused historical surface reference color image, and the evaluation can be combined with qualitative visual inspection (spatial texture clarity, spectral fidelity, color naturalness, and the presence of fusion defects) and quantitative indicators. Quantitative indicators such as the relative global dimensional synthesis error (ERGAS), spectral angle mapping (SAM), correlation coefficient (CC), root - mean - square error (RMSE), etc. can be selected. Evaluate the spectral distortion degree by comparing the fused image with the original multi - spectral image, and evaluate the spatial detail injection effect by comparing it with the original panchromatic image. Acceptable index thresholds should be set to ensure that the quality of the fused image meets the requirements of subsequent applications.

[0088] 5) Integrate with DEM for 3D visualization: Integrate the historical surface reference color image generated in the previous step with the Digital Elevation Model (DEM) data to create a 3D visualization scene of the target area in GIS software (such as ArcGIS Pro). When creating the 3D visualization scene, the historical surface reference color image should be used as the surface color source, and the DEM data should be used as the elevation source. It should be noted that the spatial range corresponding to the DEM data should correspond to the spatial range of the historical surface reference color image itself, and there is no time limit. The latest DEM data can be used. Here, the spatial resolution and elevation accuracy of the DEM used for 3D visualization will affect the 3D visualization effect and the terrain-based interpretation accuracy. In this embodiment, it is recommended to use DEM data with a spatial resolution better than 30 meters (such as ALOS PALSAR 12.5m DEM or higher-precision DEM) and higher elevation accuracy. Insufficient DEM accuracy may affect the recognition of micro-topographic features and the accuracy of terrain-based boundary delineation. Thus, the 3D visualization scene generated by 3D visualization technology can be used as the 3D reference base map for subsequent geometric correction and modification of the second national soil survey soil map.

[0089] 6) Geometric correction of the second national soil survey soil map: Use the high-precision historical surface reference color image generated in step 4 (or combined with the 3D visualization scene in step 5) as an accurate reference base map to perform geometric correction on the original second national soil survey soil map. Select homologous ground control points on the second national soil survey soil map that can be clearly identified on the historical surface reference color image (note: at this time, find homologous points on the image of the same historical period, rather than spanning 40 years), and then perform geometric correction on the second national soil survey soil map. The geometric correction algorithm for the second national soil survey soil map preferably uses the Rubber Sheet algorithm, which can locally adjust vector data according to control points to adapt to the irregular deformation of the map. The correction target is to make the geometric accuracy of the corrected soil map reach the meter level (for example, 2 - 5 meters).

[0090] 7) Content correction and improvement of the second national soil survey soil map: Overlay the second national soil survey soil map after geometric correction in step 6 on the historical surface reference color image in step 4 (or in the 3D visualization scene in step 5), and let soil experts or trained technicians correct and improve the map patch boundaries and content. Of course, an automatic algorithm can also be designed to correct and improve the map patch boundaries and content.

[0091] The content correction and improvement of the second national soil survey soil map need to follow certain correction standards. In the embodiments of the present invention, the following boundary adjustment basis and content update standard requirements can be set:

[0092] Basis for boundary adjustment: In a two-dimensional or three-dimensional environment, the main bases for adjusting the boundaries of soil patches are as follows: 1) The boundaries of geomorphic units shown on the fused high-precision historical color images (such as ridge lines, valley lines, slope / slope aspect change lines, terrace edges); 2) The obvious demarcation lines of vegetation or land use types reflected by colors and textures on the images (the correlation with soil types needs to be judged in combination with soil knowledge); 3) The topographic position characteristics (such as sunny slope / shady slope, slope top / slope shoulder / slope middle / slope foot) shown in combination with the DEM. The goal is to make the corrected soil boundaries as consistent as possible with the restored surface landscape characteristics in the 1980s.

[0093] Content update: Check whether the soil type markings of the original patches are accurate against the historical reference images, and correct the errors; identify and delineate the soil patches that were omitted in the original map or could not be expressed due to scale limitations; merge or subdivide the original patches as needed.

[0094] 8) Topological check and final result output: Conduct a topological relationship check and correction on the corrected Second National Soil Census map.

[0095] Topological relationship rules: The content of the topological check should at least include no overlap, no gap, closed boundaries between patches, and the integrity of patch attributes. The topological rule verification tool of GIS software (such as ArcGIS Pro) can be used for checking, and the discovered errors can be repaired manually or automatically.

[0096] Quality evaluation of the final result: For the corrected result of the Second National Soil Census map, the quality evaluation should include: 1) Geometric accuracy: By selecting independent check points (ICPs), measuring the coordinates on the historical reference color image, and calculating the RMSE of the corrected patch boundary position, it should meet the design requirements (such as 2 - 5 meters); 2) Topological consistency: Confirm that there are no topological errors; 3) Content rationality: Soil experts conduct a sampling inspection and evaluation on the geographical / geomorphic compliance of the patch boundaries and the accuracy of the patch attribute markings by comparing the historical reference image and the three-dimensional view.

[0097] Final result output: Output the digital result of the Second National Soil Census map that meets the standards and has been corrected with high precision and content improvement.

[0098] It can be seen that the present invention reconstructs a high-precision historical surface reference color image by fusing high- and low-resolution remote sensing images of a specific historical period, effectively solving the long-existing key technical bottleneck that it is difficult to perform high-precision geometric correction on historical maps such as the second national soil survey map due to the lack of stable homologous points across different periods caused by the long time span and significant surface changes. This historical surface reference color image establishes the correlation between today's high-precision remote sensing images and the early second national soil survey map, enabling the use of the rich information provided by the historical surface reference color image to achieve geometric correction, content correction, and improvement of the second national soil survey map patches.

[0099] Of course, although the above embodiments are described with the second national soil survey map as a typical example, it is obvious that the method of the present invention is not only applicable to the second national soil survey map, but also applicable to correcting, repairing, and improving other historical geoscience maps of the same period, such as land use maps, vegetation maps, geological maps, geomorphic maps, etc., or directly used for making high-precision thematic maps of this period, and there is no limitation in this regard.

[0100] It should be noted that the method steps shown in S1~S5 above can essentially be implemented in the form of a computer program.

[0101] Thus, based on the same inventive concept, as Figure 2 shown, the present invention also provides a computer electronic device corresponding to a method for geometric precision correction and improvement of a historical geoscience map provided in the above embodiment, which includes a memory and a processor;

[0102] The memory is used to store a computer program;

[0103] The processor is used to implement the method for geometric precision correction and improvement of the historical geoscience map as described above when executing the computer program;

[0104] In addition, when the logical instructions in the above memory are implemented in the form of software functional units and sold or used as independent products, they can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, in essence, or the part that contributes to the prior art or a part of this technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions for causing a computer device (which can be a personal computer, a server, or a network device, etc.) to execute all or part of the steps of the methods described in various embodiments of the present invention.

[0105] Thus, based on the same inventive concept, the present invention provides a computer-readable storage medium corresponding to a method for geometric precision correction and improvement of a historical geoscience map. The storage medium stores a computer program, and when the computer program is executed by a processor, it can implement the method for geometric precision correction and improvement of the historical geoscience map as described above.

[0106] Accordingly, based on the same inventive concept, the present invention provides a computer program product including computer programs / instructions, which, when executed by a processor, can implement the method for improving the geometric precise correction of historical geological maps as described above.

[0107] Specifically, in the computer-readable storage media of the above three embodiments, the stored computer program, when executed by a processor, can execute the steps of S1 to S4 described above.

[0108] It can be understood that the above storage medium may include a random access memory (RAM), and may also include a non-volatile memory (NVM), such as at least one disk memory. At the same time, the storage medium may also be various media such as a USB flash drive, a mobile hard disk, a magnetic disk, or an optical disc that can store program codes.

[0109] It can be understood that the above processor may be a general-purpose processor, including a central processing unit (CPU), a network processor (NP), etc.; it may also be a digital signal processor (DSP), an application specific integrated circuit (ASIC), a field-programmable gate array (FPGA), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components.

[0110] In addition, it should be noted that those skilled in the art can clearly understand that for the convenience and simplicity of description, the specific working process of the above-described system can refer to the corresponding process in the foregoing method embodiments, and will not be elaborated herein. In each embodiment provided in the present application, the division of steps or modules in the system and method is only a logical function division, and there may be other division methods in actual implementation. For example, multiple modules or steps may be combined or integrated together, and one module or step may also be split.

[0111] To better understand the technical effects that can be achieved by the above method of the present invention, the following shows the specific implementation manner of the method for improving the geometric precise correction of historical geological maps in the above embodiments, as well as the correction and improvement effects, through an exemplary local spatial region.

[0112] Embodiment

[0113] In this embodiment, the historical geological map targeted is the second national soil survey map, the corresponding high-resolution panchromatic image is the Conora satellite image, and the medium-resolution multispectral image is the Landsat MSS / TM image or the SPOT HRV image. The specific steps for improving the geometric precise correction of the second national soil survey map are as follows:

[0114] Step 1: Data acquisition and preprocessing

[0115] Obtain the third batch of declassified Conora satellite images with high spatial resolution (0.6 - 2 meters) during the second national soil survey period (1979 - 1986), and obtain remote sensing satellite images with medium spatial resolution (20 - 80 meters) and containing multispectral bands (3 - 6 bands) during the same period and complete data preprocessing. The specific data and preprocessing methods are as follows:

[0116] Obtaining high-resolution panchromatic images: Obtain declassified Corona KH series images (such as KH-4B, KH-9, etc., with a spatial resolution of about 0.6 - 2 meters) from the USGS EarthExplorer website. Specifically, select Declassified Data -> Declass 3 in DataSets and download.

[0117] Obtaining medium-resolution multispectral images: Obtain Landsat MSS (about 80-meter resolution, 4 bands) or Landsat TM (30-meter resolution, 6 visible / infrared bands) images (Landsat Collection 2 Level-1) during the same period from the USGS EarthExplorer; or obtain SPOT-1 HRV images (20-meter multispectral resolution, 3 bands) from the SPOT official website. If there are multiple images during the same period, images with less cloud cover and good quality can be selected.

[0118] Preprocessing of Corona image data: Visually inspect the downloaded scanned images to confirm that there are no serious scratches or stains. Image processing software (such as Photoshop, GIMP, or built-in tools in GIS software) can be used for contrast stretching or histogram equalization to enhance ground object details.

[0119] Preprocessing of Landsat / SPOT image data: For Landsat / SPOT image data, professional remote sensing image processing software (such as ENVI, ERDAS IMAGINE, ArcGIS Pro, QGIS) or online platforms (such as GoogleEarth Engine - GEE) can be used for processing. The specific preprocessing includes radiometric calibration, atmospheric correction, and cloud and noise processing, which are as follows:

[0120] Radiometric calibration: Convert DN values ​​to radiance.

[0121] Atmospheric correction: Use the surface reflectance product algorithm provided by FLAASH, QUAC, 6S model or GEE to convert radiance or apparent reflectance into surface reflectance to eliminate the effects of atmospheric scattering and absorption.

[0122] Cloud and noise processing: Use the image's built-in quality assessment band (QA band) or a dedicated cloud detection algorithm (such as Fmask) to identify and generate masks for clouds and cloud shadows. In addition, it is necessary to repair existing stripe noise or bad lines (such as using the de-striping tool).

[0123] Step 2: High-resolution panchromatic image high-precision geometric correction

[0124] In the embodiment, Conora images of the required place and time are selected and downloaded from USGS, and loaded into professional geospatial commercial software, such as ArcGIS Pro, MapGIS, SuperMap, etc., with an ESRI global image base map with a resolution higher than 1 meter without offset.

[0125] Select control points (GCPs): Manually select points with the same name on the Corona image and reference data, following the control point selection principles (even distribution, stable and clear objects) and density requirements (differentiate terrain, such as 0.5-1 / km² in plains, 1.5-3 / km² in mountainous areas, a total of 200-300+). Record the image coordinates and corresponding geographic coordinates (or reference image coordinates) of each control point;

[0126] Select the calibration model and execute: Select the cubic spline interpolation function (Spline) algorithm as the geometric transformation model. Enter the selected GCPs and execute the calibration. Accuracy assessment: Check the RMSE value after calibration to ensure that it reaches the target accuracy (such as 1-3 meters). After completing the geometric calibration of the calibration model, a quasi-orthorectified Corona image is obtained.

[0127] Step 3: Multispectral image (Landsat / SPOT) geometric calibration

[0128] In the GIS software, the quasi-orthophoto Corona image corrected in step 2 is used as the reference image, and the multispectral image (Landsat / SPOT) is used as the image to be calibrated for geometric calibration. In this embodiment, the Landsat-5 TM image is used as an exemplary image to be calibrated (the calibration process of SPOT is similar), and the specific calibration process is as follows:

[0129] Select homologous points: Select homologous points between the multi-spectral image to be calibrated and the ortho-rectified Corona image.

[0130] Perform calibration: Use the cubic spline interpolation function (Spline) algorithm as the geometric transformation model. Based on the selected homologous points, perform geometric calibration (registration) on the multi-spectral image to ensure that the calibrated multi-spectral image and the Corona image are precisely superimposed spatially.

[0131] Step 4: Generate the fused high-resolution historical surface reference color image

[0132] This step aims to combine the spatial details of the high-resolution panchromatic image (PAN, i.e., the Corona image with a spatial resolution of 1 meter) corrected in Step 2 with the spectral information of the medium-resolution multi-spectral image (MS, i.e., the Landsat-5 TM image with a spatial resolution of 30 meters) calibrated in Step 3 to generate a historical surface reference color image with both high spatial resolution (1 meter) and multi-spectral information (inheriting the bands of Landsat TM). This step can be executed on local professional remote sensing / GIS software (such as ENVI, ArcGIS Pro) or cloud platforms (such as Google Earth Engine, GEE). Considering the potentially large data volumes of Corona and Landsat images, using the computing power of cloud platforms such as GEE has an efficiency advantage for large-area processing.

[0133] In this embodiment, the Gram-Schmidt (GS) transform fusion algorithm is adopted, and radiometric normalization (taking histogram matching as an example) is performed as a preprocessing step before fusion to optimize the fusion effect, especially color fidelity. The specific fusion process is as follows:

[0134] 1) Radiometric normalization - Histogram Matching

[0135] (a) Generate the Target Intensity Image: From the original Landsat-5 TM surface reflectance image with a 30-meter resolution, calculate a single-band intensity image representing its overall brightness. In this embodiment, the arithmetic mean method is used. Bands 1 (blue), 2 (green), 3 (red), and 4 (near-infrared) of Landsat-5 TM that overlap well with the Corona spectral range (visible light - part of the near-infrared) are selected for averaging. The formula for the averaged single-band intensity image Intensity_MS is: Intensity_MS = (Band1_reflectance + Band2_reflectance + Band3_reflectance + Band4_reflectance) / 4, where Band1_reflectance, Band2_reflectance, Band3_reflectance, and Band4_reflectance represent Bands 1 (blue), 2 (green), 3 (red), and 4 (near-infrared) in the Landsat-5 TM image, respectively. This calculation is performed at the original 30-meter resolution of the Landsat image.

[0136] (b) Upsample the Target Intensity Image: Using the Bilinear Interpolation method, resample the 30-meter resolution Intensity_MS image generated in the previous step to the same spatial resolution (1 meter) as the Corona image, thereby obtaining the high-resolution target intensity image Intensity_MS_resampled.

[0137] (c) Perform histogram matching: Use the original high-resolution Corona image (PAN) corrected in step 2 as the source image (Source) and Intensity_MS_resampled as the target / reference image (Target / Reference). Perform histogram matching using image processing software (such as ENVI's Histogram_Match tool or the skimage.exposure.match_histograms function in Python). This process adjusts the pixel values ​​(grayscale) of the original Corona image so that its cumulative frequency histogram matches the cumulative frequency histogram of Intensity_MS_resampled as much as possible. The output is the radiometrically matched panchromatic image Corona_matched. This image retains the original 1-meter spatial resolution and ground object details of Corona, but its overall brightness and contrast distribution has been adjusted to be closer to the intensity characteristics of Landsat images.

[0138] 2) Gram-Schmidt (GS) fusion:

[0139] Image fusion can be implemented on local GIS / remote sensing software (such as ArcGIS Pro's "Create Pan-sharpenedRaster Dataset" tool selecting the GS method, or the corresponding module of ENVI) or cloud platforms (such as GEE).

[0140] Local software example (ArcGIS Pro): Use the Create Pan-sharpened Raster Dataset tool, input the rectified Corona image (as the panchromatic band) and the calibrated Landsat TM multispectral image (as the multispectral band), and select the Gram-Schmidt fusion method. Adjust the weight parameters as needed (usually the software has default values ​​or provides optimization options).

[0141] Cloud platform example (GEE): Upload the corrected / calibrated image to GEE Assets. Write a script to call the pan-sharpening function provided by GEE (you may need to find a function or community script that supports the GS algorithm) to perform fusion. GEE's advantage lies in processing large amounts of data and parallel computing.

[0142] Parameter optimization and quality evaluation: Visually inspect the spatial clarity, color naturalness, and presence of obvious artifacts of the fusion results, calculate quantitative indicators (such as ERGAS, SAM, CC, RMSE) and compare them with the original images, and adjust the fusion parameters (if the software allows) to obtain the best results. Ensure that the quality of the fused image meets the preset standards. Finally, output the fused high-resolution, multi-spectral historical surface benchmark color image.

[0143] Step 5: Perform 3D visualization by integrating DEM

[0144] In this embodiment, in the GIS software ArcGIS Pro that supports 3D scenes, load the historical surface reference color image generated in Step 4. At the same time, load the DEM data corresponding to this area (it is recommended to use ALOS 12.5m DEM or higher-precision data, which needs to be pre-processed to be consistent with the image coordinate system). Set the historical surface reference color image as the source of the surface color of the 3D scene, and set the DEM as the source of the elevation of the 3D scene, and create a local scene or a global scene for 3D display. Thus, a 3D reference base map can be generated through 3D visualization technology, and users can freely zoom in, rotate, and tilt the viewing angle through mouse operations for stereoscopic observation.

[0145] As Figure 3 shown, it shows the 2D reference base map (historical surface reference color image) and the 3D reference base map (3D visualization scene) generated based on the original Conora satellite image and Landsat TM multispectral image. It can be seen from this that the original Conora (black and white) satellite image has a high spatial resolution and clearer details of ground objects, but it is a black and white image, a central projection, without geographical coordinates, and the scales of different parts are inconsistent, and the whole image has distortion. The Landsat TM multispectral (color) image has a spatial resolution of 30 meters, but the details of ground objects cannot be seen clearly and there are geometric distortions. After fusing the two according to the method of this embodiment, the obtained historical surface reference color image is a color image (partial) with both high spatial resolution and multispectral information, and can clearly display different land type information and its boundary details. Further, after assisting with DEM information and establishing a 3D visualization scene, it is more convenient to display, observe, analyze, and judge soil types and ground object types from different angles and directions, which is very beneficial for the correction and high-quality improvement of soil maps.

[0146] Step 6: Geometric correction of the second national soil survey soil map

[0147] In this embodiment, after overlaying the second national soil survey map in the GIS software, the historical reference color image completed in step 4 is used as an accurate reference base map. Using a georegistration tool (such as the Georeferencing toolbar or Spatial Adjustment tool in ArcGIS Pro), homologous ground control points are selected between the second national soil survey map and the historical reference image (such as features that can be clearly identified on both, like road intersections, river confluences, reservoir corners, special geomorphic points, etc.). The Rubber sheet transformation model is selected to perform geometric correction. Accuracy verification: The independent check points (ICPs) are used to evaluate the correction accuracy, with the goal of controlling the position error of the patch boundaries of the second national soil survey map within the range of 2 - 5 meters.

[0148] As Figure 4 shown in the figure, a) and b) respectively exemplarily show the effects formed after overlaying the second national soil survey map on the two-dimensional historical reference color image in mountainous areas and plain areas. It can be seen from this that the deviation (displacement) of the boundaries of mountain soil types (features) often exceeds 60 meters, and locally reaches 120 meters. After overlaying, relying on the image fused with high spatial resolution and multi-spectral information as the base map, it is convenient to select place names and ground control points for geometric correction of the second national soil survey map, and at the same time, it is also easy to adjust the boundaries of features and verify the types of features. Similarly, the displacement (deviation) of the soil map (feature type) in the plain area is usually within 30 - 40 meters. After overlaying, relying on the satellite image fused with high spatial resolution and multi-spectral information as the base map, it is easy to adjust the boundaries of features and verify the types of features.

[0149] Step 7: Content correction and improvement of the second national soil survey map

[0150] In the GIS software, the second national soil survey map corrected in step 6 is overlaid on the historical surface reference color image in step 4 (or the 3D visualization scene in step 5) to correct and improve the content of the second national soil survey map. Specifically, it can include starting editing, boundary adjustment and addition, and content verification and correction. The methods are as follows:

[0151] Starting editing: Start an editing session for the vector data of the second national soil survey map.

[0152] Boundary adjustment and addition: Referring to the base map image (and 3D scene), according to the aforementioned boundary adjustment basis and content update standard requirements, use editing tools to accurately adjust the boundaries of the existing soil patches to make them coincide with the historical surface landscape features. For soil units omitted in the original map or not drawn due to scale limitations, new boundaries and patches are added.

[0153] Content verification and correction: Check the soil type attribute annotations of each patch. Referring to the image features and soil knowledge, correct the incorrect annotations. Merge overly fragmented or misinterpreted patches, or subdivide large patches as needed.

[0154] In addition, Figure 5 further shows the soil maps before and after correction in an exemplary scenario, as well as the comparison results of the map patches between the two. It can be seen that by using the color satellite images fused with high spatial resolution and multi-spectral information, assisted by DEM information, and performing three-dimensional stereoscopic observation, it is possible to well correct the patch boundaries of the second general soil map, verify and correct the patch content, and complete the topological verification. Figure 5 From the comparison results before and after correction, it can be seen that the blue line segment is the original soil map, and the red one is the boundary of the corrected soil map. It can be clearly shown that the spatial position boundaries of the ground objects in the soil map after spatial geometric correction generally coincide with the ground objects. In this scenario, the total number of map patches in the display area is 16. By using the fused color images and terrain digital elevation models for interpretation, eight new ground object boundaries are added, and eight new map patches are added. Among them, the three map patches corresponding to the three boundaries were missed during the original mapping at the original mapping scale, and the five ground object boundaries corresponding to the five map patches are newly added according to the current resolution; at the same time, two map patch boundaries are deleted, and one map patch is merged. It can be seen that through this spatial position correction, patch content verification and correction, the quality improvement effect of the second general soil map in this area is extremely obvious.

[0155] Step 8: Topological Check and Final Result Output

[0156] In this embodiment, further create topological rules (such as Must Not Overlap, Must Not Have Gaps, etc.) for the corrected second general soil layer in the GIS software (ArcGIS Pro), and find and repair all topological errors by running the topological verification.

[0157] Final Quality Evaluation: Conduct a final quality check, including three aspects: geometric accuracy, topological consistency, and content rationality:

[0158] Geometric Accuracy: Use ICPs to confirm again that the RMSE is within the range of 2 - 5 meters.

[0159] Topological Consistency: Confirm that there are no topological errors.

[0160] Content Rationality: Soil experts conduct a sampling inspection to evaluate the boundary rationality and attribute accuracy.

[0161] When passing the final quality check, the final result output can be carried out, and the digital results of the finally completed high-precision and content-perfect second general soil map (such as standard formats like shapefile, geodatabase feature class, etc.) can be saved.

[0162] In addition, two comparative experiments were designed in this embodiment to verify the superiority of using the high-precision historical reference color images of the present invention for the correction and improvement of the second national soil survey maps.

[0163] 1. Experimental Area and Data

[0164] The longitude and latitude range of the experimental area: 120.5315 - 120.5684 E, 28.0844 - 28.1071 N

[0165] The original data obtained for the experimental area is as follows:

[0166] 1) The original paper or scanned version of the second national soil survey map of the experimental area; 2) The high-resolution panchromatic images of the Corona KH series in November 1983 during the second national soil survey covering the area (spatial resolution of about 1 - 2 meters); 3) The Landsat TM multispectral images covering the area during the same period (spatial resolution of 30 meters); 4) The modern high-precision reference data ESRI basemap for initial geometric correction; 5) Digital elevation model (DEM), ALOS PALSAR 12.5m DEM.

[0167] The evaluation indicators used in the experiment are as follows:

[0168] Geometric correction accuracy: Using independent check points (ICPs), calculate the root mean square error (RMSE) of the corrected map patch boundaries of the second national soil survey map relative to the "true" historical surface position (using the high-precision historical color reference image constructed by the method of the present invention or the high-precision Corona image based on it as a reference).

[0169] The content improvement effect is evaluated from two dimensions: qualitative and quantitative:

[0170] Qualitative evaluation: Soil science or remote sensing interpretation experts, with the assistance of their respective reference images (or 3D views), evaluate the coincidence degree of the map patch boundaries of the corrected soil map with the landform / features, the homogeneity within the map patches, and the reliability of the map patch attribute interpretation (excellent, good, medium, poor).

[0171] Quantitative evaluation (reference): During the content correction process, count the degree of map patch boundary adjustment, the number and area of newly added / deleted / merged map patches, the number of attribute corrections, etc.

[0172] 2. Comparative Experiment Design

[0173] The following two comparative experiments were set up. Except for the different core reference images, other processing steps (such as initial image preprocessing, the transformation model used for geometric correction of the second national soil survey map - rubber sheet transformation, the operation process of content correction, topological inspection, etc.) were kept as consistent as possible with this embodiment.

[0174] Comparative Experiment 1: Using only historical medium-resolution multispectral images as the benchmark

[0175] (1) Obtain and preprocess contemporaneous Landsat TM images (radiometric calibration, atmospheric correction, etc.).

[0176] (2) Geometrically correct the Landsat TM images using modern high-precision reference data.

[0177] (3) Use the corrected Landsat TM images (30-meter resolution) as the benchmark base map.

[0178] (4) Select homologous points on this benchmark base map and geometrically correct the original Second National Soil Survey map (rubber stretching transformation).

[0179] (5) Amend and improve the content of the Second National Soil Survey map on this benchmark base map (which can be combined with DEM to construct a 3D view).

[0180] (6) Conduct topological checks and output the results, and then evaluate them.

[0181] Comparative Experiment 2: Using only historical high-resolution panchromatic images as the benchmark

[0182] (1) Obtain and preprocess contemporaneous Corona images (contrast enhancement, etc.).

[0183] (2) Geometrically correct the Corona images with high precision using modern high-precision reference data (such as Spline transformation, the same as step 2 of the present invention).

[0184] (3) Use the corrected Corona images (1-2-meter resolution, black and white) as the benchmark base map.

[0185] (4) Select homologous points on this benchmark base map and geometrically correct the original Second National Soil Survey map (rubber stretching transformation).

[0186] (5) Amend and improve the content of the Second National Soil Survey map on this benchmark base map (which can be combined with DEM to construct a 3D view).

[0187] (6) Conduct topological checks and output the results, and then evaluate them.

[0188] 3. The experimental results and analysis are as follows:

[0189] Table 1 Experimental results of different experimental groups

[0190]

[0191] As can be seen from the experimental results in Table 1, the base map obtained from Comparative Experiment 1 has a low resolution and cannot provide fine feature boundary information, resulting in insufficient geometric positioning accuracy and difficulty in accurately identifying and delineating small patches and precisely adjusting the boundaries. The base map obtained from Comparative Experiment 2 has a high spatial resolution and good geometric positioning accuracy. However, it lacks color and spectral information, making it difficult to interpret surface differences caused by vegetation, water bodies, and different soil types, seriously affecting the accurate verification and correction of patch content (especially attributes). The base map obtained in this embodiment (i.e., the historical surface reference color image) combines high spatial resolution and multi-spectral information, ensuring both meter-level geometric positioning accuracy and providing rich feature information (color, texture), greatly facilitating and improving the efficiency and quality of accurate adjustment of patch boundaries, supplementation of missing patches, and accurate verification and correction of patch attributes.

[0192] The results of the above comparative examples clearly show that using only historical medium-resolution multi-spectral images cannot meet the accuracy requirements of meter-level geometric correction for the second national soil survey maps, and the content correction effect is poor. Using only historical high-resolution panchromatic images, although high geometric positioning accuracy can be achieved, due to the lack of spectral information, the accuracy and reliability of content correction (especially attribute interpretation) are greatly reduced. The technical solution of constructing and using the fused high-precision historical color reference image proposed in the present invention is the only effective way to simultaneously achieve meter-level geometric correction accuracy and high-quality content correction and improvement. By integrating the advantages of the two historical data sources, it overcomes the limitations of a single data source and is the key to solving the problems of accurate correction and high-quality improvement of the second national soil survey maps.

[0193] Statistics of the experimental results for a large number of regions show that by applying the method of the present invention in this embodiment, the average geometric position error of the second national soil survey maps in the region can be successfully reduced from more than about 50 meters to 2 - 5 meters. In terms of content, through verification and correction against the high-precision historical reference image, the number of patches increased from 8 to 16, doubling, greatly improving the accuracy and detail of the soil maps. This proves that the method of the present invention has obvious effects on the correction and high-quality improvement of historical geological maps such as the second national soil survey maps.

[0194] The above-described embodiments are only some preferred implementation solutions of the present invention, but are not intended to limit the present invention. Those of ordinary skill in the relevant technical fields can make various changes and modifications without departing from the spirit and scope of the present invention. Therefore, all technical solutions obtained by means of equivalent replacement or equivalent transformation fall within the protection scope of the present invention.

Claims

1. A method for geometric precision correction and improvement of historical geological maps, characterized in that, Including: S1. Obtain corresponding high-resolution panchromatic images and medium-resolution multispectral images according to the spatio-temporal range corresponding to the historical geoscience map to be corrected, and perform preprocessing on them respectively; S2. Use the standard geospatial information as the first reference image, select control points on the preprocessed high-resolution panchromatic image and the first reference image, and perform geometric correction on the high-resolution panchromatic image; S3. Use the geometrically corrected high-resolution panchromatic image as the second reference image, select control points on the preprocessed medium-resolution multispectral image and the second reference image, and perform geometric correction on the medium-resolution multispectral image; S4. According to the spectral range of the high-resolution panchromatic image, extract the multi-band image within the spectral range from the geometrically corrected medium-resolution multispectral image and convert it into a single-band intensity image, and then perform histogram matching on the single-band intensity image and the geometrically corrected high-resolution panchromatic image, and generate a historical surface reference color image through image fusion; S5. Use the historical surface reference color image as a two-dimensional reference base map to perform geometric correction and patch content correction on the historical geoscience map.

2. The geometric precision correction improvement method for historical geological maps according to claim 1, wherein, After obtaining the historical surface reference color image, use the historical surface reference color image as the surface color source and the digital elevation model data within the spatio-temporal range as the elevation source, and generate a three-dimensional reference base map through three-dimensional visualization technology to perform geometric correction and patch content correction on the historical geoscience map.

3. The geometric precision correction improvement method for historical geological maps according to claim 1, wherein The historical geoscience map is a soil map, a land use map, a vegetation map, a geological map or a geomorphic map.

4. The geometric precision correction improvement method for historical geological maps according to claim 1, characterized in that The historical geoscience map is a second national soil survey soil map, the high-resolution panchromatic image is a Conora satellite image, and the medium-resolution multispectral image is a Landsat MSS image, a Landsat TM image or a SPOT HRV image.

5. The method for enhancing the geometric precise correction of historical geological maps according to claim 1, wherein The georeference data uses a topographic map with a scale of not less than 1:10,000 or an orthophoto image with a spatial resolution of not less than 1 meter.

6. The method for improving the geometric precise correction of historical geological maps according to claim 1, characterized in that, The single-band intensity image is obtained by weighted averaging of the multi-band images within the spectral range, or uses the first principal component of the multi-band images within the spectral range, or uses the first component of the Gram-Schmidt transform of the multi-band images within the spectral range.

7. The geometric precision correction improvement method for historical geological maps according to claim 1, characterized in that When performing histogram matching, it is necessary to upsample the single-band intensity image to the same resolution as the geometrically corrected high-resolution panchromatic image, and then perform histogram matching on the two, so that the cumulative histogram of the high-resolution panchromatic image approaches the upsampled single-band intensity image. Finally, the two types of images after histogram matching are fused through the Gram-Schmidt algorithm to generate a historical surface reference color image.

8. The method for improving the geometric precise correction of historical geological maps according to claim 2, characterized in that, When performing geometric correction and patch content correction on the historical geoscience map based on the two-dimensional reference base map and the three-dimensional reference base map, it is necessary to select the same-name ground feature points in the same period on the reference base map and the historical geoscience map as control points, perform geometric correction on the historical geoscience map, and then compare with the reference base map to correct and improve the soil patch boundaries and attribute contents on the geometrically corrected historical geoscience map. Finally, after topological relationship inspection and correction, the final map is output.

9. A computer-readable storage medium, characterized in that, A computer program is stored on the storage medium, and when the computer program is executed by a processor, the method for improving the geometric precise correction of historical geological maps according to any one of claims 1 to 7 is implemented.

10. A computer electronic device, characterized in that, It includes a memory and a processor; The memory is used for storing a computer program; The processor is used for implementing the method for improving the geometric precise correction of historical geological maps according to any one of claims 1 to 7 when executing the computer program.

Images