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

29.07

Browse Source
This commit is contained in:
xaver 2020-07-29 18:54:50 +02:00
parent 0d550e2a5d
commit 0d98a7842e
4 changed files with 102 additions and 189 deletions

140
jar_functions.py
View File

@ -171,117 +171,29 @@ def average(freq_all, time_all, start, stop, timespan, dm):
return mf_all, tnew_all, values_all
def iload_traces(basedir, repro='', before=0.0, after=0.0):
"""
returns:
info : metadata from stimuli.dat
key : key from stimuli.dat
time : an array for the time axis
data : the data of all traces of a single trial
"""
p = re.compile('([-+]?\d*\.\d+|\d+)\s*(\w+)')
reproid = 'RePro'
deltat = None
# open traces files:
sf = []
for trace in range(1, 1000000):
if path.isfile('%s/trace-%i.raw' % (basedir, trace)):
sf.append(open('%s/trace-%i.raw' % (basedir, trace), 'rb'))
else:
break
basecols = None
baserp = True
for info, key, dat in iload('%s\\stimuli.dat' % (basedir,)):
if deltat is None:
deltat, tunit = p.match(info[0]['sample interval%i' % 1]).groups()
deltat = float(deltat)
if tunit == 'ms':
deltat *= 0.001
if 'repro' in info[-1] or 'RePro' in info[-1]:
if not reproid in info[-1]:
reproid = 'repro'
if len(repro) > 0 and repro != info[-1][reproid] and \
basecols is None:
continue
baserp = (info[-1][reproid] == 'BaselineActivity')
duration_indices = [i for i, x in enumerate(key[2]) if x == "duration"]
else:
baserp = True
duration_indices = []
if dat.shape == (1, 1) and dat[0, 0] == 0:
warnings.warn("iload_traces: Encountered incomplete '-0' trial.")
yield info, key, array([])
continue
if baserp:
basecols = []
basekey = key
baseinfo = info
if len(dat) == 0:
for trace in range(len(sf)):
basecols.append(0)
for d in dat:
if not baserp and not basecols is None:
x = []
xl = []
for trace in range(len(sf)):
col = int(d[trace])
sf[trace].seek(basecols[trace] * 4)
buffer = sf[trace].read((col - basecols[trace]) * 4)
tmp = fromstring(buffer, float32)
x.append(tmp)
xl.append(len(tmp))
ml = min(xl)
for k in range(len(x)):
if len(x[k]) > ml:
warnings.warn("trunkated trace %d to %d" % (k, ml))
x[k] = x[k][:ml]
xtime = arange(0.0, len(x[0])) * deltat - before
yield baseinfo, basekey, xtime, asarray(x)
basecols = None
if len(repro) > 0 and repro != info[-1][reproid]:
break
durations = [d[i] for i in duration_indices if not isnan(d[i])]
if not baserp:
if len(durations) > 0:
duration = max(durations)
if duration < 0.001: # if the duration is less than 1ms
warnings.warn(
"iload_traces: Skipping one trial because its duration is <1ms and therefore it is probably rubbish")
continue
l = int(before / deltat)
r = int((duration + after) / deltat)
else:
continue
x = []
xl = []
for trace in range(len(sf)):
col = int(d[trace])
if baserp:
if col < 0:
col = 0
basecols.append(col)
continue
sf[trace].seek((col - l) * 4)
buffer = sf[trace].read((l + r) * 4)
tmp = fromstring(buffer, float32)
x.append(tmp)
xl.append(len(tmp))
if baserp:
break
ml = min(xl)
for k in range(len(x)):
if len(x[k]) > ml:
warnings.warn("trunkated trace %d to %d" % (k, ml))
x[k] = x[k][:ml]
time = arange(0.0, len(x[0])) * deltat - before
yield info, key, time, asarray(x)
for trace in range(len(sf)):
sf[trace].close()
def import_data(dataset):
print(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']
dat = []
pre_dat = []
for mt in b.multi_tags:
data = mt.retrieve_data(0, 'EOD-1')[:] # data[0]
dat.append(data)
i0 = int(mt.positions[0][0]/dt)
pre_data = eod[i0-di:i0]
pre_dat.append(pre_data)
return dat, pre_dat, dt
#nf.close()
if __name__ == '__main__':
import_data(os.path.join('JAR', '2020-07-21-ak', '2020-07-21-ak.nix'))

BIN
jar_spectogram.png Normal file
View File

Binary file not shown.
Side by Side

After

(image error) Size: 308 KiB

65
scratch.py
View File

@ -2,47 +2,42 @@ import os
import numpy as np
from IPython import embed
from jar_functions import mean_noise_cut
from matplotlib.mlab import specgram
import DataLoader as dl
datasets = [(os.path.join('D:\\jar_project\\JAR\\2020-06-22-ac\\beats-eod.dat'))]
"""
#second_try scratch
minimum = min(len(frequency[0]), len(frequency[1]))
f1 = frequency[0][:minimum]
f2 = frequency[1][:minimum]
frequency = np.array(frequency)
mean = np.mean(frequency, axis=0)
for i in range(len(minimum)):
mean(frequency[0][i], frequency[1][i])
for f in frequency:
print(np.mean(f))
#print(np.logspace(-3, 1, 10))
'''for idx, dataset in enumerate(datasets):
datapath = os.path.join(base_path, dataset)
for info, key, time, data in dl.iload_traces(datapath, repro='Beats', before=0.0, after=0.0):
'''
# mean_f = np.mean(x) for x in zip(freqeuncies1, frequencies2)
"""
base_path = 'D:\\jar_project\\JAR'
'''
g = [1, 2]
#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
]
for idx, dataset in enumerate(datasets):
datapath = os.path.join(base_path, dataset)
for info, key, time, data in dl.iload_traces(datapath, repro='Beats', before=0.0, after=0.0):
print(data[0])
h = [3, 4]
dat = np.arange(100)
z = np.array([[g], [h]])
for d in range(int(len(data)/10)):
nfft = 2
mean0 = np.mean(z, axis=0)
mean1 = np.mean(z, axis=1)
spec, freqs, times = specgram(data[0][d*10:(d+1)*10], NFFT=nfft, noverlap=nfft*0.5)
print(mean0)
print(mean1)
''''''
t = [600, 650]
#print(freqs)
#print(times)
embed()
x = 1 - np.exp(t / 11)
print(x)
a, b ,c,d = normalized_JAR(fre)
'''
print(np.logspace(-3, 1, 10))
embed()

86
step_response_specto.py
View File

@ -17,6 +17,7 @@ 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'
@ -50,44 +51,49 @@ for infodataset in datasets:
for idx, dataset in enumerate(datasets):
datapath = os.path.join(base_path, dataset)
for info, key, time, data in dl.iload_traces(datapath, repro='Beats', before=0.0, after=0.0):
#plt.plot(time, data[0]) # V(t)
#plt.show()
nfft = 2**18
spec, freqs, times = specgram(data[0], Fs=1.0/(time[1]-time[0]), detrend='mean', NFFT=nfft, noverlap=nfft*0.8)
dbspec = 10.0*np.log10(spec) # in dB
power = dbspec[:, 10]
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', extent=(times[0], times[-1], freqs[0], freqs[-1]), aspect='auto', vmin=-80, vmax=-30)
# 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
plt.plot(times, jar4, '-k')
#plt.plot(freqs, dbspec[:,100], label = '0')
#plt.legend()
plt.ylim(lim0, lim1)
plt.show()
embed()
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()