xaveros/jar_project
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

04.08

Browse Source
This commit is contained in:
xaver 2020-08-04 14:36:06 +02:00
parent 0d98a7842e
commit 47e386000f
5 changed files with 154 additions and 103 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
jar_functions.py
View File

@ -10,7 +10,7 @@ def step_response(t, a1, a2, tau1, tau2):
return r_step
def sin_response(t, f, p, A):
r_sin = A*sin(2*np.pi*t*f + p)
r_sin = A*np.sin(2*np.pi*t*f + p)
return r_sin
def parse_dataset(dataset_name):
@ -172,7 +172,6 @@ def average(freq_all, time_all, start, stop, timespan, dm):
def import_data(dataset):
print(dataset)
import nixio as nix
nf = nix.File.open(dataset, nix.FileMode.ReadOnly)
b = nf.blocks[0]
@ -193,6 +192,17 @@ def import_data(dataset):
return dat, pre_dat, dt
#nf.close()
def import_amfreq(dataset):
import nixio as nix
nf = nix.File.open(dataset, nix.FileMode.ReadOnly)
b = nf.blocks[0]
eod = b.data_arrays['EOD-1']
dt = eod.dimensions[0].sampling_interval
di = int(50.0/dt)
t = b.tags['Beats_1']
amfreq = t.metadata['RePro-Info']['settings']['amfreq']
return amfreq
if __name__ == '__main__':
import_data(os.path.join('JAR', '2020-07-21-ak', '2020-07-21-ak.nix'))

9
scratch.py
View File

@ -1,5 +1,6 @@
import os
import numpy as np
import matplotlib.pyplot as plt
from IPython import embed
from jar_functions import mean_noise_cut
from matplotlib.mlab import specgram
@ -11,7 +12,7 @@ import DataLoader as dl
for info, key, time, data in dl.iload_traces(datapath, repro='Beats', before=0.0, after=0.0):
'''
base_path = 'D:\\jar_project\\JAR'
'''base_path = 'D:\\jar_project\\JAR'
#nicht: -5Hz delta f, 19-aa, 22-ae, 22-ad (?)
datasets = [#'2020-06-19-aa', #-5Hz delta f, horrible fit
@ -38,6 +39,10 @@ for d in range(int(len(data)/10)):
#print(freqs)
#print(times)
embed()
embed()'''
g = [1.2917623576698833, -5.479055166593157, -2.689492238578325, -0.11604244418416806, -0.05353823781665627]
a = [0.2, 0.002, 0.02, 0.5, 1.0]
plt.plot(a, g)
plt.show()

0
sin_response_fit.py Normal file
View File

135
sin_response_specto.py Normal file
View File

@ -0,0 +1,135 @@
import matplotlib.pyplot as plt
import matplotlib as cm
from matplotlib.colors import ListedColormap, LinearSegmentedColormap
from matplotlib.mlab import specgram
import os
import glob
import IPython
import numpy as np
import DataLoader as dl
from IPython import embed
from scipy.optimize import curve_fit
from jar_functions import step_response
from jar_functions import sin_response
from jar_functions import parse_dataset
from jar_functions import parse_infodataset
from jar_functions import mean_traces
from jar_functions import mean_noise_cut
from jar_functions import norm_function
from jar_functions import sort_values
from jar_functions import average
from jar_functions import import_data
from jar_functions import import_amfreq
base_path = 'D:\\jar_project\\JAR\\sin\\2018lepto98'
#nicht: -5Hz delta f, 19-aa, 22-ae, 22-ad (?)
#dat = glob.glob('D:\\jar_project\\JAR\\2020*\\beats-eod.dat')
#infodat = glob.glob('D:\\jar_project\\JAR\\2020*\\info.dat')
datasets = [#'2020-06-19-aa', #-5Hz delta f, horrible fit
#'2020-06-19-ab', #-5Hz delta f, bad fit
#'2020-06-22-aa', #-5Hz delta f, bad fit
#'2020-06-22-ab', #-5Hz delta f, bad fit
#'2020-06-22-ac', #-15Hz delta f, good fit
#'2020-06-22-ad', #-15Hz delta f, horrible fit
#'2020-06-22-ae', #-15Hz delta f, horrible fit
#'2020-06-22-af', #-15Hz delta f, good fit
'2020-07-21-al', #sin
'2020-07-21-am',
'2020-07-21-ak',
#'2020-07-21-an',
#'2020-07-21-ao'
]
time_all = []
freq_all = []
amfrequencies = []
gains = []
ID = []
col = ['dimgrey', 'grey', 'darkgrey', 'silver', 'lightgrey', 'gainsboro', 'whitesmoke']
labels = zip(ID, datasets)
for infodataset in datasets:
infodataset = os.path.join(base_path, infodataset, 'info.dat')
i = parse_infodataset(infodataset)
identifier = i[0]
if not identifier[1:-2] in ID:
ID.append(identifier[1:-2])
for idx, dataset in enumerate(datasets):
datapath = os.path.join(base_path, dataset, '%s.nix' % dataset)
print(datapath)
data, pre_data, dt = import_data(datapath)
nfft = 2**17
for d, dat in enumerate(data):
amfreq = float(import_amfreq(datapath))
print(amfreq)
amfrequencies.append(amfreq)
file_name = []
file_name.append(ID)
file_name.append(amfreq)
file_name.append(d)
spec, freqs, times = specgram(dat, Fs=1/dt, detrend='mean', NFFT=nfft, noverlap=nfft*0.6)
dbspec = 10.0*np.log10(spec) # in dB
power = dbspec[:, 50]
fish_p = power[(freqs > 400) & (freqs < 1000)]
fish_f = freqs[(freqs > 400) & (freqs < 1000)]
index = np.argmax(fish_p)
eodf = fish_f[index]
eodf4 = eodf * 4
lim0 = eodf4-10
lim1 = eodf4+15
# control of plt.imshow
df = freqs[1] - freqs[0]
ix0 = int(np.floor(lim0/df)) # back to index
ix1 = int(np.ceil(lim1/df)) # back to index
spec4 = dbspec[ix0:ix1, :]
freq4 = freqs[ix0:ix1]
jar4 = freq4[np.argmax(spec4, axis=0)] # all freqs at max specs over axis 0
jar = jar4 / 4
jm = jar4 - np.mean(jar4)
cut_times = times[:len(jar4)]
#np.save('%s.npy' % file_name, jar4)
plt.plot(cut_times, jm, '-k')
cf, ct = mean_noise_cut(jar4, cut_times, n = int(round(len(jar4)/((times[-1] - times [0]) * amfreq))))
#plt.plot(ct, cf, '-k')
#plt.imshow(spec4, cmap='jet', origin='lower', extent=(times[0], times[-1], lim0, lim1), aspect='auto', vmin=-80, vmax=-10)
#np.save( , spec4)
embed()
b, a = signal.butter(4, 0.01 / 10000, 'high', analog=True)
y = signal.filtfilt(b, a, jm)
sinv, sinc = curve_fit(sin_response, cut_times, jm, [amfreq, 2, 0.5])
print('frequency, phaseshift, amplitude:', sinv)
gains.append(sinv[2])
plt.plot(cut_times, sin_response(cut_times, *sinv), label='fit: f=%f, p=%.2f, A=%.2f' % tuple(sinv))
#plt.legend()
#plt.ylim(lim0, lim1)
plt.legend()
plt.show()
#print(np.load('%s.npy' % file_name))
embed()

99
step_response_specto.py
View File

@ -1,99 +0,0 @@
import matplotlib.pyplot as plt
import matplotlib as cm
from matplotlib.colors import ListedColormap, LinearSegmentedColormap
from matplotlib.mlab import specgram
import os
import glob
import IPython
import numpy as np
import DataLoader as dl
from IPython import embed
from scipy.optimize import curve_fit
from jar_functions import parse_dataset
from jar_functions import parse_infodataset
from jar_functions import mean_traces
from jar_functions import mean_noise_cut
from jar_functions import norm_function
from jar_functions import step_response
from jar_functions import sort_values
from jar_functions import average
from jar_functions import import_data
base_path = 'D:\\jar_project\\JAR'
#nicht: -5Hz delta f, 19-aa, 22-ae, 22-ad (?)
datasets = [#'2020-06-19-aa', #-5Hz delta f, horrible fit
#'2020-06-19-ab', #-5Hz delta f, bad fit
#'2020-06-22-aa', #-5Hz delta f, bad fit
#'2020-06-22-ab', #-5Hz delta f, bad fit
#'2020-06-22-ac', #-15Hz delta f, good fit
#'2020-06-22-ad', #-15Hz delta f, horrible fit
#'2020-06-22-ae', #-15Hz delta f, horrible fit
#'2020-06-22-af', #-15Hz delta f, good fit
'2020-07-21-ak' #sin
]
#dat = glob.glob('D:\\jar_project\\JAR\\2020*\\beats-eod.dat')
#infodat = glob.glob('D:\\jar_project\\JAR\\2020*\\info.dat')
time_all = []
freq_all = []
ID = []
col = ['dimgrey', 'grey', 'darkgrey', 'silver', 'lightgrey', 'gainsboro', 'whitesmoke']
labels = zip(ID, datasets)
for infodataset in datasets:
infodataset = os.path.join(base_path, infodataset, 'info.dat')
i = parse_infodataset(infodataset)
identifier = i[0]
ID.append(identifier)
for idx, dataset in enumerate(datasets):
datapath = os.path.join(base_path, dataset, '%s.nix' % dataset)
data, pre_data, dt = import_data(datapath)
nfft = 2**15
spec, freqs, times = specgram(data[0], Fs=1/dt, detrend='mean', NFFT=nfft, noverlap=nfft*0.95)
dbspec = 10.0*np.log10(spec) # in dB
power = dbspec[:, 50]
fish_p = power[(freqs > 400) & (freqs < 1000)]
fish_f = freqs[(freqs > 400) & (freqs < 1000)]
index = np.argmax(fish_p)
eodf = fish_f[index]
eodf4 = eodf * 4
print(eodf4)
lim0 = eodf4-20
lim1 = eodf4+20
#plt.imshow(dbspec, cmap='jet', origin='lower', aspect='auto', vmin=-80, vmax=-10, extent=(times[0], times[-1], freqs[0], freqs[-1]))
# control of plt.imshow
df = freqs[1] - freqs[0]
ix = int(np.round(eodf4/df))
ix0 = int(np.round(lim0/df))
ix1 = int(np.round(lim1/df))
spec4 = dbspec[ix0:ix1, :]
freq4 = freqs[ix0:ix1]
jar4 = freq4[np.argmax(spec4, axis=0)]
jar = jar4 / 4
#np.save
plt.plot(times[:4862], jar4, '-k')
#plt.plot(times, jar)
plt.imshow(spec4, cmap='jet', origin='lower', extent=(times[0], times[-1], freq4[0], freq4[-1]), aspect='auto',
vmin=-80, vmax=-10)
#plt.plot(freqs, dbspec[:,10], label = '0')
#plt.legend()
plt.ylim(lim0, lim1)
plt.show()