import matplotlib.pyplot as plt import numpy as np from IPython import embed import matplotlib as matplotlib import math import scipy.integrate as integrate from scipy import signal from scipy.interpolate import interp1d from scipy.interpolate import CubicSpline import scipy as sp import pickle from scipy.spatial import distance from myfunctions import * import time from matplotlib import gridspec #from matplotlib_scalebar.scalebar import ScaleBar import matplotlib.mlab as ml import scipy.integrate as si import pandas as pd def fish_indirect(delta_t, amplitude, size, beat, phase_zero, sampling, show_figure = False): for a in range(len(amplitude)): for s in range(len(size)): for p in range(len(phase_zero)): for b in range(len(beat)): # sigma calculation sigma = delta_t /math.sqrt((2*math.log(10))) # width of the chirp # time calculation #if np.abs(1 / beat[b]) * 10 < 6 * delta_t: time = np.arange(-3 * delta_t, 3 * delta_t, 1 / sampling) #else: # time = np.arange(-np.abs(1/beat[b]) * 3, abs(1/beat[b]) * 3, 1/ sampling_frequency) # frequency calculation chirp = size[s] * np.exp((-(time**2)/(2*sigma**2))) frequency = beat[b] + chirp # phase calculation I1 = [] for i in range(len(time)): y1 = ((np.pi**0.5)/2)*math.erf(time[i]/sigma) - ((np.pi**0.5)/2)*math.erf(-np.inf) I1.append(y1) phase = 2 * np.pi * beat[b] * time + 2 * np.pi * size[s] * sigma * np.array(I1) + phase_zero[p] # am calculation am = amplitude[a] * np.cos(phase) # plot plt.plot(time*1000, am, label='Amplitude Modulation', color='red', linewidth=2) plt.axvline(x=-7.5, color='black', linewidth=1, linestyle='-') plt.axvline(x=7.5, color='black', linewidth=1, linestyle='-') plt.xlabel('time [ms]') plt.ylabel('AM [mV]') plt.title('beat(%5.2f' % (beat[b]) + 'Hz) ' + 'size(%5.2f' % (size[s]) + 'Hz) ' + 'phase(%5.2f' % (phase_zero[p]) + ') ' + 'amplitude(%1.2f' % (amplitude[a]) + ')') #plt.show() plt.savefig(fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude[a]) +')_size(%5.2f' % (size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_beat(%5.2f' % (beat[b]) + ')_AM_fish%5.2f' % sampling + '.png') if show_figure: plt.show() else: plt.close() def stimulus_indirect(delta_t,amplitude ,size, beat, phase_zero, show_figure = False): for a in range(len(amplitude)): for s in range(len(size)): for p in range(len(phase_zero)): for b in range(len(beat)): # sigma calculation sigma = delta_t /math.sqrt((2*math.log(10))) # width of the chirp # time assignment if np.abs(1 / beat[b]) * 2 < 3 * delta_t: time = np.arange(-1.5 * delta_t, 1.5 * delta_t, 1 / 40000) else: time = np.arange(-np.abs(1/beat[b]) * 1, abs(1/beat[b]) * 1, 1/40000) # frequency calculation chirp = size[s] * np.exp((-(time**2)/(2*sigma**2))) frequency = beat[b] + chirp # phase calculation I1 = [] for i in range(len(time)): y1 = ((np.pi**0.5)/2)*math.erf(time[i]/sigma) - ((np.pi**0.5)/2)*math.erf(-np.inf) I1.append(y1) phase = 2 * np.pi * beat[b] * time + 2 * np.pi * size[s] * sigma * np.array(I1) + phase_zero[p] # am calculation am = amplitude[a] * np.cos(phase) plt.plot(time*1000, am, label='Amplitude Modulation', color='red', linewidth=2) plt.axvline(x=-7.5, color='black', linewidth=1, linestyle='-') plt.axvline(x=7.5, color='black', linewidth=1, linestyle='-') plt.xlabel('time [ms]') plt.ylabel('AM [mV]') plt.title('beat(%5.2f' % (beat[b]) + 'Hz) ' + 'size(%5.2f' % (size[s]) + 'Hz) '+ 'phase(%5.2f' % (phase_zero[p]) + ') '+ 'amplitude(%1.2f' % (amplitude[a]) + ')') plt.savefig(fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude[a]) +')_size(%5.2f' % (size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_beat(%5.2f' % (beat[b]) + ')_AM_stimulus.png') if show_figure: plt.show() else: plt.close() def create_distance(conv_chirp,conv_beat,a): dm_wo_max = np.sqrt(np.mean(np.transpose( (np.transpose(conv_chirp) - np.mean(conv_chirp, a) - np.transpose(conv_beat) + np.mean(conv_beat, a)) ** 2), a)) range_chirp = np.max(conv_chirp,a)-np.min(conv_chirp,a)/np.sqrt(2) range_beat = np.max(conv_beat,a)-np.min(conv_beat,a)/np.sqrt(2) dm_metz = dm_wo_max/np.max(np.transpose([range_chirp,range_beat]),a) dm_wo_mean2 = np.sqrt(np.mean((conv_chirp - conv_beat) ** 2,a)) / np.max(np.transpose([np.max(conv_chirp,a)-np.min(conv_chirp,a)/np.sqrt(2),np.max(conv_beat,a)-np.min(conv_beat,a)/np.sqrt(2)]),a) #dm_wo_mean2_max = np.sqrt(np.mean((conv_chirp - conv_beat) ** 2,a)) dm_wo_mean2_max = mean_euc_d(conv_beat, conv_chirp, a) d_euclidean = euc_d(conv_beat, conv_chirp, a) #dst = distance.euclidean(conv_beat, conv_chirp) return dm_metz, dm_wo_mean2,dm_wo_max,dm_wo_mean2_max,range_chirp,range_beat, d_euclidean def snippets2(peaks_interp, maxima_interp, am_indirect, am, d, period_shift, period_distance, middle_am, eod_convolved_up, middle_conv, delta_t, sampling, deviation, time, eod_conv_downsampled): conv_chirp = eod_convolved_up[ middle_conv - int(0.5 * period_distance * (sampling)):middle_conv + int(0.5 * period_distance * (sampling))] peaks_interp_chirp = peaks_interp[ 3*deviation[d]+middle_conv - int(0.5 * period_distance * (sampling)):3 * deviation[d] + middle_conv + int(0.5 * period_distance * (sampling))] maxima_interp_chirp = peaks_interp[ 3*deviation[d]+middle_conv - int(0.5 * period_distance * (sampling)):3 * deviation[d] + middle_conv + int(0.5 * period_distance * (sampling))] am_indirect_chirp = am_indirect[ middle_conv - int(0.5 * period_distance * (sampling)):middle_conv + int(0.5 * period_distance * (sampling))] time_chirp = time[middle_conv - int(0.5 * period_distance * (sampling)):middle_conv + int( 0.5 * period_distance * (sampling))] am_fish_chirp = am[middle_am - int(0.5 * period_distance * (eod_fr)):middle_am + int( 0.5 * period_distance * (eod_fr))] conv_chirp_corr = eod_conv_downsampled[middle_am - int(0.5 * period_distance * (eod_fr)):middle_am + int( 0.5 * period_distance * (eod_fr))] shift = int(period_shift * (sampling) * 2) shift_fish = int(period_shift * (eod_fr) * 2) am_fish_beat = am[middle_am - int(0.5 * period_distance * (eod_fr)) - shift_fish:middle_am + int( 0.5 * period_distance * (eod_fr)) - shift_fish] conv_beat = eod_convolved_up[middle_conv - int(0.5 * period_distance * (sampling)) - shift:middle_conv + int( 0.5 * period_distance * (sampling)) - shift] peaks_interp_beat = peaks_interp[3 * deviation[d] + middle_conv - int(0.5 * period_distance * (sampling)) - shift:3 * deviation[d] + middle_conv + int( 0.5 * period_distance * (sampling)) - shift] maxima_interp_beat = maxima_interp[3 * deviation[d] + middle_conv - int(0.5 * period_distance * (sampling)) - shift:3 * deviation[d] + middle_conv + int( 0.5 * period_distance * (sampling)) - shift] am_indirect_beat = am_indirect[middle_conv - int(0.5 * period_distance * (sampling)) - shift:middle_conv + int( 0.5 * period_distance * (sampling)) - shift] time_beat = time[middle_conv - int(0.5 * period_distance * (sampling)) - shift: middle_conv + int( 0.5 * period_distance * (sampling)) - shift] return am_fish_beat,maxima_interp_beat,peaks_interp_beat,maxima_interp_chirp, peaks_interp_chirp,am_indirect_beat,am_indirect_chirp,time_beat,conv_beat,conv_chirp_corr,am_fish_chirp,time_chirp,conv_chirp def snippets(am_indirect,am, d,period, middle_am,eod_convolved_up,middle_conv,delta_t,sampling,deviation,time,eod_conv_downsampled): if period <= delta_t*2: conv_chirp = eod_convolved_up[ middle_conv - int(delta_t * (sampling)):middle_conv + int( delta_t * (sampling))] time_chirp = time[middle_conv - int( delta_t * (sampling)):middle_conv + int( delta_t * (sampling))] am_indirect_chirp = am_indirect[ middle_conv - int(delta_t * (sampling)):middle_conv + int(delta_t * (sampling))] am_fish_chirp = am[middle_am - int(delta_t * (eod_fr)):middle_am + int( delta_t * (eod_fr))] conv_chirp_corr = eod_conv_downsampled[middle_am - int( delta_t * (eod_fr)):middle_am + int( delta_t * (eod_fr))] if period == 0: shift = int((delta_t)*3* sampling) else: shift = int((5 * period) * sampling)*2 #if period == 0: # shift = int(((delta_t + 3 * period) * sampling)) #else: # shift = int(((np.ceil(delta_t / period) * (period) + 3 * period) * sampling)) conv_beat = eod_convolved_up[middle_conv - int(delta_t * (sampling)) - shift:middle_conv + int( delta_t * (sampling)) - shift] am_indirect_beat = am_indirect[middle_conv - int( delta_t * (sampling)) - shift:middle_conv + int( delta_t * (sampling)) - shift] time_beat = time[ middle_conv - int( delta_t * (sampling)) - shift:middle_conv + int( delta_t * (sampling)) - shift] else: conv_chirp = eod_convolved_up[ middle_conv - int( period * (sampling)):middle_conv + int( period * (sampling))] am_indirect_chirp = am_indirect[ middle_conv - int(period * (sampling)):middle_conv + int( period * (sampling))] time_chirp = time[middle_conv - int( period * (sampling)):middle_conv + int( period * (sampling))] am_fish_chirp = am[middle_am - int( period * (eod_fr)):middle_am + int( period * (eod_fr))] conv_chirp_corr = eod_conv_downsampled[middle_am - int( period * (eod_fr)):middle_am + int( period * (eod_fr))] shift = int(period * (sampling) * 2) conv_beat = eod_convolved_up[middle_conv - int( period * (sampling)) - shift:middle_conv + int( period * (sampling)) - shift] am_indirect_beat = am_indirect[middle_conv - int(period * (sampling)) - shift:middle_conv + int( period * (sampling)) - shift] time_beat = time[middle_conv - int( period * (sampling)) - shift: middle_conv + int( period * (sampling)) - shift] return am_indirect_beat,am_indirect_chirp,time_beat,conv_beat,conv_chirp_corr,am_fish_chirp,time_chirp,conv_chirp def global_maxima_array(period_fish_e,period_fish_r,eod_rectified_up): period_length = np.max([period_fish_e, np.ones(len(period_fish_e))*len(period_fish_r)],0) split_windows = np.linspace(period_length, len(eod_rectified_up), len(eod_rectified_up)/period_length) splits = np.split(eod_rectified_up, split_windows) steps = np.arange(0, len(eod_rectified_up), len(splits[0])) maxima_values = np.max(splits[0:-1], 1) maxima_index = np.argmax(splits[0:-1], 1) maxima_index = maxima_index + steps[0:-1] maxima_interp = np.interp(np.arange(0, len(eod_rectified_up), 1), maxima_index, maxima_values) return maxima_values,maxima_index, maxima_interp def global_maxima(sampling, eod_fr, period_fish_e,period_fish_r,eod_rectified_up): #period_length = max(len(period_fish_e), len(period_fish_r)) period_length = len(period_fish_r) period_length = int(np.round((1/eod_fr)*sampling)) #embed() if period_length >len(eod_rectified_up): maxima_values = np.max(eod_rectified_up) maxima_index = np.argmax(eod_rectified_up) maxima_interp = [maxima_values]*len(eod_rectified_up) else: split_windows = np.arange(period_length, len(eod_rectified_up), period_length) splits = np.split(eod_rectified_up, split_windows) steps = np.arange(0, len(eod_rectified_up), len(splits[0])) maxima_values = np.max(splits[0:-1], 1) maxima_index = np.argmax(splits[0:-1], 1) maxima_index = maxima_index + steps[0:-1] maxima_interp = np.interp(np.arange(0, len(eod_rectified_up), 1), maxima_index, maxima_values) return maxima_values,maxima_index, maxima_interp def snippets3(skretch,peaks_interp, maxima_interp, am_indirect,am_fish, d,period, middle_am,eod_convolved_up,middle_conv,delta_t,sampling,deviation,time,eod_conv_downsampled): #period_distance = 0.1 period_dist = [] if period <= delta_t*skretch: if period == 0: period_shift = delta_t*skretch period_distance = delta_t*skretch else: period_shift = np.ceil(delta_t*skretch/period)*period period_distance = np.ceil(delta_t*skretch/period)*period am_fish_beat,maxima_interp_beat,peaks_interp_beat,maxima_interp_chirp, peaks_interp_chirp,am_indirect_beat,am_indirect_chirp,time_beat,conv_beat,conv_chirp_corr,am_fish_chirp,time_chirp,conv_chirp = \ snippets2(peaks_interp, maxima_interp, am_indirect, am_fish, d, period_shift,period_distance, middle_am, eod_convolved_up, middle_conv, delta_t, sampling, deviation, time, eod_conv_downsampled) else: period_shift = period period_distance = period am_fish_beat,maxima_interp_beat,peaks_interp_beat,maxima_interp_chirp, peaks_interp_chirp,am_indirect_beat,am_indirect_chirp,time_beat,conv_beat,conv_chirp_corr,am_fish_chirp,time_chirp,conv_chirp = \ snippets2( peaks_interp, maxima_interp, am_indirect, am_fish, d, period_shift,period_distance, middle_am, eod_convolved_up, middle_conv, delta_t, sampling, deviation, time, eod_conv_downsampled) period_dist.append(period_distance) return period_dist, am_fish_beat,maxima_interp_beat,peaks_interp_beat,maxima_interp_chirp, peaks_interp_chirp,am_indirect_beat,am_indirect_chirp,time_beat,conv_beat,conv_chirp_corr,am_fish_chirp,time_chirp,conv_chirp def conv(eod_fr, sampling, cut,deviation, eod_rectified_up, eod_rectified_down): if deviation* 5 % 2: points = deviation * 5 else: points = deviation * 5 - 1 gaussian = signal.gaussian(points, std=deviation, sym=True) gaussian_normalised = (gaussian * 2) / np.sum(gaussian) length_convolved = int(len(gaussian_normalised) / 2) eod_convolved_up = np.convolve(gaussian_normalised, eod_rectified_up) eod_convolved_up = eod_convolved_up[length_convolved + cut:-length_convolved - cut] eod_convolved_down = np.convolve(gaussian_normalised, eod_rectified_down) eod_convolved_down = eod_convolved_down[length_convolved + cut:-length_convolved - cut] middle_conv = int(len(eod_convolved_up) / 2) eod_conv_downsampled = eod_convolved_up[0:-1:int(np.round(sampling / eod_fr))] return middle_conv, eod_convolved_down, eod_convolved_up,eod_conv_downsampled def find_periods(a_fe, time, eod_fr,a_fr,eod_fe,e): time_fish_r = time * 2 * np.pi * eod_fr eod_fish_r = a_fr * np.sin(time_fish_r) period_fish_r = time_fish_r[(time_fish_r <= np.mean(time_fish_r)+2 * np.pi) & (time_fish_r > np.mean(time_fish_r))] time_fish_e = time * 2 * np.pi * eod_fe[e] eod_fish_e = a_fe * np.sin(time_fish_e) period_fish_e = time_fish_e[(time_fish_e <= np.mean(time_fish_e)+ 2 * np.pi) & (time_fish_e > np.mean(time_fish_e))] return eod_fish_e,eod_fish_r,period_fish_r,period_fish_e def rectify(eod_fish_r,eod_fe_chirp): eod_rec_up = eod_fish_r + eod_fe_chirp eod_rectified_down = eod_fish_r + eod_fe_chirp eod_rec_up[eod_rec_up < 0] = 0 # rectify eod_rectified_down[eod_rectified_down > 0] = 0 # rectify return eod_rectified_down, eod_rec_up def integrate_chirp(a_fe,time,beat,phase_zero,size, sigma): I = ((np.pi ** 0.5) / 2) * sp.special.erf(time / sigma) - ((np.pi ** 0.5) / 2) * sp.special.erf(-np.inf) phase = time * 2 * np.pi * beat+ 2 * np.pi * size * sigma * I + phase_zero eod_fe_chirp = a_fe * np.sin(phase) return eod_fe_chirp def find_lm(eod_rec_up): x = signal.find_peaks(eod_rec_up) index_peaks = x[0] value_peaks = eod_rec_up[index_peaks] peaks_interp = np.interp(np.arange(0, len(eod_rec_up), 1), index_peaks, value_peaks) return index_peaks, value_peaks, peaks_interp def find_times(left_c,right_c, sampling,deviation_s): for_conv = 5 * deviation_s time = np.arange(int(np.round(left_c))-1000, int(np.round(right_c))+1000, 1) time = time[(time >left_c) &(time < right_c)] time = time/sampling #time = np.arange(-for_conv+left_c,for_conv+right_c, 1 / sampling) cut = int(np.ceil(for_conv*sampling)) time_cut = time[cut:-cut] return time, time_cut, cut def snip(left_c,right_c,e,g,sampling, deviation_s,d,eod_fr, a_fr, eod_fe,phase_zero,p, size,s, sigma,a_fe,deviation,beat_corr, chirp = True): time, time_cut, cut = find_times(left_c[g], right_c[g], sampling, deviation_s[d]) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) #embed() if chirp == True: eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) else: eod_fe_chirp = eod_fish_e eod_rec_down, eod_rec_up = rectify(eod_fish_r, eod_fe_chirp) # rectify eod_overlayed_chirp = (eod_fish_r + eod_fe_chirp)[cut:-cut] threshold_cube = (eod_rec_up) ** 3 period_length = len(period_fish_r) #embed() sampled = np.abs(eod_overlayed_chirp[0:len(eod_overlayed_chirp):period_length]) time_new = time[0:len(eod_overlayed_chirp):period_length] maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_rec_up[cut:-cut]) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_rec_up[cut:-cut]) # local maxima middle_conv, eod_conv_down, eod_conv_up, eod_conv_downsampled = conv(eod_fr,sampling, cut, deviation[d], eod_rec_up, eod_rec_down) # convolve eod_fish_both = integrate_chirp(a_fe, time, eod_fe[e] - eod_fr, phase_zero[p], size[s], sigma) am_corr_full = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p], size[s], sigma) # indirect am calculation _, time_fish, cut_f = find_times(left_c[g], right_c[g], eod_fr, deviation_s[d]) # downsampled through fish EOD am_corr_ds = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) am_df_ds = integrate_chirp(a_fe, time_fish, eod_fe[e] - eod_fr, phase_zero[p], size[s], sigma) # indirect am calculation return time_new,sampled, cut, threshold_cube , time_cut, eod_conv_up, am_corr_full, peaks_interp, maxima_interp, am_corr_ds,am_df_ds,eod_fish_both,eod_overlayed_chirp def same_lengths(maxima_interp_b, maxima_interp_c, time_b, time_c, am_fish_c, am_fish_b, conv_c, conv_b, am_indirect_b, am_indirect_c, peaks_interp_b, peaks_interp_c): diff_s = len(time_b) - len(time_c) diff_f = len(am_fish_c) - len(am_fish_b) if len(time_b) > len(time_c): time_b = time_b[0:len(time_c)] conv_b = conv_b[0:len(conv_c)] am_indirect_b = am_indirect_b[0:len(am_indirect_c)] peaks_interp_b = peaks_interp_b[0:len(peaks_interp_c)] maxima_interp_b = maxima_interp_b[0:len(maxima_interp_c)] if len(time_c) > len(time_b): time_c = time_c[0:len(time_b)] conv_c = conv_c[0:len(conv_b)] am_indirect_c = am_indirect_c[0:len(am_indirect_b)] peaks_interp_c = peaks_interp_c[0:len(peaks_interp_b)] maxima_interp_c = maxima_interp_c[0:len(maxima_interp_b)] if len(am_fish_b) > len(am_fish_c): am_fish_b = am_fish_b[0:len(am_fish_c)] if len(am_fish_c) > len(am_fish_b): am_fish_c = am_fish_c[0:len(am_fish_b)] return diff_s, diff_f, maxima_interp_b, maxima_interp_c, time_b, time_c, am_fish_c, am_fish_b, conv_c, conv_b, am_indirect_b, am_indirect_c, peaks_interp_b, peaks_interp_c def power_func(restrict = [],a_fr = 1, a_fe = 0.2, eod_fr = 500, eod_fe = [500], win = 'w2', deviation_s = [0.0015], sigma = 50, sampling = 100000, deviation_ms = [1.5], beat_corr = [0], size = [120],phase_zero = [0], delta_t = 0.014,deviation_dp = [180], show_figure = True, plot_dist = False, save = False,bef_c = -2,aft_c = -1): print('plot_eod') results = [[]]*10 if len(restrict) >0: results = [[]]*len(restrict) for d in range(len(deviation_s)): left_c, right_c, left_b, right_b, period_distance_c, period_distance_b, _, period, to_cut,exclude,consp_needed,deli,interval = period_calc( [float('Inf')]*len(beat_corr), 50, win, deviation_s[d], sampling, beat_corr, bef_c, aft_c, 'stim') for s in range(len(size)): for p in range(len(phase_zero)): for e in range(len(eod_fe)): #embed() left_b = [bef_c * sampling] * len(beat_corr) right_b = [aft_c * sampling] * len(beat_corr) time_fish,sampled, cut, cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b,am_corr_ds_b, am_df_ds_b, am_df_b,eod_overlayed_chirp = snip(left_b, right_b,e,e, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation_dp, beat_corr, chirp = False) thresholding = eod_overlayed_chirp*1 thresholding[eod_overlayed_chirp<0] = 0 squaring = eod_overlayed_chirp**2 fifth = (thresholding)**5 #embed() rectified = np.abs(eod_overlayed_chirp)#'df ds','corr','corr ds', name = ['df','conv','thresholding','rectified','squaring','cubed','fifth','global max','local max', 'samped threshold']#am_df_b,am_df_ds_b, am_corr_b, am_corr_ds_b, var = [am_df_b,conv_b, thresholding, rectified, squaring,cubed ,fifth, maxima_interp_b,peaks_interp_b, sampled] samp = [sampling,sampling,sampling, sampling, sampling,sampling,sampling,sampling, sampling, eod_fr]#eod_fr,sampling,eod_fr, #['Hilbert', 'Thresholded', 'Rectified', 'Squaring', 'Threshold cubed'] #embed() if len(restrict) == 1: var = [np.array(var)[restrict[0]]] name = [np.array(name)[restrict[0]]] samp = [np.array(samp)[restrict[0]]] elif len(restrict)>1: var = np.array(var)[restrict] name = np.array(name)[restrict] samp = np.array(samp)[restrict] #embed() for i in range(len(results)): if samp[i]>1000: nfft = int(4096 * samp[0] / 40000 * 2) else: nfft = 4096 results[i],pp,f = power_parameters(results[i], var[i], samp[i],nfft,name[i],eod_fr) results = pd.DataFrame(results) if save: results = pd.DataFrame(results) try: results.to_pickle('../data/power_simulation.pkl') np.save('../data/Ramona/power_simulation.npy', results) except: a = 0 return results def power_parameters(results, var,samp,nfft,name, eod_fr,max_all = 'max_all',max = 'max',integral = 'amp',plot = False,title = '',freq_half = 'f',freq_whole = 'ff'): pp, f = ml.psd(var - np.mean(var), Fs=samp, NFFT=nfft, noverlap=nfft / 2) # # embed() # print(nfft * 5 < len(var[i]): if plot == True: plt.figure() plt.subplot(1, 2, 1) plt.plot(var) plt.title(name) plt.subplot(1, 2, 2) plt.plot(f, pp) # plt.xlim([0,2000]) plt.show() # print(results) freq_step = np.abs(f[0] - f[1]) if type(results) != dict: #embed() results = {} results['type'+title] = name # embed() results['result_frequency_whole'+title] = list([f[np.argmax(pp)]]) results['result_frequency'+title] = list([f[np.argmax(pp[f < 0.5 * eod_fr])]]) results['result_amplitude'+title] = list([np.sqrt(si.trapz(pp, f, np.abs(f[1] - f[0])))]) results['result_amplitude_max'+title] = list([np.sqrt(np.max(pp[f < 0.5 * eod_fr]) * freq_step)]) results['result_amplitude_max_whole'+title] = list([np.sqrt(np.max(pp) * freq_step)]) else: #embed() results['result_frequency'+title].extend(list([f[np.argmax(pp[f < 0.5 * eod_fr])]])) results['result_frequency_whole'+title].extend(list([f[np.argmax(pp)]])) # embed() results['result_amplitude'+title].extend(list([np.sqrt(si.trapz(pp, f, np.abs(f[1] - f[0])))])) results['result_amplitude_max'+title].extend(list([np.sqrt(np.max(pp[f < 0.5 * eod_fr]) * freq_step)])) results['result_amplitude_max_whole'+title].extend(list([np.sqrt(np.max(pp) * freq_step)])) return results,pp,f def distances_func(disc, bef_c, aft_c, win, deviation_s, sigma, sampling, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = False, plot_dist = False, save = False): distance_metz_all = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_all = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_max_all = only_one_dict_str(d_ms, eod_fe) distance_wo_max_all = only_one_dict_str(d_ms, eod_fe) range_chirp_all = only_one_dict_str(d_ms, eod_fe) range_beat_all = only_one_dict_str(d_ms, eod_fe) corr_all = only_one_dict_str(d_ms, eod_fe) distance_metz_all_corr = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_all_in = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_max_all_in = only_one_dict_str(d_ms, eod_fe) distance_wo_max_all_in = only_one_dict_str(d_ms, eod_fe) range_chirp_all_in = only_one_dict_str(d_ms, eod_fe) range_beat_all_in = only_one_dict_str(d_ms, eod_fe) distance_metz_all_gm = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_all_gm = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_max_all_gm = only_one_dict_str(d_ms, eod_fe) distance_wo_max_all_gm = only_one_dict_str(d_ms, eod_fe) distance_metz_all_lm = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_all_lm = only_one_dict_str(d_ms, eod_fe) distance_wo_mean_max_all_lm = only_one_dict_str(d_ms, eod_fe) distance_wo_max_all_lm = only_one_dict_str(d_ms, eod_fe) range_chirp_all_lm = only_one_dict_str(d_ms, eod_fe) distance_metz_all_df_ds,distance_metz_all_corr_ds,range_beat_all_lm,range_chirp_all_gm ,range_beat_all_gm,distance_metz_all_df,distance_wo_mean_all_df,distance_wo_mean_max_all_df,distance_wo_max_all_df,range_chirp_all_df ,euc_all_in,euc_all_df,euc_all_gm,euc_all_lm,euc_all,range_beat_all_df = (only_one_dict_str(d_ms, eod_fe) for i in range(16)) period_distance_c,period_distance_b,period_shift = ({} for i in range(3)) for d in range(len(deviation)): left_c, right_c, left_b, right_b, period_distance_c[d_ms[d]], period_distance_b[d_ms[d]], period_shift[d_ms[d]], period, to_cut,exclude,consp_needed,deli,interval = period_calc( [float('Inf')]*len(beat_corr), 50, win, deviation_s[d], sampling, beat_corr, bef_c, aft_c, 'stim') save_n = win for s in range(len(size)): for p in range(len(phase_zero)): beats = eod_fe - eod_fr for e in range(len(eod_fe)): time_fish,sampled, cut, cube_b, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b,am_corr_ds_b, am_df_ds_b, am_df_b,eod_overlayed_chirp_b = snip(left_b, right_b,e,e, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation, beat_corr) time_fish,sampled,cut, cube_c, time_c, conv_c, am_corr_c, peaks_interp_c, maxima_interp_c,am_corr_ds_c, am_df_ds_c, am_df_c,eod_overlayed_chirp_c = snip(left_c, right_c,e,e, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation, beat_corr) #print(len(time_c)) #print(len(time_b)) #diff_s[e], diff_f[e],maxima_interp_b, maxima_interp_c, time_b, time_c, am_df_c, am_df_b, conv_c,conv_b, am_corr_b, am_corr_c, peaks_interp_b, peaks_interp_c = same_lengths( # maxima_interp_b, maxima_interp_c, time_b, time_c, am_df_c, am_df_b,conv_c, conv_b, am_corr_b, # am_corr_c, peaks_interp_b, peaks_interp_c) #time_b[e], conv_b[e], am_corr_b[e], peaks_interp_b[e], maxima_interp_b[e], am_df_b[e] = snip(left_c-shift[e], right_c-shift[e], e, # sampling, deviation_s, # d, eod_fr, a_fr, eod_fe, # phase_zero, p, size, s, # sigma, a_fe, deviation, # beat_corr) #time_c[e], conv_c[e], am_corr_c[e], peaks_interp_c[e], maxima_interp_c[e], am_df_c[e] = snip(left_c, right_c, e, # sampling, deviation_s, # d, eod_fr, a_fr, eod_fe, # # phase_zero, p, size, s, # sigma, a_fe, deviation, # beat_corr) distance_metz_all_df_ds[d_ms[d]][e], _,_,_,_,_,_ = create_distance(np.array(am_df_ds_c), np.array(am_df_ds_b), None) distance_metz_all_corr_ds[d_ms[d]][e], _,_,_,_,_,_ = create_distance(np.array(am_corr_ds_c), np.array(am_corr_ds_b), None) distance_metz_all[d_ms[d]][e], distance_wo_mean_all[d_ms[d]][e],distance_wo_max_all[d_ms[d]][e],distance_wo_mean_max_all[d_ms[d]][e],range_chirp_all[d_ms[d]][e],range_beat_all[d_ms[d]][e],euc_all[d_ms[d]][e] = create_distance(np.array(conv_c),np.array(conv_b),None) distance_metz_all_corr[d_ms[d]][e], distance_wo_mean_all_in[d_ms[d]][e],distance_wo_max_all_in[d_ms[d]][e],distance_wo_mean_max_all_in[d_ms[d]][e],range_chirp_all_in[d_ms[d]][e],range_beat_all_in[d_ms[d]][e],euc_all_in[d_ms[d]][e] = create_distance(np.array(am_corr_c),np.array(am_corr_b),None) distance_metz_all_gm[d_ms[d]][e], distance_wo_mean_all_gm[d_ms[d]][e],distance_wo_max_all_gm[d_ms[d]][e],distance_wo_mean_max_all_gm[d_ms[d]][e],range_chirp_all_gm[d_ms[d]][e],range_beat_all_gm[d_ms[d]][e],euc_all_gm[d_ms[d]][e] = create_distance(np.array(maxima_interp_b),np.array(maxima_interp_c),None) distance_metz_all_lm[d_ms[d]][e], distance_wo_mean_all_lm[d_ms[d]][e], distance_wo_max_all_lm[d_ms[d]][e], distance_wo_mean_max_all_lm[d_ms[d]][e], range_chirp_all_lm[d_ms[d]][e], range_beat_all_lm[d_ms[d]][e],euc_all_lm[d_ms[d]][e] = create_distance( np.array(peaks_interp_c), np.array(peaks_interp_b),None) distance_metz_all_df[d_ms[d]][e], distance_wo_mean_all_df[d_ms[d]][e], distance_wo_max_all_df[d_ms[d]][e], distance_wo_mean_max_all_df[d_ms[d]][e], range_chirp_all_df[d_ms[d]][e], range_beat_all_df[d_ms[d]][e],euc_all_df[d_ms[d]][e] = create_distance( np.array(am_df_c), np.array(am_df_b),None) #distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, range_chirp_all, range_beat_all, euc_all = create_distance( # np.array(conv_c), np.array(conv_b), None) #distance_metz_all_corr, distance_wo_mean_all_in, distance_wo_max_all_in, distance_wo_mean_max_all_in, range_chirp_all_in, range_beat_all_in, euc_all_in = create_distance( # np.array(am_corr_c), np.array(am_corr_b), None) #distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_max_all_gm, distance_wo_mean_max_all_gm, range_chirp_all_gm, range_beat_all_gm, euc_all_gm = create_distance( # np.array(maxima_interp_b), np.array(maxima_interp_c), None) #distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_max_all_lm, distance_wo_mean_max_all_lm, range_chirp_all_lm, range_beat_all_lm, euc_all_lm = create_distance( # np.array(peaks_interp_c), np.array(peaks_interp_b), None) #distance_metz_all_df, distance_wo_mean_all_df, distance_wo_max_all_df, distance_wo_mean_max_all_df, range_chirp_all_df, range_beat_all_df, euc_all_df = create_distance( # np.array(am_df_b), np.array(am_df_c), None) # correlation #corr = np.corrcoef(am_chirp, conv_chirp_corr) #corr_all.append(corr[1,0]) if save: np.save('../data/' + 'sim_dist3_' + str( d_ms[d]) + 'ms' + 'beats(%1.2f' % (beats[0]) + '_%1.2f' % (beats[-1]) + '_%1.2f' % (abs(beats[1] - beats[0])) + ')' + save_n + '.npy',[distance_metz_all_df_ds,distance_metz_all_corr_ds,period_shift,period_distance_c, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_df, range_beat_all_df, range_chirp_all_df, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm, beat_corr, distance_metz_all_df, distance_wo_mean_all_df,distance_wo_mean_max_all_df, distance_wo_max_all_df, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, d_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_corr, distance_wo_mean_all_in, eod_fe, eod_fr,distance_metz_all, distance_wo_mean_all, distance_wo_max_all,distance_wo_mean_max_all]) with open('../data/' + 'sim_dist3_' + str( d_ms[d]) + 'ms' + 'beats(%1.2f' % (beats[0]) + '_%1.2f' % (beats[-1]) + '_%1.2f' % (abs(beats[1] - beats[0])) + ')' + save_n + '.pkl', 'wb') as f: # Python 3: open(..., 'wb') pickle.dump([distance_metz_all_df_ds,distance_metz_all_corr_ds,period_shift,period_distance_c, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_df, range_beat_all_df, range_chirp_all_df, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm, beat_corr, distance_metz_all_df, distance_wo_mean_all_df, distance_wo_mean_max_all_df, distance_wo_max_all_df, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, d_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_corr, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all], f) return distance_metz_all_df_ds,distance_metz_all_corr_ds,period_shift,period_distance_c,period_distance_b, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_df, range_beat_all_df, range_chirp_all_df, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm \ , beat_corr, distance_metz_all_df, distance_wo_mean_all_df, distance_wo_mean_max_all_df, distance_wo_max_all_df, \ distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, \ distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, \ distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, \ d_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, \ distance_metz_all_corr, distance_wo_mean_all_in, eod_fe, eod_fr, \ distance_metz_all, distance_wo_mean_all, distance_wo_max_all, \ distance_wo_mean_max_all def plot_dist_interactiv3(bef_c, aft_c,d, sampling, size, phase_zero, sigma, beat_corr, delta_t, a_fr, a_fe, deviation_dp, deviation_s, d_ms, eod_fe, eod_fr, share = False, show = False): fig = plt.figure(figsize = [13,6]) plt.title('convolved (red), indirect (magenta), loc max (green), glo max (black),sampled f (orange), period (purple)') gs0 = gridspec.GridSpec(3, 2) win = ['w1', 'w2', 'w3', 'w4'] ax = {} gs00 = gridspec.GridSpecFromSubplotSpec(7, 2, subplot_spec=gs0[3]) ax[0] = plt.subplot(gs00[4]) p = 0 s = 0 for g in range(4): period_shift, period_distance_c, period_distance_b, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, range_beat_all_ds, range_chirp_all_ds, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm, beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, d_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, \ distance_wo_mean_max_all = distances_func(bef_c, aft_c, win[g], deviation_s, sigma, sampling, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation_dp, show_figure=False, plot_dist=False, save=True) gs00 = gridspec.GridSpecFromSubplotSpec(7, 2, subplot_spec=gs0[g]) #plt.title('w1') #fig = frame_subplots(6.5, 3, 'Beats [Hz]', '', 'Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (d_ms[d]) + ') '+win) #two_ylabels(1,2,1,fig, 'distance (not normalised)', 'distance (normalised)') subplots = 14 #ax = create_subpl(0,fig, subplots, 7, 2, share2 = share) plots = [distance_metz_all[d_ms[d]], euc_all[d_ms[d]],distance_metz_all_in[d_ms[d]],euc_all_in[d_ms[d]],distance_metz_all_lm[d_ms[d]], euc_all_lm[d_ms[d]], distance_metz_all_gm[d_ms[d]],euc_all_gm[d_ms[d]], distance_metz_all_ds[d_ms[d]], euc_all_ds[d_ms[d]],period_distance_c[d_ms[d]],period_distance_b[d_ms[d]],period_shift[d_ms[d]],period_shift[d_ms[d]]] labels = ['conv mean','conv sum', 'indirect mean', 'indirect sum', 'lm mean', 'lm sum', 'gm mean', 'gm sum','ds mean', 'ds sum','period lengths','','period shift',''] colors = ['red','red', 'magenta','magenta','green','green','black','black','orange','orange','darkblue','darkblue','lime','lime'] I = {} for i in range(subplots): if share: ax[i] = plt.subplot(gs00[i],sharex=ax[0]) else: ax[i] = plt.subplot(gs00[i]) I[i] = ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') ax[i].set_yticks([]) axvline_everywhere(ax, subplots, eod_fe, eod_fr) ax[0].set_title(win[g] +' Mean ') ax[1].set_title(win[g]+' Sum (Euclidean) ') plots = [distance_wo_mean_all[d_ms[d]],distance_wo_mean_all_in[d_ms[d]],distance_wo_mean_all_lm[d_ms[d]],distance_wo_mean_all_gm[d_ms[d]],distance_wo_mean_all_ds[d_ms[d]]] for i in np.arange(0,subplots-4,2): ax[i].plot(eod_fe - eod_fr, plots[int(i/2)], color=colors[i], linestyle='--') plots1 = [distance_wo_max_all[d_ms[d]],distance_wo_max_all_in[d_ms[d]],distance_wo_max_all_lm[d_ms[d]],distance_wo_max_all_gm[d_ms[d]],distance_wo_max_all_ds[d_ms[d]]] plots2 = [distance_wo_mean_max_all[d_ms[d]], distance_wo_mean_max_all_in[d_ms[d]], distance_wo_mean_max_all_lm[d_ms[d]], distance_wo_mean_max_all_gm[d_ms[d]], distance_wo_mean_max_all_ds[d_ms[d]]] for i in np.arange(0, subplots - 4, 2): ax[20+i] = ax[i].twinx() ax[20 +i].set_yticks([]) ax[20+i].plot(eod_fe - eod_fr, plots1[int(i/2)],color=colors[i],linestyle='dashed') ax[20+i].plot(eod_fe - eod_fr, plots2[int(i/2)],linestyle='-.', color=colors[i]) plt.subplots_adjust(hspace = 0.15, wspace = 0.15) #fig.legend(I,labels=labels) #fig.legend([I[1],I[3],I[5],I[7],I[9],I[11]],labels = ['convolved','indirect','lm','gm','ds','period']) if win[g] == 'w1': ax[0].set_ylabel('convolved',color = 'red', rotation = 360, labelpad = 40) ax[2].set_ylabel('perfect', color='magenta', rotation = 360, labelpad = 40) ax[4].set_ylabel('local max', color='green', rotation = 360, labelpad = 40) ax[6].set_ylabel('global max', color='black', rotation = 360, labelpad = 40) ax[8].set_ylabel('sampled f', color='orange', rotation = 360, labelpad = 40) ax[10].set_ylabel('period dist', color='darkblue', rotation = 360, labelpad = 40) ax[12].set_ylabel('period shift', color='lime', rotation=360, labelpad=40) ax[8].set_xticks([]) ax[10].set_xticks([]) if win[g] == 'w3': ax[0].set_ylabel('w1: pd and ps, mult', rotation=360, labelpad=60) ax[2].set_ylabel('w2: same pd', rotation=360, labelpad=60) ax[4].set_ylabel('w3: pd constant high f', rotation=360, labelpad=60) ax[6].set_ylabel('w4, same pd and ps', rotation=360, labelpad=60) if win[g] == 'w2': ax[8].set_xticks([]) ax[10].set_xticks([]) #ax = plt.gca() #leg = ax.get_legend() #leg.legendHandles[0].set_color('red') #leg.legendHandles[1].set_color('darkblue') def onclick(event): e = np.argmin(np.abs((eod_fe - eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t,sampling, deviation_s,d,a_fe, eod_fr, a_fr, eod_fe, e,phase_zero,p, size,s, sigma) plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) #plt.savefig( # fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( # size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( # eod_fr) + ')_beat_%1.2f' % ( # beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1]) + ')_deviation(%1.2f' % ( # d_ms[d]) + '_distances_all' + 'several' + '.png') g = 4 gs00 = gridspec.GridSpecFromSubplotSpec(4, 2, subplot_spec=gs0[g]) with open('C:/Users/alexi/Desktop/Masterarbeit/labrotation/data/consp4_1.5msbeats(-500.00_1995.00_5.00).pkl', 'rb') as f: # Python 3: open(..., 'rb') eod_fe, beat_corr, beats, integral_ds, integral_lm, integral_gm, integral_ind, integral_conv, max_values_gm, max_positions_gm, max_values_lm, max_positions_lm, max_values_conv, max_positions_conv, max_values_ind, max_positions_ind, max_values_ds, max_positions_ds = pickle.load( f) d = 1 #two_ylabels(1,2,1,fig, 'Integral beneath the chirp Consp', '') #two_ylabels(1, 2, 2, fig, '', 'Conspicousness [range 0 to 1]') subplots = 8 #ax = create_subpl(0.3,fig, subplots, 4, 2, share2 = share) plots = [integral_conv, max_values_conv,integral_lm,max_values_lm, integral_gm, max_values_gm, integral_ind, max_values_ind] labels = ['Integral conv','Max Values conv','Integral lm','Max Values lm','Integral gm','Max Values gm','Integral ind','Max Values ind'] colors = ['red','red', 'green','green','black','black','magenta','magenta','orange','orange','darkblue','darkblue'] for i in range(subplots): if share: ax[i] = plt.subplot(gs00[i], sharex=ax[0]) else: ax[i] = plt.subplot(gs00[i]) ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') axvline_everywhere(ax,subplots, eod_fe, eod_fr) #ax[0].set_ylabel('convolved',color = 'red', rotation) #ax[2].set_ylabel('loc max', color='green') #ax[4].set_ylabel('glo max', color='black') #ax[6].set_ylabel('perfect beat corr', color='magenta') ax[0].set_title('Integral') ax[1].set_title('Max Value') for i in np.arange(0,8,2): ax[i].set_ylim([0,int(20/3) * 3]) for i in np.arange(1, 8, 2): ax[i].set_ylim([0, 1.1]) #plt.legend(loc='upper left',bbox_to_anchor=(0, 2)) period = period_func(beat_corr, sampling, 'stim') def onclick2(event): e = np.argmin(np.abs((eod_fe-eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t, sampling, deviation_s, d, a_fe, eod_fr, a_fr, eod_fe, e, phase_zero, p, size, s, sigma) rang = [0.5*sampling] time_fish,sampled, cut, cube_c, time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp = snip([-rang[0] * 0.5]*len(eod_fe), [rang[0] * 0.5]*len(eod_fe), e,e,sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, 0, size, 0, sigma, a_fe, deviation, beat_corr) B_conv, B_indirect,B_downsampled,B_lm,B_gm,time = snippets4(deviation_s,a_fr, sigma,period,time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp, eod_fe, beat_corr, a_fe, eod_fr, e, size, 0, phase_zero, 0, deviation, d, rang[0], delta_t, sampling, plot = False,show = False) #fig = plt.figure(figsize = [7,6]) gaussian = gauss(10,sampling,0.02,1) pl = True sh = False if event.inaxes == ax[4] or event.inaxes == ax[5]: max_position, max_value, I,ConSpi2 = single_consp(1,eod_fe[e],eod_fr,'black',gaussian,B_gm, gm_chirp, sampling,plot = pl,show = sh) if event.inaxes == ax[2] or event.inaxes == ax[3]: max_position, max_value, I,ConSpi2 = single_consp(1,eod_fe[e],eod_fr,'green',gaussian, B_lm, lm_chirp, sampling,plot = pl,show = sh) if event.inaxes == ax[0] or event.inaxes == ax[1]: max_position, max_value, I,ConSpi2 = single_consp(1,eod_fe[e],eod_fr,'red',gaussian, B_conv, conv_chirp, sampling,plot = pl,show = sh) if event.inaxes == ax[6] or event.inaxes == ax[7]: max_position, max_value, I,ConSpi2 = single_consp(1,eod_fe[e],eod_fr,'magenta',gaussian,B_indirect, am_indirect_chirp, sampling,plot = pl,show = sh) plt.xlabel('Time [ms]') plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) if show: plt.show() def plot_dist_interactiv4(bef_c, aft_c,disc,d, sampling, size, phase_zero, sigma, beat_corr, delta_t, a_fr, a_fe, deviation_dp, deviation_s, d_ms, eod_fe, eod_fr, share = False, show = False): fig = plt.figure(figsize = [12,5]) plt.title('convolved (red), indirect (magenta), loc max (green), glo max (black),sampled f (orange), period (purple)') gs0 = gridspec.GridSpec(2, 2) win = ['w1', 'w2', 'w3', 'w4'] ax = {} gs00 = gridspec.GridSpecFromSubplotSpec(7, 2, subplot_spec=gs0[3]) ax[0] = plt.subplot(gs00[4]) p = 0 s = 0 for g in range(4): period_shift, period_distance_c, period_distance_b, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, range_beat_all_ds, range_chirp_all_ds, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm, beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, d_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, \ distance_wo_mean_max_all = distances_func(disc, bef_c, aft_c, win[g], deviation_s, sigma, sampling, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation_dp, show_figure=False, plot_dist=False, save=True) gs00 = gridspec.GridSpecFromSubplotSpec(7, 2, subplot_spec=gs0[g]) #plt.title('w1') #fig = frame_subplots(6.5, 3, 'Beats [Hz]', '', 'Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (d_ms[d]) + ') '+win) #two_ylabels(1,2,1,fig, 'distance (not normalised)', 'distance (normalised)') subplots = 14 #ax = create_subpl(0,fig, subplots, 7, 2, share2 = share) plots = [distance_metz_all[d_ms[d]], euc_all[d_ms[d]],distance_metz_all_in[d_ms[d]],euc_all_in[d_ms[d]],distance_metz_all_lm[d_ms[d]], euc_all_lm[d_ms[d]], distance_metz_all_gm[d_ms[d]],euc_all_gm[d_ms[d]], distance_metz_all_ds[d_ms[d]], euc_all_ds[d_ms[d]],period_distance_c[d_ms[d]],period_distance_b[d_ms[d]],period_shift[d_ms[d]],period_shift[d_ms[d]]] labels = ['conv mean','conv sum', 'indirect mean', 'indirect sum', 'lm mean', 'lm sum', 'gm mean', 'gm sum','ds mean', 'ds sum','period lengths','','period shift',''] colors = ['red','red', 'magenta','magenta','green','green','black','black','orange','orange','darkblue','darkblue','lime','lime'] I = {} for i in range(subplots): if share: ax[i] = plt.subplot(gs00[i],sharex=ax[0]) else: ax[i] = plt.subplot(gs00[i]) I[i] = ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') ax[i].set_yticks([]) axvline_everywhere(ax, subplots, eod_fe, eod_fr) ax[0].set_title(win[g] +' Mean ') ax[1].set_title(win[g] +' Sum (Euclidean) ') plots = [distance_wo_mean_all[d_ms[d]],distance_wo_mean_all_in[d_ms[d]],distance_wo_mean_all_lm[d_ms[d]],distance_wo_mean_all_gm[d_ms[d]],distance_wo_mean_all_ds[d_ms[d]]] for i in np.arange(0,subplots-4,2): ax[i].plot(eod_fe - eod_fr, plots[int(i/2)], color=colors[i], linestyle='--') plots1 = [distance_wo_max_all[d_ms[d]],distance_wo_max_all_in[d_ms[d]],distance_wo_max_all_lm[d_ms[d]],distance_wo_max_all_gm[d_ms[d]],distance_wo_max_all_ds[d_ms[d]]] plots2 = [distance_wo_mean_max_all[d_ms[d]], distance_wo_mean_max_all_in[d_ms[d]], distance_wo_mean_max_all_lm[d_ms[d]], distance_wo_mean_max_all_gm[d_ms[d]], distance_wo_mean_max_all_ds[d_ms[d]]] for i in np.arange(0, subplots - 4, 2): ax[20+i] = ax[i].twinx() ax[20 +i].set_yticks([]) ax[20+i].plot(eod_fe - eod_fr, plots1[int(i/2)],color=colors[i],linestyle='dashed') ax[20+i].plot(eod_fe - eod_fr, plots2[int(i/2)],linestyle='-.', color=colors[i]) plt.subplots_adjust(hspace = 0.15, wspace = 0.15) #fig.legend(I,labels=labels) #fig.legend([I[1],I[3],I[5],I[7],I[9],I[11]],labels = ['convolved','indirect','lm','gm','ds','period']) if win[g] == 'w1': ax[0].set_ylabel('convolved',color = 'red', rotation = 360, labelpad = 40) ax[2].set_ylabel('perfect', color='magenta', rotation = 360, labelpad = 40) ax[4].set_ylabel('local max', color='green', rotation = 360, labelpad = 40) ax[6].set_ylabel('global max', color='black', rotation = 360, labelpad = 40) ax[8].set_ylabel('sampled f', color='orange', rotation = 360, labelpad = 40) ax[10].set_ylabel('period dist', color='darkblue', rotation = 360, labelpad = 40) ax[12].set_ylabel('period shift', color='lime', rotation=360, labelpad=40) ax[8].set_xticks([]) ax[10].set_xticks([]) if win[g] == 'w3': ax[0].set_ylabel('w1: pd and ps, mult', rotation=360, labelpad=60) ax[2].set_ylabel('w2: same pd', rotation=360, labelpad=60) ax[4].set_ylabel('w3: pd constant high f', rotation=360, labelpad=60) ax[6].set_ylabel('w4, same pd and ps', rotation=360, labelpad=60) if win[g] == 'w2': ax[8].set_xticks([]) ax[10].set_xticks([]) #ax = plt.gca() #leg = ax.get_legend() #leg.legendHandles[0].set_color('red') #leg.legendHandles[1].set_color('darkblue') def onclick(event): e = np.argmin(np.abs((eod_fe - eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t,sampling, deviation_s,d,a_fe, eod_fr, a_fr, eod_fe, e,phase_zero,p, size,s, sigma) plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) #plt.savefig( # fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( # size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( # eod_fr) + ')_beat_%1.2f' % ( # beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1]) + ')_deviation(%1.2f' % ( # d_ms[d]) + '_distances_all' + 'several' + '.png') if show: plt.show() def plot_dist_interactiv2_more(beat_corr, deviation_dp, delta_t, deviation_s, a_fe, a_fr, sigma,disc,win,period_shift,d,d_ms, period_distance, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share = False, show = False): fig = frame_subplots(12, 5, 'Beats [Hz]', '', 'Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (d_ms[d]) + ') '+win) two_ylabels(110,1,2,1,fig, 'distance (not normalised)', 'distance (normalised)') subplots = 14 ax = create_subpl(0,fig, subplots, 7, 2, share2 = share) plots = [distance_metz_all[d_ms[d]], euc_all[d_ms[d]],distance_metz_all_in[d_ms[d]],euc_all_in[d_ms[d]],distance_metz_all_lm[d_ms[d]], euc_all_lm[d_ms[d]], distance_metz_all_gm[d_ms[d]],euc_all_gm[d_ms[d]], distance_metz_all_ds[d_ms[d]], euc_all_ds[d_ms[d]],period_distance[d_ms[d]],period_distance[d_ms[d]],period_shift[d_ms[d]],period_shift[d_ms[d]]] labels = ['conv mean','conv sum', 'indirect mean', 'indirect sum', 'lm mean', 'lm sum', 'gm mean', 'gm sum','ds mean', 'ds sum','period lengths','','period shift',''] colors = ['red','red', 'magenta','magenta','green','green','black','black','orange','orange','darkblue','darkblue','lime','lime'] I = {} for i in range(subplots): I[i] = ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') axvline_everywhere(ax, subplots, eod_fe, eod_fr) ax[0].set_title('Mean (Mean instead sum)') ax[1].set_title('Sum (Euclidean distance)') plots = [distance_wo_mean_all[d_ms[d]],distance_wo_mean_all_in[d_ms[d]],distance_wo_mean_all_lm[d_ms[d]],distance_wo_mean_all_gm[d_ms[d]],distance_wo_mean_all_ds[d_ms[d]]] for i in np.arange(0,subplots-4,2): ax[i].plot(eod_fe - eod_fr, plots[int(i/2)], color=colors[i], linestyle='--') plots1 = [distance_wo_max_all[d_ms[d]],distance_wo_max_all_in[d_ms[d]],distance_wo_max_all_lm[d_ms[d]],distance_wo_max_all_gm[d_ms[d]],distance_wo_max_all_ds[d_ms[d]]] plots2 = [distance_wo_mean_max_all[d_ms[d]], distance_wo_mean_max_all_in[d_ms[d]], distance_wo_mean_max_all_lm[d_ms[d]], distance_wo_mean_max_all_gm[d_ms[d]], distance_wo_mean_max_all_ds[d_ms[d]]] for i in np.arange(0, subplots - 4, 2): ax[20+i] = ax[i].twinx() ax[20+i].plot(eod_fe - eod_fr, plots1[int(i/2)],color=colors[i],linestyle='dashed') ax[20+i].plot(eod_fe - eod_fr, plots2[int(i/2)],linestyle='-.', color=colors[i]) #fig.legend(I,labels=labels) #fig.legend([I[1],I[3],I[5],I[7],I[9],I[11]],labels = ['convolved','indirect','lm','gm','ds','period']) ax[0].set_ylabel('convolved',color = 'red',rotation = 360, labelpad = 40) ax[2].set_ylabel('indirect', color='magenta',rotation = 360, labelpad = 40) ax[4].set_ylabel('loc max', color='green',rotation = 360, labelpad = 40) ax[6].set_ylabel('glo max', color='black',rotation = 360, labelpad = 40) ax[8].set_ylabel('sampled f', color='orange',rotation = 360, labelpad = 40) ax[10].set_ylabel('period distance', color='darkblue',rotation = 360, labelpad = 40) ax[12].set_ylabel('period shift', color='lime', rotation=360, labelpad=40) #ax = plt.gca() #leg = ax.get_legend() #leg.legendHandles[0].set_color('red') #leg.legendHandles[1].set_color('darkblue') def onclick(event): e = np.argmin(np.abs((eod_fe - eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t,sampling, deviation_s,d,a_fe, eod_fr, a_fr, eod_fe, e,phase_zero,p, size,s, sigma) plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) plt.savefig( fname= disc+'Masterarbeit/labrotation/results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( eod_fr) + ')_beat_%1.2f' % ( beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1]) + ')_deviation(%1.2f' % ( d_ms[d]) + '_distances_all' + 'several' + '.png') if show: plt.show() def plot_dist_interactiv2(df_ds,corr_ds,beat_corr, deviation_dp, delta_t, deviation_s, a_fe, a_fr, sigma, disc, win, period_shift, d, d_ms, period_distance, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, beat_corr_org, distance_metz_all_df, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_corr, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share = False, show = False): #embed() fig = frame_subplots(7, 7, 'Beats [Hz]', '', 'Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (d_ms[d]) + ') '+win) two_ylabels(130,1,2,1,fig, 'distance (normalised)', 'distance (normalised)') plots = [distance_metz_all[d_ms[d]],distance_metz_all_df[d_ms[d]], df_ds[d_ms[d]],distance_metz_all_corr[d_ms[d]],corr_ds[d_ms[d]], distance_metz_all_gm[d_ms[d]], distance_metz_all_lm[d_ms[d]], period_distance[d_ms[d]], period_shift[d_ms[d]]] subplots = len(plots) ax = create_subpl(0,fig, subplots, subplots, 1, share2 = share) #plots = [distance_metz_all[d_ms[d]], euc_all[d_ms[d]], distance_metz_all_corr[d_ms[d]], euc_all_in[d_ms[d]], distance_metz_all_lm[d_ms[d]], euc_all_lm[d_ms[d]], distance_metz_all_gm[d_ms[d]], euc_all_gm[d_ms[d]], distance_metz_all_df[d_ms[d]], euc_all_ds[d_ms[d]], period_distance[d_ms[d]], period_distance[d_ms[d]], period_shift[d_ms[d]], period_shift[d_ms[d]]] labels = ['conv mean', 'df','df ds','corr mean','corr ds', 'gm mean','lm mean', 'period lengths',''] colors = ['black','magenta','pink','orange','moccasin','red','green','darkblue','lime'] I = {} for i in range(subplots): I[i] = ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') axvline_everywhere(ax, subplots, eod_fe, eod_fr) ax[0].set_title('Mean (Mean instead sum)') #ax[1].set_title('Sum (Euclidean distance)') #plots = [distance_wo_mean_all[d_ms[d]],distance_wo_mean_all_in[d_ms[d]],distance_wo_mean_all_lm[d_ms[d]],distance_wo_mean_all_gm[d_ms[d]],distance_wo_mean_all_ds[d_ms[d]]] #for i in np.arange(0,subplots,1): # ax[i].plot(eod_fe - eod_fr, plots[int(i/2)], color=colors[i], linestyle='--') #plots1 = [distance_wo_max_all[d_ms[d]],distance_wo_max_all_in[d_ms[d]],distance_wo_max_all_lm[d_ms[d]],distance_wo_max_all_gm[d_ms[d]],distance_wo_max_all_ds[d_ms[d]]] #plots2 = [distance_wo_mean_max_all[d_ms[d]], distance_wo_mean_max_all_in[d_ms[d]], distance_wo_mean_max_all_lm[d_ms[d]], distance_wo_mean_max_all_gm[d_ms[d]], # distance_wo_mean_max_all_ds[d_ms[d]]] #for i in np.arange(0, subplots - 4, 2): # ax[20+i] = ax[i].twinx() # ax[20+i].plot(eod_fe - eod_fr, plots1[int(i/2)],color=colors[i],linestyle='dashed') # ax[20+i].plot(eod_fe - eod_fr, plots2[int(i/2)],linestyle='-.', # color=colors[i]) #fig.legend(I,labels=labels) #fig.legend([I[1],I[3],I[5],I[7],I[9],I[11]],labels = ['convolved','indirect','lm','gm','ds','period']) ax[0].set_ylabel('convolved', color='black', rotation=360, labelpad=40) ax[1].set_ylabel('DF', color='magenta', rotation=360, labelpad=40) ax[2].set_ylabel('DF downsampled', color='pink', rotation=360, labelpad=40) ax[3].set_ylabel('beat corr', color='orange',rotation = 360, labelpad = 40) ax[4].set_ylabel('beat corr downsampled', color='moccasin', rotation=360, labelpad=40) ax[5].set_ylabel('global max', color='red', rotation=360, labelpad=40) ax[6].set_ylabel('local max', color='green', rotation=360, labelpad=40) ax[7].set_ylabel('segments length', color='darkblue',rotation = 360, labelpad = 40) ax[8].set_ylabel('timewindow', color='lime', rotation=360, labelpad=40) #ax = plt.gca() #leg = ax.get_legend() #leg.legendHandles[0].set_color('red') #leg.legendHandles[1].set_color('darkblue') plt.subplots_adjust(left = 0.3) def onclick(event): e = np.argmin(np.abs((eod_fe - eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t,sampling, deviation_s,d,a_fe, eod_fr, a_fr, eod_fe, e,phase_zero,p, size,s, sigma) plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( eod_fr) + ')_beat_%1.2f' % ( beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1]) + ')_deviation(%1.2f' % ( d_ms[d]) + '_distances_all' + 'several' + '.png') if show: plt.show() def onclick_effect(beat_corr, deviation_dp, delta_t,sampling, deviation_s,d,a_fe, eod_fr, a_fr, eod_fe, e,phase_zero,p, size,s, sigma): #embed() left_c = -200 * delta_t*sampling right_c = 200 * delta_t*sampling time, time_cut, cut = find_times(left_c, right_c, sampling, deviation_s[d]) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) eod_rec_down, eod_rec_up = rectify(eod_fish_r, eod_fe_chirp) # rectify eod_overlayed_chirp = (eod_fish_r + eod_fe_chirp)[cut:-cut] maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_rec_up[cut:-cut]) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_rec_up[cut:-cut]) # local maxima middle_conv, eod_conv_down, eod_conv_up, eod_conv_downsampled = conv(eod_fr,sampling, cut, deviation_dp[d], eod_rec_up, eod_rec_down) # convolve am_corr = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p], size[s], sigma) # indire# ct am calculation am_df = integrate_chirp(a_fe, time_cut, eod_fe[e]-eod_fr, phase_zero[p], size[s], sigma) # indire# ct am calculation _, time_fish, cut_f = find_times(left_c, right_c, eod_fr, deviation_s[d]) # downsampled through fish EOD am_df_ds = integrate_chirp(a_fe, time_fish, eod_fe[e]-eod_fr, phase_zero[p], size[s], sigma) # indire# ct am calculation am_corr_ds = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_corr_ds) / 2) plt.figure(figsize=[6.5, 5]) #conv plt.plot(time_cut * 1000, eod_conv_up + 3.3, color='black',label='convolved', linewidth=2) plt.plot(time_cut * 1000, am_df + 3.4, color='magenta', label='DF', linewidth=2) plt.plot( time_fish * 1000, am_df_ds +2.8, color='pink', label='DF downsampled', linewidth=2) plt.plot(time_cut * 1000, am_corr + 2.2, color='orange', label='Beat corr', linewidth=2) plt.plot( time_fish * 1000, am_corr_ds + 1.6, color='moccasin', label='Beat corr downsampled', linewidth=2) plt.plot(time_cut * 1000, eod_conv_down - 1.3, color='black', linewidth=2) #signal plt.plot(time_cut * 1000, eod_overlayed_chirp, label='EOD both fish', color='silver', linewidth=0.5) #AM maxima plt.scatter((index_peaks - 0.5 * len(eod_rec_up[cut:-cut])) / (sampling / 1000), value_peaks, color='green', s=10, linewidth=2) plt.plot((index_peaks - 0.5 * len(eod_rec_up[cut:-cut])) / (sampling / 1000), value_peaks, color='green', label='all maxima', linewidth=2) plt.scatter((maxima_index - 0.5 * len(eod_rec_up[cut:-cut])) / (sampling / 1000), maxima_values, color='red', s=10) plt.plot((maxima_index - 0.5 * len(eod_rec_up[cut:-cut])) / (sampling / 1000), maxima_values, color='red', label='maxima on lower frequency', linewidth=2) # beat corr # df plt.subplots_adjust(left = 0.4,hspace = 0.8) #embed() plt.title('Mult' + str(int(((eod_fe[e] - eod_fr)/eod_fe[e]+1)*100)/100)) plt.legend(loc=(-0.75,0.6), ncol=1) plt.xlim([-200, 200]) def plot_cons_interactiv2(deviation_dp, deviation_s, a_fe, a_fr, sigma, delta_t,type,disc,integral_ind,p, s, d, size, amplitude, phase_zero, beats, sampling, deviation,deviation_ms, eod_fe, eod_fr, beat_corr, max_values_ds, max_positions_ds, integral_ds, max_values_lm, max_positions_lm, integral_lm, max_values_gm, max_positions_gm, integral_gm, max_values_ind, max_positions_ind, integral_corr, max_values_conv, max_positions_conv, integral_conv, share = False): fig = frame_subplots(13, 6, 'Beats [Hz]', '', 'Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (deviation_ms[d]) + ')') two_ylabels(60,1,2,1,fig, 'Integral beneath the chirp Consp', '') two_ylabels(60,1, 2, 2, fig, '', 'Conspicousness [range 0 to 1]') subplots = 8 ax = create_subpl(0.3,fig, subplots, 4, 2, share2 = share) plots = [integral_conv, max_values_conv,integral_lm,max_values_lm, integral_gm, max_values_gm, integral_ind, max_values_ind] labels = ['Integral conv','Max Values conv','Integral lm','Max Values lm','Integral gm','Max Values gm','Integral ind','Max Values ind'] colors = ['red','red', 'green','green','black','black','magenta','magenta','orange','orange','darkblue','darkblue'] for i in range(subplots): ax[i].plot(eod_fe - eod_fr, plots[i], label=labels[i], color=colors[i], linestyle='-') axvline_everywhere(ax,subplots, eod_fe, eod_fr) ax[0].set_ylabel('convolved',color = 'red') ax[2].set_ylabel('loc max', color='green') ax[4].set_ylabel('glo max', color='black') ax[6].set_ylabel('perfect beat corr', color='magenta') ax[0].set_title('Integral') ax[1].set_title('Max Value') for i in np.arange(0,8,2): ax[i].set_ylim([0,int(20/3) * 3]) for i in np.arange(1, 8, 2): ax[i].set_ylim([0, 1.1]) #plt.legend(loc='upper left',bbox_to_anchor=(0, 2)) period = period_func(beat_corr, sampling,'stim') def onclick(event): e = np.argmin(np.abs((eod_fe-eod_fr) - event.xdata)) onclick_effect(beat_corr, deviation_dp,delta_t, sampling, deviation_s, d, a_fe, eod_fr, a_fr, eod_fe, e, phase_zero, p, size, s, sigma) rang = [0.5*sampling] time_fish,sampled, cut, cube_c, time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp = snip([-rang[0] * 0.5]*len(eod_fe), [rang[0] * 0.5]*len(eod_fe), e,e,sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, 0, size, 0, sigma, a_fe, deviation, beat_corr) B_conv, B_indirect,B_downsampled,B_lm,B_gm,time = snippets4(deviation_s,a_fr, sigma,period,time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp, eod_fe, beat_corr, a_fe, eod_fr, e, size, 0, phase_zero, 0, deviation, d, rang[0], delta_t, sampling, plot = False,show = False) #fig = plt.figure(figsize = [7,6]) gaussian = gauss(10,sampling,0.02,1) pl = True sh = False if event.inaxes == ax[4] or event.inaxes == ax[5]: max_position, max_value, I,ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'black',gaussian,B_gm, gm_chirp, sampling,'stim',plot = pl,show = sh) if event.inaxes == ax[2] or event.inaxes == ax[3]: max_position, max_value, I,ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'green',gaussian, B_lm, lm_chirp, sampling,'stim',plot = pl,show = sh) if event.inaxes == ax[0] or event.inaxes == ax[1]: max_position, max_value, I,ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'red',gaussian, B_conv, conv_chirp, sampling,'stim',plot = pl,show = sh) if event.inaxes == ax[6] or event.inaxes == ax[7]: max_position, max_value, I,ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'magenta',gaussian,B_indirect, am_indirect_chirp, sampling,'stim',plot = pl,show = sh) plt.xlabel('Time [ms]') plt.show() if share == False: cid = fig.canvas.mpl_connect('button_press_event', onclick) plt.savefig( fname='../results/' + 'simulations/' +str(type[0]) +'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( eod_fr) + ')_beat_%1.2f' % ( beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1])+ ')_deviation(%1.2f' % ( deviation_ms[d]) + '_consp' + 'several' + '.png') plt.show() def plot_dist_interactiv(distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, maxima_index_all, maxima_all, am_indirect_all, am_fish_all, time_fish_all, eod_convolved_down_all, eod_convolved_up_all, value_peaks_all, p, s, d, size, amplitude, phase_zero, beats, sampling, eod_rectified_up_all, index_peaks_all, time, eod_overlayed_chirp_all, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, interact2 = False): fig = plt.figure(figsize=[13, 6]) #plt.axes(frameon=False) ax = fig.add_subplot(1, 1, 1) plt.xticks([]) plt.yticks([]) ax.set_xticks([]) ax.set_yticks([]) ax0 = ax.twinx() ax0.set_xticks([]) ax0.set_yticks([]) ax1 = fig.add_subplot(5, 1, 1) ax2 = fig.add_subplot(5, 1, 4) ax3 = fig.add_subplot(5, 1, 3) ax4 = fig.add_subplot(5, 1, 2) ax5 = fig.add_subplot(5, 1, 5) ax1.plot(eod_fe - eod_fr, distance_metz_all, label='Metzen', color='red',linestyle=':') #ax1.set_xticks([]) ax1.plot(eod_fe - eod_fr, distance_wo_mean_all, label='Metzen ohne Mean', color='red',linestyle='-') ax11 = ax1.twinx() ax11.plot(eod_fe - eod_fr, distance_wo_max_all, label='Metzen ohne Normierung', color='red',linestyle='--') ax11.plot(eod_fe - eod_fr, distance_wo_mean_max_all, label='Metzen ohne Mean und Normierung',linestyle='-.', color='red') ax2.plot(eod_fe - eod_fr, distance_metz_all_in, label='Metzen',linestyle=':',color ='magenta') #ax2.set_xticks([]) ax2.plot(eod_fe - eod_fr, distance_wo_mean_all_in, label='Metzen ohne Mean',linestyle='-',color ='magenta') ax22 = ax2.twinx() ax22.plot(eod_fe - eod_fr, distance_wo_max_all_in, label='Metzen ohne Normierung',linestyle='--',color ='magenta') ax22.plot(eod_fe - eod_fr, distance_wo_mean_max_all_in, label='Metzen ohne Mean und Normierung',linestyle='-.',color ='magenta') ax3.plot(eod_fe - eod_fr, distance_metz_all_lm, label='Metzen', color='green',linestyle=':') ax3.plot(eod_fe - eod_fr, distance_wo_mean_all_lm, label='Metzen ohne Mean', color='green',linestyle='-') ax33 = ax3.twinx() ax33.plot(eod_fe - eod_fr, distance_wo_max_all_lm, label='Metzen ohne Normierung', color='green',linestyle='--') ax33.plot(eod_fe - eod_fr, distance_wo_mean_max_all_lm, label='Metzen ohne Mean und Normierung',linestyle='-.', color='green') ax4.plot(eod_fe - eod_fr, distance_metz_all_gm, label='Metzen', color='black',linestyle=':') ax4.plot(eod_fe - eod_fr, distance_wo_mean_all_gm, label='Metzen ohne Mean', color='black',linestyle='-') ax44 = ax4.twinx() ax44.plot(eod_fe - eod_fr, distance_wo_max_all_gm, label='Metzen ohne Normierung', color='black',linestyle='--') ax44.plot(eod_fe - eod_fr, distance_wo_mean_max_all_gm, label='Metzen ohne Mean und Normierung', color='black',linestyle='-.') ax5.plot(eod_fe - eod_fr, distance_metz_all_ds, label='Metzen', color='orange',linestyle=':') ax5.plot(eod_fe - eod_fr, distance_wo_mean_all_ds, label='Metzen ohne Mean', color='orange',linestyle='-') ax55 = ax5.twinx() ax55.plot(eod_fe - eod_fr, distance_wo_max_all_ds, label='Metzen ohne Normierung', color='orange',linestyle='--') ax55.plot(eod_fe - eod_fr, distance_wo_mean_max_all_ds, label='Metzen ohne Mean und Normierung', color='orange',linestyle='-.') ax1.axvline(x=-eod_fr, color='black', linewidth=1, linestyle='-') for i in range(int(np.max(eod_fe) / np.max(eod_fr))): ax1.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') ax2.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') ax3.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') ax4.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') ax5.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') #ax2.axvline(x=-eod_fr, color='black', linewidth=1, linestyle='-') #for i in range(int(np.max(eod_fe) / np.max(eod_fr))): # ax1.legend(loc = 'upper left',bbox_to_anchor=(1.15, 1.6)) # ax11.legend(loc='upper right', bbox_to_anchor=(1.15, 1.6)) ax1.legend(loc='upper left',bbox_to_anchor=(0, 2)) ax11.legend(loc='upper right',bbox_to_anchor=(1.0, 2)) #ax2.legend(loc='upper left') #ax22.legend(loc='upper right',bbox_to_anchor=(1, 1)) #ax3.legend(loc='upper left') #ax33.legend(loc='upper right',bbox_to_anchor=(1, 1)) #ax4.legend(loc='upper left') #ax44.legend(loc='upper right',bbox_to_anchor=(1, 1)) ax1.set_title('Convolved: EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (deviation_ms[d]) + ')') ax2.set_title('Beat corr') ax3.set_title('Lobal maxima') ax4.set_title('Global maxima') ax5.set_title('Down sampled') #ax1.set_xlabel('beat [hz]') #ax1.set_ylabel('distance') #ax2.set_xlabel('beat [hz]') #ax2.set_ylabel('distance') ax5.set_xlabel('beat [hz]') #ax3.set_ylabel('distance') ax.set_ylabel('distance (not normalised)', labelpad = 40) ax0.set_ylabel('distance (normalised)', rotation = 270,labelpad = 40) plt.subplots_adjust(hspace = 0.65) def onclick(event): i = np.argmin(np.abs((eod_fe-eod_fr) - event.xdata)) plt.figure(figsize = [7,6]) plt.plot(time * 1000, eod_overlayed_chirp_all[i][3 * deviation_dp[d]:-3 * deviation_dp[d]], label='EOD both fish', color='blue', linewidth=1) #ax3.title('beat:' + str(beat) + 'Hz ' + 'rf:' + str(eod_fr) + ' ef:' + str(eod_fe[e])) plt.scatter((index_peaks_all[i] - 0.5 * len(eod_rectified_up_all[i])) / (sampling / 1000), value_peaks_all[i], color='green', s=10, linewidth=2) plt.plot((index_peaks_all[i] - 0.5 * len(eod_rectified_up_all[i])) / (sampling / 1000), value_peaks_all[i], color='green', label='all maxima', linewidth=2) plt.plot(time * 1000, eod_convolved_up_all[i], color='red', linewidth=2) plt.plot(time * 1000, eod_convolved_down_all[i], color='red', label='convolved', linewidth=2) plt.plot( time_fish_all[i][int(3 * deviation_dp[d] / factor):int(-3 * deviation_dp[d] / factor)] * 1000, am_fish_all[i][int(3 * deviation_dp[d] / factor):int(-3 * deviation_dp[d] / factor)] - 0.3, color='orange', label='indirect am - downgesampled', linewidth=2) plt.scatter((maxima_index_all[i] - 0.5 * len(eod_rectified_up_all[i])) / (sampling / 1000), maxima_all[i], color='black', s=10) plt.plot((maxima_index_all[i] - 0.5 * len(eod_rectified_up_all[i])) / (sampling / 1000), maxima_all[i], color='black', label='all maxima', linewidth=2) plt.plot(time * 1000, am_indirect_all[i] + 0.3, color='magenta', label='beat corr', linewidth=2) plt.title('beat'+str(eod_fe[i]-eod_fr)) plt.legend(loc = 'lower center', ncol = 3, bbox_to_anchor = (0.5,1.05)) plt.xlim([-400,400]) plt.show() if interact2: cid = fig.canvas.mpl_connect('button_press_event', onclick) plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % ( eod_fr) + ')_beat_%1.2f' % ( beats[0]) + '_%1.2f' % (abs(beats[1] - beats[0])) +'_%1.2f' % (beats[-1])+ ')_deviation(%1.2f' % ( deviation_ms[d]) + '_distances_all' + 'several' + '.png') plt.show() def plot_eod_phase(beat_corr, minus_bef ,plus_bef,row,col,a_fe, delta_t, a_fr, size,beat_corr_org,s,pp, amplitude, phase_zero, deviation_ms, eod_fe,deviation, eod_fr, sampling,factor, fs = [6.2992, 4], shift_phase = 0): beats = eod_fe - eod_fr #fig = plt.figure(figsize=[6.2992/3, 4]) fig = plt.figure(figsize=fs) ax = {} counter = -1 ax[0] = plt.subplot(row, col, 1) # int(np.ceil(np.sqrt(len(eod_fe)))) for p,p_nr in zip(pp,range(len(pp))): d = 0 for e in range(len(eod_fe)): counter += 1 #time, time_cut = find_times(time_range, sampling, deviation[d], 1) left_c = minus_bef * delta_t * sampling right_c = plus_bef * delta_t * sampling time, time_cut,cut = find_times(left_c, right_c, sampling, deviation[d]/(1000*sampling)) #embed() time_fish_both = time * 2 * np.pi * (eod_fr-eod_fe[e]) eod_fish_both = 0.05 * np.sin(time_fish_both) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe,time, eod_fr, a_fr, eod_fe, e) eod_fish_both = integrate_chirp(a_fe, time, eod_fe[e]-eod_fr, phase_zero[p]+shift_phase, size[s], sigma) eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified_down, eod_recitified_up = rectify(eod_fish_r, eod_fe_chirp) # rectify maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_recitified_up) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_recitified_up) # local maxima middle_conv, eod_convolved_down, eod_convolved_up,eod_conv_downsampled = conv(eod_fr,sampling, cut,deviation[d], eod_recitified_up, eod_rectified_down) # convolve left_c = -200 * delta_t * sampling right_c = 200 * delta_t * sampling _, time_fish,_ = find_times(left_c,right_c, eod_fr, deviation[d]) # downsampled through fish EOD am_fish = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_fish) / 2) print(beat_corr[e]) am_corr = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p]+shift_phase, size[s], sigma) # indirect am calculation if plus_bef < 0: green_true = False else: green_true = True plt.axvline(x=-7.5, color='black', linestyle='dotted', linewidth=1) plt.axvline(x=7.5, color='black', linestyle='dotted', linewidth=1) ax[counter] = pl_eods(row, col,eod_fish_both,cut, beats,beat_corr, peaks_interp, maxima_interp, maxima_index, maxima_values,fig, ax[0],counter, e, time_cut ,am_corr,eod_fe,eod_overlayed_chirp,deviation,d,eod_fr,sampling,value_peaks, time_fish,am_fish,factor,eod_convolved_down,index_peaks,eod_convolved_up,eod_recitified_up,green_true = green_true,add = 'simple') ax[counter].title.set_text('Beat:' + str(eod_fe[e] - eod_fr) + 'Hz'+' phase'+str(int(phase_zero[p]*100)/100)) #embed() xticks = 'off' yticks = 'off' plot_pos = col * row - col + 1 ax[counter].set_xlabel('Time [ms]', labelpad=5) ax[counter].set_ylabel('[mv]', labelpad=5) #if counter == 2: # ax[counter].set_xlabel('Time [ms]', labelpad=5) # ax[counter].set_ylabel('[mv]', labelpad=5) #else: # if xticks == 'off': # ax[counter].set_xticks([]) # if yticks == 'off': # ax[counter].set_yticks([]) #lower_left_label(e+1, col, row, 'Time [ms]', '[mv]',xticks = 'off',yticks = 'off',) #plt.xlabel('Time [ms]') plt.legend(loc=(-1, -1.3),ncol = 2) plt.subplots_adjust(wspace=.2, hspace=.7,bottom = 0.3) #embed() plt.savefig( fname='../results/' + 'thesis/ '+'phasepic'+str(p)+str(eod_fe[e])+'.png') plt.show() def plot_eod(beat_corr, minus_bef ,plus_bef,row,col,a_fe, delta_t, a_fr, size,beat_corr_org,s,p, amplitude, phase_zero, deviation_ms, eod_fe,deviation, eod_fr, sampling,factor,fs = [6.2992, 4], shift_phase = 0): beats = eod_fe - eod_fr for d in range(len(deviation)): fig = plt.figure(figsize=fs) ax = {} #shift_phase = 0 ax[d] = fig.add_subplot(row,col, 1)#int(np.ceil(np.sqrt(len(eod_fe)))) for e in range(len(eod_fe)): #time, time_cut = find_times(time_range, sampling, deviation[d], 1) left_c = minus_bef * sampling right_c = plus_bef * sampling time, time_cut,cut = find_times(left_c, right_c, sampling, deviation[d]/(1000*sampling)) #embed() time_fish_both = time * 2 * np.pi * (eod_fr-eod_fe[e]) eod_fish_both = 0.05 * np.sin(time_fish_both) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) time_test = period_fish_e * 2 * np.pi * eod_fe[e] eod_test = a_fe * np.sin(time_test) #embed() eod_fish_both = integrate_chirp(a_fe, time, eod_fe[e]-eod_fr, phase_zero[p]+shift_phase, size[s], sigma) eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified_down, eod_recitified_up = rectify(eod_fish_r, eod_fe_chirp) # rectify maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_recitified_up) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_recitified_up) # local maxima middle_conv, eod_convolved_down, eod_convolved_up,eod_conv_downsampled = conv(eod_fr,sampling, cut,deviation[d], eod_recitified_up, eod_rectified_down) # convolve left_c = -200 * delta_t * sampling right_c = 200 * delta_t * sampling _, time_fish,_ = find_times(left_c,right_c, eod_fr, deviation[d]) # downsampled through fish EOD am_fish = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_fish) / 2) print(beat_corr[e]) am_corr = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p]+shift_phase, size[s], sigma) # indirect am calculation if plus_bef < 0: green_true = False else: green_true = True plt.axvline(x=-7.5, color='black', linestyle='dotted', linewidth=1) plt.axvline(x=7.5, color='black', linestyle='dotted', linewidth=1) ax[e] = pl_eods(row, col,eod_fish_both,cut, beats,beat_corr, peaks_interp, maxima_interp, maxima_index, maxima_values,fig, ax[d],e, e, time_cut ,am_corr,eod_fe,eod_overlayed_chirp,deviation,d,eod_fr,sampling,value_peaks, time_fish,am_fish,factor,eod_convolved_down,index_peaks,eod_convolved_up,eod_recitified_up,add = False,green_true = green_true) #embed() xticks = 'off' yticks = 'off' plot_pos = col * row - col + 1 #if e+1 == plot_pos: # ax[e].set_xlabel('Time [ms]', labelpad=5) if e+1 in np.arange(row*col-col+1,row*col+1,1): ax[e].set_xlabel('Time [ms]', labelpad=5) if e + 1 in np.arange(1,row*col+1,col): ax[e].set_ylabel('[mv]', labelpad=5) ax[e].set_yticks([]) #else: # if xticks == 'off': # ax[e].set_xticks([]) # if yticks == 'off': # ax[e].set_yticks([]) #lower_left_label(e+1, col, row, 'Time [ms]', '[mv]',xticks = 'off',yticks = 'off',) #plt.xlabel('Time [ms]') plt.legend(loc=(-1, -0.7),ncol = 2) plt.subplots_adjust(wspace=.2, hspace=.4,bottom = 0.20) #embed() #scalebar = ScaleBar(50,'ms') #plt.gca().add_artist(scalebar) #embed() # plt.xlim([-200,200]) plt.savefig( fname='../highbeats_pdf/cell_simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_deviation(%1.2f' % ( deviation_ms[d]) + 'beats(%1.2f' % (beats[0]) + '_%1.2f' % (beats[-1]) + '_%1.2f' % ( abs(beats[1] - beats[0])) + ')_eods_several_range'+str(minus_bef)+'_'+str(plus_bef)+ '.png') plt.show() embed() def lower_left_label(number, ncol, nrow, x, y,xticks = 'off',yticks = 'off',labelpad = 5): plot_pos = ncol * nrow - ncol+1 if number == plot_pos: plt.xlabel(x,labelpad = labelpad) plt.ylabel(y,labelpad = labelpad) else: if xticks == 'off': plt.xticks([]) if yticks == 'off': plt.yticks([]) def plot_am(a_fe, delta_t, a_fr, size,beat_corr_org,s,p, amplitude, phase_zero, deviation_ms, eod_fe,deviation, eod_fr, sampling,factor): beats = eod_fe - eod_fr for d in range(len(deviation)): fig = plt.figure(figsize=[12, 6]) ax = {} ax = fig.add_subplot(int(np.ceil(np.sqrt(len(eod_fe)))), int(np.ceil(np.sqrt(len(eod_fe)))), 1) for e in range(len(eod_fe)): time, time_cut = find_times(time_range, sampling, deviation[d], 1) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified_down, eod_rec_up = rectify(eod_fish_r, eod_fe_chirp) # rectify maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_rec_up) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_rec_up) # local maxima middle_conv, eod_conv_down, eod_conv_up,eod_conv_downsampled = conv(eod_fr,sampling, deviation[d], eod_rec_up, eod_rectified_down) # convolve _, time_fish = find_times(time_range, eod_fr, deviation[d], 50) # downsampled through fish EOD am_fish = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_fish) / 2) am_indirect = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p], size[s], sigma) # indirect am calculation period_dist,period_dist,am_beat, maxima_interp_beat, peaks_interp_beat, maxima_interp_chirp, peaks_interp_chirp, am_indirect_beat, \ time_fish,sampled,cube, am_indirect_chirp, time_beat, conv_beat, conv_chirp_corr, am_chirp, time_chirp, conv_chirp = snippets3( 1.2, peaks_interp, maxima_interp, am_indirect, am_fish, d, period, middle_am, eod_convolved_up, middle_conv, delta_t, sampling, deviation, time, eod_conv_downsampled) sav, ax = plot_amp(eod_fr,time, ax, e, fig, 'indirect', eod_fe, time_chirp, time_beat, am_indirect_chirp, am_indirect_beat, am_indirect) plt.xlabel('Time [ms]') plt.legend(loc=(1, 0)) plt.subplots_adjust(wspace=.2, hspace=.7) # plt.xlim([-200,200]) plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_deviation(%1.2f' % ( deviation_ms[d]) + 'beats(%1.2f' % (beats[0]) + '_%1.2f' % (beats[-1]) + '_%1.2f' % ( abs(beats[1] - beats[0])) + ')_am.png') plt.show() def pl_eods(row,col,eod_fish_both, cut, beats, beat_corr, peaks_interp, maxima_interp, maxima_index, maxima, fig, ax, i, e, time, am_corr, eod_fe, eod_overlayed_chirp, deviation, d, eod_fr, sampling, value_peaks, time_fish, am_fish, factor, eod_convolved_down, index_peaks, eod_convolved_up, eod_rectified_up, add = False, lw = 1, add1 = False,share = False,green_true = True): if share == True: ax = fig.add_subplot(row, col, i + 1, sharex=ax, sharey=ax) else: ax = fig.add_subplot(row, col, i + 1) #embed() #ax = plt.subplot(1,1,1) ax.plot(time * 1000, eod_overlayed_chirp[cut:-cut], label='EOD both fish', color='silver', linewidth = lw) #ax.plot(time * 1000, eod_overlayed_chirp[3 * deviation[d]:-3 * deviation[d]], # label='EOD both fish', # color='blue', linewidth=2) if add == True: ax.title.set_text('Beat:' + str(eod_fe[e]-eod_fr) + 'Hz ' + 'rf:' + str(eod_fr) + ' ef:' + str(eod_fe[e])) elif add == 'simple': ax.title.set_text('Beat:' + str(eod_fe[e] - eod_fr) + 'Hz') else: ax.title.set_text('Beat:' + str(eod_fe[e]-eod_fr) + 'Hz, Mult:'+str(int(((eod_fe[e]-eod_fr)/eod_fr+1)*100)/100)) if eod_fe[e]-eod_fr <0: green_true = True black_true = False else: green_true = False black_true = True if green_true == True: #ax.scatter((index_peaks - 0.5 * len(eod_rectified_up[cut:-cut])) / (sampling / 1000), value_peaks, # color='green', s=10) ax.plot(time * 1000, peaks_interp[cut:-cut], color='red', label = 'all maxima',linewidth = lw) #embed() #ax.plot((index_peaks - 0.5 * len(eod_rectified_up)) / (sampling / 1000), value_peaks, # color='green', label='all maxima') #plt.axvline(x=-7.5, color='black', linestyle='dotted', linewidth=1) #plt.axvline(x=7.5, color='black', linestyle='dotted', linewidth=1) if add == True: ax.plot(time * 1000, eod_convolved_up, color='red',linewidth = lw) ax.plot(time * 1000, eod_convolved_down, color='red', label='convolved',linewidth = lw) if add1 == True: ax.plot(time_fish[int(3 * deviation[d] / factor):int(-3 * deviation[d] / factor)] * 1000, am_fish[int(3 * deviation[d] / factor):int(-3 * deviation[d] / factor)]+ 0.4, color='purple', label='indirect am - downgesampled',linewidth = lw) ax.plot(time * 1000, am_corr +1.55, color='orange', label='EOD adjusted beat', linewidth = lw) if add1 == True: ax.scatter((maxima_index - 0.5 * len(eod_rectified_up)) / (sampling / 1000), maxima, color='red', s=10) if black_true == True: ax.plot(time*1000, maxima_interp[cut:-cut], color='red', label= 'AM',linewidth = lw)#[int(3 * deviation[d]):int(-3 * deviation[d])] ax.plot(time*1000,eod_fish_both[cut:-cut]+ 2.10,color='magenta',label= 'Difference frequency', linewidth = lw) #ax.set_xlim([-(2.5*1/beat_corr)*1000,(2.5*1/beat_corr)*1000]) #ax.plot((maxima_index - 0.5 * len(eod_rectified_up)) / (sampling / 1000), maxima, # color='black', label='one maxima pro lower frequency') #ax.set_ylabel('Time [ms]') #ax.set_ylabel('[mV]') return ax def plot_amp(eod_fr,time ,ax, e, fig, sav, eod_fe,time_chirp,time_beat,conv_chirp,conv_beat, eod_convolved_up): ax = fig.add_subplot(int(np.ceil(np.sqrt(len(eod_fe)))), int(np.ceil(np.sqrt(len(eod_fe)))), e + 1,sharex = ax,sharey = ax) ax.plot(time * 1000, eod_convolved_up,color='red', linewidth=2,label = 'beat') ax.plot(time_chirp * 1000, conv_chirp, color='blue', linewidth=2, label = 'chirp') ax.plot(time_beat * 1000, conv_beat, color='green', linewidth=2,label = 'beat') ax.set_title('beat' + str(eod_fe[e]-eod_fr) + 'Hz') #ax.set_ylim([0.2, 1]) return sav,ax def plot_all_distances(d,s,p,show_figure,beats,amplitude,size,phase_zero,corr_all,deviation_ms,eod_fe,eod_fr,distance_metz_all,distance_wo_mean_all,distance_wo_max_all,distance_wo_mean_max_all,range_chirp_all,range_beat_all, version): # plt.title('deviation(%1.2f' % (deviation[d]) + ')_beat(_%1.2f' % (beats[0]) + '_%1.2f' % ( # beats[-1]) + '_%1.2f' % (abs(beats[1] - beats[0])) + ')recieverEOD(_%1.2f' % (eod_fr) + ')') fig = plt.figure(figsize=[10, 10]) ax1 = fig.add_subplot(3, 1, 1) ax2 = fig.add_subplot(3, 1, 2, sharex = ax1) ax3 = fig.add_subplot(3, 1, 3, sharex = ax1) ax1.plot(eod_fe - eod_fr, distance_metz_all, label='Metzen', color='orange') ax1.plot(eod_fe - eod_fr, distance_wo_mean_all, label='Metzen ohne Mean', color='blue') ax11 = ax1.twinx() ax11.plot(eod_fe - eod_fr, distance_wo_max_all, label='Metzen ohne Normierung', color='red') ax11.plot(eod_fe - eod_fr, distance_wo_mean_max_all, label='Metzen ohne Mean und Normierung', color='green') ax22 = ax2.twinx() ax2.plot(eod_fe - eod_fr, range_chirp_all, label='range chirps',linestyle='-') ax2.plot(eod_fe - eod_fr, range_beat_all, label='range beat',linestyle='--') ax1.axvline(x=-eod_fr, color='black', linewidth=1, linestyle='-') for i in range(int(np.max(eod_fe) / np.max(eod_fr))): ax1.axvline(x=eod_fr * i, color='black', linewidth=1, linestyle='-') # ax1.legend(loc = 'upper left',bbox_to_anchor=(1.15, 1.6)) # ax11.legend(loc='upper right', bbox_to_anchor=(1.15, 1.6)) ax1.legend(loc='upper left') ax11.legend(loc='upper right') ax1.set_title('EOD reciever ' + str(eod_fr) + '_deviation(%1.2f' % (deviation_ms[d]) + ')') ax2.legend() ax1.set_xlabel('beat [hz]') ax1.set_ylabel('distance') ax2.set_xlabel('Amplitude') ax2.set_ylabel('distance') ax3.plot(eod_fe - eod_fr, corr_all) ax3.set_xlabel('beat [hz]') ax3.set_ylabel('Correlation') plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')recieverEOD(_%1.2f' % (eod_fr) + '_%1.2f' % ( beats[-1]) + '_%1.2f' % (abs(beats[1] - beats[0])) + ')_deviation(%1.2f' % ( deviation_ms[d]) + ')_beat(_%1.2f' % (beats[0]) + '_distances_all'+version+'.png') if show_figure: plt.show() else: plt.close() def plot_distances(eod_fe,p, s,deviation_ms,d, distance,size,phase_zero,amplitude, show_figure = False): for e in range(len(eod_fe)): time, time_cut = find_times(time_range, sampling, deviation_dp[d], 1) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe,time, eod_fr, a_fr, eod_fe, e) eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified_down, eod_rec_up = rectify(eod_fish_r, eod_fe_chirp) # rectify maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_rec_up) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_rec_up) # local maxima middle_conv, eod_conv_down, eod_conv_up,eod_conv_downsampled = conv(eod_fr,sampling, deviation_dp[d], eod_rec_up, eod_rectified_down) # convolve _, time_fish = find_times(time_range, eod_fr, deviation_dp[d], 50) # downsampled through fish EOD am_fish = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_fish) / 2) am_indirect = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p], size[s], sigma) # indirect am calculation period_dist,am_beat, maxima_interp_beat, peaks_interp_beat, maxima_interp_chirp, peaks_interp_chirp, am_indirect_beat, \ cube,am_indirect_chirp, time_beat, conv_beat, conv_chirp_corr, am_chirp, time_chirp, conv_chirp = snippets3( 1.2, peaks_interp, maxima_interp, am_indirect, am_fish, d, period, middle_am, eod_convolved_up, middle_conv, delta_t, sampling, deviation_dp, time, eod_conv_downsampled) plt.plot(distance) plt.plot(conv_chirp) plt.plot(conv_beat) if show_figure: plt.show() else: plt.close() plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_deviation(%1.2f' % (deviation_ms[d])+')_beat(%5.2f' % ( beat) + ')_absolute_distance.png') def plot_snippets(time_chirp, time_beat, time_range,p, sigma, s,delta_t,period,beat,sampling, deviation,d, time, eod_convolved_up,conv_chirp,conv_beat,size,phase_zero,amplitude,show_figure = False): plt.plot(time* 1000, eod_convolved_up, label='Amplitude Modulation',color='red', linewidth=2) plt.plot(time_chirp * 1000, conv_chirp, label='Amplitude Modulation', color='blue', linewidth=2) plt.plot(time_beat * 1000, conv_beat, label='Amplitude Modulation', color='green', linewidth=2) if show_figure: plt.show() else: plt.close() plt.savefig( fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_deviation(%1.2f' % (deviation_ms[d])+'_beat(%5.2f' % ( beat) + ')_distances_snippets_single.png') def plot_subplots(eod_fe,sides,p,s,deviation,d, deviation_ms,sampling,size,phase_zero,amplitude, show_figure = False, chirp = True): # frequency calculation for e in range(len(eod_fe)): time, time_cut = find_times(time_range, sampling, deviation[d], 1) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe,time, eod_fr, a_fr, eod_fe, e) if chirp == True: eod_fe_chirp = integrate_chirp(a_fe, time, eod_fe[e], phase_zero[p], size[s], sigma) else: eod_fe_chirp = eod_fish_e eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified_down, eod_rec_up = rectify(eod_fish_r, eod_fe_chirp) # rectify maxima_values, maxima_index, maxima_interp = global_maxima(sampling, eod_fr, period_fish_e, period_fish_r, eod_rec_up) # global maxima index_peaks, value_peaks, peaks_interp = find_lm(eod_rec_up) # local maxima middle_conv, eod_conv_down, eod_conv_up,eod_conv_downsampled = conv(eod_fr,sampling, deviation[d], eod_rec_up, eod_rectified_down) # convolve _, time_fish = find_times(time_range, eod_fr, deviation[d], 50) # downsampled through fish EOD am_fish = integrate_chirp(a_fe, time_fish, beat_corr[e], phase_zero[p], size[s], sigma) middle_am = int(len(am_fish) / 2) am_indirect = integrate_chirp(a_fe, time_cut, beat_corr[e], phase_zero[p], size[s], sigma) # indirect am calculation chirp = size[s] * np.exp((-(time ** 2) / (2 * sigma ** 2))) frequency = eod_fr - eod_fe[e] + chirp fig = plt.figure(figsize=[10, 6]) ax = {} ax[0] = fig.add_subplot(3, 2, 1) for i in range(6): ax[i] = fig.add_subplot(3, 2, i + 1, sharex = ax[0]) ax[i].axvline(x=-sides, color='black', linewidth=1, linestyle='-') ax[i].axvline(x=sides, color='black', linewidth=1, linestyle='-') ax[0].plot(time * 1000, eod_fish_r[3 * deviation[d]:-3 * deviation[d]]) ax[0].title.set_text('EOD reciever fish') ax[2].plot(time * 1000, eod_fe_chirp[3 * deviation[d]:-3 * deviation[d]]) ax[2].title.set_text('EOD and chirp of emmiting fish') ax[4].plot(time * 1000, frequency) ax[4].title.set_text('Frequency') ax[1].plot(time* 1000, eod_overlayed_chirp[3 * deviation[d]:-3 * deviation[d]]) # ax[3].plot(time[3 * deviation[d]:-3 * deviation[d]] * 1000,eod_rectified[3 * deviation[d]:-3 * deviation[d]]) ax[1].title.set_text('EOD both fish + chirp') (ax[1]).plot(time * 1000, eod_convolved_up, color='red') (ax[1]).plot(time * 1000, eod_convolved_down, color='red') corr_up = np.array([eod_convolved_up - np.min(eod_convolved_up) - np.mean(eod_convolved_up)]) # calculate AM directly time_sam = np.arange(-100 * delta_t, 100 * delta_t, 1 / 500*50) I1 = [] for i in range(len(time_sam)): y1 = ((np.pi ** 0.5) / 2) * math.erf(time_sam[i] / sigma) - ((np.pi ** 0.5) / 2) * math.erf(-np.inf) I1.append(y1) phase = 2 * np.pi * beat * time_sam + 2 * np.pi * size[s] * sigma * np.array(I1) + phase_zero[p] # am calculation am = amplitude * np.cos(phase) (ax[3]).plot(time_sam[3 * deviation[d]:-3 * deviation[d]] * 1000, am[ 3 * deviation[d]: - 3 * deviation[d]], color='green') #ax33 = ax[3].twinx() #ax33.twinx().plot(time * 1000, eod_convolved_up, color='red') corr_down = np.array([eod_convolved_down - np.min(eod_convolved_down) - np.mean(eod_convolved_down)]) # (ax[5]).plot(time * 1000, eod_convolved_down[ # length_convolved + 3 * deviation[ # d]:-length_convolved - 3 *deviation[d]]-np.min(eod_convolved_down)-np.mean(eod_convolved_down), color = 'red') #time = np.arange(-5 * delta_t, 5 * delta_t, 1 / sampling) # calculate AM directly time_fish = np.arange(-100 * delta_t, 100 * delta_t, 1 / 500) I1 = [] for i in range(len(time_fish)): y1 = ((np.pi ** 0.5) / 2) * math.erf(time_fish[i] / sigma) - ((np.pi ** 0.5) / 2) * math.erf(-np.inf) I1.append(y1) phase = 2 * np.pi * beat * time_fish + 2 * np.pi * size[s] * sigma * np.array(I1) + phase_zero[p] # am calculation am = amplitude * np.cos(phase) #FIXME #(ax[5]).twinx().plot(time * 1000, eod_convolved_down, # color='red') (ax[5]).plot(time_fish[int(np.round((3 * deviation[d]) / sampling)):int( np.round((-3 * deviation[d]) / sampling))] * 1000, am[int( np.round((3 * deviation[d]) / sampling)):int(np.round((-3 * deviation[d]) / sampling))], color='green') ax[1].set_xlim([-60,60]) ax[2].set_xlim([-60, 60]) ax[3].set_xlim([-60, 60]) ax[4].set_xlim([-60, 60]) ax[5].set_xlim([-60, 60]) # eod_convolved2 = (eod_convolved[length_convolved + int(3.395 * deviation[d]): len( # eod_convolved) - length_convolved - int(3.395 * deviation[d])]) plt.subplots_adjust(wspace=.2, hspace=.7) plt.savefig(fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % ( size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_deviation(%1.2f' % (deviation_ms[d])+'_beat(%5.2f' % (beat) + ')_subplots.png') if show_figure: plt.show() else: plt.close() def stimulus_harmonics(size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = False): sampling = 40000 time = np.arange(-5 * delta_t, 5 * delta_t, 1 / sampling) for s in range(len(size)): for p in range(len(phase_zero)): for d in range(len(deviation)): save_div = deviation[d] / (sampling / 1000) if deviation[d] * 5 % 2: points = deviation[d] * 5 else: points = deviation[d] * 5 -1 gaussian = signal.gaussian(points, std=deviation[d], sym=True) gaussian_normalised = (gaussian * 100) / np.sum(gaussian) length_convolved = int(len(gaussian_normalised) / 2) #cut_convolved = np.arange(length_convolved + (int(3.395 * deviation[d])), # (len(time) - length_convolved - int(3.395 * deviation[d])), 1) fig = plt.figure(figsize=[5, 4]) ax = {} # plot only the maxima for i in range(1): ax[i] = fig.add_subplot(1, 1, i + 1) ax[i].axvline(x=-7.5, color='black', linewidth=1, linestyle='-') ax[i].axvline(x=7.5, color='black', linewidth=1, linestyle='-') for e in range(len(eod_fe)): eod_fish_r = a_fr*np.sin(time * 2 * np.pi * eod_fr) sigma = delta_t / math.sqrt((2 * math.log(10))) # width of the chirp I1 = [] for i in range(len(time)): y1 = ((np.pi ** 0.5) / 2) * math.erf(time[i] / sigma) - ((np.pi ** 0.5) / 2) * math.erf(-np.inf) I1.append(y1) phase1 = time * 2 * np.pi * eod_fe[e] + 2 * np.pi * size[s] * sigma * np.array(I1) + phase_zero[p] eod_fe_chirp = a_fe*np.sin(phase1) eod_overlayed_chirp = eod_fish_r + eod_fe_chirp eod_rectified = eod_fish_r + eod_fe_chirp amplitude = a_fe beat = eod_fe[e] - eod_fr # rectify and gaussian eod_rectified[np.where(eod_overlayed_chirp < 0)] = 0 # rectify eod_convolved = np.convolve(gaussian_normalised, eod_rectified) # deviation = 60 # 40 = (1/1000)/(1 / 40000): 40 datapoints for 1 milisecond: 40*5 = 200 # frequency calculation chirp = size[s] * np.exp((-(time**2)/(2*sigma**2))) frequency = beat + chirp #ax[0].plot(time[(length_convolved + (int(3.395 * deviation[d]))):(len(time) - length_convolved - int(3.395 * deviation[d]))]* 1000, eod_convolved[(length_convolved + (int(3.395 * deviation[d]))):(len(time) - length_convolved - int(3.395 * deviation[d]))],label=str(beat)) ax[0].plot(time[3*deviation[d]:-3*deviation[d]] * 1000, eod_convolved[length_convolved+3*deviation[d]:-length_convolved-3*deviation[d]], label=str(beat)) ax[0].title.set_text('deviation('+str(save_div)+')') ax[0].legend(loc = 'lower left') plt.subplots_adjust(wspace=.2, hspace=.7) plt.savefig(fname='../results/' + 'simulations/' + 'amplitude(%1.2f' % (amplitude) + ')_size(%5.2f' % (size[s]) + ')_phase(%1.2f' % (phase_zero[p]) + ')_beat(%5.2f' % (beat) + ')_AM_stimulus_illustration3.png') if show_figure: plt.show() else: plt.close() def snippets4(deviation_s,a_fr, sigma, period, time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp, eod_fe, beat_corr, a_fe,eod_fr, e,size, ss, phase_zero, p, deviation, d, window, delta_t, sampling, plot = False,show = False): shifts_b = 100 B_conv = [[]]*shifts_b B_indirect = [[]]*shifts_b B_downsampled = [[]]*shifts_b B_lm = [[]]*shifts_b time = [[]] * shifts_b B_gm = [[]] * shifts_b ran = np.linspace(0, 3* period[e], shifts_b) shift = np.round(window * 2 + ran) for s in range(shifts_b): time_fish,sampled,cut, time[s], B_conv[s] , B_indirect[s], B_lm[s], B_gm[s] , B_downsampled[s] = snip(-window*0.5-shift, window*0.5-shift, e,s, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, ss, sigma, a_fe, deviation, beat_corr) if plot: fig = frame_subplots(6, 6, 'Time', 'AM', str(eod_fe[e]-eod_fr)+'Hz') ax = create_subpl(0.2, fig, 5, 3, 2, share2=True) ax[0].set_title('convovled') ax[0].plot(np.transpose(B_conv[0:100:20])) ax[1].set_title('global maximum') ax[1].plot(np.transpose(B_gm[0:100:20])) ax[2].set_title('local maximum') ax[2].plot(np.transpose(B_lm[0:100:20])) ax[3].set_title('downsampled') ax[3].plot(np.transpose(B_downsampled[0:100:20])) ax[4].set_title('indirect') ax[4].plot(np.transpose(B_indirect[0:100:20])) if show: plt.show() return B_conv, B_indirect,B_downsampled,B_lm,B_gm,time def conspicousy(disc,deviation_s,sigma, type,eod_fe, beats, deviation_ms, d, sampling, beat_corr, a_fe, delta_t, eod_fr, a_fr, size, phase_zero, deviation): max_values_gm = [[]]*len(eod_fe) max_positions_gm = [[]]*len(eod_fe) max_values_lm = [[]]*len(eod_fe) max_positions_lm = [[]]*len(eod_fe) max_values_conv = [[]]*len(eod_fe) max_positions_conv = [[]]*len(eod_fe) max_values_ind = [[]]*len(eod_fe) max_positions_ind = [[]]*len(eod_fe) max_values_ds = [[]]*len(eod_fe) max_positions_ds = [[]]*len(eod_fe) integral_ds = [[]]*len(eod_fe) integral_lm= [[]]*len(eod_fe) integral_gm= [[]]*len(eod_fe) integral_ind= [[]]*len(eod_fe) integral_conv= [[]]*len(eod_fe) rang = [0.4*sampling] max = [[]]*len(eod_fe) min = [[]] * len(eod_fe) period = period_func(beat_corr,sampling, 'stim') for e in range(len(eod_fe)): time_fish,sampled,cut, time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp = snip([-rang[0] * 0.5]*len(eod_fe), [rang[0] * 0.5]*len(eod_fe), e,e,sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, 0, size, 0, sigma, a_fe, deviation, beat_corr) B_conv, B_indirect,B_downsampled,B_lm,B_gm, time = snippets4(deviation_s,a_fr, sigma,period,time_c, conv_chirp, am_indirect_chirp, lm_chirp, gm_chirp, am_chirp, eod_fe, beat_corr, a_fe, eod_fr, e, size, 0, phase_zero, 0, deviation, d, rang[0], delta_t, sampling, plot = False,show = False) gaussian = gauss(10,sampling,0.02,1) max[e] = np.max(np.concatenate(time)) min[e] = np.min(np.concatenate(time)) pl = False sh = False max_positions_gm[e], max_values_gm[e], integral_gm[e],ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'black',gaussian,B_gm, gm_chirp, sampling,'stim',plot = pl, show = sh) max_positions_lm[e], max_values_lm[e], integral_lm[e],ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'green',gaussian,B_lm, lm_chirp, sampling,'stim',plot = pl,show = sh) max_positions_conv[e], max_values_conv[e], integral_conv[e],ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'red',gaussian,B_conv, conv_chirp, sampling,'stim',plot = pl,show = sh) max_positions_ind[e], max_values_ind[e], integral_ind[e],ConSpi2 = single_consp(type,0,1,eod_fe[e]-eod_fr,eod_fr,'magenta',gaussian,B_indirect, am_indirect_chirp, sampling,'stim',plot = pl,show = sh) #max_positions_ds[e], max_values_ds[e], integral_ds[e],ConSpi2 = single_consp(gaussian,B_downsampled, am_chirp, eod_fr) with open(disc +'Masterarbeit/labrotation/data/' + 'consp4_' +str(type[0])+ str( deviation_ms[d]) + 'ms' + 'beats(%1.2f' % (beats[0]) + '_%1.2f' % (beats[-1]) + '_%1.2f' % ( abs(beats[1] - beats[0])) + ')' + '.pkl', 'wb') as f: # Python 3: open(..., 'wb') pickle.dump( [eod_fe,beat_corr,beats,integral_ds, integral_lm, integral_gm, integral_ind, integral_conv, max_values_gm, max_positions_gm, max_values_lm, max_positions_lm, max_values_conv, max_positions_conv, max_values_ind, max_positions_ind, max_values_ds, max_positions_ds], f) return eod_fe,beat_corr,beats,integral_ds, integral_lm, integral_gm, integral_ind, integral_conv, max_values_gm,max_positions_gm,max_values_lm,max_positions_lm,max_values_conv,max_positions_conv,max_values_ind,max_positions_ind,max_values_ds,max_positions_ds def print_len(time): for i in range(len(time)): print(len(time[i])) def test_periods(): plt.plot(period) plt.plot(period <= delta_t * skretch) plt.plot(period_shift) plt.plot(np.ceil(delta_t * skretch / period) * period) plt.plot(np.ceil(delta_t * skretch / period)) a = np.ones(100) b = np.ones(100)*2 C1 = [] C2 = [] C3 = [] C4 = [] euc = [[]]*len(a) mean = [[]] * len(a) for c in range(len(a)): g = np.sqrt(np.mean((a[0:c]-np.mean(a[0:c]) - b[0:c] -np.mean(b[0:c]) ) ** 2)) C1.append(g) d = np.sqrt(np.mean((a[0:c] - b[0:c]) ** 2)) C2.append(d) k = np.sqrt(np.mean((a[0:c]-np.mean(a[0:c]) - b[0:c] -np.mean(b[0:c]) )** 2)) #C3.append(k) o = np.mean(np.sqrt(((a[0:c]-np.mean(a[0:c]) - b[0:c] -np.mean(b[0:c]) ) ** 2))) k = np.sqrt(np.sum((a[0:c]-np.mean(a[0:c]) - b[0:c] -np.mean(b[0:c]) )** 2)/len(a[0:c])) C3.append(k) o = np.sum(np.sqrt(((a[0:c]-np.mean(a[0:c]) - b[0:c] -np.mean(b[0:c]) ) ** 2)))/len(a[0:c]) C4.append(o) mean[c] = dm_wo_mean2_max = mean_euc_d(a[0:c], b[0:c], None) euc[c] = d_euclidean = euc_d(a[0:c], b[0:c], None) def find_beats(start,end,step,eod_fr): eod_fe = np.arange(start, end, step) beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] return eod_fe, beat_corr, beats def find_dev(x, sampling): deviation_ms = np.array(x) deviation_s = deviation_ms/1000 deviation_dp = sampling*deviation_s deviation_dp = list(map(int, deviation_dp)) return deviation_ms, deviation_s, deviation_dp def options(run, win,type): disc = 'D:/' delta_t = 0.014 # ms interest_interval = delta_t * 1.2 bef_c = interest_interval / 2 aft_c = interest_interval / 2 sigma = delta_t / math.sqrt((2 * math.log(10))) # width of the chirp size = [120] # maximal frequency excursion during chirp / 60 or 100 here phase_zero = [0] # phase when the chirp occured (vary later) / zero at the peak o a beat cycle eod_fr = 500 # eod fish reciever a_fr = 1 # amplitude fish reciever amplitude = a_fe = 0.2 # amplitude fish emitter factor = 200 sampling = eod_fr * factor sampling_fish = 500 x = [0.5, 1.5, 2.5] deviation_ms, deviation_s, deviation_dp = find_dev(x, sampling) time_range = 200 * delta_t if run == 'run_int': start = 0 end = 2500 step = 10 eod_fe, beat_corr, beats = find_beats(start,end,step,eod_fr) eod_fe,beat_corr,beats, integral_ds, integral_lm, integral_gm, integral_ind, integral_conv, max_values_gm,max_positions_gm,max_values_lm,max_positions_lm,max_values_conv,max_positions_conv,max_values_ind,max_positions_ind,max_values_ds,max_positions_ds = conspicousy(disc,deviation_s, sigma, type,eod_fe,beats, deviation_ms, 1, sampling, beat_corr, a_fe, delta_t, eod_fr, a_fr, size, phase_zero, deviation_dp) d = 1 plot_cons_interactiv2(deviation_dp, deviation_s, a_fe, a_fr, sigma, delta_t,type,disc,integral_ind,0, 0, d, size, a_fe, phase_zero, beats, sampling,deviation_dp, deviation_ms, eod_fe, eod_fr, beat_corr, max_values_ds, max_positions_ds, integral_ds, max_values_lm, max_positions_lm, integral_lm, max_values_gm, max_positions_gm, integral_gm, max_values_ind, max_positions_ind, integral_ind, max_values_conv, max_positions_conv, integral_conv, share = True) if run == 'show_int': with open(disc +'Masterarbeit/labrotation/data/consp4_1.5msbeats(-495.00_1995.00_10.00).pkl', 'rb') as f: # Python 3: open(..., 'rb') eod_fe,beat_corr,beats, integral_ds, integral_lm, integral_gm, integral_ind, integral_conv, max_values_gm, max_positions_gm, max_values_lm, max_positions_lm, max_values_conv, max_positions_conv, max_values_ind, max_positions_ind, max_values_ds, max_positions_ds = pickle.load(f) d = 1 plot_cons_interactiv2(deviation_dp, deviation_s, a_fe, a_fr, sigma, delta_t,type,disc,integral_ind,0, 0, d, size, a_fe, phase_zero, beats, sampling,deviation_dp, deviation_ms, eod_fe, eod_fr, beat_corr, max_values_ds, max_positions_ds, integral_ds, max_values_lm, max_positions_lm, integral_lm, max_values_gm, max_positions_gm, integral_gm, max_values_ind, max_positions_ind, integral_ind, max_values_conv, max_positions_conv, integral_conv, share = False) if run == 'run_dist': start = 5 end = 2500 step = 10 eod_fe, beat_corr, beats = find_beats(start,end,step,eod_fr) period_shift,period_distance_c,period_distance_b,euc_all_in,euc_all_gm,euc_all_lm,euc_all,euc_all_ds,range_beat_all_ds, range_chirp_all_ds, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm,beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds,distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm,distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm,distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm,distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling,deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in,distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr,distance_metz_all, distance_wo_mean_all, distance_wo_max_all,\ distance_wo_mean_max_all = distances_func(disc, bef_c, aft_c, win, deviation_s, sigma, sampling, deviation_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation_dp, show_figure = True, plot_dist = False, save = True) #plot_dist_interactiv2(period_shift,0,deviation_ms,period_distance_c,euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds,beat_corr_org,distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s,size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share = False,show = False) plot_dist_interactiv2(beat_corr, deviation_dp, delta_t, deviation_s, a_fe, a_fr, sigma,disc,win,period_shift,1,deviation_ms,period_distance_c,euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds,beat_corr_org,distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s,size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share = False,show = True) if run == 'show_dist': with open('../data/sim_dist3_1.5msbeats(-500.00_1995.00_5.00)'+win+'.pkl', 'rb') as f: # Python 3: open(..., 'rb') period_distance, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, range_beat_all_ds, range_chirp_all_ds, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm, beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds, \ distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, \ distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, \ distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, \ distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling, \ deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, \ distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, \ distance_metz_all, distance_wo_mean_all, distance_wo_max_all, \ distance_wo_mean_max_all = pickle.load(f) plot_dist_interactiv2(period_distance, euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds, beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s, d, size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share=True) if run == 'run_all_dist': start = 5 end = 2500 step = 10 eod_fe, beat_corr, beats = find_beats(start, end, step, eod_fr) plot_dist_interactiv4(bef_c, aft_c,disc,1, sampling, size, phase_zero, sigma, beat_corr, delta_t, a_fr, a_fe, deviation_dp, deviation_s, deviation_ms, eod_fe, eod_fr, share=True, show=True) def plot_beat_corr(eod_fr): eod_fe = np.arange(1, 1500, 5) beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] gs0 = gridspec.GridSpec(3, 1, height_ratios=[4, 1, 1], hspace=0.7) plt.figure(figsize=(4.5, 6)) style = 'dotted' color_v = 'black' color_b = 'silver' # plt.subplot(3,1,1) plt.subplot(gs0[0]) np.max(beats) / eod_fr plt.plot(beats, beats, color=color_b) # plt.axvline(x = -250, Linestyle = style,color = color_v) # plt.axvline(x = 250, Linestyle = style,color = color_v) # plt.axvline(x = 750, Linestyle = style,color = color_v) # plt.axvline(x = 1500, Linestyle = style) # plt.subplot(3,1,2) plt.xlabel('Beats [Hz]') plt.ylabel('Difference frequency [Hz]') plt.subplot(gs0[1]) plt.plot(beats, beat_corr, color=color_b) plt.axvline(x=-250, Linestyle=style, color=color_v) plt.axvline(x=250, Linestyle=style, color=color_v) plt.axvline(x=750, Linestyle=style, color=color_v) plt.xlabel('EOD adjusted beat [Hz]') # plt.axvline(x = 1250, Linestyle = style,color = color_v) # plt.axvline(x = 1500, Linestyle = style,color = color_v) mult = np.array(beats) / eod_fr + 1 # plt.subplot(3,1,3) plt.xlabel('Beats [Hz]') plt.ylabel('EOD adj. beat [Hz]', fontsize = 10) plt.subplot(gs0[2]) plt.plot(mult, beat_corr, color=color_b) # plt.axvline(x = 0, Linestyle = style) plt.axvline(x=0.5, Linestyle=style, color=color_v) # plt.axvline(x = 1, Linestyle = style) plt.axvline(x=1.5, Linestyle=style, color=color_v) plt.axvline(x=2.5, Linestyle=style, color=color_v) plt.xlabel('EOD multiples') plt.ylabel('EOD adj. beat [Hz]', fontsize = 10) plt.tight_layout() plt.subplots_adjust(left = 0.2,top = 0.97) plt.savefig(fname='../results/thesis/' + 'beat_corr_illustration') plt.show() embed() def plot_chirp(delta_t, sampling, deviation_ms, sigma, eod_fr, eod_fe,minus_bef = -3,plus_bef = +3): plt.figure(figsize=[3, 2]) left_c = minus_bef * delta_t * sampling right_c = plus_bef * delta_t * sampling time, time_cut, _ = find_times(left_c, right_c, sampling, deviation_ms[0] / (1000 * sampling)) chirp = 60 * np.exp((-(time ** 2) / (2 * sigma ** 2))) frequency = eod_fr - eod_fe[0] + chirp plt.plot(time * 1000, frequency, color='black') plt.ylabel('Frequency [Hz]') plt.xlabel('Time [ms]') plt.subplots_adjust(left=0.3, bottom=0.3) # plt.tight_layout() plt.savefig(fname='../results/thesis/chirpfigure') plt.show() def plot_power(beats, results, win, deviation_s, sigma, sampling, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = False, plot_dist = False, save = False,bef_c = 1.1,aft_c =-0.1): plt.rcParams['figure.figsize'] = (12, 7) plt.rcParams["legend.frameon"] = False colors = ['black','blue','purple','magenta','pink','orange', 'brown', 'red', 'green','lime'] fig, ax = plt.subplots(nrows=len(results), ncols=3, sharex=True) all_max = [[]] * len(results) #embed() for i in range(len(results)): ax[i, 0].set_ylabel(results['type'][i], rotation=0, labelpad=40, color=colors[i]) ax[i, 0].plot(beats / eod_fr + 1, np.array(results['result_frequency'][i]) / eod_fr + 1, color=colors[i]) # plt.title(results['type'][i]) ax[i, 1].plot(beats / eod_fr + 1, np.array(results['result_amplitude_max'][i]), color=colors[i]) ax[i, 2].plot(beats / eod_fr + 1, np.array(results['result_amplitude_max'][i]), color=colors[i]) ax[i,2].set_ylim([1,1.6]) all_max[i] = np.max(np.array(results['result_amplitude_max'][i])) for i in range(len(results)): ax[i, 2].set_ylim([0, np.max(all_max)]) plt.subplots_adjust(left=0.25) ii, jj = np.shape(ax) ax[0, 0].set_title('Most popular frequency') ax[0, 1].set_title('Modulation depth') ax[0, 2].set_title('Modulation depth (same scale)') for i in range(ii): for j in range(jj): ax[-1, j].set_xlabel('EOD multiples') ax[i, j].spines['right'].set_visible(False) ax[i, j].spines['top'].set_visible(False) def onclick(event): #embed() eod_fe = [(event.xdata-1)*eod_fr] nfft = 4096 e = 0 #sampling = 121000 left_b = [bef_c * sampling] * len(beat_corr) right_b = [aft_c * sampling] * len(beat_corr) p = 0 s = 0 d = 0 time_fish,sampled,cut, cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp = snip( left_b, right_b, e, e, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation, beat_corr, chirp=False) #cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp = onclick_func(e,nfft,beats, results, disc, win, deviation_s, sigma, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = show_figure, plot_dist = plot_dist, save = save,bef_c = bef_c,aft_c =aft_c, sampling = sampling) nfft = 4096 results = [[]] * 1 name = ['cubed'] var = [cubed] var = [maxima_interp_b] samp = [sampling] nfft = int(4096 * samp[0] / 10000 * 2) i = 0 pp, f = ml.psd(var[i] - np.mean(var[i]), Fs=samp[i], NFFT=nfft, noverlap=nfft / 2) plt.figure() plt.subplot(1,2,1) plt.plot(f, pp) plt.xlim([0,2000]) #plt.subplot(1,3,2) #plt.plot(time_b, cubed) plt.subplot(1,2,2) plt.plot(time_b, maxima_interp_b) plt.show() cid = fig.canvas.mpl_connect('button_press_event', onclick) plt.savefig('../highbeats_pdf/simulate_power') plt.show() def onclick_func(e,nfft, beats, results, disc, win, deviation_s, sigma, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = False, plot_dist = False, save = False,bef_c = 1.1,aft_c =-0.1,sampling = 100000): left_b = [bef_c * sampling] * len(beat_corr) right_b = [aft_c * sampling] * len(beat_corr) p = 0 s = 0 time_fish,sampled,cut, cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp = snip( left_b, right_b, e, e, sampling, deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation, beat_corr, chirp=False) # time_c, conv_c, am_corr_c, peaks_interp_c, maxima_interp_c,am_corr_ds_c, am_df_ds_c, am_df_c = snip(left_c, right_c,e,e, # sampling, deviation_s, # d, eod_fr, a_fr, eod_fe, # phase_zero, p, size, s, # sigma, a_fe, deviation, # beat_corr) # embed() # plt.plot(cubed) # embed() name = ['conv', 'df', 'df ds', 'corr', 'corr ds', 'global max', 'local max', 'cubed'] var = [conv_b, am_df_b, am_df_ds_b, am_corr_b, am_corr_ds_b, maxima_interp_b, peaks_interp_b, cubed] samp = [sampling, sampling, eod_fr, sampling, eod_fr, sampling, sampling, sampling] # pp, f = ml.psd(eod_overlayed_chirp - np.mean(eod_overlayed_chirp), Fs=sampling, NFFT=nfft, noverlap=nfft / 2) return pp, f, cubed, cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp def test_different_resolutions(): delta_t = 0.014 # ms interest_interval = delta_t * 1.2 bef_c = interest_interval / 2 aft_c = interest_interval / 2 sigma = delta_t / math.sqrt((2 * math.log(10))) # width of the chirp size = [120] # maximal frequency excursion during chirp / 60 or 100 here phase_zero = [0] # phase when the chirp occured (vary later) / zero at the peak o a beat cycle phase_zero = np.arange(0,2*np.pi,2*np.pi/10) eod_fr = 500 # eod fish reciever eod_fr = 537# eod fish reciever a_fr = 1 # amplitude fish reciever amplitude = a_fe = 0.2 # amplitude fish emitter factor = 200 sampling = eod_fr * factor sampling = 100000 sampling = 121730 sampling_fish = 500 #start = 0 #end = 2500 #step = 10 start = 510 end = 3500 end = eod_fr*5 step = 500 win = 'w2' d = 1 x = [ 1.5, 2.5,0.5,] x = [ 1.5] time_range = 200 * delta_t deviation_ms, deviation_s, deviation_dp = find_dev(x, sampling) start = 5 end = 3500 step = 250 eod_fe, beat_corr, beats = find_beats(start,end,step,eod_fr) bef_c = 1.1 aft_c = 1 bef_c = 0.02 aft_c = 0.01 bef_c = -2 aft_c = -1 samp = [10000,10000,40000,11000,11234,100000, 110000,1100000] fig = plt.figure() ax = {} nr = 4 ax['second'+str(-1)] = fig.add_subplot(len(samp), nr, 0 * nr + 2) ax['third'+str(-1)] = fig.add_subplot(len(samp), nr, 0 * nr + 3) ax['forth' + str(-1)] = fig.add_subplot(len(samp), nr, 0 * nr + 4) pp_max = [] for i in range(len(samp)): eod_fe = [eod_fr*0.75] eod_fe = [0] add = 1 nfft = 4096 e = 0 # sampling = 121000 #embed() left_b = [bef_c * samp[i]] * len(beat_corr) right_b = [aft_c * samp[i]] * len(beat_corr) p = 0 s = 0 d = 0 time_fish,sampled, cut, cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp = snip( left_b, right_b, e, e, samp[i], deviation_s, d, eod_fr, a_fr, eod_fe, phase_zero, p, size, s, sigma, a_fe, deviation_dp, beat_corr, chirp=False) # cubed, time_b, conv_b, am_corr_b, peaks_interp_b, maxima_interp_b, am_corr_ds_b, am_df_ds_b, am_df_b, eod_overlayed_chirp = onclick_func(e,nfft,beats, results, disc, win, deviation_s, sigma, d_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation, show_figure = show_figure, plot_dist = plot_dist, save = save,bef_c = bef_c,aft_c =aft_c, sampling = sampling) if single == True: eod_fe = [0] add = 1 time, time_cut, cut = find_times(left_c[e], right_c[e], sampling, deviation_s[d]) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) eod_overlayed_chirp = eod_fish_r+1 nfft = int(4096*samp[i]/samp[0]*2*2) results = [[]] * 1 name = ['cubed'] var = [cubed] var = cubed time1 = time_fish time2 = time_b time1 = time_b #samp = [sampling] #i = 0 pp, f = ml.psd(var - np.mean(var), Fs=samp[i], NFFT=nfft, noverlap=nfft / 2) #plt.figure() #plt.subplot(len(samp),3,i*3+1,sharex=True, sharey=True) #ax[i] = fig.add_subplot(len(samp),3,i*3+1, sharex=ax[0]) fig.add_subplot(len(samp), 4, i * 4 + 1) plt.plot(f, pp, label = samp[i]) #if i == 0: # pp_max = (np.max(pp)) plt.title(str(len(pp))) #plt.ylim([0,pp_max]) plt.xlim([0, 2000]) # plt.subplot(len(samp), 3, i*3+2, color = 'red',sharex=True, sharey=True) ax['second'+str(i)] = fig.add_subplot(len(samp), nr, i * nr + 2,sharex=ax['second'+str(i-1)],sharey=ax['second'+str(i-1)]) #embed() plt.plot(time1, var[cut:-cut], color='red')#[cut:-cut] # plt.subplot(len(samp), 3, i * 3 + 3,sharex=True, sharey=True) #ax['third'+str(i)] = fig.add_subplot(1, 1, 1) plt.xlim([time_b[0], -1.9]) ax['third'+str(i)] = fig.add_subplot(len(samp), nr, i * nr + 3,sharey=ax['third'+str(i-1)], sharex=ax['third'+str(i-1)])#sharey=ax['third'+str(i)], sharex=ax['third'+str(i)] plt.plot(time2, maxima_interp_b,color = 'red') pp, f = ml.psd(eod_overlayed_chirp - np.mean(eod_overlayed_chirp ), Fs=samp[i], NFFT=nfft, noverlap=nfft / 2) plt.plot(time2, eod_overlayed_chirp ) plt.xlim([time_b[0],-1.9]) ax['forth' + str(i)] = fig.add_subplot(len(samp), nr, i * nr + 4, sharex=ax['forth' + str(i - 1)], sharey=ax['forth' + str(i - 1)]) plt.plot(f, pp) plt.xlim([0, 2000]) # plt.subplot(1,3,2) # plt.plot(time_b, cubed) #plt.subplot(1, 2, 2) #plt.plot(time_b, maxima_interp_b) plt.legend() plt.show() def singlemodulation(beat_corr, deviation_s, a_fe, a_fr): bef_c = 0 aft_c = 1 sampling = 100000 i = 0 fig = plt.figure() ax = {} pp_max = [] add = 1 nfft = 4096 e = 0 left_c = [bef_c * sampling] * len(beat_corr) right_c = [aft_c * sampling] * len(beat_corr) p = 0 s = 0 d = 0 eod_fe = [0] eod_fr = 200 time, time_cut, cut = find_times(left_c[e], right_c[e], sampling, deviation_s[d]) eod_fish_e, eod_fish_r, period_fish_r, period_fish_e = find_periods(a_fe, time, eod_fr, a_fr, eod_fe, e) eod_overlayed_chirp = eod_fish_r+1 nfft = int(4096*2*2) ax[1] = fig.add_subplot(1,2,1)#sharey=ax['third'+str(i)], sharex=ax['third'+str(i)] plt.plot(time, eod_overlayed_chirp ) plt.xlim([0,0.05]) ax[2] = fig.add_subplot(1,2,2) pp, f = ml.psd(eod_overlayed_chirp, Fs=sampling, NFFT=nfft, noverlap=nfft / 2) plt.plot(f, pp) plt.xlim([-50, 2000]) plt.legend() plt.savefig('../highbeats_pdf/psd_singlemodulation.pdf') plt.show() if __name__ == "__main__": delta_t = 0.014 # ms interest_interval = delta_t * 1.2 bef_c = interest_interval / 2 aft_c = interest_interval / 2 sigma = delta_t / math.sqrt((2 * math.log(10))) # width of the chirp phase_zero = [0] # phase when the chirp occured (vary later) / zero at the peak o a beat cycle phase_zero = np.arange(0,2*np.pi,2*np.pi/10) eod_fr = 500 # eod fish reciever eod_fr = 537# eod fish reciever factor = 200 sampling_fish = 500 size = [120] d = 1 x = [ 1.5, 2.5,0.5,] x = [ 1.5] sampling = 100000 time_range = 200 * delta_t deviation_ms, deviation_s, deviation_dp = find_dev(x, sampling) #sampling = 750 #sampling = 522300 bef_c = -2 aft_c = -1 deviation_ms, deviation_s, deviation_dp = find_dev(x, sampling) #minus_bef = bef_c/(delta*sampling)#minus_bef = -30 #plus_bef = aft_c/(delta*sampling) minus_bef = bef_c plus_bef = aft_c #plus_bef = -10 eod_fe = np.array([510, 560,850,1010]) eod_fr = 750 eod_fe = np.array([eod_fr + 10, eod_fr + 100,eod_fr + eod_fr/2,eod_fr + eod_fr/2+100,eod_fr*3 + 10,eod_fr*4 + 10]) a_fe = 0.2 a_fr = 1 beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] plot_eod(beat_corr, minus_bef, plus_bef,3,2,a_fe, delta_t, a_fr, size, beat_corr, 0, 0, a_fe, phase_zero, deviation_ms, eod_fe, deviation_dp, eod_fr, sampling, factor,fs = [6.2992, 5],shift_phase = (2 * np.pi) / 4) embed() start = 5 end = 3500 step = 250 eod_fe, beat_corr, beats = find_beats(start,end,step,eod_fr) singlemodulation(beat_corr, deviation_s, a_fe, a_fr) embed() delta_t = 1 results = power_func(a_fr = 1, a_fe = 0.2, eod_fr = eod_fr, eod_fe = eod_fe, win = 'w2', deviation_s = deviation_s, sigma = sigma, sampling = 100000, deviation_ms = deviation_ms, beat_corr = beat_corr, size = [120],phase_zero = [phase_zero[0]], delta_t = delta_t,deviation_dp = deviation_dp, show_figure = True, plot_dist = False, save = False,bef_c = -2,aft_c = -1) plot_power(beats, results,'w2', deviation_s, sigma, sampling, deviation_ms, beat_corr, size, [phase_zero[0]], delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation_dp, show_figure = True, plot_dist = False, save = True,bef_c = bef_c,aft_c =aft_c) embed() start = 0 end = 2500 step = 25 eod_fe, beat_corr, beats = find_beats(start,end,step,eod_fr) distance_metz_all_df_ds,distance_metz_all_corr_ds,period_shift,period_distance_c,period_distance_b,euc_all_in,euc_all_gm,euc_all_lm,euc_all,euc_all_ds,range_beat_all_ds, range_chirp_all_ds, range_chirp_all, range_beat_all, range_chirp_all_in, range_beat_all_in, range_beat_all_gm, range_chirp_all_gm, range_beat_all_lm, range_chirp_all_lm,beat_corr_org, distance_metz_all_ds, distance_wo_mean_all_ds,distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm,distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm,distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm,distance_wo_max_all_gm, p, s, d, size, phase_zero, beats, sampling,deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in,distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr,distance_metz_all, distance_wo_mean_all, distance_wo_max_all,\ distance_wo_mean_max_all = distances_func(disc, bef_c, aft_c, 'w2', deviation_s, sigma, sampling, deviation_ms, beat_corr, size, phase_zero, delta_t, a_fr, a_fe, eod_fr, eod_fe, deviation_dp, show_figure = True, plot_dist = False, save = True) plot_dist_interactiv2(distance_metz_all_df_ds,distance_metz_all_corr_ds,beat_corr, deviation_dp, delta_t, deviation_s, a_fe, a_fr, sigma,disc,win,period_shift,0,deviation_ms,period_distance_c,euc_all_in, euc_all_gm, euc_all_lm, euc_all, euc_all_ds,beat_corr_org,distance_metz_all_ds, distance_wo_mean_all_ds, distance_wo_mean_max_all_ds, distance_wo_max_all_ds, distance_metz_all_lm, distance_wo_mean_all_lm, distance_wo_mean_max_all_lm, distance_wo_max_all_lm, distance_metz_all_gm, distance_wo_mean_all_gm, distance_wo_mean_max_all_gm, distance_wo_max_all_gm, p, s,size, amplitude, phase_zero, beats, sampling, deviation_ms, distance_wo_mean_max_all_in, distance_wo_max_all_in, distance_metz_all_in, distance_wo_mean_all_in, eod_fe, eod_fr, distance_metz_all, distance_wo_mean_all, distance_wo_max_all, distance_wo_mean_max_all, share = False, show = True) embed() minus_bef = -30 plus_bef = -10 eod_fe = np.array([510, 560,850,1010]) beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] plot_eod(beat_corr, minus_bef, plus_bef,2,2,a_fe, delta_t, a_fr, size, beat_corr, 0, 0, a_fe, phase_zero, deviation_ms, eod_fe, deviation_dp, eod_fr, sampling, factor,fs = [6.2992, 5],shift_phase = (2 * np.pi) / 4) minus_bef = -4 plus_bef = +4 eod_fe = np.array([510, 580,750,1010,1250,1510]) beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] plot_eod(beat_corr, minus_bef, plus_bef,3,2,a_fe, delta_t, a_fr, size, beat_corr, 0, 0, a_fe, phase_zero, deviation_ms, eod_fe, deviation_dp, eod_fr, sampling, factor,fs = [6.2992, 7.5],shift_phase = (2 * np.pi) / 4) minus_bef = -4 plus_bef = +4 eod_fe = np.array([510]) eod_fe = np.array([510,750]) eod_fe = np.array([580])#[510] beats = eod_fe - eod_fr beat_corr = eod_fe % eod_fr beat_corr[beat_corr > eod_fr / 2] = eod_fr - beat_corr[beat_corr > eod_fr / 2] pp = [1,3,6,9]#[2,7] #embed() plot_eod_phase(beat_corr, minus_bef, plus_bef,2,2,a_fe, delta_t, a_fr, size, beat_corr, 0, pp, a_fe, phase_zero, deviation_ms, eod_fe, deviation_dp, eod_fr, sampling, factor,fs = [6.2992, 5],shift_phase = (2 * np.pi) / 4) embed()