Multi-scene SAR (Synthetic Aperture Radar) image uniform color processing method and device considering terrain radiation

Abstract

The invention discloses a multi-scene SAR image uniform color processing method and device considering terrain radiation, and belongs to the technical field of remote sensing image processing, and the method comprises the steps: obtaining an SAR image of a to-be-detected region, carrying out the preprocessing of the obtained SAR image, and extracting a sampling point; the sampling points are substituted into a sinc function model, global optimization is carried out by fitting the overall SAR image of the to-be-detected area, and a global sinc function compensation coefficient is obtained; constructing polynomial fitting for each original SAR image in the to-be-measured region by using the idea of block adjustment, and integrally calculating polynomial model correction parameters of all SAR images; and performing uniform color correction on the original SAR image in combination with the global sinc function compensation coefficient and the polynomial correction coefficient. According to the method, the problems of path error propagation and accumulation, color difference residue in an overlapping region and the like generally existing in an existing color equalization processing method are effectively avoided, and an effective way is provided for uniform color processing of the SAR image.

Description

Technical Field

[0001] The invention relates to a method and device for color-homogenizing multi-scene SAR images taking terrain radiation into consideration, and belongs to the technical field of remote sensing image processing. Background Art

[0002] Remote sensing technology is a technology for obtaining earth information from space and is one of the main sources of earth big data. Because remote sensing technology has the advantages of macroscopicity, current situation, comprehensiveness and accuracy, it can effectively help people gain a deep understanding of the state and situation of earth resources and environment. It is widely used in environmental monitoring, geological mapping, resource investigation and disaster warning, and plays an important role in the sustainable development of environment, economy and society.

[0003] The active detection and range imaging methods of the SAR system give SAR unique advantages, enabling it to detect the earth's surface at night and during the day and in any weather conditions. It has extremely high resolution and can penetrate the surface to a certain extent to obtain sufficient information. SAR imaging technology can effectively fill the gaps of optical and infrared remote sensing, especially in extreme climate conditions and complex environmental conditions. It has become a mature and powerful remote sensing technology used in a variety of applications worldwide. This technology has great development potential in the fields of topographic mapping, disaster emergency response, land and resources survey, and has a very broad application prospect. However, the limitations of the SAR system itself cause great interference between the backscattered signals between pixels, which will cause a certain deviation between the acquired data imaging and the actual situation, mainly showing uneven brightness distribution and inconsistent color difference in the image. Due to this significant difference in SAR images, the image quality is significantly reduced, which makes the image interpretation work complicated and difficult, and the scope of application and prospects are greatly limited. In addition, there are limitations in the design of SAR sensors. The coverage of single-view SAR images is limited. Obtaining large-scale image information of the entire area to be measured usually requires the splicing of multiple SAR images. Ideally, the radiation intensity of the overlapping areas of adjacent images is consistent after radiation calibration. However, due to the influence of factors such as radar wave signal attenuation, antenna pattern changes, and imaging processor gain errors, there are radiation differences between images in certain areas, especially data from different strips, which reduces the quality of SAR images, makes the images difficult to interpret and read, and hinders the subsequent work. For example, when generating digital orthophotos, if the radiation differences are not eliminated, the subsequent mosaic images will have obvious seams, poor visualization, and the quality of digital orthophoto products will be directly affected. In order to eliminate this difference between images, it is necessary to perform radiation consistency processing on SAR images so that the overall brightness of the SAR images is uniform and the contrast is consistent, so as to obtain images with good visual effects, moderate contrast, and clear texture information.

[0004] At present, the color consistency method of optical remote sensing images is mostly used to directly process the radiation consistency of SAR images. In terms of imaging mechanism, the optical system belongs to a passive central projection system, while the SAR system belongs to an active slant range projection imaging system, and the SAR image can clearly reflect the backscattering characteristics of the ground objects to microwaves. Due to the dielectric constant, roughness and other factors of the ground objects, the backscattering characteristics can change accordingly, and then the structure and properties of different ground objects can be understood. The use of optical image color consistency method to process SAR images does not fully consider the radiation principle of SAR, cannot solve the radiation difference problem generated by the SAR image itself, and is difficult to adapt to actual needs. Therefore, according to the characteristics of SAR images themselves, a color uniformity processing method suitable for SAR images is studied, which has important theoretical significance and research value. Summary of the invention

[0005] In order to solve the above problems, the present invention proposes a method and device for color uniformity processing of multi-scene SAR images taking terrain radiation into consideration.

[0006] The technical solution adopted by the present invention to solve the technical problem is:

[0007] In a first aspect, an embodiment of the present invention provides a method for color-homogenizing a multi-scene SAR image taking terrain radiation into consideration, comprising the following steps:

[0008] Step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points;

[0009] Step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the overall SAR image of the test area, and obtaining the global sinc function compensation coefficient;

[0010] Step S3, uses the idea of ​​regional block adjustment to construct polynomial fitting for each original SAR image in the area to be measured, and solves the polynomial model correction parameters of all SAR images;

[0011] Step S4, combining the global sinc function compensation coefficient and the polynomial correction coefficient to perform uniform color correction on the original SAR image.

[0012] As a possible implementation of this embodiment, the step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points, includes:

[0013] The original SAR images are geometrically corrected at the same resolution, all the original SAR images are unified in the same object projection coordinate system, and an overlapping relationship based on the unified geographical location is established;

[0014] In the object coordinate system of the SAR image, a sampling grid is drawn every 50 m in the horizontal and vertical directions;

[0015] At the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, the average brightness value of 2500 pixels is calculated as the brightness value of the sampling point;

[0016] The geometric model of the SAR image is used to transform the position of the grid point to the position on the original SAR image, and the brightness value of the sampling point is used as the brightness value of the original image position, and the incident angle information of the sampling point is extracted:

[0017] sita=(H s -H g ) / (R0+c×m c )

[0018] Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

[0019] As a possible implementation of this embodiment, step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the entire SAR image of the test area, and obtaining the global sinc function compensation coefficient, includes:

[0020] Build the sinc function model:

[0021]

[0022] Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle;

[0023] The error equation is established, the brightness and incident angle information of all original image sampling points in the test area are substituted into the error equation, and the optimal estimates of χ1, χ2, and χ3 are solved using the nonlinear least squares method:

[0024]

[0025] Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point.

[0026] Calculate the sinc fitting value G0 for the image point coordinates of any image in the test area and the incident angle of the sampling point;

[0027] The global sinc function compensation coefficient is obtained by calculating the brightness value after radiation correction:

[0028]

[0029] In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

[0030] As a possible implementation of this embodiment, step S3 uses the idea of ​​regional block adjustment to construct a polynomial fit for each original SAR image in the test area, and solves the polynomial model correction parameters of all SAR images, including:

[0031] According to the principle that the brightness of different original images at the same location should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal, the error equation for each group of sampling points is established:

[0032]

[0033] dV=a0+a1c+a2c 2 +...

[0034] In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the column coordinate of the image, and a0, a1, a2... are the polynomial coefficients of the column coordinate c;

[0035] Select all the images with overlapping relationship between two grid points in turn to construct the error equation group, that is, if there are n overlaps on a grid point, n>=2, then n-1 independent error equations can be constructed, n is the number of overlapping SAR images in the test area;

[0036] Solve the error equations of all sampling points jointly, and then use the least squares principle to calculate the optimal solution of the parameters;

[0037] The normal equations are solved using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

[0038] In a second aspect, an embodiment of the present invention provides a multi-scene SAR image color uniformity processing device taking terrain radiation into consideration, comprising:

[0039] SAR image acquisition module, used to acquire SAR images of the area to be measured, pre-process the acquired SAR images and extract sampling points;

[0040] The global optimization module is used to substitute the sampling points into the sinc function model, perform global optimization by fitting the overall SAR image of the test area, and obtain the global sinc function compensation coefficient;

[0041] The polynomial fitting module is used to construct a polynomial fitting for each original SAR image in the test area using the idea of ​​regional block adjustment, and calculate the polynomial model correction parameters of all SAR images;

[0042] The uniform color correction module is used to perform uniform color correction on the original SAR image by combining the global sinc function compensation coefficient and the polynomial correction coefficient.

[0043] As a possible implementation of this embodiment, the SAR image acquisition module includes:

[0044] The geometric correction module is used to perform geometric correction processing on the original SAR images at the same resolution, unify all the original SAR images under the same object projection coordinate system, and establish an overlapping relationship based on the unified geographical location;

[0045] The sampling grid drawing module is used to draw the sampling grid every 50m in the horizontal and vertical directions in the object coordinate system of the SAR image;

[0046] A brightness value calculation module is used to calculate the average brightness value of 2500 pixels at the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, as the brightness value of the sampling point;

[0047] The incident angle extraction module is used to convert the position of the grid point to the position on the original SAR image using the geometric model of the SAR image, and use the brightness value of the sampling point as the brightness value of the original image position, and extract the incident angle information of the sampling point:

[0048] sita=(H s -H g ) / (R0+c×m c )

[0049] Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

[0050] As a possible implementation of this embodiment, the global optimization module includes:

[0051] Function model building module, used to build the sinc function model:

[0052]

[0053] Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle;

[0054] The best estimate solution module is used to establish the error equation, substitute the brightness and incident angle information of all original image sampling points in the test area into the error equation, and use the nonlinear least squares method to solve the best estimates of χ1, χ2, and χ3:

[0055]

[0056] Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point.

[0057] The sinc fitting module is used to calculate the sinc fitting value G0 of the image point coordinates and the incident angle of the sampling point of any image in the test area;

[0058] The compensation coefficient acquisition module is used to obtain the global sinc function compensation coefficient by calculating the brightness value after radiation correction:

[0059]

[0060] In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

[0061] As a possible implementation of this embodiment, the uniform color correction module includes:

[0062] The error equation establishment module is used to establish the error equation for each group of sampling points based on the principle that the brightness of different original images at the same position should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal:

[0063]

[0064] dV=a0+a1c+a2c 2

[0065] In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the column coordinate of the image, and a0, a1, and a2 are the polynomial coefficients of the column coordinate c;

[0066] The error equation group establishment module is used to sequentially select images with two-to-two overlapping relationships on all grid points to construct the error equation group, that is, if there are n-fold coverage on a certain grid point, n>=2, then n-1 independent error equations can be constructed, and n is the number of overlapping relationships of SAR images in the test area;

[0067] The joint solution module is used to jointly solve the error equations of all sampling points and then calculate the optimal solution of the parameters using the least squares principle;

[0068] The polynomial coefficient acquisition module is used to solve the normal equation group using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

[0069] The beneficial effects of the technical solution of the embodiment of the present invention are as follows:

[0070] A method for color-homogenizing multi-scene SAR images considering terrain radiation in the technical solution of an embodiment of the present invention comprises the following steps: step S1, obtaining a SAR image of a test area, pre-processing the obtained SAR image and extracting sampling points; step S2, substituting the sampling points into a sinc function model, performing global optimization by fitting the overall SAR image of the test area, and obtaining a global sinc function compensation coefficient; step S3, using the idea of ​​regional network adjustment to construct a polynomial fitting for each original SAR image in the test area, and overall solving the polynomial model correction parameters of all SAR images; step S4, combining the global sinc function compensation coefficient and the polynomial correction coefficient to perform color-homogenizing correction on the original SAR image. The present invention optimizes all images as a whole, realizes radiation correction of all images by fitting a sinc radiation correction model, further performs polynomial fitting processing on the overlapping area, and realizes color-homogenizing processing of all images. The present invention uses the sinc function model in the SAR image uniform color processing, takes into account the influence of terrain radiation on the SAR image, and adopts a processing strategy combining global and local, effectively avoiding the common problems of path error propagation and accumulation, chromatic aberration residual in overlapping areas, etc. in the existing color equalization processing methods, and provides an effective way for the uniform color processing of SAR images. BRIEF DESCRIPTION OF THE DRAWINGS

[0071] Figure 1 is a flow chart of a method for color-homogenizing a multi-scene SAR image considering terrain radiation according to an exemplary embodiment;

[0072] Figure 2 is a structural schematic diagram of a multi-scene SAR image color uniformity processing device considering terrain radiation according to an exemplary embodiment; DETAILED DESCRIPTION

[0073] In order to more clearly illustrate the technical features of the solution of the present invention, the present invention is described in detail below through specific implementation methods and in conjunction with the accompanying drawings.

[0074] like Figure 1 As shown, a method for color-homogenizing a multi-scene SAR image considering terrain radiation is provided in an embodiment of the present invention, comprising the following steps:

[0075] Step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points;

[0076] Step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the overall SAR image of the test area, and obtaining the global sinc function compensation coefficient;

[0077] Step S3, uses the idea of ​​regional block adjustment to construct polynomial fitting for each original SAR image in the area to be measured, and overall solves the polynomial model correction parameters of all SAR images;

[0078] Step S4, combining the global sinc function compensation coefficient and the polynomial correction coefficient to perform uniform color correction on the original SAR image.

[0079] As a possible implementation of this embodiment, the step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points, includes:

[0080] The original SAR images are geometrically corrected at the same resolution, all the original SAR images are unified in the same object projection coordinate system, and an overlapping relationship based on the unified geographical location is established;

[0081] In the object coordinate system of the SAR image, a sampling grid is drawn every 50 m in the horizontal and vertical directions;

[0082] At the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, the average brightness value of 2500 pixels is calculated as the brightness value of the sampling point;

[0083] The geometric model of the SAR image is used to transform the position of the grid point to the position on the original SAR image, and the brightness value of the sampling point is used as the brightness value of the original image position, and the incident angle information of the sampling point is extracted:

[0084] sita=(Hs -H g ) / (R0+c×m c )

[0085] Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

[0086] As a possible implementation of this embodiment, step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the entire SAR image of the test area, and obtaining the global sinc function compensation coefficient, includes:

[0087] Build the sinc function model:

[0088]

[0089] Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle;

[0090] The error equation is established, the brightness and incident angle information of all original image sampling points in the test area are substituted into the error equation, and the optimal estimates of χ1, χ2, and χ3 are solved using the nonlinear least squares method:

[0091]

[0092] Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point.

[0093] Calculate the sinc fitting value G0 for the image point coordinates of any image in the test area and the incident angle of the sampling point;

[0094] The global sinc function compensation coefficient is obtained by calculating the brightness value after radiation correction:

[0095]

[0096] In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

[0097] As a possible implementation of this embodiment, step S3 uses the idea of ​​regional block adjustment to construct a polynomial fit for each original SAR image in the test area, and overall solves the polynomial model correction parameters of all SAR images, including:

[0098] According to the principle that the brightness of different original images at the same location should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal, the error equation for each group of sampling points is established:

[0099]

[0100] dV=a0+a1c+a2c 2 +...

[0101] In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the column coordinate of the image, and a0, a1, and a2 are the polynomial coefficients of the column coordinate c;

[0102] Select all the images with pairwise overlapping relationships on all grid points in turn to construct the error equation group. That is, if there are n overlaps on a grid point, n>=2, then n-1 independent error equations can be constructed.

[0103] Solve the error equations of all sampling points jointly, and then use the least squares principle to calculate the optimal solution of the parameters;

[0104] The normal equations are solved using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

[0105] like Figure 2 As shown, an embodiment of the present invention provides a multi-scene SAR image color uniformity processing device considering terrain radiation, comprising:

[0106] SAR image acquisition module, used to acquire SAR images of the area to be measured, pre-process the acquired SAR images and extract sampling points;

[0107] The global optimization module is used to substitute the sampling points into the sinc function model, perform global optimization by fitting the overall SAR image of the test area, and obtain the global sinc function compensation coefficient;

[0108] The polynomial fitting module is used to construct a polynomial fitting for each original SAR image in the test area using the idea of ​​regional block adjustment, and to calculate the polynomial model correction parameters of all SAR images as a whole;

[0109] The uniform color correction module is used to perform uniform color correction on the original SAR image by combining the global sinc function compensation coefficient and the polynomial correction coefficient.

[0110] As a possible implementation of this embodiment, the SAR image acquisition module includes:

[0111] The geometric correction module is used to perform geometric correction processing on the original SAR images at the same resolution, unify all the original SAR images under the same object projection coordinate system, and establish an overlapping relationship based on the unified geographical location;

[0112] The sampling grid drawing module is used to draw the sampling grid every 50m in the horizontal and vertical directions in the object coordinate system of the SAR image;

[0113] A brightness value calculation module is used to calculate the average brightness value of 2500 pixels at the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, as the brightness value of the sampling point;

[0114] The incident angle extraction module is used to convert the position of the grid point to the position on the original SAR image using the geometric model of the SAR image, and use the brightness value of the sampling point as the brightness value of the original image position, and extract the incident angle information of the sampling point:

[0115] sita=(H s -H g ) / (R0+c×m c )

[0116] Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

[0117] As a possible implementation of this embodiment, the global optimization module includes:

[0118] Function model building module, used to build the sinc function model:

[0119]

[0120] Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle;

[0121] The optimal estimation solution module is used to establish the error equation, substitute the brightness and incident angle information of all original image sampling points in the test area into the error equation, and use the nonlinear least squares method to solve the optimal estimation of χ1, χ2, and χ3:

[0122]

[0123] Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point.

[0124] The sinc fitting module is used to calculate the sinc fitting value G0 of the image point coordinates and the incident angle of the sampling point of any image in the test area;

[0125] The compensation coefficient acquisition module is used to obtain the global sinc function compensation coefficient by calculating the brightness value after radiation correction:

[0126]

[0127] In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

[0128] As a possible implementation of this embodiment, the uniform color correction module includes:

[0129] The error equation establishment module is used to establish the error equation for each group of sampling points based on the principle that the brightness of different original images at the same position should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal:

[0130]

[0131] dV=a0+a1c+a2c 2

[0132] In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the distance column coordinate of the image, and a0, a1, and a2 are the polynomial coefficients of the column coordinate c;

[0133] The error equation group establishment module is used to sequentially select images with two-to-two overlapping relationships on all grid points to construct the error equation group. That is, if there are n-fold coverages on a certain grid point, n>=2, then n-1 independent error equations can be constructed;

[0134] The joint solution module is used to jointly solve the error equations of all sampling points and then calculate the optimal solution of the parameters using the least squares principle;

[0135] The polynomial coefficient acquisition module is used to solve the normal equation group using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

[0136] The specific process of the present invention for performing color uniformity processing on multi-scene SAR images taking terrain radiation into consideration is as follows.

[0137] Step 1: Based on the SAR image of the area to be measured, preprocessing and fast noise-resistant extraction of sampling point brightness are performed to eliminate part of the SAR image speckle noise and reduce the amount of calculation.

[0138] Due to the large number of images and pixels in the test area, and the large brightness difference of the original SAR images, a lot of computing resources and time will be consumed in data processing, and the image data is extracted for sampling points. The steps to obtain the sampling points are as follows.

[0139] (1) Initial orthorectification of original images. The original images are independent in their respective image-space coordinate systems, and it is difficult to find corresponding overlapping relationships. Therefore, the original images are first geometrically corrected at the same resolution, and all the original images are unified in the same object-space projection coordinate system. An overlapping relationship based on a unified geographic location is established, which provides convenient conditions for obtaining points of the same name.

[0140] (2) Establish a sampling grid and obtain sampling points. It is planned to select 50m as the sampling interval. In the object coordinate system of the image, draw a sampling grid every 50m in the horizontal and vertical directions; at the intersection center of the sampling grid, use 50×50 pixels as the sampling unit, and calculate the average grayscale value of 2500 pixels as the grayscale value of this sampling point.

[0141] (3) Geometric transformation back to the original SAR image. After sampling is completed, the position of the grid point is transformed to the position on the original image using the geometric model of the image, and the extracted brightness value is used as the brightness value of the original image position. At the same time, the incident angle information of the point is extracted according to the following formula:

[0142] sita=(H s -H g ) / (R0+c×m c )

[0143] Where, sita is the incident angle of the point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

[0144] The original SAR images are orthorectified to obtain the corresponding position relationship between overlapping images. Uniform vertical and horizontal grids are divided in the object coordinate system. The pixel mean value and position information within a certain range are statistically extracted, which can eliminate the speckle noise in some SAR images and reduce the amount of calculation.

[0145] Step 2: Global sinc model optimization: Apply the global sinc function model to the uniform color processing of the SAR image of the entire test area to eliminate the total color difference from the entire test area.

[0146] The global sinc function model is an ideal function model that can fit SAR signals. The sinc function has periodicity and good spectral characteristics, which can well describe the spectral characteristics of SAR signals, thereby realizing the fitting of SAR images. The fitting method based on the sinc function model has the advantages of simplicity, high efficiency, and easy implementation, and can meet the actual needs of SAR image processing and analysis. The fitting model of the sinc function is:

[0147]

[0148] In the formula, χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle. In order to obtain accurate estimates of χ1, χ2, and χ3, an error equation is established and the nonlinear least squares method is used to solve the best estimates of χ1, χ2, and χ3. The error equation is:

[0149]

[0150] In the formula, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area; the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the average of the incident angles of the sampling points. Substitute the brightness and incident angle information of all the original image sampling points of the test area into the above equation, and obtain the best estimated values ​​of χ1, χ2, and χ3 by least squares.

[0151] After obtaining the sinc parameters, the sinc fitting value G0 can be calculated for the image point coordinates of any image in the test area and the incident angle of the point, and the brightness value after radiation correction can be calculated according to the following formula:

[0152]

[0153] In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at the image point. Through the above formula, the sinc function correction value of all images in the test area can be calculated.

[0154] The sinc function model has a good effect on the range antenna pattern fitting and can ensure the non-negativity of the antenna pattern. Substituting the sampling point information into the sinc function model for the overall SAR image fitting of the test area, the compensation parameters can be obtained globally to improve the overall color difference of the test area.

[0155] Step 3: Local polynomial model optimization is performed on the basis of global optimization to further eliminate the grayscale deviation in the overlapping area.

[0156] A set of polynomial fitting models is constructed for each original SAR image based on the sinc function corrected pixel values. The polynomial fitting model is a method of approximating data through a polynomial function. According to the principle that the brightness of different original images at the same location should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal, the error equation of each group of sampling points is established.

[0157]

[0158] dV=a0+a1c+a2c 2

[0159] In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j They are the polynomial fitting models of the i-th and j-th images respectively. c in dV represents the column coordinate of the image, and a0, a1, a2, ... are the polynomial coefficients of the column coordinate c.

[0160] Select all grid points with overlapping images to construct error equations. That is, if a grid point has n-fold coverage (n>=2), n-1 independent error equations can be constructed. Solve the error equations of all sampling points jointly, and then use the least squares principle to calculate the optimal solution of the parameters. The basic idea of ​​the least squares method is to find a set of coefficients that minimizes the square error between the predicted value and the actual value of the polynomial function, satisfying:

[0161]

[0162] In order to avoid the ill-conditioned solution of the equation system, the ridge estimation is used to solve the normal equation system to obtain the polynomial model coefficients of the distance coordinate c of each original image.

[0163] Using the idea of ​​regional block adjustment, a polynomial fitting is constructed for each original image in the measured area, and the polynomial model correction parameters of all images are calculated as a whole to improve the chromatic aberration of different images in the overlapping area.

[0164] Step 4: Combine the global sinc function compensation coefficient and the polynomial correction coefficient to achieve uniform color correction for each original image.

[0165] Combining global and local processing, the color differences between different images can be eliminated without specifying a reference image, thus avoiding error transmission between images. The global optimization strategy introduces the global sinc function model to solve and obtain compensation parameters, transforming the color difference elimination problem into a least squares optimization problem, effectively eliminating the total color difference. The local optimization strategy is a supplement to the global optimization strategy. Based on the pixel values ​​corrected by the global sinc function model, each original image is further processed using a polynomial fitting model to eliminate the color difference in the overlapping areas of the original images.

[0166] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the above embodiments, ordinary technicians in the relevant field should understand that the specific implementation methods of the present invention can still be modified or replaced by equivalents. Any modification or equivalent replacement that does not depart from the spirit and scope of the present invention should be covered within the scope of protection of the claims of the present invention.

Claims

1. A method for color-homogenizing multi-scene SAR images considering terrain radiation, characterized in that: The steps include: Step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points; Step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the overall SAR image of the test area, and obtaining the global sinc function compensation coefficient; Step S3, uses the idea of ​​regional block adjustment to construct polynomial fitting for each original SAR image in the area to be measured, and solves the polynomial model correction parameters of all SAR images; Step S4, combining the global sinc function compensation coefficient and the polynomial correction coefficient to perform uniform color correction on the original SAR image.

2. The method for color-homogenizing multi-scene SAR images considering terrain radiation according to claim 1, characterized in that: The step S1, obtaining a SAR image of the area to be measured, preprocessing the obtained SAR image and extracting sampling points, includes: The original SAR images are geometrically corrected at the same resolution, all the original SAR images are unified in the same object projection coordinate system, and an overlapping relationship based on the unified geographical location is established; In the object coordinate system of the SAR image, a sampling grid is drawn every 50 m in the horizontal and vertical directions; At the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, the average brightness value of 2500 pixels is calculated as the brightness value of the sampling point; The geometric model of the SAR image is used to transform the position of the grid point to the position on the original SAR image, and the brightness value of the sampling point is used as the brightness value of the original image position, and the incident angle information of the sampling point is extracted: sieve=(H s -H g ) / (R0+c×m c ) Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

3. The method for color-homogenizing multi-scene SAR images considering terrain radiation according to claim 1 or 2, characterized in that: The step S2, substituting the sampling points into the sinc function model, performing global optimization by fitting the overall SAR image of the test area, and obtaining the global sinc function compensation coefficient, comprises: Build the sinc function model: Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle; The error equation is established, the brightness and incident angle information of all original image sampling points in the test area are substituted into the error equation, and the optimal estimates of χ1, χ2, and χ3 are solved using the nonlinear least squares method: Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point. Calculate the sinc fitting value G0 for the image point coordinates of any image in the test area and the incident angle of the sampling point; The global sinc function compensation coefficient is obtained by calculating the brightness value after radiation correction: In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

4. The method for color-homogenizing multi-scene SAR images considering terrain radiation according to claim 3, characterized in that: Described step S3, uses the idea of ​​regional block adjustment to construct polynomial fitting for each original SAR image in the test area, and solves the polynomial model correction parameters of all SAR images, including: According to the principle that the brightness of different original images at the same location should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal, the error equation for each group of sampling points is established: dV=a0+a1c+a2c 2 +... In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the column coordinate of the image, and a0, a1, a2... are the polynomial coefficients of the column coordinate c; Select images with pairwise overlapping relationships on all grid points in turn to construct error equation groups; Solve the error equations of all sampling points jointly, and then use the least squares principle to calculate the optimal solution of the parameters; The normal equations are solved using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

5. A multi-scene SAR image color processing device considering terrain radiation, characterized in that: include: SAR image acquisition module, used to acquire SAR images of the area to be measured, pre-process the acquired SAR images and extract sampling points; The global optimization module is used to substitute the sampling points into the sinc function model, perform global optimization by fitting the overall SAR image of the test area, and obtain the global sinc function compensation coefficient; The polynomial fitting module is used to construct a polynomial fitting for each original SAR image in the test area using the idea of ​​regional block adjustment, and calculate the polynomial model correction parameters of all SAR images; The uniform color correction module is used to perform uniform color correction on the original SAR image by combining the global sinc function compensation coefficient and the polynomial correction coefficient.

6. The device for color-homogenizing multi-scene SAR images considering terrain radiation according to claim 5, characterized in that: The SAR image acquisition module includes: The geometric correction module is used to perform geometric correction processing on the original SAR images at the same resolution, unify all the original SAR images under the same object projection coordinate system, and establish an overlapping relationship based on the unified geographical location; The sampling grid drawing module is used to draw the sampling grid every 50m in the horizontal and vertical directions in the object coordinate system of the SAR image; A brightness value calculation module is used to calculate the average brightness value of 2500 pixels at the intersection center of the sampling grid, with 50×50 pixels as the sampling unit, as the brightness value of the sampling point; The incident angle extraction module is used to convert the position of the grid point to the position on the original SAR image using the geometric model of the SAR image, and use the brightness value of the sampling point as the brightness value of the original image position, and extract the incident angle information of the sampling point: sieve=(H s -H g ) / (R0+c×m c ) Where, sita is the incident angle of the sampling point; Hs is the sensor elevation; Hg is the elevation of the grid point; R0 is the initial slant distance; c is the converted image column coordinate, m c is the resolution of column coordinate c.

7. The device for color-homogenizing multi-scene SAR images taking terrain radiation into consideration according to claim 5 or 6, characterized in that: The global optimization module comprises: Function model building module, used to build the sinc function model: Wherein, θ and G are the incident angle corresponding to the image of the entire test area and the average brightness of the sampling points representing the entire test area, respectively; χ1 is the amplitude coefficient; χ2 is the incident angle scale coefficient; χ3 is the offset of the incident angle on the horizontal axis; θ is the magnitude of the incident angle; The best estimate solution module is used to establish the error equation, substitute the brightness and incident angle information of all original image sampling points in the test area into the error equation, and use the nonlinear least squares method to solve the best estimates of χ1, χ2, and χ3: Where, the initial values ​​of the fitting function parameters χ1, χ2, and χ3 are respectively the maximum pixel value of the sampling point, twice the difference between the maximum incident angle and the minimum incident angle of the sampling point, and the mean of the incident angles of the sampling points. N is the number of sampling points, and i is the serial number of the sampling point. The sinc fitting module is used to calculate the sinc fitting value G0 of the image point coordinates and the incident angle of the sampling point of any image in the test area; The compensation coefficient acquisition module is used to obtain the global sinc function compensation coefficient by calculating the brightness value after radiation correction: In the formula, Pixel i sinc is the brightness value after sinc correction; Pixel i 0 is the brightness value of the original image; G0 is the fitting value of sinc at this image point.

8. The device for color-homogenizing multi-scene SAR images considering terrain radiation according to claim 7, characterized in that: The uniform color correction module comprises: The error equation establishment module is used to establish the error equation for each group of sampling points based on the principle that the brightness of different original images at the same position should be equal after correction, that is, the sinc correction value + polynomial model correction value on different images of the same object point should be equal: dV=a0+a1c+a2c 2 In the formula, represents the grayscale value corrected by the sinc function of the i-th and j-th images respectively; dV i , dV j are the polynomial fitting models of the i-th and j-th images respectively, c in dV represents the column coordinate of the image, and a0, a1, and a2 are the polynomial coefficients of the column coordinate c; The error equation group establishment module is used to sequentially select images with two-to-two overlapping relationships on all grid points to construct the error equation group; The joint solution module is used to jointly solve the error equations of all sampling points and then calculate the optimal solution of the parameters using the least squares principle; The polynomial coefficient acquisition module is used to solve the normal equation group using ridge estimation to obtain the polynomial model coefficients of the distance coordinate c of each original image.

Images