Merge pull request #50 from jfsehuanes/csv_with_data

Csv with data

thunderfish/csvmaker.py

@@ -0,0 +1,43 @@
+import numpy as np
+
+
+def write_csv(filename, csv_header, data_matrix, format_str='%.3f'):
+
+    """
+    Write matrix to csv-file.
+
+    Parameters
+    ----------
+    filename: string
+        String with the path and filename of the csv-file to be produced.
+    csv_header: list or array with strings
+        List with the names of the columns of the csv-matrix.
+    data_matrix: n-d array
+        Matrix with the data to be converted into a csv-file.
+    format_str: string
+        Format string for writing data_matrix (decimals)
+    """
+
+    # Check if header has same row length as data_matrix
+    if len(csv_header) != len(data_matrix[0]):
+        raise ValueError('The length of the header does not match the length of the data matrix!')
+
+    with open(filename, 'wb') as fin:
+        fin.write(str.encode(','.join(csv_header)))
+        fin.write(str.encode('\n'))
+        for row_id in range(len(data_matrix)):
+
+            fin.write(str.encode(','.join([format_str % e for e in data_matrix[row_id]])))  # convert to strings!
+            fin.write(str.encode('\n'))
+    pass
+
+
+if __name__ == '__main__':
+    print("\nChecking csvmaker module ...")
+
+    header = ['fundamental frequency', 'dB']
+    data = [[1.2, 2.3], [3.4, 4.5], [5.6, 6.7]]
+    filename = 'csvmaker_testfile.csv'
+    write_csv(filename, header, data)
+
+    print('\ncsv_file created in %s' % filename)

thunderfish/harmonicgroups.py

@@ -755,6 +755,41 @@ def fundamental_freqs(group_lists):
     return fundamentals
 
 
+def fundamental_freqs_and_db(grouplist):
+
+    """
+    Extract the fundamental frequencies and their power in dB from lists of harmonic groups.
+
+    Parameters
+    ----------
+    grouplist: list of 2-D arrays or list of list of 2-D arrays
+            Lists of harmonic groups as returned by extract_fundamentals() and
+            harmonic_groups() with the element [0][0] of the harmonic groups being the fundamental frequency,
+            and element[0][1] being the corresponding power.
+
+    Returns
+    -------
+    eodf_db_matrix: 2-D array or list of 2-D arrays
+        Matrix with fundamental frequencies in first column and corresponding power in dB in second column.
+
+    """
+
+    if len(grouplist) == 0:
+        eodf_db_matrix = np.array([])
+    elif hasattr(grouplist[0][0][0], '__len__'):
+        eodf_db_matrix = []
+        for groups in grouplist:
+            f = [np.array([harmonic_group[0][0], harmonic_group[0][1]]) for harmonic_group in grouplist]
+            f[:, 1] = 10.0 * np.log10(f[:, 1])  # calculate decibel using 1 as P0
+            eodf_db_matrix.append(f)
+    else:
+        eodf_db_matrix = np.array([np.array([harmonic_group[0][0], harmonic_group[0][1]])
+                                   for harmonic_group in grouplist])
+        eodf_db_matrix[:, 1] = 10.0 * np.log10(eodf_db_matrix[:, 1])  # calculate decibel using 1 as P0
+
+    return eodf_db_matrix
+
+
 def add_psd_peak_detection_config(cfg, low_threshold=0.0, high_threshold=0.0,
                                   thresh_bins=100, noise_fac=6.0, peak_fac=0.5,
                                   max_peak_width_fac=3.5, min_peak_width=1.0):

thunderfish/thunderfish.py

@@ -19,9 +19,10 @@
 from .bestwindow import clip_amplitudes, best_window_indices
 from .checkpulse import check_pulse_width, check_pulse_psd
 from .powerspectrum import plot_decibel_psd, multi_resolution_psd
-from .harmonicgroups import harmonic_groups, harmonic_groups_args, psd_peak_detection_args
+from .harmonicgroups import harmonic_groups, harmonic_groups_args, psd_peak_detection_args, fundamental_freqs_and_db
 from .consistentfishes import consistent_fishes_psd_plot, consistent_fishes
 from .eodanalysis import eod_waveform_plot, eod_waveform
+from .csvmaker import write_csv
 
 
 def output_plot(audio_file, pulse_fish_width, pulse_fish_psd, EOD_count, median_IPI, inter_eod_intervals,
@@ -193,6 +194,8 @@ def thunderfish(audio_file, channel=0, output_folder='', verbosearg=0):
     if verbosearg is not None:
         verbose = verbosearg
 
+    outfilename = os.path.splitext(os.path.basename(sys.argv[1]))[0]
+
     # load data:
     raw_data, samplerate, unit = load_data(audio_file, channel)
     if len(raw_data) == 0:
@@ -212,37 +215,44 @@ def thunderfish(audio_file, channel=0, output_folder='', verbosearg=0):
     # pulse-type fish?
     pulse_fish_width, pta_value = check_pulse_width(data, samplerate)
 
-    # calculate powerspectrums with different frequency resolutions
+    # calculate powerspectra with different frequency resolutions:
     psd_data = multi_resolution_psd(data, samplerate, fresolution=[0.5, 2 * 0.5, 4 * 0.5])
 
-    # find the fishes in the different powerspectrums:
+    # find the fishes in the different powerspectra:
     fishlists = []
     for i in range(len(psd_data)):
         h_kwargs = psd_peak_detection_args(cfg)
         h_kwargs.update(harmonic_groups_args(cfg))
         fishlist = harmonic_groups(psd_data[i][1], psd_data[i][0], verbose, **h_kwargs)[0]
         fishlists.append(fishlist)
 
-    # find the psd_type
+    # find the psd_type:
     pulse_fish_psd, proportion = check_pulse_psd(psd_data[0][0], psd_data[0][1])
 
     # filter the different fishlists to get a fishlist with consistent fishes:
     if not pulse_fish_width and not pulse_fish_psd:
         filtered_fishlist = consistent_fishes(fishlists)
+
+        # write csv file with main EODF and corresponding power in dB of detected fishes:
+        csv_matrix = fundamental_freqs_and_db(filtered_fishlist)
+        csv_name = outfilename + '-eodfs.csv'
+        header = ['fundamental frequency (Hz)', 'power (dB)']
+        write_csv(csv_name, header, csv_matrix)
     else:
         filtered_fishlist = []
 
     # analyse eod waveform:
     mean_eod, std_eod, time, eod_times = eod_waveform(data, samplerate, th_factor=0.6)
+    header = ['time (ms)', 'mean', 'std']
+    write_csv(outfilename + '-eodwaveform.csv', header, np.column_stack((1000.0*time, mean_eod, std_eod)))
+
     period = np.mean(np.diff(eod_times))
 
-    # inter-peal interval
+    # analyze inter-peak intervals:
     inter_peak_intervals = np.diff(eod_times) # in sec
-
     lower_perc, upper_perc = np.percentile(inter_peak_intervals, [1, 100-1])
     inter_eod_intervals = inter_peak_intervals[(inter_peak_intervals > lower_perc) &
                                                (inter_peak_intervals < upper_perc)]
-
     median_IPI = np.median(inter_eod_intervals)
     std_IPI = np.std(inter_eod_intervals, ddof=1)
 
@@ -267,4 +277,4 @@ def main():
 
 if __name__ == '__main__':
     main()
-
+    print('\nThanks for using thunderfish! Your files have been analyzed!')