Clear-Air Propagation Modeling using Parabolic Equation Method

50 M. GRÁBNER, V. KVIČERA, CLEAR-AIR PROPAGATION MODELING USING PARABOLIC EQUATION METHOD Clear-Air Propagation Modeling Using Parabolic Equation Method Martin GRÁBNER, Václav KVIČERA Lab. of Microwave Links, TESTCOM, Hvožďanská 3, 148 01 Praha 4, Czech Republic [email protected], [email protected] Abstract. Propagation of radio waves under clear-air conditions is affected by the distribution of atmospheric refractivity between the transmitter and the receiver. The measurement of refractivity was carried out on the TV Tower Prague to access evolution of a refractivity profile. In this paper, the parabolic equation method is used in modeling propagation of microwaves when using the measured data. This paper briefly describes the method and shows some practical results of simulation of microwave propagation using real vertical profiles of atmospheric refractivity. Starting with the ray equation and taking into account only vertical variations of refractivity, the ray bending can be expressed by: h n nd d)cos(1 ⋅ Φ −= ρ (1) where ρ is the radius of curvature of the ray, Φ is the elevation angle of the ray , h is the height over the Earth’s surface and dn/dh is the vertical gradient of refractivity. The value of refractivity n is close to one, and for near horizontal rays, it can be approximated by: h n d d1 −≈ ρ (2) Keywords When the ray path radius ρ is equal to the Earth’s radius R, then substituting to (2) and solving, we obtain: Clear-air propagation, refractivity profile, parabolic equation method, microwave communication. )(constant)( hm R h hn ==+ (3) where we introduced the modified refractivity m. When m(h) is used instead of n(h), the propagation is considered as it was over a flat Earth. For reference, the quantity M = (m-1)·106, which is called the refractive modulus, is often used in clear-air propagation studies. M is equivalent to radio refractivity N [1] and is usually used for depicting ducting layers in a vertical profile. 1. Introduction Spatial distribution of refractivity index of the air affects the propagation of electromagnetic waves in the troposphere. The path of the microwave ray is being bent, and the bending changes according to fluctuations of the vertical gradient of refractivity. Recently, a refractivity vertical profile was measured at the TV Tower Prague to obtain time evolution and statistics of the refractivity gradient [1]. 3. Parabolic Equation Method The parabolic equation method [2], [3] is used to compute the received signal level at a certain distance from the transmitter under various spatial distributions of refractivity. In this paper, we use the measured time series of the refractivity gradient as input data for a simulation of time series of the received signal level. The parabolic equation (PE) approximation is used to solve the Maxwell’s equations by means of a simplifying computation of the resulting partial differential wave equation. In this section, we outline the derivation of the method and describe a numerical technique for its solution. For a precise derivation see [4]. 2. Influence of Refraction 3.1 Paraxial Approximation The influence of refraction on the propagation is important mainly on horizontal or near horizontal paths due to the horizontal stratification of refractivity in the atmosphere. Now we will derive a parameter useful for characterization of the refractivity vertical profile. Let’s start with 3D scalar wave equation for an electric or magnetic field component ψ: 0 222 =+∇ ψψ nk (4) RADIOENGINEERING, VOL. 12, NO. 4, DECEMBER 2003 51 where k = 2π/λ is the wave number in the vacuum and n(r,θ,φ) 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, ψ(r,θ,φ) = ψ(r,θ), and express the wave equation in cylindrical coordinates: ∫ +∞ ∞− − ==≡ dzezxuzxupxUU jpz ),()],([F),( (12) From (11), we obtain: ),( 2 )1(),( 222 pxU jk mkp x pxU        −+− −= ∂ ∂, (13) 0),( 122 2 2 2 2=+ ∂ ∂ + ∂ ∂ + ∂ ∂ ψ ψψψ zxmk xx xz (5) )2/)1(())2/(( 22 ),( −− ⋅= mkjxkpjx eepxU , (14) where: Rzzxnzxm /),(),( += (6) and we get the formula for step-by-step solution: ( ) ),( ),( )2/)1(())2/(( 22 pxUee pxxU mkxjkpxj ⋅⋅ = ∆ + −∆∆− . (15) is the modified refractive index which takes account of the Earth’s radius R, and where x = rθ 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 following substitution: The field in the next layer u(x+∆x,z) is computed using the field in the previous layer u(x,z): x e zxuzx jkx ),(),( = ψ (7) ]),([F ),( ))2/((1)2/)1(( 22 kpxjmkxj epxUe zxxu ∆−−−∆ ⋅ = ∆ + . (16) in equation (5) to obtain: 0 )2( 1 12 2 22 2 2 2 2 =        +−+ ∂ ∂ + ∂ ∂ + ∂ ∂u kx mk x u jk x u z u (8) Fourier transformation is applied in z-direction and the variable p represents the “spatial frequency” (wave number) of this direction: p = kz = k sin(ξ) and ξ is the angle of propagation. The assumption that m is constant is not fulfilled however equation (16) is used regardless of that. It can be shown that the resulting error is proportional to ∆x and to horizontal and vertical gradients of refractivity. In practice, the solution converges as ∆x is reduced, and the value of ∆x can be of several hundred wavelengths. Finally, we make parabolic (paraxial) approximation. The field u(x,z) depends little on z, because main dependence of ψ(x,z) is covered in the exp(jkx) factor in (7). Then we assume: x u k x u ∂ ∂ << ∂ ∂2 2 2 (9) 3.3 Numerical Implementation and remove the 1/(2kx)2 term since kx >> 1 when the field is calculated far enough from a source. We obtain the following parabolic equation: 0)1),((2 22 2 2 =−+ ∂ ∂ + ∂ ∂uzxmk x u jk z u. (10) A Fast Fourier Transform (FFT) algorithm is used to implement iteration (16). It is necessary to restrict the z range (altitude). To avoid reflections from the upper boundary, the field is smoothly attenuated above a certain altitude value z < zmax. For a perfectly conducting surface, the boundary condition is u(x,0) = 0 for horizontal polarization and ∂u/∂z = 0 for vertical polarization. The conditions are fulfilled if u(x,z) is an odd or even function of z respectively. The implementation of impedance boundary conditions for imperfectly conducting surfaces is described in [4]. An elliptic wave equation is simplified to a parabolic equation when near horizontal propagation is assumed. This equation can be solved by step-by-step (iterative) methods more efficiently. The section, which follows, describes the Fourier split-step method for solution of (10). The transform size is chosen to prevent aliasing due to sampling the field (Nyquist’s theorem). The frequency and maximum angle of propagation determine the minimum transform size Nz according to: 3.2 Split-Step Fourier Method Initially, let’s assume the modified refractivity m is constant. Then we can apply Fourier transform on the equation (10) to get: z z Nz∆ =max , (17) )sin( maxmax ξ π π ⋅ =≤∆ kp z (18) 0)1(2 222 =−+ ∂ ∂ +− Umk x U jkUp (11) where ξmax is the largest propagation angle assumed. where Fourier transform is defined as: 52 M. GRÁBNER, V. KVIČERA, CLEAR-AIR PROPAGATION MODELING USING PARABOLIC EQUATION METHOD The initial field u(0,z) is calculated from the radiation pattern envelope of the transmitting antenna FRPE(p) using inverse FFT. Again, p = k sin(ξ) and the pattern is limited to angles below ξmax. The height of the transmitting antenna and the elevation angle can be easily taken into account using Fourier shift theorems. 4. Simulation of Propagation We will introduce the simulation of propagation over a flat Earth first, in order to better understand the relationship between the propagation and refractivity profile. Then we will show more practical results taking into account a real terrain profile. Fig. 2. Received signal level distribution – with earth curvature (the same parameters as Fig. 1) -250 -200 -150 -100 -50 0 0 10 20 30 40 50 60 70 80 90 100 relative received signal level (dB) range x (km) PE method free space loss Omnidirectional antenna is considered here, and the received signal level (RSL) is expressed as a relative value to the transmitting power. Calculations were made by GNU Octave (Matlab-like code). 4.1 Propagation over a Flat Earth First, standard atmosphere is assumed with a refractivity gradient G = dN/dh = -40 N-units/km. Figures 1 and 2 show the space distribution of received signal level which may be expected. The vertical axis z is the altitude and the horizontal axis x is the range (distance). The transmitting antenna is located on the left side in the height of 30 meters above a perfectly conducting ground and the electromagnetic wave propagates to the right. The frequency is 6.2 GHz and the polarization is horizontal. Fig. 3. Received signal level in the height of the transmitting antenna along the horizontal range 0 200 400 290 310 330 M (-) height (m) Figure 3 shows the range dependence of the received signal level calculated using the PE method and using free space attenuation. We can clearly see the shadowing of the field over the horizon. On the other hand, Figure 4 shows the resulting field in the presence of a ducting layer depicted in Figure 5. Between the heights of 120 and 190 meters, the refractivity gradient is G = -400 N-units/km. We can obtain the gradient of refractive modulus dM/dh = G + 157 (N-units/km) from the equation (3). Hence, there is a negative slope of the modulus M and a greater curvature of the wave than the earth’s curvature R. The wave is guided over the horizon. Fig. 4. Received signal level distribution – presence of ducting layer. Coordinates as in Figs. 1, 2 Fig. 5. Example of vertical profile of the refractive modulus in the presence of a ducting layer relative RSL: -80dB…….…-140dB….…….-200dB 4.2 Propagation over an Irregular Terrain In this section, a simulation of propagation above a real terrain profile from the west of Bohemia is presented. The transmitting site A is 516 meters above the sea level, site B is 788 meters a. s. l. The span of the link is 41.3 km and the link operates in the 4 GHz band. There is one critical obstacle along the path of the microwave beam that potentially corrupts the first Fresnel zone. We make a Fig. 1. Spatial distribution of received signal level over a flattened earth. Parameters: TX_height=30m, f=6.2GHz. Coordinates: xmin=0, xmax=100km, zmin=0, zmax=400m RADIOENGINEERING, VOL. 12, NO. 4, DECEMBER 2003 53 4.3 Time Evolution of Refractivity simulation at the frequency of 6.2 GHz because measured data exist from the receiving site for this frequency. In the course of refractivity measurement campaign at the TV Tower Prague [1], the time series of refractivity and of its gradient were obtained. An example of measured refractivity gradient is shown in Figure 9. Hereafter, we will assume the same series of the gradient appear at the transmitting site A. A number of PE simulations are performed with the vertical profile changing, and the time series of possible received signal level at the site B is calculated. Figure 8 explains the assumed structure of the vertical profile at the point of transmitter (A). This structure is assumed to be range independent. Figure 10 shows the resulting received signal in the height of 50 meters at the point B at the frequency 4 GHz. Figure 6 shows the coverage diagram calculated by the PE method, and Figure 7 compares the resulted field with the measured data for standard propagation conditions. The height of the trees growing at the highest path obstacle is estimated to fit the measured and simulated data. Then we calculate the field for the frequency 4 GHz and for standard and adverse propagation conditions. The adverse case is described by the effective Earth’s radius factor [1] (sometimes called refractivity coefficient) k = 0.75 that corresponds with the gradient G = 50 N-units/km. ground layer elevated layer 126m 191m M standard atmosphere height relative RSL: -80dB……….-140dB….…….-200dB Fig. 8. Two gradients (ground and elevated) measured at the TV Tower Prague are used for simulation of propagation over the terrain profile A - B. The heights are related to the transmitting site A. -400 -300 -200 -100 0 100 200 0:00 4:00 8:00 12:00 16:00 20:00 0:00 Time (hour) dN/dh (1/km) ground layer elevated layer Fig. 6. Received signal level distribution and terrain profile of microwave link A - B (flattened earth). Parameters: TX_height=54m, f=6.2GHz, coordinates: xmin=0, xmax=41.3km, zmin=360m, zmax=1121m (above sea level). 0 20 40 60 80 100 120 140 160 180 200 -155 -150 -145 -140 -135 height (m) received signal level (dB) PE method 6.2 GHz measured 6.2 GHz PE method 4.0 GHz PE method 4.0 GHz k=0.75 Fig. 9. Evolution of refractivity gradients on 9th August 2002 The resulted signal might be surprising. We would expect the fade at about 8:45 o’clock due to the large positive gradient in the elevated layer. However experiments show the strong dependence of the received signal on the relative height of the ducting layer to the height of the transmitting antenna and to the path obstacle position. It is clear that reliable vertical and horizontal distributions of refractivity are needed to use this approach successfully. In the example above, we assumed range independent distribution, which is rather unrealistic. In the vertical profile, an extent to higher altitudes is needed, too. Fig. 7. Vertical profiles of relative received signal level at the receiving microwave site B. 54 M. GRÁBNER, V. KVIČERA, CLEAR-AIR PROPAGATION MODELING USING PARABOLIC EQUATION METHOD A promising solution is to estimate the spatial distribution of the gradient from meteorological ground station data using meteorological objective analysis methods [5]. Larger extent of altitudes of radiosounding measurements [6] is suitable even for modeling propagation in radar applications. meteorological techniques to obtain real spatial distribution of atmospheric refractivity should be studied. Acknowledgements The research work described in this paper was supported by the Czech Science Foundation under the Grant No. 102/ 01/0029. -145 -140 -135 -130 0:00 4:00 8:00 12:00 16:00 20:00 0:00 Time (hour) RSL (dB) References [1] GRÁBNER, M., KVIČERA, V. Refractive Index Measurement at TV Tower Prague. Radioengineering, 2003, vol. 12, no. 1, p. 5 - 7. [2] DOCKERY, G. D. Modeling Electromagnetic Wave Propagation in the Troposphere using the Parabolic Equation. IEEE Transaction on Antennas and Propagation, 1988, vol. 3, no. 10, p. 1464 - 1470. [3] BARRIOS, A. E. A Terrain Parabolic Equation Model for Propagation in the Troposphere. IEEE Transaction on Antennas and Propagation, 1994, vol. 42, no. 1, p. 90 - 98. [4] KUTTLER, J. R., DOCKERY G. D. Theoretical Description of the Parabolic Approximation/Fourier Split Step Method of Representing Electromagnetic Propagation in the Troposphere. Radio Science, 1991, vol. 26, no. 2, p. 381 - 393. Fig. 10. Simulated received signal level at the receiving site B in the height of 50 meters (free space loss 136.8 dB) [5] SOKOL, Z., ŠTEKL, J. 3-D Mesoscale Objective Analysis of Selected Elements from SYNOP and SYRED Reports. Meteorol. Z., 1994, N.F.3, H.4, p. 242-246. 5. Conclusion [6] ŘEZÁČOVÁ, D., FIŠER, O., RAMÓN SÁEZ, L. Statistics of the Radio Refractivity Derived from Radiosounding Data. Radioengineering, 2003, vol. 12, no. 4, pp. 84-88. In this paper, the parabolic equation method for radio wave propagation modeling under clear air conditions was introduced and simulation examples were presented. The PE method is capable to take account of vertical and horizontal distribution of atmospheric refractivity, of real terrain profile and of electrical properties of the ground. About Authors... Martin GRÁBNER was born in Prague in 1976. He graduated from the Czech Technical University in Prague in 2000. He has been working at the Laboratory of Microwave Links in TESTCOM. In the literature, the importance of refractivity gradient in the first one hundred meters over the ground is usually mentioned for microwave communications. The simulations made and described above point out that a larger vertical extent should be taken into account, depending on relevant terrain profile. Václav KVIČERA was born in Poděbrady in 1948. He received MSc. and PhD. degree from the Czech Technical University in Prague. He is a Senior Researcher at the Laboratory of Microwave Links in TESTCOM. His research is concentrated on radio wave propagation above 1 GHz on fixed links. In our simulations, refractivity data derived from onepoint measurements were extended to the entire microwave path.. This is said to be quite pessimistic and unreliable approach. In the future, the possibility of utilization of