Water body hyperspectral abnormal region spectral information correction method

Abstract

The invention belongs to the field of water quality remote sensing monitoring, and particularly relates to a water body hyperspectral abnormal region spectral information correction method which comprises the following steps: arranging a hyperspectral scanning imaging monitoring system; obtaining hyperspectral data; extracting a flare area; calculating the spectral difference of the flare region relative to the normal region; calculating the correlation between the reflectivity and the turbidity; selecting a characteristic wave band based on the spectral difference, the correlation and the signal-to-noise ratio, and carrying out nonlinear regression on the visible light wave band based on the spectral difference of the characteristic wave bands of the flare region and the normal region; s8, calculating a normalized distance coefficient based on the position from the point to the center of the flare; and based on the distance coefficient, a physical method and a geostatistical interpolation method are integrated to carry out weighted resampling. According to the method, samples at the positions corresponding to the flares of the coastal water body and the center of the water body are collected and the turbidity of the samples is analyzed, so that the universality of the flares correction method to water bodies with different water depths is improved, and a geostatistical interpolation method is comprehensively applied to improve the correction effect.

Description

Technical Field

[0001] The invention belongs to the field of water quality remote sensing monitoring, and in particular relates to a method for correcting spectral information of a water body hyperspectral abnormal region. Background Art

[0002] The hyperspectral scanning imaging monitoring system is carried on an unmanned aerial vehicle (UAV) carrier or deployed on a high tower on the shore to obtain rich two-dimensional spatial hyperspectral information of the target area, and converts the hyperspectral features into multi-parameter information of the water ecological environment with the help of intelligent inversion algorithms.

[0003] During the hyperspectral scanning imaging monitoring process, the flare area produced by the specular reflection of sunlight is likely to interfere with the hyperspectral water quality inversion results. The most common processing method is to perform remote sensing image correction on the flare based on the assumption that the water-off radiation in the near-infrared band is close to zero. However, this method has certain limitations. For shallow or highly turbid water bodies, the bottom reflected radiation in the near-infrared band is difficult to ignore. Summary of the invention

[0004] In order to solve the problems existing in the above-mentioned prior art, the present invention provides a method for correcting spectral information of abnormal hyperspectral areas of water bodies, which can improve the interference of flare areas produced by mirror reflection of sunlight on hyperspectral water quality inversion results, and help improve the accuracy of hyperspectral monitoring.

[0005] The specific technical solution adopted by the present invention is:

[0006] The method for correcting spectral information of abnormal hyperspectral regions of water bodies comprises the following steps:

[0007] S1. Deployment of hyperspectral scanning imaging monitoring system;

[0008] S2. Acquire hyperspectral data;

[0009] S3. extracting the flare area;

[0010] S4. Calculate the spectrum difference between the flare area and the normal area;

[0011] S5. Calculate the correlation between reflectance and turbidity;

[0012] S6. Select characteristic bands based on spectral differences, correlation, and signal-to-noise ratio;

[0013] S7. Perform nonlinear regression on the visible light band based on the spectral difference of the characteristic bands between the flare area and the normal area;

[0014] S8. Calculate a normalized distance coefficient based on the position of the point to the center of the flare;

[0015] S9. Based on the distance coefficient, weighted sampling is performed by integrating physical methods and geostatistical interpolation methods.

[0016] In the step S1, the fan-shaped research area is determined according to the range of the water body of interest, the installation height and inclination angle are determined, and the deployment of the hyperspectral scanning imaging monitoring system is completed.

[0017] In the step S2, scanning is performed by a hyperspectral scanning imaging monitoring system; water samples are collected at positions corresponding to the flare phenomenon in the nearshore water body and the center of the water body within the camera field of view within the scanning range, and the collection time, geographical location and turbidity are recorded.

[0018] In the step S3, the flare area in the scanning result is extracted, and the most spectrally similar adjacent pixels are identified within the water body through a seed region growing algorithm, and the geographic coordinate system of the sampling position is used as the seed coordinate system of the seed region growing algorithm to minimize the spectral Euclidean distance as a measure of regional spectral similarity; a threshold is used to stop regional growth, starting from a lower threshold and increasing in steps of 1; if a numerical explosion is encountered, the last threshold before the numerical explosion is used and the process is repeated with a step of 0.1 until another numerical explosion is caused.

[0019] In the step S4, within the detectable near-infrared band of 700-1000 nm, the spectrum difference D of the flare area extracted in a certain band relative to the normal area is calculated. glint (λ) is defined as:

[0020] D glint (λ)=∑ x，y，t |R glint (x, y, t, λ)-R normal (x, y, t, λ)|

[0021] Among them, R glint (x, y, t, λ), R normal (x, y, t, λ) represent the reflectivity of the flare area and the normal area adjacent to the flare at wavelength λ, t represents time, and x and y represent the two dimensions of the geographic coordinate system.

[0022] In step S5, the correlation between reflectivity and turbidity is:

[0023]

[0024] In step S6, the characteristic band λ is calculated min The method is,

[0025]

[0026] Among them, P T(λ)represents the correlation coefficient between the reflectivity and turbidity of the flare area at wavelength λ, T x，y，t , They represent the turbidity and the average turbidity on the two-dimensional plane (x, y) at time t, respectively, and SNR is the signal-to-noise ratio.

[0027] In the step S7, nonlinear regression is performed on the visible light band in the flare region to correct the reflectivity, and the significant influence of different aerosol types in the multi-temporal and spatial backgrounds on the correction process of the flare region is eliminated;

[0028]

[0029] Among them, R cor，i is the corrected reflectivity of the ith band, R i is the reflectivity of the i-th band before correction, is the flare area λ min The reflectivity, is the pixel reflectivity of the normal area near the flare, is the nonlinear regression parameter, R air is the atmospheric correction factor.

[0030] In the step S8, in each flare area The smallest point P is taken as the center, the distance from the point to the center P is normalized, and the distance coefficient w of the farthest point is recorded as 1 as a parameter for resampling the reflectivity;

[0031] In step S9, the reflectivity resampling calculation method is:

[0032] R final，i (x, y, t) = wR geo，i (x,y,t)+(1-w)R cor，i (x, y, t)

[0033] Where R final，i (x, y, t) is the reflectivity of the ith band finally obtained, R geo，i (x, y, t) is the reflectivity obtained by geostatistical interpolation based on the reflectivity of each band in the non-flare area.

[0034] The beneficial effects of the present invention are:

[0035] The present invention collects samples of flare corresponding positions in nearshore water bodies and water body centers and analyzes their turbidity to increase the versatility of the flare correction method for water bodies of different depths; selects characteristic bands with a small correlation with turbidity for image correction in the flare area to reduce the influence of bottom reflected radiation in turbid water bodies. Considering that the correction effect of physical methods at the edge of the flare area is not ideal, a geostatistical interpolation method is also comprehensively applied to improve the correction effect. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1 It is a schematic diagram of the process of the present invention; DETAILED DESCRIPTION

[0037] The present invention will be further described below in conjunction with the accompanying drawings and specific embodiments:

[0038] Specific implementation examples Figure 1 As shown, the present invention is a method for correcting spectral information of abnormal areas of water body hyperspectral, comprising the following steps:

[0039] S1. Deployment of hyperspectral scanning imaging monitoring system: Establish a monitoring system to collect hyperspectral data of water bodies;

[0040] S2. Obtain hyperspectral data: collect reflectance data of water bodies at different wavelengths;

[0041] S3. Extract flare area: Identify the spectral anomaly area caused by water surface reflection;

[0042] S4. Calculate the spectral difference between the flare area and the normal area: quantify the spectral difference between the abnormal area and the surrounding normal area;

[0043] S5. Calculate the correlation between reflectance and turbidity: Establish the correlation between spectral reflectance and water turbidity for subsequent correction;

[0044] S6. Select characteristic bands based on spectral differences, correlation, and signal-to-noise ratio: Select characteristic bands that distinguish abnormal and normal areas;

[0045] S7. Perform nonlinear regression on the visible light band based on the spectral difference of the characteristic bands in the flare area and the normal area: adjust the spectral data of the abnormal area to make it closer to the data of the normal area;

[0046] S8. Calculate a normalized distance coefficient based on the position of the point to the center of the flare;

[0047] S9. Based on the distance coefficient, weighted sampling is performed by integrating physical methods and geostatistical interpolation methods.

[0048] Furthermore, in the step S1, the fan-shaped research area is determined according to the range of the water body of interest, the installation height and inclination are determined, the hyperspectral host, light source, turntable, camera, etc. are arranged on a high tower or column, and the circuit is connected and then powered on for detection to complete the arrangement of the hyperspectral scanning imaging monitoring system.

[0049] Furthermore, in the step S2, scanning is performed by a hyperspectral scanning imaging monitoring system; at the same time, water samples are collected from the positions corresponding to the flare phenomenon in the nearshore water body and the center of the water body within the camera field of view within the scanning range, and the time and geographical location of the collection are recorded. The samples are sent to the laboratory for chemical analysis of their turbidity for subsequent selection of characteristic bands, and the water depth difference between the nearshore water body and the center of the water body is used to increase the versatility of the flare correction method for different water depths, and a preprocessing method is used to achieve lens correction, reflectivity correction, and atmospheric correction.

[0050] Furthermore, in the step S3, the flare area in the scanning result is extracted, and the most spectrally similar adjacent pixels are identified within the water body through a seed region growing algorithm, and the geographic coordinate system of the sampling position is used as the seed coordinate system of the seed region growing algorithm to minimize the spectral Euclidean distance as a measure of regional spectral similarity; a threshold is used to stop regional growth, starting from a lower threshold and increasing in steps of 1; if a numerical explosion is encountered, the last threshold before the numerical explosion is used and the process is repeated with a step of 0.1 until another numerical explosion occurs.

[0051] Among them, the initial threshold usually starts from a lower threshold, such as 0 or a value close to 0, to ensure that pixels similar to the seed point can be captured, and the threshold is gradually increased in steps of 1. The growth process of the seed region is re-executed after each increase. If the result suddenly becomes very large or out of control during the increase of the threshold, it is a numerical explosion. The last threshold before the numerical explosion is used, and the threshold is continued to be adjusted in steps of 0.1 until the numerical explosion occurs again. The threshold before this numerical explosion is a more appropriate threshold, which can be used for image segmentation operations of the seed region growing algorithm.

[0052] Furthermore, in step S4, within the detectable near-infrared band of 700-1000 nm, the spectrum difference D of the flare area extracted in a certain band relative to the normal area is calculated. glint (λ) is defined as:

[0053] D glint (λ)=∑ x，y，t |R glint (x, y, t, λ)-R normal (x, y, t, λ)|

[0054] Among them, R glint (x, y, t, λ), R normal (x, y, t, λ) represent the reflectivity of the flare area and the normal area adjacent to the flare at wavelength λ, t represents time, and x and y represent the two dimensions of the geographic coordinate system.

[0055] When determining the characteristic band for flare area correction, the following conditions are considered comprehensively:

[0056] 1. Assume that at the selected characteristic wavelength, the spectrum of the flare area is significantly different from that of the adjacent normal area, and there is a high degree of distinction between the flare area and the non-flare area.

[0057] 2. Assuming that the correlation coefficient between reflectance and turbidity is small at the selected characteristic wavelength, the model is more versatile for turbid water bodies.

[0058] 3. Assuming that at the selected characteristic wavelength, the reflectivity is low, the impact of bottom reflected radiation is even less significant.

[0059] 4. Assume that the signal-to-noise ratio is large at the selected characteristic wavelength to avoid introducing too much noise. Based on the above assumption, in step S5, the correlation between reflectivity and turbidity is:

[0060]

[0061] In step S6, the characteristic band λ is calculated min The method is,

[0062]

[0063] Among them, P T(λ) represents the correlation coefficient between the reflectivity and turbidity of the flare area at wavelength λ, T x，y，t , They represent the turbidity and the average turbidity on the two-dimensional plane (x, y) at time t, respectively, and SNR is the signal-to-noise ratio.

[0064] In the step S7, nonlinear regression is performed on the visible light band in the flare region to correct the reflectivity, and the significant influence of different aerosol types in the multi-temporal and spatial backgrounds on the correction process of the flare region is eliminated;

[0065]

[0066] Among them, R cor，i is the corrected reflectivity of the ith band, R i is the reflectivity of the i-th band before correction, is the flare area λ min The reflectivity, is the pixel reflectivity of the normal area near the flare, is the nonlinear regression parameter, R air is the atmospheric correction factor.

[0067] Considering that the correction effect of physical method at the edge of flare area is not ideal, the reflectivity of flare area is resampled. In step S8, the reflectivity of each flare area is resampled. The smallest point P is taken as the center, the distance from the point to the center P is normalized, and the distance coefficient w of the farthest point is recorded as 1 as a parameter for resampling the reflectivity;

[0068] In step S9, the reflectivity resampling calculation method is:

[0069] R final，i (x, y, t) = wR geo，i (x,y,t)+(1-w)R cor，i (x, y, t)

[0070] Where R final，i (x, y, t) is the reflectivity of the ith band finally obtained, R geo，i (x, y, t) is the reflectivity obtained by geostatistical interpolation based on the reflectivity of each band in the non-flare area.

Claims

1. A method for correcting spectral information of abnormal areas of water body hyperspectral, characterized in that The steps include: S1. Deployment of hyperspectral scanning imaging monitoring system; S2. Acquire hyperspectral data; S3. extracting the flare area; S4. Calculate the spectrum difference between the flare area and the normal area; S5. Calculate the correlation between reflectance and turbidity; S6. Select characteristic bands based on spectral differences, correlation, and signal-to-noise ratio; S7. Perform nonlinear regression on the visible light band based on the spectral difference of the characteristic bands between the flare area and the normal area; S8. Calculate a normalized distance coefficient based on the position of the point to the center of the flare; S9. Based on the distance coefficient, weighted sampling is performed by integrating physical methods and geostatistical interpolation methods.

2. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S1, the fan-shaped research area is determined according to the range of the water body of interest, the installation height and inclination angle are determined, and the deployment of the hyperspectral scanning imaging monitoring system is completed.

3. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S2, scanning is performed by a hyperspectral scanning imaging monitoring system; water samples are collected at positions corresponding to the flare phenomenon in the nearshore water body and the center of the water body within the camera field of view within the scanning range, and the collection time, geographical location and turbidity are recorded.

4. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S3, the flare area in the scanning result is extracted, and the most spectrally similar adjacent pixels are identified within the water body by using a seed region growing algorithm, and the geographic coordinate system of the sampling position is used as the seed coordinate system of the seed region growing algorithm to minimize the spectral Euclidean distance as a measure of regional spectral similarity; Use a threshold to stop region growth, starting from a lower threshold and increasing in steps of 1; if a numerical explosion occurs, use the last threshold before the numerical explosion and repeat the process in steps of 0.1 until another numerical explosion occurs.

5. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S4, within the detectable near-infrared band of 700-1000 nm, the spectrum difference D of the flare area extracted in a certain band relative to the normal area is calculated. glint (λ) is defined as: D glint (λ)=Σ x,y,t |R glint (x,y,t,λ)-R normal (x,y,t,λ)| Among them, R glint (x, y, t, λ), R normal (x, y, t, λ) represent the reflectivity of the flare area and the normal area adjacent to the flare at wavelength λ, t represents time, and x and y represent the two dimensions of the geographic coordinate system.

6. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In step S5, the correlation between reflectivity and turbidity is: In step S6, the characteristic band λ is calculated min The method is, Among them, P T (λ) represents the correlation coefficient between the reflectivity and turbidity of the flare area at wavelength λ, T x，y，t , They represent the turbidity and the average turbidity on the two-dimensional plane (x, y) at time t, respectively, and SNR is the signal-to-noise ratio.

7. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S7, nonlinear regression is performed on the visible light band in the flare region to correct the reflectivity, and the significant influence of different aerosol types in the multi-temporal and spatial backgrounds on the correction process of the flare region is eliminated; Among them, R cor，i is the corrected reflectivity of the ith band, R i is the reflectivity of the i-th band before correction, is the flare area λ min The reflectivity, is the pixel reflectivity of the normal area near the flare, is the nonlinear regression parameter, R air is the atmospheric correction factor.

8. The method for correcting spectral information of abnormal areas of water body hyperspectral according to claim 1, characterized in that: In the step S8, in each flare area The smallest point P is the center, and the point to the center P The distance is normalized, and the distance coefficient w of the farthest point is recorded as 1, which is used as a parameter for resampling the reflectivity; In step S9, the reflectivity resampling calculation method is: R final，i (x,y,t)=wR geo,i (x,y,t)+(1-w)R cor，i (x,y,t) Where R final，i (x, y, t) is the reflectivity of the ith band finally obtained, R geo，i (x, y, t) is the reflectivity obtained by geostatistical interpolation based on the reflectivity of each band in the non-flare area.

Images