Compare commits
7 commits
fff11e4648
...
58d5188c11
| Author | SHA1 | Date | |
|---|---|---|---|
| 58d5188c11 | |||
| ae4e8ec66a | |||
| 129f3cb7e7 | |||
| 1626446c72 | |||
| f334161b70 | |||
| cd9ced7b95 | |||
|
|
0da36eb56e |
2 changed files with 76 additions and 27 deletions
|
|
@ -27,7 +27,7 @@ fp_out = args.outfile
|
||||||
current_time = timeit.default_timer()
|
current_time = timeit.default_timer()
|
||||||
filename, rest = fp.split(".")
|
filename, rest = fp.split(".")
|
||||||
#filename = input("Filename (without extension): ")
|
#filename = input("Filename (without extension): ")
|
||||||
file = open(filename+'.AHA')
|
file = open(filename+'.AHA', encoding='utf-8', errors='ignore')
|
||||||
content = file.readlines()
|
content = file.readlines()
|
||||||
file.close()
|
file.close()
|
||||||
|
|
||||||
|
|
@ -75,8 +75,11 @@ for i in range(len(EDB_Array)):
|
||||||
# EDBnpf_array[i][j-7] = float(parts[j+1])
|
# EDBnpf_array[i][j-7] = float(parts[j+1])
|
||||||
#else:
|
#else:
|
||||||
# EDBnp_array[i][j] = int(parts[j+1])
|
# EDBnp_array[i][j] = int(parts[j+1])
|
||||||
|
try:
|
||||||
if not(int(parts[3])>= 24) and float(parts[5])/10000 > min_value and float(parts[5])/10000 < max_value:
|
if not(int(parts[3])>= 24) and float(parts[5])/10000 > min_value and float(parts[5])/10000 < max_value:
|
||||||
len_array[int(parts[2])][int(parts[3])] += 1
|
len_array[int(parts[2])][int(parts[3])] += 1
|
||||||
|
except:
|
||||||
|
pass
|
||||||
#EDB_File.close()
|
#EDB_File.close()
|
||||||
print("Element counting completed", timeit.default_timer()-current_time)
|
print("Element counting completed", timeit.default_timer()-current_time)
|
||||||
current_time = timeit.default_timer()
|
current_time = timeit.default_timer()
|
||||||
|
|
|
||||||
86
seth_hk.py
86
seth_hk.py
|
|
@ -23,6 +23,8 @@ fp = args.file
|
||||||
plot_flag = args.plot
|
plot_flag = args.plot
|
||||||
live = args.live
|
live = args.live
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
global time
|
global time
|
||||||
global data
|
global data
|
||||||
time = []
|
time = []
|
||||||
|
|
@ -64,8 +66,12 @@ def counter(l):
|
||||||
lower_count = (math.fsum(sl0)-sl0[3]-sl0[11]-sl0[20]-sl0[11])/20
|
lower_count = (math.fsum(sl0)-sl0[3]-sl0[11]-sl0[20]-sl0[11])/20
|
||||||
if lower_count <= 0:
|
if lower_count <= 0:
|
||||||
lower_count = 'nan'
|
lower_count = 'nan'
|
||||||
if upper_count == 'nan' or lower_count == 'nan':
|
if upper_count == 'nan' and lower_count == 'nan':
|
||||||
return [upper_count, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count] # divide by 12 to get counts per second
|
return [upper_count, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count] # divide by 12 to get counts per second
|
||||||
|
elif upper_count == 'nan':
|
||||||
|
return [upper_count, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count/12] # divide by 12 to get counts per second
|
||||||
|
elif lower_count == 'nan':
|
||||||
|
return [upper_count/12, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count] # divide by 12 to get counts per second
|
||||||
else:
|
else:
|
||||||
return [upper_count/12, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count/12] # divide by 12 to get counts per second
|
return [upper_count/12, BGO_upper/12, HETA/12, HETB/12, BGO_lower/12, lower_count/12] # divide by 12 to get counts per second
|
||||||
|
|
||||||
|
|
@ -136,14 +142,6 @@ def calibrate_vectors(acc_raw, mag_raw):
|
||||||
"""
|
"""
|
||||||
mag_sens = 6842 # LSB/gauss
|
mag_sens = 6842 # LSB/gauss
|
||||||
acc_sens = 16000 # 1000 LSB/g and it is only a 12 bit number so divide by 16 (1hex bit) equivalent to shifting 4 bits
|
acc_sens = 16000 # 1000 LSB/g and it is only a 12 bit number so divide by 16 (1hex bit) equivalent to shifting 4 bits
|
||||||
mag_raw[0] = mag_raw[0] + 700 #Values for 360-29e (guessed)
|
|
||||||
mag_raw[1] = mag_raw[1] +4000 # For tvw-36e use -500 and +2500
|
|
||||||
if mag_raw[0] < 0:
|
|
||||||
mag_raw[0] = mag_raw[0]+65535
|
|
||||||
if mag_raw[0] > 65535:
|
|
||||||
mag_raw[0] = mag_raw[0]-65535
|
|
||||||
if mag_raw[1] > 65535:
|
|
||||||
mag_raw[1] = mag_raw[1]-65535
|
|
||||||
mag_calib = [1.02305, 1.13811, 1.16184, 0.0265714, -0.188856, -0.177294] # fit coefficients from BA thesis
|
mag_calib = [1.02305, 1.13811, 1.16184, 0.0265714, -0.188856, -0.177294] # fit coefficients from BA thesis
|
||||||
acc_calib = [0.972148, 0.968608, 0.974061, -0.00109779, -0.0154984, 0.0135722]
|
acc_calib = [0.972148, 0.968608, 0.974061, -0.00109779, -0.0154984, 0.0135722]
|
||||||
acc_raw = [struct.unpack('<h', struct.pack('<H', int(num)))[0] for num in acc_raw] # turns int into hex, turns hex into signed int
|
acc_raw = [struct.unpack('<h', struct.pack('<H', int(num)))[0] for num in acc_raw] # turns int into hex, turns hex into signed int
|
||||||
|
|
@ -178,16 +176,14 @@ def calculate_heading(mag, phi, theta):
|
||||||
thus, the heading is simply the angle between x and north. Translate this to compass.
|
thus, the heading is simply the angle between x and north. Translate this to compass.
|
||||||
"""
|
"""
|
||||||
declination = 0 # Kiruna
|
declination = 0 # Kiruna
|
||||||
#phi = 0
|
|
||||||
theta = 0 # Set zero because it works...
|
|
||||||
mag_world = np.matmul(Rz(-phi), np.matmul(Ry(theta-np.pi), np.matmul(Rz(phi), mag)))
|
mag_world = np.matmul(Rz(-phi), np.matmul(Ry(theta-np.pi), np.matmul(Rz(phi), mag)))
|
||||||
mag_world_psi = np.arctan2(mag_world[1], mag_world[0]) + declination # angle mag vector in wf, x axis in wf + delta because TN is delta left of the MN (=mag vector)
|
mag_world_psi = np.arctan2(mag_world[1], mag_world[0]) + declination # angle mag vector in wf, x axis in wf + delta because TN is delta left of the MN (=mag vector)
|
||||||
if mag_world_psi <= 0:
|
if mag_world_psi <= 0:
|
||||||
mag_world_psi = 2*np.pi + mag_world_psi
|
#mag_world_psi = 2*np.pi + mag_world_psi
|
||||||
#mag_world_psi = abs(mag_world_psi) # because compass counts clock-wise and atan2 counter CW
|
mag_world_psi = abs(mag_world_psi) # because compass counts clock-wise and atan2 counter CW
|
||||||
else:
|
else:
|
||||||
pass
|
#pass
|
||||||
#mag_world_psi = 2*np.pi - mag_world_psi
|
mag_world_psi = 2*np.pi - mag_world_psi
|
||||||
heading = mag_world_psi * 180/np.pi # heading in degrees (0deg = True North)
|
heading = mag_world_psi * 180/np.pi # heading in degrees (0deg = True North)
|
||||||
return heading
|
return heading
|
||||||
|
|
||||||
|
|
@ -206,14 +202,20 @@ def altitude(pressure):
|
||||||
'''
|
'''
|
||||||
Uses the barometric formula to solve for the height
|
Uses the barometric formula to solve for the height
|
||||||
'''
|
'''
|
||||||
|
if type(pressure) is int or type(pressure) is float:
|
||||||
H = 8400 # scale height for T=288K
|
H = 8400 # scale height for T=288K
|
||||||
return -math.log(pressure/1013)*H
|
return -math.log(pressure/1013)*H
|
||||||
|
else:
|
||||||
|
return 'nan'
|
||||||
|
|
||||||
def main(input):
|
def main(input):
|
||||||
i = -1 # HK comes in 12s intervals i ended by C64 lines: ... C64 P H P HD HD P C64 P H P HD HD P C64 ...
|
i = -1 # HK comes in 12s intervals i ended by C64 lines: ... C64 P H P HD HD P C64 P H P HD HD P C64 ...
|
||||||
global time
|
global time
|
||||||
global data
|
global data
|
||||||
global live
|
global live
|
||||||
|
time_wrong_begin = False
|
||||||
|
time_wrong_middle = False
|
||||||
|
time_changed = False
|
||||||
dorn_temps = [] # [Tadc0, Tpa0, Tbgo01, Tbgo02, Tpwr, Tadc1, Tpa1, Tbgo11, Tbgo12, Text] slice 0 and slice1
|
dorn_temps = [] # [Tadc0, Tpa0, Tbgo01, Tbgo02, Tpwr, Tadc1, Tpa1, Tbgo11, Tbgo12, Text] slice 0 and slice1
|
||||||
temp = []
|
temp = []
|
||||||
pres = []
|
pres = []
|
||||||
|
|
@ -277,18 +279,49 @@ def main(input):
|
||||||
|
|
||||||
if l[:2] == 'H ': # One every 12s
|
if l[:2] == 'H ': # One every 12s
|
||||||
ts = int(l.split()[1])
|
ts = int(l.split()[1])
|
||||||
if fp and fp.split('/')[-1] == "2025-08-14-longtime-23e.AHA" and ts > 1800000000: # timestamps wrong in this file
|
#if fp and fp.split('/')[-1] == "2025-08-14-longtime-23e.AHA" and ts > 1800000000: # timestamps wrong in this file
|
||||||
ts -= 1963615744
|
# ts -= 1963615744
|
||||||
if i < 2000 and ts >= 3718824668-1963615744:
|
# if i < 2000 and ts >= 3718824668-1963615744:
|
||||||
ts -= 12*3600
|
# ts -= 12*3600
|
||||||
if fp and fp.split('/')[-1] == '2025-08-08-seth-17e.AHA' and ts > 1800000000:
|
if fp and fp.split('/')[-1] == '2025-08-08-seth-17e.AHA' and ts > 1800000000:
|
||||||
ts -= 2138848512
|
ts -= 2138848512
|
||||||
if i < 500 and ts > 3893525631-2138848512:
|
if i < 500 and ts > 3893525631-2138848512:
|
||||||
ts -= 16404
|
ts -= 16404
|
||||||
|
|
||||||
|
#if live and (ts > 1800000000 or ts < 1700000000) and not time_wrong_begin and time == []:
|
||||||
|
# time_wrong_begin = True
|
||||||
|
#if live and (ts < 1800000000 and ts > 1700000000) and time_wrong_begin:
|
||||||
|
# for i in range(len(time)):
|
||||||
|
# time[i] += ts - time[-1] -12
|
||||||
|
# time_wrong_begin = False
|
||||||
|
# time_changed = True
|
||||||
|
# if live and (ts > 1800000000 or ts < 1700000000) and not time_wrong_middle and not time == []:
|
||||||
|
# time_wrong_middle == len(time)
|
||||||
|
# if live and (ts < 1800000000 and ts > 1700000000) and time_wrong_middle:
|
||||||
|
# for i in range(i, len(time)):
|
||||||
|
# time[i] += ts - time[-1] -12
|
||||||
|
# time_wrong_middle = False
|
||||||
|
# time_changed = True
|
||||||
|
#if live and len(time) and abs(ts - time[-1]) > 120:
|
||||||
|
#print(ts-time[-1])
|
||||||
|
#print(time, "time (pre)")
|
||||||
|
#print(len(time))
|
||||||
|
#print(list(data), "data")
|
||||||
|
#print(len(time)==len(list(data)), len(time), len(data))
|
||||||
|
# for i in range(len(time)-1):
|
||||||
|
# data[time[i] + ts - time[-1] -12] = data.pop(time[i])
|
||||||
|
# time[i] += ts - time[-1] -12
|
||||||
|
#print(i)
|
||||||
|
#time[i] += ts - time[-1] -12
|
||||||
|
# time_changed = True
|
||||||
|
#print(time, "time (post)")
|
||||||
|
#print(len(time))
|
||||||
if len(time) > i: # i is zero-indexed sometimes there are >1 H lines in one interval, use the last
|
if len(time) > i: # i is zero-indexed sometimes there are >1 H lines in one interval, use the last
|
||||||
time[-1] = ts
|
time[-1] = ts
|
||||||
else:
|
else:
|
||||||
time.append(ts) # [time0, time1, ...] times of the intervals
|
time.append(ts) # [time0, time1, ...] times of the intervals
|
||||||
|
if time_changed:
|
||||||
|
print(time) #debugging code
|
||||||
|
|
||||||
if l[:6] == 'HDORN ': # two lines, sl0 and sl1 per interval
|
if l[:6] == 'HDORN ': # two lines, sl0 and sl1 per interval
|
||||||
sl = l.split()[1]
|
sl = l.split()[1]
|
||||||
|
|
@ -302,8 +335,14 @@ def main(input):
|
||||||
|
|
||||||
if l[:4] == 'C64 ':
|
if l[:4] == 'C64 ':
|
||||||
heading, roll, pitch = calculate_ads(acc_raw, mag_raw)
|
heading, roll, pitch = calculate_ads(acc_raw, mag_raw)
|
||||||
|
if pres:
|
||||||
mean_pres = math.fsum(pres)/len(pres)
|
mean_pres = math.fsum(pres)/len(pres)
|
||||||
|
else:
|
||||||
|
mean_pres = 'nan'
|
||||||
|
if temp:
|
||||||
mean_temp = math.fsum(temp)/len(temp)
|
mean_temp = math.fsum(temp)/len(temp)
|
||||||
|
else:
|
||||||
|
mean_temp = 'nan'
|
||||||
if len(dorn_temps) < 10:
|
if len(dorn_temps) < 10:
|
||||||
dorn_temps = ['nan' for i in range(10)]
|
dorn_temps = ['nan' for i in range(10)]
|
||||||
data[time[-1]] = [mean_pres, mean_temp, counter(l), dorn_temps, [heading, roll, pitch]]
|
data[time[-1]] = [mean_pres, mean_temp, counter(l), dorn_temps, [heading, roll, pitch]]
|
||||||
|
|
@ -314,6 +353,13 @@ def main(input):
|
||||||
# Plotting
|
# Plotting
|
||||||
if plot_flag and live: # only plot if -pl option is True
|
if plot_flag and live: # only plot if -pl option is True
|
||||||
### Date formatting
|
### Date formatting
|
||||||
|
if time_changed:
|
||||||
|
#print(x, "x (pre)")
|
||||||
|
#print(data)
|
||||||
|
for i in range(len(time)-1):
|
||||||
|
x[i] = datetime.datetime.fromtimestamp(list(data)[i], datetime.timezone.utc)
|
||||||
|
time_changed = False
|
||||||
|
#print(x, "x (post)")
|
||||||
x.append(datetime.datetime.fromtimestamp(list(data)[-1], datetime.timezone.utc))
|
x.append(datetime.datetime.fromtimestamp(list(data)[-1], datetime.timezone.utc))
|
||||||
fmtdate = md.DateFormatter('%d-%H:%M:%S')
|
fmtdate = md.DateFormatter('%d-%H:%M:%S')
|
||||||
ax1.xaxis.set_major_formatter(fmtdate)
|
ax1.xaxis.set_major_formatter(fmtdate)
|
||||||
|
|
@ -383,7 +429,7 @@ def main(input):
|
||||||
|
|
||||||
|
|
||||||
if fp:
|
if fp:
|
||||||
with open(fp) as file:
|
with open(fp, encoding='utf-8', errors='ignore') as file:
|
||||||
input = file
|
input = file
|
||||||
main(input)
|
main(input)
|
||||||
if plot_flag and not live:
|
if plot_flag and not live:
|
||||||
|
|
|
||||||
Loading…
Add table
Add a link
Reference in a new issue