From 8f5c2f65e6aeaf622109e06caf7e780a137cb6d4 Mon Sep 17 00:00:00 2001
From: mbergmann <maximilian.bergmann@student.uni-tuebingen.de>
Date: Thu, 24 Oct 2024 09:17:25 +0200
Subject: [PATCH] finished basic tuning curve

---
 code/tuning_curve_max.py | 70 +++++++++++++++++++++++++++++-----------
 1 file changed, 52 insertions(+), 18 deletions(-)

diff --git a/code/tuning_curve_max.py b/code/tuning_curve_max.py
index bea6610..64e669e 100644
--- a/code/tuning_curve_max.py
+++ b/code/tuning_curve_max.py
@@ -10,11 +10,9 @@ import useful_functions as f
 
 
 
-# variables
-delta = 2.5 # radius for peak detection
 
 # all files we want to use
-files = glob.glob("../data/2024-10-16-af*.nix")
+files = glob.glob("../data/2024-10-*.nix")
 
 # get only the good and fair filepaths
 new_files = f.remove_poor(files)
@@ -22,6 +20,9 @@ new_files = f.remove_poor(files)
 
 # loop over all the good files
 for file in new_files:
+    
+    contrast_frequencies = []
+    contrast_powers = []
     # load a file
     dataset = rlx.Dataset(file)
     # extract sams
@@ -30,14 +31,40 @@ for file in new_files:
     stim_frequencies = np.zeros(len(sams))
     peak_powers = np.zeros_like(stim_frequencies)
     # loop over all sams
-    for i, sam in enumerate(sams):
-        # get sam frequency and stimuli
-        avg_dur, _, _, _, _, _, stim_frequency = f.sam_data(sam)
-        print(avg_dur)
-        if np.isnan(avg_dur):
+    # dictionary for the contrasts
+    contrast_sams = {20 : [],
+                     10 : [],
+                     5 : []}
+    # loop over all sams
+    for sam in sams:
+        # get the contrast
+        avg_dur, contrast, _, _, _, _, _ = f.sam_data(sam)
+        # check for valid trails
+        if np.isnan(contrast):
+            continue
+        elif sam.stimulus_count < 3: #aborted trials
+            continue
+        elif avg_dur < 1.7:
             continue
-        # use this to change lists basically and add the contrast somewhere
         else:
+            contrast = int(contrast) # get integer of contrast
+            # sort them accordingly
+            if contrast == 20:
+                contrast_sams[20].append(sam)
+            if contrast == 10:
+                contrast_sams[10].append(sam)
+            if contrast == 5:
+                contrast_sams[5].append(sam)
+            else:
+                continue
+    # loop over the contrasts
+    for key in contrast_sams:
+        stim_frequencies = np.zeros(len(contrast_sams[key]))
+        peak_powers = np.zeros_like(stim_frequencies)
+        
+        for i, sam in enumerate(contrast_sams[key]):
+            # get stimulus frequency and stimuli
+            _, _, _, _, _, _, stim_frequency = f.sam_data(sam)
             stimuli = sam.stimuli
             # lists for the power spectra
             frequencies = []
@@ -52,20 +79,27 @@ for file in new_files:
             #average over the stimuli
             sam_frequency = np.mean(frequencies, axis = 0)
             sam_power = np.mean(powers, axis = 0)
-            # detect and validate peaks
+            # detect peaks
             integral, surroundings, peak_power = f.calculate_integral(sam_frequency, 
                                                           sam_power, stim_frequency)
-            valid = f.valid_integrals(integral, surroundings, stim_frequency)
-            #if there is a peak get the power in the peak powers
-            if valid == True:
-                peak_powers[i] = peak_power
+            
+            peak_powers[i] = peak_power
             # add the current stimulus frequency
             stim_frequencies[i] = stim_frequency
+        
+        # replae zeros with NaN
+        peak_powers = np.where(peak_powers == 0, np.nan, peak_powers)
+        
+        contrast_frequencies.append(stim_frequencies)
+        contrast_powers.append(peak_powers)
             
-    # replae zeros with NaN
-    peak_powers = np.where(peak_powers == 0, np.nan, peak_powers)
-    
-plt.plot(stim_frequencies, peak_powers)
+    fig, ax = plt.subplots(layout = 'constrained')
+    ax.plot(contrast_frequencies[0], contrast_powers[0])
+    ax.plot(contrast_frequencies[1], contrast_powers[1])
+    ax.plot(contrast_frequencies[2], contrast_powers[2])
+    ax.set_xlabel('stimulus frequency [Hz]')
+    ax.set_ylabel(r' power [$\frac{\mathrm{mV^2}}{\mathrm{Hz}}$]')
+    ax.set_title(f"{file}")