[AGIPD] [PC] Summary notebook for PC processing

Description

This MR provides summary notebook for processing AGIPD PC constants, visualising all 16 modules, similar to summary for FF constants.

How Has This Been Tested?

Tested with runs 61-68 in /gpfs/exfel/exp/SPB/202331/p900376/raw/

To Do

• Comparison plots to previous version of constants.

Relevant Documents (optional)

Processing report: PC_SPB_DET_AGIPD1M1_mc202_05mhz_gs0_900376_r61_r68.pdf

Types of changes

• New feature (non-breaking change which adds functionality)

Checklist:

Reviewers

@ahmedk @jsztuk

notebooks/AGIPD/Chracterize_AGIPD_Gain_PC_Summary.ipynb

 %% Cell type:markdown id: tags:
 
 # Pulse Capacitor Characterisation Summary
 
 This notebook is used as a dependency notebook for a pulse capacitor characterisation to provide summary for all modules of the AGIPD.
 
 %% Cell type:code id: tags:
 
 ``` python
 in_folder = "/gpfs/exfel/exp/SPB/202331/p900376/raw/" # path to input data, required
 out_folder = "/gpfs/exfel/exp/SPB/202331/p900376/scratch/PC_test/202cells0.5MHz_gs0_20clk/" # path to output to, required
 metadata_folder = ""  # Directory containing calibration_metadata.yml when run by xfel-calibrate
-runs = [] # runs to use, required, range allowed
+runs = range(61,69) # runs to use, required, range allowed
 
 karabo_id_control = "SPB_IRU_AGIPD1M1"  # karabo-id for the control device e.g. "MID_EXP_AGIPD1M1", or "SPB_IRU_AGIPD1M1"
 karabo_id = "SPB_DET_AGIPD1M-1"
 ctrl_source_template = '{}/MDL/FPGA_COMP' # path to control information
 
 creation_time = "" # To overwrite the measured creation_time. Required Format: YYYY-MM-DD HR:MN:SC e.g. "2022-06-28 13:00:00"
 creation_date_offset = "00:00:00" # add an offset to creation date, e.g. to get different constants
+cal_db_timeout = 3000000 # timeout on caldb requests"
+cal_db_interface = "tcp://max-exfl-cal001:8015#8045"
 
 bias_voltage = -1  # detector bias voltage, negative values for auto-detection.
 mem_cells = -1  # number of memory cells used, negative values for auto-detection.
 acq_rate = -1.  # the detector acquisition rate, negative values for auto-detection.
 gain_setting = -1  # gain setting can have value 0 or 1, negative values for auto-detection.
 integration_time = -1  # integration time, negative values for auto-detection.
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 import warnings
 warnings.filterwarnings('ignore')
 
 import os
 from dateutil import parser
 from datetime import timedelta
 
 import h5py
 import matplotlib.pyplot as plt
 import matplotlib.gridspec as gridspec
 import matplotlib.cm as cm
 import numpy as np
 import tabulate
+import multiprocessing
 
 from cal_tools.agipdlib import AgipdCtrl
 from cal_tools.ana_tools import get_range
 from cal_tools.tools import (
     calcat_creation_time,
-    module_index_to_qm
+    module_index_to_qm,
+    get_from_db,
 )
 
 from extra_data import RunDirectory
 from extra_geom import AGIPD_1MGeometry, AGIPD_500K2GGeometry
 from IPython.display import Latex, display
 from XFELDetAna.plotting.simpleplot import simplePlot
 
+from iCalibrationDB import Conditions, Constants
+
 %matplotlib inline
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Evaluate creation time
 creation_time = calcat_creation_time(in_folder, runs[0], creation_time)
 offset = parser.parse(creation_date_offset)
 delta = timedelta(hours=offset.hour, minutes=offset.minute, seconds=offset.second)
 creation_time += delta
 print(f"Creation time: {creation_time}\n")
 
 # Get operation conditions
 ctrl_source = ctrl_source_template.format(karabo_id_control)
 run_folder = f'{in_folder}/r{runs[0]:04d}/'
 raw_dc = RunDirectory(run_folder)
 
 # Read operating conditions from AGIPD00 files
 instrument_src_mod = [
     s for s in list(raw_dc.all_sources) if "0CH" in s][0]
 ctrl_src = [
     s for s in list(raw_dc.all_sources) if ctrl_source in s][0]
 
 agipd_cond = AgipdCtrl(
     run_dc=raw_dc,
     image_src=instrument_src_mod,
     ctrl_src=ctrl_src,
     raise_error=False,  # to be able to process very old data without mosetting value
 )
 if mem_cells == -1:
     mem_cells = agipd_cond.get_num_cells()
 if mem_cells is None:
     raise ValueError(f"No raw images found in {run_folder}")
 if acq_rate == -1.:
     acq_rate = agipd_cond.get_acq_rate()
 if gain_setting == -1:
     gain_setting = agipd_cond.get_gain_setting(creation_time)
 if bias_voltage == -1:
     bias_voltage = agipd_cond.get_bias_voltage(karabo_id_control)
 if integration_time == -1:
     integration_time = agipd_cond.get_integration_time()
 
 # Evaluate detector instance for mapping
 instrument = karabo_id.split("_")[0]
 if instrument == "HED":
     nmods = 8
 else:
     nmods = 16
 
 print(f"Using {creation_time} as creation time\n")
 print(f"Operating conditions are:\n• Bias voltage: {bias_voltage}\n• Memory cells: {mem_cells}\n"
       f"• Acquisition rate: {acq_rate}\n• Gain setting: {gain_setting}\n• Integration time: {integration_time}\n"
 )
 ```
 
+%% Output
 
+    Reading creation_date from input files metadata `INDEX/timestamp`
+    Creation time: 2023-08-21 11:02:29.330613+00:00
+    
+    Using 2023-08-21 11:02:29.330613+00:00 as creation time
+    
+    Operating conditions are:
+    • Bias voltage: 299
+    • Memory cells: 202
+    • Acquisition rate: 0.5
+    • Gain setting: 0
+    • Integration time: 20
+    
 
 %% Cell type:code id: tags:
 
 ``` python
 keys = ['ml', 'bl', 'mh', 'bh', 'BadPixelsPC']
 
 fit_data = {}
 labels = {'m': 'Slope (m)',
           'b': 'Intercept (b)',
           'H-M': "High-to-Medium Ratio",
           'mask': 'Fraction of bad pixels over cells' }
 
 modules = []
 karabo_da = []
 for mod in range(nmods):
     qm = module_index_to_qm(mod)
     fit_data[mod] = {}
     constants_file = f'{out_folder}/agipd_pc_store_{"_".join([str(run) for run in runs])}_{mod}_{mod}.h5'
 
     if os.path.exists(constants_file):
         print(f'Data available for module {qm}')
         with h5py.File(constants_file, 'r') as hf:
             for key in keys:
                 fit_data[mod][key]= hf[f'/{qm}/{key}/0/data'][()].swapaxes(1,2)
 
         modules.append(mod)
         karabo_da.append(f"AGIPD{mod:02d}")
     else:
         print(f"No fit data available for module {qm}")
 ```
 
+%% Output
 
+    Data available for module Q1M1
+    Data available for module Q1M2
+    Data available for module Q1M3
+    Data available for module Q1M4
+    Data available for module Q2M1
+    Data available for module Q2M2
+    Data available for module Q2M3
+    Data available for module Q2M4
+    Data available for module Q3M1
+    Data available for module Q3M2
+    Data available for module Q3M3
+    Data available for module Q3M4
+    Data available for module Q4M1
+    Data available for module Q4M2
+    Data available for module Q4M3
+    Data available for module Q4M4
 
 %% Cell type:code id: tags:
 
 ``` python
-#Define AGIPD geometry
-if instrument == "HED":
-    geom = AGIPD_500K2GGeometry.from_origin()
-else:
-    geom = AGIPD_1MGeometry.from_quad_positions(quad_pos=[
-        (-525, 625),
-        (-550, -10),
-        (520, -160),
-        (542.5, 475),
-    ])
+condition = Conditions.Dark.AGIPD(
+    memory_cells=mem_cells,
+    bias_voltage=bias_voltage,
+    acquisition_rate=acq_rate,
+    gain_setting=gain_setting,
+    integration_time=integration_time
+)
 ```
 
+%% Cell type:code id: tags:
 
+``` python
+def retrieve_old_PC(mod):
+    dconst = getattr(Constants.AGIPD, 'SlopesPC')()
+    old_PC = get_from_db(karabo_id=karabo_id,
+            karabo_da=karabo_da[mod],
+            constant=dconst,
+            condition=condition,
+            empty_constant=None,
+            cal_db_interface=cal_db_interface,
+            creation_time=creation_time-timedelta(seconds=1) if creation_time else None,
+            strategy="pdu_prior_in_time",
+            verbosity=1,
+            timeout=cal_db_timeout)
+    return old_PC
+
+with multiprocessing.Pool(processes=len(modules)) as pool:
+    old_PC_consts = pool.map(retrieve_old_PC, range(len(modules)))
 ```
 
+%% Output
 
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
+    Retrieved SlopesPC with creation time: 2023-01-11 16:48:58+01:00
 
+%% Cell type:code id: tags:
 
+``` python
+old_PC_consts = {}
+for i, mod in enumerate(modules):
+    dconst = getattr(Constants.AGIPD, 'SlopesPC')()
+    old_PC = get_from_db(karabo_id=karabo_id,
+            karabo_da=karabo_da[i],
+            constant=dconst,
+            condition=condition,
+            empty_constant=None,
+            cal_db_interface=cal_db_interface,
+            creation_time=creation_time-timedelta(seconds=1) if creation_time else None,
+            strategy="pdu_prior_in_time",
+            verbosity=1,
+            timeout=cal_db_timeout)
+    old_PC_consts[mod] = old_PC
 ```
 
+%% Cell type:code id: tags:
 
+``` python
+old_PC_consts[0][0].shape
 ```
 
+%% Output
 
+    (11, 202, 128, 512)
 
+%% Cell type:code id: tags:
 
+``` python
+old_PC = np.full(())
+for mod in range(len(modules)):
+
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Create the arrays that will be used for figures.
 # Each array correponds to the data for all processed modules.
 
 pixel_range = [0,0,512,128]
 const_data = {}
+const_order = [0, 1, 3, 4]
 
-for key in keys:
+for (key, c) in zip(keys, const_order):
     const_data[key] = np.full((nmods, mem_cells, 512, 128), np.nan)
+    old_const[key] = np.full((nmods, mem_cells, 512, 128), np.nan)
     for i in range(nmods):
         if key in fit_data[i]:
             const_data[key][i,:,pixel_range[0]:pixel_range[2],
                                pixel_range[1]:pixel_range[3]] = fit_data[i][key]
+            old_const[key][i,:,pixel_range[0]:pixel_range[2],
+                               pixel_range[1]:pixel_range[3]] =
 ```
 
+%% Cell type:code id: tags:
 
+``` python
+const_data.keys()
+#DB const order [0, 1, 3, 4]
 ```
 
+%% Output
 
+    dict_keys(['ml', 'bl', 'mh', 'bh', 'BadPixelsPC'])
 
+%% Cell type:code id: tags:
 
+``` python
+#Define AGIPD geometry
+if instrument == "HED":
+    geom = AGIPD_500K2GGeometry.from_origin()
+else:
+    geom = AGIPD_1MGeometry.from_quad_positions(quad_pos=[
+        (-525, 625),
+        (-550, -10),
+        (520, -160),
+        (542.5, 475),
+    ])
 ```
 
 %% Cell type:markdown id: tags:
 
 ## Summary across pixels ##
 
 %% Cell type:code id: tags:
 
 ``` python
 gain_data = {'HG': {},
              'MG': {},
              'BP': {}
             }
 
 for key in const_data.keys():
     if key is 'ml' or key is 'bl':
         gain_data['HG'][key] = const_data[key]
     if key is 'mh' or key is 'bh':
         gain_data['MG'][key] = const_data[key]
     if key is 'BadPixelsPC':
         gain_data['BP'][key] = const_data[key]
 
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 g_label = ['High gain', 'Medium gain', 'Bad pixels PC']
 for idx, g in enumerate(gain_data.keys()):
     gs = gridspec.GridSpec(1, 2)
     fig = plt.figure(figsize=(15, 7))
     plt.suptitle(f'{g_label[idx]}', fontsize=16)
 
     for i, key in enumerate(gain_data[g].keys()):
         if key is 'BadPixelsPC':
             data = np.nanmean(gain_data[g][key]>0, axis=1)
             vmin, vmax = (0,1)
             ax = fig.add_subplot(gs[0, :])
 
         else:
             data = np.nanmean(gain_data[g][key], axis=1)
             vmin, vmax = get_range(data, 5)
             ax = fig.add_subplot(gs[0, i])
         geom.plot_data_fast(data,
                                 vmin=vmin, vmax=vmax, ax=ax, cmap="jet", figsize=(13,7),
                                 colorbar={'shrink': 1,
                                           'pad': 0.04,
                                           'fraction': 0.1
                                          })
         colorbar = ax.images[0].colorbar
         ax.set_title(key, fontsize=14)
         ax.set_xlabel('Columns', fontsize=12)
         ax.set_ylabel('Rows', fontsize=12)
 
     plt.show()
 ```
 
 %% Cell type:markdown id: tags:
 
 ## Summary across cells ##
 
 Good pixels only.
 
 %% Cell type:code id: tags:
 
 ``` python
 ratio = gain_data['HG']['ml'] / gain_data['MG']['mh']
 
 fig = plt.figure(figsize=(7, 7))
 ax = fig.add_subplot(111)
 data = np.nanmean(ratio, axis=1)
 vmin, vmax = get_range(data, 5)
 ax = geom.plot_data_fast(data,
                                 vmin=vmin, vmax=vmax, ax=ax, cmap="jet", figsize=(6,7),
                                 colorbar={'shrink': 1,
                                           'pad': 0.04,
                                           'fraction': 0.1
                                          })
 colorbar = ax.images[0].colorbar
 colorbar.set_label('HG slope / MG slope', fontsize=12)
 ax.set_title('High/Medium Gain Slope Ratio', fontsize=14)
 ax.set_xlabel('Columns', fontsize=12)
 ax.set_ylabel('Rows', fontsize=12)
 
 plt.show()
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 for idx, g in enumerate(gain_data.keys()):
 
     for key in gain_data[g].keys():
         data = np.copy(gain_data[g][key])
         if key=='BadPixelsPC':
             data = data>0
         else:
             data[gain_data['BP']['BadPixelsPC']>0] = np.nan
 
         d = []
         for i in range(nmods):
             d.append({'x': np.arange(data[i].shape[0]),
                       'y': np.nanmean(data[i], axis=(1,2)),
                       'drawstyle': 'steps-pre',
                       'label': f'{i}',
                       'linewidth': 2,
                       'linestyle': '--' if i>7 else '-'
                       })
 
         fig = plt.figure(figsize=(15, 6))
         plt.suptitle(f'{g_label[idx]} - {key}', fontsize=16)
         ax = fig.add_subplot(111)
 
         _ = simplePlot(d, xrange=(-12, 510),
                             x_label='Memory Cell ID',
                             y_label=key,
                             use_axis=ax,
                             legend='top-left-frame-ncol8',)
         ylim = ax.get_ylim()
         ax.set_ylim(ylim[0], ylim[1] + np.abs(ylim[1]-ylim[0])*0.2)
         ax.grid()
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 d = []
 for i in range(nmods):
     d.append({'x': np.arange(ratio[i].shape[0]),
               'y': np.nanmean(ratio[i], axis=(1,2)),
               'drawstyle': 'steps-pre',
               'label': f'{i}',
               'linewidth': 2,
               'linestyle': '--' if i>7 else '-'
               })
 
 fig = plt.figure(figsize=(15, 6))
 plt.suptitle('High/Medium Gain Slope Ratio', fontsize=16)
 ax = fig.add_subplot(111)
 
 _ = simplePlot(d, xrange=(-12, 510),
                     x_label='Memory Cell ID',
                     y_label=key,
                     use_axis=ax,
                     legend='top-left-frame-ncol8',)
 ylim = ax.get_ylim()
 ax.set_ylim(ylim[0], ylim[1] + np.abs(ylim[1]-ylim[0])*0.2)
 ax.grid()
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 table = []
 for mod in modules:
 
     table.append((mod,
                   f"{np.nanmean(ratio[mod]):0.1f} +- {np.nanstd(ratio[mod]):0.2f}",
                   f"{np.nanmean(gain_data['BP']['BadPixelsPC'][mod]>0)*100:0.1f} ({np.nansum(gain_data['BP']['BadPixelsPC'][mod]>0)})"
                 ))
 
 all_HM = []
 for mod in modules:
     all_HM.extend(ratio[mod])
 all_HM = np.array(all_HM)
 
 all_MSK = np.array([list(msk) for msk in gain_data['BP']['BadPixelsPC']])
 
 table.append(('overall',
               f"{np.nanmean(all_HM):0.1f} +- {np.nanstd(all_HM):0.2f}",
               f"{np.nanmean(all_MSK>0)*100:0.1f} ({np.nansum(all_MSK>0)})"
             ))
 
 md = display(Latex(tabulate.tabulate(table, tablefmt='latex',
                                      headers=["Module", "Ratio of high/medium gain slope","Bad pixels [%(Count)]"])))
 ```
 
 %% Cell type:markdown id: tags:
 
 ## Summary high-medium gain ratio (good pixels only) + histograms
 
 %% Cell type:code id: tags:
 
 ``` python
 if instrument == 'HED':
     n = 8
 else:
     n = 16
 
 colors = cm.rainbow(np.linspace(0, 1, n))
 
 gs = gridspec.GridSpec(1, 2)
 fig = plt.figure(figsize=(17, 8))
 
 
 ratio[gain_data['BP']['BadPixelsPC'] > 0] = np.nan
 data = np.nanmean(ratio, axis=1)
 vmin, vmax = get_range(data, 5)
 ax = fig.add_subplot(gs[0, 0])
 geom.plot_data_fast(data,
                         vmin=vmin, vmax=vmax, ax=ax, cmap="jet", figsize=(13,7),
                         colorbar={'shrink': 1,
                                   'pad': 0.04,
                                   'fraction': 0.1
                                  })
 colorbar = ax.images[0].colorbar
 colorbar.set_label('High/Medium Gain Ratio', fontsize=12)
 ax.set_xlabel('Columns', fontsize=12)
 ax.set_ylabel('Rows', fontsize=12)
 
 ax = fig.add_subplot(gs[0,1])
 for mod in modules:
     h, e = np.histogram(ratio[mod].flatten(), bins=100, range=(vmin, vmax))
     ax.plot(e[:-1], h, color=colors[mod],linewidth=2, label=f'{mod}', alpha=0.8)
     ax.set_xlabel('High/Medium Gain Ratio', fontsize=12)
     ax.set_ylabel('Counts', fontsize=12)
 ax.grid()
 ax.legend()
 plt.show()
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 ```