CN112348872B

CN112348872B - Modeling method of infrared accumulated cloud based on physical process

Info

Publication number: CN112348872B
Application number: CN202011044260.0A
Authority: CN
Inventors: 黄曦; 雷越; 孙浩; 马向超; 吴鑫; 刘德连; 张建奇; 单仝
Original assignee: Xidian University
Current assignee: Xidian University
Priority date: 2020-09-28
Filing date: 2020-09-28
Publication date: 2023-10-03
Anticipated expiration: 2040-09-28
Also published as: CN112348872A

Abstract

According to the modeling method of the infrared accumulation cloud based on the physical process, which is provided by the embodiment of the invention, the processes of water vapor transmission, phase change characteristics and the like existing in the grids are calculated through a three-dimensional grid model of the accumulation cloud, so as to obtain the distribution of the characteristics of cloud liquid drop density, water vapor density, temperature and the like, then, a multiple scattering simplified model is adopted, the illumination phenomenon existing in a visible light model of the accumulation cloud is classified and calculated according to the illumination direction and the sight direction, after a scattering phase function and a basic scattering theory are determined, the radiation energy of multiple scattering is calculated, the spontaneous radiation energy of the accumulation cloud in an infrared scene is calculated, so that the total energy is obtained, the radiation energy of the visible light model is corrected, the calculation process of spontaneous radiation under different wave bands is designed, so as to obtain the accumulation cloud infrared illumination model, and then, the illumination of the cloud is rendered in real time. Therefore, the cloud accumulation of the infrared band and the non-infrared band is considered, the reality and the instantaneity of the cloud accumulation rendering are higher, and the accuracy of the cloud accumulation rendering is improved.

Description

Modeling method of infrared accumulated cloud based on physical process

Technical Field

The invention belongs to the field of infrared waves, and particularly relates to a modeling method of infrared accumulation clouds based on a physical process.

Background

The infrared scene simulation is an effective means for analyzing and researching the infrared scene, and relies on accurate three-dimensional modeling and strict physical process simulation, so that the infrared virtual scene has strong sense of reality, and can meet many military requirements and civil requirements. The infrared air scene simulation needs to analyze targets and scenes, especially some common air scenes including clouds, sunlight, rain, snow and the like. Clouds are typical air scenes and often occur in flight simulation and outdoor scene simulation, but because their appearance of optical characteristics, shapes, etc. is related to altitude, climate, weather, sunlight, viewing angle, etc., the simulation of clouds is quite complex and difficult. In the prior art, a modeling method of stacking spherical particles in a three-dimensional space is adopted to obtain a cloud accumulation shape with layering sense, and the result is compared with the result of an improved noise modeling method, so that the advantages and disadvantages of the two methods are analyzed on the basis of ensuring the rendering efficiency. Because the three-dimensional model in the existing cloud simulation system mostly uses a traditional modeling method, the physical process in the cloud is lack of accurate description, and in the aspect of illumination models, a simpler forward scattering model is adopted, and various illumination phenomena in the cloud, such as secondary scattering and multiple scattering effects, are adopted. On the whole, the existing cloud simulation system is difficult to achieve both reality and real-time performance after comprehensive three-dimensional modeling and illumination modeling, is mostly applicable to the field of visible light, and lacks an illumination model for infrared band adaptation.

Disclosure of Invention

In order to solve the problems in the prior art, the invention provides a modeling method of infrared accumulation cloud based on a physical process. The technical problems to be solved by the invention are realized by the following technical scheme:

the embodiment of the invention provides a modeling method of infrared accumulation cloud based on a physical process, which comprises the following steps:

s1, acquiring a three-dimensional grid model of the accumulated cloud;

wherein each grid in the three-dimensional grid model stores first variable information of a plurality of particles, the particles representing clouds at respective locations within the grid, the variable information comprising: drop density, temperature, and velocity;

s2, determining cloud accumulation after diffusion according to the density, temperature and speed of the liquid drops based on the density, temperature and speed of the liquid drops;

s3, determining a scattering phase function based on an included angle between an incident direction and a scattering direction of illumination irradiation of the accumulated clouds;

s4, determining an effective phase function of the size distribution of the liquid drops based on the width of the size distribution of the liquid drops in the cloud and the characteristic radius of the size distribution of the liquid drops;

s5, aiming at the cloud accumulation in different wave bands, calculating an average scattering efficiency factor according to the scattering efficiency factors in different wave bands;

S6, averaging scattering phase functions of the clouds in different wave bands, and determining an average scattering function;

s7, calculating first radiation energy of light scattered for multiple times in the cloud according to the sight direction and the light direction of the colorant after the light is incident into the cloud based on the average scattering function and the average scattering efficiency factor;

s8, aiming at the cloud accumulation in the far infrared band and the middle infrared band, respectively calculating the spontaneous radiation energy corresponding to the cloud accumulation;

s9, summing the spontaneous radiation energy and the first radiation energy, and determining the total energy of the cloud in the sight line direction;

s10, multiplying the total energy by the maximum pixel value of each channel in the three primary color channels of the image respectively to obtain the pixel value of the accumulated cloud in the sight line direction;

and S11, rendering the diffused clouds at the positions in the direction of the line of sight according to the pixel values, and repeatedly executing S8 to S11 until each position of the diffused clouds is rendered.

According to the embodiment of the invention, through a three-dimensional grid model of the cloud, processes of water vapor transmission, phase change characteristics and the like existing in the grid are calculated to obtain the distribution of the characteristics of the cloud such as droplet density, water vapor density, temperature and the like, then a multiple scattering simplified model is used for classifying and calculating illumination phenomena existing in a visible light model of the cloud according to illumination directions and sight directions, after a scattering phase function and a basic scattering theory are determined, radiation energy of multiple scattering is calculated, spontaneous radiation energy of the cloud in an infrared scene is calculated, so that total energy is obtained, the radiation energy of the visible light model is corrected, the calculation process of spontaneous radiation under different wave bands is designed to obtain a cloud infrared illumination model, and then illumination of the cloud is rendered in real time. Therefore, the cloud accumulation of the infrared band and the non-infrared band is considered, the reality and the instantaneity of the cloud accumulation rendering are higher, and the accuracy of the cloud accumulation rendering is improved.

The present invention will be described in further detail with reference to the accompanying drawings and examples.

Drawings

FIG. 1 is a schematic flow chart of a modeling method of infrared accumulation clouds based on a physical process according to an embodiment of the present invention;

FIG. 2 is a schematic diagram of the present invention for realizing integration of the direction of the line of sight and the direction of the light, and circularly processing spontaneous scattering calculation and scattering calculation of the infrared band, thereby obtaining radiant energy;

FIG. 3 is a detailed flow chart of the invention for realizing integration of the line of sight direction and the light direction, and circularly processing spontaneous scattering calculation and scattering calculation of the infrared band, thereby obtaining radiation energy;

FIG. 4 is a schematic diagram of a three-dimensional mesh model provided by an embodiment of the present invention;

FIG. 5 is a flow chart of a semi-Lagrangian method for solving transport effects of velocity, temperature, and density due to a fluid velocity field provided by an embodiment of the present invention;

FIG. 6 is a flow chart of cloud computing diffusion provided by an embodiment of the present invention;

FIG. 7 is a schematic diagram of an initial cloud computing provided by an embodiment of the present invention;

fig. 8 is an illumination schematic diagram of the cloud computing provided by the embodiment of the invention;

FIG. 9 is a schematic diagram of calculating total energy from sampling points set in the direction of the line of sight and the direction of the light provided by the embodiment of the invention;

Fig. 10 is an infrared illumination simulation result of cloud accumulation in an infrared scene according to an embodiment of the present invention.

Detailed Description

The present invention will be described in further detail with reference to specific examples, but embodiments of the present invention are not limited thereto.

Example 1

As shown in fig. 1, the modeling method of infrared accumulation cloud based on physical process provided by the invention comprises the following steps:

s1, acquiring a three-dimensional grid model of the accumulated cloud;

the total energy needs to be quantized, so that the quantized total energy value is located in a range of 0-1.

Referring to fig. 2, after the three-dimensional model is obtained in fig. 2, integration of the line of sight direction and the light direction is needed, spontaneous scattering calculation and scattering calculation of the infrared band are circularly processed, so as to obtain radiant energy, and the radiant energy is rendered, and formulas and principles of each step are described in more detail later.

Example two

As an alternative embodiment of the present invention, the step of S2 includes:

s21, setting initial cloud accumulation;

s22, initializing the droplet density, temperature and speed of particles in each grid in the initial cloud;

S23, moving each particle in the grid according to the initialized speed;

s24, calculating the external force of each moved particle based on the wind power, buoyancy and vortex suppression force of the particle;

s25, calculating the droplet density and the temperature of the moved particles based on the current atmospheric pressure of the accumulated cloud; 2.26

S26, simplifying the acquired fluid motion equation based on the external force to obtain a simplified fluid motion equation;

s27, solving the simplified fluid motion equation to obtain the speed of the particles in the current unit time;

and S28, replacing the initialized speed with the speed in the current unit time, and repeatedly executing the steps from S23 to S27 until the iteration times are reached, so as to obtain the diffused cloud.

Example III

As an alternative embodiment of the present invention, step S26 includes:

acquiring a Navie-Stokes equation as a fluid motion equation;

setting the minimum unit of the intermediate velocity field and the time division, and simplifying the fluid motion equation.

Example IV

As an optional embodiment of the present invention, the step of S4 includes:

acquiring a scattering phase function;

determining different droplet size distributions based on the total particle number in the droplets, the width of the droplet size distribution, the gamma distribution function, the characteristic radius of the size distribution, and the effective radius corresponding to the size distribution;

And weighting and summing scattering phase functions of different droplet sizes based on the different droplet size distributions to obtain an effective phase function of the droplet size distribution.

Example five

As an alternative embodiment of the present invention, the step of S8 includes:

aiming at the cloud in the far infrared band, calculating the spontaneous radiation energy corresponding to the cloud by using a Planck formula;

for the cloud of the mid-infrared band, the spontaneous emission energy is calculated based on the set emissivity of the radiation.

Example six

As an optional embodiment of the present invention, before the step S8, the modeling method further includes:

the emissivity of light in the far infrared band and the emissivity of radiation in the mid infrared band are determined according to the law of light absorptivity and energy conservation.

Example seven

As an alternative embodiment of the present invention, the step of S7 includes:

s71, setting light parameters; the light parameters include: the light projection direction, the illumination direction, the light projection step length and the maximum number of light projections;

s72, determining a sight line direction sampling point along the sight line direction according to the light projection step length;

s73, determining second variable information at the current sampling point of the sight line direction;

S74, calculating the transmittance in the sight line direction based on the second variable information;

s75, determining second variable information at the current sampling point of the light direction according to the light projection step length along the light direction of the current sampling point;

s76, calculating the transmittance of the current sampling point in the light direction based on the second variable information;

s77, accumulating the transmittance;

s78, repeating the steps from S75 to S76 until the light transmittance of all sampling points in the light direction is calculated;

s79, calculating the primary scattering rate and the secondary scattering rate of the light subjected to scattering along the light ray direction;

s80, calculating second radiation energy based on the primary scattering rate and the secondary scattering rate;

s81, accumulating the second radiation energy respectively;

s82, taking the next sampling point along the sight line direction as the current sampling point, and repeatedly executing the steps from S73 to S82 until the sampling point along the sight line direction is polled;

s83, calculating the radiation proportion after multiple scattering along the light ray direction;

s84, calculating the first radiation energy based on the radiation proportion and the accumulated second radiation energy.

Referring to fig. 3, the entire flow of polling the gaze direction sampling point and the light direction sampling point is shown in fig. 3, and is a flow of embodiment seven, and the specific process is the same as the implementation process of embodiment seven.

Example eight

As an alternative embodiment of the present invention, the step of S78 includes:

determining an extinction cross-section based on an effective radius of the droplet and a population density of the droplet;

and determining the light transmittance based on the extinction interface.

Example nine

As an alternative embodiment of the present invention, the step of S79 includes:

determining a projection function of the light passing through the cloud based on the transmittance, the sight direction and the light direction when the light passes through the cloud layer;

determining a primary scattering function of light passing through the cloud;

determining a secondary scattering function of light passing through the cloud;

calculating a primary scattering rate using the primary scattering effect function and calculating a secondary scattering rate using the secondary scattering effect function.

Examples ten

As an alternative embodiment of the present invention, the step of S83 includes:

when the cloud is in a thermal equilibrium state and is in multiple scattering, determining an energy conservation equation of the transmittance, the primary scattering proportion, the secondary scattering proportion and the multiple scattering proportion of the light;

determining the radiation proportion of the individual bins based on the energy conservation equation when the multiple scattering is uniform;

when the multiple scattering is non-directional, the radiation proportion of different directions on a single surface source is determined based on the energy conservation equation.

The following is a detailed description of the basic principles in the embodiments of the present invention in the steps of the first to tenth embodiments.

The meaning of the parameters is described first.

_u Indicating the speed, t indicating the time, ρ indicating the density, p indicating the pressure, f indicating the force, u ^* Represents an intermediate velocity field, θ represents potential temperature, q represents particle density, T represents Kelvin temperature, pi represents Exner function, p ₀ Represents standard atmospheric pressure, R _d Representing the gas constant of dry air, c _p Represents the constant pressure specific heat capacity of dry air, c _v Represents the volumetric specific heat capacity of dry air, f _wind Indicating the wind force, f _buo Representing buoyancy effect, f _vc Represents the vortex suppression force, theta _v Representing the temperature of the virtual potential, θ _v0 Represents the reference virtual temperature, q _c Represents the drop density, q _v Represents the water vapor density, ω represents the rotation of the fluid velocity field, l represents the mesh unit length, q _vs Represents saturated vapor density, z represents altitude, Γ represents vertical air temperature reduction rate, L represents latent heat of vaporization of water, θ _b The intermediate potential temperature, p (θ) the scattering phase function, k the wave number, csca the scattering cross section, X the vector scattering amplitude, N (r) the droplet size distribution, γ the width of the size distribution, Γ (γ) the gamma distribution function, N ₀ Represents the total particle number, r, in the water drop _n Characteristic radius r representing size distribution _e An effective radius, P (θ) an effective phase function, τ (d) a transmittance,indicating the direction of vision>The light direction, Q the scattering efficiency factor, sigma the extinction cross section, kappa the extinction coefficient, S ₀ Indicating transmission effect, S ₁ Representing the primary scattering effect, S ₂ Representing the effect of secondary scattering, < >>The total transmission ratio is indicated by the ratio,represents the total primary scattering ratio, +.>Representing the total secondary scattering ratio,/>Representing the total multiple scattering ratio, ">Representing the average effective phase function, M _λbb Representing spectral radiant emittance, c ₁ Representing a first radiation constant, c ₂ Representing the second radiation constant, α represents the absorption, epsilon represents the emissivity, and G represents the quantization gray.

In the first embodiment, the three-dimensional model contains three-dimensional information, and contains information of three directions of XYZ axes. In order to use and access this information conveniently and quickly, a suitable model needs to be designed. According to the embodiment of the invention, the information is stored by adopting a three-dimensional grid, the information such as density, temperature and speed in the cloud is stored, and the result is stored in a three-dimensional map, so that the use and access in the GPU computing process are facilitated.

As shown in fig. 4, the three-dimensional model in this chapter is built under a left-hand coordinate system, with XYZ axes forward as shown. The cube in the figure is a three-dimensional grid for storing cloud model information, the grid is divided into a plurality of uniform grids, 4 x 4 grids in the figure, and the grid fineness is far greater in the actual calculation process. Each uniform cell after being divided is called a voxel (voxel), and each voxel stores cloud information of the current position. Cloud generation in a grid, the volume of the cloud is typically smaller than the volume of the grid in order to control the shape of the cloud and observe changes in the cloud. When the meshing is fine enough, voxels within the cloud volume can represent particles in the cloud, storing the droplet information and the vapor information of the current position.

The coordinate system in fig. 4 is also a coordinate system of the simulated virtual scene according to the embodiment of the present invention, and the observation camera is generally located on the axis in the positive direction of the X-axis, and the observation direction (-1, 0). Some azimuthal descriptions of embodiments of the invention are relative to a cloud grid, which is at the origin of the coordinate axes. Wherein "front" refers to an axial positive direction with respect to the cloud grid, "up" refers to a Y-axial positive direction with respect to the cloud grid, "right" refers to a Z-axial positive direction with respect to the cloud grid, and "rear", "lower" and "left" respectively correspond thereto.

The implementation of the two steps of the embodiment and the parameters required in the process are described in detail below.

(1) Incompressible conditions

The actual fluid has different degrees of compressibility, and in aerodynamics, whether the density change of the gas is negligible or not is determined according to the speed of the gas flow. Generally, low velocity aerodynamic methods are applicable when the gas flow rate is below Mach 0.3, where the gas may be considered incompressible; when the gas flow rate is higher than mach 0.3, high-speed aerodynamic methods are applicable, where the gas should be considered to be compressible.

(2) The viscous forces are negligible

Reynolds number (Reynolds number) is a dimensionless number that can be used to characterize fluid flow, and when the Reynolds number is less than or equal to 1, i.e., the fluid viscosity is small and the relative sliding speed is not large, the viscous stress is small outside the boundary layer of the flowing object, and the viscous term in the equation can be ignored, so that the fluid can be regarded as ideal fluid. With water, air, etc. under natural conditions, the viscous force can be ignored because the movement speed is small and the density is low.

The fluid motion equation and the simplified process in the second embodiment of the present invention are described in detail.

(1) Equation of fluid motion

The Navier-Stokes equation (Navier-Stokes equation, N-S equation for short) is a motion equation describing conservation of momentum of viscous incompressible fluids. Because the expansion speed of the fluid is far lower than the sonic speed, the viscous force term of the fluid can be ignored, and thus the N-S equation can be simplified into an Euler equation, and the form is as follows:

at the same time, the continuous conditions are as follows:

where u is the speed, ρ is the density, t is the time, p is the pressure, f is the effect of the force, equation (2-1) is the conservation of fluid momentum, the right term is the motion update of the fluid velocity field itself, the second term is the acceleration caused by the fluid pressure gradient, and the last term is the acceleration caused by the action of external force, typically including gravity, wind force, etc. The formula (2-2) is a continuity equation, which shows that the fluid velocity field has no divergence, and the mass conservation is ensured.

It is not practical to solve the euler equation completely under real-time conditions, so some approximation is required. The embodiment of the invention divides the formula (2-1) into two parts, and the first step is to calculate an intermediate velocity field u ^* And assuming that the minimum unit of time division is Δt, there are:

for equation (2-3), the first term on the right of the equal sign is the transport result of the fluid velocity along with the fluid motion, and can be approximated by using the semi-lagrangian method, and the same method can be used to update other constants, such as temperature and water density (here also referred to as water vapor and liquid droplets contained in the voxel):

Where θ is the potential temperature of the fluid and q represents the density of the particles (gaseous or liquid). The potential temperature is constant at an adiabatic approximation of the change in height, which can be calculated:

in the formula (2-6), T is temperature, the unit is K, pi is Exner function, and the temperature can be regarded as dimensionless pressure of gas. P is p ₀ Is the standard atmospheric pressure at sea level, and p is the pressure at the height. R is R _d Is the gas constant of dry air (287 Jkg) ^-1 K ^-1 )，c _p And c _v The constant pressure specific heat capacity and the volume specific heat capacity of dry air are respectively.

(2) Water vapor transport effect

The process of solving the transport effect of speed, temperature and density due to the fluid speed field by the semi-Lagrangian method mainly comprises five steps, as shown in FIG. 5:

taking the temperature as an example, firstly obtaining the temperature and the speed v of the position pos to be updated, and then rewinding the position pos to half step length to calculate the intermediate position mdpos=pos-v 0.5 x dt, and obtaining the speed v of the intermediate position mdpos _mid One step back to calculate the final position finalpos=mdpos-v _mid * dt, finally updating the pos position with the temperature of the finalpos positionTemperature. In this way, a temperature-dependent fluid velocity field transport is accomplished.

(3) Force application and correction

After the speed is conveyed along with the speed field, the acceleration f caused by the external force is considered, and the external force generally comprises the following three items, namely, wind force f _wind Buoyancy f _buo Vortex suppression force f _vc Thus, there are formulae (2-8).

f＝f _vc +f _buo +f _wind (2-8)

For wind force f _wind The simulation is typically independent with respect to the cloud, so that the setting or taking of typical wind speeds under simulated conditions can be done directly. Buoyancy f _buo The calculation of (1) comprises gravity and is actually the result of stress analysis on the vertical direction of the air mass in a voxel. The higher the water vapor ratio in the air mass is, the larger the buoyancy is; the higher the drop ratio, the greater the gravity and the less the buoyancy. In addition, the calculation of the buoyancy force is also related to the temperature and the pressure of the environment. The specific calculation is as follows:

wherein g is the gravity acceleration rate,is a unit vector vertically upwards, q _c Is the density of liquid drop, theta _v0 Is a reference virtual temperature and typically takes a value between 290K and 300K. θ _v The virtual potential temperature represents the potential temperature of dry air under the same density and pressure conditions, and the calculation method is as follows:

θ _v ＝θ(1+0.61q _v ) (2-10)

wherein q is _v Is the density of water vapor.

Equation (2-10) is a continuous equation, and only discrete grid calculation can be adopted in the actual calculation process, namely, a discrete difference equation is used to approach the solution of the continuous equation. Deviation of the discrete difference equation from the continuous equation in the result toAnd coarser grid division of the calculation force requirement determines that details of the calculation result on a small scale are difficult to embody, such as the rotation effect of fluid. To solve this problem, a vortex suppression force f is introduced _vc These fluid details are presented using a series of calculations based on the velocity field. The vortex suppression force can be calculated by first calculating the rotation ω of the fluid velocity field as follows:

recalculating the rotation field _ω Is unitized by the gradient direction of (a):

finally, the vortex suppression force is calculated:

f _vc ＝ηl(N×ω) (2-13)

in the formula (2-13), l is the unit length of the grid, eta is the vortex suppression coefficient, the strength of the vortex suppression force is determined, namely the small-scale action compensation effect in the process of the cloud accumulation motion is determined, and the value is generally between 0.1 and 10, and is determined according to the specific situation.

(4) Pressure field solution

After correction of the calculated force we get an intermediate velocity field and this velocity field satisfies the equation after ignoring the effect of the density field gradient:

the right side of the above equation is known, so that the poisson equation is a poisson equation, the pressure p can be solved by using an iteration method, the speed is corrected by using the pressure, and a final speed field can be obtained:

thus, the motion characteristics of the fluid at the unit time are calculated.

(5) Phase change characteristics of particles

The phase change of the cloud refers to the state change of water occurring in the cloud, i.e. the conversion process between droplets, water vapor, ice crystals. In the embodiment of the invention, the phase change is assumed to be generated by taking voxels in the cloud as units, namely, the state of water molecules in each voxel is changed. Typically, the vapor density exceeds the saturation point at the current conditions, and the vapor is converted to droplets, which in turn form a cloud that is seen later. The saturation point of the water vapor density in the air is related to the temperature T and the atmospheric pressure p of the current position, and can be calculated by the following method:

The pressure in equations (2-16) requires a pressure correction taking into account the current position:

in the formula (2-17), z is the altitude of the current position, p ₀ For standard atmospheric pressure, 100kPa, T may be taken ₀ For the surface temperature, 280k-310k may be taken, Γ 10k/km, indicating a decrease in air temperature of 10k per 1km increase.

For clouds, conservation of water during phase change must be considered, and in clouds, the following continuous conditions are satisfied:

q in the formula (2-18) _v For water vapor density, q _c Is the drop density. C is a constant, i.e. the total amount of water molecules undergoing phase change is constant.

When the water vapor density is larger than the saturation point, phase change occurs, and the water vapor density and the liquid drop density change:

Δq _v ＝-C＝min(q _vs -q _vb ,q _cb ) (2-19)

q _v ＝q _vb +Δq _v (2-20)

q _c ＝q _cb -Δq _v (2-21)

q in the formula (2-19) _vb 、q _cb Refers to the moisture density and drop density within the voxel in the previous state.

Although the adiabatic approximation is a valid approximation of the movement of air with unsaturated water vapor, it cannot be assumed that the potential temperature of saturated air is constant. If the density of the water vapor in the voxel exceeds the saturation point, the water vapor is condensed and releases latent heat, so that the voxel air mass is heated. If the latent heat and cooling due to condensation and evaporation are the only non-adiabatic changes, then the first law of thermodynamics holds that:

l in the formula (2-22) is latent heat of vaporization of water, and 2.501J/kg is taken at 0 ℃. The above equation illustrates that the change in potential temperature is related to the phase change at the current location in addition to the change in transport with the velocity field. So after considering the transport variation, there is the following formula:

θ in the formula (2-23) _b Representing the intermediate potential temperature after the potential temperature is delivered with the velocity field. Thus, the change characteristics of the particle phase change are calculated.

The density field of the cloud which is rendered in a unit time (which is expressed as a frame in the simulation process) is calculated, and besides the motion, the phase change, the heat energy conversion and the like of each voxel are calculated, then the partial differential equation is also required to be calculated in an iterative mode to calculate the pressure, and the pressure distribution is used for correcting the speed field. Such complex calculations require a considerable number of calculations, but the real-time requirements of the simulation place limits on the calculation time per frame. In order to realize real-time simulation (the frame frequency reaches at least 30), the embodiment of the invention transfers all simulation calculation to the GPU for carrying out, and optimizes repeated calculation by utilizing the high parallelism of the GPU, thereby improving the frame frequency to meet the real-time simulation requirement while ensuring high authenticity. The calculation process of each voxel in the cloud is shown in fig. 3 in a frame time, fig. 6 shows the process of implementing S2 in the second embodiment, in which, in the initial stage, the initial cloud is shown in fig. 7, the white area represents the initial cloud, the box represents a grid, and the particles of each grid need to be diffused to obtain the diffused cloud.

The following details the process of solar ray irradiation cloud emission according to embodiments one to ten.

As shown in fig. 8, the illumination schematic diagram of the cloud is that the sunlight is different from the sight line direction, and the sunlight enters the view point through scattering. For one point P in the cloud, consider the ray of the point P into the viewpoint. Firstly, light enters a P point through part of cloud layers, and after being scattered by the P point, primary scattered light is emitted along the direction of sight; second, the illumination passes through part of the cloud to other points, e.g. P ₁ 、P ₂ The illumination enters the sight after secondary and tertiary scattering; finally, diffuse reflection of the cloud surface, self-radiation, and the like enter the line of sight. In addition, since the end point of the line of sight is different for different points, sky background radiation needs to be considered, for example, sky radiation transmitted in the line of sight direction needs to be considered in fig. 8. When the radiation simulation is in the sky, the sky radiation may be processed with a background radiation metering algorithm. When the viewing direction is the same as the illumination direction, it is necessary to consider solar radiation transmitted in the viewing direction.

(a) Scattering phenomena in clouding

(1) Scattering calculation method

The scattering effect in the cloud is calculated, first requiring the calculation of a scattering phase function. The scattering phase function describes the angular distribution of energy after the incident energy is scattered by the particles. Let the angle θ between the incident direction and the scattering direction be the scattering angle, introduce a dimensionless parameter to describe the variation of the scattering energy with the scattering angle, called the phase function p (θ). For the phase function, there are normalization conditions:

The Mie scattering phase function (hereinafter referred to as Mie function) can be obtained by solving maxwell's equations, but the process is too complex and will not be described in detail herein. The scattering phase function p of the final spherical particle at a single wavelength can be expressed as:

p＝2π(|S ₁ | ² +|S ₂ | ² )/(k ² Csca) (3-2)

where k is the wave number, csca= ≡ _4π |X|/k ² dΩ is a scattering cross section, X is a vector scattering amplitude, and S1 and S2 are complex amplitudes, which can be obtained by solving the bessel equation, and are not described in detail in the embodiments of the present invention.

The cloud in the low altitude is formed by mixing water drops with different sizes, the water drops are assumed to be spherical, and the size of the water drops has great influence on the optical characteristics of the cloud at the moment: for a certain amount of water, the smaller the drop, the more opaque the cloud. Approximately half the size is reduced and the extinction efficiency is doubled, the Mie scattering phase function also varying with droplet size. To analyze this, it is first necessary to model the droplet size. The Droplet Size Distribution (DSD) model in the cloud approximates the modified gamma distribution. The definition is as follows:

wherein N is ₀ ＝300cm ^-3 Gamma is the width of the size distribution, typically taking gamma=2, for the total number of particles in the water droplet. Γ (γ) is a gamma distribution function, r _n Is the characteristic radius of the size distribution. r is (r) _e ＝(γ+2)r _n For the effective radius corresponding to this distribution, the embodiment of the present invention takes r _e =7μm. After the parameters are determined, the droplet size distribution can be obtained.

The effective phase function P (θ) corresponding to the droplet size distribution can be obtained by summing the weighted sums of the phase functions of different droplet sizes, the summation formula being as follows:

/>

and finally obtaining the corresponding effective phase function under the specified wavelength.

Because the whole calculation flow of the effective phase function is quite complex, the effective phase function is difficult to embed into real-time calculation, a pre-calculation method is mostly adopted, and the effective phase function is stored in tools such as textures and the like for direct value taking in real-time calculation, or calculation can be simplified through high-order function simulation, so that the real-time performance of the simulation process is ensured.

After determining the droplet size model and the phase function, a specific scattering process can be considered. Assuming that the line of sight direction is vertically downward, the sunlight direction is vertically downward transmitted, and the transmission refers to the process that light directly penetrates through the cloud layer without being scattered by particles and finally enters the view point. When spontaneous radiation of the cloud layer is ignored, only the attenuation effect of the cloud layer on light rays is considered in the process, and the scattering effect of particles is not considered, so that the method can be finally expressed as:

where τ (d) is the transmittance of light passing through the cloud layer when the direction of the line of sightIs +. >Transmission need not be considered at the same time.

Assuming that the individual voxels within the cloud are uniformly distributed, the permeability of the voxels may then be determined by the effective radius r of the drop _e And the particle number density N of the droplets. Using these two parameters, the extinction cross-section can be defined:

where Q is the extinction efficiency factor and is a function of wavelength. The extinction coefficient is obtained from the extinction cross section:

κ＝Nσ (3-7)

the transmittance of the cloud can be obtained by the method:

where x is the specified ray path and the visible transmittance τ is the integral over that path. Considering all voxels through which the path passes, a voxel-wise additive transmission result is obtained, where l represents the corresponding length of the voxel grid. From the frequency point of view, this is the attenuation result of a beam of light after it has passed a path; from a probability point of view, this is the probability that a ray will pass through a path without encountering a particle, i.e. the probability that no scattering will occur through a path. The probability of scattering of the hit particles in this path x is thus 1- τ, for which the derivative function:

s(x)＝κe ^-∫κdx (3-9)

it is understood that the probability density of hitting a particle at a certain position on path x, so s (x) is the probability density of hitting just at the end of a free path of length x.

(2) Primary scattering process

Primary scattering refers to the process by which light is scattered by only one particle and then transmitted through the cloud layer to finally enter the viewpoint. Neglecting spontaneous emission of the cloud layer, the attenuation effect of the cloud layer on light rays and the primary scattering effect of particles need to be considered in the process, and finally the method can be expressed as:

S ₁ ＝s(d _l )P(θ _vl )τ(d _v ) (3-10)

the first term on the right of the equal sign is the probability that a ray hits a given particle (e.g., particle P) directly through the cloud; the second term is a phase function representing the probability of scattering of the particle in a given direction (e.g., line of sight direction v) at the current angle; the third term is the probability that the ray leaves the particle and directly enters the viewpoint without collision.

(3) Secondary scattering process

The secondary scattering refers to the process that light rays are scattered by two particles and finally enter the view point through the cloud layer. Neglecting spontaneous emission of the cloud layer, the attenuation effect of the cloud layer on light rays and the primary and secondary scattering effects of particles need to be considered in the process, and finally the method can be expressed as:

for the secondary scattering process, one determined scattering path contains two scattering particles: and first and second scattering particles. The light rays hit the first scattering particles through the cloud layer, then hit the second scattering particles, and finally transmit out of the cloud layer and finally enter the view point. As with the primary scattering process, scattering of the exit particles, i.e. the secondary scattering particles, is considered. The scattered light particles scatter all of the primary scattered light hitting the particles, e.g., for point P, P1 and P2, and other particles act as primary scattered particles of P, providing primary scattered light thereto.

In fact, all particles in the cloud can be used as primary scattering particles, and the angle and brightness of only primary scattering rays are distinguished. However, in the real-time simulation process, due to the excessive calculation cost, the scattering of each emergent particle cannot be considered by all possible primary scattering particles, so that some approximation algorithms are needed for simplification.

As shown in the first right term of the formula (3-11), light rays with different angles and different distances are integrated, and because the integration is too complex, simplified calculation is needed, and the second and third formulas are approximate calculation. The first to third formulas are based on the assumption that the radiation weights of different paths under the same scattering angle are the same, and the second to third formulas replace the extinction coefficient of the primary scattering particles with the extinction coefficient of the secondary scattering particles to simplify calculation. The first simplified calculation is mainly based on the fact that the larger the path difference is, the lower the transmittance is and the smaller the weight is; the second assumption is that the calculation is based mainly on the fact that the particle number distribution is not very different when the distance between the liquid drops in the cloud is relatively close.

After the approximate integration, the result is independent of the position of the primary scattering point, and the final calculation result is only related to the position of the secondary scattering point, similar to the primary scattering calculation. Therefore, the scattering effect of the scattering point as a secondary scattering point can be calculated while the primary scattering is calculated, thereby avoiding multiple calculations and reducing the calculation pressure.

(4) Multiple scatter calculation

Multiple scattering refers to the process that light rays pass through a cloud layer and finally enter a view point after being scattered by three or more particles. Neglecting the spontaneous emission of the cloud, the light attenuation effect of the cloud and the scattering effect of the particles three times and more need to be considered in the process.

Multiple scattering in a cloud can be approximated by a similar and quadratic scattering method, but neither the efficiency of the calculation nor the accuracy of the calculation results is satisfactory. This is mainly due to the following reasons: 1. multiple scattering effects theoretically involve 3 to infinite scattering process calculations, and even for a single determined scattering exit point, scattering effects cannot be calculated one by one, and finally approximation algorithms still have to be used to ignore or simplify a portion of the calculations. Even so, the result obtained is affected by the approximation algorithm. 2. The approximate calculation of the secondary scattering always comprises the characteristics of directivity and symmetry, wherein the directivity refers to the approximate result of the secondary scattering and the approximate result of the primary scattering, and most of emergent energy is contained in a small angle of the incident direction; symmetry means that the approximate result of the secondary scattering is symmetrical, i.e. the scattering results are the same at both positive and negative scattering angles. For multiple scattering, the directivity of scattering is obviously reduced due to the multiple scattering of light, and as the scattering times are increased, the directivity of the final scattering completely disappears, and the scattering energy distribution at different angles is the same. Therefore, if the approximation algorithm of the secondary scattering is adopted to process the secondary scattering, the calculation result also contains strong directivity and is quite different from the actual result.

For the multiple scattering calculation result, the following assumption is used for the approximation process: 1. neglecting the influence of scattering times on scattering results, and uniformly considering scattering effects above three scattering levels. 2. Neglecting the effect of the incident light direction on the multiple scattering effect, it is assumed that the multiple scattering effect is uniform. 3. The directionality of the scattering effect is ignored and approximated by a diffuse reflection model.

When the cloud is in a thermal equilibrium state, based on assumption 1 and the law of conservation of energy, the following formula can be obtained:

in the middle ofRepresenting the transmittance of light, +.>Represents the total primary scattering, secondary scattering ratio, < ->Representing the multiple scattering ratio, there are

Thus, the total multiple scattering ratio is obtained, including

Based on assumption 2, the radiation ratio of the single bin A can be obtained, including

Based on assumption 3, the radiation proportion of different directions on a single surface source can be obtained, and the radiation proportion is that

Wherein θ is the angle between the observation direction and the normal vector of the patch at the exit position.

When the influence of the cloud accumulation shape is ignored, the simplest uniform sphere model is used for approximate calculation, and a rapid calculation method can be obtained. The only variables in this approach are the direction of the incident light, the size of the cloud, and the average particle count of the cloud. In a simulation scene, the direction of sunlight is generally fixed, and when the preset shape of cloud to be simulated is determined, the final stable shape is also fixed, so that the three variables can be preset and pre-calculated, the simulation efficiency is accelerated, and the real-time requirement is met.

(5) Sky radiation

Because the sky radiation effect is far smaller than the sunlight illumination effect, the scattering effect after the sky radiation enters the cloud is not considered, and the sky radiation is treated as background radiation. As shown in fig. 9, the patch where the intersection point of the observation direction and the sky is located emits radiant energy as a main body of sky radiation in the corresponding line-of-sight direction, and enters the viewpoint through penetration of the clouds, thereby finally affecting the observation effect. Therefore, the energy of sky radiation is related to the position, temperature, time period, attenuation of energy through an atmospheric path, attenuation of energy through clouding, and the like of the small patch, and the sky radiation is only related to the position of the small patch and the attenuation of energy through clouding under the condition that the temperature and the sun position are kept unchanged in the simulation time period under the condition that the atmospheric influence is ignored. The position of the small patch, especially the zenith angle of the view point corresponding to the small patch position, determines the initial energy of sky radiation, the energy corresponding to the position is stored in the texture after the completion of the pretreatment method, and the energy is directly sampled on the texture when the radiation energy is calculated; attenuation through the cloud can be represented using the transmittance of the cloud, which is related to the shape, density distribution of the cloud, and can be obtained in later calculations.

Further details of the self-radiation process in an infrared scene are required below.

(a) Infrared illumination model

There are many differences in simulation between the cloud accumulation in the infrared scene and in the visible scene, the most important of which is the processing of spontaneous radiation of the cloud, and secondly the scattering calculation problem caused by the different wave bands.

(1) Scattering calculation correction

The scattering calculation of the infrared band requires two corrections based on the calculation of the visible band: 1. scattering efficiency factor variation due to band difference; 2. the infrared band pass needs to take the band width into account and correct the scattering phase function accordingly.

For the first point, the average value of the efficiency factors in the band is substituted when the formula calculation is performed, and when the difference between the efficiency factors in the band is not large, the efficiency factor corresponding to the intermediate wavelength can be used as the average value. At a wavelength of 4 μm, the refractive index of the droplet was 1.351, and the effective radius of the particle was 7 μm, which gave a extinction efficiency factor of about 2.6.

For the second point, a scattering phase function for each wavelength within the corresponding band needs to be calculated. After the scattering phase function at a single wavelength is calculated according to the formula (3-4), an average phase function corresponding to the whole wave band is calculated.

Equation (4-1) gives accurate results of the effective phase function in the corresponding band, but is complex, and in fact, in a shorter band (e.g., mid-infrared band 3-5 μm), the change in the effective phase function is negligible, so that the effective phase function in the entire band can be calculated using the effective phase function corresponding to the mid-wavelength in the band as an average value.

(2) Spontaneous emission correction

In the visible light scene, the rendering of the virtual scene is based on the operation of RGB three colors, namely, the radiation brightness in three wave bands of visible light is a calculation object, but according to the Planckian radiation law, the spontaneous radiation energy of cloud in the three wave bands is extremely small, and compared with the illumination radiation of sunlight and the diffuse reflection effect caused by illumination, the virtual scene is completely negligible. In contrast, in the infrared band, especially in the mid-infrared band (3-5 μm) and the far-infrared band (8-12 μm), the spontaneous emission energy of the cloud is concentrated, and compared with the scattering effect caused by illumination, the spontaneous emission effect of the cloud in the infrared scene simulation needs to be considered.

From the aspects of reality and real-time of simulation, the traditional method for processing the spontaneous radiation effect of cloud in infrared scene simulation generally adopts the following two calculation models: black body approximation model and gray body approximation model.

The blackbody approximation model is used for solving the cloud infrared simulation problem in the far infrared band (8-12 mu m), and in the band, the solar radiation intensity is smaller, the illumination effect has smaller influence on the cloud simulation, and the infrared radiation of the cloud is mainly spontaneous radiation. Therefore, when only the sun is considered as the only light source and the sky is considered as the simulation background, the cloud can be regarded as a black body, and the spontaneous radiation energy at the moment can be directly calculated through the Planck formula (4-2).

Wherein M is _λbb For spectral radiant emittance, unit W/(m) ² μm), λ is the wavelength, in μm, T is the thermodynamic temperature, in K, h is the Planck constant, c is the speed of light, c ₁ Is the firstRadiation constant, c ₂ For a second radiation constant, c ₁ ＝2πhc ³ ≈3.7415×10 ⁸ W·μm ⁴ /m ² ，c ₂ ＝hc/K _B ≈1.43879×10 ⁴ μm·K。

The gray body approximation model is used for solving the cloud infrared simulation problem under the middle infrared band (3-5 mu m), solar radiation is dominant under the band, and the radiation effect of the cloud is mainly scattering solar radiation and spontaneous radiation is secondary. Therefore, when only the sun is considered as the only light source and the sky is considered as the simulation background, the cloud can be regarded as a gray body, and the emissivity of radiation is generally between 0.7 and 0.9.

Both models are approximate assumptions based on measured data, and although relatively simple, the simulation results are close to the real infrared image. And referring to the two simple approximate models, the self-emission correction is carried out on the illumination model proposed in the third chapter so as to ensure the accuracy of the illumination model. First, whether in the mid-infrared band or the far-infrared band, the absorption of infrared radiation by the cloud needs to be considered, and the absorption rate is assumed to be α, so according to the law of conservation of energy, the equation (3-12) should become:

For blackbody and gray body absorption, as known from kirchhoff's law,

α＝ε (4-4)

i.e. the absorptivity is equal to the emissivity, in particular the emissivity of a black body is equal to 1.

Based on the gray body approximation model, the absorption rate of the mid-infrared band is 0.7-0.9, the current multiple scattering energy proportion can be calculated by substituting the energy proportion into the energy proportion (4-3), and the multiple scattering proportion is less than 10% according to the calculation result. In addition, the multiple scattering of the embodiment of the invention adopts a treatment method similar to diffuse reflection, so that the multiple scattering effect can be completely ignored, and the energy is treated by a spontaneous radiation method, so that the method has the following advantages of

In the far infrared band, the extinction efficiency factor is smaller, the absorption is stronger, and the multiple scattering effect is completely negligible. Physically, the strong absorption in the cloud causes the scattering times of solar radiation in the cloud to be reduced, and the scattering intensity to be weakened, so that the detection is difficult. In addition, for primary scattering and secondary scattering effects in the far infrared band, since the phase function, scattering efficiency factor and solar radiation are smaller at the corresponding wavelengths, the relative spontaneous radiation is also negligible within the error allowable range, and thus the cloud exhibits blackbody radiation characteristics in the far infrared band.

Given the calculation method of spontaneous emission, it is necessary to determine the relevant simulation scheme, which is included in the previous simulation method. Similar to the transmission effect of light scattered on particles, spontaneous radiation starts from the particles and propagates along the observation path, and finally enters the viewpoint through absorption and scattering of the path, resulting in an observation result. As shown in fig. 9, P ₂ The particles generate spontaneous radiation, and the radiation effect is along P ₂ A enters the view point, also P ₁ The particles at P generate spontaneous radiation, and the radiation effect is along P ₁ A. PA enters the viewpoint, except for the transmittance difference caused by the path length difference. From this characteristic of the spontaneous emission, a simulation scheme can be determined, and the spontaneous emission characteristics are carried out in a rendering shader determined by a three-dimensional model, similarly to the illumination model, and when the integration is carried out by the observation path, the results of the spontaneous emission integration are accumulated, and finally the total spontaneous emission under the path is output.

After calculation of the transmitted, scattered, spontaneous radiant energy, the total energy in this path is obtained, as in fig. 2. Unlike visible light simulation, which directly uses RGB three colors for calculation, infrared simulation requires a certain quantization scale. And quantifying the result of radiation output within the determined scale, and outputting the final quantified result, so that the simulation experiment can be quantitatively analyzed. Assuming that the radiation brightness output of an object is L _output Upper limit of quantization L for simulation scene _max Lower limit L _min Obtaining the final quantization result GThere is

The quantization result G obtained by the formulas (4-6) is generally represented by gray scale, and the obtained radiation brightness output is L in the common 0-255 gray scale scene _output Corresponding gray g.255.

Before display, firstly determining quantization scale and temperature setting of simulation scene, wherein the scene quantization scale is radiation brightness scale, the simulation wave band is middle infrared wave band (3-5 μm), and the setting lower limit is black body radiation brightness corresponding to 270K, namely 0.6064W/(sr.m) ² ) The upper limit is 330K, namely 5.3518W/(sr.m) ² ). The ground scene temperature is 300K, the height of the bottom of the cloud is 1000m, so that the bottom temperature of the grid of the cloud is 290K, and a water vapor generation area is arranged in the center of the grid. Setting the sunlight illumination direction as (1, -1, 0), observing the sunlight illumination direction (-1, 0), and taking the solar surface temperature as 5500K to be used as black body approximation treatment.

Three illumination effects should exist in the infrared illumination simulation, namely primary scattering, secondary scattering and spontaneous radiation effects. Unlike visible light illumination simulations, since spontaneous radiation of the cloud in the infrared band is not negligible, the radiant energy concentrated in primary, secondary and multiple scattering is weakened by the presence of spontaneous radiation. The energy of the scattering decay is absorbed by the particles and output by spontaneous emission. Scene rendering is performed in combination with the simulation scheme mentioned in the first section, and the results are shown in terms of strong directional illumination (primary scattering, secondary scattering) and spontaneous emission.

Fig. 10 shows simulation results of infrared light, which are cloud-deposited in an infrared scene, from top to bottom, which are the combined effect of simulated light, strong directional light (primary scattering, secondary scattering) and spontaneous emission effect. As is apparent from the figure, compared with the dominant position of sunlight in the visible light simulation, in the infrared simulation scene, the proportion of strong directional illumination (primary scattering and secondary scattering) is obviously reduced, but the strong directional illumination is still the main body of illumination energy, and the illumination effect generated in the current simulation view angle has extremely strong directivity due to the illumination direction of sunlight; some of the energy in the simulation results comes from the spontaneous emission, and the spontaneous emission changes due to the change of the three-dimensional model of the cloud, but the energy of the spontaneous emission is scattered and absorbed as well, which results in a change of the spontaneous emission results with time, that is, a significant change with the divergence of the cloud.

Analyzing simulation results of stable states, wherein the highest gray 255 appears at the upper edge of the cloud, mainly because strong forward scattering occurs, and scattered light contains a large amount of solar energy in a small angle, so that the radiation brightness is large; average gray level 35 in spontaneous radiation effect corresponds to radiation brightness 1.3508W/(sr.m) ² ) The corresponding blackbody temperature 287K is slightly lower than the atmospheric temperature of the height of the cloud body, and meets the actual measurement result of the cloud storage emissivity.

As can be seen from fig. 10, the spontaneous emission effect is not uniform in the multiple scattering effect, but rather, is significantly changed due to the change of the three-dimensional model of the cloud. This is mainly because in the simulation scheme, the calculation object of the spontaneous radiation is cloud particles in the cloud, so the spontaneous radiation starts from the particles, and enters the view point after traveling a path, and the effects of absorption, scattering and the like still occur on the path, which leads to obvious dynamics of the simulation result.

In addition, the simulation result in the figure has an obvious characteristic that the spontaneous radiation effect of a part of the region in the center of the cloud is stronger than that of other regions in the simulation process. This is mainly because the central area is the location of the water vapor supply, so two consequences can result: 1. the phase change process generated in the central area is stronger than that generated in other areas, and the temperature is increased due to the heat generated by the water vapor condensation; 2. the concentration in the central region is higher than in the other regions and therefore the transmittance is relatively low. Since this result only affects the spontaneous emission, it is not obvious in the visible simulation but rather in the infrared simulation. This phenomenon can be effectively reduced by reducing the initial concentration of water vapor and reducing the water vapor supply area during simulation, but has an effect on the final simulation result of the three-dimensional model.

Furthermore, the terms "first," "second," and the like, are used for descriptive purposes only and are not to be construed as indicating or implying a relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defining "a first" or "a second" may explicitly or implicitly include one or more such feature. In the description of the present invention, the meaning of "a plurality" is two or more, unless explicitly defined otherwise.

In the description of the present specification, a description referring to terms "one embodiment," "some embodiments," "examples," "specific examples," or "some examples," etc., means that a particular feature, structure, material, or characteristic described in connection with the embodiment or example is included in at least one embodiment or example of the present invention. In this specification, schematic representations of the above terms are not necessarily directed to the same embodiment or example. Furthermore, the particular features, structures, materials, or characteristics described may be combined in any suitable manner in any one or more embodiments or examples. Further, one skilled in the art can engage and combine the different embodiments or examples described in this specification.

The foregoing is a further detailed description of the invention in connection with the preferred embodiments, and it is not intended that the invention be limited to the specific embodiments described. It will be apparent to those skilled in the art that several simple deductions or substitutions may be made without departing from the spirit of the invention, and these should be considered to be within the scope of the invention.

Claims

1. The modeling method of the infrared accumulation cloud based on the physical process is characterized by comprising the following steps of:

s1, acquiring a three-dimensional grid model of a cloud, wherein each grid in the three-dimensional grid model stores first variable information of a plurality of particles, the particles represent clouds at all positions in the grid, and the variable information comprises: drop density, temperature, and velocity;

s11, rendering the diffused clouds at the positions in the direction of the line of sight according to the pixel values, and repeatedly executing S8 to S11 until each position of the diffused clouds is rendered;

the step of S7 includes:

s77, accumulating the transmittance;

s81, accumulating the second radiation energy respectively;

2. Modeling method in accordance with claim 1, characterized in that the step of S2 comprises:

s21, setting initial cloud accumulation;

s23, moving each particle in the grid according to the initialized speed;

s25, calculating the droplet density and the temperature of the moved particles based on the current atmospheric pressure of the accumulated cloud;

3. The modeling method of claim 2, wherein the step of S26 includes:

Acquiring a Navie-Stokes equation as a fluid motion equation;

4. Modeling method in accordance with claim 1, characterized in that the step of S4 comprises:

acquiring a scattering phase function;

5. The modeling method of claim 4, wherein the step of S8 includes:

6. The modeling method of claim 5, wherein prior to the step of S8, the modeling method further comprises:

7. The modeling method of claim 1, wherein the step of S78 includes:

based on the extinction cross section, light transmittance is determined.

8. The modeling method of claim 1, wherein the step of S79 includes:

determining a primary scattering function of light passing through the cloud;

determining a secondary scattering function of light passing through the cloud;

9. The modeling method of claim 1, wherein the step of S83 includes: