Path loss predictions in the presence of buildings on flat terrain a 3D vector parabolic equation approach.pdf

1716

IEEE TRANSACTIONS ON ANTENNAS AND PROPAGATION, VOL. 51, NO. 8, AUGUST 2003

Path Loss Predictions in the Presence of Buildings on Flat Terrain: A 3-D Vector Parabolic Equation Approach Ramakrishna Janaswamy, Fellow, IEEE

Abstract—Starting from a parabolic approximation to the Helmholtz equation, a three–dimensional (3-D) vector parabolic equation technique for calculating path loss in an urban environment is presented. The buildings are assumed to be polygonal in cross section with vertical sides and flat rooftops and the terrain is assumed to be flat. Both buildings and ground are allowed to be lossy and present impedance-type boundary condition to the electromagnetic field. Vector fields are represented in terms of the two components of Hertzian potentials and depolarization of the fields is automatically included in the formulation. A split-step algorithm is presented for marching the aperture fields along the range. Boundary conditions on the building surfaces are treated by using a local Fourier representation of the aperture fields. Several test cases are considered to check the boundary treatment used in the technique as well as to validate the overall approach. Comparison is shown with uniform theory of diffraction (UTD), exact solutions, as well as with measurements. Index Terms—Outdoor models, split-step algorithm, three– dimensional (3-D) propagation, vector parabolic equation.

I. INTRODUCTION ITE-specific methods have been considered in the past as an alternative to empirical models for predicting path loss in outdoor wireless environments. Some of the most widely used empirical models, such as the Hata model or the COST 231 model [8], are based on curve-fitting measurement results, and yield quick and reasonable statistical answers. These models do not take into account the details of the local terrain topography or building database in the prediction step. Most site-specific (or deterministic) models, on the other hand, attempt to solve an approximate version of the problem of electromagnetic wave propagation over terrain or in the presence of vertical buildings or combinations thereof. These methods yield results that are more representative of the local environment than the empirical ones, but are much slower than the latter. Still there is an interest in developing fast deterministic methods for better network planning and antenna siting. Among the site-specific methods, ray-based methods [1] are by far the most popular ones for outdoor urban environments, although some investigators have tried integral-equation-based methods [13], [7], [3] or finite-difference-type methods [9] over irregular terrain with some success. An entirely different paradigm has also emerged

S

Manuscript received April 25, 2002; revised August 7, 2002. This work was supported by the U.S. Army Research Office, Research Triangle Park, NC. The author is with the Department of Electrical and Computer Engineering, University of Massachusetts, Amherst, MA 01003-9292 USA (e-mail: [email protected]). Digital Object Identifier 10.1109/TAP.2003.815415

recently—rather than working with the Helmholtz equation (elliptic type) that leads to bidirectional waves and global coupling of fields, some formulations have considered the parabolic approximation to the Helmholtz equation, which assumes forward scattering at the outset [16], [10]. Under the assumption of forward scattering, the fields are only coupled along the principal direction, but not in the opposite direction. The problem of single or multiple knife-edge diffraction under the Kirchhoff’s approximation is one such example. Despite its seemingly flagrant assumption, the parabolic equation yields very useful results for long-distance propagation problems, as observed by several investigators (see [16] for a comprehensive reference list). Most of the formulations of parabolic equation are two-dimensional (2-D) in nature, in that they assume propagation to take place in a 2-D horizontal or vertical plane. The 2-D PE model presented in [10] had an average error of 10 dB and a standard deviation of 8–10 dB compared with field measurements. It is generally believed that a full 3-D formulation that includes both vertically and laterally propagating waves would not only reduce the mean error but also lower the standard deviation of error to around 6 dB [15], [14]. This is particularly true for small ranges [18]. Preliminary 3-D work with the scalar parabolic equation was reported in [23], where reflection and diffraction by one isolated building were considered. However, its implementation to realistic scenarios with multiple buildings is neither straightforward not feasible. A useful 3-D implementation of the scalar parabolic equation over irregular terrain with gentle slopes was recently considered in [25]. Finite-difference technique was used to discretize the parabolic equation, and the Crank–Nicholson scheme was used to march in range. As only scalar Helmholtz equation was used as the starting point in both of the above 3-D formulations, depolarization of waves was completely ignored. A natural extension of the 3-D scalar formulation to vector problems is to consider three scalar equations, one for each component of the electric field. Such an approach was taken by Zaporozhets [24], who also used Padé approximation and finite differences to discretize the various differential operators. Results were only shown for perfectly conducting buildings. The three scalar functions must be coupled at all points in space owing to the divergence-free condition of the electric field in a source-free region. This will effectively reduce the number of independent scalar functions to two. In [24], the divergence-free condition for electric field was only applied near an object and was indeed treated as a separate boundary condition. The number of unknown scalar functions outside the object was still regarded as three. It is not clear whether such an assumption has any theoretical justification.

0018-926X/03$17.00 © 2003 IEEE

JANASWAMY: PATH LOSS PREDICTIONS IN THE PRESENCE OF BUILDINGS ON FLAT TERRAIN

In this paper, we consider the vector problem of wave propagation in an urban area comprised of vertical buildings on a flat terrain. The buildings are assumed to be polygonal in cross section (horizontal cut), have flat rooftops, and finite height. A nonflat rooftop is, however, still within the scope of the present formulation. The building sides, its rooftop, and the terrain are made lossy to represent realistic conditions. Impedance boundary condition is used on the building walls and ground to uncouple the interior problem from the exterior one. The total fields are expressed in terms of the two components of Hertzian potentials, which become coupled at the building boundaries, thus, depolarization is included in the present formulation. This overall approach taken here is in contrast to the approach taken in [24], where three independent scalar functions are used to represent the vector fields. A principal direction ( axis) is first chosen along which only forward propagation is assumed. Waves are, however, bidirectional along the other two spatial axes. The principal axis is called the range axis; for a highly built-up area, this axis is entirely arbitrary, but remains fixed throughout the field computation once it is chosen. Under the assumption of forward propagation, the scalar Helmholtz equation for each potential component is approximated by a wide-angle, scalar parabolic equation. The split-step Fourier algorithm is used to march the potentials from one aperture plane) to the next along the range axis. When plane ( an aperture plane crosses a building, interior and boundary conditions are enforced within that region of the aperture plane occupied by the building. The process is repeated at several locations along the range until the building is exited by the marching aperture plane. Because of the finite widths and heights of the buildings and because of the lossy nature of the walls, the two potentials can get coupled at building surfaces. , The overall field is computed within a region is the maximum dimension in the lateral (or ) where is the maximum height direction about the transmitter, is the maximum range from the above the ground, and transmitter. Comparison is shown for test problems vis-a-vis ray methods and measurements to establish the validity of the overall numerical procedure. II. FORMULATION Fig. 1 shows a transmitter (marked T) and a receiver (marked R) located in a multipath environment. Waves arrive at the receiver via reflections and diffractions off of building sides, edges, rooftops, and ground. In contrast to the 2-D model developed in [10], which included only vertical plane rays, the model presented here will include also the laterally propagating waves (shown as dashed lines in Fig. 1) and is expected to improve the predictions of the field strength. The main source of the backward wave is a possible building directly behind the receiver, shown as a dash-dot line in Fig. 1. Even though the uniaxial parabolic equation used in this paper ignores the backward wave, its first-order contribution to the overall field can be determined separately by ray methods, for example, once the forward propagated aperture field is known at the range preceeding the receiver. In fact, the Walfish–Ikegami empirical model [1] is based on such reasoning. Reference [10] also discusses such an approach for a 2-D parabolic equation model.

1717

Fig. 1. Typical rays in an urban multipath environment.

In the ensuing theory, the transmitting antenna is assumed to be an extended, vertically polarized source, although a dual analysis would yield results for a horizontally polarized time convention, where is the source. An radian frequency of the wave at the frequency , is assumed for the sources and fields and is suppressed throughout. The wavenumber and intrinsic impedance in free space are denoted and , respectively. The free-space wavelength is deby noted by . The electric current density of the distributed source is assumed to be , where located at is its current moment and is the unit delta function. Expressing the fields in terms of two -directed potentials, it is easy to see from Maxwell’s equations that [4] (1) (2) and are the scalar potentials arising from the elecwhere tric and magnetic currents, respectively. The superscripts “ ” and “ ” identify the sources (electric or magnetic) that give rise and to the potentials and fields. The total fields are decomposed in terms of a TE mode and a TM mode with constituent fields (3) (4) (5) (6) The potentials are related to the sources through the inhomogeneous Helmholtz equations (7) (8) is the Laplacian operator and the magnetic current where is assumed to be zero here. The TE and TM modes propa-

1718


gate uncoupled on flat ground, but become coupled in the presence of lossy and finite-height buildings. In this work, we use impedance boundary conditions to represent the lossy nature of building surfaces and ground. If is the unit outward normal on a boundary with normalized surface impedance , then the impedance boundary condition states that [19] (9) For the most general case, the surface impedance will depend , where upon the material constitutive parameters is the dielectric constant and is the conductivity of the medium, and the incidence angle of the incoming wave. How, the ever, if the material is dense enough transmitted wave will be approximately normal in the medium for normally propaand it is reasonable to take the value of gating waves. The surface impedance for normally propagating plane waves is [2] (10) We will assume this approximation to hold on both the ground and building surfaces, but with different values of material parameters. (If propagation takes place only over a flat ground, it may be more appropriate to take the impedance to be that of a , , grazing wave [11].) On the flat ground with and , (9) gives (11) Expressing the various fields in terms of the potentials using reduce (3)–(6), the boundary conditions on the ground at to (12) and its dual (13) The boundary conditions on the rooftop of a building will take replaced by the same forms as (12) and (13) with

where and are the material constants for the rooftop. , On the vertical sides of a building, assuming , , , and expressing the various fields in terms of the potentials, we have from (9) that

On horizontal edges formed by the intersection of the rooftops and vertical sides of buildings, we combine the boundary conditions on the two faces. Of course, such a combination of boundary conditions will be approximate and ignores the edge condition that the true near fields have to satisfy. It may be noted that the satisfaction of edge condition by the near fields is not necessary to effectively predict the far-zone diffracted field. Diffraction by knife-edge under Kirchhoff’s approximation is a simple example. Hence, we regard the above approximation as no more severe than the impedance boundary condition itself, which will also become invalid near edges. If the buildings and ground are treated as perfect conductors (PEC), then , , and the boundary conditions become ground and building rooftops

building sidewalls

(16)

(17)

and on the horCombining these two we get izontal edges for a PEC building. It is seen from (16) and (17) that the TE and TM modes propagate uncoupled when the buildings and ground are treated as perfect conductors. Even for perfect conductors, polarization coupling will arise if the buildings have sloping rooftops or nonvertical sides. But such geometries are not considered in this paper. For the lossy case, mode coupling will occur even for flat rooftops and vertical sides. A. Parabolic Approximation Equations (7) and (8) together with the boundary conditions (12)–(15) and a radiation boundary condition at infinity constitute the exact boundary value problem. However, their numerical solution for long-range propagation problems requires excessive computational resources, which current technology cannot provide. This is due to the occurrence of full matrices that are the result of global coupling of fields. Reasonable answers for path loss can still be obtained by assuming forward propagation at the outset [10]. To do this, a principal direction is first chosen along which only forwardly propagating waves are permitted. Let this axis be the axis and be known as the range axis. In the transverse directions (i.e., along and ) waves are still allowed to propagate in both directions. In a source-free region, , with a second-order the Helmholtz equation derivative along the axis, is replaced with the wide-angle parabolic equation (PE) (18) which has a first-order derivative along the

(14) and its dual

(15)

axis, where (19)

is the transverse Laplacian operator. Derivation of (18) is simply based on factorizing the Helmholtz operator into two factors and retaining only the forward propagating operator [12]. Equation (18) models waves that travel in a forward cone about the axis and will be assumed to satisfy and is causal in . Both


1719

Fig. 2. Plan view of aperture plane marching over buildings. The dashed lines show the outline of the buildings, whereas the solid vertical lines show the intersections with the aperture plane.

this equation outside the source region. In the numerical solution of the parabolic equation by the split-step Fourier technique, the is determined in terms of the aperture aperture field at field at by decomposing the field in a spectrum of plane waves [20]. In the spectral representation, both incoming and outgoing plane waves along the -direction are chosen to automatically satisfy the impedance boundary condition on the ground. An important point to recognize is that when a wave field modeled by a uniaxial parabolic equation encounters a vertical obstacle in range, it only affects a local change to the aperture field. This is due to the causal nature of the forward-propagating fields along the range that prevents information about a future (in range) obstacle to flow backward. Boundary conditions on the building surface are then imposed on the field as the aperture plane crosses it. Due to the finite range-step size, they are, however, only applied at discrete values of the range as shown in Fig. 2. Between two adjacent aperture planes, we assume that the field marches as if it were propagating over flat ground. The situation is not unlike finite-differencing of differential operators on a discrete grid. The next subsection discusses the treatment of object boundary conditions in more detail. We assume that the frequency is high enough so that the surface wave on the ground could be neglected. In the far zone, the over flat ground could be written in terms potential at of the potential at as [11]

(20) where

is the wavenumber along the range axis for a plane wave travdirection, and is the 2-D Fourier transeling in the form of the aperture field that satifies the impedance boundary condition over ground

(21) These equations follow from extending the results for the axisymmetric case presented in [11] to a more general case with is transverse (i.e., -) variation. The propagator due to the accumulation of phase by a plane wave traveling in the . The quantity is positive -direction over a distance simply the plane-wave reflection coefficient over an impedance ground and is TM mode TE mode.

(22)

and with the Equations (20) and (21) are used for both chosen from (22). As with all paraappropriate value of bolic equation formulations, they will continue to be used near the source even though they have been derived for the far zone. The error due to this assumption will be most severe for plane waves propagating at very large angles near the source, but will decrease as the waves become shallow with respect to the axis. It can also be shown that, under this approximation, the initial field is [11] (23)

1720

Fig. 3.


Boundary treatment on an aperture plane intersecting a building.

Because the magnetic current is assumed to be zero here, the is simply initial condition for (24)

B. Treatment of Object Boundary Conditions A heuristic treatment of boundary conditions that captures the essence of the reflection and diffraction by building’s sides and rooftops is described in the following. Fig. 2 shows the discrete representation of buildings in our methodology. Boundary conditions are enforced on a building as often as it is crossed by the marching aperture plane. Let the intersection of an aperture plane with the building be termed as a strip. In this discrete representation, the space between adjacent intersections is taken to be empty. The number of strips depends on the range step chosen. For example, the first building in Fig. 2 gets intersected four times. It is to be noted that buildings are not replaced with physical strips having edges but only with mathematical strips. , whose right-hand Consider an intersection at corner is shown in Fig. 3 as dashed lines. Also shown in the figure are the aperture mesh points used in the numerical computations. The potentials on an aperture plane just before it hits the , where is an arbitrarily small number) building (i.e., at are computed using (20) and (21). Let this field be termed as the original field. The goal is to determine the actual field taking . The presence of into account the presence of a building at an obstacle will, in general, alter the aperture field not just on the building outline, but also in a band in the exterior region conlayers around it. The first layer is the building sisting of, say, th layer we postulate outline itself. On and beyond the the original field. The fields in this -layered region should be determined from the original field propagating over a distance along the range, the boundary condition on the obstacle, and layers, where . In known fields in the the following, in the interest of simplicity and speed, we adopt . The the lowest order approximation obtained by using

final justification of such an approach will be evident from the numerical results presented later. Because the field interior to the building is zero, we null it at all mesh points inside the strip. These points are indicated by open ovals in the figure. Exterior to the strip we maintain the original field. These points are indicated by crosses in the figure. On the boundary of the strip we assume that the potentials satisfy the boundary conditions given in (12)–(15). Equation (18) is used to express in terms of and . Hence, only transverse derivatives are involved in the boundary condition on the aperture field. The potentials at these boundary points (indicated by solid ovals in the figure) must be determined numerically from those at the surrounding points. In the split-step algorithm, the mesh points in the and directions are dictated by the highest angle of the wave (about the horizontal axis) that will be modeled. Large angles of propagation require a finer mesh. To include wide angles (up to about , where is the 45 ), the spacing required will be around free-space wavelength. With such large mesh sizes, the usual interpolating schemes, such as finite differences, tend to be inaccurate. We will, instead, consider spectral interpolation (or Fourier interpolation) by assuming periodicity of potentials in a 2-D region containing the boundary points. This is reasonable to do because the potential at a boundary point will depend mostly on its immediate neighboring points, and periodicity of the potentials in a remote region is not expected to influence at it much. Consider first the determination of the potential the boundary points on a rooftop at height . Let denote the vertically displaced axis and be the vertical mesh size. In the following, we assume that the rooftop height so that for an integer is rounded off to within . We express the electric potential as an -point Fourier series [17] in the vertical direction (25)

(26) where the prime after the summation indicates that the terms are multiplied by , , and the subscript denotes known quantities. The primed summation is used in the forward series so that the highest wavenumber is treated symmetrically during the analysis. The goal here is to given , . Using the representadetermine tion (25) into the boundary condition at directly yields

(27)


Equation (27) and a dual version of it are used to determine and at every point of the aperture plane lying on the rooftop. To determine the potentials on the vertical sides of the building, a 2-D Fourier representation is needed. Consider a points about the vertical subregion consisting of side as shown by the shaded region in Fig. 3. We assume that the width of this region is such that it only includes one side (viz., the right side, in this illustration) of the building. Furthermore, it must be high enough to include the rooftop of the building. Let us redefine the origin of the coordinate correspond to the system on the aperture plane so that vertical side where the potentials are being determined, and be the increment in the lateral direction. Note that the let unit normal on this face can have both and components. We use the 2-D Fourier representation in the shaded region ,

(28)

1721

and

where the superscript denotes transpose

(31) , , , , , and are The entries of the matrices , the number of equagiven in the Appendix. As tions in (31) is less than or equal to the number of unknowns. The systems of equations in (31) is then solved using pseudoinand verse or least squares techniques. After it is solved for , an inverse fast Fourier transform (FFT) will yield and . Only the first values pertaining to the boundary are finally used. C. Treatment of Aperture Truncation

(29) where

and

is the unknown potential. It is to be noted that sum-

in the definition of because the pomation starts at tential has been assumed to be zero below the ground (i.e., for ). The null potential for yields (30)

. Inserting (28) and A similar equation is obtained for its dual into the boundary conditions (14) and (15) for and into (12)–(15) for and evaluating the various derivatives, we arrive at two coupled sets of linear equations in the unknown vectors

The various Fourier integrals in (20) and (21) can be computed efficiently by means of 2-D FFTs. This will necessarily truncate the aperture in both the vertical and lateral directions. Spurious reflections from the fictitious boundary can arise if a simple rectangular truncation is employed. By carefully designing smoothing windows in the spatial and wavenumber domains it is possible to minimize these unwanted reflections. In the literature, these are referred to as sponge absorbers [12] and tantamount to using a lossy material near the periphery of the aperture. A smooth and gradually increasing conductivity profile in the lossy layer will attenuate the plane waves without causing undue reflections. These sponge type of absorbers have been successfully used in the PE community for well over 25 years. In this paper, we make the last quarter of space occupied by the aperture gradually lossy and use a Hanning type of is the aperture size in the -direcwindow. For example, if tion, the spatial window along is defined as (32) D. Propagation Factor For reference, we will compare the field in an urban environment with the field in free space. If the antenna source function at height is expressed in terms of a baseband funcby means of , then the tion far-zone potential in free space is (33) where (34) (35)

1722


and

(36) The -component of the free-space electric field in the far zone is (37) is defined as . In most of the Propagation factor numerical results shown below we use a Gaussian source with a linear phase taper of the form (38) is the standard deviation of aperture distribution and where is the elevation angle of main beam, with being the horizontal. The width of the main lobe in the elevation plane can be controlled by varying . For example, to produce [11]. an elevation beam width of 15 we choose III. NUMERICAL RESULTS Using the theory presented in Section II, we have computed the normalized field for several test geometries and compared the results with those obtained by using series solution, or uniform theory of diffraction (UTD) or with measurements. The test geometry for series solution is a lossy or perfectly conducting circular cylinder. The geometries considered for UTD comparisons were an infinitely tall, perfectly conducting, offset circular cylinder, an infinitely tall, inclined strip, a finite-height strip, a series of finite-height strips, a finite-height circular cylinder, and a finite-height square cylinder. The geometries involved in the comparisons with measurements include a centered, lossy rectangular block and an offset, lossy rectangular building. Some of the geometries chosen are general enough without having symmetries with respect to the source and axis planes that they provide the same challenging conditions that one encounters in practice. The important parameters which will affect the numerical , , and accuracy the PE solution are the grid sizes along , , and directions, respectively, and the 2-D FFT size . The increments and required on the aperture plane are dictated by the maximum propagation angle being is the maximum propagation angle with respect modeled. If , then, we have to the positive axis and if (39) Equation (39) can be used to predict the maximum propagation , angle for a given mesh size. For example, if and the maximum propagation angle modeled then without aliasing is 20.7 . and , it can be Using shown that the range step size over flat ground satisfies (40) , are chosen large enough to prevent spurious The sizes reflections off of the aperture boundaries from reaching the , useful computational domain. As an example, with

Fig. 4. Normalized field over flat ground versus lateral displacement showing the angular performance of wide-angle PE.

, we get , , and . At a frequency of 1 GHz, the maximum range will step over flat ground is about 362 m. In urban areas, be governed more by the building density and its representation rather than (40). Nevertheless, (40) gives the maximum range step one can ever use using the split-step algorithm. We would first like to illustrate the dependence of the wide ). angle capability of the algorithm on the mesh size ( , Fig. 4 shows the normalized field over flat ground at a distance . The numerical of 1 km from a 15 source located at parameters chosen in the PE computations are shown in the inset. The PE solution was computed out to a distance of 1 km and the results are shown for a horizontal cut. Comparison is shown with a two-ray model that includes a direct ray and a ground-reflected ray. The variable on the bottom abscissa is the lateral displacement, whereas the variable on the top abscissa is the equivalent propagation angle as computed from (34). As stated previously, the maximum propagation angle . This is treated without aliasing is 20.7 for clearly seen from the close agreement of the PE result with the . Beyond about 21 , aliasing two-ray model until about completely destroys the PE solution and only smaller mesh sizewill remedy the solution. A similar observation was also seen in the vertical cut. We now consider 2–D obstacles placed asymmetrically with respect to the source to provide rich conditions for lateral reflection and diffraction. The purpose of the following two test cases is to bring out the effectiveness (or otherwise) of the boundary treatment discussed in Section II-B. Both geometries will have vertical surfaces that are not parallel to the axis planes. The first test case we consider is that of propagation around a perfectly conducting circular cylinder with infinite dimension along the axis. No ground plane is assumed in this case. A line source is located at and radiates a vertically polarized field toward , as shown a circular cylinder of radius and center at


1723

Fig. 5. Direct and reflected rays from a line source exciting a PEC circular cylinder.

Fig. 7. A line source exciting an inclined conducting strip. The various shadow boundaries are shown as dashed lines.

Fig. 6. 2–D propagation around a laterally displaced circular cylinder.

in Fig. 5. All fields are invariant with respect to the axis and plane. The field is determined propagation takes place in the and compared with the exact series solution. at The purpose of this example is to demonstrate the effectiveness of the boundary treatment of solid objects. Boundary conditions 0.5 m. Other paon the cylinder were enforced at every , rameters are as shown in the inset of Fig. 6. At 10 m incident shadow boundaries (ISB-1 and ISB-2) exist at 30 m. Beyond this region, the field is comprised of and a direct ray and a ray reflected from the circular boundary of the cylinder. It is seen that even though the agreement with the series solution is not perfect, the PE solution predicts the field variation in the shadow zone and the illuminated zone very well. The agreement at large values of can be significantly improved . However, there is a limit to the improvement by decreasing will necessitate more marching steps and hence as a small

more frequent application of the spatial and wavenumber windowing. The latter have a tendency to distort the calculations at high angles. The next test case we consider is an inclined conducting strip as shown in Fig. 7. Once again, this geometry is chosen to show the effectiveness of the boundary treatment. Note in this case that the marching aperture plane intersects the obstacle only at one point, and, consequently, the field is modified at that one point only under our boundary treatment. Various incident and reflection shadow boundaries (ISB and RSB) are indicated in the figure. In Fig. 8, the normalized field is shown at a dis, and the various parameters used in tance the computation are shown in the inset. Comparison is shown with calculations performed by UTD. It is seen that the PE solution very nicely follows the UTD solution in all of the lit and shadow regions. Both of the preceding examples demonstrate that reflections by solid obstacles are adequately modeled even though the boundary conditions are applied only at discrete locations on them. The next comparison we show is for 3-D propagation around a finite screen. A perfectly conducting screen having a height and width is located at as shown in Fig. 9. at zero range The transmitting antenna is located at height and the ground is assumed to be perfectly conducting. The field is computed at a distance from the transmitter. We will show propagation factor comparisons with UTD in both the horizontal and vertical cuts. Fig. 10 shows the propagation factor (PF) as of a function of lateral displacement for a receiver height 2.1 m. The other parameters used in the numerical computations are shown in the figure inset. Also shown in the figure is

1724


Fig. 8. Reflected and diffracted fields of a line source exciting an offset and inclined conducting strip.

Fig. 10.

Propagation factor versus lateral displacement around a finite screen.

Fig. 9. A finite conducting screen between a transmitter and receiver. Fig. 11. Propagation factor versus vertical displacement behind a finite screen.

the result due to a 2-D parabolic equation, run for different radials about the transmitter. The 2-D model assumes propagation to take place in the meridian planes. If a radial line intersects the screen, then an infinite-width screen perpendicular to the radial line is erected at the intersecting point. In the shadow region, waves will reach the receiver only via diffraction from the top edge. In the 3-D model, there will be diffraction from the 500 m, the incident side and top edges of the screen. At 100 m. It is seen that the 3-D shadow boundary exists at PE results agrees very well with the UTD results both in the lit and shadow regions. In contrast, the 2-D model only predicts the average field variation and misses all of the interference lobes. Fig. 11 shows the field variation in the vertical cut as a func. Fine fluctuations, due to tion of the receiver height at the constructive and destructive interference of laterally propagating waves, are seen in both the UTD and PE results. The is believed to large spike seen in the UTD result at

be an aberration. (We have utilized numerical electromagnetics code-basic scattering code (NEC-BSC) to compute the UTD results.) Once again the 3-D PE results are in good agreement with the UTD results both in the lit and shadow regions. We have also considered multiple screens of varying widths and a perfectlyconducting square cylinder and compared the PE results with UTD results. The agreement obtained there was similar to the ones shown in Figs. 10 and 11. The next two examples we consider bring out the effectiveness of the boundary treatment in calculating the cross-polar component. The first example is that of a circular conducting placed over perfect ground building of radius and height along the axis from the transmitter. The at a distance of and components of the electric field are calculated at 500 m and for a –directed Gaussian source placed at . Fig. 12 shows the comparison of the transverse vari-


Fig. 12. Copolar and cross-polar components of the electric field behind a perfectly conducting circular building.

ation of the fields with UTD. All fields have been normalized to in free space. In this example, the cross-polar component reitself and no is excited. The circular cylinder sults from times. Other pawas intersected by the marching plane rameters are shown in the figure inset. It is seen that the agreeand . ment with UTD is within a few decibels both for It is well known that an infinite dielectric cylinder produces depolarization of waves when the incident wave is oblique to the cylinder axis [22]. The next example we consider is that of oblique scattering by a cylindrical building having a surface impedance . A Gaussian source producing a main beam at el30 and located at is used for this evation angle , a radius purpose. The cylindrical building has its axis at , and a height . The copolar field and the depolarized (resulting from the excitation of ) are calculated at field . The elevation beam width of the Gaussian source and the building height are chosen such that there is essentially no field incident on the building rooftop. Consequently, propagation is almost entirely off of the building side walls. This is done so that the exact series solution of a infinite cylinder could be used for comparison. The Green’s function for a point source and an infinite impedance cylinder was derived in terms of cylindrical harmonics and integrated over the Gaussian aperture to produce the exact field. The received fields are calculated in the direction of the main beam at and normalized to produced in free space. Fig. 13 shows 50.85 m and the comparison of the results for 200 m. Other parameters are shown in the figure inset. In the PE computations, the ground impedance is taken to be the same . Considering the fact that the depolarized field is about as 40 dB below the copolar field, the agreement seen in the figure is excellent over the wide range of azimuth angles considered (0–26.6 ). The example further testifies the effectiveness of the boundary treatment of Section II-B. We will next a show comparison with measurement for a lossy building-shaped obstacle. Van Dooren and Herben [21]

1725

Fig. 13. Depolarization due to oblique incidence on an impedance cylinder of finite height.

Fig. 14. Field strength measurement behind a lossy block-shaped obstacle, w 30 cm, d 7.7 cm, H 50 cm, x 68.9 cm, x 25.2 cm, H 46 cm, H 40 cm, and 114.8 cm.

=

=

=

= =

=

=

=

performed scaled model measurements and measured the field strength along a circular arc behind a lossy, block-shaped obstacle at a frequency of 50 GHz. The geometry of the measurement setup is shown in Fig. 14. A vertically polarized signal was generated using open-ended WR-90 waveguide. A comparison of propagation factor is shown in Fig. 15 as a function of the azimuth angle . Measured data have been read off the plots given in [21]. Both copolar, , and cross-polar, , components were computed using the 3-D PE algorithm described in the previous sections. The constitutive parameters used for the building material and ground are indicated in the figure text. In the PE comsuch that there were putations, we had a variable sections before the block, sections inside the block, and sections beyond the block. For the boundary treatment . The results were rather insensitive to the we chose over the range . It is seen that the choice of PE results match very closely with measurements having a mean error of only 0.5 dB and a standard deviation of error of 4.6 dB. Most of the contribution to the standard deviation comes from the deep shadow region where the fields are already attenuated 30 dB. Furthermore, it is seen for this case that the cross-polar component is about 30 dB below the copolar component.

1726


Fig. 17. Normalized field behind a rectangular building. Fig. 15. Propagation factor for copolar and cross-polar components behind lossy obstacle.

Fig. 16.

Site geometry for outdoor field strength measurements.

The final comparison we show is that of propagation behind a true building. The building we consider has vertical faces that are not parallel to any axis plane. Haslett [6] made outdoor measurement of field strength behind an isolated rectangular building, whose plan view is shown in Fig. 16. The receiver was moved horizontally through the entire shadow region, while keeping all the other parameters fixed. A vertically polarized transmitter was used at a frequency of 11.2 GHz. Fig. 17 shows the normalized field as a function of the lateral displacement. Four sections were chosen before and behind the building and six sections were chosen within the building. For the boundary , and . The results treatment we chose were insensitive to the choice of these parameters over the range , . The other parameters are as indicated in the figure inset. It is seen that agreement with measurement is very good and is indeed better than the agreement shown in [6] with UTD. The mean error in prediction with PE is 0.7 dB and the standard deviation of error is 3.7 dB. Once again, most of the contribution to the standard deviation comes from the deep shadow region where the PE result has deeper nulls. If only 2-D rooftop diffractions were considered, as in the 2-D PE model

considered in Fig. 10, the predicted normalized field strength would be largely independent of at approximately 28 dB and lead to highly erroneous results. A word about the computational resources required for the implementation of a full 3-D PE algorithm is fitting. Assuming range steps are used, the field will be calculated within a . At a frequency of 1 GHz volume , the wavelength is 0.3 m. Assuming that , , , 5 m, the field will be calculated in a volume 500 1228 307 m . For complex arithmetic the storage requirements per aperture plane per bytes of RAM. In the above example, potential are we would need 134 MBytes per potential component. The CPU time required would depend on the level of programming language. We used MATLAB to code the entire algorithm on a 1.7-GHz Pentium IV PC. It took about 84 s per march and a and total time of 85 min for all marches with . Most of the time was spent in computing the four 2-D FFTs required per march in the split-step algorithm. The overhead required in the CPU time for treating buildings was a very small fraction of the total time even when buildings were present during 60% of the marches. IV. CONCLUSION Starting from a parabolic approximation to the Helmholtz equation, we have presented a 3-D vector parabolic equation technique for calculating path loss in an urban environment. The buildings were assumed to be polygonal in cross section with vertical sides and flat rooftops and the terrain was assumed to be flat. Both buildings and ground were allowed to be lossy and present impedance type boundary condition to the electromagnetic field. Vector fields were represented by using two components of Hertzian potentials and depolarization of the fields was automatically included in the formulation. A split-step algorithm was presented for marching the aperture fields along the range. Boundary conditions on the building surfaces were treated by using a local Fourier representation of the aperture fields.


Several test cases were considered to check the boundary treatment used in the technique as well as to validate the overall approach. Reflection and diffraction from offset circular cylinder, offset strip, and finite cylinder as well as oblique scattering from lossy cylinder were considered to demonstrate the effectiveness of the boundary treatment. The effect of the mesh increments on the accuracy of the solution was discussed. Favorable comparison with UTD and measurement was shown for 3-D diffraction and reflection by a finite screen, a lossy rectangular block, and an asymmetric lossy building. The mean error relative to measurement was shown to be less than 1 dB and the standard deviation was shown to be less than 5 dB. Computational resources required to implement the algorithm was discussed. The method is very suitable for multiple buildings with no great addition in the computational resources required. It is believed that path loss predictions within a 1-km layout at the cellular frequencies could be made on a PC within reasonable times. As a byproduct of the split-step Fourier algorithm, where the aperture fields are expressed in terms of a spectrum of plane waves, the method is also capable of predicting the angular arrival of waves in the forward hemisphere. These will be extremely useful in the design of smart antennas [8]. It would be highly desirable to compare the PE results to field measurements in an environment having multiple buildings. We wish to initiate some measurement campaigns in the future to address this issue. More rigorous boundary treatments than the heuristic ones considered in the paper are certainly worth pursuing. The uniaxial technique presented in this paper ignores backscattering. The main source of the backward component that will be seen in a highly built-up area will be reflection from a building possibly behind the receiver. A first-order contribution of the backward wave can always be handled separately starting from the aperture field determined by the present method at a range step preceeding the receiver. A parabolic or a ray-based method may be used equation with respect to for this purpose. The 3-D parabolic equation used in the paper is wide enough to handle waves traveling 45 with respect to the axis. The angular capability is limited by the mesh increments chosen in the aperture plane and numerical dispersion arising from the structure of the mesh. A simple Cartesian mesh has been chosen in this paper. More advanced mesh structures, such as those proposed in [5], would be certainly worthwhile. Another area of future work is the incorporation of more advanced absorbing boundary conditions at the upper boundary. Even though well-designed spatial windows have been successfully used to minimize reflections from the truncation of the domain, more robust schemes will certainly be welcome. APPENDIX

1727

(43) (44) (45) (46) (47) (48) (49) (50) (51)

(52) (53) (54)

(55) (56) (57) the th row of the matrices element of the vectors and

, , , are given by

, and the th

(58)

(59) (60)

(61) (62)

With

(63) (41) (42)

(64)

1728


(65) (66) (67) (68) (69) (70) where

denotes element-by-element array multiplication. REFERENCES

[1] H. L. Bertoni, Radio Propagation for Modern Wireless Systems. Upper Saddle River, NJ: Prentice Hall PTR, 2000. [2] L. M. Brekhovskikh, Waves in Layered Media. San Diego, CA: Academic, 1980. [3] C. Brennan and P. J. Cullen, “Tabulated interaction method for UHF terrain propagation problems,” IEEE Trans. Antennas Propagat., vol. 46, pp. 738–739, May 1998. [4] W. C. Chew, Waves and Fields in Inhomogeneous Media. Piscataway, NJ: IEEE Press, 1995. [5] F. Collino and P. Joly, “Splitting of operators, alternate directions, and paraxial approximations for the three-dimensional wave equation,” SIAM J. Sci. Comp., vol. 16, no. 5, pp. 1019–1048, Sept. 1995. [6] C. J. Haslett, “Modeling and measurements of the diffraction of microwaves by buildings,” IEE Proc. Microwawe Antennas Propag., vol. 141, no. 5, pp. 397–401, Oct. 1994. [7] J. T. Hviid, J. B. Andersen, J. Toftgard, and J. Bojer, “Terrain based propagation model for rural area–An integral equation approach,” IEEE Trans. Antennas Propagat., vol. 43, pp. 41–46, Jan. 1995. [8] R. Janaswamy, Radiowave Propagation & Smart Antennas for Wireless Communications. Boston, MA: Kluwer, 2000. , “A fast finite difference method for propagation predictions over [9] irregular, inhomogeneous terrain,” IEEE Trans. Antennas Propagat., vol. 42, pp. 1257–1267, Sept. 1994. [10] R. Janaswamy and J. B. Andersen, “Path loss predictions in urban areas with irregular terrain topography,” Wireless Personal Commun., vol. 12, no. 3, pp. 255–268, Mar. 2000. Also see the addendum in Wireless Personal Commun., vol. 14, no. 3, pp. 303-304, Sept.2000.. [11] R. Janaswamy, “Radio wave propagation over a nonconstant immittance plane,” Radio Sci., vol. 36, no. 3, pp. 387–405, May-June 2001. [12] F. B. Jensen, W. A. Kuperman, M. B. Porter, and H. Schmidt, Computational Ocean Acoustics. New York: Amer. Institute Phys., 1994. [13] J. T. Johnson et al., “A method of moments model for VHF propagation,” IEEE Trans Antennas Propagat., vol. 45, pp. 115–125, Jan. 1997. [14] S.-C. Kim et al., “Radio propagation measurement and prediction using three-dimensional ray tracing in urban environments at 908 MHz and 1.9 GHz,” IEEE Trans. Veh. Technol., vol. 48, pp. 931–946, May 1999. [15] T. Kurner, D. J. Cichon, and W. Wiesbek, “Concepts and results for 3D digital terrain based wave propagation models–An overview,” IEEE J. Select. Areas Commun., vol. 11, pp. 1002–1012, July 1993. [16] M. F. Levy, Parabolic Equation Methods for Electromagnetic Wave Propagation. London, U.K.: IEE Press, 2000.

[17] A. Papoulis, The Fourier Integral and its Applications. New York: McGraw-Hill, 1962. [18] K. Rizk, R. Valenzuela, S. Fortune, D. Chizhik, and F. Gardiol, “Lateral, Full 3D and Vertical Plane Propagation in Microcells and Small Cells,” European Cooperation in the Field of Scientific and Technical Research, Switzerland, COST 259 Tech. Doc., TD(98)-47, 1998. [19] T. B. A. Senior and J. L. Volakis, Approximate Boundary Conditions in Electromagnetics. Piscataway, NJ: IEEE Press, 1995. [20] F. Tappert, “The parabolic equation method,” in Wave Propagation in Underwater Acoustics, J. B. Keller and J. S. Papadakis, Eds. New York: Springer-Verlag, 1977, pp. 224–287. [21] G. A. J. van Dooren and M. H. A. J. Herben, “Field strength prediction behind lossy dielectric objects by using UTD,” Electron. Lett., vol. 29, no. 11, pp. 1016–1018, May 1993. [22] J. R. Wait, “Scattering of a plane wave from a circular dielectric cylinder at oblique incidence,” Canad. J. Phys., vol. 33, no. 1, pp. 189–195, Jan. 1955. [23] A. A. Zaporozhets and M. F. Levy, “Modeling of radiowave propagation in urban environment with parabolic equation method,” Electron. Lett., vol. 32, no. 17, pp. 1615–1616, Aug. 1996. [24] A. A. Zaporozhets, “Application of the vector parabolic equation method to urban propagation problems,” Proc. Inst. Elec. Eng., Pt. H, vol. 146, pp. 253–256, 1999. [25] C. A. Zelley and C. C. Constantinou, “A three-dimensional parabolic equation applied to VHF/UHF propagation over irregular terrain,” IEEE Trans. Antennas Propagat., vol. 47, pp. 1586–1596, Oct. 1999.

Ramakrishna Janaswamy (S’82–M’83–SM’93– F’03) received the Bachelor’s degree in 1981 from REC Warangal, India, the Master’s degree in 1983 from IIT Kharagpur, India, both in electronics and communications engineering, and the Ph.D. degree in electrical engineering in 1986 from the University of Massachusetts, Amherst. From August 1986 to May 1987, he was an Assistant Professor of Electrical Engineering at Wilkes University, Wilkes Barre, PA. From August 1987 to August 2001, he was on the faculty of the Department of Electrical and Computer Engineering, Naval Postgraduate School, Monterey, CA. In September 2001, he joined the Department of Electrical and Computer Engineering, University of Massachusetts, Amherst, where is a Full Professor. He was a Visiting Researcher at the Center for PersonKommunikation, Aalborg, Denmark from September 1997 to June 1998. His research interests include deterministic and stochastic radio wave propagation modeling, antenna analysis and design, and applied electromagnetics. He is the author of the book Radiowave Propagation and Smart Antennas for Wireless Communications (Boston, MA: Kluwer, 2000) and a contributing author in Handbook of Antennas in Wireless Communications, L. Godara (Ed.) (Boca Raton, FL: CRC, 2001). Dr. Janaswamy was the recipient of the R. W. P. King Prize Paper Award of the IEEE TRANSACTIONS ON ANTENNAS AND PROPAGATION in 1995 and a recipient of the IEEE 3rd Millennium Medal from the Santa Clara Valley Section in 2000. He also received Certificates of Recognition for Research Contributions from ONR/ASEE in 1995 and from the Naval Postgraduate School, Monterey, CA in 1991. He is an elected member of U.S. National Committee of International Union of Radio Science, Commissions B and F.

Path loss predictions in the presence of buildings on flat terrain a 3D vector parabolic equation approach.pdf

Recommend Documents