-
Notifications
You must be signed in to change notification settings - Fork 0
/
LISA_filter.py
158 lines (133 loc) · 4.35 KB
/
LISA_filter.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
#!/usr/bin/env python3
import sys
import numpy as np
import matplotlib.pyplot as plt
from scipy import fftpack
from datetime import datetime, timedelta
import time
import tqdm
#load data
fnamedat='datafiles/LISA_full_25_gf.txt'
cols=2
#fnametimes='datafiles/g1_huge_file_times.dat'
if len(sys.argv) == 2:
fnamedat = sys.argv[1]# they can override the file name
#fnametimes = sys.argv[2]
elif len(sys.argv) == 3:
fnamedat = sys.argv[1]
#fnametimes = sys.argv[2]
cols = sys.argv[3]
else:
sys.stderr.write(f'no arguments given')
data=np.loadtxt(fname=fnamedat, usecols=range(cols))
file = open(f'datafiles/g2_alldat_neg1001.txt',"a")
file.truncate(0)
file.close()
#splitting data into individual colums
time=np.array(data[0:, 0])
combined=np.array(data[0:, 1])
#imaginary=np.array(data[0:, 2])
#data parse check
np.set_printoptions(suppress=True, precision=30)
#print(time, combined, imaginary)
#time_stuff=np.loadtxt(fname=fnametimes)
def main():
for i in range(combined.size):
if combined[i] > 1.5*10**-8:
#print(i, combined[i])
combined[i]=-1001
elif combined[i] < -1.5*10**-8:
#print(i, combined[i])
combined[i]=-1001
file = open(f'datafiles/g2_alldat_neg1001.txt',"a")
for i in range(combined.size):
file.write(str(time[i]))
file.write("\t")
file.write(str(combined[i]))
file.write("\n")
file.close()
freq=time
#freq=freq+gps
#
# plt.rcParams["figure.figsize"] = (16,11)
#
# figure, axis = plt.subplots(2,1)
#
# axis[0].plot(freq, combined)
# axis[0].set_title("TimeSeriesDomain: with Gaps Filled")
#
# #scale/unit of signal of time is currently unknown
# plt.setp(axis[0], xlabel="time(gps)")
# plt.setp(axis[0], ylabel="signal")
# ################33 LEGENDS NOT WORKING
# axis[0].legend()
#
# axis[1].plot(freq, combined)
# axis[1].set_title("TimeSeriesDomain: with Gaps filled")
#
# #scale/unit of signal of time is currently unknown
# plt.setp(axis[1], xlabel="time(gps)")
# plt.setp(axis[1], ylabel="signal")
#
# #saving plot
# #ftypes=['jpg', 'svg']
# ftypes=['png']
# saveplot('plots/LISA_neg_1001', ftypes)
#plt.show()
def get_data(fname):
#load data
data=np.loadtxt(fname, usecols=range(3))
#splitting column data into individual arrays
time=np.array(data[0:, 0])
combined=np.array(data[0:, 1])
#imaginary=np.array(data[0:, 2])
#data parse check
np.set_printoptions(suppress=True, precision=30)
#print(time, combined, imaginary)
return time, combined
def get_UTC_datetime(gps):
utc = datetime(1980, 1, 6) + timedelta(seconds=gps - (37-18))#apparently leap seconds between gps and utc team need to be calculated
print(utc)
return utc
def saveplot(title, filetypes):
for ftype in filetypes:
filename=f'{title}.{ftype}'
#print(f'saving file {filename}')
#plt.savefig(filename)
def saving_original_data_to_dat_file(time, combined):
#creating 3 column data
newdata=np.column_stack((time, combined))
#making data printable
originaldata=np.array2string(newdata, precision=30,suppress_small=True)
#checking data
print(newdata)
#data file and appending data in a way so that its human-readable, not binary and its in 3 column format
file = open('gaps_plot.txt',"a")
for i in range(combined.size):
file.write(float(time[i]))
file.write("\t")
file.write(str(combined[i]))
# file.write("\t")
#file.write(str(combined.imag[i]))
file.write("\n")
file.close()
def plot_original(times, sig, subp, ylab):
plt.subplot(subp)
plt.ylabel(ylab)
#not the best way to make title and labels, will improve later
#plt.title(ylab+" ", fontsize=13, ha="right")
#plt.ylabel("signal")
#plt.xlabel("time")
plt.plot(times, sig, label=ylab)
def plot_fft(freqs, sigfft, subp, ylab):
plt.subplot(subp)
plt.ylabel(ylab)
#not the best way to make the title and labels, will improve later
#plt.title(ylab+" ", fontsize=13, ha="right")
#plt.ylabel("signal strength")
#plt.xlabel("time")
markerline, stemlines, baseline = plt.stem(freqs, np.abs(sigfft), '-.')
plt.setp(stemlines, 'linewidth', 0.2)
# plt.stem(freqs, np.abs(sigfft))
if __name__ == '__main__':
main()