Construction of Rayleigh Scattering Look-up Table and Search Method for Environmental Disaster Reduction Hyperspectral Satellite
By calculating the radiation parameters of each band of HJ-2A/B satellite and the numerical solution of the ocean-atmospheric coupled vector radiation transmission model, a high-precision Rayleigh scattering lookup table was generated, which solved the problem of insufficient atmospheric correction accuracy of ocean remote sensing in the existing technology, and realized high-precision Rayleigh scattering radiation brightness calculation.
Patent Information
- Application Number
- CN202310111681.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-02-14
- Publication Date
- 2025-06-13
- Estimated Expiration
- 2043-02-14
AI Technical Summary
It is difficult to establish a high-precision Rayleigh scattering lookup table suitable for environmental disaster reduction AB satellite hyperspectral imager, resulting in insufficient accuracy of ocean remote sensing atmospheric correction, affecting the performance of quantitative remote sensing monitoring.
By calculating the equivalent atmospheric top solar incident irradiance, atmospheric molecule Rayleigh scattering optical thickness and ozone ratio absorption coefficient of each band of HJ-2A/B satellite, the radiation transmission equation is numerically solved with the ocean-atmospheric coupled vector radiation transmission model, a high-precision Rayleigh scattering lookup table is generated, and the atmospheric molecule Rayleigh scattering optical thickness is revised according to the actual atmospheric pressure.
A high-precision numerical solution to the ocean-atmospheric vector radiation transmission equation is achieved, meeting the atmospheric correction calculation requirements of 0.5%, and the Rayleigh scattering luminance is quickly and accurately calculated, filling the gap in the atmospheric correction Rayleigh scattering lookup table of the HJ-2A/B satellite hyperspectral imager.
Smart Images

Figure CN116256316B_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the field of ocean remote sensing technology, and particularly relates to a method for constructing and searching a Rayleigh scattering look-up table for an environmental disaster reduction hyperspectral satellite. Background Art
[0002] Ocean satellite remote sensing realizes the quantitative remote sensing inversion of ocean water color elements (chlorophyll, suspended sediment, yellow substance) through the water-leaving radiation in the visible light band received by satellite sensors, and is applied to the monitoring of the ecological environment of ocean and lake waters, the monitoring of ocean disasters, and the maintenance of ocean rights. Typically, when the hyperspectral imager of the Environmental Disaster Reduction Satellite 2A / B (hereinafter referred to as: HJ2A / B satellite) monitors the ocean environment, more than 80% of the total radiation received comes from atmospheric scattering radiation and sea surface reflection radiation. The process of subtracting atmospheric scattering and sea surface reflection radiation from the total radiation received by the hyperspectral imager is called atmospheric correction. Generally, to achieve a water-leaving radiation inversion accuracy of 5% after atmospheric correction, the Rayleigh scattering calculation accuracy is required to reach 1%. Usually, the Rayleigh scattering calculation adopts the single-scattering approximation calculation method and the numerical solution of the atmospheric vector radiative transfer equation method. To achieve a Rayleigh scattering calculation accuracy of 1%, the atmospheric vector radiative transfer equation must be accurately solved. The Rayleigh scattering calculation accuracy plays a key role in the accuracy of ocean water color remote sensing atmospheric correction. Therefore, in order to achieve high-precision atmospheric correction of the HJ2A / B satellite hyperspectral imager, a high-precision Rayleigh look-up table applicable to the HJ2A / B satellite hyperspectral imager needs to be established.
[0003] However, due to the lack of a high-precision Rayleigh look-up table, a business-oriented ocean remote sensing atmospheric correction algorithm has not been established in the process of applying the image data of the HJ2A / B satellite hyperspectral imager in ocean remote sensing, which seriously affects the quantitative remote sensing monitoring application efficiency of HJ2A / B in the ecological environment of ocean and lake waters. There are mainly three methods for calculating atmospheric Rayleigh scattering, namely the single-scattering approximation method, the multi-scattering approximation method, and the numerical solution of the atmospheric vector radiative transfer equation method. The calculation error of the single-scattering approximation method is less than 5% when the solar zenith angle is less than 40° and the remote sensor observation zenith angle is less than 60°. However, when the solar zenith angle is greater than 40° or the remote sensor observation zenith angle is greater than 60°, the calculation error increases rapidly. Therefore, it is not applicable to the high-precision ocean quantitative remote sensing of HJ-2A / B. The multi-scattering method refers to the Rayleigh approximation calculation method that considers the multiple scattering of atmospheric molecules and ignores the scattering polarization characteristics. Typically, it is the Rayleigh scattering calculation model in the 6S model. When the solar zenith angle is less than 40° and the remote sensor observation zenith angle is less than 60°, the accuracy of the 6S Rayleigh scattering approximation calculation model is worse than that of the single-scattering approximation method, and the relative error is generally greater than 5%. When the remote sensor observation zenith angle is greater than 60°, the 6S Rayleigh scattering approximation calculation model is significantly less error-prone than the single Rayleigh scattering approximation. Generally speaking, the 6S Rayleigh scattering approximation calculation model is not applicable to the atmospheric correction of the HJ2A / B satellite hyperspectral imager with very high accuracy requirements.
[0004] To achieve a Rayleigh scattering calculation accuracy of 1%, it is necessary to numerically solve the complex atmospheric vector radiative transfer equation. Since the atmospheric correction of the HJ-2A / B satellite hyperspectral imager is carried out pixel by pixel, it is impossible to solve the vector radiative transfer equation for each pixel. Therefore, an accurate Rayleigh scattering calculation method that can ensure both accuracy and speed is required. Currently, the lookup table method is adopted, that is, the ocean-atmosphere vector radiative transfer equation is numerically solved in advance to generate accurate Rayleigh scattering values under various HJ-2A / B satellite observation geometric conditions, and they are saved in the lookup table file. During the actual application process, only the accurate Rayleigh scattering value under the geometric condition closest to the actual situation needs to be found and linear interpolation is performed. However, the lookup tables currently applied in ocean remote sensing satellites are generated for specific remote sensors of ocean satellites and the bands are multispectral. In addition, due to the problem of the air pressure correction method in the current Rayleigh scattering lookup table of traditional ocean color satellite remote sensors, it cannot be directly applied to the HJ-2A / B satellite hyperspectral imager. Moreover, due to the weak anisotropy of atmospheric molecules, the scattered radiation in the direction of 90° of the incident radiation is not fully polarized but has a certain "depolarization", and its scattering matrix needs to be corrected. Therefore, it is necessary to establish a more accurate Rayleigh scattering lookup table specifically for the band characteristics of the HJ2A / B satellite hyperspectral imager. Summary of the Invention
[0005] The present invention proposes a method for constructing and looking up a Rayleigh scattering lookup table applicable to environmental disaster reduction hyperspectral satellites, which fully considers the band characteristics of the HJ2A / B instrument, the anisotropy of atmospheric molecules, and the influence of actual air pressure and standard air pressure, and can effectively overcome the influence of the anisotropy band difference of atmospheric molecules and the air pressure correction on the calculation accuracy of the Rayleigh scattering radiation vector under the equivalent Rayleigh scattering optical thickness of each band of the HJ-2A / B satellite hyperspectral imager.
[0006] The concept of the present invention is as follows: Calculate the solar irradiance incident on the top of the atmosphere, the Rayleigh scattering optical thickness of atmospheric molecules under standard atmospheric pressure, and the specific absorption coefficient of ozone per unit volume concentration according to the spectral response function of the HJ-2A / B satellite hyperspectral imager, and calculate the Rayleigh scattering radiation vector under the equivalent Rayleigh scattering optical thickness of each band under different wind speed and solar zenith angle conditions. Then generate a high-precision Rayleigh scattering lookup table for the HJ-2A / B satellite hyperspectral imager. When actually performing ocean atmosphere correction on the satellite data of the HJ-2A / B satellite hyperspectral imager, revise the Rayleigh scattering optical thickness of atmospheric molecules according to the actual atmospheric pressure at the imaging time, and then look up the Rayleigh scattering value in the Rayleigh scattering lookup table according to the observation geometry and solar zenith angle at the imaging time.
[0007] To achieve the above object, the technical solution adopted by the present invention is as follows:
[0008] A method for constructing and searching a Rayleigh scattering lookup table for an environmental disaster reduction hyperspectral satellite, which is characterized in that it includes the following steps:
[0009] Step 1: According to the spectral response function of the hyperspectral imager of the HJ-2A / B satellite, calculate the equivalent solar irradiance F at the top of the atmosphere for each band of the HJ-2A / B satellite 0 (i), the Rayleigh scattering optical thickness τ of atmospheric molecules r (i) and the specific absorption coefficient A of ozone per unit volume concentration OZ (i);
[0010] Step 2: Based on the anisotropic band difference and depolarization effect of atmospheric molecules, calculate the spherical atmospheric molecular scattering phase matrix P(Θ) of the HJ-2A / B satellite hyperspectral imager;
[0011] Step 3: Combining the parameter values obtained in Step 1 and Step 2, numerically solve the radiative transfer equation based on the ocean-atmosphere coupled vector radiative transfer model, and calculate the Rayleigh scattering radiation vector (I, Q, U, V) received by the satellite at the top of the atmosphere under different atmospheric molecular optical thicknesses τ r (λ), solar zenith angle θ 0 , observation zenith angle θ v , and wind speed W conditions, and construct a Rayleigh scattering lookup table applicable to the environmental disaster reduction satellite imager; T where I is the total radiation intensity, Q is the linearly polarized radiation intensity in the 45° direction, U is the linearly polarized radiation intensity in the vertical direction, and V is the circularly polarized light intensity; E
[0012] and E l and E r are the electric vectors of the parallel and vertical reference planes respectively; and are the conjugate matrices of E l and E r respectively;
[0013] Step 4: Calculate the equivalent Rayleigh scattering optical thickness τ of atmospheric molecules according to the Rayleigh scattering optical thickness of the atmosphere under the actual atmospheric pressure r (i), and complete the search calculation of Rayleigh scattering applicable to the environmental disaster reduction satellite imager under different atmospheric conditions based on the Rayleigh scattering lookup table constructed in Step 3.
[0014] Furthermore, in Step 1, the equivalent solar irradiance F at the top of the atmosphere for each band of the HJ-2A / B satellite 0 (i), the Rayleigh scattering optical thickness τ of atmospheric molecules r (i) and the specific absorption coefficient A of ozone per unit volume concentration OZ (i), the calculation formulas are as follows:
[0015]
[0016] Among them, \(i\) is the band number of the hyperspectral imager on HJ-2A / B satellite, and \(1\leq i\leq100\); \(\lambda\) i is the wavelength corresponding to the \(i\)-th band; \(\lambda\) i1 and \(\lambda\) i2 are the left and right end wavelengths corresponding to the spectral responsivity of 1% in the \(i\)-th band respectively, and \(\lambda\) i1 \(\leq\lambda\) i \(\leq\lambda\) i2 ;
[0017] F 0 (\(\lambda\) i ) is the solar irradiance at the top of the atmosphere corresponding to the \(i\)-th band;
[0018] S(\(\lambda\) i ) is the spectral responsivity corresponding to the \(i\)-th band;
[0019] \(\tau\) r (\(\lambda\) i ) is the Rayleigh scattering optical depth of atmospheric molecules under standard atmospheric pressure corresponding to the \(i\)-th band;
[0020] A OZ (\(\lambda\) i ) is the ozone specific absorption coefficient per unit volume concentration corresponding to the \(i\)-th band.
[0021] Furthermore, in step 2, the spherical atmospheric molecular scattering phase matrix \(P(\Theta)\) of the HJ-2A / B satellite hyperspectral imager is as follows:
[0022]
[0023] Among them, \(\rho(\lambda\) i ) is the depolarization factor corresponding to the wavelength \(\lambda\) in the \(i\)-th band, which takes into account the band difference of atmospheric molecular anisotropy;
[0024] The calculation formula of the depolarization factor \(\rho(\lambda\) i ) is as follows:
[0025]
[0026] Furthermore, step 3 is specifically as follows:
[0027] 3.1. Set the input parameters of the ocean-atmosphere coupled vector radiative transfer model;
[0028] The ocean-atmosphere coupled vector radiative transfer model includes the coupling of atmospheric stratification, sea-air interface stratification, and ocean stratification;
[0029] 3.1.1. Set the atmosphere to be divided into N layers, each layer consisting of pure atmospheric molecules. The input parameters include the optical thickness τ of atmospheric molecules r (λ), the scattering phase matrix P(Θ), and the single-scattering albedo ω of the atmosphere, where N≥10, and the optical thickness τ of atmospheric molecules r (λ) is numerically input as 0.002 - 0.6 with a step size of 0.002; ω = 1.0;
[0030] 3.1.2. Set the air-sea interface layer to be the probability density function of the distribution of the wavelet surface constructed based on statistical results, where the input parameter is the wind speed W, and the wind speed ranges from 0 to 24 m / s with an interval of 4 m / s;
[0031] The probability density function p(e n ) of the wavelet surface is:
[0032]
[0033] where e n is the unit normal vector of the wavelet surface; μ n is the cosine of the angle between the wavelet surface normal and the Z-axis; n is the number of randomly distributed wavelet surfaces on the rough sea surface; σ is the variance of the probability density of the wavelet surface distribution; σ 2 = 0.003 + 0.00512W;
[0034] 3.1.3. Set the ocean to be divided into M layers. The inherent optical property parameters of the seawater layer include the absorption coefficient a water , the scattering phase matrix, and the single-scattering albedo ω water of the ocean, where M≥1, a water ≥10, ω water = 1, and the ocean depth is infinitely deep;
[0035] 3.2. Conduct vector radiative transfer simulation to calculate the Rayleigh scattering radiation vector (I, Q, U, V) r under different optical thicknesses τ of atmospheric molecules 0 (λ), solar zenith angle θ v , observation zenith angle θ T , and wind speed W, and establish a high-precision Rayleigh scattering look-up table for the HJ-2A / B satellite hyperspectral imager.
[0036] Furthermore, in step 3, the Rayleigh scattering radiation vector (I, Q, U, V) T is calculated from the Stokes vector S, and the Stokes vector S is:
[0037]
[0038] where E l and Er They are the electric vectors of the parallel and perpendicular reference planes respectively; and are E l and E r 's conjugate matrix respectively.
[0039] Furthermore, step 4 is specifically as follows:
[0040] 4.1. Revise the Rayleigh optical depth of atmospheric molecules based on the atmospheric correction algorithm:
[0041]
[0042] 4.2. Calculate the Rayleigh scattering optical depth of the atmosphere under the actual atmospheric pressure:
[0043]
[0044] where P r and P 0 are the actual sea surface atmospheric pressure and the standard atmospheric pressure respectively;
[0045] 4.3. Calculate the equivalent Rayleigh scattering optical depth τ r (i) of the atmosphere under the actual atmospheric pressure:
[0046]
[0047] 4.4. At the transit time of the HJ-2A / B satellite hyperspectral imager, the solar zenith angle θ 0 and the observation zenith angle θ v of each pixel. Based on the on-site wind speed data W, then look up and calculate the Rayleigh scattering radiance value of atmospheric molecules in the Rayleigh scattering look-up table constructed in step 3, and realize the high-precision look-up and calculation of the Rayleigh scattering radiation vector (I, Q, U, V) T of the HJ-2A / B satellite hyperspectral imager image.
[0048] Furthermore, in step 3.2, the solar zenith angle θ 0 takes values from 0 to 88°, with a step interval of 2°;
[0049] the observation zenith angle θ v takes the zeros of the 100th-order Legendre polynomial in the interval (0, 1), with an interval of approximately 1.8°;
[0050] The Legendre polynomial P m (x) is:
[0051]
[0052] where m represents the order, and x is the zero of the Legendre polynomial in the interval (0, 1).
[0053] Compared with the prior art, the beneficial technical effects of the present invention are as follows:
[0054] 1. The Rayleigh scattering look-up table construction and look-up method for the environmental disaster reduction hyperspectral satellite provided by the present invention fully considers the instrument band characteristics of the HJ-2A / B satellite hyperspectral imager and the spectral characteristics of the atmospheric molecular depolarization factor. By using the spectral response function to correct the optical thickness of atmospheric molecular Rayleigh scattering, a numerical solution of the ocean-atmosphere vector radiative transfer equation with high precision can be obtained, which can not only meet the atmospheric correction calculation requirement of 0.5%, but also achieve fast and accurate calculation of Rayleigh scattering radiance.
[0055] 2. Compared with the traditional look-up table of ocean color multispectral satellite, the look-up table of ocean color multispectral satellite usually calculates the Rayleigh scattering radiance using the standard atmosphere, and then corrects the Rayleigh scattering radiance according to the ratio of the actual sea surface atmospheric pressure to the standard atmospheric pressure. When there is a large difference between the actual sea surface atmospheric pressure and the standard atmospheric pressure, the calculation error is relatively large. The present invention fully considers factors such as the influence of the actual sea surface atmospheric pressure on the Rayleigh optical thickness of atmospheric molecules, and the calculation error is less than 0.5%. Accurate Rayleigh scattering radiance in any visible and near-infrared bands can be obtained, which can fill the blank of the Rayleigh scattering look-up table for atmospheric correction of hyperspectral remote sensing data of the HJ-2A / B satellite hyperspectral imager.
[0056] 3. The present invention uses the numerical solution of the ocean-atmosphere vector radiative transfer equation to establish an accurate Rayleigh look-up table applicable to the HJ-2A / B satellite hyperspectral imager, which has higher algorithm accuracy compared with the approximate solutions of models such as 6S, the single-scattering approximation method, and the multi-scattering approximation method. In addition, for other hyperspectral remote sensing payloads, the present invention only needs to modify the optical thickness of atmospheric Rayleigh scattering according to the spectral characteristics of other hyperspectral satellite payloads to calculate the Rayleigh scattering radiance of other hyperspectral satellite payloads, and has good adaptability and practicability. Description of the Drawings
[0057] Figure 1 It is a flow chart of the Rayleigh scattering look-up table construction and look-up method for the environmental disaster reduction hyperspectral satellite of the present invention;
[0058] Figure 2It is the spatial geometric distribution map of Rayleigh scattering radiance of HJ-2A / B hyperspectral satellite; among them, (a), (d) and (g) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 20°, (b), (e) and (h) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 40°, and (c), (f) and (i) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 60°;
[0059] Figure 3 It is the spatial geometric distribution map of the error of Rayleigh scattering radiance of HJ-2A / B hyperspectral satellite obtained by using the construction method of the present invention; among them, (a), (d) and (g) are respectively the spatial geometric distribution maps of the error of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 20°, (b), (e) and (h) are respectively the spatial geometric distribution maps of the error of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 40°, and (c), (f) and (i) are respectively the spatial geometric distribution maps of the error of Rayleigh scattering radiance at 488nm, 531nm, and 670nm when the solar zenith angle is 60°;
[0060] Figure 4 It is the graph of the variation of the error of Rayleigh scattering radiance of HJ-2A / B hyperspectral satellite with the wavelength band obtained by using the construction method of the present invention; among them, (a) to (d) are respectively the graphs of the variation of the error of Rayleigh scattering radiance with the wavelength band when the azimuth angles are 0°, 45°, 90°, 135°, 180° under the solar zenith angles of 0°, 20°, 40° and 60°. Specific embodiments
[0061] To make the objectives, advantages and features of the present invention clearer, the following further details a Rayleigh scattering look-up table construction and look-up method for an environmental disaster reduction satellite imager proposed by the present invention in conjunction with the accompanying drawings and specific embodiments. Those skilled in the art should understand that these embodiments are only used to explain the technical principle of the present invention, and the purpose is not to limit the protection scope of the present invention.
[0062] Typically, for the hyperspectral imager of the Environmental Disaster Reduction 2A / B satellite, more than 80% of the total radiation received during ocean environment monitoring comes from atmospheric scattered radiation and sea surface reflected radiation. Generally, the process of subtracting atmospheric scattered radiation and sea surface reflected radiation from the total radiation received by the hyperspectral imager is called "atmospheric correction". Generally, to achieve a 5% accuracy in retrieving the water-leaving radiance after atmospheric correction, the calculation accuracy of Rayleigh scattering is required to reach 1%. Usually, the Rayleigh scattering calculation adopts the single-scattering approximation calculation method and the numerical solution method of the ocean-atmosphere vector radiative transfer equation. To achieve a 1% accuracy in Rayleigh scattering calculation, the ocean-atmosphere vector radiative transfer equation must be accurately solved.
[0063] The present invention proposes a method for constructing and looking up the Rayleigh scattering lookup table of the Environmental Disaster Reduction satellite imager, which calculates the Rayleigh scattering radiance at the top of the atmosphere required for atmospheric correction of the HJ-2A / B satellite hyperspectral imager using the ocean-atmosphere coupled vector radiative transfer model. Considering factors such as the instrument band characteristics of the HJ-2A / B satellite hyperspectral imager, the spectral characteristics of the atmospheric molecular depolarization factor, and the influence of the actual sea surface atmospheric pressure on the Rayleigh optical thickness of atmospheric molecules, the spectral response function and the actual sea surface atmospheric pressure are used to correct the Rayleigh scattering optical thickness of atmospheric molecules, obtaining a numerical solution of the high-precision ocean-atmosphere vector radiative transfer equation, thereby establishing a Rayleigh scattering lookup table applicable to the Environmental Disaster Reduction satellite imager.
[0064] As Figure 1 shown, the method for constructing and looking up the Rayleigh scattering lookup table of the Environmental Disaster Reduction satellite imager proposed in this embodiment specifically includes the following steps:
[0065] Step 1: According to the spectral response function of the HJ-2A / B satellite hyperspectral imager, calculate the equivalent solar irradiance F at the top of the atmosphere for each band of the HJ-2A / B satellite 0 (i), the Rayleigh scattering optical thickness τ r (i) of atmospheric molecules and the specific absorption coefficient A OZ (i) of ozone per unit volume concentration, and the calculation formulas are as follows:
[0066]
[0067] Among them, i is the band number of the HJ-2A / B satellite hyperspectral imager, 1 ≤ i ≤ 100; λ i is the wavelength corresponding to the i-th band;
[0068] F 0 (λ i ) is the solar irradiance at the top of the atmosphere corresponding to the wavelength λ of the i-th band;
[0069] S(λ i ) is the spectral responsivity corresponding to the wavelength λ of the i-th band;
[0070] τ r (λ i ) is the Rayleigh scattering optical depth of atmospheric molecules at the wavelength λ of the i-th band under standard atmospheric pressure;
[0071] A OZ (λ i ) is the specific absorption coefficient of ozone per unit volume concentration at the wavelength λ of the i-th band;
[0072] λ i1 and λ i2 are the left and right end wavelengths corresponding to 1% spectral responsivity of the i-th band, respectively, and λ i1 ≤λ i ≤λ i2 ;
[0073] Step 2: Calculate the spherical atmospheric molecular scattering phase matrix P(Θ) of the HJ-2A / B satellite hyperspectral imager based on the anisotropic band differences of atmospheric molecules;
[0074] Based on the anisotropic band differences of atmospheric molecules, when revising the anisotropic band differences of atmospheric molecules, the spherical atmospheric molecular scattering phase matrix is as follows:
[0075]
[0076] where Θ is the scattering angle, 0° ≤ Θ ≤ 180°; P(Θ) is the Rayleigh scattering phase matrix of atmospheric molecules in each band of the HJ-2A / B satellite hyperspectral imager (the band λ is omitted here).
[0077] For actual atmospheric molecules, due to the anisotropy of molecules, the depolarization effect needs to be considered in the atmospheric molecular scattering. According to Young's conclusion, the depolarization factor ρ is 0.0279. Considering the depolarization factor ρ of each band of the HJ-2A / B satellite hyperspectral imager with atmospheric molecular anisotropy, the spherical atmospheric molecular Rayleigh scattering phase matrix P(Θ) after considering the depolarization effect is:
[0078]
[0079] In addition, in the process of revising the anisotropic band differences of atmospheric molecules, only considering the depolarization factor ρ still has limitations, and the variation of the depolarization factor ρ with the band needs to be considered. For short waves, the depolarization factor ρ varies more significantly with the band. Therefore, the calculation formula for a more accurate depolarization coefficient is as follows:
[0080]
[0081] Step 3: Numerically solve the radiative transfer equation based on the ocean-atmosphere coupled vector radiative transfer model PCOART to calculate and obtain different atmospheric molecular optical depths τ r(λ), solar zenith angle θ 0 , observation zenith angle θ v , the Rayleigh scattering radiation vector (I, Q, U, V) received by the satellite at the top of the atmosphere under the condition of wind speed W T The value of this radiation vector can be represented by the Stokes vector:
[0082]
[0083] Among them, I is the total radiation intensity, Q is the linearly polarized radiation intensity in the 45° direction, U is the linearly polarized radiation intensity in the vertical direction, and V is the circularly polarized light intensity, which is usually ignored in the sea - air system. E l and E r are the electric vectors of the parallel and vertical reference planes respectively; and are the conjugate matrices of E l and E r respectively.
[0084] The numerical value of this Rayleigh scattering radiation vector (I, Q, U, V) T is the numerical solution of the radiative transfer equation for each band of the HJ - 2A / B satellite hyperspectral imager under different observation geometric conditions stored in the lookup table, and then saves and constructs a high - precision Rayleigh scattering lookup table for the HJ - 2A / B satellite hyperspectral imager.
[0085] The PCOART model uses a matrix algorithm to describe the vector radiative transfer process in the ocean - atmosphere medium, and considers the effects of multiple scattering, refraction at the wind - generated rough air - water interface, and bottom - boundary reflection on the vector radiation light field of the sea - air system. This model meets the calculation accuracy requirements verified by 7 standard problems proposed by international ocean optical remote sensing. The PCOART model considers the influence of the earth's curvature on the solution of the radiative transfer equation for the sea - air system. It still has a high solution accuracy under large solar zenith angles and observation zenith angles, while most other sea - air system radiative transfer models are based on the parallel - plane hypothesis theory and do not consider the influence of the earth's curvature, resulting in poor accuracy of the radiative transfer equation for the sea - air system under large solar zenith angles and observation zenith angles and unable to meet the accuracy requirements of high - spectral satellite quantitative remote sensing.
[0086] The PCOART model sets the coupling of the atmosphere, sea - air interface, and ocean stratification in the sea - air system radiative transfer model. The following specifically introduces the optical characteristic parameters of the atmosphere, sea - air interface, and ocean stratification of the ocean - atmosphere coupled vector radiative transfer model PCOART.
[0087] In the ocean - atmosphere coupled vector radiative transfer model PCOART, the atmosphere is stratified into N layers (N≥10), and each layer is composed of pure atmospheric molecules. The input parameters include the optical thickness τ of atmospheric molecules r(λ), scattering phase matrix P(Θ), single scattering albedo ω. Atmospheric molecular optical thickness τ r (λ) is numerically input as 0.002 - 0.6 with a step size of 0.002; the scattering phase matrix P(Θ) is calculated according to step 2, and the single atmospheric scattering rate is 1.0.
[0088] In the ocean - atmosphere coupled vector radiative transfer model PCOART, the ocean is stratified into M layers (M ≥ 1), and the inherent optical property parameters of the seawater layer include the absorption coefficient a water , scattering phase matrix, and single ocean scattering albedo ω water . The ocean stratification is pure absorption without scattering, a water ≥10; the single ocean scattering albedo ω water is 1; the scattering phase matrix P(Θ) is calculated according to step 2, and the ocean stratification water depth is infinitely deep.
[0089] In the ocean - atmosphere coupled vector radiative transfer model PCOART, the air - sea interface is based on the distribution probability density function of the small - wave surface constructed by Cox and Munk based on statistical results (Cox and Munk, 1954), where the input parameter is the wind speed W, and the wind speed value ranges from 0 to 24 m / s with an interval of 4 m / s. The distribution probability density function of the small - wave surface is as follows:
[0090]
[0091] where, e n is the outer normal direction vector of the small - wave surface; μ n is the cosine of the angle between the small - wave surface normal and the Z - axis; n is the number of randomly distributed small - wave surfaces on the rough sea surface; σ is the variance of the distribution probability density of the small - wave surface; the relationship between σ and the sea surface wind speed W is as follows:
[0092] σ 2 = 0.003 + 0.00512W
[0093] The ocean - atmosphere coupled vector radiative transfer model PCOART also needs to input the solar irradiance F 0 (i) at the top of the atmosphere according to the spectral response characteristics of each band of the HJ - 2A / B satellite, and calculate F 0 (i) according to step 1. The solar zenith angle θ 0 takes values from 0 to 88° with a step interval of 2°. The observation zenith angle θ v of the hyperspectral imager on the HJ - 2A / B satellite takes the zeros of the 100 - order Legendre polynomial in the interval (0, 1), with an interval of about 1.8°. The Legendre polynomial P m (x) is defined as:
[0094]
[0095] Among them, m represents the order, and x is the zero point of the Legendre polynomial in the interval (0, 1).
[0096] Based on the above input parameters, perform vector radiative transfer simulation, and calculate and obtain the Rayleigh scattering radiation vectors (I, Q, U, V) under different atmospheric molecular optical thicknesses τ r (λ), solar zenith angle θ 0 , observation zenith angle θ v , and wind speed W, and establish a high-precision Rayleigh scattering look-up table for the HJ-2A / B satellite hyperspectral imager. T
[0097] Step 4: Revise the atmospheric molecular Rayleigh optical thickness model with reference to the latest atmospheric correction algorithm adopted by the OBPG of the US Marine Biology Expert Group. The calculation formula for the atmospheric molecular Rayleigh optical thickness at each wavelength under standard atmospheric pressure before revision is as follows:
[0098]
[0099] The calculation formula for the atmospheric molecular Rayleigh optical thickness at each wavelength under the revised standard atmospheric pressure is as follows:
[0100]
[0101] Step 5: Considering the difference between the atmospheric pressure at the actual imaging time of the HJ-2A / B satellite hyperspectral imager and the standard atmospheric pressure, the traditional water color satellite look-up table usually calculates the Rayleigh scattering radiance using the standard atmosphere, and then corrects the Rayleigh scattering radiance according to the ratio of the actual sea surface atmospheric pressure to the standard atmospheric pressure. When there is a large difference between the actual sea surface atmospheric pressure and the standard atmospheric pressure, the calculation error is large. Therefore, use the actual atmospheric pressure to revise the atmospheric Rayleigh scattering optical thickness:
[0102]
[0103] Among them, τ r (λ i ) is the atmospheric Rayleigh optical thickness under standard atmospheric pressure, and P r and P 0 are the actual sea surface atmospheric pressure and the standard atmospheric pressure respectively.
[0104] Step 6: Look up the atmospheric Rayleigh scattering radiance value in the high-precision Rayleigh scattering look-up table established in Step 3 according to the atmospheric conditions at the imaging time of the HJ2A / B satellite hyperspectral imager. Since the atmospheric conditions vary significantly in space and time when the HJ2A / B satellite hyperspectral imager images in different regions of the world, the most significant one is the influence of atmospheric pressure on atmospheric molecular Rayleigh scattering.
[0105] First, calculate the Rayleigh optical thickness τ of atmospheric molecules at each wavelength under standard atmospheric pressure according to step 4 r (λ i ); then, according to the atmospheric pressure data at the satellite overpass time and step 5, correct the optical thickness τ r (λ i ) of atmospheric molecules for air pressure to obtain the Rayleigh optical thickness τ(λ i ) of atmospheric molecules at each wavelength under the actual sea surface atmospheric pressure; subsequently, calculate the equivalent Rayleigh scattering optical thickness τ r (i) of each band of the HJ-2A / B satellite hyperspectral imager according to step 1
[0106] Secondly, determine the solar zenith angle θ 0 at the overpass time of the HJ-2A / B satellite hyperspectral imager and the observation zenith angle θ v of each pixel, based on the on-site wind speed data W
[0107] Finally, look up and calculate the Rayleigh scattering radiance value of atmospheric molecules in the high-precision Rayleigh scattering look-up table of the HJ-2A / B satellite hyperspectral imager constructed in step 3 to achieve high-precision look-up calculation of the Rayleigh scattering radiation vector (I, Q, U, V) T of the HJ-2A / B satellite hyperspectral imager image
[0108] By using the construction and look-up method of the present invention, a spatial geometric distribution map of Rayleigh scattering radiance, an error spatial geometric distribution map, and a graph of the variation of Rayleigh scattering radiance error with bands of the HJ-2A / B satellite hyperspectral imager are obtained, and the results of the effect display are as follows
[0109] Figure 2 Shown is a spatial geometric distribution map of Rayleigh scattering radiance of the HJ-2A / B satellite hyperspectral imager Figure 2 (a), (d), (g) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488 nm, 531 nm, and 670 nm when the solar zenith angle is 20°, (b), (e), (h) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488 nm, 531 nm, and 670 nm when the solar zenith angle is 40°, and (c), (f), (i) are respectively the spatial geometric distribution maps of Rayleigh scattering radiance at 488 nm, 531 nm, and 670 nm when the solar zenith angle is 60°
[0110] Figure 3 Shown is a spatial geometric distribution map of Rayleigh scattering radiance error of the HJ-2A / B satellite hyperspectral imager Figure 3(a), (d), and (g) are the spatial geometric distribution diagrams of Rayleigh scattering radiance errors at 488 nm, 531 nm, and 670 nm respectively when the solar zenith angle is 20°. (b), (e), and (h) are the spatial geometric distribution diagrams of Rayleigh scattering radiance errors at 488 nm, 531 nm, and 670 nm respectively when the solar zenith angle is 40°. (c), (f), and (i) are the spatial geometric distribution diagrams of Rayleigh scattering radiance errors at 488 nm, 531 nm, and 670 nm respectively when the solar zenith angle is 20°.
[0111] It can be easily observed from the figure that, taking 488 nm, 531 nm, and 670 nm as examples, the Rayleigh scattering radiance error of the HJ-2A / B satellite hyperspectral imager is less than 0.5%, and good calculation accuracy can still be obtained under large sensor observation geometries. In addition, as the solar zenith angle gradually increases (from 20° to 60°), the relative error of the Rayleigh scattering look-up table of the HJ-2A / B satellite hyperspectral imager constructed by the present invention is still less than 0.5%.
[0112] Figure 4 It is a diagram showing the variation of Rayleigh scattering radiance error of the HJ-2A / B satellite hyperspectral imager with the wavelength band. Figure 4 (a) is the diagram showing the variation of Rayleigh scattering radiance error with the wavelength band when the solar zenith angle is 0° and the azimuth angles are 0°, 45°, 90°, 135°, and 180°. (b) is the diagram showing the variation of Rayleigh scattering radiance error with the wavelength band when the solar zenith angle is 20° and the azimuth angles are 0°, 45°, 90°, 135°, and 180°. (c) is the diagram showing the variation of Rayleigh scattering radiance error with the wavelength band when the solar zenith angle is 40° and the azimuth angles are 0°, 45°, 90°, 135°, and 180°. (d) is the diagram showing the variation of Rayleigh scattering radiance error with the wavelength band when the solar zenith angle is 60° and the azimuth angles are 0°, 45°, 90°, 135°, and 180°. Under different solar zenith angles, the relative error of each wavelength band of the Rayleigh scattering look-up table of the HJ-2A / B satellite hyperspectral imager constructed by the present invention is better than 0.2%, meeting the requirements of the HJ-2A / B satellite hyperspectral imager for operational ocean remote sensing applications.
[0113] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit them. Although the present invention has been described in detail with reference to the foregoing embodiments, those of ordinary skill in the art should understand that they can still modify the technical solutions described in the foregoing embodiments, or perform equivalent replacements on some or all of the technical features; and these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the present invention.
Claims
1. Construction and lookup method of Rayleigh scattering lookup table for environmental disaster reduction hyperspectral satellite, Characterized in that, It includes the following steps: Step 1: Calculate the solar irradiance F at the top of the atmosphere equivalent to each band of the HJ-2A / B satellite according to the spectral response function of the hyperspectral imager on the HJ-2A / B satellite 0 (i) The Rayleigh scattering optical depth τ of atmospheric molecules r (i) and the specific absorption coefficient A of ozone per unit volume concentration OZ (i); Step 2: Calculate the spherical atmospheric molecular scattering phase matrix P(Θ) of the HJ-2A / B satellite hyperspectral imager based on the anisotropic band difference and depolarization effect of atmospheric molecules; Step 3: Combine the parameter values obtained in Step 1 and Step 2, and numerically solve the radiative transfer equation based on the ocean-atmosphere coupled vector radiative transfer model to calculate the Rayleigh scattering radiative vector (I, Q, U, V) received by the satellite at the top of the atmosphere under different optical thicknesses τ r (λ) of atmospheric molecules, solar zenith angle θ 0 , observation zenith angle θ v , and wind speed W, and construct a Rayleigh scattering look-up table applicable to the environmental disaster reduction satellite imager; T Wherein, I is the total radiation intensity, Q is the linearly polarized radiation intensity in the 45° direction, U is the linearly polarized radiation intensity in the vertical direction, and V is the circularly polarized light intensity; E l and E r are the electric vectors of the parallel and perpendicular reference planes respectively; and are the conjugate matrices of E l and E r respectively; specifically: 3.
1. Set the input parameters of the ocean-atmosphere coupled vector radiative transfer model; The ocean-atmosphere coupled vector radiative transfer model includes the coupling of atmospheric stratification, sea-air interface stratification and ocean stratification; 3.1.
1. Set the atmosphere to be stratified into N layers, each layer consisting of pure atmospheric molecules. The input parameters include the optical thickness τ of atmospheric molecules r (λ), the scattering phase matrix P(Θ), and the single-scattering albedo ω of the atmosphere, where N ≥ 10, and the optical thickness τ of atmospheric molecules r (λ) is numerically input as 0.002 - 0.6 with a step size of 0.002; ω = 1.0; 3.1.
2. Set the sea-air interface stratification as the probability density function of the distribution of the wavelet surface constructed based on the statistical results, where the input parameter is the wind speed W, and the wind speed value ranges from 0 to 24 m / s, with an interval of 4 m / s; The distribution probability density function p(e n ) is as follows: Among them, e n is the outer normal direction vector of the wavelet surface; μ n is the cosine of the angle between the wavelet surface normal and the Z-axis; n is the number of wavelet surfaces randomly distributed on the rough sea surface; σ is the variance of the probability density of the wavelet surface distribution; σ 2 = 0.003 + 0.00512W; 3.1.
3. Set the ocean stratification to M layers. The inherent optical property parameters of the seawater layer include the absorption coefficient a water , the scattering phase matrix, and the single-scattering albedo ω water , where M ≥ 1, a water ≥ 10, ω water = 1, and the ocean depth is infinitely deep; 3.
2. Conduct vector radiative transfer simulation to calculate the Rayleigh scattering radiation vectors (I, Q, U, V) under different optical thicknesses τ r (λ) of atmospheric molecules, solar zenith angles θ 0 , observation zenith angles θ v , and wind speeds W, and establish a high-precision Rayleigh scattering look-up table for the HJ-2A / B satellite hyperspectral imager; T Step 4: Calculate the equivalent Rayleigh scattering optical depth τ of atmospheric molecules based on the atmospheric Rayleigh scattering optical depth under actual atmospheric pressure. r (i) Based on the Rayleigh scattering lookup table constructed in Step 3, complete the lookup calculation of Rayleigh scattering applicable to the environmental disaster reduction satellite imager under different atmospheric conditions.
2. The construction and lookup method of the Rayleigh scattering lookup table for the environmental disaster reduction hyperspectral satellite according to claim 1, Characterized in that: In Step 1, the solar irradiance F incident on the top of the atmosphere equivalent to each band of the HJ-2A / B satellite 0 (i), the Rayleigh scattering optical depth τ of atmospheric molecules r (i) and the specific absorption coefficient A of ozone per unit volume concentration OZ (i), the calculation formula is as follows: where \(i\) is the band number of the hyperspectral imager on HJ-2A / B satellite, and \(1\leq i\leq100\); \(\lambda\) i is the wavelength corresponding to the \(i\)-th band; \(\lambda\) i1 and \(\lambda\) i2 are the left and right end wavelengths corresponding to 1% spectral responsivity of the \(i\)-th band, respectively, and \(\lambda\) i1 \(\leq\lambda\) i \(\leq\lambda\) i2 ; F 0 (λ i ) is the solar irradiance at the top of the atmosphere corresponding to the i-th band; S(λ i ) is the spectral responsivity corresponding to the i-th band; τ r (λ i ) is the Rayleigh scattering optical depth of atmospheric molecules at the standard atmospheric pressure corresponding to the i-th band; A OZ (λ i ) is the specific absorption coefficient of ozone corresponding to the concentration per unit volume in the i-th band.
3. The construction and lookup method of the Rayleigh scattering lookup table for the environmental disaster reduction hyperspectral satellite according to claim 2, Characterized in that: In step 2, the spherical atmospheric molecular scattering phase matrix P(Θ) of the HJ-2A / B satellite hyperspectral imager is: where ρ(λ i ) is the depolarization factor corresponding to the wavelength λ in the i-th band, taking into account the anisotropic band differences of atmospheric molecules; The partial recession factor ρ(λ i ) is calculated as follows:
4. The construction and lookup method of the Rayleigh scattering lookup table for the environmental disaster reduction hyperspectral satellite according to claim 3, Characterized in that: In step 3, the Rayleigh scattering radiation vector (I, Q, U, V) T is calculated from the Stokes vector S, and the Stokes vector S is as follows: where, E l and E r are the electric vectors parallel and perpendicular to the reference plane, respectively; and are the conjugate matrices of E l and E r respectively.
5. The construction and lookup method of the Rayleigh scattering lookup table for the environmental disaster reduction hyperspectral satellite according to any one of claims 1-4, Characterized in that, Step 4 is specifically: 4.
1. Revise the Rayleigh optical depth of atmospheric molecules based on the atmospheric correction algorithm: 4.
2. Calculate the Rayleigh scattering optical depth of the atmosphere under the actual atmospheric pressure: Where, P r and P 0 are the actual sea surface atmospheric pressure and the standard atmospheric pressure, respectively; 4.
3. Calculate the equivalent atmospheric molecular Rayleigh scattering optical thickness τ under the actual atmospheric pressure r (i): 4.4 Solar zenith angle θ of HJ-2A / B satellite hyperspectral imager at transit time 0 and the observation zenith angle θ of each pixel v , based on the on-site wind speed data W, and then look up and calculate the Rayleigh scattering radiance value of atmospheric molecules in the Rayleigh scattering look-up table constructed in step 3, so as to achieve high-precision look-up calculation of the atmospheric Rayleigh scattering radiation vector (I, Q, U, V) of the HJ-2A / B satellite hyperspectral imager image T of.
6. The construction and lookup method of the Rayleigh scattering lookup table for the environmental disaster reduction hyperspectral satellite according to claim 1, Characterized in that: In step 3.2, the solar zenith angle θ 0 takes values from 0 to 88°, with a step interval of 2°; The observed zenith angle θ v Take the zeros of the Legendre polynomial of order 100 in the interval (0, 1) with an interval of 1.8°; The Legendre polynomial P m (x) is as follows: Where, m represents the order, and x is the zero point of the Legendre polynomial in the interval (0,1).
Citation Information
Patent Citations
Self-adaptive waveband selection method based on hyperspectral water body reservoir
CN111912799A
Characterizing tropospheric boundary layer thermodynamic and refractivity profiles utilizing selected waveband infrared observations
US20210041299A1