A multi-source satellite remote sensing image automatic rectification and quality inspection method and system

Abstract

The application discloses a kind of multi-source satellite remote sensing image automation rectification and quality inspection method and system, it is related to remote sensing image processing technical field, the method includes: obtaining remote sensing image data and accessory parameter, digital elevation model data and reference orthographic image data;Rational polynomial coefficient geometric positioning model is constructed, through image matching collection same name control point, utilize digital elevation model to calculate topographic feature parameter and construct prior distortion intensity index and spatial heterogeneity weight matrix, through regularized weighted adjustment solution affine transformation correction parameter completes orthographic rectification;After rectification image is carried out quality inspection;Evaluation result is constructed as the basic probability assignment of three evidence body, conflict coefficient is calculated and is fused decision, and quality grade is output;According to quality grade, unqualified image triggers re-correction processing, and qualified image is cut and result output.The application improves the rectification precision and reliability of quality inspection of mountainous area image.

Description

Technical Field

[0001] This invention relates to the field of remote sensing image processing technology, specifically to a method and system for automated correction and quality inspection of multi-source satellite remote sensing images. Background Technology

[0002] With the continuous increase in the number of remote sensing satellites, multi-source remote sensing image data is being used more and more widely in fields such as natural resource surveys, ecological monitoring, and land mapping. Image data from different satellite platforms and different sensor types have significant differences in spatial resolution, radiometric characteristics, and geometric positioning accuracy, posing a huge challenge to automated correction and quality inspection.

[0003] Currently, orthorectification of remote sensing images mainly employs the regional network adjustment method based on the RPC model. However, existing methods have the following shortcomings: traditional regional network adjustment treats all control points equally, failing to fully consider the varying impacts of terrain undulations on correction accuracy, leading to difficulties in effectively suppressing local streaking in mountainous areas; current quality inspections primarily rely on planar position errors, averaging all checkpoints equally, making it impossible to distinguish whether errors occur in flat farmland or steep slopes; in radiometric quality assessment, information entropy cannot distinguish between clear edges and noise textures, the average gradient is sensitive to illumination changes, and the comparability of images from different time phases is poor.

[0004] Therefore, in order to address the above problems, there is an urgent need for a comprehensive quality assessment method that can achieve regularized adjustment for terrain perception, structural fidelity evaluation based on visual perception mechanisms, and multi-evidence conflict decision-making. Summary of the Invention

[0005] The technical problem to be solved by the present invention is to provide an automated correction and quality inspection method and system for multi-source satellite remote sensing images, which addresses the shortcomings of the existing technology.

[0006] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0007] An automated correction and quality inspection method for multi-source satellite remote sensing images includes the following steps:

[0008] Step S1: Acquire remote sensing image data and its associated parameters, and simultaneously acquire digital elevation model data and reference orthophoto data. Perform format parsing and specification verification on the remote sensing image data, and extract image parameters.

[0009] Step S2: Construct the original geometric positioning model of the image based on the rational polynomial coefficient model, collect corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, construct the prior distortion intensity index and spatial heterogeneity weight matrix based on the terrain feature parameters, solve the affine transformation correction parameters through regularized weighted adjustment, and complete the orthorectification of panchromatic and multispectral images respectively.

[0010] Step S3: Fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image;

[0011] Step S4: Perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension, and semantic usability dimension.

[0012] Step S5: Construct the evaluation results of geometric accuracy dimension, radiometric quality dimension and semantic usability dimension into basic probability values ​​for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image.

[0013] Step S6: Based on the quality level, mark the error exceeding the limit area for unqualified images and trigger re-correction processing; for qualified and above images, perform cropping, metadata generation and output according to standard framing rules.

[0014] Furthermore, in step S1, the remote sensing image data is multi-source remote sensing image data from different satellite platforms, including satellite image data with different sensor types, different spatial resolutions, different imaging phases, and different coordinate systems; the image parameters include sensor type, spatial resolution, number of bands, imaging phase, satellite side tilt angle, and flight azimuth angle; the auxiliary parameters include at least rational polynomial coefficient parameters, metadata, quick view, and combined chart.

[0015] Furthermore, step S2 specifically includes the following steps:

[0016] Step S2.1: Construct the original geometric positioning model of the image based on the rational polynomial coefficient model, and use the reference orthophoto to perform image matching and automatically collect corresponding control points.

[0017] Step S2.2: Calculate the terrain feature parameters of each control point based on the digital elevation model, and construct the prior distortion intensity index based on the terrain feature parameters, satellite side sway angle and satellite flight azimuth angle;

[0018] Step S2.3: Construct a spatial heterogeneity weight matrix based on the prior distortion intensity index, and obtain the affine transformation correction parameters of the rational polynomial coefficient model by solving the regularized weighted least squares objective function.

[0019] Further, in step S2.2, the terrain feature parameters include slope angle and aspect angle; the method for constructing the prior distortion intensity index is as follows: calculate the square cosine of the angle between the aspect angle and the satellite flight azimuth angle, and the prior distortion intensity index is equal to the tangent of the satellite side swing angle multiplied by the slope magnitude and then multiplied by the product of a directional coupling coefficient and the square cosine of the angle.

[0020] Further, in step S2.3, the adaptive weight of the k-th control point in the spatial heterogeneity weight matrix is ​​calculated as follows: first, the ratio of the prior distortion intensity index to the arithmetic mean of the prior distortion intensity indices of all control points is calculated, the negative of this ratio is taken as the exponent of the natural constant e, and then multiplied by the ratio of the control point matching correlation coefficient to the preset confidence threshold.

[0021] The regularized weighted least squares objective function includes two terms: the first term is the sum of the products of the adaptive weights of all control points and the squares of their image-side residual vectors; the second term is the sum of the regularization penalty coefficient multiplied by the squares of the image-side coordinate correction difference vectors of all adjacent control point pairs; the image-side residual vector is the difference between the actual measured image-side coordinates and the theoretical image-side coordinates after six-parameter affine correction.

[0022] Furthermore, step S3 specifically includes the following steps:

[0023] Step S3.1: Perform histogram matching on the panchromatic image and the multispectral image to ensure that the grayscale distribution of the two is consistent.

[0024] Step S3.2: Convert the multispectral image from RGB space to IHS space and extract the luminance component;

[0025] Step S3.3: Replace the luminance component with the histogram-matched panchromatic image to obtain a new luminance component;

[0026] Step S3.4: Perform an inverse IHS transformation on the new luminance component and the original chrominance and saturation components to transform them back to RGB space and generate a fused digital orthophoto.

[0027] Furthermore, step S4 specifically includes the following steps:

[0028] Step S4.1: Evaluate the geometric accuracy dimension based on the prior distortion weighted geometric deviation index. The geometric deviation index is calculated as follows: For each independent high-precision check point, the weight is the prior distortion intensity index of the check point divided by the reciprocal of the arithmetic mean of the prior distortion intensity indices of all control points plus one. The weight is then multiplied by the planar distance between the measured coordinates of the check point on the corrected image and the true coordinates on the high-precision reference orthophoto. Finally, the sum of the weighted products of all check points is divided by the sum of the weights of all check points.

[0029] Step S4.2: Evaluate the radiometric quality dimension based on the fidelity index. The fidelity index is calculated as follows: For each pixel position in the image, a multi-directional Log-Gabor filter is used to calculate the local amplitude and phase deviation in each direction. The phase consistency value of the pixel is equal to the sum of the products of the local amplitude in each direction and the positive part of the phase deviation in that direction, divided by the sum of the local amplitudes in all directions plus a very small positive number; the fidelity index is equal to the arithmetic mean of the phase consistency values ​​of all effective pixels in the image.

[0030] Step S4.3: Extract invalid coverage areas from the image using a semantic segmentation model and calculate the effective coverage rate as an evaluation dimension of semantic availability. The invalid coverage areas include at least clouds, cloud shadows, water bodies, and ice and snow. The effective coverage rate is equal to the total number of effective pixels divided by the total number of pixels in the image.

[0031] Furthermore, step S5 specifically includes the following steps:

[0032] Step S5.1: Map the evaluation results of the geometric precision dimension, radiation quality dimension, and semantic usability dimension to the basic probability assignments of three independent evidence bodies, namely geometric evidence, radiation evidence, and semantic evidence, respectively, through Gaussian membership functions.

[0033] Step S5.2: Calculate the conflict coefficient between the three pieces of evidence. The method for calculating the conflict coefficient is as follows: traverse the quality levels supported by each of the geometric evidence, radiation evidence and semantic evidence, and accumulate the product of the three basic probability assignments corresponding to each combination where the intersection of the three evidence judgment results is empty.

[0034] Step S5.3: Calculate the overall trust level after fusion using the conflict allocation combination rule, and output the quality level of the image based on the overall trust level.

[0035] Further, in step S5.3, the conflict allocation combination rule is as follows: for each possible quality level, the fusion confidence level is equal to the sum of the products of the basic probability assignments of all combinations of the three evidence judgment results, which are exactly the intersection of the three evidence judgment results, and divided by the denominator; the denominator is equal to a product of the conflict allocation coefficient and the conflict coefficient minus the conflict coefficient; if the confidence level corresponding to the quality level with the highest fusion confidence level is greater than 0.6, then the quality level is selected as the final quality level, otherwise it is judged as suspicious and awaits manual review.

[0036] An automated correction and quality inspection system for multi-source satellite remote sensing images, used to implement any one of the methods for automated correction and quality inspection of multi-source satellite remote sensing images, includes:

[0037] The data acquisition module is used to acquire remote sensing image data and its associated parameters, perform format parsing and specification verification on the remote sensing image data, and extract image parameters.

[0038] The image correction module is used to construct the original geometric positioning model of the image based on the rational polynomial coefficient model, acquire corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, and then construct the prior distortion intensity index and spatial heterogeneity weight matrix. The affine transformation correction parameters are solved through regularized weighted adjustment to complete the orthorectification of panchromatic and multispectral images respectively.

[0039] The fusion processing module is used to fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image.

[0040] The quality inspection module is used to perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension and semantic usability dimension.

[0041] The comprehensive scoring module is used to construct the evaluation results of geometric accuracy, radiometric quality and semantic usability into basic probability assignments for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image.

[0042] The feedback optimization module is used to mark the error-exceeding areas of unqualified images and trigger re-correction processing based on the quality level.

[0043] The output module is used to crop qualified and above images according to standard framing rules, generate metadata, and output the results.

[0044] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0045] 1. This invention models the relative direction of slope magnitude, slope aspect, satellite side swing angle, and flight azimuth angle by coupling the model, so that the distortion penalty is maximized when the slope aspect is consistent with the flight direction and minimized when it is perpendicular. Based on this, a spatial heterogeneity weight matrix is ​​constructed, and the weight of control points with high distortion intensity is actively reduced, avoiding the drawback of traditional equal weight adjustment that sacrifices overall accuracy to fit high distortion areas, and effectively suppressing the streaking distortion of mountain images.

[0046] 2. This invention introduces a fidelity index to replace the traditional information entropy for radiation quality evaluation. It is naturally insensitive to changes in illumination and contrast, can accurately distinguish between clear edges and noisy textures, and is more in line with the perception mechanism of the human visual system.

[0047] 3. This invention uses the DS evidence theory for comprehensive quality assessment. By calculating the conflict coefficient between the three evidence bodies, the geometric, radial, and semantic indicators are mutually constrained. When the evidence is highly contradictory, it is automatically downgraded or judged as suspicious / awaiting manual review, thus avoiding the drawback of the weighted summation method where a high score can cover up other defects. Attached Figure Description

[0048] Other features, objects, and advantages of the invention will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings:

[0049] Figure 1 This is a flowchart illustrating an embodiment of the present invention;

[0050] Figure 2 This is a system schematic diagram according to an embodiment of the present invention;

[0051] Figure 3 This is a comparison image of the road deformation before and after correction according to an embodiment of the present invention. Detailed Implementation

[0052] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be described in detail below with reference to the accompanying drawings and specific embodiments.

[0053] like Figure 1 As shown, an automated correction and quality inspection method for multi-source satellite remote sensing images includes the following steps:

[0054] Step S1: Acquire remote sensing image data and its associated parameters, and simultaneously acquire digital elevation model data and reference orthophoto data. Perform format parsing and specification verification on the remote sensing image data, and extract image parameters.

[0055] Step S2: Construct the original geometric positioning model of the image based on the rational polynomial coefficient model, collect corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, construct the prior distortion intensity index and spatial heterogeneity weight matrix based on the terrain feature parameters, solve the affine transformation correction parameters through regularized weighted adjustment, and complete the orthorectification of panchromatic and multispectral images respectively.

[0056] Step S3: Fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image;

[0057] Step S4: Perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension, and semantic usability dimension.

[0058] Step S5: Construct the evaluation results of geometric accuracy dimension, radiometric quality dimension and semantic usability dimension into basic probability values ​​for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image.

[0059] Step S6: Based on the quality level, mark the error exceeding the limit area for unqualified images and trigger re-correction processing; for qualified and above images, perform cropping, metadata generation and output according to standard framing rules.

[0060] Acquire remote sensing image data and its associated parameter files from different satellite platforms, including RPC parameter files, metadata files, quick views, and stitching charts. The RPC parameter files contain rational polynomial coefficients for subsequent geometric positioning model construction. The metadata files contain information such as the image's sensor type, spatial resolution, number of bands, imaging time phase, satellite side tilt angle, and flight azimuth angle. The quick views provide an overview of the data quality. The stitching charts are used to determine the image coverage and its edge-to-edge relationship with neighboring scenes.

[0061] Simultaneously acquire digital elevation model data covering the same area and existing high-precision reference orthophoto data. The DEM data is used for terrain projection difference correction in the subsequent correction process, and its accuracy directly affects the elevation correction effect of orthophoto correction. The reference DOM is used as the reference base map for automatic control point acquisition, and its planar accuracy should be better than the target accuracy of the image to be corrected.

[0062] The acquired raw image data is parsed and its integrity is verified to ensure that the data is readable and undamaged. The sensor type is extracted to determine the correction strategy, the spatial resolution is extracted to determine the subsequent cropping and framing specifications, the number of bands is extracted to determine the fusion scheme, and the side sway angle and flight azimuth angle are extracted for subsequent calculation of the a priori distortion intensity index.

[0063] The original geometric positioning model of the image is constructed using RPC parameters. The RPC model is a mathematical model widely used in satellite image geometric processing. It does not depend on specific sensor physical parameters, but establishes a mathematical mapping relationship between image-side coordinates and object-side coordinates through rational polynomial functions. Its coefficients are directly read from the RPC parameter file. This mapping relationship describes the correspondence between the three-dimensional coordinates of ground points and the image pixel coordinates. Based on this model, the scale-invariant feature transformation algorithm or the correlation coefficient method is used to automatically match the image to be corrected with the reference orthophoto to obtain corresponding control points. A total of K control points are obtained.

[0064] Traditional adjustment methods treat all control points equally, without considering the differential impact of terrain undulation on correction accuracy. In mountainous areas, the projection difference of steep terrain is much greater than that of flat areas. Equal weight adjustment makes it difficult to effectively suppress local ripple deformation, often requiring manual point addition and redoing. This invention uses DEM data to pre-calculate the terrain features at each control point, and then constructs adaptive weights, enabling the adjustment process to perceive terrain differences.

[0065] For the kth control point, the slope angle and aspect angle of that point are calculated using the digital elevation model, with the following formula:

[0066]

[0067] in, This represents the slope angle of the k-th control point, reflecting the steepness of the terrain at that location. A larger value indicates a steeper terrain. This represents the elevation value in the digital elevation model. This represents the X-direction coordinate in the Gauss-Kruger projection coordinate system, i.e., eastward. This represents the Y-direction coordinate in the Gauss-Kruger projection coordinate system, i.e., the north direction. This represents the coordinates of the kth control point in the digital elevation model. Rate of elevation change along the X direction at a given location This represents the coordinates of the kth control point in the digital elevation model. The rate of change of elevation along the Y direction, This represents the k-th control point. This represents the slope angle of the k-th control point, reflecting the orientation of the terrain slope.

[0068] The formula for calculating the prior distortion intensity index is as follows:

[0069]

[0070] in, This represents the prior distortion intensity index of the k-th control point. The larger the value, the greater the projection difference caused by the coupling of terrain undulation and lateral observation at that ground point, and the easier it is to produce streaking distortion during correction. This represents the tilt angle of the satellite sensor, directly read from the image metadata file. Its tangent value characterizes the projection amplification effect caused by tilted observation; the larger the tilt angle, the larger this value. The value is a non-negative angle, i.e., the absolute value of the tilt angle. If the metadata record is a signed angle, its absolute value is used in the calculation. This represents the directional coupling coefficient, which is dimensionless and ranges from 0.5 to 1.5, with 1.0 being preferred. This coefficient maximizes the distortion penalty when the slope is aligned with the flight direction and minimizes it when it is perpendicular to the flight direction. Indicates the satellite's azimuth angle;

[0071] Based on the prior distortion intensity index, the adaptive weight of the k-th control point in the adjustment objective function is defined, and the specific formula is as follows:

[0072]

[0073] in, This represents the arithmetic mean of the prior distortion intensity indices of all control points on the current image. This represents the adaptive weight of the k-th control point. This represents the correlation coefficient of the k-th control point in image matching, ranging from 0 to 1. It is directly output by the image matching algorithm and reflects the matching reliability of this control point. The closer the value is to 1, the higher the matching quality. This represents a preset confidence threshold, ranging from 0.60 to 0.70, with 0.65 being preferred.

[0074] The objective function for regularized adjustment is constructed as follows:

[0075]

[0076] in, This represents the objective function of regularized adjustment. This indicates the total number of control points. This represents the residual vector of the k-th control point in the image-side coordinate system, which is the difference between the actual measured image-side coordinates and the theoretical image-side coordinates after six-parameter affine correction. It reflects the geometric positioning deviation of the RPC model at this control point. This represents the regularization penalty coefficient, ranging from 0.01 to 0.05. This represents the set of control point indices adjacent to the k-th control point. This represents the difference vector between adjacent control points k and j in the image-side coordinate correction, i.e., the difference between the correction of the k-th control point and the correction of the j-th control point. It is constrained by the regularization penalty coefficient so that spatially adjacent control points have similar correction values, thereby ensuring a smooth transition of correction parameters in the local area.

[0077] For the k-th control point, its object coordinates are: The initial image-side coordinates obtained by substituting them into the RPC model projection calculation are: After applying a six-parameter affine correction, the corrected theoretical image coordinates are... The calculation formula is:

[0078]

[0079] in, Indicates the amount of translation in the row direction. Indicates the column direction translation amount. The first linear coefficient representing the row direction. The second linear coefficient, representing the row direction, corresponds to the linear corrections for the theoretical x-axis and y-axis of the RPC, respectively. Indicates the first linear coefficient in the column direction. The second linear coefficient, representing the column direction, corresponds to the linear corrections for the theoretical horizontal and vertical axes of the RPC, respectively.

[0080] residual vector For measured image-side coordinates Compared with the corrected theoretical image coordinates The difference is used to substitute the residual vector into the objective function. , The equation becomes a quadratic function with respect to 6 parameters. Taking the partial derivatives with respect to each parameter and setting them to zero yields a system of 6 × 6 linear equations. The Cholesky decomposition method is used to solve this system when the parameter variation is less than 10. -6 When convergence is determined, the above correction process is performed on the panchromatic image and the multispectral image respectively to complete the orthorectification of both.

[0081] The orthorectified panchromatic image and the multispectral image are fused. First, histogram matching is performed on the panchromatic and multispectral images to ensure that their grayscale distributions are consistent, thus eliminating brightness differences caused by different imaging methods and dynamic ranges. Then, the multispectral image is converted from RGB space to IHS space, and the luminance component I, chrominance component H, and saturation component S are extracted. The original luminance component I is replaced by the histogram-matched panchromatic image to obtain a new luminance component, while the chrominance component H and saturation component S remain unchanged. Finally, the new luminance component is combined with H and S through inverse IHS transformation and transformed back to RGB space to generate the fused digital orthorectified image.

[0082] To perform quality inspection on the fused digital orthophoto, traditional quality inspection uses the Planar Position Error (RMSE) as a geometric accuracy evaluation index. However, RMSE is an equal-weighted average of all checkpoints, which cannot distinguish whether the error occurs in flat farmland or steep slopes. This invention proposes a geometric accuracy evaluation based on a priori distortion-weighted geometric deviation index. It uses independent high-precision checkpoints and a priori distortion intensity index for weighting, with the specific formula as follows:

[0083]

[0084] in, Indicates the geometric deviation index. This represents the total number of independent checkpoints. This represents the planar coordinates of the t-th checkpoint measured on the corrected image. This represents the true planar coordinates of the corresponding point on the high-precision reference DOM for the t-th checkpoint. This represents the prior distortion intensity index at the t-th checkpoint;

[0085] For radiation quality assessment, traditional information entropy cannot distinguish between clear noise and clear edges. This invention introduces a fidelity index to quantify the degree of preservation of image texture details. For any pixel location in the image, the specific formula for the phase consistency value is as follows:

[0086]

[0087] in, This represents the phase consistency value at any pixel position x in the image. This represents the direction index of the Log-Gabor filter in the frequency domain. This represents the local amplitude of the pixel under filtering in the nth direction, and signifies the energy intensity in that direction. This represents the phase deviation in the nth direction, reflecting the degree of agreement between the actual phase and the expected symmetrical phase. To indicate a positive number divisible by zero, take 10. -4 This is used to avoid the denominator being zero. The positive part of the phase deviation ranges from 0 to 1. It takes 1 when the actual phase is exactly the same as the expected phase and 0 when they are orthogonal, thus ensuring that the phase consistency result has a clear range.

[0088] The specific formula for phase deviation is:

[0089]

[0090] in, This represents the actual phase in the nth direction. This represents the desired phase in the nth direction;

[0091] The specific formula for calculating the fidelity index of the effective region in the entire image is as follows:

[0092]

[0093] in, This represents the fidelity index; the closer the value is to 1, the clearer and sharper the structural features such as edges and corners in the image. This represents the total number of valid pixels in the image after deducting invalid areas such as clouds, cloud shadows, water bodies, and ice and snow.

[0094] The semantic availability dimension uses a deep learning semantic segmentation model to extract invalid coverage areas such as clouds, cloud shadows, water bodies, and ice and snow in the image, and calculates the effective coverage rate, which is the ratio of the total number of effective pixels after deducting invalid areas such as clouds, cloud shadows, water bodies, and ice and snow to the total number of pixels in the image.

[0095] The evaluation results of the three dimensions of geometric accuracy, radiometric quality, and semantic usability are mapped to the basic probability assignments of three independent evidence bodies using Gaussian membership functions, resulting in... , , Recognition framework The basic probability assignment for each piece of evidence is the sum of the probability assignments at the four levels, which is 1.

[0096] The Gaussian membership function maps each index value to a basic probability assignment of four quality levels. The Gaussian membership parameters for each level are shown in Table 1.

[0097] Table 1

[0098]

[0099] For a certain quality level, after calculating the Gaussian membership value of that level, the Gaussian membership values ​​of the current evidence body at the four levels are summed respectively as the normalized denominator. The probability assignment for each level is the Gaussian membership value of that level divided by the normalized denominator. After processing, the sum of the probability assignments of the four levels of the same evidence body is 1, which satisfies the basic attribute requirements of the basic probability assignment.

[0100] The conflict coefficient between the three pieces of evidence is calculated using the following formula:

[0101]

[0102] in, This represents the conflict coefficient among the three pieces of evidence; a larger value indicates a greater discrepancy in the judgments based on the three pieces of evidence. Representing geometric evidence The supported quality levels belong to a subset of Indicating radiation evidence The supported quality levels belong to a subset of Semantic evidence The supported quality levels belong to A subset of the three pieces of evidence, when the intersection is an empty set, indicates that the judgment results of the three pieces of evidence have no common level, and no single level is supported by all three pieces of evidence simultaneously. This represents the basic probability assignment for a geometric evidence determination result of grade A. This represents the basic probability assignment indicating that the radiation evidence assessment result is Level B. This represents the basic probability assignment for a semantic evidence determination result of level C.

[0103] The overall trust level after fusion is calculated using the conflict allocation combination rule. The specific formula is as follows:

[0104]

[0105] in, Indicates overall trust level. Let R represent a specific quality level to be determined after fusion, belonging to one of {Excellent, Good, Satisfactory, Unsatisfactory}. The summation of the numerators iterates through all combinations satisfying A∩B∩C=R, meaning the intersection of the three pieces of evidence equals R. This indicates that the three pieces of evidence jointly support a quality level R. The sum of the probability products of these combinations is then calculated. The conflict allocation coefficient, ranging from 0.2 to 0.4, is used to allocate a portion of the probability to the unknown state when there is a high degree of conflict, preventing blind trust in conflicting evidence. When the conflict coefficient is 0, all evidence is completely consistent, the denominator is 1, and the fusion trust level degenerates into the sum of the trust levels supported by all evidence. When the conflict coefficient is 1, all evidence is completely conflicting, the numerator is 0, and therefore the fusion trust level for each quality level is 0. The excess probability is allocated to the unknown state, which represents a situation where the quality level cannot be reliably determined. The probability is equal to 1 minus the sum of the fusion trust levels of the four known levels. When the conflict coefficient is large, the judgment results of the evidence are highly divergent, and the probability of the unknown state also increases accordingly.

[0106] Finally, the maximum value of the overall trust score among the four quality levels after merging is calculated. If the maximum value of the overall trust score is greater than 0.6, the corresponding level is selected as the final quality level; if it is less than or equal to 0.6 or the conflict coefficient is lower, the remaining level is selected. If the value is greater than 0.8, it is considered suspicious / requires manual review.

[0107] Based on the assessed quality level, images are processed. For images with a quality level of unqualified, questionable, or requiring manual review, the point-by-point geometric deviation is calculated for each independent high-precision checkpoint t. When the error exceeds twice the allowable mean square error of the corresponding terrain category, the area where the checkpoint is located is marked as an error-over-limit area. The allowable error for flat and hilly areas is usually 2.5m, and the allowable error for mountainous and high mountain areas is usually 3.75m. The area where the corresponding checkpoint is located is marked as an error-over-limit area, triggering a re-correction process. When step S2 is re-executed, the weight coefficient of the smoothing constraint term in this area is increased, and the regularization penalty coefficient is increased to 2 to 3 times the original value to suppress local excessive distortion.

[0108] For images with a quality level of qualified or above, cropping is performed according to standard map sheet division rules. Images with a resolution of 0.5 meters or better are cropped using the 1:25000 standard map sheet division rule, that is, the images are divided according to the map sheet division rules of a 1:25000 scale topographic map. Images with a resolution of 1 meter or better but less than 0.5 meters are cropped using the 1:50000 standard map sheet division rule, that is, the images are divided according to the map sheet division rules of a 1:50000 scale topographic map. After cropping, an MDB format metadata file is generated, and the final output is packaged and exported.

[0109] like Figure 2 As shown, an automated correction and quality inspection system for multi-source satellite remote sensing images is used to implement any of the methods for automated correction and quality inspection of multi-source satellite remote sensing images, comprising:

[0110] The data acquisition module is used to acquire remote sensing image data and its associated parameters, perform format parsing and specification verification on the remote sensing image data, and extract image parameters.

[0111] The image correction module is used to construct the original geometric positioning model of the image based on the rational polynomial coefficient model, acquire corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, and then construct the prior distortion intensity index and spatial heterogeneity weight matrix. The affine transformation correction parameters are solved through regularized weighted adjustment to complete the orthorectification of panchromatic and multispectral images respectively.

[0112] The fusion processing module is used to fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image.

[0113] The quality inspection module is used to perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension and semantic usability dimension.

[0114] The comprehensive scoring module is used to construct the evaluation results of geometric accuracy, radiometric quality and semantic usability into basic probability assignments for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image.

[0115] The feedback optimization module is used to mark the error-exceeding areas of unqualified images and trigger re-correction processing based on the quality level.

[0116] The output module is used to crop qualified and above images according to standard framing rules, generate metadata, and output the results.

[0117] Data transfer between modules follows this process: After data acquisition, the data acquisition module transmits standardized image data, RPC parameters, DEM data, and reference DOM to the image correction module; after orthorectification, the image correction module transmits the corrected panchromatic and multispectral images to the fusion processing module; the fusion processing module generates a fused image and transmits it to the quality inspection module; after quality evaluation across three dimensions, the quality inspection module transmits the geometric deviation index, fidelity index, and effective coverage values ​​to the comprehensive scoring module; after evidence fusion decision-making, the comprehensive scoring module transmits the quality level judgment results to the feedback optimization module and the output module respectively; when re-correction is triggered, the feedback optimization module sends the marked error-exceeding area information back to the image correction module, forming a closed loop; the output module completes the final output after receiving a qualified or higher judgment result. Modules communicate via standardized data interfaces, and logs are automatically recorded and breakpoint resumption is supported when a single module experiences an execution anomaly.

[0118] like Figure 3 As shown, this invention addresses the road deformation problem in local areas after orthorectification. The left image is the image before correction, and the right image is the image after correction. This invention automatically identifies deformed areas using a priori distortion intensity index. When the priori distortion intensity index at a control point exceeds a certain multiple of the mean, the system automatically reduces the weight of that control point and increases the regularization penalty coefficient for that region, strengthening smoothing constraints to eliminate deformation. A comparison before and after correction is provided. Figure 3 As shown, road deformation was effectively eliminated and the geometric integrity of surrounding features was maintained.

[0119] The examples described herein are merely preferred embodiments of the invention and are not intended to limit the concept and scope of the invention. Any modifications and improvements made by those skilled in the art to the technical solutions of the invention without departing from the design concept of the invention should fall within the protection scope of the invention.

Claims

1. A method for automated correction and quality inspection of multi-source satellite remote sensing images, characterized in that, Includes the following steps: Step S1: Acquire remote sensing image data and its associated parameters, and simultaneously acquire digital elevation model data and reference orthophoto data. Perform format parsing and specification verification on the remote sensing image data, and extract image parameters. Step S2: Construct the original geometric positioning model of the image based on the rational polynomial coefficient model, collect corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, construct the prior distortion intensity index and spatial heterogeneity weight matrix based on the terrain feature parameters, solve the affine transformation correction parameters through regularized weighted adjustment, and complete the orthorectification of panchromatic and multispectral images respectively. Step S3: Fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image; Step S4: Perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension, and semantic usability dimension. Step S5: Construct the evaluation results of geometric accuracy dimension, radiometric quality dimension and semantic usability dimension into basic probability values ​​for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image. Step S6: Based on the quality level, mark the error exceeding the limit area for unqualified images and trigger re-correction processing; for qualified and above images, perform cropping, metadata generation and output according to standard framing rules.

2. The method according to claim 1, characterized in that, In step S1, the remote sensing image data is multi-source remote sensing image data from different satellite platforms, including satellite image data with different sensor types, spatial resolutions, imaging phases, and coordinate systems; the image parameters include sensor type, spatial resolution, number of bands, imaging phase, satellite side tilt angle, and flight azimuth angle; the auxiliary parameters include at least rational polynomial coefficient parameters, metadata, quick view, and combined chart.

3. The method according to claim 2, characterized in that, Step S2 specifically includes the following steps: Step S2.1: Construct the original geometric positioning model of the image based on the rational polynomial coefficient model, and use the reference orthophoto to perform image matching and automatically collect corresponding control points. Step S2.2: Calculate the terrain feature parameters of each control point based on the digital elevation model, and construct the prior distortion intensity index based on the terrain feature parameters, satellite side sway angle and satellite flight azimuth angle; Step S2.3: Construct a spatial heterogeneity weight matrix based on the prior distortion intensity index, and obtain the affine transformation correction parameters of the rational polynomial coefficient model by solving the regularized weighted least squares objective function.

4. The method according to claim 3, characterized in that, In step S2.2, the terrain feature parameters include slope angle and aspect angle; the method for constructing the prior distortion intensity index is as follows: calculate the square cosine of the angle between the aspect angle and the satellite flight azimuth angle, and the prior distortion intensity index is equal to the tangent of the satellite side angle multiplied by the slope magnitude and then multiplied by the product of a directional coupling coefficient and the square cosine of the angle.

5. The method according to claim 4, characterized in that, In step S2.3, the adaptive weight of the kth control point in the spatial heterogeneity weight matrix is ​​calculated as follows: first, the ratio of the prior distortion intensity index to the arithmetic mean of the prior distortion intensity indices of all control points is calculated, the negative of this ratio is taken as the exponent of the natural constant e, and then multiplied by the ratio of the control point matching correlation coefficient to the preset confidence threshold. The regularized weighted least squares objective function includes two terms: the first term is the sum of the products of the adaptive weights of all control points and the squares of their image-side residual vectors; the second term is the sum of the regularization penalty coefficient multiplied by the squares of the image-side coordinate correction difference vectors of all adjacent control point pairs; the image-side residual vector is the difference between the actual measured image-side coordinates and the theoretical image-side coordinates after six-parameter affine correction.

6. The method according to claim 5, characterized in that, Step S3 specifically includes the following steps: Step S3.1: Perform histogram matching on the panchromatic image and the multispectral image to ensure that the grayscale distribution of the two is consistent. Step S3.2: Convert the multispectral image from RGB space to IHS space and extract the luminance component; Step S3.3: Replace the luminance component with the histogram-matched panchromatic image to obtain a new luminance component; Step S3.4: Perform an inverse IHS transformation on the new luminance component and the original chrominance and saturation components to transform them back to RGB space and generate a fused digital orthophoto.

7. The method according to claim 6, characterized in that, Step S4 specifically includes the following steps: Step S4.1: Evaluate the geometric accuracy dimension based on the prior distortion weighted geometric deviation index. The geometric deviation index is calculated as follows: For each independent high-precision check point, the weight is the prior distortion intensity index of the check point divided by the reciprocal of the arithmetic mean of the prior distortion intensity indices of all control points plus one. The weight is then multiplied by the planar distance between the measured coordinates of the check point on the corrected image and the true coordinates on the high-precision reference orthophoto. Finally, the sum of the weighted products of all check points is divided by the sum of the weights of all check points. Step S4.2: Evaluate the radiometric quality dimension based on the fidelity index. The fidelity index is calculated as follows: For each pixel position in the image, a multi-directional Log-Gabor filter is used to calculate the local amplitude and phase deviation in each direction. The phase consistency value of the pixel is equal to the sum of the products of the local amplitude in each direction and the positive part of the phase deviation in that direction, divided by the sum of the local amplitudes in all directions plus a very small positive number; the fidelity index is equal to the arithmetic mean of the phase consistency values ​​of all effective pixels in the image. Step S4.3: Extract invalid coverage areas from the image using a semantic segmentation model and calculate the effective coverage rate as an evaluation dimension of semantic availability. The invalid coverage areas include at least clouds, cloud shadows, water bodies, and ice and snow. The effective coverage rate is equal to the total number of effective pixels divided by the total number of pixels in the image.

8. The method according to claim 7, characterized in that, Step S5 specifically includes the following steps: Step S5.1: Map the evaluation results of the geometric precision dimension, radiation quality dimension, and semantic usability dimension to the basic probability assignments of three independent evidence bodies, namely geometric evidence, radiation evidence, and semantic evidence, respectively, through Gaussian membership functions. Step S5.2: Calculate the conflict coefficient between the three pieces of evidence. The method for calculating the conflict coefficient is as follows: traverse the quality levels supported by each of the geometric evidence, radiation evidence and semantic evidence, and accumulate the product of the three basic probability assignments corresponding to each combination where the intersection of the three evidence judgment results is empty. Step S5.3: Calculate the overall trust level after fusion using the conflict allocation combination rule, and output the quality level of the image based on the overall trust level.

9. The method according to claim 8, characterized in that, In step S5.3, the conflict allocation combination rule is as follows: for each possible quality level, the fusion trust level is equal to the sum of the products of the basic probability assignments of all combinations of the three evidence judgment results, which are exactly the intersection of the three evidence judgment results and the quality level; the denominator is equal to one minus the conflict coefficient plus the product of the conflict allocation coefficient and the conflict coefficient; if the trust level corresponding to the quality level with the highest fusion trust level is greater than 0.6, then the quality level is selected as the final quality level; otherwise, it is judged as suspicious and awaits manual review.

10. An automated correction and quality inspection system for multi-source satellite remote sensing images, used to implement the automated correction and quality inspection method for multi-source satellite remote sensing images as described in any one of claims 1-9, characterized in that, include: The data acquisition module is used to acquire remote sensing image data and its associated parameters, perform format parsing and specification verification on the remote sensing image data, and extract image parameters. The image correction module is used to construct the original geometric positioning model of the image based on the rational polynomial coefficient model, acquire corresponding control points through image matching, calculate the terrain feature parameters of each control point using the digital elevation model, and then construct the prior distortion intensity index and spatial heterogeneity weight matrix. The affine transformation correction parameters are solved through regularized weighted adjustment to complete the orthorectification of panchromatic and multispectral images respectively. The fusion processing module is used to fuse the orthorectified panchromatic image with the multispectral image to generate a fused digital orthorectified image. The quality inspection module is used to perform multi-dimensional quality inspection on the fused digital orthophoto, including geometric accuracy dimension, radiometric quality dimension and semantic usability dimension. The comprehensive scoring module is used to construct the evaluation results of geometric accuracy, radiometric quality and semantic usability into basic probability assignments for three independent evidence bodies, calculate the conflict coefficient between evidence and make fusion decisions, and output the quality level of the image. The feedback optimization module is used to mark the error-exceeding areas of unqualified images and trigger re-correction processing based on the quality level. The output module is used to crop qualified and above images according to standard framing rules, generate metadata, and output the results.

Images