forked from Amisto/butterfly
-
Notifications
You must be signed in to change notification settings - Fork 2
/
2_hilbertize.py
executable file
·96 lines (84 loc) · 3.09 KB
/
2_hilbertize.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
import png
def final(file, val, pref):
max_data = max(map(max, val))
min_data = min(map(min, val))
data_color = [[int(255*(x-min_data)/(max_data-min_data)) for x in l] for l in val]
with open(file+pref, 'wb') as fpng:
w = png.Writer(len(val[0]), len(val), greyscale=True)
w.write(fpng, data_color)
dir = 'data/baseline/Sensor{}' # function to create final files
import matplotlib.pyplot as plt # directory for saving
import os
import re
import math
import json
import csv
import numpy as np
from scipy.signal import hilbert
import matplotlib
matplotlib.use('Agg')
# Readin config
with open('res/config.json') as jf:
configuration = json.load(jf)
window = configuration['WINDOW']
cnt = 0
cnt_dirs = 0
file_names = []
for i in range(0, 32):
file_name = os.listdir('data/baseline/Sensor{}/'.format(i))
file_names.append(file_name)
raw_data=[]
tpl = r'^raw\d+.csv$'
for dr in file_names:
for i in dr:
if re.match(tpl, i):
raw_data.append(i)
else:
continue
for one_file in raw_data:
# for saving files by dirs
with open(dir.format(cnt) + '/'+ one_file) as fi:
with open(dir.format(cnt)+"/hilb.csv", 'w') as fo:
vals_init = [[float(x) for x in l.split()] for l in fi.readlines()]
vals = list(zip(*vals_init))
# Narrowband filtering + Hilbert transformation
res_fft = []
maxf = 0
for i, row in enumerate(vals):
fft = np.fft.rfft(row)
# this is to check the spectrum and find the carrying freauency
# near-zero freqencies somehow are exremely large
plt.plot(np.abs(fft)[10:])
maxf += np.abs(fft)[10:].argmax()
# plt.clf()
mid = maxf/len(vals)
for i, row in enumerate(vals):
fft = np.fft.rfft(row)
for j, freq in enumerate(fft):
if not (mid - window/2.0 <= j <= mid + window/2.0):
fft[j] = 0
res_fft.append(np.abs(hilbert(np.fft.irfft(fft))))
print(mid)
# saving spectrum
plt.savefig(dir.format(cnt) +'/spectrum' + ".png" )
plt.cla()
res = list(zip(*res_fft))
max_data = max(map(max, res))
min_data = min(map(min, res))
for r in res:
stri = ""
for t in r:
stri += str(round((t - min_data)/(max_data-min_data), 2))
stri += " "
print(stri, file=fo)
with open(dir.format(cnt) + "/2hilb.csv", 'w') as fu:
for r in vals_init:
stri = ""
for t in r:
stri += str(round(t, 2))
stri += " "
print(stri, file=fu)
final(file=dir.format(cnt), val=vals_init, pref="/prehilb.png")
# writing final graphics
final(file=dir.format(cnt), val=res, pref="/2hilb.png")
cnt += 1 # plussing to copy next files in another folder