Problem 3¶

Part a: Using two Cornu spiral terms: $$|E(x,z)| = E_0 \frac{1}{\sqrt{2}} \left|Fr\left(\sqrt{\frac{2}{\lambda z}}(h_2 - x)\right) - Fr\left(\sqrt{\frac{2}{\lambda z}}(h_1 - x)\right)\right| $$

Take $h_1 = 5\lambda, ~~~ h_2 = -5\lambda$.

In [241]:

from scipy.special import fresnel

def fr(x):
    '''
    Return complex Fresnel function evaluated at x (needs to be array)
    '''
    fsc = fresnel(x) # tuples (real, imag)
    f = zeros(len(fsc[0]), dtype=complex)
    for i in range(len(f)):
        re = fsc[0][i]
        im = fsc[1][i]
        f[i] = re + 1j*im
        if x[i] == inf:
            f[i] = 0.5 + 1j*0.5
    return f

def plot_beam(z):
    '''
    Plot the beam shape at the given distance z
    '''
    lam = 1.0
    h_1 = -5*lam
    h_2 = 5*lam
    N = 1000
    x = linspace(-20*lam, 20*lam, N)
    Fsc = 1/sqrt(2)*abs(fr(sqrt(2/(lam*z))*(h_1-x)) - fr(sqrt(2/(lam*z))*(h_2-x)))
    plot(x,Fsc,'k-')
    xlabel('Transverse distance [m]')
    ylabel('Beam Amplitude')

In [242]:

figure(figsize=(7,4))
zvec = [1,10,100,1000]
for i in range(len(zvec)):
    subplot(2,2,i+1)
    plot_beam(zvec[i])
    title('z=%i m' % zvec[i])
tight_layout()

Part b

In [88]:

def interactive_beam(log10z):
    plt.figure()
    plot_beam(10**log10z)
    plt.show()

In [111]:

w = interactive(interactive_beam, log10z=(0.,3.))
display(w)

In [94]:

interactive_beam(1.7)

Based on my qualitative experiments with different values of z, I found that $10^{1.7} = 50$ m seems to be the transition to a far field beam pattern. I based this off of qualitatively deciding when the shape looked more like a gradual falloff than a torch-light.

In [89]:

rayleigh_distance = 2*(10*lam)**2/lam
print 'Rayleigh distance is %.1f m' % rayleigh_distance

Rayleigh distance is 200.0 m

My transition distance of 50 m is shorter than the Rayleigh distance of 200 m, but they are on the same order of magnitude. It seems that the concept of Rayleigh distance is applicable in this case.

Problem 4¶

Part a

In [104]:

x1 = 0
z1 = 30
x2 = 1000
z2 = 40
xe = [500,500]
ze = [0,40]
figure(figsize=(6,2))
plot([x1,x2],[z1,z2],'ko')
plot(xe,ze,'k-',linewidth=5)
xlim((-100,1100))
ylim((0,50))
xlabel('x [m]')
ylabel('z [m]')

Out[104]:

<matplotlib.text.Text at 0x7fc0d215b7d0>

From Problem 2, the edge-wave attenuation is given by:

$$ \frac{\sqrt{\lambda}}{\sqrt{r}} \frac{1}{2 \pi \sin \theta} $$

The following approximations are reasonable for this problem:

$$ r \approx 500 $$$$ sin \theta \approx \tan \theta = \frac{10}{500} $$

Plugging in, and taking -20*log10 to get dB attenuation, we get:

In [107]:

lam = 1./3
r = 500
sinth = 10./500.
G = sqrt(lam)/sqrt(r)/(2*pi*sinth)
print 'Attenuation is %.2f dB' % (-20*log10(G))

Attenuation is 13.75 dB

Part b

Let's check the two requirements specified in Problem 2 for the edge-wave approximation:

$\sqrt{\frac{2}{\lambda z}} |h-x| > 1 $. See calculation below
$\sin{\theta} << 1$. Yes, easily.

In [109]:

sqrt(2/(lam*500))*10

Out[109]:

1.0954451150103321

Since this value is close to 1, our point is very near to the parabola edge, as discussed in class. Thus, the edge-wave approximation is suspect. I would use the full theory in this case.

Problem 5¶

Part a: $\hat{p} = -\sin \theta_i \hat{x} + cos \theta_i \hat{y}$. This ensures that $\hat{E} \times \hat{H} = \hat{k}$, where we take $\hat{H}$ in the $+\hat{z}$ direction.

Part b: We start with the following FT pair, which ignores the delta function and thus implicitly assumes that $k_y$ has some negative imaginary part:

$$ u(y) \leftrightarrow \frac{1}{j k_y} $$

Using the modulation property of the Fourier transform, we obtain the FT of $E_y(0,y)$:

$$ E_y(0,y) = E_0 \cos\theta_i\, e^{-jk \sin \theta_i y} u(y) \leftrightarrow \frac{E_0 \cos \theta_i}{j (k_y - k \sin \theta_i)} $$

So we write the plane wave spectrum as

$$ P(s) = \frac{E_0 c_i}{jk(s - s_i)}$$

with $s = \frac{k_y}{y}$, $s_i = \sin \theta_i$, and $c_i = \cos \theta_i$.

As discussed above, the region of convergence is $Im(s) < 0$

Part c: Interpreting $s$ as a direction cosine, $\sin\theta$, we arrive at the edge-wave expression, which is valid near (but not too near) to the shadow edge. In this region, $c_i \approx 1$. From the course notes (i.e., directly from Huygen's):

$$E(x,y) = \frac{e^{-jkr}}{\sqrt{r}} \sqrt{\frac{j}{\lambda}} P(s) [-s,c,0] = \frac{E_0 e^{-jkr}}{\sqrt{r}} \sqrt{\frac{j}{\lambda}}\frac{1}{jk(s-s_i)}$$

Part d: First we have to reorganize the above expression to relate it to the unobstructed wave. Taking magnitudes, we get an expression that compares nicely to the diffraction coefficient in Problem 2:

$$\frac{|E|}{E_0} =\frac{1}{k|s-s_i|\sqrt{r \lambda}} = \frac{1}{2\pi |s-s_i|} \sqrt{\frac{\lambda}{r}}$$

In [119]:

lam = 1./3
s = 0 # x = 0 ==> s = 0
si = sin(pi/6) # given in problem
r = 30 
gain = 1.0/(2*pi*abs(s-si)) * sqrt(lam/r)
print 'Attenuation is %.2f dB' % (-20*log10(gain))

Attenuation is 29.49 dB

Problem 7¶

In [437]:

def F(x,sig2,f):
    ''' Return ground wave factor in dB '''
    lam = 3.0e8/f
    epsr = 15.0
    eps0 = 8.85e-12
    eps2 = epsr*eps0
    x0 = lam/pi*(epsr + sig2/(1j*2*pi*f*eps0))
    
    nx = x/x0

    frterm2sc = fresnel(sqrt(2.0/pi * nx))
    frterm2 = frterm2sc[1] - 1j*frterm2sc[0]
    frterm1sc = fresnel(1e20)
    frterm1 = (frterm1sc[1] - 1j*frterm1sc[0])*ones_like(frterm2)
    expterm = exp(1j*nx)
    sqrtterm = 1j*sqrt(2*pi*nx)
    F = 1-sqrtterm*expterm*(frterm1 - frterm2)
    FdB = 20*log10(abs(F))
    
    return FdB

def bisection(f,x0,x1, xtol=1e-3):
    '''
    Solve f(x) == 0 for x using the bisection method.
    x0 and x1 are two values that straddle the solution.
    x0 < x1.
    '''
    f0 = f(x0)
    f1 = f(x1)
    
    while(abs(x1-x0) > xtol):

        xn = sqrt(x0*x1) # bisect equally in log domain
        fn = f(xn)
        if sign(fn)==sign(f0):
            x0 = xn
            f0 = fn
        else:
            x1 = xn
            f1 = fn
        if sign(f0)==sign(f1):
            raise Exception('Something went horribly wrong')
    return x1

In [457]:

fvec = logspace(5,8,100)
xvec0 = zeros(len(fvec))
xvec1 = zeros(len(fvec))

for i in range(len(fvec)):
    f = fvec[i]
    
    ######## sig = 1e-3 #######
    # Define function that equals 0 at -40dB point
    myfunc = lambda x: F(x, 1e-3, f)+40
    # Solve myfunc(x)==0
    xvec0[i] = bisection(myfunc, 1e-5, 1e10)
    
    ######## sig = 1e-2 #######
    # Define function that equals 0 at -40dB point
    myfunc = lambda x: F(x, 1e-2, f)+40
    # Solve myfunc(x)==0
    xvec1[i] = bisection(myfunc, 1e-5, 1e10)

In [459]:

loglog(fvec,xvec0/1e3,'k-',linewidth=2,label='$\sigma = 10^{-3}$')
loglog(fvec,xvec1/1e3,'k--',linewidth=2,label='$\sigma = 10^{-2}$')
ylim((1e-1,1e3))
legend(loc='best')
xlabel('f (MHz)')
ylabel('x (km) at -40 dB')

Out[459]:

<matplotlib.text.Text at 0x7fc0be6d0590>

Problem 8¶

Part a: By inspection of the Collins 1985 plot, at 1 MHz the -40 dB point occurs at $p \approx 15$. At 2MHz, it is $p \approx 33$.

To convert this to actual distance: $x = p |x_0| = p \frac{\lambda}{\pi} \left|\epsilon_r + \frac{\sigma_2}{j \omega \epsilon_0}\right|$, which is calculated below.

In [469]:

epsr = 15
eps0 = 8.85e-12
sig2 = 1e-2
xkm = lambda((p,f)): p * 3e5/f/pi*abs(epsr + sig2/(1j*2*pi*f*eps0))
print 'At 1MHz, the -40 dB distance is %.0f km' % xkm((15, 1e6))
print 'At 2MHz, the -40 dB distance is %.0f km' % xkm((33, 2e6))

At 1MHz, the -40 dB distance is 258 km
At 2MHz, the -40 dB distance is 144 km

Part b: The -40 dB distance is larger at 1 MHz. This is relevant to AM radio, since it means that the lower-frequency stations will transmit farther than the higher-frequency ones.

Problem 3¶

Problem 4¶

Problem 5¶

Problem 7¶

Problem 8¶

Convert to pdf or slides¶