You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
258 lines
8.2 KiB
Python
258 lines
8.2 KiB
Python
import datetime
|
|
import time
|
|
import math
|
|
|
|
# Brett Miller 9th October 2015.
|
|
# This program reads Harmonics for Fort Denison as provided by BOM on the 9th October 2015
|
|
#
|
|
# NOT FOR DISTRIBUTION
|
|
# the predictions from this code have been cross checked against some daily predictions but
|
|
# not against all.
|
|
# USE AT YOUR OWN CAUTION. NO GAURENTEE MADE THAT THESE PREDICTIONS ARE CORRECT.
|
|
#
|
|
# The x-tide method for calculating tides is:-
|
|
# z(t,y,site)= MWL(site) + sum[c=1,N](nodefactor(c,y).amp(c,site).cos( speed(c).t + eqarg(c,y) - kappa(c,site) ))
|
|
# where
|
|
# t is the time since the beginning of the year
|
|
# y is the year
|
|
# site is the site of interest
|
|
# N is the number of consituents being included
|
|
# nodefactor is an amplitude adjustment for each constituent for each year
|
|
# eqarg is a phase adjustment for each site for each year
|
|
# amp is the amplitude of each constituent for the site
|
|
# kappa is the phase adjustment for each constituent for the site
|
|
#
|
|
# the harmonics included below in the python code were provided directly to WRL from BOM.
|
|
|
|
harm = {}
|
|
speed = 0
|
|
amp = 1
|
|
lag = 2
|
|
|
|
|
|
def readlinenocomm(F):
|
|
l = F.readline()
|
|
while l[0] == '#':
|
|
l = F.readline()
|
|
return l
|
|
|
|
|
|
#
|
|
H = open("HARMONIC", "r")
|
|
l = readlinenocomm(H)
|
|
numconstits = int(l)
|
|
|
|
print("Number of Constituents = ", numconstits)
|
|
|
|
# Read in the constituent speeds
|
|
|
|
constitspeed = {}
|
|
i = 0
|
|
while i < numconstits:
|
|
l = readlinenocomm(H)
|
|
ll = l.split()
|
|
constitspeed[ll[0]] = float(ll[1])
|
|
i = i + 1
|
|
|
|
l = readlinenocomm(H)
|
|
sy = int(l)
|
|
l = readlinenocomm(H)
|
|
numyears = int(l)
|
|
|
|
# Read in the equalibrium arguements
|
|
|
|
equalibrium_arg = {}
|
|
i = 0
|
|
while i < numconstits:
|
|
l = readlinenocomm(H)
|
|
ll = l.split()
|
|
constit = ll[0]
|
|
v = []
|
|
j = 0
|
|
while j <= numyears / 10:
|
|
l = readlinenocomm(H)
|
|
ll = l.split()
|
|
for lll in ll:
|
|
v.append(float(lll))
|
|
j = j + 1
|
|
equalibrium_arg[constit] = v
|
|
i = i + 1
|
|
|
|
l = readlinenocomm(H)
|
|
if l[0:5] != '*END*':
|
|
print("Error, expected end after eq args. Got: ", l)
|
|
|
|
# Read in the node factors
|
|
|
|
l = readlinenocomm(H)
|
|
numyears = int(l)
|
|
|
|
node_factor = {}
|
|
i = 0
|
|
while i < numconstits:
|
|
l = readlinenocomm(H)
|
|
ll = l.split()
|
|
constit = ll[0]
|
|
v = []
|
|
j = 0
|
|
while j <= numyears / 10:
|
|
l = readlinenocomm(H)
|
|
ll = l.split()
|
|
for lll in ll:
|
|
v.append(float(lll))
|
|
j = j + 1
|
|
node_factor[constit] = v
|
|
i = i + 1
|
|
|
|
l = readlinenocomm(H)
|
|
if l[0:5] != '*END*':
|
|
print("Error, expected end after node factors. Got: ", l)
|
|
|
|
harm['sa'] = [0.0410686, 0.04, 44.3373]
|
|
harm['ssa'] = [0.0821373, 0.0233, 140.8642]
|
|
harm['Mm'] = [0.5443747, 0.0022, 163.3925]
|
|
harm['MSF'] = [1.0158958, 0.0019, 73.4003]
|
|
harm['MF'] = [1.0980331, 0.0025, 203.8471]
|
|
harm['2Q1'] = [12.8542862, 0.0036, 9.8664]
|
|
harm['SIGMA1'] = [12.9271398, 0.0043, 16.6066]
|
|
harm['Q1'] = [13.3986609, 0.023, 52.3723]
|
|
harm['RHO1'] = [13.4715145, 0.0043, 55.9113]
|
|
harm['O1'] = [13.9430356, 0.0963, 79.6431]
|
|
harm['MP1'] = [14.0251729, 0.001, 45.9175]
|
|
harm['M1'] = [14.4920521, 0.0044, 99.5325]
|
|
harm['CHI1'] = [14.5695476, 0.0013, 102.4847]
|
|
harm['PI1'] = [14.9178647, 0.0028, 125.9361]
|
|
harm['P1'] = [14.9589314, 0.0439, 115.4037]
|
|
harm['S1'] = [15, 0.0033, 129.8967]
|
|
harm['K1'] = [15.0410686, 0.148, 119.6599]
|
|
harm['PSI1'] = [15.0821353, 0.002, 11.0025]
|
|
harm['PHI1'] = [15.1232059, 0.0015, 103.8567]
|
|
harm['THETA1'] = [15.5125897, 0.0017, 137.3796]
|
|
harm['J1'] = [15.5854433, 0.0096, 142.7051]
|
|
harm['SO1'] = [16.0569644, 0.0017, 173.2765]
|
|
harm['OO1'] = [16.1391017, 0.0062, 180.0375]
|
|
harm['OQ2'] = [27.3416964, 0.0006, 75.4143]
|
|
harm['MNS2'] = [27.4238337, 0.0061, 193.1893]
|
|
harm['2N2'] = [27.8953548, 0.0205, 199.6233]
|
|
harm['MEU2'] = [27.9682084, 0.021, 212.8684]
|
|
harm['N2'] = [28.4397295, 0.1127, 223.5144]
|
|
harm['NEU2'] = [28.5125831, 0.0206, 223.5377]
|
|
harm['OP2'] = [28.9019669, 0.0024, 128.4082]
|
|
harm['M2'] = [28.9841042, 0.5022, 237.061]
|
|
harm['MKS2'] = [29.0662415, 0.0014, 313.3572]
|
|
harm['LAMDA2'] = [29.4556253, 0.0058, 224.4572]
|
|
harm['L2'] = [29.5284789, 0.0135, 243.1576]
|
|
harm['T2'] = [29.9589333, 0.0073, 281.2802]
|
|
harm['S2'] = [30, 0.1247, 260.7625]
|
|
harm['R2'] = [30.0410667, 0.0022, 206.1633]
|
|
harm['K2'] = [30.0821373, 0.0373, 250.3219]
|
|
harm['MSN2'] = [30.5443747, 0.0016, 124.387]
|
|
harm['KJ2'] = [30.626512, 0.0023, 62.9849]
|
|
harm['2SM2'] = [31.0158958, 0.0019, 145.4123]
|
|
harm['MO3'] = [42.9271398, 0.0004, 148.654]
|
|
harm['M3'] = [43.4761563, 0.0023, 327.537]
|
|
harm['SO3'] = [43.9430356, 0.0002, 307.6915]
|
|
harm['MK3'] = [44.0251729, 0.0005, 248.2678]
|
|
harm['SK3'] = [45.0410686, 0.0016, 113.2852]
|
|
harm['MN4'] = [57.4238337, 0.001, 90.6269]
|
|
harm['M4'] = [57.9682084, 0.0031, 114.2894]
|
|
harm['SN4'] = [58.4397295, 0.0003, 103.5135]
|
|
harm['MS4'] = [58.9841042, 0.0014, 183.7788]
|
|
harm['MK4'] = [59.0662415, 0.0004, 196.3685]
|
|
harm['S4'] = [60, 0.0006, 299.9646]
|
|
harm['SK4'] = [60.0821373, 0.0002, 45.1703]
|
|
harm['2MN6'] = [86.407938, 0.0009, 6.8059]
|
|
harm['M6'] = [86.9523127, 0.002, 52.2368]
|
|
harm['MSN6'] = [87.4238337, 0.0007, 64.4301]
|
|
harm['2MS6'] = [87.9682084, 0.0025, 114.6867]
|
|
harm['2MK6'] = [88.0503457, 0.0007, 131.372]
|
|
harm['2SM6'] = [88.9841042, 0.0003, 157.8974]
|
|
harm['MSK6'] = [89.0662415, 0.0003, 204.7339]
|
|
harm['2MN2S2'] = [26.4079379, 0.0003, 1.3443]
|
|
harm['3M(SK)2'] = [26.8701753, 0.0005, 3.6145]
|
|
harm['3M2S2'] = [26.9523127, 0.0008, 78.7665]
|
|
harm['MNK2S2'] = [27.505971, 0.0001, 225.556]
|
|
harm['SNK2'] = [28.3575922, 0.0003, 218.9001]
|
|
harm['2SK2'] = [29.9178627, 0.0011, 321.0943]
|
|
harm['MQ3'] = [42.3827651, 0.0004, 43.8537]
|
|
harm['2MP3'] = [43.009277, 0.0001, 134.0684]
|
|
harm['2MQ3'] = [44.5695475, 0.0001, 329.7593]
|
|
harm['3MK4'] = [56.8701754, 0, 39.731]
|
|
harm['3MS4'] = [56.9523127, 0.0001, 97.9154]
|
|
harm['2MSK4'] = [57.8860711, 0.0002, 283.9222]
|
|
harm['3MK5'] = [71.9112441, 0.0002, 182.3575]
|
|
harm['M5'] = [72.4602605, 0, 24.2117]
|
|
harm['3MO5'] = [73.0092771, 0.0005, 264.7303]
|
|
harm['2(MN)S6'] = [84.8476674, 0, 129.7631]
|
|
harm['3MNS6'] = [85.3920422, 0.0001, 260.1927]
|
|
harm['4MK6'] = [85.8542795, 0.0001, 332.2482]
|
|
harm['4MS6'] = [85.9364168, 0.0001, 249.9087]
|
|
harm['2MSNK6'] = [86.3258006, 0, 280.3885]
|
|
harm['2MNU6'] = [86.4807915, 0.0003, 28.4884]
|
|
harm['3MSK6'] = [86.8701754, 0.0001, 222.0912]
|
|
harm['4MN6'] = [87.4966873, 0.0003, 292.3403]
|
|
harm['3MSN6'] = [88.5125832, 0.0006, 338.3166]
|
|
harm['MKL6'] = [88.5947204, 0.0001, 176.1921]
|
|
harm['2MN8'] = [114.8476674, 0.0001, 219.9735]
|
|
harm['3MN8'] = [115.3920422, 0.0002, 245.3801]
|
|
harm['M8'] = [115.9364168, 0.0002, 284.0175]
|
|
harm['2MSN8'] = [116.4079379, 0.0002, 203.8838]
|
|
harm['3MS8'] = [116.9523127, 0.0001, 280.466]
|
|
harm['3MK8'] = [117.03445, 0.0001, 286.3325]
|
|
harm['MSNK8'] = [117.505971, 0, 219.0404]
|
|
harm['2MS8'] = [117.9682084, 0.0002, 188.3974]
|
|
harm['2MSK8'] = [118.0503457, 0.0001, 195.8549]
|
|
harm['4MS10'] = [145.9364168, 0.0002, 293.7057]
|
|
harm['3M2S10'] = [146.9523127, 0.0002, 307.5005]
|
|
harm['4MSN12'] = [174.3761463, 0.0001, 110.9933]
|
|
harm['5MS12'] = [174.920521, 0, 180.6601]
|
|
harm['4M2S12'] = [175.9364168, 0, 108.791]
|
|
harm['SP3'] = [44.9589314, 0.0019, 112.6443]
|
|
harm['S3'] = [45, 0.0012, 23.4633]
|
|
harm['MA2'] = [28.9430375, 0.0039, 117.2716]
|
|
harm['Ma2'] = [29.0251709, 0.0027, 150.6304]
|
|
harm['MSV2'] = [30.4715211, 0, 98.718]
|
|
harm['SKM2'] = [31.0980331, 0.0008, 156.4092]
|
|
harm['2MNS4'] = [56.407938, 0.0001, 26.7109]
|
|
harm['MV4'] = [57.4966873, 0.0003, 94.1175]
|
|
harm['3MN4'] = [58.5125832, 0.0003, 339.4208]
|
|
harm['2MSN4'] = [59.5284789, 0.0001, 264.4342]
|
|
harm['NA2'] = [28.3986628, 0.0004, 16.9727]
|
|
harm['MSO5'] = [72.9271398, 0.0001, 247.7639]
|
|
harm['MSK5'] = [74.0251729, 0.0002, 189.3859]
|
|
|
|
MSL = 0.9952
|
|
|
|
# Now do the maths
|
|
|
|
start_time = datetime.datetime(2019, 11, 29, 0, 0, 0)
|
|
end_time = datetime.datetime(2019, 12, 5, 0, 0, 0)
|
|
|
|
time_step = datetime.timedelta(minutes=15)
|
|
|
|
x = []
|
|
y = []
|
|
F = open("tide.out", "w")
|
|
t = start_time
|
|
while (t < end_time):
|
|
offset = t - datetime.datetime(t.year, 1, 1, 0, 0, 0)
|
|
yr = t.year - sy
|
|
t1 = (offset.days) * 24.0 + offset.seconds / 3600.0
|
|
z = MSL
|
|
for constit in harm.keys():
|
|
if constit in node_factor:
|
|
nf = node_factor[constit][yr]
|
|
ea = equalibrium_arg[constit][yr]
|
|
else:
|
|
nf = 1.0
|
|
ea = 0.0
|
|
z = z + nf * harm[constit][amp] * math.cos(
|
|
(harm[constit][speed] * t1 + ea - harm[constit][lag]) *
|
|
(2 * 3.1415927 / 360.0))
|
|
|
|
# print t,z
|
|
F.write("%s,%s\n" % (t, z))
|
|
x.append(t)
|
|
y.append(z)
|
|
t = t + time_step
|