Atmospheric Refraction and Propagation in Lower Troposphere (Electromagnetic Waves) Part 2

Refractivity statistics

As already mentioned, the physical processes in troposphere are complex enough to allow only statistical description of spatial and temporal characteristics of atmospheric refractivity. Nevertheless the statistics of important refractivity parameters such as an average vertical gradient are extremely useful in practical design of terrestrial radio paths when the long term statistics of the received signal have to be estimated, see (Rec. ITU-R P.530-12, 2009).

Average vertical gradient of refractivity

The prevailing vertical gradient of refractivity can be regarded as the single most important characteristics of atmospheric refractivity. According to (6), it is related to the effective Earth radius discussed above and it specifically determines the influence of terrain obstacles on terrestrial radio propagation paths. The examples of measured vertical profiles presented in the previous section show that the near-ground refractivity profile evolution is complex enough to not be described by only a single value of the gradient. The question arises what should be considered as a prevailing vertical gradient at a particular time. The gradient value is usually obtained from the refractivity difference at fixed heights, e.g. at 0 and 65 meters above the ground (Rec. ITU-R P.453-9, 2009). If more accurate data is available, the prevailing vertical gradient of refractivity can be calculated using a linear regression approach. Two year data (2008-2009) of measured vertical profiles were analysed by means of linear regression of refractivity in the heights (0 – 120 m) and the statistics of the vertical gradient so obtained were calculated. The results are in Fig. 7a where the annual cumulative distribution functions of the gradient are depicted. The quantiles provided by ITU-R datasets are also shown for comparison. It is clear that extreme gradients are less probable in reality than predicted by ITU-R. Linear regression tends to filter out the extreme gradients (otherwise obtained from two-point measurements) which do not fully represent the vertical distribution as a whole.

Fig. 7. Annual cumulative distributions of the vertical gradient of atmospheric refractivity obtained in 2008, 2009 (a), cumulative distribution obtained from the whole season (2 years) and fitted model (b).

Taking into account the importance of the gradient statistics for the design of terrestrial radio path, it seems desirable to have a suitable model. Several models of the gradient statistics were proposed, see (Brussaard, 1996), that can be fitted to measured data. Since they are often discontinuous in the probability density, they can be thought to be little unnatural. One can see in Fig. 7b where the two-year cumulative distribution is shown that the distribution consists of three parts: the part around the standard (median) gradient and two other parts – tails. Therefore the following model of the probability density f(x) and of the cumulative distribution function F(x) is proposed:

where the pi, fi and Oi are the relative probabilities, the mean values and the standard deviations of the Gaussian distributions forming the three parts of the whole distribution. Fitted model parameters (see Fig. 7b) are summarized in Table 2.

i	pi		Oi
1	0.086	-128.0	75.1
2	0.793	-46.1	11.8
3	0.121	-99.6	24.8

Table 2. Vertical refractivity gradient distribution parameters

Ducting layers

Although the ducting layers appearing in the first several tens or hundreds meters above the ground have significant impact on the propagation of EM waves on nearly horizontal paths, surprisingly little is known about their occurrence probabilities or about their spatial/temporal properties (Ikegami et al., 1966). This is true especially in the lowest troposphere where the usual radio-sounding data suffers from insufficient spatial and also time resolution. In the following, the parameters of ducting layers observed during the experiment are analysed by means of the modified Webster duct model.

An analytic approach to the modelling of refractivity profiles was proposed in (Webster, 1982). The refractivity profile with the height h (m) was to be approximated by the formula similar to the following modified model:

where the refractivity No (N-units), the gradient GN (N-units/m), the duct depth dN (N-units), the duct height h0 (m) and the duct width dh (m) are model parameters. A hyperbolic tangent is used in (10) instead of arctangent in the original Webster model because the "tanh" function converges faster to a constant value for increasing arguments than the "arctan" does. As a consequence, there is a sharper transition between the layer and the ambient gradient in the modified model and so the duct width values dh are more clearly recognizable in profiles. Figure 8 shows the meaning of the model parameters by an example where the modified refractivity profile is also included. It is seen from (7) and (10) that the model for modified refractivity profiles differs only in the value of the gradient: G = Gn + 0.157 (N-units/m).

Fig. 8. Duct model parameter definition with the values of parameters: N0 = 300 N-units, Gn = -40 N-units/km, dN = -20 N-units, ho = 80 m, dh = 40 m.

The above model was fitted to the refractivity profiles measured in between May and November 2010. More than 3 • 105 profiles were analysed and related model parameters were obtained. Figure 9 shows two examples of 1-hour measured data and fitted models. Significant dynamics is clearly seen in the evolving elevated ducting layers. It is also clear from the examples in Fig. 9 that the model is not able to capture all the fine details of measured profiles but it serves very well to describe the most important features relevant for radio propagation studies. Sometimes, the part or the whole ducting layer is located above the measurement range and so it is out of reach of modelling despite its effect on the propagation might be serious. This should be kept in mind while studying the statistical results presented below.

Fig. 9. The examples of time evolution of elevated ducting layers observed on the 1st of August 2010 at 00:00-00:50 (a) and on the 14th of July 2010 at 22:00-22:50 (b), measured data with points, fitted profiles with lines.

Figure 10 shows the empirical cumulative distributions of duct model parameters obtained from the fitting procedure. The medians (50% of time) of duct parameters can be read as N0 = 320 N-units, G = 116 N-units/km, dN = -2.2 N-units, h0 = 61 m, dh = 73 m. The probability distributions of N0 and G are almost symmetric around the median. On the other hand, the depth dN and width dh distributions are clearly asymmetric showing that the smaller negative values of the depth and the smaller values of width are observed more frequently. Almost linear cumulative distribution of the duct height h0 between 50 and 100 m above the ground suggests that there is no preferred duct height here.

Figure 10 shows the empirical cumulative distributions of duct model parameters obtained from the fitting procedure. The medians (50% of time) of duct parameters can be read as N0 = 320 N-units, G = 116 N-units/km, dN = -2.2 N-units, h0 = 61 m, dh = 73 m. The probability distributions of N0 and G are almost symmetric around the median. On the other hand, the depth dN and width dh distributions are clearly asymmetric showing that the smaller negative values of the depth and the smaller values of width are observed more frequently. Almost linear cumulative distribution of the duct height h0 between 50 and 100 m above the ground suggests that there is no preferred duct height here.

Fig. 10. The cumulative distribution functions of duct parameters obtained from measured profiles of atmospheric refractivity at Podebrady, 05/2010 – 11/2010.

Important interrelations between duct parameters are revealed by empirical joint probability density functions (PDF) presented in Fig. 11 – 15. The 2D maps show the logarithm of joint PDFs of all combinations of 5 parameters of the duct model (10). In these plots, dark areas mean the high probability values and light areas mean the low probability values. It is generally observed that there are certain preferred areas in the parameter space where the combinations of duct parameters usually fall in. For example, it is seen in Fig. 13a that the absolute value of the negative duct depth is likely to increase with the increasing gradient G. On the other hand, there are empty areas in the parameter space where the combinations of parameters are not likely to appear. One may find this information helpful when analysing terrestrial propagation using random ducts generated by the Monte Carlo method.

Fig. 11. The logarithm of the joint probability density function of duct parameters, obtained from measured profiles of atmospheric refractivity at Podebrady, 05/2010 – 11/2010.

Fig. 13. The logarithm of the joint probability density function of duct parameters, obtained from measured profiles of atmospheric refractivity at Podebrady, 05/2010 – 11/2010.

Fig. 14. The logarithm of the joint probability density function of duct parameters, obtained from measured profiles of atmospheric refractivity at Podebrady, 05/2010 – 11/2010.

Fig. 15. The logarithm of the joint probability density function of duct parameters, obtained from measured profiles of atmospheric refractivity at Podebrady, 05/2010 – 11/2010.

Modelling of EM waves in the troposphere

Several numerical methods have been used in order to assess the effects of atmospheric refractivity on the propagation of electromagnetic waves in the troposphere. They can be roughly divided into two categories – ray tracing methods based on geometrical optics and full-wave methods. The ray tracing methods numerically solve the ray equation (5) in order to get the ray trajectories of the electromagnetic wave within inhomogeneous refractivity medium. The ray tracing provides a useful qualitative insight into refraction phenomena such as bending of electromagnetic waves. Its utilization for quantitative modelling is limited to conditions where the electromagnetic waves of sufficiently large frequency may be approximated by rays. Geometrical optics description is known to fail at focal points and caustics where the full-wave methods provide more accurate results.

The full-wave numerical methods solve the wave equation that is a partial differential equation. Among time domain techniques, finite difference time domain (FDTD) based approaches were proposed (Akleman & Sevgi, 2000) that implement sliding rectangular window where 2D FDTD algorithm is applied. Nevertheless, tropospheric propagation simulation in frequency domain is more often. In particular , there is a computationally efficient approach based on the paraxial approximation of Helmholtz wave equation, so called Parabolic Equation Method (PEM), which is the most often used full-wave method in tropospheric propagation.

Split step parabolic equation method

We start the brief summary of PEM (Levy, 2000) with the scalar wave equation for an electric or magnetic field component y:

where k = 2n/X is the wave number in the vacuum and n(r,d,q>) is the refractive index. Spherical coordinates with the origin at the center of the Earth are used here. Further, we assume the azimuthal symmetry of the field, y(r,d,q>) = y(r,d), and express the wave equation in cylindrical coordinates:

where:

is the modified refractive index which takes account of the Earth’s radius R and where x = rd is a horizontal range and z = r – R refers to an altitude over the Earth’s surface. We are interested in the variations of the field on scales larger than a wavelength. For near horizontal propagation we can separate "phase" and "amplitude" functions by the substitution of:

in equation (12) to obtain:

Paraxial approximation is made now. The field u(x,z) depends only little on z, because main dependence of f(x,z) is covered in the exp(jkx) factor in (14). Then it is assumed that:

and the 1/(2kx)2 term can be removed from (15) since kx >> 1 when the field is calculated far enough from a source. We obtain the following parabolic equation:

An elliptic wave equation is therefore simplified to a parabolic equation where near horizontal propagation is assumed. This equation can be solved by the efficient iterative methods such as the Fourier split-step method. Let us assume the modified refractivity m is constant. Then we can apply Fourier transform on the equation (17) to get:

where Fourier transform is defined as:

From (18), we obtain:

and we get the formula for step-by-step solution:

The field in the next layer u(x+Ax,z) is computed using the field in the previous layer u(x,z):

Fourier transformation is applied in z-direction and the variable p represents the "spatial frequency" (wave number) of this direction: p = kz = ksin(§) and § is the angle of propagation.

The assumption that m is constant is not fulfilled, but equation (23) is used anyway. The resulting error is proportional to Ax and to horizontal and vertical gradients of refractivity. In practice, the value of Ax can be of several hundred wavelengths.

Application example and comparison with measured data

The parabolic equation method outlined above has been applied frequently to investigate the propagation characteristics on terrestrial (and also on Earth – space) paths under the influence of different refractivity conditions (Barrios, 1992, 1994; Levy, 2000) including the ducting layers described in the section 4.2. Users agree the method gives reliable results provided all the relevant details of terrain profile and of refractivity distribution are known and modelled correctly. This is however not always the case in practice. It is believed that the modelling results have to be compared with real world data whenever possible in order to validate the method under different propagation conditions and to know more about the expected errors due to incomplete knowledge of propagation medium.

Let us illustrate the particular example of conditions where the parabolic equation method performs successfully regardless the fact that refractivity profile along the propagation path is only roughly estimated. Figures 16a and 16b show the results of PEM propagation simulation performed using refractivity gradients measured during the 4th of November, 2008 at the receiver site. Sub-refractive conditions that occurred early morning caused a significant diffraction fading of more than 20 dB on the two lowest paths see Fig. 16b. On the other hand, the higher paths (receiving antennas located at 90 m and above) were not affected by diffraction effects.

Fig. 16. Spatial distribution of received power loss during sub-refractive condition on the path TV Tower Prague – Podebrady calculated by PEM (a), received signal levels measured in 5 receivers located in different heights and received signal levels modelled by PEM using time dependent vertical gradient of refractivity (b).

The results shown in Fig. 16b confirm that a very good agreement between PEM simulation and measurement can be achieved if the diffraction fading due to sub-refractive conditions (see time about 2:00) is the most important effect influencing the received power. It suggests that sub-refractive gradients are likely to be approximately the same along the whole propagation path and the approximation of horizontally independent refractivity, which is usually applied in PEM, is reasonable in this case. On the other hand, similar conclusion cannot be reached when multipath propagation occurs because only slight change in a refractivity profile along the propagation path may vary the received power distribution profoundly. These facts have to be kept in mind when the simulation results are interpreted.

Conclusion

Some results of the ongoing studies focussed on the propagation impairments of the atmospheric refractivity in the lowest troposphere were presented. Concurrent measurements of the vertical distribution of atmospheric refractivity together with the multi-receiver microwave propagation experiment were described. A new statistical model of vertical refractivity gradient was introduced. The unique joint statistics of ducting layers parameters were presented. The application of parabolic equation method was demonstrated on the example of a diffraction fading event. Simulated and measured time series were compared. A good agreement between simulation and measured data has been witnessed.

Future works in the area of the atmospheric refractivity related propagation effects should, for example, investigate the relations between the time evolution of duct parameters and multipath propagation characteristics, which is the area where only little is known at this moment.