1038.sagews

#Index of Refraction (Exp Model): n = 1+10^-6*a*exp(-b*h)
#Where: n is Index of Refraction. n is also dimensionless.
#a is a dimensionless constant for a given set of atmospheric conditions. 
#b (per m) is a constant for a given set of atmospheric conditions.
#h is height in meters above mean sea level.

def expmodel(h,a,b): #Index of Refraction (Exponential Model) defined.
    return float(1+10^-6*a*exp(-b*h))

import scipy #Transfers files from scipy (open source library of algorithms and mathematical tools from Python programming language) into this worksheet.

def expindex(h0,hn,ht,a,b): #Function creates array of index of refraction values for expmodel.
    
    height = scipy.arange(h0,hn,ht) #Creates array range (arange) of height values. 
    n=[] #Stores index of refraction values.
    
    for h in height:
        n0 = expmodel(h,a,b)
        n.append(n0) #Attaches index of refraction values at end of variable n.
    return height,n

expindex(0,5,0.1,340,1.4090) #h=0,a=340,b=1.4090(per m)*10^-4. Parameters are for curve fit that would be obtained using data from 5 km and higher.

(array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1. , 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. , 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3. , 3.1, 3.2, 3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4. , 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9]), [1.00034, 1.0002953158960903, 1.0002565043484812, 1.0002227935633023, 1.0001935131787933, 1.0001680809347078, 1.0001459910936734, 1.0001268043842628, 1.0001101392657845, 1.0000956643410879, 1.0000830917665067, 1.0000721715278931, 1.0000626864689175, 1.000054447972768, 1.0000472922113772, 1.0000410768875911, 1.0000356784054929, 1.0000309894126151, 1.0000269166651641, 1.000023379173802, 1.0000203065931212, 1.0000176378227768, 1.0000153197924659, 1.0000133064065881, 1.0000115576276039, 1.0000100386798605, 1.0000087193580545, 1.0000075734265801, 1.0000065780978145, 1.0000057135789726, 1.0000049626785124, 1.0000043104642702, 1.000003743966525, 1.0000032519200861, 1.0000028245402772, 1.0000024533283618, 1.0000021309025398, 1.0000018508511559, 1.000001607605199, 1.0000013963275585, 1.0000012128168358, 1.0000010534237962, 1.0000009149788009, 1.0000007947287779, 1.0000006902824741, 1.0000005995629042, 1.0000005207660478, 1.0000004523249768, 1.0000003928786936, 1.0000003412450691])

h,n=expindex(0,5,0.1,340,1.4090) #Variable h,n is assigned height and index of refraction values.

import numpy as np #Transfer files from fundamental library numpy (N-dimensional array manipulation) used with Python.
import matplotlib.pyplot as plt #Transfer files from matplotlib.pyplot (plotting library) used with Python language. 
#Index of Refraction is always greater than 1.  Index of Refraction approaches 1 from above.

x = n
y = h

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x, y)
ax.set_xlim(1.00000031, 1.0003)
ax.set_ylim(0, 5)
ax.set_xlabel('Index of Refraction (n)')
ax.set_ylabel('Height (m)*10^4 above mean sea level')
ax.set_title('Index of Refraction')
ax.grid()

plt.savefig('myfig.png')

def expmodeldot(h,a,b): #First Derivative of Exponential Model defined.
    return float(1+10^-6*a*-b*exp(-b*h))

import scipy #Transfers files from scipy (open source library of algorithms and mathematical tools from Python programming language) into this worksheet.

def expindexdot(h0,hn,ht,a,b): #Function creates array of index of refraction values dependent on height.
    
    height = scipy.arange(h0,hn,ht) #Creates array range (arange) of height values. 
    n=[] #Stores index of refraction values.
    
    for h in height:
        n0 = expmodeldot(h,a,b)
        n.append(n0) #Attaches index of refraction values at end of variable n.
    return height,n

expindexdot(0,5,0.1,340,1.4090) #h=0,a=340,b=1.4090(per m)*10^-4. Parameters are for curve fit that would be obtained using data from 5 km and higher.

(array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1. , 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. , 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3. , 3.1, 3.2, 3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4. , 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9]), [0.99952094000000002, 0.99958389990240881, 0.99963858537299, 0.99968608386930702, 0.99972733993108009, 0.99976317396299663, 0.99979429854901425, 0.99982133262257378, 0.99984481377450962, 0.99986520894340714, 0.99988292370099219, 0.99989831031719867, 0.9999116747652953, 0.99992328280637, 0.99993336527416943, 0.99994212266538418, 0.99994972912666069, 0.99995633591762512, 0.99996207441878382, 0.99996705874411296, 0.99997138801029228, 0.99997514830770751, 0.99997841441241553, 0.99998125127311721, 0.99998371530270602, 0.99998585550007668, 0.99998771442450118, 0.99998932904194882, 0.99999073146017925, 0.99999194956722759, 0.9999930075859762, 0.99999392655584329, 0.99999472475116624, 0.99999541804459868, 0.99999602022274947, 0.99999654326033827, 0.99999699755832128, 0.99999739215072114, 0.99999773488427457, 0.99999803257446995, 0.9999982911410783, 0.99999851572587128, 0.99999871079486957, 0.99999888022715189, 0.99999902739199398, 0.99999915521586813, 0.99999926624063851, 0.99999936267410772, 0.99999944643392069, 0.99999951918569774])

h,n=expindexdot(0,5,0.1,340,1.4090) #Variable h,n is assigned height and index of refraction values.

import numpy as np #Transfer files from fundamental library numpy (N-dimensional array manipulation) used with Python.
import matplotlib.pyplot as plt #Transfer files from matplotlib.pyplot (plotting library) used with Python language. 
#Derivative of Index of Refraction approaches 0 from above.

x = n
y = h

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x, y)
ax.set_xlim(0.9999994, 0.9994)
ax.set_ylim(0, 5)
ax.set_xlabel('Index of Refraction (n)')
ax.set_ylabel('Height (m)*10^4 above mean sea level')
ax.set_title('First Derivative of Index of Refraction')
ax.grid()

plt.savefig('myfig.png')

def SWmodel(P,T,e): #Refractivity (Smith-Weintraub Model) defined.
    return float(77.6*(P/T)+3.73*10^5*(e/T^2))

SWmodel(0.002,4,0.001) #Function call to SWmodel with predetermined assignments.

23.351299999999998

#Thayer Model

def Zwinverse(e,T,t):
    return float(1+1650*(e/T^3)*(1-0.01317*t+1.75*10^-4*t^2+1.44*10^-6*t^3))

Zwinverse(0.001,4,0.002)

1.0257805709399221

def Thmodel(K1,Pd,T,Zdinverse,K2,e,K3,Zwinverse): #Refractivity (Thayer Model) defined.
    return float(K1*(Pd/T)*Zdinverse+(K2*(e/T)+K3*(e/T^2))*Zwinverse)

Thmodel(4,0.007,4,3,7,0.001,1,1.03) #Function call to Thmodel with predetermined assignments.

0.022866875000000002