Short-Wave Solar Radiation in the Earth’s Atmosphere Part 3 doc

52 Theoretical Base of Solar Irradiance and Radiance Calculation in the Earth Atmosphere

squares ξ

↓

2 (z) and ξ

↑

2 (z) with zeroth initial values and write together with

(2.14) the following:

ξ

↓

2 (z) := ξ

↓

2 (z)+[ξ↓(z)]2 , ξ

↑

2 (z) := ξ

↑

2 (z)+[ξ↑(z)]2 . (2.15)

Using the known expression for variance D(ξ) = M(ξ2) − M2(ξ), where M(. . .)

is expectation, we obtain:

D(ξ↓) = 1

K ξ

↓

2 (z) −  1

K ξ

↓

1

2

, D(ξ↑) = 1

K ξ

↑

2 (z) −  1

K ξ

↑

1

2

. (2.16)

The behavior of the distribution of random valuesξ↓(z) and ξ↑(z) is unknown.

However, the distribution of its expectations according to the central limit the￾orem tends to the normal distribution as K → ∞. Hence, desired irradiances

(2.13), which are also considered as random values, have the distributions

asymptotically close to the normal distribution. It is known that the normal

distribution is characterized with the expectation and the variance expressed

by (2.16). Forthe standard deviation (SD) (s(. . .) = √D(. . .)) of the irradiances

in accordance with the study by Marchuk et al. (1980) with taking into account

the known rule for the variances addition the following is obtained:

s(F↓(z)) = F0µ0



D(ξ↓)|K , s(F↑(z)) = F0µ0



D(ξ↑)|K . (2.17)

As follows from (2.17), the increasing of the number of trajectories K leads

to the minimization of the standard deviation (SD), i. e. of the random error

of the irradiances calculation. Evaluating the SD with (2.15)–(2.17) is of great

practical interest because it allows accomplishment of the calculations with

the accuracy fixed in advance. Actually, the calculation of the SD gives the

possibility of estimating the necessary number of photon trajectories and as

soon as the SD is less than the fixed value, the simulating is finished.

The above-considered scheme of the simulating of photon trajectories is

called “direct modeling” (Kargin 1984) as it directly reflects our implication

concerning photon motion throughout the atmosphere. However, direct mod￾eling is not enough for accelerating the calculation according to the algorithm

of the Monte-Carlo method or for the radiance calculation (Kargin 1984). Con￾sider two approaches to increase the calculation effectiveness that we have

applied. It is possible to find detailed descriptions of other approaches in the

books by Kargin (1984), and Marchuk et al. (1980).

The basis of optimizing the calculation with the Monte-Carlo method is an

idea of decreasing the spread in the values written to the counters. Then the

variance expressed by (2.16) decreases too and fewer trajectories are necessary

for reaching the fixed accuracy according to (2.17).

Assume that the photon could be divided into parts (as it is a mathematical

object and not a real quantum here). Then a part of the photon equal to

1 − ω0(τ

) is absorbed at every interaction with the atmosphere and the rest

ω0(τ

) is scattered and, then, continues the motion. During the interaction with

the surface these parts are equal to 1 − A and to A (A is the surface albedo)

Monte-Carlo Method for Solar Irradiance and Radiance Calculation 53

correspondingly. We specify the value w
called the weight of a photon (Kargin

1984; Marchuk et al. 1980), which it is possible to formally consider as a fourth

coordinate. Assume value w
= 1 in the beginning of every trajectory and

while writing to the counters, (2.12) will be assigned not unity but value w

.

Then the simulation of the interaction with the atmosphere is reducing to the

assignment w
:= w

ω0(τ

) at every step, and the simulation of the interaction

with the surface is reducing to the assignment w
:= w

A. Now the photon

trajectory can’t break (the surviving part of the photon always remains), the

break of the trajectory occurs only when the photon is outgoing from the

atmosphere top. Usually for not driving the photon with too small weight

within the atmosphere parameter of the trajectory break W is introduced: the

trajectory is broken if w
< W. It is suitable to evaluate value W based on the

accuracy needed for the calculation: W = sδ, where s is the minimal (over all

altitudes z for the downward and upward irradiances) needed relative error of

the calculation; δ is the small value (we have used δ = 10−2). This approach of

the photon “dividing” is known under the unsuccessful name “the analytical

averaging of the absorption” (Kargin 1984) (the words “analytical averaging”

are associated with a certain approximation, which is not used in reality).

Consider a photon at the beginning of the trajectory at the top of the atmo￾sphere. In this case, before the simulation of the first free path (τ
= 0, µ
= µ0,

w
= 1) using Beer’s Law (1.42) it is possible to account direct radiation, i. e.

radiation reaching level τ(z) without interaction with the atmosphere. For that

it is necessary to write to all counters ξ↓(z) the value depending on z instead

unity:

ψ = w
exp

−τ − τ

µ



, (2.18)

and further writing to the counters is not implemented for the first free path

(direct radiation). This approach is easy to extend to other parts of the trajec￾tory: before the writing of the free path to the counter, which the photon can

reach (ξ↓(z), for µ
> 0 and τ ≥ τ

, or ξ↑(z), for µ
< 0 and τ ≤ τ

) value ψ

calculated with (2.18) is writing and the further photon flight through the

counters is not registering. Note that as it has been shown above the exponent

in (2.18) is a probability of the photon started from level τ
to reach level τ.

This general approach of writing to the counter the probability of the photon

to reach the counter is called “a local estimation” (Kargin 1984; Marchuk et al.

1980).

The analysis of the above-described algorithm of the irradiances calculation

indicates that the irradiances are not depending on photon azimuth ϕ

. Actu￾ally, calculated only in two cases with (2.10) and (2.11), azimuth ϕ
does not

influence other coordinates and hence, the values written to the counters. Thus,

the “photon azimuth” coordinate is excessive in the task and it could be ex￾cluded for accelerating the calculations (but only in this task of the irradiances

calculations above the orthotropic surface).

Consider the second of the problems described above: the problem of ra￾diance I(z, µ, ϕ) calculation. It is obvious that the procedures either of the

54 Theoretical Base of Solar Irradiance and Radiance Calculation in the Earth Atmosphere

simulation of photon trajectories or of the calculation of the expectation and

variance are depending on the desired value, and hence they wouldn’t change.

The difference is concerning the procedure of writing the values to the coun￾ters. Encircle the cone with small solid angle ∆Ω(µ, ϕ) around direction (µ, ϕ).

We will be writing to the counter all photons, which have reached level z and

have come to the cone for radiance I according to the equation analogous

to (2.12). Moreover, in this case value 1||µ| has to be written to the counter

instead of unity in the case of the irradiance calculation to satisfy the link of

the radiance and irradiance (1.4). Pass further from the above-described (but

not realized) scheme of the direct modeling of the radiances to the schemes

of the weight modeling and local estimation. Let the photon have coordinates

(τ

, µ

, ϕ

). According to the definition of the phase function as a density of the

probability of the scattering (Sect. 1.2), the probability of the photon coming

to solid angle ∆Ω(µ, ϕ) after scattering at level τ
is equal to the integral of the

phase function over the angle intervals defined by (1.17) (i. e. ∆Ω and scatter￾ing angle (µ

, ϕ

)(µ, ϕ)) with taking into account normalizing factor 1|4π. Let

value ∆Ω decrease toward zero. Then we are revealing that the density of the

probability of the photon to reach direction (µ, ϕ) coincides with the value of

the phase function for argument χ
= cos( (µ

, ϕ

)(µ, ϕ)), which is computed

with (1.46). This probability is necessary to multiply by factor ψ defined with

(2.18), i. e. by the probability of the photon to reach level τ(z). Finally, the local

estimation for the radiance is obtained according to the results of the books by

(Kargin 1984, Marchuk et al. 1980).

ψ = w

4π



µ





x(τ

, χ

) exp

−τ − τ

µ



χ
= µµ
+



(1 − µ2)(1 − µ
2) cos(ϕ − ϕ

) .

(2.19)

Thus, the considered algorithm of the radiance computation according to the

Monte-Carlo method differs from the irradiance computation algorithm just

with the other equation for the local estimation (2.19) instead of (2.18) and with

other equations for the counters: for radiance over single trajectory ξ(z, µ, ϕ),

for expectation ξ1(z, µ, ϕ) and for the square of the expectation ξ2(z, µ, ϕ). Both

algorithms (for radiance andirradiance) could be carried out on computer with

one computer code. It is pointed out that the condition of the clear atmosphere

(the small optical thickness) has not been assumed so the Monte-Carlo method

algorithms can be also applied for the cloudy atmosphere.

In conclusion, illustrate that the considered algorithms actually correspond

to the solution of the equation of radiative transfer (1.47).

The desired radiation characteristic (radiance, irradiance) could be written

in the operator form according to expressions of the radiance through the

source function (1.52), and as per the link of the irradiance and the radiance

(1.4):

ΨB =

Ψ(u)B(u)du , (2.20)