From b18119106ff8023c06edbc7970eb3659ec2bbf59 Mon Sep 17 00:00:00 2001
From: Martin Killenberg <martin.killenberg@desy.de>
Date: Thu, 26 Oct 2023 16:11:06 +0200
Subject: [PATCH] feat: post analysis for RF cables

---
 Python_script/analysis.py | 70 +++++++++++++++++++++++++++++++++------
 1 file changed, 59 insertions(+), 11 deletions(-)

diff --git a/Python_script/analysis.py b/Python_script/analysis.py
index 2fa4aaa..94bbf77 100644
--- a/Python_script/analysis.py
+++ b/Python_script/analysis.py
@@ -3,7 +3,7 @@ import matplotlib.pyplot as plt
 import numpy as np
 from matplotlib import gridspec
 
-def extract_stable_data(datafile, measurement_set, reference_signal_names):
+def extract_stable_data(datafile, measurement_set, reference_signal_names, extra_signal_names):
     datapoint = {}
     # df is a pandas data frame
     df = pd.read_csv(datafile)
@@ -29,23 +29,32 @@ def extract_stable_data(datafile, measurement_set, reference_signal_names):
     datapoint['humidity_mean'] = humidities.mean()
     datapoint['humidity_var'] = humidities.var()
 
+    for extra_signal_name in extra_signal_names:
+        extra_signal_values = df.loc[(df['EQUILIBRIUM_INDICATOR'] == 31) & (df['SET_NAME'] == measurement_set), extra_signal_name]
+        datapoint[extra_signal_name] = extra_signal_values.mean()
+
     return datapoint
 
 
 # sweep_type is either 'temperature' or 'humidity'
-def plot_sweep(temperatures, humidities, basename, sweep_type, measurement_sets, reference_signal_names, normalise = [True, False]):
+def plot_sweep(temperatures, humidities, basename, sweep_type, measurement_sets, reference_signal_names,
+               analysis_config = {'type': 'rf_cable', 'normalise': [True, False], 'cable_length': 10,
+                                  'extra_signal_names' : ['RF_FREQUENCY']}):
     set_data = {}
     derivatives = {}
     for measurement_set in measurement_sets:
         set_data[measurement_set] = {'signal0_means': [], 'signal0_vars': [], 'signal1_means': [], 'signal1_vars': [],
                                      'x_data': []}
+        for external_signal_name in analysis_config['extra_signal_names']:
+            set_data[measurement_set][external_signal_name] = []
         derivatives[measurement_set] = {'signal0_deltas': [], 'signal1_deltas': [], 'x_deltas': []}
-
+        
     for temp, hum in zip(temperatures, humidities):
         datafile = basename+'_'+str(temp)+'deg_'+str(hum)+'rh.csv'
         print(datafile)
         for measurement_set in measurement_sets:
-            datapoint = extract_stable_data(datafile, measurement_set, reference_signal_names)
+            datapoint = extract_stable_data(datafile, measurement_set, reference_signal_names,
+                                            analysis_config['extra_signal_names'])
             if datapoint is None:
                 continue
 
@@ -62,11 +71,14 @@ def plot_sweep(temperatures, humidities, basename, sweep_type, measurement_sets,
             data['signal1_means'].append(datapoint['signal1_mean'])
             data['signal1_vars'].append(datapoint['signal1_var'])
 
+            for external_signal_name in analysis_config['extra_signal_names']:
+                data[external_signal_name].append(datapoint[external_signal_name])
+
     for measurement_set in measurement_sets:
         data = set_data[measurement_set]
-        if normalise[0]:
+        if analysis_config['normalise'][0]:
             data['signal0_means'] -= data['signal0_means'][0]
-        if normalise[1]:
+        if analysis_config['normalise'][1]:
             data['signal1_means'] -= data['signal1_means'][0]
 
     fig = plt.figure(figsize=(10, 15))
@@ -78,9 +90,9 @@ def plot_sweep(temperatures, humidities, basename, sweep_type, measurement_sets,
     fig.text(0.1, 0.9, upper_text_block)
     
     ax1 = fig.add_subplot(gs[1, 0])
-    ax2 = fig.add_subplot(gs[2, 0],sharex=ax1)
+    ax2 = fig.add_subplot(gs[2, 0], sharex=ax1)
     ax4 = fig.add_subplot(gs[4, 0])
-    ax5 = fig.add_subplot(gs[5, 0],sharex=ax4)
+    ax5 = fig.add_subplot(gs[5, 0], sharex=ax4)
     
     ax1.tick_params(bottom=True, top=True, direction='in')
     ax2.tick_params(bottom=True, top=True, direction='inout')
@@ -132,12 +144,48 @@ def plot_sweep(temperatures, humidities, basename, sweep_type, measurement_sets,
                                                 (data['x_data'][i+1]+data['x_data'][i]))
             deriv_data['x_deltas'].append((data['x_data'][i+1]+data['x_data'][i])/2)
 
-        ax4.scatter(derivatives[measurement_set]['x_deltas'], derivatives[measurement_set]['signal0_deltas'],marker='+', )
-        ax5.scatter(derivatives[measurement_set]['x_deltas'], derivatives[measurement_set]['signal1_deltas'],marker='+')
-
     ax4.set_ylabel('$\\Delta$ '+reference_signal_names[0]+' '+denominator_name)
     ax5.set_ylabel('$\\Delta$ '+reference_signal_names[1]+' '+denominator_name)
 
+    #######################################################################################################################
+    # RF cable post analysis: Normalise to cable length and convert phase to time
+    #######################################################################################################################
+    
+    if analysis_config['type'] == 'rf_cable':
+        for measurement_set in measurement_sets:
+            data = set_data[measurement_set]
+            deriv_data = derivatives[measurement_set]
+            phase_to_time = 1e15 / (analysis_config['cable_length'] * data['RF_FREQUENCY'][0] * 360.)
+            deriv_data['signal0_deltas'][:] = [x * phase_to_time for x in deriv_data['signal0_deltas']]
+            deriv_data['signal1_deltas'][:] = [x / analysis_config['cable_length'] for x in deriv_data['signal1_deltas']]
+
+    ax4.set_ylabel('$\\Delta$t [fs/m/K]')
+    ax5.set_ylabel('$\\Delta$A [dB/m/K]')
+
+
+    #######################################################################################################################
+    # Closeout
+    #######################################################################################################################
+    #manual scaling if y axis. Auto is broken with scatter.
+    ax4_mins = []
+    ax4_maxs = []
+    ax5_mins = []
+    ax5_maxs = []
+
+    for measurement_set in measurement_sets:
+        ax4.scatter(derivatives[measurement_set]['x_deltas'], derivatives[measurement_set]['signal0_deltas'], marker='+')
+        ax5.scatter(derivatives[measurement_set]['x_deltas'], derivatives[measurement_set]['signal1_deltas'], marker='+')
+        ax4_mins.append(min(derivatives[measurement_set]['signal0_deltas']))
+        ax4_maxs.append(max(derivatives[measurement_set]['signal0_deltas']))
+        ax5_mins.append(min(derivatives[measurement_set]['signal1_deltas']))
+        ax5_maxs.append(max(derivatives[measurement_set]['signal1_deltas']))
+
+    # set auto y range which does not work properly for small values in scatter
+    dy4 = (max(ax4_maxs) - min(ax4_mins))*0.1
+    ax4.set_ylim(min(ax4_mins)-dy4, max(ax4_maxs)+dy4)
+    dy5 = (max(ax5_maxs) - min(ax5_mins))*0.1
+    ax5.set_ylim(min(ax5_mins)-dy5, max(ax5_maxs)+dy5)
+
     fig.tight_layout()  # otherwise the legend is clipped
 
     fig.savefig(basename+measurement_set+'_analysis.pdf')
-- 
GitLab