[eodanalysis] improved waveform and spectrum plots

thunderfish/eodanalysis.py

@@ -333,7 +333,8 @@ def analyze_pulse(eod, eod_times, thresh_fac=0.01, min_dist=50.0e-6,
     peak_list = np.array([prange[np.argmax(np.abs(meod[prange,1]))]
                  for prange in np.split(peak_l, np.where(np.diff(meod[peak_l,1]>0.0)!=0)[0]+1)])
     # remove multiple peaks that are too close and too small:
-    ridx = [(k, k+1) for k in np.where((np.diff(meod[peak_list,0]) < min_dist) & (meod[peak_list[:-1],1]<0.1))]
+    ridx = [(k, k+1) for k in np.where(np.diff(meod[peak_list,0]) < min_dist)]
+    #  & (meod[peak_list[:-1],1]<0.1)
     peak_list = np.delete(peak_list, ridx)
     # find P1:
     p1i = np.where(peak_list == max_idx)[0][0]
@@ -368,6 +369,9 @@ def analyze_pulse(eod, eod_times, thresh_fac=0.01, min_dist=50.0e-6,
     props['max-amplitude'] = max_ampl
     props['min-amplitude'] = min_ampl
     props['p-p-amplitude'] = ppampl
+    props['peakfrequency'] = freqs[np.argmax(power)]
+    props['lowfreqattenuation5'] = decibel(np.mean(power[freqs<5.0])/np.max(power))
+    props['lowfreqattenuation50'] = decibel(np.mean(power[freqs<50.0])/np.max(power))
     props['flipped'] = flipped
     props['n'] = len(eod_times)
 
@@ -427,16 +431,20 @@ def eod_waveform_plot(eod_waveform, peaks, ax, unit=None,
     if peaks is not None and len(peaks)>0:
         maxa = np.max(peaks[:,2])
         for p in peaks:
-            ax.scatter(1000.0*p[1], p[2], s=80,
+            ax.scatter(1000.0*p[1], p[2], s=80, clip_on=False,
                        c=mkwargs['color'], edgecolors=mkwargs['color'])
             label = 'P%d' % p[0]
             if p[0] != 1:
                 label += '(%.0f%% @ %.0fus)' % (100.0*p[3], 1.0e6*p[1])
             va = 'bottom' if p[2] > 0.0 else 'top'
+            y = np.sign(p[2])*0.02*maxa
+            if p[0] == 1 or p[0] == 2:
+                va = 'bottom'
+                y = 0.0
             if p[1] >= 0.0:
-                ax.text(1000.0*p[1]+0.1, p[2]+np.sign(p[2])*0.05*maxa, label, ha='left', va=va)
+                ax.text(1000.0*p[1]+0.1, p[2]+y, label, ha='left', va=va)
             else:
-                ax.text(1000.0*p[1]-0.1, p[2]+np.sign(p[2])*0.05*maxa, label, ha='right', va=va)
+                ax.text(1000.0*p[1]-0.1, p[2]+y, label, ha='right', va=va)
     ax.set_xlabel('Time [msec]')
     if unit is not None and len(unit)>0:
         ax.set_ylabel('Amplitude [%s]' % unit)

thunderfish/thunderfish.py

@@ -82,11 +82,11 @@ def output_plot(base_name, pulse_fish, inter_eod_intervals,
     """
 
     fig = plt.figure(facecolor='white', figsize=(14., 10.))
-    ax1 = fig.add_axes([0.02, 0.9, 0.96, 0.1])  # title
-    ax2 = fig.add_axes([0.075, 0.05, 0.9, 0.1]) # whole trace
-    ax3 = fig.add_axes([0.075, 0.6, 0.7, 0.3])  # psd
-    ax4 = fig.add_axes([0.075, 0.2, 0.4, 0.3])  # mean eod
-    ax5 = fig.add_axes([0.575, 0.2, 0.4, 0.3])  # inter EOD histogram
+    ax1 = fig.add_axes([0.02, 0.9, 0.96, 0.1])   # title
+    ax2 = fig.add_axes([0.075, 0.06, 0.9, 0.09]) # whole trace
+    ax3 = fig.add_axes([0.075, 0.6, 0.7, 0.3])   # psd
+    ax4 = fig.add_axes([0.075, 0.2, 0.4, 0.3])   # mean eod
+    ax5 = fig.add_axes([0.575, 0.2, 0.4, 0.3])   # pusle spectrum
 
     # plot title
     wavetitle = ""
@@ -168,12 +168,23 @@ def output_plot(base_name, pulse_fish, inter_eod_intervals,
 
     if not usedax5 and len(eod_props) > 0 and eod_props[0]['type'] == 'pulse':
         usedax5 = True
+        box = mpatches.Rectangle((1,-60), 49, 60, linewidth=0, facecolor='black', alpha=0.1)
+        ax5.add_patch(box)
+        att = eod_props[0]['lowfreqattenuation50']
+        ax5.text(10.0, att+1.0, '%.0f dB' % att, ha='left', va='bottom')
+        box = mpatches.Rectangle((1,-60), 4, 60, linewidth=0, facecolor='black', alpha=0.1)
+        ax5.add_patch(box)
+        att = eod_props[0]['lowfreqattenuation5']
+        ax5.text(4.0, att+1.0, '%.0f dB' % att, ha='right', va='bottom')
         db = decibel(spec_data[0][:,1])
         smax = np.nanmax(db)
         ax5.plot(spec_data[0][:,0], db - smax, 'b', lw=3)
+        peakfreq = eod_props[0]['peakfrequency']
+        ax5.scatter([peakfreq], [0.0], c='b', edgecolors='b', s=80)
+        ax5.text(peakfreq*1.2, 1.0, '%.0f Hz' % peakfreq, va='bottom')
         ax5.set_xlim(1.0, 10000.0)
         ax5.set_xscale('log')
-        ax5.set_ylim(-60.0, 1.0)
+        ax5.set_ylim(-60.0, 2.0)
         ax5.set_xlabel('Frequency [Hz]')
         ax5.set_ylabel('Power [dB]')
         ax5.set_title('Single-pulse spectrum', fontsize=14, y=1.05)