A method and system for simulating Rayleigh scattering vector radiation transmission on non-uniform underlying surface

Abstract

The present invention provides a method and system for simulating the radiation transmission of Rayleigh scattering vectors on a non-uniform underlying surface, simulating the radiation transmission process of photons, and distinguishing the different contributions of photons to the top-of-atmosphere radiance after reaching the atmosphere, the land surface, and the water surface. The top-of-atmosphere radiance received by the satellite under preset conditions is calculated by simulating a preset number of photons, and accurate correction of the Rayleigh scattering radiance is achieved during the atmospheric correction process.

Description

Technical Field

[0001] The invention belongs to the field of remote sensing, and in particular relates to a method and system for simulating the transmission of Rayleigh scattering vector radiation of a non-uniform underlying surface. Background Art

[0002] Water color remote sensing can measure the spectral information reflected by the water surface through satellite, aerial or ground sensors, and then study the biochemical characteristics and dynamic changes of water bodies. It has wide applications in water quality monitoring, environmental protection, resource management and ecological research.

[0003] In order to obtain the true reflection characteristics of water bodies, it is necessary to perform atmospheric correction on remote sensing images to remove the influence of the atmosphere on the light signals received by the sensor, that is, to remove the effects of atmospheric scattering and absorption, so that the spectral information obtained from remote sensing images can be as close as possible to the reflection spectrum of the water body itself.

[0004] Accurate atmospheric correction is the basis for quantitative inversion of water body parameters by remote sensing, and the first step in atmospheric correction of water color remote sensing is to deduct the Rayleigh scattering radiance of atmospheric molecules. Currently, Rayleigh scattering correction usually uses a lookup table method based on radiation transfer simulation. Therefore, accurately simulating the Rayleigh scattering radiation transfer process is the prerequisite and basis for accurate atmospheric correction.

[0005] The most commonly used SEADAS accurate Rayleigh scattering lookup table in atmospheric correction of water color remote sensing is established by solving the atmospheric vector radiation transfer equation by the successive scattering method. The lookup table is targeted at a specific remote sensor. In addition, the universal ocean water color remote sensing accurate Rayleigh scattering lookup table generated based on the PCOART model can be applied to the accurate Rayleigh scattering calculation of all water color remote sensors. However, the above models all assume that the underlying surface is a uniform water body, which is not applicable to coastal or inland water bodies with non-uniform underlying surfaces, especially narrow and long river water bodies with land on both sides, because these water bodies are easily affected by the radiation signals of the surrounding land pixels, resulting in a proximity effect.

[0006] Currently, there is a lack of water-land interactive vector radiation transfer models that can autonomously set complex non-uniform underlying surfaces. For the application of water color remote sensing in complex water bodies such as rivers, the accuracy of its Rayleigh scattering correction cannot be guaranteed, which limits the subsequent atmospheric correction and inversion accuracy of water color elements. Summary of the invention

[0007] The present invention provides a method and system for simulating Rayleigh scattering vector radiation transmission on a non-uniform underlying surface, so as to solve the problem in the prior art that there is a lack of a water-land interactive vector radiation transmission model capable of autonomously setting a complex non-uniform underlying surface.

[0008] In order to solve the above technical problems, the embodiments of the present invention disclose the following technical solutions:

[0009] One aspect of the present invention provides a method for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface, which is applied to atmospheric correction of satellite remote sensing images. The method comprises:

[0010] S100, simulating the motion state of a satellite emitting photons to a preset observation point in an atmosphere-surface system according to preset conditions, wherein the preset conditions at least include the position relationship of land / water body and the height of the atmosphere at the observation point, as well as the initial weight, position, motion direction and polarization state of the photons;

[0011] S200, randomly sampling the free path of the photon, and determining the position state of the photon after the photon runs according to the current movement direction and free path of the photon, wherein the position state includes the interior of the atmosphere, the land surface, and the water surface;

[0012] S300, updating the polarization state and weight of the photon using a preset updating method corresponding to the position state;

[0013] S400, calculating the contribution value of the photon to the top-of-atmosphere radiance using a corresponding preset calculation method according to the position state of the photon;

[0014] S500, determining whether the current weight of the photon is lower than a preset threshold, if so, executing step S600; if not, executing step S700;

[0015] S600, determining that the photon has disappeared, and judging whether the number of photons currently simulated has reached a preset total number, if so, executing step S800; if not, re-executing step S100 and subsequent steps to simulate the motion state of the next photon emitted by the satellite;

[0016] S700, updating the movement direction of the photon; re-executing S200 and subsequent steps to simulate the next movement trajectory of the photon until the photon disappears;

[0017] S800, calculating the average value of all photon contribution values ​​as the top-of-atmosphere radiance received by the satellite under the preset conditions, so as to be used for correcting the Rayleigh scattering radiance in the atmospheric correction process.

[0018] Optionally, the step S100 includes:

[0019] Step S101, taking the preset observation point as the origin, the normal direction perpendicular to the surface where the observation point is located as the z-axis, the projection direction of the sun on the xoy plane as the positive direction of the x-axis, and the y-axis determined by the right-hand screw rule, to construct a spatial coordinate system for the movement of photons;

[0020] Step S102, obtaining geometric conditions for sun and satellite observation and atmospheric height, wherein the geometric conditions at least include an observed solar zenith angle, an observed zenith angle, and a relative azimuth angle;

[0021] Step S103, determining the initial movement direction and position of the photon according to the geometric conditions, comprising:

[0022] The satellite's sensor is located at the top of the atmosphere, and the initial direction of motion of the photon is calculated according to the following formula: and location :

[0023]

[0024]

[0025] in, is the observation zenith angle, is the relative azimuth, is the atmospheric height;

[0026] Step S104: Set the initial weight of the photon to 1, and set the polarization state matrix to the unit matrix.

[0027] Optionally, the step S200 includes:

[0028] S201, randomly sampling to obtain the free path of photons;

[0029] S202. After the photon runs along the free path, the position of the photon is calculated according to the following formula:

[0030]

[0031] in, is the free path of the photon;

[0032] S203, determining whether the photon is located in the atmosphere according to the position of the photon, if so, determining that the current position state of the photon is inside the atmosphere; if not, executing step S204;

[0033] S204. Determine whether the photon reaches the water surface or the land surface after running the free path according to the land / water position relationship, and determine that the current position state of the photon is the water surface when the photon reaches the water surface, and determine that the current position state of the photon is the land surface when the photon reaches the land surface.

[0034] Optionally, the step S300 includes:

[0035] When the position state of the photon is inside the atmosphere, the polarization state and weight of the photon are updated according to the Rayleigh scattering phase function;

[0036] When the position state of the photon is the land surface, the weight of the photon is updated according to the Lambertian reflection matrix, and the polarization state is set to unbiased;

[0037] When the position state of the photon is the surface of the water body, the polarization state and weight of the photon are updated according to the mirror reflection matrix.

[0038] Optionally, the step S400 includes:

[0039] When the photon's position state is inside the atmosphere, the photon contribution value is calculated based on the Rayleigh scattering phase function and the scattering angle;

[0040] When the photon's position state is the land surface, the photon contribution value is calculated based on the surface albedo and reflection angle;

[0041] When the position state of the photon is the water surface, the photon contribution value is calculated based on the Rayleigh scattering phase function and the specular reflection matrix.

[0042] Optionally, the step S700 includes:

[0043] When the position state of the photon is inside the atmosphere, it is determined that the photon is scattered, and a random sampling method and a rejection-acceptance method are used to obtain a scattering angle, a random sampling method is used to obtain a scattering azimuth, and the next movement direction of the photon is calculated based on the scattering angle and the scattering azimuth;

[0044] When the position state of the photon is the land surface, it is determined that the photon has Lambertian reflection, and the reflection angle and azimuth of the photon are obtained by a random sampling method, and the next movement direction of the photon is calculated based on the reflection angle and azimuth;

[0045] When the position state of the photon is the surface of the water body, it is determined that the photon undergoes mirror reflection, the reflection angle of the photon is obtained according to the incident angle of the mirror reflection, the azimuth of the photon is obtained by a random sampling method, and the next movement direction of the photon is calculated based on the reflection angle and the azimuth.

[0046] Optionally, when the photon's position state is inside the atmosphere, the photon's polarization state and weight are updated as follows:

[0047] Photon weight According to Rayleigh scattering single scattering albedo Attenuation, preset , photon polarization state matrix According to the scattering matrix Make changes:

[0048]

[0049]

[0050] in, is the photon weight, is the polarization state matrix, Calculated by the following formula:

[0051]

[0052]

[0053]

[0054] in, is the Rayleigh scattering phase matrix, is the rotation matrix, is the preset depolarization factor, are the scattering angle and the reference plane rotation angle when scattering occurs, which are calculated by the meridian plane method. The angles between the incident meridian plane, the exit meridian plane and the scattering plane are and .

[0055] Optionally, when the position state of the photon is the land surface, the polarization state and weight of the photon are updated as follows:

[0056] The weight of the photons is based on the preset land albedo Attenuation, photon polarization state matrix According to the reflection matrix Make the corresponding changes:

[0057]

[0058]

[0059] in, is the land albedo, A is The value of Represents the pre-acquired cosine value of the incident angle.

[0060] Optionally, when the position state of the photon is the surface of the water body, the polarization state and weight of the photon are updated as follows:

[0061] The weight of the photon is based on the preset specular albedo Attenuation, photon polarization state matrix According to the reflection matrix Make changes:

[0062]

[0063]

[0064]

[0065] in,

[0066]

[0067] in, Represents the preset incident angle between the incident light and the water surface normal vector, Represents the refractive index of water, with a value of 1.34.

[0068] Another aspect of the present invention provides a non-uniform underlying surface Rayleigh scattering vector radiation transmission simulation system, which is applied to the atmospheric correction of satellite remote sensing images. The system comprises:

[0069] The photon simulation module is configured to simulate the motion state of a satellite emitting photons to a preset observation point in the atmosphere-surface system according to preset conditions, wherein the preset conditions at least include the position relationship of the land / water body and the atmospheric height at the observation point, as well as the initial weight, position, motion direction and polarization state of the photons;

[0070] A position determination module is configured to randomly sample the free path of the photon, and determine the position state of the photon after the photon runs according to the current movement direction and free path of the photon, wherein the position state includes the interior of the atmosphere, the surface of the land, and the surface of the water body;

[0071] An updating module is configured to update the polarization state and weight of the photon using a preset updating method corresponding to the position state;

[0072] A contribution value calculation module is configured to calculate the contribution value of the photon to the top-of-atmosphere radiance using a corresponding preset calculation method according to the position state of the photon;

[0073] The judgment and radiance determination module is configured to determine whether the current weight of the photon is lower than a preset threshold,

[0074] When the current weight of the photon is lower than the preset threshold, it is determined that the photon has died, and it is judged whether the number of photons currently simulated has reached the preset total number; when the number of photons currently simulated has not reached the preset total number, the photon simulation module is used to simulate the motion state of the next photon emitted by the satellite; when the number of photons currently simulated has reached the preset total number, the average value of all photon contribution values ​​is calculated as the top-of-atmosphere radiance received by the satellite under the preset conditions, so as to be used for the correction of Rayleigh scattering radiance in the atmospheric correction process;

[0075] When the current weight of the photon is not lower than a preset threshold, the movement direction of the photon is updated, and the position determination module is used to simulate the next movement trajectory of the photon until the photon disappears.

[0076] The method and system for simulating the radiation transmission of Rayleigh scattering vectors on non-uniform underlying surfaces disclosed in the present invention simulate the radiation transmission process of photons and distinguish the different effects of photons on the top of the atmosphere radiance after reaching the atmosphere, the land surface and the water surface. By simulating a preset number of photons, the top of the atmosphere radiance received by the satellite under preset conditions is counted, and accurate correction of the Rayleigh scattering radiance is achieved during the atmospheric correction process.

[0077] This Summary is provided to introduce a selection of concepts in a simplified form that are further described below in the Detailed Description. This Summary is not intended to identify key features or essential features of the disclosure, nor is it intended to limit the scope of the disclosure. BRIEF DESCRIPTION OF THE DRAWINGS

[0078] The above and other objects, features and advantages of the present disclosure will become more apparent through a more detailed description of exemplary embodiments of the present disclosure in conjunction with the accompanying drawings, wherein like reference numerals generally represent like components throughout the exemplary embodiments of the present disclosure.

[0079] Figure 1 A schematic diagram of a flow chart of a method for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface disclosed in an embodiment of the present invention;

[0080] Figure 2 An implementation disclosed in an embodiment of the present invention Figure 1 Schematic diagram of the process of step S100;

[0081] Figure 3 An implementation disclosed in an embodiment of the present invention Figure 1 Schematic diagram of the process of step S200;

[0082] Figure 4 The present invention is a schematic structural diagram of a system for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface disclosed in an embodiment of the present invention. DETAILED DESCRIPTION

[0083] Embodiments of the present disclosure will be described in more detail below with reference to the accompanying drawings. Although embodiments of the present disclosure are shown in the accompanying drawings, it should be understood that the present disclosure can be implemented in various forms and should not be limited by the embodiments set forth herein. On the contrary, these embodiments are provided to make the present disclosure more thorough and complete, and to fully convey the scope of the present disclosure to those skilled in the art.

[0084] As used herein, the term "including" and its variations mean open inclusion, i.e., "including but not limited to". Unless otherwise stated, the term "or" means "and / or". The term "based on" means "based at least in part on". The terms "an example embodiment" and "an embodiment" mean "at least one example embodiment". The term "another embodiment" means "at least one additional embodiment". The terms "first", "second", etc. may refer to different or the same objects. Other explicit and implicit definitions may also be included below.

[0085] Figure 1 The present invention discloses a method for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface, which is applied to atmospheric correction of satellite remote sensing images. Figure 1 As shown, the method comprises the following steps:

[0086] S100, simulating the motion state of a satellite emitting photons to a preset observation point in an atmosphere-surface system according to preset conditions.

[0087] The entire simulation process is initialized according to preset conditions, which at least include the land / water position relationship and the atmospheric height at the observation point, as well as the initial weight, initial position, movement direction and polarization state of the photon. Among them, the land / water position relationship is the distribution of land and water on the earth's surface, which can be set by the user; the initial position of the photon is the position of the satellite sensor. In the disclosed embodiment of the present invention, the satellite sensor is located at the top of the atmosphere.

[0088] According to the above preset conditions, the motion state of a photon is simulated when it is emitted from a predetermined satellite to a preset observation point in the atmosphere-surface system. And after determining that the photon has disappeared in the process of radiation transmission, the motion state of another photon is re-simulated in the above manner until the number of simulated photons reaches a preset total number, for example, 1000.

[0089] In one embodiment disclosed in the present invention, Figure 2 As shown, the simulation process can be initialized in the following ways:

[0090] Step S101, taking the observation point as the origin, the normal direction perpendicular to the surface where the observation point is located as the z-axis, the projection direction of the sun on the xoy plane as the positive direction of the x-axis, and the y-axis determined by the right-hand screw rule, to construct a spatial coordinate system for the movement of photons.

[0091] The following is a description of the spatial coordinate system:

[0092] Origin (o): The location of the observation point.

[0093] z-axis: perpendicular to the normal direction of the surface where the observation point is located, the positive direction is pointing to the upper atmosphere.

[0094] x-axis: the projection direction of the sun on the xoy plane, the projection direction pointing to the sun is the positive direction.

[0095] y-axis: Perpendicular to the x-axis and z-axis, the positive direction of the y-axis can be determined according to the right-hand screw rule.

[0096] Step S102: Obtain the geometric conditions for solar and satellite observations and the altitude of the atmosphere.

[0097] The geometric conditions at least include the observed solar zenith angle, the observed zenith angle and the relative azimuth angle.

[0098] Step S103: determine the initial movement direction and position of the photon according to geometric conditions.

[0099] In the embodiments disclosed in the present invention, the initial movement direction and position of the photon can be determined in the following manner:

[0100] The satellite's sensor is located at the top of the atmosphere, and the initial direction of motion of the photon is calculated according to the following formula: and initial position :

[0101]

[0102]

[0103] in, is the observation zenith angle, is the relative azimuth, is the atmospheric height;

[0104] Step S104: Set the initial weight of the photon to 1, and set the polarization state matrix to the unit matrix.

[0105] The weight of a photon represents the energy carried by the photon and is used to track the energy changes of the photon during transmission. Setting the initial weight to 1 means that each photon is assigned a standard, unit energy when it is emitted. As the photon propagates in the atmosphere and undergoes scattering or reflection events, its weight will change, representing the increase, decrease or loss of the photon's energy.

[0106] S200, randomly sampling the free path of the photon, and determining the position state of the photon after running according to the current moving direction of the photon and the free path of the photon.

[0107] In the embodiments disclosed in the present invention, the position states include the interior of the atmosphere, the surface of land, and the surface of a water body, wherein the interior of the atmosphere represents that the photon is located in the atmosphere; the surface of land represents that the photon reaches the surface of land; and the surface of a water body represents that the photon reaches the surface of the water body.

[0108] In one embodiment disclosed in the present invention, Figure 3 As shown, step S200 may be implemented in the following manner:

[0109] S201, randomly sampling to obtain the free path of photons.

[0110] The free path is the distance a photon travels when propagating in the atmosphere without scattering or absorption, and is related to the physical properties of the atmosphere (such as density and scattering coefficient).

[0111] The photon free path can be obtained by random sampling in the following way:

[0112] 1. Generate a random number R uniformly distributed between (0,1);

[0113] 2. Use the following formula to calculate the free path:

[0114]

[0115] in, is the free path of the photon; β is the preset extinction coefficient of the atmosphere, which indicates the probability of scattering or absorption occurring per unit distance.

[0116] S202. After the photon runs along the free path, calculate the position of the photon.

[0117] The photon position is calculated as follows:

[0118]

[0119] in, represents the free path of the photon.

[0120] S203. Determine whether the photon is in the atmosphere according to the position of the photon.

[0121] According to the z value of the photon in the spatial coordinate system, that is, the height of the photon from the ground, it is judged whether the photon is running in the atmosphere or has reached the surface of the earth. If the z value is greater than 0, it means that the photon is running in the atmosphere. If the z value is less than or equal to 0, it means that the photon has reached the surface of the earth.

[0122] If the photon is in the atmosphere, it is determined that the current position state of the photon is inside the atmosphere; if not, it is determined that the photon has reached the surface of the earth, and step S204 is executed.

[0123] S204. Determine whether the photon reaches the water surface or the land surface after running the free path according to the land / water position relationship, and determine that the current position state of the photon is the water surface when the photon reaches the water surface, and determine that the current position state of the photon is the land surface when the photon reaches the land surface.

[0124] If the photon has reached the surface of the earth, it can be determined whether the photon is on the land surface or the water surface based on the position of the photon on the surface.

[0125] In one embodiment disclosed in the present invention, the land / water position relationship function can be preset , judged in the following way:

[0126]

[0127] in, is the land / water position relationship, is the position of the photon corresponding to Axis and The coordinate values ​​of the axis.

[0128] In another embodiment disclosed in the present invention, the coordinates of the land and water body contours can be obtained in advance to construct the land / water body position relationship. and The coordinate value of the photon determines whether the photon reaches the land surface or the water surface.

[0129] In the disclosed embodiment of the present invention, the land / water position relationship can be preset according to the actual needs of the user.

[0130] For example, suppose the underlying surface is a river-land interaction surface, and the river width is km, the land / water position relationship can be expressed as:

[0131]

[0132] in, is the position of the photon corresponding to The coordinate values ​​of the axis.

[0133] S300, updating the polarization state and weight of the photon using a preset updating method corresponding to the position state.

[0134] Different location states correspond to different update methods.

[0135] In one embodiment disclosed in the present invention, step S300 may be completed in the following manner:

[0136] When the position state of the photon is inside the atmosphere, scattering occurs, and the polarization state and weight of the photon are updated according to the Rayleigh scattering phase function.

[0137] Photon weight According to Rayleigh scattering single scattering albedo Attenuation, preset , photon polarization state matrix According to the scattering matrix Make changes:

[0138]

[0139]

[0140] in, is the photon weight, is the polarization state matrix, Calculated by the following formula:

[0141]

[0142]

[0143]

[0144] in, is the Rayleigh scattering phase matrix, is the rotation matrix, is the preset depolarization factor, and i are the scattering angle and the reference plane rotation angle when scattering occurs, respectively, which are calculated by the meridian plane method.

[0145] The incident light before the photon is scattered, and the outgoing light after the photon is scattered. The plane formed by the two is the scattering surface, and the angle between the two is the scattering angle The incident meridian plane is the plane formed by the incident light and the z-axis, and the exit meridian plane is the plane formed by the exit light and the z-axis. The angles between the two meridian planes and the scattering plane are and .

[0146] If the photon reaches the surface of the earth, corresponding reflection occurs depending on the surface type (land / water).

[0147] When the position state of the photon is the land surface, Lambertian reflection occurs, the weight of the photon is updated according to the Lambertian reflection matrix, and the polarization state is set to unbiased.

[0148] The weight of the photons is based on the preset land albedo Attenuation, photon polarization state matrix According to the reflection matrix Make the corresponding changes:

[0149]

[0150]

[0151] in, is the land albedo, A is The value of Represents the pre-acquired cosine value of the incident angle.

[0152] When the position state of the photon is the surface of the water body, mirror reflection occurs, and the polarization state and weight of the photon are updated according to the mirror reflection matrix.

[0153] The weight of the photon is based on the preset specular albedo Attenuation, photon polarization state matrix According to the reflection matrix Make changes:

[0154]

[0155]

[0156]

[0157] in,

[0158]

[0159] in, Represents the preset incident angle between the incident light and the water surface normal vector, Represents the refractive index of water, with a value of 1.34.

[0160] S400, according to the position state of the photon, using a corresponding preset calculation method to calculate the contribution value of the photon to the top-of-atmosphere radiance.

[0161] Different position states correspond to different calculation methods.

[0162] In one embodiment disclosed in the present invention, step S400 may be implemented in the following manner:

[0163] When the photon is in the atmosphere, scattering occurs, and the part scattered toward the sun contributes to the radiance at the top of the atmosphere. The photon contribution is calculated based on the Rayleigh scattering phase function and the scattering angle.

[0164] When the photon's position state is the land surface, the part of the Lambertian reflection toward the sun contributes to the top-of-atmosphere radiance, and the photon contribution value is calculated based on the surface albedo and reflection angle;

[0165] When the photon's position state is the water surface, the part scattered toward the water surface and then reflected by the mirror toward the sun contributes to the top-of-atmosphere radiance. The photon contribution value is calculated based on the Rayleigh scattering phase function and the mirror reflection matrix.

[0166] S500: Determine whether the current weight of the photon is lower than a preset threshold.

[0167] If the current weight of the photon is lower than the preset threshold, execute step S600;

[0168] If the current weight of the photon is not lower than the preset threshold, it is necessary to determine the direction of movement of the photon during the next movement, and execute step S700.

[0169] S600: Determine whether the photon has disappeared, and judge whether the number of simulated photons reaches a preset total number.

[0170] If the number of photons currently simulated reaches the preset total number, step S800 is executed.

[0171] If the number of photons currently simulated does not reach the preset total number, step S100 and subsequent steps are re-executed to re-simulate the motion state of the next photon emitted by the satellite.

[0172] S700: Update the moving direction of the photon.

[0173] In one embodiment disclosed in the present invention, step S700 may be implemented in the following manner:

[0174] When the position state of the photon is inside the atmosphere, it is determined that the photon is scattered, and the scattering angle is obtained by using the random sampling method and the rejection-acceptance method, the scattering azimuth is obtained by using the random sampling method, and the next movement direction of the photon is calculated based on the scattering angle and the scattering azimuth.

[0175] The following is a further explanation of the above content:

[0176] First, generate a random number uniformly distributed in [-1, 1] As , then, randomly sample from [0, 1] using a uniform distribution If satisfied ,but As new Perform subsequent calculations; otherwise, resample until the conditions are met. Finally, the random sampling method is used to determine the scattering azimuth , is a random number uniformly distributed between [0, 1]. The next direction of motion of the photon is determined based on the scattering angle and scattering azimuth.

[0177] When the position state of the photon is the land surface, it is determined that the photon undergoes Lambertian reflection, and a random sampling method is used to obtain the reflection angle and azimuth of the photon, and the next movement direction of the photon is calculated based on the reflection angle and azimuth.

[0178] When the position state of the photon is the surface of the water body, it is determined that the photon undergoes mirror reflection, the reflection angle of the photon is obtained according to the incident angle of the mirror reflection, the azimuth of the photon is obtained by a random sampling method, and the next movement direction of the photon is calculated based on the reflection angle and the azimuth.

[0179] After determining the next movement direction of the photon, S200 and subsequent steps are re-executed to simulate the next movement trajectory of the photon until the photon disappears.

[0180] S800, calculating the average value of all photon contribution values ​​as the top-of-atmosphere radiance received by the satellite under preset conditions, so as to be used for correcting the Rayleigh scattering radiance in the atmospheric correction process.

[0181] Atmospheric correction is the process of removing atmospheric effects from satellite observation data to obtain the true value of surface reflectance or radiance. The average value of the contribution of all simulated photons is calculated as the top-of-atmosphere radiance received by the satellite under the aforementioned preset conditions, and the Rayleigh scattering radiance in the atmospheric correction process is corrected based on the top-of-atmosphere radiance.

[0182] Figure 4 This is a schematic diagram of the structure of a non-uniform underlying surface Rayleigh scattering vector radiation transmission simulation system disclosed in an embodiment of the present invention, which is applied to atmospheric correction. Figure 4 As shown, the system includes the following modules:

[0183] The photon simulation module 1 is configured to simulate the trajectory of a satellite emitting photons to a preset observation point in the atmosphere-surface system according to preset conditions. The preset conditions include at least the land / water position relationship and the atmospheric height at the observation point, as well as the initial weight, position, movement direction and polarization state of the photons.

[0184] The position determination module 2 is configured to randomly sample the free path of the photon and determine the position state of the photon after running according to the current movement direction and free path of the photon. The position state includes the interior of the atmosphere, the land surface and the water surface.

[0185] The updating module 3 is configured to update the polarization state and weight of the photon using a preset updating method corresponding to the position state.

[0186] The contribution value calculation module 4 is configured to calculate the contribution value of the photon to the top-of-atmosphere radiance according to the position state of the photon using a corresponding preset calculation method.

[0187] The judging and radiance determination module 5 is configured to judge whether the current weight of the photon is lower than a preset threshold value,

[0188] When the current weight of a photon is lower than a preset threshold, it is determined that the photon has disappeared, and it is determined whether the number of photons currently simulated has reached the preset total number; when the number of photons currently simulated has not reached the preset total number, a photon simulation module is used to simulate the motion state of the next photon emitted by the satellite; when the number of photons currently simulated has reached the preset total number, the average value of all photon contribution values ​​is calculated as the top-of-atmosphere radiance received by the satellite under preset conditions, so as to be used for the correction of Rayleigh scattering radiance during the atmospheric correction process;

[0189] When the current weight of the photon is not lower than the preset threshold, the direction of movement of the photon is updated, and the position determination module is used to simulate the next movement trajectory of the photon until the photon disappears.

[0190] The embodiments of the present disclosure have been described above, and the above description is exemplary, not exhaustive, and is not limited to the disclosed embodiments. Many modifications and variations will be apparent to those of ordinary skill in the art without departing from the scope and spirit of the described embodiments. The selection of terms used herein is intended to best explain the principles of the embodiments, practical applications, or technical improvements to the technology in the market, or to enable other persons of ordinary skill in the art to understand the embodiments disclosed herein.

Claims

1. A method for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface, applied to atmospheric correction of satellite remote sensing images, characterized in that: The radiation transfer simulation method comprises: S100, simulating the motion state of a satellite emitting photons to a preset observation point in an atmosphere-surface system according to preset conditions, wherein the preset conditions at least include the position relationship of land / water body and the height of the atmosphere at the observation point, as well as the initial weight, position, motion direction and polarization state of the photons; S200, randomly sampling the free path of the photon, and determining the position state of the photon after the photon runs according to the current movement direction and free path of the photon, wherein the position state includes the interior of the atmosphere, the land surface, and the water surface; S300, updating the polarization state and weight of the photon using a preset updating method corresponding to the position state, including: When the position state of the photon is inside the atmosphere, the polarization state and weight of the photon are updated according to the Rayleigh scattering phase function; the polarization state and weight of the photon are updated in the following manner: The weight of the photon is attenuated according to the Rayleigh scattering single scattering albedo ω0, which is preset to ω0 = 1. The photon polarization state matrix pmc is based on the scattering matrix F s Make changes: weight=weight*ω0# pmc=pmc*F s # Among them, weight is the photon weight, pmc is the polarization state matrix, F s Calculated by the following formula: F s =R(π-i2)P(Θ)R(-i1),# Where P(θ) is the Rayleigh scattering phase matrix, R(i) is the rotation matrix, ρ is the preset depolarization factor, θ and i are the scattering angle and the reference plane rotation angle when scattering occurs, respectively, which are calculated by the meridian plane method, and the angles between the incident meridian plane, the exit meridian plane and the scattering plane are i1 and i2 respectively; When the position state of the photon is the land surface, the weight of the photon is updated according to the Lambertian reflection matrix, and the polarization state is set to unbiased; the polarization state and weight of the photon are updated as follows: The weight of the photons is based on the preset land albedo Attenuation is performed, and the photon polarization state matrix pmc is based on the reflection matrix F rl Make the corresponding changes: in, is the land albedo, A is , μ represents the pre-acquired cosine value of the incident angle; When the position state of the photon is the surface of the water body, the polarization state and weight of the photon are updated according to the mirror reflection matrix; the polarization state and weight of the photon are updated in the following manner: The weight of the photon is based on the preset specular albedo Attenuation is performed, and the photon polarization state matrix pmc is based on the reflection matrix F rw Make changes: in, Among them, α i Represents the preset incident angle between the incident light and the water surface normal vector, n w Represents the refractive index of water, with a value of 1.34; S400, according to the position state of the photon, using a corresponding preset calculation method to calculate the contribution value of the photon to the top-of-atmosphere radiance, including: When the photon's position state is inside the atmosphere, the photon contribution value is calculated based on the Rayleigh scattering phase function and the scattering angle; When the photon's position state is the land surface, the photon contribution value is calculated based on the surface albedo and reflection angle; When the position state of the photon is the surface of the water body, the photon contribution value is calculated based on the Rayleigh scattering phase function and the specular reflection matrix; S500, determining whether the current weight of the photon is lower than a preset threshold, if so, executing step S600; if not, executing step S700; S600, determining that the photon has disappeared, and judging whether the number of photons currently simulated has reached a preset total number, if so, executing step S800; if not, re-executing step S100 and subsequent steps to simulate the motion state of the next photon emitted by the satellite; S700, updating the movement direction of the photon; re-executing S200 and subsequent steps to simulate the next movement trajectory of the photon until the photon disappears; S800, calculating the average value of all photon contribution values ​​as the top-of-atmosphere radiance received by the satellite under the preset conditions, so as to be used for correcting the Rayleigh scattering radiance in the atmospheric correction process.

2. The radiation transfer simulation method according to claim 1, characterized in that: The step S100 includes: Step S101, taking the preset observation point as the origin, the normal direction perpendicular to the surface where the observation point is located as the z-axis, the projection direction of the sun on the xoy plane as the positive direction of the x-axis, and the y-axis determined by the right-hand screw rule, to construct a spatial coordinate system for the movement of photons; Step S102, obtaining geometric conditions for sun and satellite observation and atmospheric height, wherein the geometric conditions at least include an observed solar zenith angle, an observed zenith angle, and a relative azimuth angle; Step S103, determining the initial movement direction and position of the photon according to the geometric conditions, comprising: The satellite's sensor is located at the top of the atmosphere, and the initial direction of motion of the photon m = (a, b, c) and position p = (x, y, z) are calculated according to the following formula: Where θ is the observation zenith angle, is the relative azimuth, height is the altitude in the atmosphere; Step S104: Set the initial weight of the photon to 1, and set the polarization state matrix to the unit matrix.

3. The radiation transmission simulation method according to claim 2, characterized in that: The step S200 includes: S201, randomly sampling to obtain the free path of photons; S202, after the photon runs according to the free path, calculating the position of the photon, including: Where l is the free path of the photon; S203, determining whether the photon is located in the atmosphere according to the position of the photon, if so, determining that the current position state of the photon is inside the atmosphere; if not, executing step S204; S204. Determine whether the photon reaches the water surface or the land surface after running the free path according to the land / water position relationship, and determine that the current position state of the photon is the water surface when the photon reaches the water surface, and determine that the current position state of the photon is the land surface when the photon reaches the land surface.

4. The radiation transmission simulation method according to claim 1, characterized in that: The step S700 includes: When the position state of the photon is inside the atmosphere, it is determined that the photon is scattered, and a random sampling method and a rejection-acceptance method are used to obtain a scattering angle, a random sampling method is used to obtain a scattering azimuth, and the next movement direction of the photon is calculated based on the scattering angle and the scattering azimuth; When the position state of the photon is the land surface, it is determined that the photon has Lambertian reflection, and the reflection angle and azimuth of the photon are obtained by a random sampling method, and the next movement direction of the photon is calculated based on the reflection angle and azimuth; When the position state of the photon is the surface of the water body, it is determined that the photon undergoes mirror reflection, the reflection angle of the photon is obtained according to the incident angle of the mirror reflection, the azimuth of the photon is obtained by a random sampling method, and the next movement direction of the photon is calculated based on the reflection angle and the azimuth.

5. A system for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface, applied to any one of the methods for simulating the transmission of Rayleigh scattering vector radiation on a non-uniform underlying surface as claimed in claims 1 to 4, characterized in that: The radiation transfer simulation system comprises: The photon simulation module is configured to simulate the motion state of a satellite emitting photons to a preset observation point in the atmosphere-surface system according to preset conditions, wherein the preset conditions at least include the position relationship of the land / water body and the atmospheric height at the observation point, as well as the initial weight, position, motion direction and polarization state of the photons; A position determination module is configured to randomly sample the free path of the photon, and determine the position state of the photon after the photon runs according to the current movement direction and free path of the photon, wherein the position state includes the interior of the atmosphere, the surface of the land, and the surface of the water body; An updating module is configured to update the polarization state and weight of the photon using a preset updating method corresponding to the position state; A contribution value calculation module is configured to calculate the contribution value of the photon to the top-of-atmosphere radiance using a corresponding preset calculation method according to the position state of the photon; The judgment and radiance determination module is configured to determine whether the current weight of the photon is lower than a preset threshold, When the current weight of the photon is lower than the preset threshold, it is determined that the photon has died, and it is judged whether the number of photons currently simulated has reached the preset total number; when the number of photons currently simulated has not reached the preset total number, the photon simulation module is used to simulate the motion state of the next photon emitted by the satellite; when the number of photons currently simulated has reached the preset total number, the average value of all photon contribution values ​​is calculated as the top-of-atmosphere radiance received by the satellite under the preset conditions, so as to be used for the correction of Rayleigh scattering radiance in the atmospheric correction process; When the current weight of the photon is not lower than a preset threshold, the movement direction of the photon is updated, and the position determination module is used to simulate the next movement trajectory of the photon until the photon disappears.

Images