School of Mathematics and Systems Engineering

Reports from MSI - Rapporter från MSI

October

2009

Methods for Path loss Prediction

Cem Akkaşlı

MSI Report 09067

Växjö University ISSN 1650-2647

SE-351 95 VÄXJÖ ISRN VXU/MSI/ED/E/--09067/--SE

Abstract

Large scale path loss modeling plays a fundamental role in designing both fixed and mobile

radio systems. Predicting the radio coverage area of a system is not done in a standard

manner. Wireless systems are expensive systems. Therefore, before setting up a system one

has to choose a proper method depending on the channel environment, frequency band and

the desired radio coverage range. Path loss prediction plays a crucial role in link budget

analysis and in the cell coverage prediction of mobile radio systems. Especially in urban

areas, increasing numbers of subscribers brings forth the need for more base stations and

channels. To obtain high efficiency from the frequency reuse concept in modern cellular

systems one has to eliminate the interference at the cell boundaries. Determining the cell size

properly is done by using an accurate path loss prediction method. Starting from the radio

propagation phenomena and basic path loss models this thesis aims at describing various

accurate path loss prediction methods used both in rural and urban environments. The

Walfisch-Bertoni and Hata models, which are both used for UHF propagation in urban areas,

were chosen for a detailed comparison. The comparison shows that the Walfisch-Bertoni

model, which involves more parameters, agrees with the Hata model for the overall path loss.

Keywords: path loss, prediction, wave propagation, rural, urban, Hata, Walfisch-Bertoni.

Email: [emailprotected]

mailto:[emailprotected]

Acknowledgments

I would like to express my sincere thanks to my supervisor Prof. Sven-Erik Sandström,

Växjö University for his support and helpful suggestions for this thesis work. I also would

like to express my special thanks to my family and friends.

Table of Contents

1. INTRODUCTION .................................................................................................................................. 1

2. THEORETICAL BACKGROUND ....................................................................................................... 2

2.1 RADIATED AND RECEIVED POWER ................................................................................................... 2

2.1.1 Radiated Power .................................................................................................................................... 2

2.1.2 Radiation Resistance and Received Power .......................................................................................... 6

2.1.3 Friis Transmission Equation ................................................................................................................ 7

2.2 PROPAGATION MODELING .................................................................................................................. 10

2.2.1 Overview of Channel Modeling ........................................................................................................ 10

2.2.2 Path loss Models due to Propagation Mechanisms ............................................................................ 14

2.2.2.a Path loss due to reflection and the Two Ray model ................................................................ 14

2.2.2.b Path loss due to diffraction ..................................................................................................... 19

3. PROPAGATION MODELS ................................................................................................................. 28

3.1 PROPAGATION MODELS FOR RURAL AREAS ............................................................................... 28

3.1.1 Deterministic Multiple Edge Diffraction Models ............................................................................ 28

3.1.2 Approximate Multiple Edge Diffraction Models ............................................................................ 30

3.1.2.a The Bullington method ........................................................................................................... 30

3.1.2.b The Epstein Petersen method .................................................................................................. 31

3.1.2.c The Japanese method ............................................................................................................... 32

3.1.2.d The Deygout method............................................................................................................... 33

3.1.2.e The Giovanelli method ............................................................................................................ 33

3.1.3 The Slope UTD method ................................................................................................................... 35

3.1.4 The Integral Equation approach ...................................................................................................... 42

3.1.5 The Parabolic Equation method ...................................................................................................... 48

3.2 PROPAGATION MODELS FOR URBAN AREAS .............................................................................. 52 3.2.1 The Okumura model ....................................................................................................................... 52

3.2.2 The Hata model ............................................................................................................................... 54

3.2.3 The Walfisch - Bertoni model ......................................................................................................... 55

4. SIMULATION ...................................................................................................................................... 63

4.1 The Hata model........................................................................................................................................ 63

4.2 The Walfisch - Bertoni model .................................................................................................................. 64

4.3 Comparison of the Hata model and the Walfisch - Bertoni model .......................................................... 66

APPENDIX A ............................................................................................................................................... 69

APPENDIX B ................................................................................................................................................ 79

APPENDIX C ............................................................................................................................................... 81

APPENDIX D ............................................................................................................................................... 85

REFERENCES ............................................................................................................................................. 89

1. INTRODUCTION

In radio propagation channels, spatial and temporal variations of signal levels are usually

observed on three main scales. Signal variation over small areas (fast fading), variations of the

small area average (shadowing), and variations over very large distances (path loss). Path loss

prediction plays a crucial role in determining transmitter-receiver distances in mobile systems

[1]. This thesis aims at describing various accurate path loss models that are used in rural and

urban areas.

In Chapter 2, the theoretical origins of the propagation phenomena and the received power

concept are introduced. The difference between the range dependent path loss and fast/slow

fading, as well as the basic path loss models are introduced. Basic models are also simulated

by using MATLAB.

Chapter 3 focuses on describing various accurate models for rural and urban areas,

respectively. For the rural case, deterministic and approximate multiple edge diffraction

models, as well as models with other approaches, are introduced. For the urban case the

Walfisch-Bertoni and the Okumura-Hata models are chosen for a detailed comparison. The

Walfisch-Bertoni model is described in detail via MATLAB simulations. A detailed

description of the repeated Kirchoff-Integral which is used in this model is given in Appendix

A. The standard Hata model is chosen as a reference for the results of the Walfisch-Bertoni

model.

In Chapter 4, the proposed models are compared for various parameters and simulated in

MATLAB.

1

2. THEORETICAL BACKGROUND

2.1 RADIATED AND RECEIVED POWER

2.1.1 Radiated Power

From the Maxwell equations for a hom*ogeneous, isotropic, linear and lossless dielectric

medium, the magnetic vector potential A can be obtained as [2, 3]:

' ,4

jkR

ve dv

R

J

A (2.1.1)

where,

| | ,R r r' (2.1.2)

is the distance from the source point to the field point. The free space wave number is,

0 0 .k (2.1.3)

Assuming a Hertzian dipole source one has the current as,

( ) Re{ } .j ti t Ie (2.1.4)

Figure 2.1.1 – The Hertzian dipole.

Considering Figure 2.1.1 and replacing 'dvJ in equation (2.1.1) with ˆ 'Idzz it can be written,

ˆ' .4

jkR

c

Ie dz

R

A z (2.1.5)

Assuming that the current is constant over the infinitesimal length l of the dipole ( )R r , and

that the point of observation is far away, one has,

2

ˆ .4

jkrIle

r

A z (2.1.6)

This expression suggests that the wave is propagating radially, in the direction of r̂ , with the

phase constant k. The amplitude of the wave is inversely proportional to the distance.

Since,

B H A , (2.1.7)

and in spherical coordinates,

ˆˆˆ (cos sin ) , z r θ (2.1.8)

one can write,

ˆˆ(cos sin ) ,4

jkrIle

r

A r θ (2.1.9)

and obtain the following expression for the magnetic field,

1 1 ˆ[ ] sin 1 4

jkrjkIle

r jkr

H A . (2.1.10)

For the electric field far away from the dipole, Maxwell’s equation ( 0J ) yields,

1[ ]

j E H , (2.1.11)

which can be written,

2 2 2

1 1 1 ˆˆcos 1 sin 1 ,2 4

jkr jkrIl jk Ile e

r jkr r jkr k r

E r θ

(2.1.12)

where the intrinsic impedance of free space is introduced as,

.

(2.1.13)

From the above equations it is understood that the electric and magnetic fields vary with the

distance r. If the kr product in the equations is much smaller than 1, i.e. kr<<1, the terms with

21

r and 31

r are predominant and the term

jkre goes to unity. After approximations the

following expression can be written for the magnetic and electric field:

3

2

sin ˆ .4

Il

r

H (2.1.14)

For the electric field some extra approximations can also be made for kr<<1. The term

1( 1)

jkr can be approximated by

1

jkr. Furthermore, substituting k

in equation (2.1.12)

with 1 , the electric field can be written as,

3 3

2cos sin ˆˆ ,4

Il

j r r

E r θ (2.1.15)

where the term with 31

r is the electric field intensity produced by a static electric dipole, and

hence called the electrostatic field component. Similarly the term with 21

r for H, represents

the magnetic field intensity produced by the very short filament of current, and is called as the

induction field component. Since r << 1/k this zone is called the near field. The time

dependent instantaneous vectors are,

( ) Re{ } andj tt e e E ( ) Re{ } .j tt e h H (2.1.16)

For the near field the instantaneous Poynting vector, which corresponds to the vector power

density (W/m2), is,

( ) ( ) ( ) Re } Re{ }j t j tt t t e e w e h {E H . (2.1.17)

One may write a vector depending on (x,y,z,t) in phasor form as a vector depending on (x,y,z)

independent of time.

By writing the real magnetic and electric field as,

*1( ) [ [ ]

2

j t j tt e e h H H , (2.1.18)

*1( ) [ [ ]

2

j t j tt e e e E E , (2.1.19)

the instantaneous Poynting vector can then be rewritten as follows:

2 *1( ) ( ) ( ) Re{[ ] [ ]} .

2

j tt t t e w e h E H E H

(2.1.20)

The time averaged vector power density then becomes,

*1( ) Re{ } ,

2t w E H (2.1.21)

4

which for the near field eventually yields,

2 22

2 5

1 | |ˆ( ) Re{ sin } 0 .

2 16r

I lt j

r

w r (2.1.22)

In equation (2.1.22), the -j indicates that the near zone has capacitive behavior, which means

that the dominant field is purely reactive and hence has zero average power.

On the other hand, assuming that the observation point is far away from the dipole (kr>>1),

the terms 21

r and 31

r get extremely small. The far field components then become,

sin ˆ ,4

jkrjkIle

r

H

(2.1.23)

sin ˆ .4

jkrjk Ile

r

E θ

(2.1.24)

From this, it is seen that the far field is a spherical wave with H and E fields that are

perpendicular, transverse and propagating in the r direction. In this case the medium’s

intrinsic impedance can be written as,

,E

H

(2.1.25)

and from equation (2.1.13) it is found that,

120 377 . (2.1.26)

For the far field time averaged vector power density one has,

* 21 1ˆ( ) Re{ } | |

2 2t E

w E H r , (2.1.27)

2 22 2 2

2 2

| |ˆ( ) sin (W/m ) .

32

I lt k

r

w r (2.1.28)

Equation (2.1.28) states that the power density in the far field is purely real and directed

radially outwards. Thus, this is also called the radiated power per unit area. The total

radiated power then becomes,

P ( ) .rad

s

t w ds (2.1.29)

Substituting equation (2.1.28) into (2.1.29) and integrating over a large sphere one has,

5

22 2 2

3 2 2 2

2

0 0

| |P sin | |

32 12rad

k I ld k I l

. (2.1.30)

In free space ( 120 ) one has,

2

2 2P 40 | |rad

lI

. (2.1.31)

2.1.2 Radiation Resistance and Received Power

One may seek a relation between the current and the total radiated power from the dipole.

Since the far field is purely real, the far field can be linked to a resistance radR called

radiation resistance,

21P | |

2rad radI R . (2.1.32)

Thus,

22

3rad

lR

, (2.1.33)

and in free space,

2

280rad

lR

. (2.1.34)

sV

Transmitter E

r

Receiver

radR

LossR

LZ

AX

I

(90 )

0V

Figure 2.1.2 – Voltage induced at the receiver antenna.

AZ

6

The scenario in Figure 2.1.2 assumes that the antennas are lossless Hertzian dipoles

( 0)LossR polarized in z direction. The electric field at distance r is,

sin

4

jkrjk IlE e

r

. (2.1.35)

The electric field impinges on the receiving antenna at an angle of (90 ) and induces a

voltage proportional to vertical component of the electric field and the length of the receiver

antenna l. So the induced voltage is,

0 0 sinV E l , (2.1.36)

where 0E from equation (2.1.24) is,

0 sin .4

jk IlE

r

(2.1.37)

To find out the average power delivered to the load in Figure 2.1.2 one can assume that the

load is matched to the antenna. That is to say that the antenna impedance *

A LZ Z . The total

impedance is 2 radR which yields maximum power to the load according to Jakobi’s law. The

power delivered to the load then becomes,

2 2 2 2

0 0

1 1P [ / 2 ] sin

2 8r rad rad

rad

V R R E lR

. (2.1.38)

2.1.3 Friis Transmission Equation

In antenna theory, the receiving antenna’s power capturing capability

( )

rer

PA

t w

, (2.1.39)

is defined as the ratio of average power received by the antenna’s load, to the time average

power density at the antenna, and is called effective area. Earlier the average power density

for the far-field was introduced as,

2

0| |( ) .

2

Et

w (2.1.40)

From the equations (2.1.38), (2.1.39) and (2.1.40) the effective area of the receiving antenna

then can be written as,

2 2sin .4

er

rad

A lR

(2.1.41)

7

By using equation (2.1.33) one obtains,

2 221.5sin

4 4er rA G

, (2.1.42)

where,

21.5sinrG , (2.1.43)

is the directive gain of the Hertzian dipole. Equation (2.1.42) shows that the receiving

antenna’s effective area is independent of its length and inversely proportional to the square

of the carrier frequency. At this point one can realize that the term frequency dependent

propagation loss is not the effect of wave propagation but the receiving antenna itself. The

average power density in terms of radiated power, transmitter gain tG and the distance r can

be written as,

2

P( )

4

rad tGt

r w . (2.1.44)

Considering,

rP ( ) ert A w , (2.1.45)

the equations (2.1.42), (2.1.44) and (2.1.45) yield the following formula:

2

rP P ( )4

rad t rG Gr

, (2.1.46)

which is called the Friis transmission formula and gives a relation between the power radiated

by the transmitting antenna and the power received by the receiving antenna. The Path loss

for the free space in dB then can be written as follows:

2

2 2

r

P( ) 10log 10log[ ] .

P (4 )

rad t rG GPL dB

r

(2.1.47)

The far-field (Fraunhofer) distance fR depends on the maximum linear dimension D of the

transmitter antenna [2],

22 .f

DR

(2.1.48)

For the distance fR to be in the far-field zone it should also satisfy and f fR D R .

Another practical choice is to determine a received power reference point 0R , which is

chosen smaller than any distance practically used in the mobile propagation, and also satisfies

8

0 fR R . Applying Friis formula provides a relation between power rP ( )R at any arbitrary

distance and power r 0P ( )R at the reference point,

2

0r r 0P ( ) P ( )

RR R

R

, (2.1.49)

where 0R R emphasizes the inverse square law. The power level in dBm is defined as the

received power with respect to 1 milli-watt [4]:

r 0 0r

P ( )P ( ) 10log 20log( )

0.001dBm

R RR

R

. (2.1.50)

In practice, 0R

could be 1 m for indoor environment and 100 m-100 km for outdoor

environment, for frequencies in the range 1-2 GHz.

9

2.2 PROPAGATION MODELING

2.2.1 Overview of Channel Modeling

Traditionally, mobile radio propagation models are categorized in two groups based on the

fading phenomena. The models that predict the overall average of the received signal strength

at a distance from the transmitter are called large-scale propagation models [4]. In general the

amount of damping is then called the path loss. Path loss is crucial in determining link

budgets, cell sizes and reuse distances (frequency planning). In this kind of modeling, the

mean signal strengths for arbitrary transmitter-receiver separation and large distances are

predicted. Therefore, large-scale propagation models are useful for estimation of radio

coverage areas, since they provide a general overview over long T-R distances. Basically, as

the mobile moves away from the transmitter, the local mean of the received signal decreases

gradually with slight variations. Here the local average of the received signal power is the

issue.

Large scale models account for reflectors, but not dense reflecting and scattering

environments. The local average received power is computed by averaging signal

measurements over a measurement track of 5 to 40 moving radially away from the base

station [4]. Measurement tracks for personal communications is generally from 1 m up to 10

meters long. Large-scale propagation curves are plotted based on these measurements. One

important thing to point out is that the dominant mechanism in large scale models is

reflection. Generally path loss models assume that the path loss is the same at a given T-R

distance without including the middle scale shadowing effects [5].

Shadowing (slow-fading, log-normal fading) is caused by obstacles between transmitter and

receiver. These obstacles attenuate the propagating signal power by absorbing, reflecting,

diffracting or scattering. Shadowing occurs over tens or some hundred meters and is also seen

as a part of large scale fading, which in fact can be also categorized as middle-scale fading.

Variation of the signal due to path loss occurs over very large distances, 100 m up to 1 km.

However in the case of shadowing the variation occurs over distances proportional to the

length of the obstructing object, typically in the order of building dimensions, 10 m up to 100

meters. Especially in cellular system design, when modeling the coverage area of a cell, path

loss and a shadowing model is combined. Figure 2.1.3 depicts contours of constant received

power due to path loss plus random or average shadowing.

10

Figure 2.1.3 – Combined path loss and shadowing.

The second group is called small-scale models that focus on the instantaneous behavior of the

received signal. These models are called small scale because here the issue is the variations of

instantaneous signal strength over very short time intervals or very short distances with

respect to the wave length. Moving the mobile the distance of one wavelength may cause a

variation in the received signal strength of up to 40 dB. Unlike the predictions in large scale

models, the variation of the signal is not a gradual decrease. The speed of motion is also a

parameter that effects the variation in the signal level. Especially in urban and heavily

populated areas the received signal is the vector sum of the scattered, reflected and diffracted

signals which arrive at the receiving point from various directions with various propagation

mechanisms. These propagation mechanisms cause rapid and dramatic variations for the

resultant received signal as the mobile moves a few wavelengths.

It is also important to point out that the dominant propagation mechanism in small scale

fading is not reflection but scattering. Due to the fact that the phases of the received signals

are random, the sum of all these different components behaves like noise. In this case the

Rayleigh fading model can be used. To overcome the fading, MIMO technology may be used

to achieve diversity by combining signals. The signal attenuation can be either Rayleigh

distributed or Rician distributed, depending on whether there is line of sight (LOS) or not. If

there are a lot of scattered components but no LOS then the attenuation coefficient can be

effectively modeled as Rayleigh distributed. If there is LOS then one can apply the Rician

distribution. The small scale fading model is a short term model and the severe signal

fluctuations usually happen around a slowly varying mean. Figure 2.1.4 illustrates the

received power variations with distance for both types of fading.

Path loss + random shadowing

Path loss + average shadowing

R

11

Figure 2.1.4 – Large scale versus small scale fading over distance.

Figure 2.1.5 illustrates the received signal strength variations in large and small scale models

over time. Red and green curves represent the means in large and small scale fading,

respectively. One can imagine that a person is holding a power-meter and taking signal

strength measurements while moving along. On top in Figure 2.1.5 large scale fading is

plotted while small scale fading is plotted at the bottom. Both fadings have slowly varying

means. But one can also see that even the mean variations in the small scale fading are faster

than in the large scale fading.

Figure 2.1.5 – Large scale versus small scale fading over time.

Figure 2.1.6 illustrates the same scenario and depicts the relation between the received signal

strength and the distance from the transmitter. If one takes a small segment from the path loss

curve and magnifies it, as in the middle figure, one finds a slow-fading (shadowing) curve.

However if one again magnifies a small section from the slow-fading curve and analyzes it in

detail one finds a fast-fading (short-term fading) curve. So fading can be subdivided into two

additional types of fading which of course have different mathematical descriptions.

Received Power (dBm)

T - R distance (meters)

Received Power strength (dBm)

Time (µs)

12

Figure 2.1.6 – Path loss vs. Slow-fading vs. Fast-fading.

As the mobile user moves away from the transmitter, measurements with a power-meter

probably yield a graphic similar to Figure 2.1.4 with a decreasing trend with fluctuations

around the mean. If the same user repeats this experiment he will get a similar but not the

same curve. Figure 2.1.7 also shows the received power variations of the overall signal and its

components as the receiver moves away from the transmitter. The horizontal axis corresponds

to the T-R distance and the perpendicular axis corresponds to the power in dB.

Figure 2.1.7 – Components of Overall Signal as a function of distance.

Fast-fading

Slow-fading

Path loss

Received Signal Strength (dB)

13

2.2.2 Path loss Models due to Propagation Mechanisms

2.2.2.a Path loss due to reflection and the Two Ray model

Reflection of an electromagnetic wave occurs when it impinges upon an object with different

electrical properties and very large dimensions compared to the wavelength [4]. The wave

impinging upon a new medium is partially transmitted into the second medium and partially

reflected back to the first medium. If the second medium is perfectly dielectric there is no

energy loss during reflection. If one assumes a second medium which is a perfect conductor,

all the incident energy is reflected back into the first medium without any loss. When the first

medium is free space and the permeabilities of two media are 1 2 , the electric field

intensity of the reflected and the transmitted waves are related by the simplified Reflection

Coefficient [4, 5]:

| | jr

i

Ee

E

, (2.2.1)

and,

2

2

2

2

sin cos for vertical polarization

sin cos

sin cos for horizontal polarization

sin cos

r i r i

v

r i r i

i r i

h

i r i

(2.2.2)

where,

i is the grazing angle as shown in Figure 2.2.1 and r is the complex dielectric constant of

the ground.

Reflection causes multipath effects. It may also be used (passive metallic reflectors) for

covering areas which are normally not covered. To know the electrical properties of the

reflector surface is important when it comes to modeling.

Materials in rooms that cause reflections are made out of dielectrics (walls, tables etc.). There

are also nearly perfect conductors such as metal frames in the windows that form perfect

reflectors.

At this point it is important to note that the statistical models begin to fail above 10 GHz.

Above 10 GHz one has to search for a deterministic model in terms of rays.

14

Figure 2.2.1 – Geometry of the two ray ground reflection model.

The two ray ground reflection model is a simple but useful model for predicting large scale

signal strength over long intervals. Ground reflection occurs when one normally has a LOS.

This model is accurate for predicting the large scale path loss for mobile radio systems with

large T-R distances (without any obstruction) and tall transmitter masts (exceeding 50 m) [4].

The surface wave is usually negligible in two ray models, since the surface wave extends

about 1 λ above ground and gives no contribution to the received field strength.

Figure 2.2.1 shows the geometry of the two ray ground reflection model in the case of

horizontal polarization [4, 5]. To find out the total field strength at the reception point R one

must first seek a formula for the phase difference between the direct field dE and the

ground reflected field rE . Since obviously depends on the difference in path lengths the

geometry in Figure 2.2.1, and imaging, yields:

2 2 2| | ( ) 1 ( )t rt r

h hTR d h h d

d

, (2.2.3)

2 2 2| | | | | | | ' | ( ) 1 ( )t rt r

h hTOR TO OR TR d h h d

d

, (2.2.4)

where d = | '' ' |T R .

Additionally one can obtain the Taylor series of 1 x for small x as,

21 11 ...

2 81 xx x

(2.2.5)

For x<<1 one has,

iE

dE

rE

T

R

'R

th

rh

O

i

2 1

d

''T

'T

15

1/2( )1

1 12

x x . (2.2.6)

In this case,

2( )t rh hx

d

or 2( )t rh h

xd

and t rd h h

. (2.2.7)

The difference in length then can be written as,

2 21 1 21 ( ) (1 ( ) )

2 2

t r t r t rh h h h h hd d

d d d

. (2.2.8)

The phase difference then becomes [6],

2 4 t rh hd

d

. (2.2.9)

Here is the phase delay due to the path length difference. But there is also another phase

delay which occurs at the point of reflection O. The reflection coefficient,

| | je , (2.2.10)

changes the amplitude and the phase of the impinging wave on reflection.

In the case of vertical polarization, for very large distances, the vertical vector components

dE and iE becomes almost the same. It can be seen from the geometry that,

1sin| ' |

d

TR and 2sin

| |

d

TR . (2.2.11)

Assuming 0i for very large distances one can see that | | | |TR TOR and 1 2sin sin

which results in equal vertical components of the electric fields in the z direction. In the case

of horizontal polarization, the incident and reflected electric field are also assumed to remain

the same, independent of the grazing angle.

The equality of the angles i and 0 can be derived by using boundary conditions from

Maxwell’s equations in a dielectric. In the case of vertical and horizontal polarization, when

0 ( =180 )i , the reflected wave rE is assumed to be 180° out of phase but equal in

magnitude with iE . The reflection then gives a phase shift | | 1 1j je e when

0i , for both vertical and horizontal polarization, in the case of grazing incidence.

16

Total electric field strength at the point R is,

(1 ) (1 )j j

tot d r d dE E E E e E e , (2.2.12)

(1 cos sin )tot dE E j . (2.2.13)

For the amplitude one can write,

2 2 2 2| | | | {(1 cos ) (sin ) } ,tot dE E

(2.2.14)

which simplifies to,

2 2 2| | 4 | | sin2

tot dE E

. (2.2.15)

Thus, the received signal can be written as,

2| | 2 | | sin( )t r

tot d

h hE E

d

. (2.2.16)

The average power density for the direct path was introduced as,

2

2

| |( )

4

rms t tE PGt

d w , (2.2.17)

so for the direct wave,

24d

PtGtE

d

, (2.2.18)

and in free space,

30 t t

d

PGE

d

. (2.2.19)

The resultant field can then be written as,

30 2| | 2 sin( )t t t r

tot

PG h hE

d d

. (2.2.20)

17

Figure 2.2.2 – Received field strengths in the Two Ray model.

Figure 2.2.2 was created in MATLAB and shows the change in total and direct field strengths

at the receiver due to distance. The chosen parameters are: transmitter height 60mth ,

receiver height 2mrh , transmitted power 5tP W , and the transmitter gain 1,5tG and an

ideal reflection coefficient 1 . The first graph depicts the scenario for 90 MHz, while in

the second graph the frequency is 900 MHz. Many observations can be made in equation

(2.2.20) and the corresponding Figure 2.2.2. It is obvious that the heights of the antennas, the

T-R distance and the wavelength determine the maxima and minima of the field strength. One

can also say that as the frequency is increased, more rapid variations occur in the resultant

field due to distance. The maximum received field strength caused by the interference is twice

the field strength of the direct wave. The dips in the graphs actually represent zeros ( dB ),

and cannot be shown clearly in Figure 2.2.2. For ,t rd h h one has,

2 2sin( )t r t rh h h h

d d

, (2.2.21)

and the field becomes,

2| | 4 30 t r

tot t t

h hE PG

d

. (2.2.22)

18

At this point an expression for the total average received power can be derived as follows:

2 22

2 22 4

r

(4 ) 30| |

P120 4

t rt t

rer

h hPG

E GdA

, (2.2.23)

which simplifies to the following equation:

2 2

r 4P t r

t t r

h hPG G

d . (2.2.24)

Equation (2.2.24) shows that for very large distances the resultant field strength falls faster

than the direct field. Independent of frequency, the received power decays as 4

1

d [4]. This is

good for avoiding interference in cellular systems. This simplistic but realistic model can hold

good for mobile application at GSM frequencies. The model is applicable in LOS scenarios

but LOS is seldom found as one moves very large distances away from the transmitter. Since

the antenna gains are fixed, the way to adjust the signal is usually to change the antenna

heights.

2.2.2.b Path loss due to diffraction

Geometrical optics doesn’t include the field in the shadow of an obstruction, since it is based

on rays. However, to some agree the E-M waves can propagate into the geometrical shadow

depending on the wavelength, shape, height and location of the obstruction. Especially in the

case of deterministic modeling diffraction models should be used.

A knife-edge geometry, depicted in Figure 2.2.3, is a special simplistic case of diffraction

where there is a single obstacle which is assumed to be extremely sharp so that no significant

reflections can occur. It has infinite width along z-axis (perpendicular to the page). Assuming

diffraction without any transverse effects one may reduce the scenario to 2-D. It is also

assumed that the obstacle in Figure 2.2.3 absorbs all the upcoming waves below its tip [4].

Figure 2.2.3 – Knife-edge geometry for a single obstacle.

h T

'h 1 'd

1d

R 2 'd

1d

1d

2d

19

According to the knife-edge diffraction geometry in Figure 2.2.3 the path difference can be

found as follows where 1 2' ', 'h d d

[4, 5]:

2 2 2 2

1 2 1 2' ' ' ' ( ' ')d h d h d d , (2.2.25)

which can also be written as,

2 2

1 2 1 2

1 2

' '' 1 ( ) ' 1 ( ) ( ' ')

' '

h hd d d d

d d , (2.2.26)

and the approximate path length is,

2 2

1 2 1 2

1 2

1 ' 1 ''(1 ( ) ) '(1 ( ) ) ( ' ')

2 ' 2 '

h hd d d d

d d

, (2.2.27)

2

1 2

1 2

' ( ' ') .

2 ' '

h d d

d d

(2.2.28)

Since in practice 1 2,h d d , the equations can be approximated in terms of distances parallel

to and normal to the ground, one has, 1 1 2 2', ' and 'h h d d d d . Hence the phase difference

becomes,

2

1 2

1 2

2 ( )h d d

d d

. (2.2.29)

Assuming that , and , the pitch angle, are very small,

1 2

1 2

' 'tan tan '( )

' '

d dh

d d

, (2.2.30)

so that approximately, in radians,

1 2

1 2

( )d d

hd d

. (2.2.31)

It is important to mention the Fresnel ellipsoids since obstructions that penetrate the first

Fresnel zone play a significant role in diffraction loss. If the propagating waves are considered

to bend at every point between T and R then the same path length forms an ellipsoid. A

Fresnel ellipsoid is a surface where the excessive path length is constant and has values,

2n

, (2.2.32)

20

where n is an integer. Figure 2.2.4 depicts Fresnel ellipsoids where ' ''TOR TO R TO R form

the nth

ellipsoid and T and R are the focal points of the ellipsoids. The nth

Fresnel zone is

usually defined with its maximum radius,

1 2

1 2

n

n d dr

d d

, (2.2.33)

which actually is the length of its minor half axis as represented by | '' |O M in Figure 2.2.4.

Figure 2.2.4 – Fresnel Ellipsoids.

Since the phase difference is 2

, the paths with 4 2

yield a phase difference

90 and the resultant phasor is affected destructively and becomes smaller than the LOS

phasor. The paths with 4

yield a phase difference 90 and the resultant phasor is

affected constructively and becomes greater than the LOS phasor. Figure 2.2.5 shows this,

Figure 2.2.5 – Destructive and constructive interferences due to the difference in path lengths.

resultant wave resultant wave

LOS wave LOS wave

secondary wave secondary wave

T

O

R

'O

''O

M

21

Often the Fresnel-Kirchhoff diffraction parameter is used as a normalized phase

difference,

1 2 1 2

1 2 1 2

2( ) 2

( )

d d d dh

d d d d

, (2.2.34)

or simply,

2v

. (2.2.35)

By applying Huygens principle to radio waves, one can predict the received signal in the

shadow region as the phasor sum of the electric field components of all secondary wavelets

from the obstacles [7]. Since Huygens principle assumes all these incoming and outgoing

wavelets to be spherical, one should use Fresnel diffraction.

It is important to emphasize that in Fraunhofer diffraction the incoming and outgoing wave-

fronts are assumed to be planar and the propagation paths after obstruction are assumed to be

parallel. So the phase differences between each consecutive element of the wave front

increases linearly which results in constant phase differences between consecutive elements of

the wave front. Also the paths taken by the waves are considered to be the same for very large

distances, again because of the parallelism, which results in equal amplitudes at a point of

observation. Thus the phasors of these elements form an arc of a circle at a point of

observation.

In Fresnel diffraction the outgoing wave-fronts are considered to be composed of tiny

spherical wavelets and the propagating wavelets’ rays from the obstruction axis to the point of

observation are not parallel. The phase differences do not increase linearly between two

consecutive elements.

Figure 2.2.6 – Huygen’s principle and Fresnel diffraction geometry.

P x

'y

P 1h

2h

3h

d

3d

h

y

22

In Figure 2.2.6, at point ,P the phase differences of the wavelets from the elements at

1 2 3, , ,..., nh h h h , where 2 1 3 2 1... n nh h h h h h , is directly related to the path

differences,

2

1( )n nd h h d d , (2.2.36)

which can be written in the form,

21 ( 1 ( ) 1)nn

h hd d

d

. (2.2.37)

Binomial approximation for nh d path yields,

2 2

1 11 ( )1 1

2 2

n nn

h h h hd d

d d

, (2.2.38)

and the phase difference for the nth

element is approximately,

2

1( )nn

h h

d

. (2.2.39)

This shows that the linear increase in height causes a quadratic increase in phase difference.

If one marks the consecutive points with a constant phase increment ' C , i.e. the points

1 2' , ' ... 'nh h h where 1 2 3 2 1' ' ' ' ' ... ' 'n n , one can write the following

expression:

2 2

1 1 1 1' ' ' [( ' ) ( ' ) ]2

n n n nh h h h Cd

. (2.2.40)

To satisfy the equation (2.2.40), one can see that 2 1 3 2 1' ' ' ' ... ' 'n nh h h h h h . This

means that for a linear step in phase, the step in height must decrease with n. Here the height

difference 1' 'n nh h is proportional to the radius and the effect of each wavelet at the

observation point. Hence, the contribution of each wavelet reduces for a linear increase in

phase difference [7]. At a point of observation, phasors of these elements form a curve called

the Cornu spiral.

In Figure 2.2.6 the incident plane wave can be defined as,

0( , ) jkx

incE x y E e . (2.2.41)

Applying the Kirchhoff integral and Huygens’s principle to the geometry in Figure 2.2.6 for

the region 0x , the total electric can be written as [8, 16],

23

(2.2.42)

where,

2

x

. (2.2.43)

and y is the axis of the point P in Figure 2.2.6. The integral can be written in terms of

Fresnel Integrals with the cosine integral,

2

( ) cos( )2

v

C v u du

, (2.2.44)

and the sine integral,

2

( ) sin( )2

v

S v u du

, (2.2.45)

and the Fresnel-Kirchhoff parameter,

'v y . (2.2.46)

Thus,

2( /2)

( ) ( )

v

j uC v jS v e du . (2.2.47)

The asymptotic expansions for 1v yield,

21 1( ) cos( )

2 2C v v

v

, (2.2.48)

21 1( ) sin( )

2 2S v v

v

. (2.2.49)

The limit,

1 1lim[ ( ) ( )]

2 2vC v jS v j

, (2.2.50)

yields,

2( /2)

1(1 )

2

j ue du j

, (2.2.51)

2( /4)

( /2)0( , ) ,2

j kxj u

tot

y

E eE x y e du

24

and the expression for the total field can be rewritten as [8],

( /4)

0 1 1( , ) [ ( ) ( )]

2 22

j kx

tot

E eE x y j C y jS y

. (2.2.52)

To illustrate the amplitude changes due to height for a fixed distance x from the edge, the

absolute value of the conjugate is of interest

1 1[ ( ) ( )]

2 2L j C y jS y . (2.2.53)

Figure 2.2.7 is a plot of the Cornu Spiral ( ) ( ) ( )F v C y jS y in MATLAB.

Figure 2.2.7 – The Cornu Spiral and the corresponding field strengths.

In Figure 2.2.7, for any value of v y , the distance L from the point ( )F v to (0.5 0.5)j

is proportional to the resultant signal. This is because the Cornu spiral in this case represents

the vector sum of the sources from the tip ' 0y to +∞. The length of the blue line from

( 0.5 0.5)j to (0.5 0.5)j represents the condition of LOS where v i.e. y in

Figure 2.2.6. The blue line from the origin to (0.5 0.5)j represents the condition where

0v which in Figure 2.2.6 corresponds to the point at the shadow boundary 0y . The green

line represents the first point where the field strength is equal to the LOS condition 1tot

inc

E

E

over the shadow boundary where 0 and v<0y . The red line represents the point P in Figure

2.2.6, somewhere below the shadow boundary where 0, 1y v .

25

The total field relative to the incident field is,

4 1 1[ ( ) ( )]

2 22

j

tot

inc

E ej C y jS y

E

, (2.2.54)

and using the identity,

41 1

2 2

j

e j

, (2.2.55)

it becomes,

[1 ( ) ( )] [ ( ) ( )]

2

tot

inc

E C v S v j S v C v

E

. (2.2.56)

Figure 2.2.8 represents the plots of the absolute values of the total field strengths relative to

the direct field tot

inc

E

E both in linear and decibel scales due to the parameter v in MATLAB.

Large positive values of v correspond to deep shadow.

Figure 2.2.8 – Relative field strength tot

inc

E

E

due to the parameter v in linear and dB scales.

The terms with Fresnel integrals in equation (2.2.56) makes the computation difficult.

Therefore, when predicting the knife-edge diffraction loss, equation (2.2.56) can be

approximated by the expression given by Recommendation ITU [8] for 0.7v ,

2

10 1020log ( ) 6.9 20log [ ( 0.1) 1 0.1] [dB]tot

inc

EJ v v v

E . (2.2.57)

26

Approximations for tot

inc

E

E in dB are also given by Lee [4],

10

0.95

10

( ) 0 -1

( ) 20log (0.5 0.62 ) -1 0

( ) 20log (0.5 ) v

J v v

J v v v

J v e

2

10

10

0 1

( ) 20log (0.4 0.1184 (0.38 0.1 ) 1 2.4

0.225( ) 20log ( ) >2.4

v

J v v v

J v vv

Figure 2.2.8 shows the diffraction loss due to height for a fixed frequency. However if the

frequency is increased for the same scenario the ratio tot

inc

E

E decreases more quickly. Figure

2.2.9 shows the diffraction losses for 1 GHz and 4 GHz at a point of observation x=100 m and

-15<y<15 meters. According to the figure, a higher frequency produces a more distinct

shadow.

Figure 2.2.9 - Relative field strength tot

inc

E

E

for 1 GHz and 4 GHz.

The first Fresnel zone is crucial for diffraction loss. If the frequency is increased for an

obstructed path, due to equation (2.2.33), the first zone gets narrower and the diffraction loss

may be reduced. However, as can be understood from Figure 2.2.9 and considering the

geometry in Figure 2.2.6, at higher frequency, for a fixed observation point x, moving the

point in -y direction or increasing the height of obstruction causes diffraction loss to happen

more quickly. The related parameters are,

1 2

1 2

2 2( )d dv y y h

d dx

. (2.2.59)

(2.2.58)

27

3. PROPAGATION MODELS

3.1 PROPAGATION MODELS FOR RURAL AREAS

3.1.1 Deterministic Multiple Edge Diffraction Models

When predicting the path loss between two fixed stations for large distances, the path profile

between the stations is often reduced to single knife edges since the wavelength is short

compared to the size of obstacles such as hills. The total path loss is then the free space loss

plus the predicted obstruction loss of the two-dimensional multiple edge diffraction model.

For multiple knife-edge geometry the waves incident on the edges will not be plane after the

first edge and the Fresnel integral should be evaluated consecutively for each edge to obtain a

deterministic loss prediction. Hence, for n edges an n-dimensional Fresnel integral must be

evaluated. This is very difficult to implement. A method with extreme accuracy for the

attenuation of electromagnetic waves by multiple edges is given by L. E. Vogler. Starting

from Furutsu’s generalized residue series formulation for electromagnetic waves over a

sequence of smooth cylindrical obstacles, Vogler derived an attenuation function for multiple

knife edges using Fresnel-Kirchhoff theory [12]. In his derivation he replaces the residue

series in Furutsu’s method with a contour integral since the series converges slowly in the

case of a grazing incidence. Vogler obtains an N-dimensional integral NI by considering the

cylindrical surfaces with zero that correspond to N multiple knife-edge geometry and

describes NC as the spreading loss factor over the propagation path in his formula [9]. The

attenuation NA , which is the ratio of total loss relative to the free space loss for N edges (see

Figure 3.1.1), is given as follows [9, 10]:

/2(dB) NN

N N NA C e I , (3.1.1)

where,

1

2

2 3

1

2

...

( )

n TN N

m m

m

d d d dC

d d

, (3.1.2)

2

1

1 1

2

1... ...

N

m

m

N N

f x

N N

x x

I e dx dx

, (3.1.3)

for 2N ; and 1NC for 1N where 1 1...T nd d d and,

1

1 2 1 1

1

( , ,..., ) ( )( ) ,N

N m m m m m

m

f x x x x x

(3.1.4)

28

for 2N ; and 0f for 1N .

2

1

N

N m

,

(3.1.5)

2

1 1 2( )( )

m mm

m m m m

d d

d d d d

, (3.1.6)

1 1

1

2( )

m mm m m

m m

jkd d r

d d

,

(3.1.7)

where 1r is the radius of the first Fresnel zone.

Figure 3.1.1 depicts the geometry in Vogler’s method where the given 1 3, , n are taken

positive since 1 2 3 4 1, > and > n nh h h h h h and 2 is taken negative since 2 3h h . nd is the

distance between ( 1) and nth thn edges.

Figure 3.1.1 – Multiple edge geometry for Vogler’s integral method.

Vogler transforms the multiple integral NI in equation (3.1.3) into an infinite summation

1

m

m

I

which can be computed numerically. The simplified solution for the attenuation NA

becomes,

/2

1

,NN

N N m

m

A C e I

(3.1.8)

where,

0 1

1

1 0 1 1 1 0

2 ( , ) (2, , ) ,m

m mm

m

m

I I m m C m m

(3.1.9)

1

3 n

2

1d

1

2d

1

3d

1

1nd

1

0h

1h

2h

3h

nh

1nh

29

( )!( , , ) ( , ) ( 1, , ) ,

( )!

jj i

N L N L

i

k iC N L j k I k i C N L i j

j i

(3.1.10)

with the repeated integrals of the error function,

21 1( , ) ( ) ,

!

n xI n x e dxn

(3.1.11)

and,

2

2 3 3 1 3 1 2( 1, , ) ( )! ( , ) ( , ) ,Nm

N N N N N N N NC N m m m I m I m

(3.1.12)

where the notations are given [9],

1 2

, ,

2 2, N 4

1, m , m 0 when k N-1

N L N L N L

N k

i m j m k m

L N

m

3.1.2 Approximate Multiple Edge Diffraction Models

Because of its computational rigor, Vogler’s method is usually used when comparing the

accuracy of approximate methods. Commonly used approximate methods are the Bullington,

Epstein-Petersen, Japanese, Deygout, and Giovanelli methods [10].

3.1.2.a The Bullington method

The terrain profile must be reduced to knife edges before the Bullington method is

implemented. To compute the diffraction loss, an equivalent knife edge is determined by

reducing all knife edges to a single knife edge. To determine the equivalent knife edge, two

lines are drawn joining transmitter and receiver to with respective dominant edges which has

the greatest angle of elevation. The intersection point of these two lines is the top of the

equivalent knife edge and the diffraction loss is calculated as if it were the only obstacle [7].

In Figure 3.1.2 the dashed red lines depicts the Bullington geometry.

(3.1.13)

30

Figure 3.1.2 – Geometry for the Bullington, Epstein Petersen and Deygout methods.

3.1.2.b The Epstein-Petersen method

Here, attenuation due to each obstacle is calculated and summed to find the total diffraction

loss. For instance, in Figure 3.1.2 the geometries of each triangle 12, 123T and 23R are taken

as a simple knife-edge geometry and diffraction losses are computed individually and

summed. For each case the obstruction lengths above the lines T2, 13 and 2R can be found by

using simple geometrical laws. This method results in large errors when the obstacles are

close [11]. In this case the Millington correction is added to the original loss,

10' 20log (cosec )L , (3.1.14)

where is the spacing parameter for the thn and ( 1)thn edges. The spacing parameter is

found from the formula [11],

1/2

1 1 2

1 1 2

( )( )cosec .

( )

n n n n

n n n n

d d d d

d d d d

(3.1.15)

The Epstein-Petersen method calculates the diffraction loss for N successive knife edges and

can be written as [12],

1

( ) | ( ) |N

epstein i

i

L N L v

, (3.1.16)

where,

2

1

3

T R

1d

1

2d

1

3d

1

3d

1

31

2

21

( )2

i

j t

i

v

jL v e dt

, (3.1.17)

is the diffraction loss for the thi knife edge. Since Vogler’s simplified attenuation function is,

vogler

1( ) ( ) | |

2

N

N m

m

L N C I

, (3.1.18)

then for m=0 the Epstein-Petersen loss can be described in terms of Vogler’s loss as,

vogler 0( ) |( )

m

epstein

N

L NL N

C

, (3.1.19)

which shows that Epstein-Petersen method turns out to be a first order approximation of

Vogler’s method and underestimates the NC spreading loss factor.

3.1.2.c The Japanese method

From the top of each obstruction a horizon ray is drawn which also passes from the top of the

previous obstruction. The intersection point of this horizon ray with the plane of the terminal

is then seen as the effective source. For instance, in Figure 3.1.3 the diffraction losses from

each obstruction are calculated due to the geometries of 12, T'23 and T''3RT . The sum of

these gives the total diffraction loss [11]. Actually this method is equivalent to the Epstein-

Petersen method plus the Millington correction given in equation (3.1.14).

Figure 3.1.3 – Geometry for the Japanese method.

''T

2

1

3

R 'T

T

32

3.1.2.d The Deygout method

First the v parameter is calculated for each obstacle as if there were no other diffracting

obstacles. Secondly, the dominant edge which has the maximum v parameter is determined

and the diffraction loss is calculated as if it were the only edge. The dominant edge now

becomes the terminal point of two sections divided by it. Then the process is repeated

recursively by finding out the maximum v parameter and determining loss until all edges are

considered. The total diffraction loss is then the sum of all these losses. For the geometry in

Figure 3.1.2 the total excess loss L in single edge diffraction terms can be computed by the

following expressions [13]:

1 2 3(dB)totalL L L L (3.1.20)

where

20log | |i iL A

(3.1.21)

2 21i

i

x

iA e e dx

(3.1.22)

where is calculated as in Vogler's method.i

If the edges are too many or too close, the Deygout method overestimates the path loss. In

these cases Causebrook introduces a correction which reduces the overestimated path loss as

follows:

1 2 1 2 1 3 1 3 3(6 )cos (6 )cos ,totalL L L L L L L L

(3.1.23)

where,

1 3 41

1 2 2 3 4

( )cos ,

( )( )

d d d

d d d d d

4 1 23

1 2 3 3 4

( )cos .

( )( )

d d d

d d d d d

3.1.2.e The Giovanelli method

Figure 3.1.4 relates to this method. As in the Deygout method firstly the dominant edge is

determined due to the v parameter. Edge M represents the dominant edge. Then a reference

point R’ is found by projecting a line on to the RR’’ plane which starts from M and passes

through the adjacent edge N. The loss due to the v parameter for the TMR’ geometry is then

(3.1.24)

33

calculated by obtaining the excessive effective height [10, 13]. The effective height for edge

M which is the excess height above TR’ is,

11 1

1 2 3

'' ,d H

h hd d d

(3.1.25)

where,

2 3 H h md and 2 1

2

.h h

md

(3.1.26)

So the loss from edge M can be written as a function of the TMR’ geometry parameters as,

1 2 3 1( , , '') .ML f d d d h h

(3.1.27)

After computing the loss from the dominant edge the dominant edge now becomes the

terminal point for two sections divided by it, as in the Deygout method. The same process is

applied to the remaining secondary edge N. The effective height for edge N which is the

excess height above MR is,

3 12 2

2 3

' .d h

h hd d

(3.1.28)

So the loss from edge N can be written as a function of the MNR parameters as,

2 3 2( , , ') .NL f d d h h

(3.1.29)

The total loss can then be written as,

.total M NL L L

(3.1.30)

For n edges this process is repeated until there are no edges left to be considered.

Figure 3.1.4 – Geometry for the Giovanelli method.

''R M N

'R

1 ''h 1h 2 'h

R 2h

T

1d

1

2d

1

3d

1

34

3.1.3 The Slope UTD method

Unlike the methods introduced previously this method uses the Uniform Theory of

Diffraction (UTD). The Slope UTD method estimates the propagation loss with much better

accuracy than the approximate multiple knife edge methods and with much shorter

computation time than Vogler’s method. The integral equation and parabolic equation

methods for terrain propagation uses the actual terrain profile but they do not provide an

optimum solution for a point to point prediction. The Helmholtz-Kirchoff theory of waves has

many computational difficulties and especially in the case with a large number of obstacles it

gets too complex to apply [14, 15, 16].

Geometrical optics (GO) is a useful method since it views propagating radio waves as light

rays for very small wavelengths, but neglects diffraction and predicts zero fields in the

shadow regions. A high frequency technique, the Geometrical theory of diffraction (GTD)

developed by J.B. Keller, is an extension of geometrical optics, basically by inclusion of

additional diffracted rays. The GTD method can accommodate 3D environments and can

predict diffracted rays very rapidly for complicated geometries. But because the GTD method

is a ray-optical theory it is valid only when outside the transition region which is the region

within the GO incidence and reflection shadow boundaries. That is because around shadow

boundaries the field has rapid spatial variations. For a perfectly conducting edge Keller's GTD

results in discontinuity problems in the vicinity of incidence and reflection shadow

boundaries.

The UTD method developed by Kouyoumjian and Pathak is a development of GTD for

perfectly conducting wedges which overcomes GTD’s transition region limitations. By

identifying a region around the shadow boundary, performing an asymptotic analysis, and

multiplying the diffraction coefficient by a transition function ( )F X , which involves a

Fresnel integral, they created a ray based theory valid at all spatial locations. It is proven that

UTD provides a more accurate prediction than GTD. Luebbers also developed this theory for

wedges with different dielectric materials by modifying the UTD coefficients, reducing the

errors in the case of grazing angles. The UTD method can also be extended to multi-shaped

canonical objects such as wedges, edges and cylinders. The UTD generally fails when the

transition regions or GO caustics overlap. For such regions the Physical Theory of Diffraction

(PTD) may be used for more accurate calculations.

However, UTD and GTD both fail when the incident field is not a ray optical field. For

example in the case of double wedge diffraction the wave impinging upon the second wedge

is already a diffracted field from the first wedge. One then cannot just apply UTD theory for

this wedge since it is illuminated by the transition region field, not the incident field from the

transmitter.

Therefore, the Slope UTD approach was introduced by Andersen as an extension of UTD for

the cases where the incident field has a rapid spatial variation i.e. has non-zero gradient and

even zero amplitude [10, 17]. The slope UTD uses a high order term asymptotic expansion

35

and includes higher-order diffraction unlike classical UTD. Hence, just after diffraction, the

total created field at a point is composed of a first order diffracted field which is proportional

to the amplitude of the incident field and a slope diffracted field which is proportional to the

derivative of the incident field.

Figure 3.1.5 depicts a slope UTD diffraction geometry with a wedge, where 'TT is the

incident shadow boundary and 'R RT T is the reflection boundary [17].

Figure 3.1.5 – UTD diffraction geometry.

According to the slope UTD theory the total diffracted field strength behind a single

absorbing half-plane is [14]:

( ) ( ) ( )s jksid i

EE E D d A s e

n

, (3.1.31)

where iE is the incident field and,

' , (3.1.32)

is the angle between the incident and the diffracted ray, i.e. the angle above the shadow

boundary.

/42( ) [2 cos ( / 2)]

2 2 cos( / 2)

jeD F kL

k

, (3.1.33)

is the amplitude diffraction coefficient.

( )( )

sA s

s s s

, (3.1.34)

is the spreading factor and,

1 ( )( ) ,s D

djk

(3.1.35)

T

'

RT

R

'RT

's s

'T

n

36

is the slope diffraction coefficient. The UTD transition function is given by,

2

( ) 2 .jx ju

x

F x j xe e du

(3.1.36)

By using the derivative of ( )F x ,

( )'( ) [ ( ) 1] ,

2

F xF x j F x

x (3.1.37)

the slope diffraction coefficient sd then can be rewritten as,

/4

( ) cos( / 2)(1 ( ))2

js

s

ed L F x

k

(3.1.38)

where the L and sL parameters are the distance factors calculated according to some

continuity equations that vary for each ray due to geometry.

Figure 3.1.5 – Slope UTD geometry for 2 wedges.

Figure 3.1.5 depicts a multiple-diffraction scenario that is suitable for the slope UTD method.

One can describe the field impinging on the left face of the first wedge 1W in terms of the

radiated field strength for 01 ms

which is defined as

00 1|sE . That is to say the incident field

at the vertex point of 1W then can be described as,

00 1

0 1

|( ) ,

sEE W

s

(3.1.39)

since the distance from the transmitter T to the vertex point of the first edge is 0s . The fields

from the first wedge due to classical UTD method can then be described [14]. Multiplying the

0s 1s

1s

1

2 n

1,1P 1,2P

1s 2s

T 1 '

2 ' 0E 1E

2,1P 2s

1W 2W

R 2s

2E

37

diffraction coefficients iD and spreading factors iA for each wedge the fields can be written

as follows:

1 2 0 1 1 2 1 1( ) ( ) ( ) ( )E W E W D W A s

(3.1.40)

is the field computed at 2W which is diffracted by the vertex of 1W .

2 1 2 2 2 2( ) ( ) ( ) ( )E R E W D R A s

(3.1.41)

is the field computed at the receiver R which is diffracted by the vertex of 2W .

Once 0E is diffracted, 1E

impinges on 2W as an already diffracted field since 2W

is in the

shadow. In the classical UTD method the slope diffraction term is neglected. However, in this

double diffraction scenario, better accuracy is obtained, if slope diffraction is taken into

account. Therefore the equations can be rewritten as follows for the slope UTD case as

follows:

1 2 0 1 1 2 1 1( ) ( ) ( ) ( )E W E W D W A s , (3.1.42)

is the field computed at 2W which is diffracted at the vertex of 1W . There is no slope term in

equation (3.1.42) as in equation (3.1.40), since the wave diffracted at the vertex of 1W is not

an already diffracted wave before impinging on 1W .

1 22 1 2 2 2 2 2

( )( ) [ ( ) ( ) ( )] ( ) ,sE W

E R E W D R d R A sn

(3.1.43)

is the field computed at the receiver R which is diffracted at the vertex of 2W . 1 2( )

E W

n

is

the directional derivative of 1E in the normal direction of the incident shadow boundary of

1E

and 2 ( )sd R

is the slope diffraction coefficient which is derived from the amplitude

diffraction coefficient as,

1 ( )( )s D

djk

where ' . (3.1.44)

Since,

/42( ) [2 cos ( / 2)]

2 2 cos( / 2)

jeD F kL

k

, (3.1.45)

and the derivative of the transition function is,

38

( )'( ) | ( ) 1 |

2

F xF x j F x

x , (3.1.46)

the slope diffraction coefficient is written as,

/4

sin( / 2)[1 ( )] .2

js

s

ed L F x

k

(3.1.47)

The directional (normal) derivative of the incident field iE

n

can also be defined in

cylindrical coordinates for the 2D case as follows:

1 ˆ ˆˆ ˆ( )i i i is s

s

E E E Ez

n z

, (3.1.48)

which reduces to,

1 1i i i

s s

E E E

n s

. (3.1.49)

Therefore, the second diffracted field in slope UTD can be rewritten as:

1 22 1 2 2 2 2 2

1 1

1 ( )( ) [ ( ) ( ) ( )] ( )sE W

E R E W D R d R A ss

, (3.1.50)

which then also can be written as,

0 1 1 22 1 2 2 2 1 1 2 2

0 1 1

( ) ( )( ) [ ( ) ( ) ( ) ( )] ( ) .sE W D W

E R E W D R d R A s A ss s

(3.1.51)

This expression indicates that the received field at R, which is generated at the second wedge

2W , can be written as a compound of two different waves, an amplitude wave and a slope

wave.

For computing the distance factors L and sL , the points on the shadow boundaries for each

successive edge are determined so that their distances are the same as the distances between

the wedges, as shown in Figure 3.1.5. Therefore, according to the method, one can write the

following expression for the field at point 1,1P which is diffracted at 1W ,

1 1,1 0 1 1 1,1 1 1( ) ( ) ( ) ( ) .E P E W D P A s (3.1.52)

At this point a known fact from knife-edge diffraction theory, stating that the diffracted field

on the shadow boundary is half of the incident field, can be employed. That is to say, if there

39

was no diffraction, the incident field on 1,1P would be twice that of the current situation. Thus,

the continuity equation is,

0 1,1 0 1 1 1,1 1 1( ) 2 ( ) ( ) ( )E P E W D P A s , (3.1.53)

and for Andersen introduces,

( ) / 2 ,D L (3.1.54)

so the equation becomes,

0 1,1 0 1 1,1 1 1( ) ( ) ( ) .E P E W L A s (3.1.55)

The amplitudes of 0 1( )E W and 0 1,1( )E P

are inversely proportional to the distances from the

transmitters which are consecutively 0s and 0 1( )s s . Therefore,

0 1 0 1( ) 01,1

0 1 0 1 0 1

1 1 ,

( )

jk s s jks jksse e L e

s s s s s s

(3.1.56)

with the distance factor,

0 11,1

0 1

s sL

s s

, (3.1.57)

and similarly for the next point,

0 1 21,2

0 1 2

( ) .

s s sL

s s s

(3.1.58)

But for the point 2,1P there should be two distance factors, since the field before diffraction

from 2W is an already diffracted wave. Therefore 2,1L for the amplitude wave and 2,1

sL for the

slope of the diffracted wave are to be calculated. 2,1L is calculated in the same way from the

continuity equation,

1 2,1 2 2,1 1 2 2 2,1 2 2( ) 2 ( ) 2 ( ) ( ) ( )E P E P E W D P A s , (3.1.59)

which becomes,

1 2,1 1 2 2,1 2 2( ) ( ) ( ) .E P E W L A s

(3.1.60)

Since 1 2,1( )E P and 1 2( )E W are already known, and by using the previously determined 1,1L

and 1,2L , 2,1L can then be calculated. For 2,1

sL one needs the continuity equation for the slope

of the field incident upon 2W ,

40

2

1 2,1 2 2,11 22 2

( ) ( )( ) 2 ( ) ,

E P d PE WA s

n n n

(3.1.61)

indicating that the diffracted field’s slope for point 2,1P without the existence of 2W is twice

the slope-diffracted field’s slope for point 2,1P when diffracted from 2W . In terms of the

distance factor,

1 2,1 3/21 2

2,1 2 2

( ) ( ) ( ) .sE P E W

L A sn n

(3.1.62)

The same scenario can be extended for three wedges as in Figure 3.1.6. The previous scenario

is modified by replacing the receiver R with another wedge 3W with the same height. And

the new receiver 'R is placed again in the shadow of the closest edge.

Figure 3.1.6 – Slope UTD geometry for 3 wedges.

In this case the received field at the receiver 'R can be written as 3( ')E R . Since 2( )E R is the

sum of two waves created at 2W and computed at R in the previous scenario, it is an already

diffracted field when it impinges on 3W . It will create a field strength 3( ')E R at the point 'R ,

which is now composed of 3 2( ') ( )D R E R and

23

( ) ( ')sE Rd R

n

, which requires the data of

normal and slope diffraction coefficients of the third wedge computed at 'R . Hence, the new

received field at receiver will be,

23 3 2 3

( )( ') ( ') ( ) ( ') .sE R

E R D R E R d Rn

(3.1.63)

From the first term of the expression one has one amplitude wave and one slope wave and

from the second part of the equation one has two slope waves. The second part of the

expression requires the second derivative of the diffraction coefficient ( )sd

n

which can be

found in related literature. This scenario can be extended to many wedges since the computing

process is a recursive algorithm. If there are reflections to the receiver they should be

considered as extra, additional rays.

2 0s

1s 1

n

0E 1E

2 ' 1 ' 2s

T

3, 3W E 3

2E

3s 2W

1W 'R

41

A comparison of these models is implemented by Tzaras, detailed in [7]. Figure 3.1.7 is a

comparison of the methods and includes the deterministic Vogler’s method. Implementing the

Vogler method is a huge undertaking, therefore in this thesis the data of Tzaras are used.

Tzaras tested the methods in various parts of the UK for frequencies between 40 MHz and

900 MHz. The results in Figure 3.1.7 show a cumulative distribution function of the

diffraction loss over the paths studied on one of the sites, containing 1390 records. The y-axis

gives the probability that the diffraction loss is less than the corresponding value on the x-

axis. For example, approximately 30% of the measured paths have a diffraction loss of less

than 10 dB.

Figure 3.1.7 – Comparison of described models.

3.1.4 The Integral Equation approach

The previous methods reduce the actual terrain profile to canonical objects. For accuracy

Vogler’s method can be used but it gets too complex as the number of canonical objects

increases. Hata presents a simple parametric model, but his model lacks dynamic variations of

the average signal due to shadowing and reflection.

Without reducing the actual surface geometry, Hufford derived an integral equation by

employing the Green’s theorem. But his integral equation method starts to experience

numerical instabilities over 100 MHz and in the UHF band and above these instabilities result

in significant loss of accuracy. By evaluating Hufford’s approach, Hviid presents another

method which can overcome the instabilities at high frequencies. Hviid derives his method by

using some important assumptions. The surface is assumed to be a perfect magnetic conductor

with a reflection coefficient of -1. To avoid problems in integration the surface has to be

42

smooth which means that the method excludes steep slopes and cliffs. Since the terrain profile

is two-dimensional side scattering effects are also neglected. Back scattering which happens

beyond and before the surface is also neglected because the method is developed to predict

the mean average field. Back scattering is negligible for slow fading and is only taken into

account for fast fading which is not an issue of interest for this method [8, 17, 18].

Since the surface is assumed to be a perfect magnetic conductor (PMC) the boundary

conditions are,

ˆ 0 n E and ˆ 0 . n H (3.1.64)

The free space spherical Green’s function is defined as

( , ') ,jkRe

GR

r r

(3.1.65)

which satisfies the hom*ogenous scalar Helmholtz equation (see Appendix A),

2 2( ) 0k G , (3.1.66)

with a source point at ' r r . R r r' is the distance from the source point to the

observation point as shown in Figure 3.1.8. By using the second Green’s theorem an integral

equation is obtained,

ˆ ˆ ˆ ˆ( ) ( ) ( ' ) ' '4

iS

TT Gds

n E r n E r n n E , (3.1.67)

where 2T when r is on the surface S and 1T otherwise.

The electric field at the point r can be written as a sum of the incident field ( )iE r and the

total scattered field ( )scatE r ,

ˆ ˆ ˆ( ) ( ) ( ) .i scatT T n E r n E r n E r

(3.1.68)

In terms of magnetic currents,

ˆ M n E , (3.1.69)

the integral equation becomes,

( ) ( ) ( ) ,i scatT T M r M r M r

(3.1.70)

ˆ( ) ( ) ' ' ,4

iS

TT Gds

M r M r n M

(3.1.71)

2

1ˆ ˆˆ ˆ( ) ( ) (( ) ( ) ) ',4

jkR

iS

T jkRT e ds

R

M r M r n R M n M R

(3.1.72)

43

where ˆ = 0n M . Hviid simplifies this expression by assuming the surface to be constant

transverse to the direction of propagation. The integration is reduced to x-z pane for a

vertically polarized source at 0r . The function z = c(x, y) = c(x) can be used to describe the

2-D surface. Figure 3.1.8 shows the geometry for 2-D scattering occurring at a point 'r where the observation point is at r and the transmitter is at

0r . The scattering points are

treated as source points for the observation point.

Figure 3.1.8 - Source point 0r , scattering point 'r , and observation point r in x-z plane.

For each scattering source point r' , there are other contributions ( ' 0y ). This is shown in

Figure 3.1.9.

Figure 3.1.9 – Scattering geometry for a scattering point on 'y -axis.

z

n̂

r 'r

0r x

ˆ 'n

1R

R

E

z

y

x

'x

'y

1( ',0)xR

1( ', ')x yR

( ', ')x yR

(x',0)R

'r

0r

r

E

n̂

44

Taking the x-z plane as reference, one may find relations between the path lengths ( ', ')R x y

and corresponding phases ( ', ')x y along 'y .

The path difference for a scattering point on 'y -axis with respect to ( ',0)R x then can be

written as,

( ', ') ( ', ') ( ',0) .R x y R x y R x

(3.1.73)

The phase difference becomes,

( ', ') ( ', ')x y k R x y

(3.1.74)

From Figure 3.1.9 the geometry yields,

22 2

2

'( ', ', ') ( ',0) ' ( ',0) 1

( ',0)

yR x y z R x y R x

R x

(3.1.75)

Assuming ' ( ',0)y R x and r' r one obtains,

2'( ', ') .

2 ( ',0)

yR x y

R x

(3.1.76)

Similarly for the path length difference between the source and the scattering point,

2

1

1

'( ', ')

2 ( ',0)

yR x y

R x 1' ( ',0) .y R x

(3.1.77)

For a point ( ', ')x y , the total difference in path length with respect to the point ( ',0)x then

becomes,

211

1

( ',0) ( ',0)' .

2 ( ',0) ( ',0)tot

R x R xR R R y

R x R x

(3.1.78)

For a point on the path, one can compute the scattering integral along the 'y -axis, considering

only the phase variation. This is because the amplitude variation is negligible compared to the

phase variation. The reduced integral is given by [18],

( ( ') )

2

1 1 ( ', ')ˆˆ( ) ' ( ') ( ', ') '.4 ( ', ')

totjk R y R

scat

jkR x yM x dl M x x y e dy

R x y

n R

In equation (3.1.79), scatM gives the magnetic current induced by the scattering surface at the

observation point. The integral with respect to 'y corresponds to the total contribution from a

point 'x . One can note that ˆˆ ( ', ') 0x y n R when the observation and scattering points are

on a planar surface and no scattering occurs (specular reflection). The equation

(3.1.79)

45

ˆˆ ( ', ') 0x y n R also holds for r' r on smooth surface. Therefore, Hviid approximates the

term 2

1 ( ', ')

( ', ')

jkR x y

R x y

by

jk

R and hence neglects the near field term assuming constant

amplitude in the 'y direction. The integral then is reduced to,

21

1 2

'21 ˆˆ( ) ' ( ') ( ', ') ' .

4

k R Rj y

R RjkR

scat

jkM x dl M x x y e e dy

R

n R

(3.1.80)

The integral with respect to 'y is a Fresnel integral and can be reduced to,

1 4

1

1 ˆˆ( ) ( ') ( ', ') ' .4

jjkR

scat

jk R RM x M x x y e e dl

R R R

n R

(3.1.81)

Equation (3.1.81) is the integral along the path in x-z plane. For the total field at the

observation point one obtains the integral equation,

1 4

1

ˆˆ( ) ( ) ( ') ' ,4

jjkR

i

T jk R RM x TM x M x e e dl

R R R

n R

(3.1.82)

where,

ˆ( ) ( ),i ix x M n E

(3.1.83)

is the magnetic current induced by the incident field at the observation point.

Numerically,

this integral is replaced by a finite summation over surface elements along the x-axis. The

increments mx in the summation are smaller than 2

. The result is a linear system for nM

at a point mx n x from the source point,

1

,

( , ) 4

n

n i n m m

m

TM TM M f n m x

, (3.1.84)

where 0, 1m N and 0, 1n N and,

( )1 4

1

ˆˆ( , )j kR

m

m

jk R R lf n m e

R R R x

n R , (3.1.85)

where ml is a line element along the surface and the distances are,

2 2

1 ( ) ,m m TR x z h

(3.1.86)

2 2( ) ( ) ,N m R mR x x h z

(3.1.87)

46

where Th and Rh represent the heights of the source point and the observation point,

respectively. Th , Rh and Nx (the end point) are the fixed parameters, and mx and mz are

numerical parameters for finding nM at an observation point mx N x . For each increment,

the distances 1R and R are updated and employed in equation (3.1.85) to obtain ( , )f n m . The

parameters are defined using Figure 3.1.10.

Figure 3.1.10 – Geometry for the 2-D integral equation method.

Smoothness of the surface is important for the approximations. Since the integral equation in

Hviid’s approach has no singularity for smooth surfaces, suitable assumptions yield a simple

numerical solution without divergences which makes the method a good alternative to UTD.

On the other hand, since there is a linear system, the conditioning of the system might

deteriorate if the system is large.

More recent work by Carin and Li, which aims at validating ray tracing models for complex

environments, uses the Multilevel Fast Multipole algorithm based on integral equations [19].

Especially in large and complex environments the algorithm seems to be more accurate and

applicable than Hviid’s method.

Th Rh

mz

mz

1Nz

mx Nx 1Nx

mN x

47

3.1.5 The Parabolic Equation method

The Parabolic Equation Method is used for long range forward propagation problems and

often for irregular terrains with buildings sprawling over them. The method takes reflection,

diffraction and forward scattering into consideration. The PE method was developed for

problems where the propagation has a paraxial direction as shown in Figure 3.1.11. By

calculating the dominant initial field distribution 0( , )u x z , the next distribution ( , )nu x x z

is found by propagating it in x. The distribution is computed by marching forward step by step

( )x in Cartesian coordinates. Each step in the calculation depends on the previous step [8,

10]. Assuming the time dependence ,j te Helmholtz equation for a field component is

given by,

2 2 2

0 0k n , (3.1.88)

where,

kn

k , (3.1.89)

is the refractive index and,

2k

, (3.1.90)

is the wave-number. 0k is the free space wave-number.

Since the earth is spherical, a spherical coordinate system would seem appropriate for

propagation in the troposphere, but one can simplify the problem by considering the earth as

flat and use cylindrical coordinates ( , , )r z . Another assumption for this case is that the field

is constant in the azimuth. One can then write the Helmholtz wave equation in cylindrical

coordinates in terms of range r and height z as,

2 22 2

02 2

10 k n

r r r z

(3.1.91)

Figure 3.1.11 illustrates paraxial propagation with a propagation angle in cylindrical

coordinates.

48

Figure 3.1.11- Paraxial propagation with a narrow beam-width in cylindrical coordinates.

To find a solution in the paraxial direction, one may write the field component in terms of an

envelope function ( , )u r z which is assumed to be slowly varying in ( , )r z , and a Hankel

function (2)

0 0( )H k r which represents the rapid oscillations in range,

(2)

0 0( , ) ( , ) ( )r z u r z H k r . (3.1.92)

In the far field 0 1k r and the asymptotic expression of the Hankel function can be used,

0( )(2) 40 0

2( )

j k r

H k r ek r

. (3.1.93)

Substituting the modified uH into the Helmholtz equation, neglecting small terms,

2 22 2

02 2

( )0

uH uH k n uH

r z

, (3.1.94)

and taking the second derivatives yields,

2 2

( )0r r

zz

u H uHu H uHk n

r

, (3.1.95)

2 2

0( ) 0rr r r r r rr zzu H u H u H uH u H uHk n , (3.1.96)

2 2

0( ) 2 ( ) 0rr zz r r rru u H u H u H Hk n . (3.1.97)

Neglecting the small terms yields,

2 22 2

0 02 22 ( 1) 0

u u ujk k n u

r z r

. (3.1.98)

For the 2-D problem, one can replace r with x in equation (3.1.98),

z

x r

49

2 22 2

0 02 22 ( 1) 0 .

u u ujk k n u

x x z

(3.1.99)

This equation describes the propagation in both the +x and –x directions. By factorizing the

equation (3.1.99),

(3.1.100)

2 22 2

0 0 0 02 2 2 2

0 0

1 11 ( 1) 1 ( 1) 0,jk jk n jk jk n u

x k z x k z

and discarding the backward propagation factor, one has the expression for the forward

propagation,

22

0 0 2 2

11 ( 1) 0.jk jk n u

x k z

(3.1.101)

Defining a differential operator,

22

2 2

1( , ) ( 1)P x z n

k z

, (3.1.102)

one can write,

0 1 1 ( , )u

jk P x z ux

, (3.1.103)

and by using Taylor expansion,

( , )1 ( , ) 1

2

P x zP x z | | 1P . (3.1.104)

The narrow angle ( 15 15 ) parabolic equation is then obtained as,

22 2

02

1( 1) .

2

u uk n u

x jk z

(3.1.105)

This partial differential equation can be solved by the split-step Fourier Transform method.

By defining,

2( ( , ) ) ,z

ua x z bD u

x

(3.1.106)

where,

50

20( , ) ( ( , ) 1)2

ka x z n x z

j and

1 ,

2b

jk

(3.1.107)

one transforms ( , )u x z to ( , )zU x k by using the Fourier transform and its derivatives,

( )z zD u jk U . (3.1.108)

Here the domain of the transform, zk , is the transverse wave-number which can be written,

sinzk k , (3.1.109)

where is the propagation angle relative to the horizontal. The equation then becomes,

2

z

Ua k b U

x

, (3.1.110)

where a is constant if a small variation is assumed for the refractive index. The solution to

this ordinary differential equation yields,

2( )( , ) ( , ) za k b x

z zU x x k U x k e , (3.1.111)

and by taking the inverse Fourier transform, the solution for ( , )u x x z in terms of Fourier

and inverse Fourier transforms is,

2

( , ) ( , )zbk xa xu x x z e e u x z (3.1.112)

The solution can be marched forward numerically along the x-axis, for a given value of z.

Figure 3.1.12 - Solution domain of the Parabolic Equation method.

51

3.2 PROPAGATION MODELS FOR URBAN AREAS

3.2.1 The Okumura model

By using vertical omni-directional antennas for both transmitter and mobile, Okumura

obtained extensive measurement data for median attenuation relative to free space muA for

different distances and frequencies. Okumura made measurements in Tokyo city and therefore

his model was basically applicable, and most accurate, for dense urban areas with quasi-

smooth terrain. With less accuracy the model can be employed also for other environments

such as suburban and rural areas by adding an area correction factor AREAG . The model is

applicable for distances from 1 km up to100 km , a frequency range from 150 MHz up to

1920 MHz , and base station mast heights from 30 mup to1000 m [4].

Based on the measurements, Okumura came up with curves for median attenuation which are

functions of distance and frequency [20]. Total mean path loss due to the model is described

in the following formula:

50 ( , ) ( ) ( )f mu te re AREAL L A f d G h G h G , (3.2.1)

where fL is the free space loss, muA is the median attenuation relative to free space, ( )teG h is

the base station antenna height gain factor, ( )reG h is the mobile antenna height gain factor,

and AREAG is the correction factor gain due to environment. ( , )muA f d

and AREAG are found

from the plots in Figure 3.2.1. The antenna height gain factors are given:

10

10

10

( ) 20log 30 m 1000 m ,200

( ) 20log 3 m ,200

( ) 20log 3 m 10 m .200

tete te

tere re

tere re

hG h h

hG h h

hG h h

(3.2.2)

52

Figure 3.2.1 – Okumura curves.

As an example the dashed red indicators in Figure 3.2.1 show how to read muA and AREAG

values from the Okumura curves at 1000 MHz frequency, 80 km away from the base station

in a suburban area.

53

3.2.2 The Hata model

The Hata model is one of the most common models in designing real systems. Many new

models are still using it as a reference model because of its simplicity and accuracy. The Hata

model uses four parameters for estimating the path loss: Carrier frequency cf in MHz ,

distance R in km, base station antenna height bh in m, and mobile antenna height mh in m.

The model is only valid for the parameter ranges 150 MHz 1500 MHzcf , 1 kmR ,

30 200bh and 1 10mh [4].

Hata empirically approximated Okumura’s graphical path loss data into formulas in the

following forms [21]:

10

10

10

Urban Areas log ( ) ,

Suburban Areas log ( ) ,

Open Areas log ( ) .

dB

dB

dB

L A B R

L A B R C

L A B R D

where,

10 10

10

2

10

2

10 10

2

10

69.55 26.16log ( ) 13.82log ( ) ,

44.9 6.55log ( ) ,

2(log ( / 28)) 5.4 ,

4.78(log ( )) 18.33log ( ) 40.94 ,

3.2(log (11.75 )) 4.97 large citi

c b

b

c

c c

m

A f h

B h

C f

D f f

h

2

10

10 10

es, f 300 MHz ,

8.29(log (1.54 )) 1.1 large cities, f 300 MHz ,

(1.1log ( ) 0.7) (1.56log ( ) 0.8) small to medium cities .

c

m c

c m c

h

f h f

In general, curves or formulas that are fit to a set of measurements are not valid outside the

range of the measurements. Since the model is less flexible comparing to deterministic

models and invalid for the range less than 1 km.

(3.2.3)

(3.2.4)

54

3.2.3 The Walfisch - Bertoni model

In the case of dense urban areas the model approximates rows of buildings as absorbing

screens with uniform height and separation distances. Being semi-deterministic the model has

a formulation for the path-loss that includes more parameters than the Hata model. The model

reduces the path loss to three components, free space loss, loss along the buildings and loss

down at the street level. For calculating the loss along buildings the model starts by using the

Repeated Kirchhoff Integral for uniform parallel screens which is an application of scalar

diffraction theory with some approximations [1]. For a derivation of this integral see

Appendix A. In the case when there are uniform screens with uniform heights ( 0)nh , the

field at the plane of screen 1ny can be written in terms of the incident wave coming towards

the first screen. This can be done by expressing the field at each consecutive screen in terms

of the field at the previous one. For the second screen the field can be written in terms of the

field above the first screen as [1],

22 1( )4

22 1 1

(2 , ) ( , ) ,

j jk y yjkd

de e

E d y E d y e dyd

(3.2.5)

and for the third screen one has,

23 2( )4

23 2 2

(3 , ) (2 , ) ,

j jk y yjkd

de e

E d y E d y e dyd

(3.2.6)

and so forth and so on until the last screen 1Ny .

(3.2.7)

2

1

1

42

1 1 2 3 1 1/2

0 0 0 0

(( 1) , ) ( , ) .

N

n n

n

jN kjkNd j y y

d

N NN

e eE N d y dy dy dy dy E d y e

d

With a parameter,

,2

n n

jkv y

d (3.2.8)

and the identity,

2 ,j

j e

(3.2.9)

substituting equation (3.2.9) into equation (3.2.8) yields,

55

4 ,j

n n

ddy e dv

(3.2.10)

and the result can be written as [1], (3.2.11)

21

1

1 1 2 3 1/2

0 0 0 0

(( 1) , ) ( , ) .

N

n n

n

jkNd v v

N NN

eE N d y dv dv dv dv E d y e

The incident plane wave in 2D coordinates can be written as,

( cos sin )

1 0 0( , ) .j jk d yE d y E e E e k r

(3.2.12)

For small angles and unit amplitude 1 sin

1( , ) jkyjkdE d y e e , and by Taylor series expansion,

1 1

1( , ) ( sin ) .

!

jkd m

m

E d y e jkym

(3.2.13)

Introducing a dimensionless parameter,

sin

dg

where 1 1

2

jkv y

d , (3.2.14)

yields,

1 1

1( , ) (2 ) .

!

jkd m

m

E d y e g j vm

(3.2.15)

At the screen 1, 1( )N Nx y , for the point where 1 0Ny ( 1 0Nv ) and 1 ( 1)Nx N d , the

expression for 1( ,0)NE x becomes,

(3.2.16)1

2 21 1

1 2

( 1) 2 ( ) 2

1 2 3 1/200 0 0 0

1(( 1) ,0) (2 ) .

!

N N

n n n

n n

jk N d v v v vm

NNm

eE N d dv dv dv dv g j v e

m

This result can be written in terms of the Boersma function which is in general defined as [1],

12 2

1 1

1 2

2 ( ) 2

, 1 2 3 1/2

0 0 0 0

1( ) .

N N

n n n

n n

v v v vm

N m NNI dv dv dv dv v e

(3.2.17)

56

Thus for 1 ,

( 1)

,

1(( 1) ,0) (2 ) (1) .

!

jk N d m

N m

m

E N d e g j Im

(3.2.18)

Boersma evaluates , ( )N mI using a recursion relation for 2m as follows:

1, 1

, , 21 11

( )( 1) 1( ) ( ) .

2( 1) 2 ( 1)

Nn m

N m N m

n

IN mI I

N N N n

(3.2.19)

For 1 the initial values are given as:

(3.2.20)

0,

1 for m 0(1)

0 for m 0mI

; ,0

(1 2)(1)

!

NNI

N ;

1

,1

1 (1 2)(1) ,

2 !

Nn

N

n

In N n

where (1 2)n is the Pochhammer symbol which is defined as,

0 1( ) 1 ; ( ) ; ( ) ( 1)...( 1) .np p p p p p p n

(3.2.21)

For propagation parallel to the rooftops, with the angle of incidence 0 , the parameter g

vanishes due to equation (3.2.14) and,

( 1)

,0(( 1) ,0) (1) .jk N d

NE N d e I (3.2.22)

Using the asymptotic approximation for ,0

(1 2)(1)

!

NNI

N one can finally obtain the following

result:

( 1) ( 1)

1

(1 2) 1( ,0) .

! 1

jk N d jk N dNNE x e e

N N

(3.2.23)

This expression implies that the amplitude of the field decreases monotonically as 1N

.

57

Figure 3.2.2 – Field settling to a constant value after initial variation.

Figure 3.2.2, created in MATLAB, relates to an urban scenario where 500 consecutive

buildings along a street are approximated with uniform parallel screens with equal spacing.

As parameters, the spacing is 250d and the maximum number of m (the upper limit in

equation (3.2.18)) is max 200m . The base station or the transmitter is assumed to be

mounted above the rooftop level before the first screen. Since the polarized incident wave is a

plane wave with a small angle of incidence one can obtain the curves in Figure 3.2.2 for

various angles. The curves show the change in electric field magnitude 1| ( ,0) |NE x as a

function of N , the number of screens. The red curve for the angle 0 drops one half just

after the first edge and continues to decrease almost monotonically as predicted before in

equation (3.2.23). As the angle is increased to 0,5 and 1 degrees the fields decrease less

dramatically and settle as shown in the magenta and blue curves. The dashed blue curve

shows the behavior of the field when the frequency is doubled for the 1 case. Moreover, in

the case of larger angles of incidence, 2 and 2,5 , the electric field even starts to increase

just after the first screen and settles at a value larger than a unity as shown in the green and

black curves. In these cases the incident plane wave is reinforced at every screen top by the

diffracted waves from the preceding screens.

In Figure 3.2.2, the black points indicate the settling points where the oscillation with

decreasing amplitude starts to settle about a certain value as N increases. In Figure 3.2.3, the

screen FN N causes the first penetration into the first Fresnel zone. After this screen all the

screens penetrate the first Fresnel zone until the last screen. These are the screens that cause

the settling of the field about a certain value in the plots. The number of screens penetrating

the first Fresnel zone is 0N . One can find an equation for 0N in terms of the known data d ,

58

and . Figure 3.2.3 depicts this scenario where the first red line segment from left is the

vertical distance from the blue screen’s top to the axis of the first Fresnel zone.

Figure 3.2.3 – First Fresnel zone along the screens.

From the blue screen FN each screen has a vertical distance from the top of the screen to the

axis of Fresnel zone that is smaller than the Fw of the ellipsoid for that particular point on the

axis. Similarly the green line segment is longer than Fw . From the geometry the vertical

distance for each screen can be written as the product tans . That is to say that screens

satisfying the condition tanFw s do not penetrate the Fresnel zone, while screens

satisfying tanFw s do. To find 0N

the condition tanFw s

should be applied to the

critical screen depicted in blue which first penetrates the Fresnel zone. In general, for the thn

Fresnel zone:

.F

sDw n

s D

(3.2.24)

For the first Fresnel zone one can write,

,

1F

sw

s

D

(3.2.25)

and assuming that the plane wave source is at infinity ( D ), the approximate result is,

.Fw s

(3.2.26)

Secondly, another approximation for the blue screen case can be made as 0s N d since the

angle of incidence is small. Applying the condition tanFw s with the approximation one

obtains,

D

s

Fw

0N d

1N 1N

1n

FN

59

0 0 tanN d N d , (3.2.27)

and therefore, for small angles,

2

0 2

1.

(tan )N

d g

(3.2.28)

Since 0N basically defines the settling value of the field, and since

0N for small angles

can

be written in terms of g , one can define the magnitude of the settled value of the field in

terms of g . In Figure 3.2.4, which is plotted in MATLAB, the black plot in the left subplot

depicts the dependency of | |settlingE with respect to g in logarithmic scales and the blue line is

a polynomial fitting curve given by,

2 3( ) 3.502 3.327 0.962P g g g g for 1 ,g (3.2.29)

which can be approximated as,

0.9

( ) 0.10.03

gP g

for 0.015 0.4 .g (3.2.30)

The right subplot shows this polynomial fitting curve and its straight line approximation in

red dashed line.

Figure 3.2.4 – Settling due to the g parameter, P(g) , and a line of approximation.

( )P g can be used to define the path loss due to the propagation over all the screens. The path

loss along the rooftops then can be written in terms of the field reduction as,

60

2

0.9

2

sin

( ) 0.1 .0.03

rooftops

d

PL P g

(3.2.31)

Here sin can be written in terms of the base station or transmitter height Th , the building

height BH , and the distance R as,

sin .T Bh H

R

(3.2.32)

Therefore the equation becomes,

1.8 0.9

0.01 .0.03

T Brooftops

h H dPL

R

(3.2.33)

To find the reduction of the field power diffracted from the last screen (building) down to the

street level, Ikegami’s approach can be used. It states that the spatial average power over a

distance is approximately the sum of the individual ray powers. Another fact is that the plane

waves incident on the building edges act as a source of cylindrical waves for the mobile level.

One can then come up with an estimate by using that the power density of cylindrical waves is

inversely proportional to the path distance i from the building edge to the mobile, and

proportional to the square of the diffraction coefficient for an absorbing wedge 2| ( ) |iD . The

diffraction coefficients differ with boundary conditions but for small angles they approach the

same dependence. For simplicity, using the Felsen coefficients,

2

2 1 1 1| ( ) | ,

2 | | 2 | |i

i i

Dk

(3.2.34)

the reduction from the rooftop down the street level can be written as follows:

2 2

2

2 1 1 1 2 2

1 1 1 | | 1 1 .

2 | | 2 | | 2 | | 2 | |downPL

k k

(3.2.35)

Considering the Figure 3.2.5, with the assumptions 2 13 and 2| | 0.1 , the rays have

nearly equal amplitudes [1]. Thus,

2

1

2 2

1 1 1 1 1

1 1 1 1

| | 2 | | 2 ( )down

b m

PLk k H h

. (3.2.36)

2

61

Figure 3.2.5 – Geometry for finding the loss down at the street level.

Finally the total path loss can be written as, (3.2.37)

2 1.8 0.9 2.12 1 1

2 4 2 3.8

5.51 ( )( ) .

4 2 ( ) 32 ( )

T Btotal

b m B m

h H dPL P g

R H h H h R

This equation can be expressed in decibels as follows: (3.2.38)

0.9

1

289.5 10log 21log 18log( ) 38log ,

( )total M T B k

B m

dPL f h H R

H h

where Mf is the frequency in megahertz and kR is the distance in kilometers.

L M

1

62

4. SIMULATION

4.1 The Hata model

Figure 4.1.1 shows the model for different frequencies in various environments:

Figure 4.1.1 – The path loss for different frequencies in different environments.

Figure 4.1.2 shows the model for different transmitter heights in various environments:

Figure 4.1.2 – The path loss for different transmitter heights in various environments.

63

4.2 The Walfisch - Bertoni model

Based on the theory the following figures show the predicted path loss graphs as a function of

distance when the propagation occurs along buildings with uniform height and spacing.

Simulations were plotted for various frequencies f, transmitter heightstH , building separation

distances dH and building heights bH .

Figure 4.2.1 indicates that the path loss increases (more shadowing) with increasing

frequency:

Figure 4.2.1 – Path loss due to distance and frequency.

Figure 4.2.2 indicates that the path loss decreases (less shadowing) by increasing the

transmitter heights:

Figure 4.2.2 – Path loss due to distance and transmitter height.

64

Figure 4.2.3 indicates that the path loss is smaller (less shadowing) for larger building

separation distances:

Figure 4.2.3 – Path loss due to distance and building separation distance.

Figure 4.2.4 indicates that the path loss is larger (more shadowing) for higher buildings:

Figure 4.2.4 – Path loss due to distance and building heights.

65

4.3 Comparison of the Hata model and the Walfisch - Bertoni model

A detailed comparison of the proposed models was obtained for four cases where in each case

four parameters are fixed and one particular parameter has two values. The distance range for

each case was taken as 1 km 10 kmR since the Hata model is valid in the range of

1 km 10 kmHATAR . The building heights bH were chosen to be smaller than the

transmitter heights tH since this is a must for the Walfisch-Bertoni model to be valid.

In the first case mobile antenna height, transmitter height, building heights and building

separation distances were and the path loss for the Hata and the Walfisch-Bertoni models

were plotted for two different carrier frequencies, 800 MHz and 1500 MHz . Figure 4.3.1

shows the results with the fixed parameters:

Mobile antenna height 1.5 mmh

Transmitter height 30 mtH

Building heights 10 mbH

Building separation distances 50 md

Figure 4.3.1 – Path loss from the Hata model and the Walfisch-Bertoni model for two

frequencies.

66

In the second case frequency, mobile antenna height, building heights and building separation

distances were fixed and the path loss for the Hata and the Walfisch-Bertoni models were

plotted for two different transmitter heights, 30 m and.60 m. Figure 4.3.2 shows the results

with the fixed parameters:

Carrier Frequency 900 MHzf

Mobile antenna height 1.5 mmh

Building heights 10 mbH

Building separation distances 50 md

Figure 4.3.2 – Path loss from the Hata model and the Walfisch-Bertoni model for two heights.

In the third case frequency, mobile antenna height, transmitter height and building heights

were fixed and the path loss for the Hata and the Walfisch-Bertoni models were plotted for

two different building separation distances, 20 m and 50 m . Figure 4.3.3 shows the results

with the fixed parameters:

Carrier Frequency 900 MHzf

Transmitter height 30 mtH

Mobile antenna height 1.5 mmh

Building heights 10 mbH

67

Figure 4.3.3 – Path loss from the Hata model and the Walfisch-Bertoni model for two

building separations.

In the fourth case frequency, mobile antenna height, transmitter height and building

separation distances were and the path loss for the Hata and the Walfisch-Bertoni models

were plotted for two different building heights, 10 m and.15 m. Figure 4.3.4 shows the results

with0 the fixed parameters:

Carrier Frequency 900 MHzf

Transmitter height 30 mtH

Mobile antenna height 1.5 mmh

Building separation distances 50 md

Figure 4.3.4 – Path loss from the Hata model and the Walfisch-Bertoni model for two

building heights.

68

APPENDIX A

For diffraction problems, and especially when the wavelength is significantly shorter than the

linear dimension of the obstacle, the influence of polarization is not as important as in the case

of reflection. Hence, one may treat the diffraction by using scalar wave equations which

simplifies derivations by neglecting the vector nature of waves. This Appendix, starting from

the Green’s Theorem, the repeated Kirchhoff integral for parallel screens is obtained in a step

by step manner. In the case of uniform parallel screens, this integral can be seen as the starting

point of the Walfisch-Bertoni model.

Derivation of Green’s Theorem

For a vector field F the divergence theorem states that,

V

ˆ ' ' ' .S

dv ds F F n

(A.1)

The vector field F can be written in terms of two scalar fields V and U as,

.V U F (A.2)

Using the identity

( )f f f A A A , (A.3)

it can be written,

2 ,V U V U F (A.4)

ˆ ˆ .'

UV U V

n

F n' n (A.5)

Substituting equations (A.4) and (A.5) into equation (A.1) yields,

2

V V

( ) ' ' ''

S

UV U dv V U dv V ds

n

. (A.6)

Interchanging U and V and applying the same steps yields,

2

V V

( ) ' ' ''

S

VU V dv U V dv U ds

n

. (A.7)

Subtracting equation (A.6) from equation (A.7) produces Green’s theorem:

2 2

V

( ) ' ( ) '' '

S

U VV U U V dv V U ds

n n

. (A.8)

69

Derivation of the scalar Helmholtz Equation

From Maxwell’s equations the wave equation for the electric field is,

22

0 0 02

( )t t

E JE . (A.9)

In general,

22

2 2

1( ) ( , ) ( , )t tc t

E r p r , (A.10)

where ( , )tp r is a prescribed source. This equation is obeyed by each scalar Cartesian

component of the electric field, , x y zE E and E . In the case of ( , , )x xE E x y z for a charge

and current free region in free space the wave equation becomes,

22

2 2

1( ) 0xec t

. (A.11)

where ( , , , )x xe e x y z t is the real part of the time-dependent scalar field and

can be written in

terms of spatial and time-harmonic components as,

Re{ }j t

x xe E e . (A.12)

Substituting j t

xE e in the wave equation where,

22 2

2

j t j t

x xE e j E et

, (A.13)

the equation becomes,

22

20 .j t j t

x xE e E ec

(A.14)

Dividing by the factor j te , and with the wave-number kc

, the Helmholtz equation for

the scalar time-independent component of the electric field is obtained as,

2 2( ) 0xk E . (A.15)

70

The Kirchhoff Integral Theorem

Helmholtz equation for a scalar function with j te dependency and spatial component U can

be written as

2 2( ) 0k U . (A.16)

Introducing a Green’s function ( , ')G r r which satisfies the Helmholtz equation with a source

at 'r ,

2 2( ) ( , ') 4 ( ')k G r r r r , (A.17)

and substituting equations (A.16) and (A.17) into Green’s theorem in equation (A.8) yields,

2 2

V

( ) ' ( ) '' '

S

U GG U U G dv G U ds

n n

, (A.18)

2 2

V

+4 ( ') ' ( ) '' '

S

U GUGk UGk U dv G U ds

n n

r r , (A.19)

1( ) ( ) '

4 ' 'S

U GU G U ds

n n

r . (A.20)

where r is the observation point inside V.

According to Green’s theorem the function G must be continuous inside the chosen volume

(see Figure A.1). An enclosed volume V (region II) can be chosen between two surfaces, S

and 1S . In Figure A.1, where r is the observation point, r' is the source point on the surface

of integration and ˆ 'n is the outward normal to the surfaces. The inner surface S and the outer

surface 1S can be chosen as arbitrary smooth surfaces. One assumes that the sources of

radiation are inside S [22]. Here S interacts with the fields generated in region I, allows the

modified field to pass through it, and cause a diffraction pattern in region II (volume V).

The surface of integration is the composite surface 1totS S S . From equation A.20 one has,

1( ) ( ) '

4 ' 'totS

U GU G U ds

n n

r , (A.21)

which can also be written in composite form,

71

1

1 1( ) ( ) ' ( ) '

4 ' ' 4 ' 'S S

U G U GU G U ds G U ds

n n n n

r . (A.22)

Since the fields from the surfaces are propagating outward, Green’s function has the form,

( , ')jkRe

GR

r r , (A.23)

where | ' |R r r .

Figure A.1 – Geometry for the Kirchhoff Integral Theorem.

In equation (A.22) one can assume 1S to be at infinity and the integration over 1S can be

written,

1

1( ) '

4 ' 'S

U GI G U ds

n n

. (A.24)

Considering,

( , ') 1ˆˆ '

jkRG ejk

n R R

r rn' R , (A.25)

r

'r

R

y

x

z

S

1S

I

II

n̂'

n̂'

R

r

'r

72

and substituting the surface element 2' 'ds r d in equation (A.24),

2

1

1 1ˆˆ ( U) '4 '

S

UI G jk r d

n R

n' R

.

(A.26)

Since 1S tends to infinity, considering Figure A.1, one may realize R and 'r become parallel

and 'R r . The integral can be written in the form,

1 1 1ˆ ˆˆ ˆ' ' ' ) ' ' ' 4 ' 4 '

S S

Ur G r jkU d r G r U d

n r

n R n R

(A.27)

In equation A.27, when 'r , the factors 'r G and 'r U

are bounded and 1

0'r .

Therefore the second surface integral vanishes. For the first surface integral, considering the

approximations at infinity ( ˆ ˆ ˆ' ' R r n ) and,

'

'

jkreU

r

,

(A.28)

the second surface integral also vanishes with the radiation condition,

'lim ' 0 .

'r

Ur jkU

r

(A.29)

Hence, the integration over 1S when R vanishes,

1

1lim ( ) ' 0 ,

4 ' 'S S

U GG U ds

n n

(A.30)

and equation (A.22) yields the Kirchhoff integral theorem as follows,

1( ) '

4 ' '

jkR jkR

S

U e eU U ds

n R n R

r . (A.31)

As a result, the Kirchhoff integral theorem relates the value of any scalar wave function U , at

any point r outside an arbitrary closed surface S , to the value of the wave function at the

surface. The theorem plays an important role in the scalar theory of diffraction and the

vectorial case is treated in a similar fashion.

73

Applying the Kirchhoff Integral Theorem to Diffraction

Figure A.2 – Geometry for the Kirchhoff integral theorem applied to diffraction.

For practical calculations the Kirchhoff integral can be approximated according to the

geometry of diffraction and the impinging wave. Figure A.2 depicts a scenario where

diffraction of a plane wave occurs at a slit in an opaque, plane, absorbing screen which does

not let the field pass through except at the slit. The arbitrary closed surface S is chosen as an

infinite half sphere to the right of the screen with observation point ( , , )x y zr . Since S is

infinite the contributions of the surface elements in Kirchhoff integral vanishes except for the

area of the slit. The integration then is limited to the slit area A . An approximation for the

Kirchhoff integral can be obtained by assuming short wavelengths. That is to say, 1

kR

lets us neglect the derivative of 1

R compared to the derivative of jkRe in equation A.31,

1( , , ) '

4 '

jkR

A

ex y z jk ds

n R

. (A.32)

At 0x one has a plane wave,

cosjkx

i e , (A.33)

where,

R

n̂

( , , )x y z

y

y

z

x

A

74

' 'n x

, (A.34)

so that,

cos and cos' ' ' '

i ix i i

R R xjk jk

n x n x R

(A.35)

Substituting these in equation (A.32) yields,

( , , ) cos cos (0, ', ') ' '4

jkR

i

A

jk ex y z y z dy dz

R

. (A.36)

In the extreme case of no screen, the slit A extends from the point (0,0,0) to infinity in the

yz plane. The equation then can be written as:

( , , ) cos cos (0, ', ') ' '4

jkR

i

jkex y z y z dy dz

R

. (A.37)

Approximation of plane wave diffraction by a half screen in 2D

Following the approximated equation (A.37), the approximate value for the field at a point

0x can be written in terms of the incident plane wave,

jkx

incE E e , (A.38)

on the aperture A at 0x ,

0(0, , )incE y z E . (A.39)

Since the incident wave propagates along the x-axis ( cos 1 ), the z component of the

diffracted field is,

0( , , ) 1 cos ' '4

jkR

z

jk eE x y z E dy dz

R

. (A.40)

For an absorbing half screen the contribution from 0y vanishes.

75

( , , ) ' 1 cos '4

jkR

z

jk eE x y z E dy dz

R

, (A.41)

where,

2 2 2( ') ( ')R x y y z

(A.42)

To reduce the problem to 2D geometry the distance R can be replaced by 2 ( ')x y y in

the denominator because the primary contribution from 'z to the integral lies in the phase.

The exponent causes considerable variation, so an approximation is needed,

22 2 '

'2

zR z

| ' |z , (A.43)

so that,

2'

2(1 cos )(1 cos ) ' '

jkzjkRjke

dz e e dzR

. (A.44)

Deformation of the contour and a new variable,

4'2

j ku z e

, (A.45)

yields,

42

'j

dz e duk

. (A.46)

The integration in 'z becomes,

24

1 cos 2(1 cos ) '

jkRj

jk uedz e e e du

R k

, (A.47)

where,

2

.ue du

(A.48)

The received field at the point ( , ,0)x y becomes,

76

4

1 cos( , ,0) '

2 2

j

jke kE x y E e dy

, (A.49)

where,

2 2( ')x y y and cos x

. (A.50)

The geometry is given in Figure A.3:

Figure A.3 – Plane wave diffraction by a half screen.

Repeated Kirchhoff Integral for Parallel Screens

Figure A.4 – Parallel screens with varying height.

(0, ',0)y y

y

z

x

( , ,0)x y

E ny 1ny y

n

n

reference planex

1 1( , )n nx y

1n

1n

2n 3n n 1n

77

Considering equation (A.49) and Figure A.4, the field 1 1( , ,0)n nE x y on the plane of 1ny can

be found from the field at the screen ny ,

4

1 1

1 cos( , ,0) ( , )

2 2n

n

j

jknn n n n n

h n

e kE x y E x y e dy

. (A.51)

For small angles where cos cos 2n n the equation becomes:

4

1 1( , ) ( , ) n

n

jjk

n n n n n

h n

e eE x y E x y dy

, (A.52)

where nh is the lower limit of the ny integration which corresponds to the height difference

between the top of the building and the reference plane.

78

APPENDIX B

In Appendix B the MATLAB codes of the plots in Chapter 2 is presented in order:

Code #1

TWO RAY GROUND REFLECTION MODEL

ht = 60;hr = 2; lambda1= 0.33; lambda2 = 3.3;

pt = 5; %watt

gt = 1.5;d = 1:.01:1000;

e0 = (sqrt(30*pt*gt))./d; % V/m

k1 = 2*pi*(ht*hr)./(lambda1.*(d));k2 = 2*pi*(ht*hr)./(lambda2.*(d));

e1 = 2*e0.*sin(k1);e2 = 2*e0.*sin(k2);

subplot(2,1,2);

plot(log(d),20*log10(e1/0.001),'-b'); hold on;

plot(log(d),20*log10(e0/0.001),'-r');

title('--- Two ray model for 900 MHz ---')

xlabel('distance in log (m)');ylabel('Field strength (dB)')

legend('resultant wave', 'direct wave')

subplot(2,1,1);

plot(log(d),20*log10(e2/0.001),'-b');hold on;

plot(log(d),20*log10(e0/0.001),'-r');

title('--- Two ray model for 90 MHz ---')

xlabel('distance in log (m)');ylabel('Field strength (dB)')

legend('resultant wave', 'direct wave')

Code #2

CORNU SPIRAL

v = -5:0.02:5;

C = mfun('FresnelC', v);

S = mfun('FresnelS', v);

figure,plot(C,S,'k','LineWidth',1.7)

xlabel('C(v) = Re F(v)');ylabel('S(v) = Im F(v)')

title('Cornu spiral and resultant phasor comparisons');hold on;

plot([-0.5,0.5],[-0.5,0.5],'-b','LineWidth',1.5);hold on;

plot([0.7799,0.5],[0.4383,0.5],'-r','LineWidth',1.5);hold on;

plot([-0.7228,0.5],[-0.2493,0.5],'-g','LineWidth',1.5)

axis([-0.8 0.8 -0.8 0.8])

axis square

grid on

79

Code #3

KNIFE EDGE DIFFRACTION

v=-6:0.02:6;

C=mfun('FresnelC',v);

S=mfun('FresnelS',v);

Erel = 0.5*((1-C-S)+1j*(C-S));

subplot(1,2,1);

plot(v,abs(Erel),'k','LineWidth',1.5)

grid on;

title('Knife-edge model')

ylabel('Resultant field strength rel. to direct field (linear ratio)')

xlabel('v parameter')

axis square;

subplot(1,2,2);

plot(v,20*log10(abs(Erel)),'-b','LineWidth',1.5)

grid on;

title('Knife-edge model')

ylabel('Resultant field strengt relative to direct field (dB)')

xlabel('v parameter')

axis square;

Code #4

DIFFRACTION LOSS FOR TWO FREQUENCIES

y=-15:0.01:15;

v=-2*y/sqrt(30);

C=mfun('FresnelC',v);

S=mfun('FresnelS',v);

C2=mfun('FresnelC',2*v);

S2=mfun('FresnelS',2*v);

Erel1 = 0.5*((1-C-S)+1j*(C-S));

Erel2 = 0.5*((1-C2-S2)+1j*(C2-S2));

plot(y,abs(Erel1),'k','LineWidth',1.5);hold on;

plot(y,abs(Erel2),'-b','LineWidth',1.4);grid on;

title('Knife-edge model');

ylabel('|E_t_o_t / E_i_n_c|')

xlabel('y (meters)');legend('1 GHz', '4 GHz')

80

APPENDIX C

In Appendix C the MATLAB codes of the plots in Chapter 3 is presented in order:

Code #1

SETTLING BEHAVIOUR OF PLANE WAVE ALONG SCREENS

maxN=500; %maximum number of screens

M=200; % M is maximum m

ratio=250; %ratio = d/lambda

g =[(sin(0)*sqrt(ratio)) %parameter g = sin(delta)*sqrt(d/lambda)

(sin(pi/360)*sqrt(ratio))

(sin(pi/180)*sqrt(ratio))

(sin(pi/180)*sqrt(2*ratio))

(sin(pi/90)*sqrt(ratio))

(sin(pi/72)*sqrt(ratio))];

U=1/(2*sqrt(pi));I1=zeros(maxN,M);I1(1,1)=U; I10=ones(maxN+1,1);

%creating auxiliary matrixes%

n=0;

while n < maxN

n=n+1;

I10(n+1)=I10(n)*(1/2+n-1)/n;

end

N=1;

for N=2:maxN,

Sigm=0;

for n=0:N-1,

Sigm=Sigm+I10(n+1)/sqrt(N-n);

end

I1(N,1)=U*Sigm;

end

I1(1,2)=1/4;

for N=2:maxN,

Sigm=0;

for n=1:N-1,

Sigm=Sigm+I1(n,1)/sqrt(N-n);

end

I1(N,2)=N/2*I10(N+1)+U*Sigm;

end

I1(:,2)=I1(:,2)/2;

for m=3:M

I1(1,m)=1/2 * I1(1,m-2)/m;

for N=2:maxN

Sigm=0;

for n=1:N-1

Sigm=Sigm+I1(n,m-1)/sqrt(N-n);

end

81

I1(N,m)=N/2*I1(N,m-2)/m + U*Sigm/m;

end

end

% E field %

U=2*sqrt(1j*pi).*g(:);

absE=zeros(length(g),maxN);

for N=1:maxN,

Sigm=I10(N+1);

for m=1:M,

Sigm=Sigm+I1(N,m).*U.^m;

end

absE(:,N)=abs(Sigm); % matrix of magnitude of E

end

%plotting%

end1=fix(10./(g(1).^2));

if end1>maxN, end1=maxN; end

semilogx(absE(1,1:end1),'.r');hold on

end2=fix(10./(g(2).^2));

if end2>maxN, end2=maxN; end

semilogx(absE(2,1:end2),'.m');hold on

end3=fix(10./(g(3).^2));

if end3>maxN, end3=maxN; end

semilogx(absE(3,1:end3),'.b');hold on

end4=fix(10./(g(4).^2));

if end4>maxN, end4=maxN; end

semilogx(absE(4,1:end4),'-.b');hold on

end5=fix(10./(g(5).^2));

if end5>maxN, end5=maxN; end

semilogx(absE(5,1:end5),'.g');hold on

end6=fix(10./(g(6).^2));

if end6>maxN, end6=maxN; end

semilogx(absE(6,1:end6),'.k'); hold on;

legend('\delta=0^{o} , f','\delta=0.5^{o}, f','\delta=1^{o} , f',...

'\delta=1^{o}, 2f','\delta=2^{o} , f','\delta=2.5^{o}, f')

%determining the settling points%

x=fix(1./(g(2).^2)); y=(absE(2,x)); plot(x,y,'-ok','MarkerFaceColor','k');

x=fix(1./(g(3).^2)); y=(absE(3,x));plot(x,y,'-ok','MarkerFaceColor','k');

x=fix(1./(g(4).^2));y=(absE(4,x));plot(x,y,'-ok','MarkerFaceColor','k');

x=fix(1./(g(5).^2));y=(absE(5,x));plot(x,y,'-ok','MarkerFaceColor','k');

x=fix(1./(g(6).^2));y=(absE(6,x));plot(x,y,'-ok','MarkerFaceColor','k');

axis([1 1000 0 1.5]);xlabel('N number of screens');ylabel('| E (x_{N+1}, 0) |')

title('Settling behavior of the Plane Wave along the screens')

82

Code #2

SETTLING FIELD WRT G PARAMETER WITH CURVE AND LINE FITTINGS

maxN=500; %maximum number of screens

M=200; % M is maximum m

ratio=250; %ratio = d/lambda

g=0:0.005:1;

U=1/(2*sqrt(pi));I1=zeros(maxN,M);I1(1,1)=U; I10=ones(maxN+1,1);

%creati1ng auxiliary matrixes%

n=0;

while n < maxN

n=n+1;

I10(n+1)=I10(n)*(1/2+n-1)/n;

end

N=1;

for N=2:maxN,

Sigm=0;

for n=0:N-1,

Sigm=Sigm+I10(n+1)/sqrt(N-n);

end

I1(N,1)=U*Sigm;

end

I1(1,2)=1/4;

for N=2:maxN,

Sigm=0;

for n=1:N-1,

Sigm=Sigm+I1(n,1)/sqrt(N-n);

end

I1(N,2)=N/2*I10(N+1)+U*Sigm;

end

I1(:,2)=I1(:,2)/2;

for m=3:M

I1(1,m)=1/2 * I1(1,m-2)/m;

for N=2:maxN

Sigm=0;

for n=1:N-1

Sigm=Sigm+I1(n,m-1)/sqrt(N-n);

end

I1(N,m)=N/2*I1(N,m-2)/m + U*Sigm/m;

end

end

% E field %

U=2*sqrt(1j*pi).*g(:);

83

absE=zeros(length(g),maxN);

for N=1:maxN,

Sigm=I10(N+1);

for m=1:M,

Sigm=Sigm+I1(N,m).*U.^m;

end

absE(:,N)=abs(Sigm); % matrix of magnitude of E

end

%Determinig settling values wrt g parameter and approximations

x=0:0.005:1;

t=3.502*x-3.327*x.^2+0.962*x.^3;

b=0.1*(x/0.03).^0.9;

subplot(1,2,1);

plot(log(x),log(t),'LineWidth',1.24);

for k=1:length(g),

stop=fix(10./(g(k).^2)); % maximum N for valid approximation

if stop>maxN, stop=maxN; end

hold on

subplot(1,2,1);

plot(log(g(k)),log(absE(k,stop)),'-ok');

legend('P(g)','|E_{settling}|',2)

end

xlabel('g')

title('|E_{settling}| & its polynomial fitting P(g)' );axis square

subplot(1,2,2);hold on;

plot(log(x),log(t),'LineWidth',1.24); hold on;

subplot(1,2,2);

plot(log(x),log(b),'--r','LineWidth',1.24);axis square

xlabel('g');title('P(g) and its line approximation')

legend('P(g)','Line approx',2)

84

APPENDIX D

In Appendix D the MATLAB codes of the plots in Chapter 4 is presented in order:

Code #1

PLOTS FOR THE HATA MODEL

r=1:0.01:10;

%Hata Frequency dependence

m=1.5; t=30; f1=800; f2=1000; f3 = 1200;

amurb=3.2*(log10(11.75*m)).^2-4.97;%hata large urban

LHurbf1 = 69.55+26.2*log10(f1)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r);

LHurbf2 = 69.55+26.2*log10(f2)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r);

LHurbf3 = 69.55+26.2*log10(f3)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r);

amsurbf1=2*(log10(f1/28)).^2+5.4;amsurbf2=2*(log10(f2/28)).^2+5.4;amsurbf3=2*(log10(f

3/28)).^2+5.4;%hata suburban

LHsubrf1 = 69.55+26.2*log10(f1)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amsurbf1;

LHsubrf2 = 69.55+26.2*log10(f2)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amsurbf2;

LHsubrf3 = 69.55+26.2*log10(f3)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amsurbf3;

amrurf1=4.78*(log10(f1)).^2-18.33*log10(f1)+40.94;

amrurf2=4.78*(log10(f2)).^2-18.33*log10(f2)+40.94;

amrurf3=4.78*(log10(f3)).^2-18.33*log10(f3)+40.94;

LHrurf1 = 69.55+26.2*log10(f1)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amrurf1;

LHrurf2 = 69.55+26.2*log10(f2)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amrurf2;

LHrurf3 = 69.55+26.2*log10(f3)-13.82*log10(t)-amurb+(44.9-6.55*log10(t))*log10(r)-

amrurf3;

figure(1);

plot(r,LHurbf1,'r',r,LHurbf2,'--r',r,LHurbf3,':r');hold on;

plot(r,LHsubrf1,'b',r,LHsubrf2,'--b',r,LHsubrf3,':b');hold on;

plot(r,LHrurf1,'g',r,LHrurf2,'--g',r,LHrurf3,':g');hold on;

legend('large urban f=800 MHz','large urban f=1000 MHz','large urban f=1200 MHz',...

'suburban f=800 MHz','suburban f=1000 MHz','suburban f=1200 MHz','rural f=800

MHz','rural f=1000 MHz','rural f=1200 MHz',4);

grid on;xlabel('d [km]');ylabel('L [dB]');title('Hata Model for different frequencies in different

environments')

%Hata height dependence

m=1.5; t1=30; t2=100; t3=200; f=100;

amurb=3.2*(log10(11.75*m)).^2-4.97;%hata large urban

LHurbt1 = 69.55+26.2*log10(f)-13.82*log10(t1)-amurb+(44.9-6.55*log10(t1))*log10(r);

85

LHurbt2 = 69.55+26.2*log10(f)-13.82*log10(t2)-amurb+(44.9-6.55*log10(t2))*log10(r);

LHurbt3 = 69.55+26.2*log10(f)-13.82*log10(t3)-amurb+(44.9-6.55*log10(t3))*log10(r);

amsurb=2*(log10(f/28)).^2+5.4;%hata suburban

LHsubrt1 = 69.55+26.2*log10(f)-13.82*log10(t1)-amurb+(44.9-6.55*log10(t1))*log10(r)-

amsurb;

LHsubrt2 = 69.55+26.2*log10(f)-13.82*log10(t2)-amurb+(44.9-6.55*log10(t2))*log10(r)-

amsurb;

LHsubrt3 = 69.55+26.2*log10(f)-13.82*log10(t3)-amurb+(44.9-6.55*log10(t3))*log10(r)-

amsurb;

amrur=4.78*(log10(f)).^2-18.33*log10(f)+40.94;

LHrurt1 = 69.55+26.2*log10(f)-13.82*log10(t1)-amurb+(44.9-6.55*log10(t1))*log10(r)-

amrur;

LHrurt2 = 69.55+26.2*log10(f)-13.82*log10(t2)-amurb+(44.9-6.55*log10(t2))*log10(r)-

amrur;

LHrurt3 = 69.55+26.2*log10(f)-13.82*log10(t3)-amurb+(44.9-6.55*log10(t3))*log10(r)-

amrur;

figure(2);

plot(r,LHurbt1,'r',r,LHurbt2,'--r',r,LHurbt3,':r');hold on;

plot(r,LHsubrt1,'b',r,LHsubrt2,'--b',r,LHsubrt3,':b');hold on;

plot(r,LHrurt1,'g',r,LHrurt2,'--g',r,LHrurt3,':g');hold on;

legend('large urban Ht=30 m','large urban Ht=100 m','large urban Ht=200 m',...

'suburban Ht=30 m','suburban Ht=100 m','suburban Ht=200 m','rural Ht=30 m','rural Ht=100

m','rural Ht=200 m',4);

grid on;xlabel('d [km]');ylabel('L [dB]');title('Hata Model for different transmitter heights in

different environments')

Code #2

PLOTS FOR THE WALFISCH-BERTONI MODEL

r=0:0.01:6;

m=1.7; d=40; b=15; t=35; f1=800; f2=1000; f3=1200;

rho=sqrt((d/2).^2+(b-m).^2);

L1f1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f1)-18*log10(t-b)+38*log10(r);

L2f2 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f2)-18*log10(t-b)+38*log10(r);

L3f3 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f3)-18*log10(t-b)+38*log10(r);

figure(1);plot(r,L1f1,'r',r,L2f2,'k',r,L3f3,'b');

legend('f_1= 800 MHz','f_2= 1000 MHz','f_3= 1200 MHz',4); xlabel('distance [km]');

ylabel('Loss [dB]')

title('Frequency dependence (Hm=1.7 d=40 Hb=15 Ht=35)')

t1=30; t2=35; t3=40; f=1000;

L1t1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t1-b)+38*log10(r);

L2t2 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t2-b)+38*log10(r);

86

L3t3 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t3-b)+38*log10(r);

figure(2);plot(r,L1t1,'r',r,L2t2,'k',r,L3t3,'b');

legend('Ht_1= 30m','Ht_2= 35m','Ht_3= 40m',4); xlabel('distance [km]'); ylabel('Loss [dB]')

title('Transmitter height dependence (f=1000 Hm=1.7 d=40 Hb=15)')

d1=10; d2=40; d3=70;

L1d1 = 89.5-10*log10((rho*d1.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

L1d2 = 89.5-10*log10((rho*d2.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

L1d3 = 89.5-10*log10((rho*d3.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

figure(3);plot(r,L1d1,'r',r,L1d2,'k',r,L1d3,'b');

legend('d_1= 10m','d_2= 40m','d_3= 70m',4); xlabel('distance [km]'); ylabel('Loss [dB]')

title('Building separation dependence (f=1000 Hm=1.7 Hb=15 Ht=35)')

b1=10; b2=15; b3=20;

L1b1 = 89.5-10*log10((rho*d.^0.9)/((b1-m).^2))+21*log10(f)-18*log10(t-b1)+38*log10(r);

L1b2 = 89.5-10*log10((rho*d.^0.9)/((b2-m).^2))+21*log10(f)-18*log10(t-b2)+38*log10(r);

L1b3 = 89.5-10*log10((rho*d.^0.9)/((b3-m).^2))+21*log10(f)-18*log10(t-b3)+38*log10(r);

figure(4);plot(r,L1b1,'r',r,L1b2,'k',r,L1b3,'b')

legend('Hb_1= 10m','Hb_2= 15m','Hb_3= 20m',4); xlabel('distance [km]'); ylabel('Loss [dB]')

title('Building height dependence (f=1000 Hm=1.7 d=40 Ht=35)')

Code #3

COMPARISON OF THE HATA AND THE WALFISCH BERTONI MODELS

r=1:0.5:10;

m=1.5; t=30; f1=800; f2=1500; f = 900;

aml=3.2*(log10(11.75*m)).^2-4.97;%hata large urban area

LHf1 = 69.55+26.2*log10(f1)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LHf2 = 69.55+26.2*log10(f2)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

d=50; b=10;

rho=sqrt((d/2).^2+(b-m).^2); %bertoni

LBf1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f1)-18*log10(t-b)+38*log10(r);

LBf2 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f2)-18*log10(t-b)+38*log10(r);

%different frequencies

figure(1);subplot(121); plot(r,LHf1,'r',r,LBf1,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('f = 800 MHz');legend('Hata','Bertoni',4);

subplot(122); plot(r,LHf2,'r',r,LBf2,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('f = 1500 MHz');legend('Hata','Bertoni',4);

t2=60;%different transmitter heights

LHt1 = 69.55+26.2*log10(f)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LHt2 = 69.55+26.2*log10(f)-13.82*log10(t2)-aml+(44.9-6.55*log10(t2))*log10(r);

LBt1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

LBt2 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t2-b)+38*log10(r);

figure(2);subplot(121); plot(r,LHt1,'r',r,LBt1,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('Ht = 30 m');legend('Hata','Bertoni',4);

subplot(122); plot(r,LHt2,'r',r,LBt2,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('Ht = 60 m');legend('Hata','Bertoni',4);

87

b2=15;%different building heights

rho2=sqrt((d/2).^2+(b2-m).^2);

LHb1 = 69.55+26.2*log10(f)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LHb2 = 69.55+26.2*log10(f)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LBb1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

LBb2 = 89.5-10*log10((rho2*d.^0.9)/((b2-m).^2))+21*log10(f)-18*log10(t-b2)+38*log10(r);

figure(3);subplot(121); plot(r,LHb1,'r',r,LBb1,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('Hb = 10 m');legend('Hata','Bertoni',4);

subplot(122); plot(r,LHb2,'r',r,LBb2,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('Hb = 15 m');legend('Hata','Bertoni',4);

d2=20;%different building separation distances

LHd1 = 69.55+26.2*log10(f)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LHd2 = 69.55+26.2*log10(f)-13.82*log10(t)-aml+(44.9-6.55*log10(t))*log10(r);

LBd1 = 89.5-10*log10((rho*d.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

LBd2 = 89.5-10*log10((rho*d2.^0.9)/((b-m).^2))+21*log10(f)-18*log10(t-b)+38*log10(r);

figure(4);subplot(121); plot(r,LHd2,'r',r,LBd2,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('d = 20 m');legend('Hata','Bertoni',4);

subplot(122); plot(r,LHd1,'r',r,LBd1,'b');grid on;xlabel('distance [km]');

ylabel('L [dB]');title('d = 50 m');legend('Hata','Bertoni',4);

88

REFERENCES

[1] Bertoni H.L., Radio Propagation for Modern Wireless Systems, Prentice-Hall PTR,

Upper Saddle River, NJ, 2000.

[2] Balanis C.A., Antenna Theory: Analysis and Design, 3rd Ed., Wiley, New York, 2005.

[3] Guru B.S. and Hiziroglu H.R., Electromagnetic Field Theory Fundamentals, PWS,

New York, 1998.

[4] Rappaport T.S., Wireless Communications: Principles and Practice, Prentice Hall PTR,

New Jersey, 1996.

[5] Goldsmith A., Wireless Communications, Cambridge University Press, New York,

2005.

[6] Parsons J.D., The Mobile Radio Propagation Channel, 2nd Ed., John Wiley & Sons,

2000.

[7] Haslett C., Essentials of Radio Wave Propagation, Cambridge University Press, New

York, 2008.

[8] Barclay L.W. and Hewitt M.T., Propagation of Radiowaves, 2nd Ed., IEE Press,

London, 1996.

[9] Vogler L.E., “The attenuation of Electromagnetic waves by multiple knife-edge

diffraction,” NTIA Report, USA, 1981.

[10] Saunders S.R., Antennas And Propagation For Wireless Communication Systems, 2nd

Ed., John Wiley & Sons, Chichester, 1999

[11] Blaunstein N., Radio Propagation in Cellular Networks, Artech House, London, 2000.

[12] Eliades D.E., Sharples P.A., Mehler, M.J. and Howell R., “Comparison of different

methods for modelling terrain diffraction loss at SHF,” IEEE AP 2, 377-381,1989.

[13] Tzaras C. and Saunders S.R., “Comparison of multiple-diffraction models for digital

broadcasting coverage prediction,” IEEE Transactions on Broadcasting 46, 221-226, 2000.

[14] Koutitas G. and Tzaras C., “A Slope UTD Solution for a Cascade of Multishaped

Canonical Objects,” IEEE AP 54, 2969-2976, 2006.

[15] Andersen J.B., “UTD Multiple-Edge Transition Zone Diffraction,” IEEE AP 45, 1093-

1097, 1997.

[16] McNamara D.A., Pistorius C.W.I. and Malherbe J.A.G., Introduction to the Uniform

Geometrical Theory of Diffraction, Artech House, Norwood, MA, 1990.

89

[17] Vaughan R. and Andersen J.B., Channels, Propagation and Antennas for Mobile

Communications, IEE Press, London, 2002.

[18] Hviid T., Andersen J.B., Toftgard J. and Bojer J., “Terrain-Based Propagation Model

for Rural Area - An Integral Equation Approach,” IEEE AP 43, 41-46, 1995.

[19] Li L. and Carin L., “Multi-level fast multipole calibration of ray models with

application to wireless propagation,” IEEE AP 52, 2794-2799, 2004.

[20] Okumura Y., Ohmori E., Kawano T., and f*ckuda K., “Field strength and its variability

in VHF and UHF land mobile radio service,” Review of the Electrical Communications

Laboratories (Tokyo) 16, 825-873, 1968.

[21] Hata M., “Empirical Formula for Propagation loss in Land Mobile Radio Services,”

IEEE VT 29, 317–325, 1980.

[22] Jackson J.D., Classical Electrodynamics, 3rd Ed., John Wiley & Sons, New York,

1999.

90

Matematiska och systemtekniska institutionen

SE-351 95 Växjö

Tel. +46 (0)470 70 80 00, fax +46 (0)470 840 04

http://www.vxu.se/msi/

91