Satellite image geometric correction method based on low-frequency resampling

Abstract

This application proposes a method for geometric correction of satellite images based on low-frequency resampling. The method includes: obtaining the basic data of an existing L1-level synthetic aperture radar (SAR) image product, calculating the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data; calculating the rotation angle of the observation area of ​​the L1-level SAR image within the minimum bounding rectangle based on the UTM projection range, and dividing the observation beam range into equal parts based on the basic data, calculating the corresponding scaling ratio for each equal-division interval; performing low-frequency resampling on the pixels between each equal-division interval according to the corresponding scaling ratio, and rotating the resampled image according to the rotation angle to obtain a geometrically corrected L2-level SAR image product. This method greatly reduces the time required for SAR image geometric correction, reduces computational complexity, and improves the convenience and applicability of satellite image geometric correction.

Description

Technical Field

[0001] The present application relates to the technical field of satellite image processing, and in particular to a satellite image geometric correction method based on low-frequency resampling. Background Art

[0002] With the rapid development of aerospace technology, the launch cost of small satellites has been reduced year by year, and a large number of small satellites have been successfully equipped with artificial intelligence (AI) real-time processors, which can realize the analysis, imaging, radiation correction, geometric correction, target detection and ground object classification of Synthetic Aperture Radar (SAR) satellite in-orbit data.

[0003] In related technologies, traditional SAR image geometric correction methods usually require the establishment of rational polynomial coefficient models (RPCs) based on anchored control points and digital elevation models (DEMs) on the ground for correction. Although this type of correction method can correct the range image distortion (for example, compression or stretching) of L1-level SAR images, the computational complexity of this type of correction method is high, the correction time is long, and the need for preset ground control points consumes a lot of computing resources. Therefore, the above-mentioned SAR image geometric correction method is not suitable for deployment on real-time AI processors on in-orbit SAR satellites with limited computing resources. Summary of the Invention

[0004] The present application aims to solve one of the technical problems in the related art at least to a certain extent.

[0005] To this end, the first purpose of this application is to propose a satellite image geometric correction method based on low-frequency resampling, which can avoid the calculation of RPC parameters, quickly implement SAR image geometric correction, greatly reduce the time required for SAR image geometric correction, reduce the computational complexity, and improve the convenience and applicability of satellite image geometric correction.

[0006] The second object of this application is to propose a satellite image geometric correction system based on low-frequency resampling;

[0007] A third object of the present application is to provide a non-transitory computer-readable storage medium.

[0008] To achieve the above objectives, a first embodiment of the present application is to propose a satellite image geometric correction method based on low-frequency resampling, the method comprising the following steps:

[0009] Obtaining basic data of an existing L1-level synthetic aperture radar (SAR) image product, and calculating a Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data;

[0010] Calculating the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle according to the UTM projection range, and dividing the observation beam range into equal parts based on the basic data, and calculating the scaling ratio corresponding to each equal-division interval;

[0011] The pixel points between each of the equal partitions are low-frequency resampled according to the corresponding scaling ratio, and the resampled image is rotated according to the rotation angle to obtain a geometrically corrected L2-level SAR image product.

[0012] Optionally, in one embodiment of the present application, the calculating the Universal Transverse Mercator UTM projection range corresponding to the L1-level SAR image based on the basic data includes: extracting the longitude and latitude coordinates of the four corner points of the observation area from the basic data, and calculating the longitude and latitude range of the minimum enclosing rectangle of the observation area based on the longitude and latitude coordinates of the four corner points; mapping the longitude and latitude range of the minimum enclosing rectangle to the UTM projection coordinates of the four boundary points of the minimum enclosing rectangle; and calculating the actual ground length in the distance direction and the azimuth direction respectively based on the UTM projection coordinates of the four boundary points.

[0013] Optionally, in one embodiment of the present application, the calculation of the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle includes: determining the angle between the satellite orbit corresponding to the L1-level SAR image and the prime meridian; when the angle is positive, calculating a first difference based on the UTM projection coordinates of the left boundary point and the UTM projection coordinates of the lower boundary point, and calculating the rotation angle based on the first difference and the actual ground length in the distance direction; when the angle is negative, calculating a second difference based on the UTM projection coordinates of the lower boundary point and the UTM projection coordinates of the right boundary point, and calculating the rotation angle based on the second difference and the actual ground length in the distance direction.

[0014] Optionally, in one embodiment of the present application, the observation beam range is divided into equal parts based on the basic data, and the scaling ratio corresponding to each equal-division interval is calculated, including: calculating the observation beam range according to the near-point observation incident angle and the far-point observation incident angle extracted from the basic data, and dividing the observation beam range into equal parts according to a preset equal-division value; calculating the ground distance set corresponding to each equal-division interval in combination with the SAR satellite earth observation model and the basic data; calculating the ratio set of the distance to the total pixel corresponding to each equal-division interval based on the ground distance set; calculating the ground distance of the central interval between each equal-division interval, and calculating the scaling ratio set corresponding to each equal-division interval based on the ground distance set and the ground distance of the central interval.

[0015] Optionally, in one embodiment of the present application, the SAR satellite earth observation model and the basic data are combined to calculate the ground distance set corresponding to each equally divided interval, including: calculating the perigee electromagnetic wave transmission distance and the apogee electromagnetic wave transmission distance in the SAR satellite earth observation model based on the basic data and the angle between each equally divided interval; calculating the ground distance corresponding to the current equally divided interval based on the perigee electromagnetic wave transmission distance, the apogee electromagnetic wave transmission distance, the angle between the current equally divided interval, and the angle of the previous interval adjacent to the current equally divided interval.

[0016] Optionally, in one embodiment of the present application, the calculation of the ratio set of distances to total pixels corresponding to each equally divided interval based on the ground distance set includes: calculating the sum of all ground distances in the ground distance set; dividing the ground distance corresponding to each of the equally divided intervals by the sum to obtain the total distance pixels corresponding to each of the equally divided intervals; the calculation of the scaling ratio set corresponding to each of the equally divided intervals includes: dividing the ground distance corresponding to each of the equally divided intervals by the ground distance of the central interval to obtain the scaling ratio corresponding to each of the equally divided intervals.

[0017] Optionally, in one embodiment of the present application, the pixel points between each of the equal partitions are low-frequency resampled according to the corresponding scaling ratio, including: calculating the pixel point set of the equal partition interval based on the ratio set of the distance direction total pixels, the number of azimuth direction points between each equal partition, and the number of distance direction pixel points between each equal partition; combining the pixel point set of the equal partition interval and the scaling ratio set, stretching the pixels of the near point interval and compressing the pixels of the far point interval.

[0018] To achieve the above objectives, the second embodiment of the present application further proposes a satellite image geometric correction system based on low-frequency resampling, comprising the following modules:

[0019] A calculation module is used to obtain basic data of existing L1-level synthetic aperture radar (SAR) image products and calculate the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data;

[0020] an equal division module, configured to calculate the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle according to the UTM projection range, and to divide the observation beam range into equal parts based on the basic data, and to calculate the scaling ratio corresponding to each equal division interval;

[0021] The correction module is used to perform low-frequency resampling on the pixel points between each of the equal partitions according to the corresponding scaling ratio, and rotate the resampled image according to the rotation angle to obtain a geometrically corrected L2-level SAR image product.

[0022] In order to implement the above-mentioned embodiments, the third aspect embodiment of the present application further proposes a non-temporary computer-readable storage medium on which a computer program is stored. When the computer program is executed by a processor, the satellite image geometric correction method based on low-frequency resampling in the above-mentioned first aspect embodiment is implemented.

[0023] The technical solution provided by the embodiments of the present application brings at least the following beneficial effects: the present application first calculates the UTM projection range corresponding to the current image based on the basic information in the L1-level SAR image product; then calculates the rotation angle of the observation area in the minimum circumscribed rectangle based on the UTM projection range; then divides the interval of the effective observation beam range into N equal parts based on the basic information, and calculates the scaling ratio of the observation pixels corresponding to each equal-divided interval compared to the central observation point; finally, resamples the pixels between each equal-divided interval according to the scaling ratio, rotates the resampled image according to the rotation angle, and obtains the corrected output result. Thus, the present application can avoid steps such as calculating RPC parameters and quickly implement SAR image geometric correction. On the basis of ensuring the accuracy of SAR image geometric correction, it greatly reduces the time required for SAR image geometric correction, reduces the computational complexity of geometric correction and the consumed computing resources, and can be applied to the on-orbit data processing tasks of SAR satellite AI real-time processors, thereby improving the convenience and applicability of satellite image geometric correction.

[0024] Additional aspects and advantages of the present invention will be set forth in part in the description which follows and, in part, will be obvious from the description which follows, or may be learned through practice of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0025] The above and / or additional aspects and advantages of the present application will become apparent and easily understood from the following description of the embodiments in conjunction with the accompanying drawings, in which:

[0026] Figure 1This is a flowchart of a satellite image geometric correction method based on low-frequency resampling proposed in an embodiment of the present application;

[0027] Figure 2 A flowchart of a method for calculating a UTM projection range proposed in an embodiment of the present application;

[0028] Figure 3 A flowchart of a method for equally dividing a beam range and calculating a scaling ratio proposed in an embodiment of the present application;

[0029] Figure 4 A schematic diagram of a SAR satellite earth observation model proposed in an embodiment of the present application;

[0030] Figure 5 A schematic diagram of a geometric correction process for an L1-level SAR image proposed in an embodiment of the present application;

[0031] Figure 6 A comparative schematic diagram of geometric correction of satellite L1 data proposed in an embodiment of the present application;

[0032] Figure 7 A comparative schematic diagram of another method of geometrically correcting satellite L1 data proposed in an embodiment of the present application;

[0033] Figure 8 This is a structural diagram of a satellite image geometric correction system based on low-frequency resampling proposed in an embodiment of the present application. DETAILED DESCRIPTION

[0034] The following describes embodiments of the present invention in detail, examples of which are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are intended to be used to explain the present invention, and are not to be construed as limiting the present invention.

[0035] It should be noted that this application proposes a satellite image geometric correction method based on low-frequency resampling for the on-orbit geometric correction task of SAR images. This task is an important part of the on-orbit data processing of SAR satellites. Under normal circumstances, the on-orbit data processing process of SAR satellites includes the following steps: 1) The ground transmitting station systematically plans the SAR satellite mission and uploads the mission instructions; 2) After the SAR satellite moves to the mission area, the SAR satellite payload is turned on and the target area is observed; 3) The SAR satellite antenna receives the echo data and stores the original echo data in the AI ​​real-time processor of the SAR satellite through the data transmission module; 4) The AI ​​real-time processor performs the data analysis task and obtains the L0-level data analysis result based on the original echo data; 5) The AI ​​real-time processor executes the on-orbit imaging instruction and generates the L1-level SAR image product (without geometric correction) based on the L0-level data; 6) The AI ​​real-time processor executes the on-orbit radiation correction and geometric correction instructions and generates the L2-level corrected SAR image based on the L1-level imaging product; 7) The AI ​​real-time processor performs downstream tasks such as target detection and ground object classification and generates the L3-level target-level product based on the L2-level corrected SAR image. This application mainly focuses on step 6 in the above-mentioned SAR satellite on-orbit data processing process.

[0036] However, the satellite image geometric correction methods in related technologies have high computational complexity, take a long time to correct, and require preset ground control points, which consumes high computing resources and is not suitable for deployment on on-orbit SAR satellite real-time AI processors with limited computing resources. To this end, this application proposes a satellite image geometric correction method based on low-frequency resampling. This method is based on a variety of basic information in L1-level SAR imaging products, and implements the strategy of maintaining azimuth, stretching the distance to the near point, and compressing the distance to the far point, thereby realizing rapid geometric correction of L1-level SAR images. This method is suitable for on-orbit data processing tasks of SAR satellite AI real-time processors, and greatly reduces the geometric correction time of SAR images without significantly reducing the geometric correction accuracy.

[0037] A satellite image geometric correction method and system based on low-frequency resampling according to an embodiment of the present invention will be described below with reference to the accompanying drawings.

[0038] Figure 1 This is a flow chart of a satellite image geometric correction method based on low-frequency resampling proposed in an embodiment of the present application, such as Figure 1 As shown, the method includes the following steps:

[0039] Step S101: obtain basic data of an existing L1-level synthetic aperture radar (SAR) image product, and calculate the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data.

[0040] Specifically, this application first obtains the uncorrected L1-level SAR image product generated in step 5 of the SAR satellite on-orbit data processing flow, determines the basic information of the L1-level SAR imaging product, and then calculates the Universal Transverse Mercator (UTM) projection range corresponding to the current image based on this basic information.

[0041] Among them, the basic data of the acquired L1-level SAR image products contain various types of data, such as, but not limited to: the latitude and longitude range of the target area of ​​satellite image acquisition, the center point resolution in the range and azimuth directions, the number of sampling points in the range direction before geometric correction, the near-point observation incident angle, the far-point observation incident angle, the data width, the satellite orbit altitude, and the average radius of the earth, etc.

[0042] In order to more clearly illustrate the specific implementation process of calculating the UTM projection range of this application, a calculation method proposed in an embodiment of this application is exemplified below. Figure 2 This is a flow chart of a method for calculating the UTM projection range proposed in an embodiment of the present application, such as Figure 2 As shown, the method includes the following steps:

[0043] Step S201 : extracting the longitude and latitude coordinates of the four corner points of the observation area from the basic data, and calculating the longitude and latitude range of the minimum circumscribed rectangle of the observation area based on the longitude and latitude coordinates of the four corner points.

[0044] Specifically, the coordinate ranges of the longitude and latitude of the four corner points of the observation area corresponding to the current image are extracted from the basic data of the L1 SAR image. The longitude and latitude coordinate values ​​of the four corner points: top left (TL), bottom left (BL), top right (TR), and bottom right (BR) are respectively recorded as (TL lon ,TL lat ),(BL lon ,BL lat ),(TR lon ,TR lat ),(BR lon ,BR lat ).

[0045] Then, the latitude and longitude coordinates of the four corner points are used as input data. Based on the latitude and longitude coordinates of the four corner points (TL, BL, TR, BR), the latitude and longitude range of the minimum bounding rectangle of the current SAR image observation area is calculated. As a possible implementation method, the latitude and longitude range can be calculated using the following formula:

[0046] (llon ,l lat )=min(TL lon ,BL lon ,TR lon ,BR lon )

[0047] (R lon ,R lat )=max(TL lon ,BL lon ,TR lon ,BR lon )

[0048] (T lon ,T lat )=max(TL lat ,BL lat ,TR lat ,BR lat )

[0049] (B lon ,B lat )=min(TL lat ,BL lat ,TR lat ,BR lat )

[0050] Among them, (L lon ,L lat ), (R lon ,R lat ), (T lon ,T lat ), (B lon ,B lat ) represent the longitude and latitude of the left, right, upper and lower boundary points of the minimum circumscribed rectangle of the current observation area, min(·) represents the minimum value of the current longitude / longitude point set and its corresponding longitude / longitude value, and max(·) represents the maximum value of the current longitude / longitude point set and its corresponding longitude / latitude value.

[0051] Step S202 : Mapping the latitude and longitude range of the minimum bounding rectangle into the UTM projection coordinates of the four boundary points of the minimum bounding rectangle.

[0052] Specifically, the latitude and longitude range of the minimum bounding rectangle of the observation area of ​​the current SAR image obtained in the previous step is mapped to UTM projection coordinates. As a possible implementation method, the UTM projection coordinates of the four boundary points can be calculated using the following formula:

[0053] (L utmx ,L utmy )=UTM(L lon ,L lat)

[0054] (R utmx ,R utmy )=UTM(R lon ,R lat )

[0055] (T utmx ,T utmy )=UTM(T lon ,T lat )

[0056] (B utmx ,B utmy )=UTM(B lon ,B kat )

[0057] Among them, UTM(·) means converting the longitude and latitude coordinates into UTM projection coordinates, (L utmx ,L utmy ), (R utmx ,R utmy ), (T utmx ,T utmy ), (B utmx ,B utmy ) represent the UTM projection coordinate values ​​of the boundary points (left, right, top, and bottom) of the minimum bounding rectangle of the current image observation area.

[0058] Step S203 : Calculate the actual ground lengths in the distance direction and the azimuth direction based on the UTM projection coordinates of the four boundary points.

[0059] Specifically, based on the four boundary points (L utmx ,L utmy ), (R utmx ,R utmy ), (T utmx ,T utmy ), (B utmx ,B utmy ) is used to calculate the actual ground lengths in the range direction (RD) and azimuth direction (AD) respectively.

[0060] As a possible implementation, the calculation can be performed in different ways based on the difference in the angle between the satellite orbit corresponding to the image and the prime meridian. In the first example, when the angle between the satellite orbit corresponding to the current image and the prime meridian is positive, the actual ground length in the range direction is calculated using the following formula:

[0061]

[0062] The actual ground length in azimuth is calculated using the following formula:

[0063]

[0064] In the second example, when the angle between the satellite orbit corresponding to the current image and the prime meridian is negative, the actual ground length in the range direction is calculated using the following formula:

[0065]

[0066] The actual ground length in azimuth is calculated using the following formula:

[0067]

[0068] Therefore, the calculation result of the UTM projection range output in this step includes the coordinates of the four boundary points of UTM: (L utmx ,L utmy ), (R utmx ,R utmy ), (T utmx ,T utmy ), (B utmx ,B utmy ), as well as the actual distance length RD and the actual azimuth length AD.

[0069] Step S102 : calculating the rotation angle of the observation area of ​​the L1 SAR image in the minimum bounding rectangle according to the UTM projection range, dividing the observation beam range into equal parts based on the basic data, and calculating the scaling ratio corresponding to each equal division.

[0070] Specifically, the rotation angle of the minimum bounding rectangle of the observation area of ​​the current L1 SAR image is calculated first. The input data of this process are the coordinates of the four UTM boundary points and the range length RD calculated in step S101.

[0071] In one embodiment of the present application, calculating the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle includes: determining the angle between the satellite orbit corresponding to the L1-level SAR image and the prime meridian; when the angle is positive, calculating a first difference based on the UTM projection coordinates of the left boundary point and the UTM projection coordinates of the lower boundary point, and calculating the rotation angle based on the first difference and the actual ground length in the range direction; when the angle is negative, calculating a second difference based on the UTM projection coordinates of the lower boundary point and the UTM projection coordinates of the right boundary point, and calculating the rotation angle based on the second difference and the actual ground length in the range direction.

[0072] Specifically, this embodiment calculates the rotation angle of the actual observation area in the minimum bounding rectangle based on the four boundary points and the distance direction. In the first case, when the angle between the satellite orbit corresponding to the current image and the prime meridian is positive, the rotation angle is calculated using the following formula:

[0073] RA=sin -1 ((L utmy -B utmy ) / RD)

[0074] In the second case, when the angle between the satellite orbit corresponding to the current image and the prime meridian is negative, the rotation angle is calculated using the following formula:

[0075] RA=sin -1 ((B utmy -R utmy ) / RD)

[0076] Among them, RA represents the rotation angle of the actual observation area in the minimum bounding rectangle.

[0077] Furthermore, the effective observation beam range corresponding to the observation area of ​​the L1-level SAR image is divided into N equal intervals, and the scaling ratio of the observation pixel corresponding to each divided interval compared with the central observation point is calculated, where N is the observation beam angle division coefficient.

[0078] It is understood that the present application divides the observation beam range into N equal parts, which is equivalent to dividing the angle between the near-point observation line and the far-point observation line in the SAR satellite Earth observation model into equal parts. Relevant data can be extracted from the basic data and substituted into the SAR satellite Earth observation model for calculation.

[0079] In order to more clearly illustrate the specific implementation process of calculating the UTM projection range of this application, a method proposed in an embodiment of this application is exemplified below. Figure 3 This is a flow chart of a method for equally dividing a beam range and calculating a scaling ratio proposed in an embodiment of the present application, such as Figure 3 As shown, the method includes the following steps:

[0080] Step S301 : Calculate the observation beam range according to the near point observation incident angle and the far point observation incident angle extracted from the basic data, and divide the observation beam range into equal parts according to a preset equal division value.

[0081] The input data in the equal division and calculation method of this embodiment include: the near point incident angle of the current image observation extracted from the basic data and the far-point angle of incidence And the preset observation beam angle division coefficient N.

[0082] Specifically, based on the near point incident angle and the far-point angle of incidence Calculate observation beam range Then divide it into N equal parts to get the set of zoom areas

[0083] As a possible implementation method, the observation beam range is first calculated using the following formula:

[0084]

[0085] The observation beam range is then divided equally using the following formula:

[0086]

[0087] Where k = 1, 2,…, N.

[0088] It should be noted that the elements in the zoom area set are obtained for each equally divided interval. This application can flexibly adjust the geometric correction accuracy by adjusting the value of the observation beam angle division coefficient N in combination with actual factors such as correction accuracy requirements and computing resources. Among them, the larger the value of N (not exceeding the actual number of distance sampling points), the higher the geometric correction accuracy, and the greater the computational complexity. Conversely, the smaller the value of N (not less than 100), the lower the geometric correction accuracy, and the greater the computational complexity.

[0089] Step S302 : Calculate the ground distance set corresponding to each equally divided interval by combining the SAR satellite earth observation model and basic data.

[0090] Among them, the SAR satellite earth observation model proposed in this application is as follows Figure 4 As shown, the observation near point A corresponds to the near point incident angle, and the observation near point B corresponds to the far point incident angle.

[0091] Specifically, based on Figure 4 The SAR satellite earth observation model shown in the figure uses the relevant basic data to substitute into the calculation formula. According to the relative near-point incidence angle and the relative far-point incidence angle in each zoom interval, the actual ground distance corresponding to each interval can be calculated, and then the ground distance set corresponding to each equally divided interval {D1, D2, D3, ..., D N}.

[0092] In this embodiment, the SAR satellite earth observation model and basic data are combined to calculate the ground distance set corresponding to each equally divided interval, including: calculating the perigee electromagnetic wave transmission distance and the apogee electromagnetic wave transmission distance in the SAR satellite earth observation model based on the basic data and the angles between each equally divided interval; and calculating the ground distance corresponding to the current equally divided interval based on the perigee electromagnetic wave transmission distance, the apogee electromagnetic wave transmission distance, the angle between the current equally divided interval, and the angle of the previous interval adjacent to the current equally divided interval.

[0093] Specifically, for the kth equally divided interval, the actual corresponding ground distance of the interval can be calculated by the following formula:

[0094]

[0095] Where k = 1, 2, ..., N. The various parameters in the above formula can be calculated using the following formula:

[0096]

[0097] Among them, SA means Figure 4 The mid-perigee electromagnetic wave transmission distance, SB represents Figure 4 The distance of electromagnetic wave transmission at mid-to-far point, expressed as ∠SOA Figure 1 The midpoint angle, ∠SOB, represents Figure 1 The midpoint angle, RE represents the average radius of the earth, SH represents Figure 4 Medium satellite orbit altitude. Indicates the angle between the current equal partitions, represents the angle between the previous equal partitions,

[0098] It is understandable that for the zoom area set For each equally divided interval in , the corresponding actual ground distance can be calculated by the above method, thus obtaining the ground distance set {D1, D2, D3, ..., D N}.

[0099] Step S303: Calculate the ratio set of the distance to total pixels corresponding to each equally divided interval based on the ground distance set.

[0100] Specifically, based on The corresponding ground distance set {D1,D2,D3,…,D N}, calculate the ratio set of distance to total pixels before geometric correction corresponding to each equal interval {R1, R2, R3, ..., R N}.

[0101] In this embodiment, based on the ground distance set, the ratio set of the total distance pixels corresponding to each equally divided interval is calculated, including: calculating the sum of all ground distances in the ground distance set; dividing the ground distance corresponding to each equally divided interval by the sum to obtain the total distance pixels corresponding to each equally divided interval.

[0102] Specifically, the total distance pixels corresponding to the kth equally divided interval can be calculated by the following formula:

[0103]

[0104] Where k = 1, 2,…, N.

[0105] Step S304 , calculating the ground distances of the center intervals between the equal sections, and calculating the zoom ratio set corresponding to the equal sections based on the ground distance set and the ground distances of the center intervals.

[0106] Specifically, based on the ground distance between the center intervals of each equal partition, the actual corresponding scaling ratio set {Z1, Z2, Z3, ..., Z N}.

[0107] In this embodiment, calculating the set of scaling ratios corresponding to each equally divided interval includes: dividing the ground distance corresponding to each equally divided interval by the ground distance of the center interval to obtain the scaling ratio corresponding to each equally divided interval.

[0108] Specifically, the scaling ratio corresponding to the kth equally divided interval can be calculated by the following formula:

[0109]

[0110] Among them, D mid Indicates the ground distance corresponding to the center interval in a set of equally divided intervals. The center interval can be the interval at the most central position among multiple equally divided intervals.

[0111] Therefore, the output of this step includes: the ground distance set {D1, D2, D3, ..., D N}、The ratio set of distance to total pixels {R1,R2,R3,…,R N} and the scaling ratio set {Z1,Z2,Z3,…,Z N}.

[0112] In step S103 , the pixels between each equal partition are low-frequency resampled according to the corresponding scaling ratio, and the resampled image is rotated according to the rotation angle to obtain a geometrically corrected L2-level SAR image product.

[0113] Specifically, the pixels between each equal partition are subjected to low-frequency resampling according to the corresponding scaling ratio calculated in step S102 to achieve fast geometric correction. The resampling may include operations such as stretching or compressing the pixels.

[0114] The input data of this step include the L1-level SAR image product to be corrected and the ground distance set {D1, D2, D3, ..., D N} and the ratio set of distance to total pixels {R1,R2,R3,…,R N}.

[0115] In one embodiment of the present application, the pixel points between each equal partition are low-frequency resampled according to the corresponding scaling ratio, including: calculating the pixel point set of the equal partition interval based on the ratio set of the total pixels in the distance direction, the number of points in the azimuth direction between each equal partition, and the number of pixel points in the distance direction between each equal partition; combining the pixel point set of the equal partition interval and the scaling ratio set, stretching the pixels in the near point interval and compressing the pixels in the far point interval.

[0116] Specifically, in this embodiment, when performing low-frequency resampling fast geometric correction, in the first step, based on the ratio set {R1, R2, R3, ..., R N}, select the matrix data between each equal partition for distance correction, and obtain the pixel point set of the equal partition interval before geometric correction as {RM1,RM2,RM3,…,RM N}.

[0117] in, P represents the number of azimuth points between the kth equal partitions (this value is fixed), Q k Indicates the distance between the equal partitions in pixels. is a parameter used to represent the dimension, Q k It can be calculated by the following formula:

[0118] Q k =R k ·RDN

[0119] Among them, RDN represents the total number of pixels in the range direction of the L1-level SAR image before geometric correction.

[0120] The second step is based on the equally divided interval pixel set {RM1,RM2,RM3,…,RM N} and the scaling ratio set {Z1,Z2,Z3,…,Z N}, stretch the pixels in the near point interval and compress the pixels in the far point interval. The specific implementation process of stretching or compressing the relevant pixels is as follows: Figure 5 As shown in Figure 2, the geometrically corrected image can be obtained.

[0121] Furthermore, after obtaining the geometrically corrected image, the method further includes: a third step of performing a vertical mirror flip on the corrected image according to the data storage rule, and then rotating the flipped image by the rotation angle RA calculated in step S102, that is, rotating the image by the rotation angle of the actual observation area in the minimum circumscribed rectangle, to obtain I sAR .

[0122] Finally, write the latitude and longitude ranges corresponding to the four corner points of the observation area boundary into I SAR The final L2-level SAR image product IG with coordinate information is obtained from the data SAR .

[0123] In order to more clearly and intuitively describe the advantages of the satellite image geometric correction method based on low-frequency resampling of the present application compared with other correction methods, the following is an example of geometric correction of specific satellite image data in actual application.

[0124] In one embodiment of the present application, taking the L1 level data of the GF3 satellite as an example, the two types of data of the satellite are geometrically corrected by three methods: geometric correction based on RPC in related technologies, geometric correction based on RPC and DEM, and geometric correction based on low-frequency resampling method proposed in this application. The geometric correction results of the two types of data are shown as follows: Figure 6 and Figure 7 shown.

[0125] Specifically, the test was completed by inputting L1A-level data of 24354*24518 pixels (the resolution of the data was 5m, the width was 50km, and dual polarization), and the test results were shown in Table 1.

[0126] Table 1 Comparative analysis of geometric correction methods

[0127]

[0128] From Table 1 and Figure 6 and Figure 7 It can be seen that the correction method of the present application can greatly reduce the time required for geometric correction of SAR images without significantly reducing the geometric correction accuracy and meeting the accuracy requirements of SAR image geometric correction, thereby achieving rapid correction and having low computational complexity, thereby improving the efficiency of geometric correction of SAR images.

[0129] In summary, the satellite image geometric correction method based on low-frequency resampling in the embodiment of the present application first calculates the UTM projection range corresponding to the current image based on the basic information in the L1-level SAR image product; then calculates the rotation angle of the observation area in the minimum circumscribed rectangle based on the UTM projection range; then divides the interval of the effective observation beam range into N equal parts based on the basic information, and calculates the scaling ratio of the observation pixels corresponding to each equal-divided interval compared to the central observation point; finally, resamples the pixels between each equal-divided interval according to the scaling ratio, rotates the resampled image according to the rotation angle, and obtains the corrected output result. Therefore, this method can avoid steps such as calculating RPC parameters, quickly implement SAR image geometric correction, and greatly reduce the time required for SAR image geometric correction while ensuring the accuracy of SAR image geometric correction, reduce the computational complexity of geometric correction and the consumed computing resources, and can be applied to the on-orbit data processing tasks of SAR satellite AI real-time processors, thereby improving the convenience and applicability of satellite image geometric correction.

[0130] In order to implement the above embodiment, the present application also proposes a satellite image geometric correction system based on low-frequency resampling. Figure 8 This is a structural diagram of a satellite image geometric correction system based on low-frequency resampling proposed in an embodiment of the present application, such as Figure 8 As shown, the system includes: a calculation module 100, an equal division module 200 and a correction module 300.

[0131] The calculation module 100 is used to obtain basic data of an existing L1-level synthetic aperture radar (SAR) image product and calculate the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data.

[0132] The equal division module 200 is used to calculate the rotation angle of the observation area of ​​the L1 SAR image in the minimum bounding rectangle according to the UTM projection range, and to divide the observation beam range into equal parts based on the basic data, and calculate the scaling ratio corresponding to each equal division interval.

[0133] The correction module 300 is used to perform low-frequency resampling on the pixels between each equal partition according to the corresponding scaling ratio, and rotate the resampled image according to the rotation angle to obtain a geometrically corrected L2-level SAR image product.

[0134] Optionally, in one embodiment of the present application, the calculation module 100 is specifically used to: extract the latitude and longitude coordinates of the four corner points of the observation area from the basic data, and calculate the latitude and longitude range of the minimum enclosing rectangle of the observation area based on the latitude and longitude coordinates of the four corner points; map the latitude and longitude range of the minimum enclosing rectangle to the UTM projection coordinates of the four boundary points of the minimum enclosing rectangle; and calculate the actual ground length in the distance direction and the azimuth direction based on the UTM projection coordinates of the four boundary points.

[0135] Optionally, in one embodiment of the present application, the equal division module 200 is specifically used to: determine the angle between the satellite orbit corresponding to the L1-level SAR image and the prime meridian; when the angle is positive, calculate a first difference based on the UTM projection coordinates of the left boundary point and the UTM projection coordinates of the lower boundary point, and calculate the rotation angle based on the first difference and the actual ground length in the distance direction; when the angle is negative, calculate a second difference based on the UTM projection coordinates of the lower boundary point and the UTM projection coordinates of the right boundary point, and calculate the rotation angle based on the second difference and the actual ground length in the distance direction.

[0136] Optionally, in one embodiment of the present application, the equal division module 200 is specifically used to: calculate the observation beam range based on the near-point observation incident angle and the far-point observation incident angle extracted from the basic data, and divide the observation beam range into equal parts according to a preset equal division value; calculate the ground distance set corresponding to each equal division interval in combination with the SAR satellite earth observation model and basic data; calculate the ratio set of the distance to the total pixel corresponding to each equal division interval based on the ground distance set; calculate the ground distance of the central interval between each equal division area, and calculate the zoom ratio set corresponding to each equal division interval based on the ground distance set and the ground distance of the central interval.

[0137] Optionally, in one embodiment of the present application, the equal division module 200 is specifically used to: calculate the perigee electromagnetic wave transmission distance and the apogee electromagnetic wave transmission distance in the SAR satellite earth observation model based on basic data and the angle between each equal division; calculate the ground distance corresponding to the current equal division interval according to the perigee electromagnetic wave transmission distance, the apogee electromagnetic wave transmission distance, the angle between the current equal divisions, and the angle of the previous interval adjacent to the current equal division interval.

[0138] Optionally, in one embodiment of the present application, the equal division module 200 is specifically used to: calculate the sum of all ground distances in the ground distance set; divide the ground distance corresponding to each equal division interval by the sum to obtain the total distance pixels corresponding to each equal division interval; divide the ground distance corresponding to each equal division interval by the ground distance of the central interval to obtain the zoom ratio corresponding to each equal division interval.

[0139] Optionally, in one embodiment of the present application, the correction module 300 is specifically used to: calculate the pixel point set of the equally divided interval based on the ratio set of the total pixels in the distance direction, the number of points in the azimuth direction between each equally divided area, and the number of pixel points in the distance direction between each equally divided area; combine the pixel point set of the equally divided interval and the scaling ratio set to stretch the pixels in the near point interval and compress the pixels in the far point interval.

[0140] It should be noted that the above explanation of the embodiment of the satellite image geometric correction method based on low-frequency resampling is also applicable to the system of this embodiment, and will not be repeated here.

[0141] To sum up, the satellite image geometric correction system based on low-frequency resampling in the embodiment of the present application can avoid steps such as calculating RPC parameters and quickly implement SAR image geometric correction. On the basis of ensuring the accuracy of SAR image geometric correction, it greatly reduces the time required for SAR image geometric correction, reduces the computational complexity of geometric correction and the consumed computing resources, and can be applied to the on-orbit data processing tasks of SAR satellite AI real-time processors, thereby improving the convenience and applicability of satellite image geometric correction.

[0142] In order to implement the above embodiments, the present application also proposes a non-temporary computer-readable storage medium on which a computer program is stored. When the computer program is executed by a processor, it implements the satellite image geometric correction method based on low-frequency resampling as described in any of the above embodiments.

[0143] In the description of this specification, the description with reference to the terms "one embodiment", "some embodiments", "example", "specific example", or "some examples" means that the specific features, structures, materials or characteristics described in conjunction with the embodiment or example are included in at least one embodiment or example of the present application. In this specification, the schematic representations of the above terms do not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials or characteristics described can be combined in any one or more embodiments or examples in a suitable manner. In addition, those skilled in the art can combine and combine different embodiments or examples described in this specification and features of different embodiments or examples without contradiction.

[0144] Furthermore, the terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the number of the technical features being referred to. Thus, a feature defined as "first" or "second" may explicitly or implicitly include at least one of such features. Throughout the description of this application, "plurality" means at least two, for example, two, three, etc., unless otherwise specifically defined.

[0145] Any process or method description in a flowchart or otherwise described herein may be understood to represent a module, segment or portion of code comprising one or more executable instructions for implementing the steps of a custom logical function or process, and the scope of the preferred embodiments of the present application includes alternative implementations in which functions may be performed out of the order shown or discussed, including performing functions in a substantially simultaneous manner or in the reverse order depending on the functions involved, which should be understood by those skilled in the art to which the embodiments of the present application belong.

[0146] The logic and / or steps represented in the flowcharts or otherwise described herein, for example, can be considered as a sequenced list of executable instructions for implementing the logical functions, and can be embodied in any computer-readable medium for use by, or in conjunction with, an instruction execution system, apparatus, or device (e.g., a computer-based system, a system including a processor, or other system that can fetch and execute instructions from an instruction execution system, apparatus, or device). For purposes of this specification, a "computer-readable medium" can be any device that can contain, store, communicate, propagate, or transport a program for use by, or in conjunction with, an instruction execution system, apparatus, or device. More specific examples (a non-exhaustive list) of computer-readable media include the following: an electrical connection with one or more wires (electronic devices), a portable computer disk cartridge (magnetic device), random access memory (RAM), read-only memory (ROM), erasable and programmable read-only memory (EPROM or flash memory), fiber optic devices, and a portable compact disc read-only memory (CDROM). Furthermore, the computer-readable medium may even be paper or other suitable medium on which the program is printed, since the program may be obtained electronically, for example, by optically scanning the paper or other medium and then editing, interpreting or processing it in another suitable manner if necessary, and then storing it in a computer memory.

[0147] It should be understood that various parts of the present application can be implemented using hardware, software, firmware, or a combination thereof. In the above embodiments, multiple steps or methods can be implemented using software or firmware stored in a memory and executed by a suitable instruction execution system. For example, if implemented using hardware, as in another embodiment, any one of the following technologies known in the art or a combination thereof can be used to implement: a discrete logic circuit having a logic gate circuit for implementing a logic function on a data signal, an application-specific integrated circuit having a suitable combination of logic gate circuits, a programmable gate array (PGA), a field programmable gate array (FPGA), etc.

[0148] Those skilled in the art will understand that all or part of the steps in the method of the above embodiment can be completed by instructing related hardware through a program, and the program can be stored in a computer-readable storage medium. When the program is executed, it includes one or a combination of the steps of the method embodiment.

[0149] In addition, the functional units in the various embodiments of the present application may be integrated into a processing module, or each unit may exist physically separately, or two or more units may be integrated into a module. The above-mentioned integrated module may be implemented in the form of hardware or in the form of a software functional module. If the integrated module is implemented in the form of a software functional module and sold or used as an independent product, it may also be stored in a computer-readable storage medium.

[0150] The storage medium mentioned above may be a read-only memory, a magnetic disk, or an optical disk, etc. Although the embodiments of the present application have been shown and described above, it is understood that the above embodiments are exemplary and should not be construed as limiting the present application. Persons skilled in the art may make changes, modifications, substitutions, and variations to the above embodiments within the scope of the present application.

Claims

1. A satellite image geometric correction method based on low-frequency resampling, characterized in that: The following steps are involved: Obtaining basic data of an existing L1-level synthetic aperture radar (SAR) image product, and calculating a Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data; Calculate the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle according to the UTM projection range, divide the observation beam range into equal parts based on the basic data, and calculate the scaling ratio corresponding to each equal-division interval, wherein dividing the observation beam range into equal parts based on the basic data and calculating the scaling ratio corresponding to each equal-division interval includes: calculating the observation beam range according to the near-point observation incident angle and the far-point observation incident angle extracted from the basic data, and dividing the observation beam range into equal parts according to a preset equal-division value; calculating the ground distance set corresponding to each equal-division interval in combination with the SAR satellite earth observation model and the basic data; calculating the ratio set of the distance to the total pixels corresponding to each equal-division interval based on the ground distance set; calculating the ground distance of the central interval between each equal-division interval, and calculating the scaling ratio set corresponding to each equal-division interval based on the ground distance set and the ground distance of the central interval; The pixel points between each of the equal partitions are low-frequency resampled according to the corresponding scaling ratio, and the resampled image is rotated according to the rotation angle to obtain a geometrically corrected L2-level SAR image product, wherein the low-frequency resampling of the pixel points between each of the equal partitions according to the corresponding scaling ratio includes: calculating a pixel point set in the equal partition interval based on the ratio set of the total pixels in the range direction, the number of points in the azimuth direction between each equal partition, and the number of pixel points in the range direction between each equal partition; and combining the pixel point set in the equal partition interval and the scaling ratio set, stretching the pixels in the near point interval and compressing the pixels in the far point interval.

2. The method according to claim 1, characterized in that Calculating the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image according to the basic data includes: Extracting the longitude and latitude coordinates of the four corner points of the observation area from the basic data, and calculating the longitude and latitude range of the minimum circumscribed rectangle of the observation area based on the longitude and latitude coordinates of the four corner points; Mapping the latitude and longitude range of the minimum bounding rectangle into the UTM projection coordinates of the four boundary points of the minimum bounding rectangle; Based on the UTM projection coordinates of the four boundary points, the actual ground lengths in the range direction and the azimuth direction are calculated respectively.

3. The method according to claim 2, characterized in that The calculating of the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle includes: Determining the angle between the satellite orbit corresponding to the L1-level SAR image and the prime meridian; When the angle is positive, a first difference is calculated based on the UTM projection coordinates of the left boundary point and the UTM projection coordinates of the lower boundary point, and the rotation angle is calculated based on the first difference and the actual ground length in the distance direction; When the angle is negative, a second difference is calculated based on the UTM projection coordinates of the lower boundary point and the UTM projection coordinates of the right boundary point, and the rotation angle is calculated based on the second difference and the actual ground length in the distance direction.

4. The method according to claim 1, wherein The method of combining the SAR satellite earth observation model and the basic data to calculate the ground distance set corresponding to each equally divided interval includes: Calculating the perigee electromagnetic wave transmission distance and the apogee electromagnetic wave transmission distance in the SAR satellite earth observation model based on the basic data and the angles between the equal partitions; The ground distance corresponding to the current equally divided interval is calculated based on the perigee electromagnetic wave transmission distance, the apogee electromagnetic wave transmission distance, the angle between the current equally divided intervals, and the angle of the previous interval adjacent to the current equally divided interval.

5. The method according to claim 1, wherein Calculating a set of ratios of distances to total pixels corresponding to each equally divided interval based on the ground distance set includes: Calculating the sum of all ground distances in the ground distance set; Divide the ground distance corresponding to each of the equally divided intervals by the sum value to obtain the total distance pixels corresponding to each of the equally divided intervals; The calculating of the set of scaling ratios corresponding to each equally divided interval includes: The ground distance corresponding to each of the equally divided intervals is divided by the ground distance of the central interval to obtain the zoom ratio corresponding to each of the equally divided intervals.

6. A satellite image geometric correction system based on low-frequency resampling, characterized in that: include: A calculation module is used to obtain basic data of existing L1-level synthetic aperture radar (SAR) image products and calculate the Universal Transverse Mercator (UTM) projection range corresponding to the L1-level SAR image based on the basic data; An equal division module is used to calculate the rotation angle of the observation area of ​​the L1-level SAR image in the minimum circumscribed rectangle according to the UTM projection range, and to divide the observation beam range into equal parts based on the basic data, and to calculate the scaling ratio corresponding to each equal division interval, wherein the equal division module is specifically used to: calculate the observation beam range according to the near point observation incident angle and the far point observation incident angle extracted from the basic data, and to divide the observation beam range into equal parts according to a preset equal division value; calculate the ground distance set corresponding to each equal division interval in combination with the SAR satellite earth observation model and the basic data; calculate the ratio set of the distance to the total pixels corresponding to each equal division interval according to the ground distance set; calculate the ground distance of the central interval between each equal division area, and calculate the scaling ratio set corresponding to each equal division interval according to the ground distance set and the ground distance of the central interval; A correction module is configured to perform low-frequency resampling of the pixel points between each of the equal partitions according to the corresponding scaling ratio, and to rotate the resampled image according to the rotation angle to obtain a geometrically corrected L2-level SAR image product, wherein the correction module is specifically configured to: calculate a set of pixel points in the equal partition interval based on the set of ratios of the total pixels in the range direction, the number of points in the azimuth direction between each equal partition, and the number of pixel points in the range direction between each equal partition; and to stretch the pixels in the near point interval and compress the pixels in the far point interval in combination with the set of pixel points in the equal partition interval and the scaling ratio set.

7. The system according to claim 6, characterized in that The computing module is specifically configured to: Extracting the longitude and latitude coordinates of the four corner points of the observation area from the basic data, and calculating the longitude and latitude range of the minimum circumscribed rectangle of the observation area based on the longitude and latitude coordinates of the four corner points; Mapping the latitude and longitude range of the minimum bounding rectangle into the UTM projection coordinates of the four boundary points of the minimum bounding rectangle; The actual ground lengths in the distance direction and the azimuth direction are respectively calculated based on the UTM projection coordinates of the four boundary points.

8. A non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the satellite image geometric correction method based on low-frequency resampling according to any one of claims 1 to 5 is implemented.

Images