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

14.11

Browse Source
This commit is contained in:
xaver
2020-11-14 19:03:48 +01:00
parent ebf365dbe0
commit 94a1697019
13 changed files with 815 additions and 211 deletions
Download Patch File Download Diff File Expand all files Collapse all files

22
eigenmannia_code/eigenmannia_jar_stacked.py
View File

@@ -12,12 +12,12 @@ from jar_functions import get_time_zeros
from jar_functions import import_data_eigen
from scipy.signal import savgol_filter
plt.rcParams.update({'font.size': 18})
plt.rcParams.update({'font.size': 12})
base_path = 'D:\\jar_project\\JAR\\eigenmannia\\deltaf'
#2015eigen8 no nix files
identifier = [#'2013eigen13',
identifier = ['2013eigen13',
'2015eigen16','2015eigen17', '2015eigen19', '2020eigen22','2020eigen32']
@@ -32,8 +32,8 @@ for ID in identifier:
delta_f, duration = parse_stimuli_dat(stimuli_dat)
dur = int(duration[0][0:2])
print(delta_f)
if delta_f != [4.0]:
continue
#if delta_f != [-4.0]:
# continue
data, pre_data, dt = import_data_eigen(datapath)
#hstack concatenate: 'glue' pre_data and data
@@ -74,24 +74,24 @@ for ID in identifier:
if i > 15 and i < 55:
j.append(jar4[idx])
r = np.median(j) - np.median(b)
r = (np.median(j) - np.median(b)) / 4 # divided by 4 cause of data at 4th harmonic, therefore response 4 times higher
print('response:', r)
deltaf.append(delta_f[0])
response.append(r)
plt.figure(figsize = (14,8))
plt.figure(figsize = (8.27,11.69/2))
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, 'k', label = 'peak detection trace', lw = 2)
plt.hlines(y=lim0 + 5, xmin=10, xmax=70, lw=4, color='yellow', label='stimulus duration')
plt.plot(cut_time_jar, jar4, color = 'k', label = 'peak detection trace', lw = 2)
plt.hlines(y=lim0 + 5, xmin=0, xmax=10, lw=4, color='red', label='pause')
plt.hlines(y=lim0 + 5, xmin=10, xmax=70, lw=4, color='gold', label='stimulus duration')
plt.title('spectogram %s, deltaf: %sHz' %tuple(ID_delta_f))
plt.xlim(times[0],times[-1])
#embed()
#plt.xticks((times[0], 10, 20, 30, 40, 50, 60, times[-1]), [0, 10, 20, 30 ,40, 50, 60, 70])
plt.xticks((times[0], 10, 20, 30, 40, 50, 60, times[-1]), [0, 10, 20, 30 ,40, 50, 60, 70])
plt.xlabel('time [s]')
plt.ylabel('frequency [Hz]')
plt.legend(loc = 'best')
plt.show()
#plt.show()
delta_f_ID = [str(delta_f[0]).split('.')[0], ID]
plt.close()
@@ -99,4 +99,4 @@ for ID in identifier:
res_df = sorted(zip(deltaf,response))
#np.save('res_df_%s_new' %ID, res_df)
np.save('res_df_%s_new' %ID, res_df)

7
eigenmannia_code/eigenmannia_jar_subplot.py
View File

@@ -12,7 +12,7 @@ from jar_functions import get_time_zeros
from jar_functions import import_data_eigen
from scipy.signal import savgol_filter
#plt.rcParams.update({'font.size': 18})
plt.rcParams.update({'font.size': 10})
base_path = 'D:\\jar_project\\JAR\\eigenmannia\\deltaf'
@@ -107,6 +107,7 @@ ax0.set_xlim(times[0],times[-1])
ax0.set_ylabel('frequency [Hz]')
ax0.axes.xaxis.set_ticklabels([])
ax0.set_title('∆F -2 Hz')
plt.xticks((1.7, 10, 20, 30, 40, 50, 60, times[-1]))
ax1 = fig.add_subplot(222)
ax1.imshow(specs[2], cmap='jet', origin='lower', extent=(times[0], times[-1], sub_lim0[2], sub_lim1[2]), aspect='auto', vmin=-80, vmax=-10)
@@ -116,6 +117,7 @@ ax1.axes.xaxis.set_ticklabels([])
#ax1.axes.yaxis.set_ticklabels([])
ax1.set_title('∆F 2 Hz')
ax1.get_shared_y_axes().join(ax0, ax1)
plt.xticks((1.7, 10, 20, 30, 40, 50, 60, times[-1]))
ax2 = fig.add_subplot(223)
ax2.imshow(specs[1], cmap='jet', origin='lower', extent=(times[0], times[-1], sub_lim0[1], sub_lim1[1]), aspect='auto', vmin=-80, vmax=-10)
@@ -124,6 +126,7 @@ ax2.set_xlim(times[0],times[-1])
ax2.set_ylabel('frequency [Hz]')
ax2.set_xlabel('time [s]')
ax2.set_title('∆F -10 Hz')
plt.xticks((1.7, 10, 20, 30, 40, 50, 60, times[-1]), [0, 10, 20, 30 ,40, 50, 60, 70])
ax3 = fig.add_subplot(224)
ax3.imshow(specs[3], cmap='jet', origin='lower', extent=(times[0], times[-1], sub_lim0[3], sub_lim1[3]), aspect='auto', vmin=-80, vmax=-10)
@@ -132,7 +135,7 @@ ax3.set_xlim(times[0],times[-1])
ax3.set_xlabel('time [s]')
#ax3.axes.yaxis.set_ticklabels([])
ax3.set_title('∆F 10 Hz')
plt.xticks((1.7, 10, 20, 30, 40, 50, 60, times[-1]), [0, 10, 20, 30 ,40, 50, 60, 70])
plt.subplots(sharex = True, sharey = True)
plt.show()

128
eigenmannia_code/figure_eigen_gain_plot.py Normal file
View File

@@ -0,0 +1,128 @@
import matplotlib.pyplot as plt
import numpy as np
import pylab
from IPython import embed
from scipy.optimize import curve_fit
from matplotlib.mlab import specgram
import os
from jar_functions import gain_curve_fit
plt.rcParams.update({'font.size': 10})
identifier = ['2015eigen8',
'2015eigen15',
'2015eigen16',
'2015eigen17',
'2015eigen19'
]
amfs = []
gains = []
maxgains = []
mingains = []
taus = []
f_cs = []
predicts = []
for ID in identifier:
predict = []
print(ID)
amf = np.load('eigen_amf_%s.npy' % ID)
amfs.append(amf)
gain = np.load('eigen_gain_%s.npy' % ID)
gains.append(gain)
print(np.max(gain))
sinv, sinc = curve_fit(gain_curve_fit, amf, gain, [2, 3])
# print('tau:', sinv[0])
taus.append(sinv[0])
f_cutoff = abs(1 / (2 * np.pi * sinv[0]))
print('f_cutoff:', f_cutoff)
f_cs.append(f_cutoff)
print('min gain:', np.min(gain))
print('max gain:', np.max(gain))
maxgains.append(np.max(gain))
mingains.append(np.min(gain))
# predict of gain
for f in amf:
G = np.max(gain) / np.sqrt(1 + (2 * ((np.pi * f * sinv[0]) ** 2)))
predict.append(G)
predicts.append(predict)
print('max of absolute max gain:', np.max(maxgains))
print('min of absolute max gain:', np.min(maxgains))
print('max of absolute min gain:', np.max(mingains))
print('min of absolute min gain:', np.min(mingains))
print('absolute max f_c:', np.max(f_cs))
print('absolute min f_c:', np.min(f_cs))
sort = sorted(zip(f_cs, identifier))
print(sort)
# order of plotting: 2018lepto1, 2018lepto5, 2018lepto76, 2018lepto98, 2019lepto24, 2020lepto06
# order of f_c: 2019lepto24, 2020lepto06, 2018lepto98, 2018lepto76, 2018lepto1, 2018lepto5
fig = plt.figure(figsize=(8.27, 11.69))
# ax0 = plt.subplot2grid(shape=(3,4), loc=(0,0), colspan = 2)
# ax1 = plt.subplot2grid((3,4), (0,2), colspan = 2)
# ax2 = plt.subplot2grid((3,4), (1,0), colspan = 2)
# ax3 = plt.subplot2grid((3,4), (1,2), colspan = 2)
# ax4 = plt.subplot2grid((3,4), (2,0), colspan = 2)
ax0 = fig.add_subplot(321)
fig.text(0.05, 0.5, 'gain [Hz/(mV/cm)]', ha='center', va='center', rotation='vertical')
fig.text(0.5, 0.04, 'envelope frequency [Hz]', ha='center', va='center')
ax0.set_xlim(0.0007, 1.5)
ax0.set_ylim(0.001, 10)
ax0.plot(amfs[1], gains[1], 'o', label='gain')
ax0.plot(amfs[1], predicts[1], label='fit')
ax0.axvline(x=f_cs[1], ymin=0, ymax=5, ls='-', alpha=0.5, label='cutoff frequency')
ax0.set_xscale('log')
ax0.set_yscale('log')
ax0.axes.xaxis.set_ticklabels([])
ax1 = fig.add_subplot(322)
ax1.set_xlim(0.0007, 1.5)
ax1.set_ylim(0.001, 10)
ax1.plot(amfs[0], gains[0], 'o', label='gain')
ax1.plot(amfs[0], predicts[0], label='fit')
ax1.axvline(x=f_cs[0], ymin=0, ymax=5, ls='-', alpha=0.5, label='cutoff frequency')
ax1.set_xscale('log')
ax1.set_yscale('log')
ax1.axes.yaxis.set_ticklabels([])
ax1.axes.xaxis.set_ticklabels([])
ax2 = fig.add_subplot(323)
ax2.set_xlim(0.0007, 1.5)
ax2.set_ylim(0.001, 10)
ax2.plot(amfs[4], gains[4], 'o', label='gain')
ax2.plot(amfs[4], predicts[4], label='fit')
ax2.axvline(x=f_cs[4], ymin=0, ymax=5, ls='-', alpha=0.5, label='cutoff frequency')
ax2.set_xscale('log')
ax2.set_yscale('log')
ax2.axes.xaxis.set_ticklabels([])
ax3 = fig.add_subplot(324)
ax3.set_xlim(0.0007, 1.5)
ax3.set_ylim(0.001, 10)
ax3.plot(amfs[2], gains[2], 'o', label='gain')
ax3.plot(amfs[2], predicts[2], label='fit')
ax3.axvline(x=f_cs[2], ymin=0, ymax=5, ls='-', alpha=0.5, label='cutoff frequency')
ax3.set_xscale('log')
ax3.set_yscale('log')
ax3.axes.yaxis.set_ticklabels([])
ax4 = fig.add_subplot(325)
ax4.set_xlim(0.0007, 1.5)
ax4.set_ylim(0.001, 10)
ax4.plot(amfs[3], gains[3], 'o', label='gain')
ax4.plot(amfs[3], predicts[3], label='fit')
ax4.axvline(x=f_cs[3], ymin=0, ymax=5, ls='-', alpha=0.5, label='cutoff frequency')
ax4.set_xscale('log')
ax4.set_yscale('log')
# plt.legend(loc = 'lower left')
plt.show()
# np.save('f_c', f_c)
# np.save('tau', tau)

168
eigenmannia_code/figure_eigen_jar_plot.py Normal file
View File

@@ -0,0 +1,168 @@
import matplotlib.pyplot as plt
import numpy as np
import pylab
from IPython import embed
from scipy.optimize import curve_fit
from scipy.optimize import curve_fit
from matplotlib.mlab import specgram
import os
from jar_functions import import_data
from jar_functions import import_amfreq
from jar_functions import sin_response
from jar_functions import mean_noise_cut
from jar_functions import gain_curve_fit
#plt.rcParams.update({'font.size': 10})
def take_second(elem): # function for taking the names out of files
return elem[1]
identifier = ['2015eigen8',
'2015eigen15',
'2015eigen16',
'2015eigen17',
'2015eigen19'
]
for ident in identifier:
times = []
jars = []
jms = []
amfreq = []
times1 = []
jars1 = []
jms1 = []
amfreq1 = []
amf = [0.001, 0.002, 0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1]
data = sorted(np.load('eigen_%s files.npy' %ident), key = take_second) # list with filenames in it
for i, d in enumerate(data):
dd = list(d)
if dd[1] == '1' or dd[1] == '0.2' or dd[1] == '0.05' or dd[1] == '0.01' or dd[1] == '0.005' or dd[1] == '0.001':
jar = np.load('eigen_%s.npy' %dd) # load data for every file name
jm = jar - np.mean(jar) # low-pass filtering by subtracting mean
time = np.load('eigen_%s time.npy' %dd) # time file
dt = time[1] - time[0]
n = int(1/float(d[1])/dt)
cutf = mean_noise_cut(jm, n = n)
cutt = time
if dd[1] == '0.001':
amfreq1.append(dd[1])
jars1.append(jm - cutf)
jms1.append(jm)
times1.append(time)
if dd[1] not in amfreq:
print(dd)
amfreq.append(dd[1])
jars.append(jm - cutf)
jms.append(jm)
times.append(time)
else:
#print('1:', dd)
amfreq1.append(dd[1])
jars1.append(jm - cutf)
jms1.append(jm)
times1.append(time)
if len(jars) != 6:
continue
ssample = 100000
fig = plt.figure(figsize=(8.27, 11.69))
fig.suptitle('%s' % ident)
fig.text(0.06, 0.5, 'fish frequency [Hz]', ha='center', va='center', rotation='vertical', color='C0')
fig.text(0.97, 0.5, 'stimulus amplitude [mV/cm]', ha='center', va='center', rotation='vertical', color='red')
fig.text(0.5, 0.04, 'time [s]', ha='center', va='center')
ax0 = fig.add_subplot(611)
print('absolute frequency shift 0.001Hz:', np.max(jars[0]) - np.min(jars[0]))
ax0.plot(times[0], jars[0], zorder=20)
# ax0.set_zorder(1)
ax0.set_ylim(-12, 12)
lower0 = 0
upper0 = 2000
x0 = np.linspace(lower0, upper0, sample)
y0 = (sin_response(np.linspace(lower0, upper0, sample), 0.001, np.pi / 2, .35) + 0.5)
ax0_0 = ax0.twinx()
ax0_0.set_ylim(-0.2, 1.2)
ax0_0.plot(x0, y0, color='red', zorder=1, alpha=0.5)
# ax0_0.set_zorder(2)
ax1 = fig.add_subplot(612)
print('absolute frequency shift 0.005 Hz:', np.max(jars[1]) - np.min(jars[1]))
ax1.plot(times[1], jars[1])
ax1.set_ylim(-12, 12)
lower1 = 0
upper1 = 400
x1 = np.linspace(lower1, upper1, sample)
y1 = (sin_response(np.linspace(lower1, upper1, sample), 0.005, np.pi / 2, .35) + 0.5)
ax1_0 = ax1.twinx()
ax1_0.set_ylim(-0.2, 1.2)
ax1_0.plot(x1, y1, color='red', alpha=0.5)
ax2 = fig.add_subplot(613)
print('absolute frequency shift 0.01 Hz:', np.max(jars[2]) - np.min(jars[2]))
ax2.plot(times[2], jars[2])
ax2.set_ylim(-12, 12)
lower2 = 0
upper2 = 400
x2 = np.linspace(lower2, upper2, sample)
y2 = (sin_response(np.linspace(lower2, upper2, sample), 0.01, np.pi / 2, 0.35) + 0.5)
ax2_0 = ax2.twinx()
ax2_0.set_ylim(-0.2, 1.2)
ax2_0.plot(x2, y2, color='red', alpha=0.5)
ax3 = fig.add_subplot(614)
print('absolute frequency shift 0.02 Hz:', np.max(jars[3]) - np.min(jars[3]))
ax3.plot(times[3], jars[3])
ax3.set_ylim(-12, 12)
lower3 = 0
upper3 = 200
x3 = np.linspace(lower3, upper3, sample)
y3 = (sin_response(np.linspace(lower3, upper3, sample), 0.05, np.pi / 2, 0.35) + 0.5)
ax3_0 = ax3.twinx()
ax3_0.set_ylim(-0.2, 1.2)
ax3_0.plot(x3, y3, color='red', alpha=0.5)
ax4 = fig.add_subplot(615)
print('absolute frequency shift 0.5 Hz:', np.max(jars[4]) - np.min(jars[4]))
ax4.plot(times[4], jars[4])
ax4.set_ylim(-12, 12)
lower4 = 0
upper4 = 200
x4 = np.linspace(lower4, upper4, sample)
y4 = (sin_response(np.linspace(lower4, upper4, sample), 0.2, np.pi / 2, 0.35) + 0.5)
ax4_0 = ax4.twinx()
ax4_0.set_ylim(-0.2, 1.2)
ax4_0.plot(x4, y4, color='red', alpha=0.5)
ax5 = fig.add_subplot(616)
print('absolute frequency shift 1 Hz:', np.max(jars[5]) - np.min(jars[5]))
ax5.plot(times[5], jars[5])
ax5.set_ylim(-12, 12)
lower5 = 0
upper5 = 200
x5 = np.linspace(lower5, upper5, sample)
y5 = (sin_response(np.linspace(lower5, upper5, sample), 1, np.pi / 2, 0.35) + 0.5)
ax5_0 = ax5.twinx()
ax5_0.plot(x5, y5, color='red', lw=0.5, alpha=0.5)
ax5_0.set_ylim(-0.2, 1.2)
plt.subplots_adjust(left=0.125,
bottom=0.1,
right=0.9,
top=0.9,
wspace=0.1,
hspace=0.35)
plt.show()