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

05.10

Browse Source
This commit is contained in:
xaver 2020-10-05 15:25:49 +02:00
parent 90d9b19d9a
commit d231c75806
7 changed files with 112 additions and 188 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
apteronotus_code/figure_apteronotus_gaincurve_cutofff_tau.py
View File

@ -7,7 +7,7 @@ from matplotlib.mlab import specgram
import os
from jar_functions import gain_curve_fit
identifier = ['2018lepto1', '2018lepto4', '2018lepto5', '2018lepto76']
identifier = ['2020lepto06']
tau = []
f_c = []
@ -42,7 +42,7 @@ for ID in identifier:
ax.set_ylabel('gain [Hz/(mV/cm)]')
ax.set_xlabel('envelope_frequency [Hz]')
ax.set_title('gaincurve %s' %ID)
plt.legend()
plt.legend(loc = 'lower left')
plt.show()

111
apteronotus_code/figure_apteronotus_jar_filter_fit.py
View File

@ -17,7 +17,7 @@ from jar_functions import gain_curve_fit
def take_second(elem): # function for taking the names out of files
return elem[1]
identifier = ['2018lepto1']
identifier = ['2018lepto98']
for ident in identifier:
predict = []
@ -35,16 +35,16 @@ for ident in identifier:
currf = None
idxlist = []
data = sorted(np.load('%s files.npy' %ident), key = take_second) # list with filenames in it
data = sorted(np.load('5Hz_%s files.npy' %ident), key = take_second) # list with filenames in it
for i, d in enumerate(data):
dd = list(d)
if dd[1] == '0.005':
jar = np.load('%s.npy' %dd) # load data for every file name
if dd[1] == '0.5':
jar = np.load('5Hz_%s.npy' %dd) # load data for every file name
jm = jar - np.mean(jar) # low-pass filtering by subtracting mean
print(dd)
time = np.load('%s time.npy' %dd) # time file
time = np.load('5Hz_%s time.npy' %dd) # time file
dt = time[1] - time[0]
n = int(1/float(d[1])/dt)
@ -54,10 +54,11 @@ for ident in identifier:
sinv, sinc = curve_fit(sin_response, time, jm - cutf, [float(d[1]), 2, 0.5]) # fitting
print('frequency, phaseshift, amplitude:', sinv)
p = sinv[1]
A = np.sqrt(sinv[2] ** 2)
f = float(d[1])
if sinv[2] < 0:
A = sinv[2]
if A < 0:
p = p + np.pi
A = -A
f = float(d[1])
phaseshift.append(p)
gain.append(A)
if f not in amfreq:
@ -66,7 +67,7 @@ for ident in identifier:
# jar trace
plt.plot(time, jar, color = 'C0')
#plt.hlines(y=np.min(jar) - 2, xmin=0, xmax=400, lw=2.5, color='r', label='stimulus duration')
plt.title('JAR trace 2018lepto1, AM-frequency:%sHz' % float(d[1]))
plt.title('JAR trace 2018lepto98, AM-frequency:%sHz, deltaf = -5Hz' % float(d[1]))
plt.xlabel('time[s]')
plt.ylabel('frequency[Hz]')
plt.show()
@ -81,7 +82,7 @@ for ident in identifier:
# filter by running average
plt.plot(time, jm, color = 'C0', label = 'JAR: subtracted by mean')
plt.plot(time, jm - cutf, color = 'darkorange', label = 'JAR: subtracted by mean and step response')
plt.title('JAR trace spectogram 2018lepto1: subtraction of mean and step response')
plt.title('JAR trace spectogram 2018lepto98: subtraction of mean and step response, deltaf = -5Hz')
plt.xlabel('time[s]')
plt.ylabel('frequency[Hz]')
plt.legend()
@ -89,96 +90,12 @@ for ident in identifier:
# jar trace and fit
plt.plot(time, jm - cutf, color = 'darkorange', label = 'JAR: subtracted by mean and step response')
phase_gain = [(((sinv[1] % (2 * np.pi)) * 360) / (2 * np.pi)), sinv[2]]
phase_gain = [(((p % (2 * np.pi)) * 360) / (2 * np.pi)), A]
plt.plot(time, sin_response(time, *sinv), color = 'limegreen', label='fit: phaseshift=%.2f°, gain=%.2f[Hz/(mV/cm)]' % tuple(phase_gain))
plt.title('JAR trace spectogram 2018lepto1 with fit')
plt.title('JAR trace spectogram 2018lepto98 with fit, deltaf = -5Hz')
plt.xlabel('time[s]')
plt.ylabel('frequency[Hz]')
plt.legend()
plt.show()
# root mean square
RMS = np.sqrt(np.mean(((jm - cutf) - sin_response(cutt, sinv[0], sinv[1], sinv[2]))**2))
thresh = A / np.sqrt(2)
# mean over same amfreqs for phase and gain
if currf is None or currf == d[1]:
currf = d[1]
idxlist.append(i)
else: # currf != f
meanf = [] # lists to make mean of
meanp = []
meanrms = []
meanthresh = []
for x in idxlist:
meanf.append(gain[x])
meanp.append(phaseshift[x])
meanrms.append(RMS)
meanthresh.append(thresh)
meanedf = np.mean(meanf)
meanedp = np.mean(meanp)
meanedrms = np.mean(meanrms)
meanedthresh = np.mean(meanthresh)
mgain.append(meanedf)
mphaseshift.append(meanedp)
rootmeansquare.append(meanedrms)
threshold.append(meanedthresh)
currf = d[1] # set back for next loop
idxlist = [i]
meanf = []
meanp = []
meanrms = []
meanthresh = []
for y in idxlist:
meanf.append(gain[y])
meanp.append(phaseshift[y])
meanrms.append(RMS)
meanthresh.append(thresh)
meanedf = np.mean(meanf)
meanedp = np.mean(meanp)
meanedrms = np.mean(meanrms)
meanedthresh = np.mean(meanthresh)
mgain.append(meanedf)
mphaseshift.append(meanedp)
rootmeansquare.append(meanedrms)
threshold.append(meanedthresh)
# as arrays
mgain_arr = np.array(mgain)
mphaseshift_arr = np.array(mphaseshift)
amfreq_arr = np.array(amfreq)
rootmeansquare_arr = np.array(rootmeansquare)
threshold_arr = np.array(threshold)
# condition needed to be fulfilled: RMS < threshold or RMS < mean(RMS)
idx_arr = (rootmeansquare_arr < threshold_arr) | (rootmeansquare_arr < np.mean(rootmeansquare_arr))
fig = plt.figure()
ax0 = fig.add_subplot(2, 1, 1)
ax0.plot(amfreq_arr[idx_arr], mgain_arr[idx_arr], 'o')
ax0.set_yscale('log')
ax0.set_xscale('log')
ax0.set_title('%s' % data[0][0])
ax0.set_ylabel('gain [Hz/(mV/cm)]')
ax0.set_xlabel('envelope_frequency [Hz]')
#plt.savefig('%s gain' % data[0][0])
ax1 = fig.add_subplot(2, 1, 2, sharex = ax0)
ax1.plot(amfreq, threshold, 'o-', label = 'threshold', color = 'b')
ax1.set_xscale('log')
ax1.plot(amfreq, rootmeansquare, 'o-', label = 'RMS', color ='orange')
ax1.set_xscale('log')
ax1.set_xlabel('envelope_frequency [Hz]')
ax1.set_ylabel('RMS [Hz]')
plt.legend()
pylab.show()
#np.save('phaseshift_%s' % ident, mphaseshift_arr[idx_arr])
#np.save('gain_%s' %ident, mgain_arr[idx_arr])
#np.save('amf_%s' %ident, amfreq_arr[idx_arr])
embed()
embed()

10
apteronotus_code/figure_apteronotus_rms_gaincurve.py
View File

@ -17,7 +17,7 @@ from jar_functions import gain_curve_fit
def take_second(elem): # function for taking the names out of files
return elem[1]
identifier = ['2018lepto1']
identifier = ['2019lepto03']
for ident in identifier:
predict = []
@ -35,16 +35,16 @@ for ident in identifier:
currf = None
idxlist = []
data = sorted(np.load('%s files.npy' %ident), key = take_second) # list with filenames in it
data = sorted(np.load('5Hz_%s files.npy' %ident), key = take_second) # list with filenames in it
for i, d in enumerate(data):
dd = list(d)
jar = np.load('%s.npy' %dd) # load data for every file name
jar = np.load('5Hz_%s.npy' %dd) # load data for every file name
jm = jar - np.mean(jar) # low-pass filtering by subtracting mean
print(dd)
time = np.load('%s time.npy' %dd) # time file
time = np.load('5Hz_%s time.npy' %dd) # time file
dt = time[1] - time[0]
n = int(1/float(d[1])/dt)
@ -127,7 +127,7 @@ for ident in identifier:
ax0.plot(amfreq_arr[idx_arr], mgain_arr[idx_arr], 'o')
ax0.set_yscale('log')
ax0.set_xscale('log')
ax0.set_title('gaincurve 2018lepto1')
ax0.set_title('gaincurve 2019lepto03, deltaf = -5Hz')
ax0.set_ylabel('gain [Hz/(mV/cm)]')
ax0.set_xlabel('envelope_frequency [Hz]')
#plt.savefig('%s gain' % data[0][0])

40
apteronotus_code/sin_all_normal.py
View File

@ -7,17 +7,17 @@ from jar_functions import gain_curve_fit
from jar_functions import avgNestedLists
identifier = ['2018lepto1',
'2018lepto4',
'2018lepto5',
'2018lepto76',
identifier = [#'2018lepto1',
#'2018lepto4',
#'2018lepto5',
#'2018lepto76',
'2018lepto98',
'2019lepto03',
'2019lepto24',
'2019lepto27',
'2019lepto30',
'2020lepto04',
'2020lepto06',
#'2019lepto03',
#'2019lepto24',
#'2019lepto27',
#'2019lepto30',
#'2020lepto04',
#'2020lepto06',
'2020lepto16',
'2020lepto19',
'2020lepto20'
@ -28,14 +28,14 @@ amf = [0.001, 0.002, 0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1]
all = []
for ident in identifier:
data = np.load('gain_%s.npy' %ident)
data = np.load('5Hz_gain_%s.npy' %ident)
all.append(data)
av = avgNestedLists(all)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(amf, av, 'o')
ax.plot(amf, av, 'o', c = 'C0', label = 'gain')
#plt.show()
@ -45,8 +45,8 @@ fit = []
fit_amf = []
for ID in identifier:
print(ID)
amf = np.load('amf_%s.npy' %ID)
gain = np.load('gain_%s.npy' %ID)
amf = np.load('5Hz_amf_%s.npy' %ID)
gain = np.load('5Hz_gain_%s.npy' %ID)
sinv, sinc = curve_fit(gain_curve_fit, amf, gain)
#print('tau:', sinv[0])
@ -56,17 +56,19 @@ for ID in identifier:
f_c.append(f_cutoff)
fit.append(gain_curve_fit(amf, *sinv))
fit_amf.append(amf)
col = plt.cm.magma(np.linspace(0,0.8,len(fit)))
for ff ,f in enumerate(fit):
ax.plot(fit_amf[ff], fit[ff], c = col[ff])
ax.axvline(x=f_c[ff], ymin=0, ymax=5, alpha=0.8, c = col[ff]) # colors_uniform[ff])
#for ff ,f in enumerate(fit):
# ax.plot(fit_amf[ff], fit[ff])
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_title('gaincurve_average_allfish')
ax.set_title('gain average all fish, deltaf: -5Hz')
ax.set_ylabel('gain [Hz/(mV/cm)]')
ax.set_xlabel('envelope_frequency [Hz]')
ax.set_ylim(0.0008, )
ax.plot(f_c, np.full(len(identifier), 0.0015), 'o', label = 'cutoff frequencies')
ax.legend()
#ax.plot(f_c, np.full(len(identifier), 0.0015), 'o', alpha = 0.5, c = 'darkorange', label = 'cutoff frequencies')
ax.legend(loc = 'center left')
plt.show()

63
apteronotus_code/sin_all_uniform.py
View File

@ -9,19 +9,19 @@ import matplotlib as mpl
from matplotlib import cm
identifier_uniform = ['2018lepto1',
# '2018lepto4',
# '2018lepto5',
#'2018lepto76',
#'2018lepto4',
'2018lepto5',
'2018lepto76',
'2018lepto98',
# '2019lepto03',
#'2019lepto03',
'2019lepto24',
#'2019lepto27',
# '2019lepto30',
'2020lepto04',
# '2020lepto06',
# '2020lepto16',
'2020lepto19',
# '2020lepto20'
#'2019lepto30',
#'2020lepto04',
'2020lepto06',
#'2020lepto16',
#'2020lepto19',
#'2020lepto20'
]
identifier = ['2018lepto1',
'2018lepto4',
@ -32,10 +32,10 @@ identifier = ['2018lepto1',
'2019lepto24',
'2019lepto27',
'2019lepto30',
'2020lepto04',
#'2020lepto04',
'2020lepto06',
'2020lepto16',
'2020lepto19',
#'2020lepto19',
'2020lepto20'
]
@ -56,15 +56,15 @@ new_av = avgNestedLists(new_all)
fig = plt.figure()
ax = fig.add_subplot(111)
#ax.plot(amf, av, 'o', color = 'orange', label = 'normal')
ax.plot(amf, new_av, 'o', label = 'uniformed')
"""
ax.plot(amf, av, 'o', color = 'darkorange', label = 'normal')
ax.plot(amf, new_av, 'o', label = 'uniform')
tau = []
f_c = []
fit = []
fit_amf = []
for ID in identifier:
#print(ID)
print(ID)
amf = np.load('amf_%s.npy' %ID)
gain = np.load('gain_%s.npy' %ID)
@ -72,47 +72,46 @@ for ID in identifier:
#print('tau:', sinv[0])
tau.append(sinv[0])
f_cutoff = abs(1 / (2*np.pi*sinv[0]))
#print('f_cutoff:', f_cutoff)
print('f_cutoff:', f_cutoff)
f_c.append(f_cutoff)
fit.append(gain_curve_fit(amf, *sinv))
fit_amf.append(amf)
"""
# uniformed: 2018lepto5, 2018lepto1, 2018lepto76, 2018lepto98, 2020lepto06, 2019lepto24, 2020lepto4
tau_uniform = []
f_c_uniform = []
fit_uniform = []
fit_amf_uniform = []
for ID in identifier_uniform:
#print(ID)
print(ID)
amf = np.load('amf_%s.npy' %ID)
gain = np.load('gain_%s.npy' %ID)
sinv, sinc = curve_fit(gain_curve_fit, amf, gain)
#print('tau:', sinv[0])
print('tau:', sinv[0])
tau_uniform.append(sinv[0])
f_cutoff = abs(1 / (2*np.pi*sinv[0]))
#print('f_cutoff:', f_cutoff)
print('f_cutoff:', f_cutoff)
f_c_uniform.append(f_cutoff)
fit_uniform.append(gain_curve_fit(amf, *sinv))
fit_amf_uniform.append(amf)
colors_uniform = plt.cm.flag(np.linspace(0.2,0.8,len(fit_uniform)))
#colors = plt.cm.flag(np.linspace(0.2,0.8,len(fit)))
# for ff ,f in enumerate(fit):
# ax.plot(fit_amf[ff], fit[ff],color = colors[ff])
# ax.axvline(x=f_c[ff], ymin=0, ymax=5, ls = '-', alpha = 0.5, color= colors[ff])#colors_uniform[ff])
colors = plt.cm.flag(np.linspace(0,1,len(fit_uniform)))
for ff, f in enumerate(fit_uniform):
ax.plot(fit_amf_uniform[ff], fit_uniform[ff], color = colorss[ff]) #colors_uniform[ff])
ax.axvline(x=f_c_uniform[ff], ymin=0, ymax=5, ls = '-', alpha = 0.5, color= colorss[ff])#colors_uniform[ff])
#for ff ,f in enumerate(fit_uniform):
# ax.plot(fit_amf_uniform[ff], fit_uniform[ff],color = colorss[ff])
# ax.axvline(x=f_c_uniform[ff], ymin=0, ymax=5, ls = '-', alpha = 0.5, color= colorss[ff])#colors_uniform[ff])
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_title('gaincurve_average_allfish')
ax.set_title('gain average all fish uniform')
ax.set_ylabel('gain [Hz/(mV/cm)]')
ax.set_xlabel('envelope_frequency [Hz]')
ax.set_ylim(0.0008, )
ax.legend()
ax.plot(f_c, np.full(len(identifier), 0.0015), 'o', alpha = 0.5, c = 'darkorange', label = 'normal cutoff frequencies')
ax.plot(f_c_uniform, np.full(len(identifier_uniform), 0.002), 'o', alpha = 0.5, c = 'C0', label = 'uniform cutoff frequencies')
ax.legend(loc = 'center left')
plt.show()
embed()

26
eigenmannia_code/eigenmannia_jar_stacked.py
View File

@ -14,7 +14,8 @@ from scipy.signal import savgol_filter
base_path = 'D:\\jar_project\\JAR\\eigenmannia\\deltaf'
identifier = ['2013eigen13','2015eigen16', '2015eigen17', '2015eigen19', '2020eigen22','2020eigen32']
#2015eigen8 no nix files
identifier = ['2015eigen16', '2013eigen13','2015eigen17', '2015eigen19', '2020eigen22','2020eigen32']
response = []
deltaf = []
@ -47,8 +48,8 @@ for ID in identifier:
eodf = fish_f[index]
eodf4 = eodf * 4
lim0 = eodf4 - 50
lim1 = eodf4 + 50
lim0 = eodf4 - 40
lim1 = eodf4 + 40
df = freqs[1] - freqs[0]
ix0 = int(np.floor(lim0/df)) # back to index
@ -58,10 +59,17 @@ for ID in identifier:
jar4 = freq4[np.argmax(spec4, axis=0)] # all freqs at max specs over axis 0
cut_time_jar = times[:len(jar4)]
#plt.imshow(spec4, cmap='jet', origin='lower', extent=(times[0], times[-1], lim0, lim1), aspect='auto', vmin=-80, vmax=-10)
#plt.plot(cut_time_jar, jar4)
ID_delta_f = [ID, str(delta_f[0]).split('.')[0]]
plt.imshow(spec4, cmap='jet', origin='lower', extent=(times[0] - 10, times[-1] - 10, lim0, lim1), aspect='auto', vmin=-80, vmax=-10)
plt.plot((cut_time_jar - 10), jar4, 'k', label = 'jar trace', lw = 2)
plt.hlines(y=lim0 + 5, xmin=0, xmax=60, lw=2.5, color='gold', label='stimulus duration')
plt.title('spectogram %s, deltaf: %sHz' %tuple(ID_delta_f))
plt.xlim(right=times[-1] - 10)
plt.legend()
#plt.show()
delta_f_ID = [str(delta_f[0]).split('.')[0], ID]
plt.savefig('%sHz_specgram_jar_%s' %tuple(delta_f_ID))
plt.close()
b = []
for idx, i in enumerate(times):
@ -79,8 +87,4 @@ for ID in identifier:
res_df = sorted(zip(deltaf,response))
np.save('res_df_%s_new' %ID, res_df)
# problem: rohdaten(data, pre_data) lassen sich auf grund ihrer 1D-array struktur nicht savgol filtern
# diese bekomm ich nur über specgram in form von freq / time auftragen, was nicht mehr savgol gefiltert werden kann
# jedoch könnte ich trotzdem einfach aus jar4 response herauslesen wobei dies dann weniger gefiltert wäre
#np.save('res_df_%s_new' %ID, res_df)

46
notes
View File

@ -1,10 +1,19 @@
machen:
- phaseshift (sin_all) nochmal, nicht richtiges Dings verwendet ( sinv[2]/p) und dann auch in phaseshift
- abbildung erstellen mit custom cutoff frequencies über ganzen bereich (0.001Hz-1Hz) um hoffentlich zu zeigen dass dabei lineare
Gerade entsteht, vergleichen mit uniformen Bereich bei Daten bei dem es sich auch eher linear verhält um zu zeigen auf was wir hinaus wollen
- filter: zieht mean von einer amfreq periode ab wodurch alles was nicht damit in Verbindung steht herausfiltert,
auch JAR. problematisch wird dies eher wenn JAR-Anstieg schneller abläuft als eine amfreq periode
- wenn fit nicht funktioniert einfach weglassen, wenn sättigung vorhanden nochmal anschauen
+ figures:
apteronotus: fundament by tims bachelor thesis, important that apteronotus only shifts his frequency up (as eigenmannia doesnt --> natalies measurements)
+ spectogram
+ jar trace out of specgram
+ filtering of jar trace: mean noise cut --> subtracting jar response over whole stimulus
+ fit and jar trace --> gain and phaseshift
!!! + this for different am-frequencies and delta f (-15/-5Hz) --> compare gain for them
+ this for different am-frequencies and delta f (-15/-5Hz) --> compare gain for them
+ gain curve for one or more single fish
+ fit of gain curve for cutoff frequency and tau
+ gain curve for all fish taken together
@ -15,37 +24,30 @@
fig_apt_gaincurve_cutoff_tau,
sin_all_normal (without single gaincurves),
sin_all_uniform (with gaincurves for 5Hz)
!!! discard 2019lepto03 -5Hz --> RMS always over threshold
!!! discard 2020lepto04 gain fit doesnt work
!!! discard 2020lepto19 gain fit doesnt work
eigenmannia:
+ deltaf / response: -2Hz different, show it
+ spectogram
+ direct to fit and jar trace --> gain and phaseshift DURCH SIN RESPONSE SPEC JAGEN!
+ direct to fit and jar trace
+ gain curve for one or more single fish
+ gain curve for all fish taken together
+ gain curve for all fish taken together !!! gains without filtering by RMS
- (step response eigen)
--> fig_eig_specgram
fig_eig_jar_filter_fit,
fig_eig_rms_gaincurve,
fig_eig_gaincurve_cutoff_tau,
sin_eig_normal (without single gaincurves),
sin_eig_uniform (with gaincurves for 5Hz)
!!! discard 2015eigen8 --> RMS always over threshold
??? mean_noise_cut_eigen?
fish properties:
+ parameters
+ cutoff frequency - dominance score
+ eigenmannia deltaf response over all fish mean
+ phaseshift_all: wenn negativer gain in fit --> +pi rechnen, dann modulo
- plot_eigenmannia_jar(compare res_df_%s / res_df_%s_new)
- eigenmannia_jar:
- specgram auch zeigen, vorallem was auch die ausreißer bei -2 Hz betreffen
- fish_properties:
- hauptsächlich auf f_c und tau konzentrieren, vor allem auch beides auftragen, gewicht/größe noch nehmen
- step_response eigen: absolute response
- Q10 Wert aus Formel von Jan auf base_frequenz rechnen (adjust-eodf in jar_functions)
- sin_all_uniform - sin_all_normal (also 5Hz, let away 0.001Hz?, gain_fit): fit als spur reinlegen damit klar wird aus was gerade besteht
long term:
- extra datei mit script drin um fertige daten darzustellen, den fit-code nur zur datenverarbeitung verwenden
- darstellung: specgram --> rausgezogene jarspur darüber --> filterung --> fit und daten zusammen dargestellt, das ganze für verschiedene frequenzen
- unterschiedliche nffts auf anderem rechner laufen lassen evtl um unterschiede zu sehen
- phase in degree: phase % (2pi) - modulo 2pi
offen: phaseshift, not working fit for some fish
(
- mit zu hohem RMS rauskicken: evtl nur ein trace rauskicken wenn nur da RMS zu hoch
- 2019lepto27/30: 27 - 0.05Hz (7-27-af, erste dat mit len(dat)=1), 30 - 0.001Hz (7-30-ah mit 0.005 anstatt gewollten 0.001Hz --> fehlt)
)