In [1]:

    

%matplotlib inline


    

In [2]:

    

import matplotlib.pyplot as plt
import numpy as np

from pyMPM import MPM


    

In [3]:

    

plt.style.use('seaborn')


    

In [4]:

    

fig, ax = plt.subplots(figsize=(14,6))

f = np.arange(1,100)     # Vector of frequencies in GHz
T = 15                # Air temperature in degree Celcius
P = 1013              # Air pressiure in mbar
RH_vec = [0, 50, 100] # List of relative humidity values used for plotting

for RH in RH_vec:
    A = MPM(f, P, T, RH, 0, 0, 0, 'att')
    ax.semilogy(f,A, label='RH ' + str(RH) + '%')

ax.set_xlabel('f in GHz')
ax.set_ylabel('A in db/km')
plt.legend(loc=2);


    

In [5]:

    

f_vec = [15,22.235,30, 80, 110]# Vector of frequencies in GHz
T = 30                # Air temperature in degree Celcius
P = 1013              # Air pressiure in mbar
RH_vec = np.arange(100)  # Vector of relative humidity values


    

In [6]:

    

A_dict = {}
for f in f_vec:
    A = [MPM(f, P, T, RH, 0, 0, 0, 'att') for RH in RH_vec]
    A_dict[str(f) + 'GHz'] = A


    

In [7]:

    

fig, ax = plt.subplots(figsize=(14,6))

for f in f_vec:
    key = str(f) + 'GHz'
    ax.plot(RH_vec, A_dict[key], label=key)
plt.legend(loc=2)
ax.set_xlabel('Relative humidity in %')
ax.set_ylabel('Attenuation in db/km');


    

In [8]:

    

def sat_water_vap_pres(T):
    '''
    sat_water_vap_pres(T)
    
    T     = Temperature in ∞C
    e_sat = Saturated water vapor pressure in HPa (=mbar)
    
    Calculates the saturated water vapor pressure with
    the approximation by Bolton (1980)
    '''
    return 6.112 * np.exp(17.67 * T / (T + 243.5))


def abs_hum(RH,T):
    '''
    Calculate absolute humidity from temperature and relative humidity
    
    abs_hum(rh,t)
    
    abs_hum = absolute humidity on g/m^3
    RH      = relative humidity in percent
    T       = temperature in degree C
    '''

    R = 8.3144621e-2 # gas constant in hPa*m^2/(mol*K)
    M = 18.02       # molar weight in g/mol
    TK = T + 273.15
    
    abs_h = RH / 100.0 * sat_water_vap_pres(T) * M / R / TK
    return abs_h


    

In [9]:

    

fig, ax = plt.subplots(figsize=(14,6))

abs_hum_vec = abs_hum(RH_vec,T)

for f in f_vec:
    key = str(f) + 'GHz'
    ax.plot(abs_hum_vec, A_dict[key], label=key)
plt.legend(loc=2)
ax.set_xlabel('Absolute humidity in g/m$^3$')
ax.set_ylabel('Attenuation in db/km');


    

In [10]:

    

def fit_linear(xd,yd):
    # determine best fit line
    par = np.polyfit(xd, yd, 1, full=True)
    slope=par[0][0]
    intercept=par[0][1]
    return slope[0], intercept[0]


    

In [11]:

    

fig, ax = plt.subplots(figsize=(14,6))

for f in f_vec:
    key = str(f) + 'GHz'
    ax.plot(abs_hum_vec, A_dict[key], marker='o', linestyle='none', label=key)
    slope, intercept = fit_linear(abs_hum_vec, A_dict[key])
    ax.plot(abs_hum_vec, intercept + slope*abs_hum_vec, 'k--')
plt.legend(loc=2)
ax.set_xlabel('Absolute humidity in g/m$^3$')
ax.set_ylabel('Attenuation in db/km');


    

In [12]:

    

T_vec = [0,20,40]
f_vec = range(1,100,1)

P = 1013                 # Air pressiure in mbar
RH_vec = np.arange(100)  # Vector of relative humidity values

fit_param_T = {}
for T in T_vec:
    abs_hum_vec = abs_hum(RH_vec, T)
    slope_vec = np.zeros_like(f_vec, dtype=float)
    intercept_vec = np.zeros_like(f_vec, dtype=float)
    for i_f,f in enumerate(f_vec):
        A = [MPM(f, P, T, RH, 0, 0, 0, 'att') for RH in RH_vec]
        slope_vec[i_f], intercept_vec[i_f] = fit_linear(abs_hum_vec, A)
    fit_param_T[str(T)] = {'slope': slope_vec,
                           'intercept': intercept_vec}


    

In [13]:

    

fig, ax = plt.subplots(2,1, sharex=True, figsize=(14,8))

for T in T_vec:
    ax[0].plot(f_vec, fit_param_T[str(T)]['slope'], label=str(T) + '$^oC$')
    ax[1].plot(f_vec, fit_param_T[str(T)]['intercept'], label=str(T) + '$^oC$')
fig.suptitle('Linear fit of absolute humidty vs. attenuation', fontsize=14)
ax[0].legend(loc=2)
ax[1].legend(loc=2)
ax[0].set_ylabel('slope')
ax[1].set_ylabel('intercept');


    

In [ ]: