Skip to content
Snippets Groups Projects

[AGIPD] [PC] Summary notebook for PC processing

Merged [AGIPD] [PC] Summary notebook for PC processing
feat/agipd-pc-summary into master
Merged Vratko Rovensky requested to merge feat/agipd-pc-summary into master
Viewing commit 74473f9d
Prev Next
Show latest version
1 file
+ 9
− 7
Compare changes
  • Side-by-side
  • Inline
notebooks/AGIPD/Chracterize_AGIPD_Gain_PC_Summary.ipynb
+ 9
− 7
%% 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
in_folder = "" # path to input data, required
out_folder = "" # path to output to, required
metadata_folder = "" # Directory containing calibration_metadata.yml when run by xfel-calibrate
runs = [] # 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"
db_output = False
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,
get_from_db,
get_pdu_from_db,
get_report,
send_to_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"
)
```
%% Cell type:code id: tags:
``` python
# ml - high gain slope
# bl - high gain intercept
# devl - absolute relative deviation from linearity for high gain
# mh - medium gain slope
# bh - medium gain intercept
# oh, ch, ah - parameters of hook function fit to medium gain (only if requested)
# devh - absolute relative deviation from linearity for linear part of medium gain
keys = ['ml', 'bl', 'mh', 'bh', 'BadPixelsPC']
keys_file = ["ml", "bl", "devl", "mh", "bh", "oh", "ch", "ah", "devh"]
fit_data = {}
bad_pixels = {}
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}')
print(f'Trying to find: {constants_file}')
print(f'Data available for module {qm}\n')
with h5py.File(constants_file, 'r') as hf:
bad_pixels[mod] = hf[f'/{qm}/BadPixelsPC/0/data'][()]
for key in keys_file:
fit_data[mod][key]= hf[f'/{qm}/{key}/0/data'][()]
modules.append(mod)
karabo_da.append(f"AGIPD{mod:02d}")
else:
print(f"No fit data available for module {qm}")
print(f"No fit data available for module {qm}\n")
```
%% Cell type:code id: tags:
``` python
def slope_dict_to_arr(d):
key_to_index = {
"ml": 0,
"bl": 1,
"devl": 2,
"mh": 3,
"bh": 4,
"oh": 5,
"ch": 6,
"ah": 7,
"devh": 8,
"tresh": 9
}
arr = np.zeros([11]+list(d["ml"].shape), np.float32)
for key, item in d.items():
if key not in key_to_index:
continue
arr[key_to_index[key],...] = item
return arr
```
%% Cell type:code id: tags:
``` python
# set the operating condition
condition = Conditions.Dark.AGIPD(memory_cells=mem_cells,
bias_voltage=bias_voltage,
acquisition_rate=acq_rate,
gain_setting=gain_setting,
integration_time=integration_time)
db_modules = get_pdu_from_db(karabo_id, karabo_da, Constants.AGIPD.SlopesPC(),
condition, cal_db_interface,
snapshot_at=creation_time)
```
%% Cell type:code id: tags:
``` python
proposal = list(filter(None, in_folder.strip('/').split('/')))[-2]
file_loc = proposal + ' ' + ' '.join(list(map(str,runs)))
report = get_report(metadata_folder)
```
%% Cell type:code id: tags:
``` python
if not db_output:
print('Injection to DB not requested.')
```
%% Cell type:code id: tags:
``` python
md = None
if db_output:
for mod, pdu in zip(modules, db_modules):
for const in ["SlopesPC", "BadPixelsPC"]:
dbconst = getattr(Constants.AGIPD, const)()
if const == "SlopesPC":
dbconst.data = slope_dict_to_arr(fit_data[mod])
else:
dbconst.data = bad_pixels[mod]
md = send_to_db(pdu, karabo_id, dbconst, condition,
file_loc, report, cal_db_interface,
creation_time=creation_time,
variant=1)
print("Constants parameter conditions are:\n")
print(f"• memory_cells: {mem_cells}\n• bias_voltage: {bias_voltage}\n"
f"• acquisition_rate: {acq_rate}\n• gain_setting: {gain_setting}\n"
f"• integration_time: {integration_time}\n"
f"• creation_time: {md.calibration_constant_version.begin_at if md is not None else creation_time}\n")
```
%% Cell type:code id: tags:
``` python
# Remove keys which won't be used for comparison plots and add BP to the rest of data
for mod in modules:
fit_data[mod]['BadPixelsPC'] = bad_pixels[mod]
for key in keys_file:
if key not in keys:
del fit_data[mod][key]
for key in keys:
fit_data[mod][key] = fit_data[mod][key].swapaxes(1,2)
```
%% 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)))
```
%% 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 = {}
old_const = {}
const_order = [0, 1, 3, 4]
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 modules:
for cnt, i in enumerate(modules):
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]] = old_PC_consts[i][0][c].swapaxes(1,2)
if old_PC_consts[0][0]:
old_const[key][i,:,pixel_range[0]:pixel_range[2],
pixel_range[1]:pixel_range[3]] = old_PC_consts[cnt][0][c].swapaxes(1,2)
const_data['BadPixelsPC'] = np.full((nmods, mem_cells, 512, 128), np.nan)
for i in modules:
const_data['BadPixelsPC'][i,:,pixel_range[0]:pixel_range[2],
pixel_range[1]:pixel_range[3]] = fit_data[i]['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': {}
}
old_gain_data = {'HG': {},
'MG': {}
}
for key in const_data.keys():
if key is 'ml' or key is 'bl':
gain_data['HG'][key] = const_data[key]
old_gain_data['HG'][key] = old_const[key]
if key is 'mh' or key is 'bh':
gain_data['MG'][key] = const_data[key]
old_gain_data['MG'][key] = old_const[key]
if key is 'BadPixelsPC':
gain_data['BP'][key] = const_data[key]
```
%% Cell type:code id: tags:
``` python
def plot_definition(data, g, key):
titel = ['Difference to previous', 'Percentage difference']
gs = gridspec.GridSpec(1, 2)
fig = plt.figure(figsize=(15, 7))
plt.suptitle(f'{g}', fontsize=16)
for pos in range(0,2):
vmin, vmax = get_range(data[pos], 2)
vmax = max(vmax, abs(vmin))
ax = fig.add_subplot(gs[0, pos])
geom.plot_data_fast(data[pos],
vmin=vmin, vmax=vmax, ax=ax, cmap="RdBu", figsize=(13,7),
colorbar={'shrink': 1,
'pad': 0.04,
'fraction': 0.1
})
colorbar = ax.images[0].colorbar
if pos == 1:
colorbar.set_label('%')
ax.set_title(f"{titel[pos]}: {key}", fontsize=14)
ax.set_xlabel("Columns", fontsize=12)
ax.set_ylabel("Rows", fontsize=12)
def plot_diff_consts(old_const, new_const, g, ratio=False):
if ratio:
old_data = old_const['HG']['ml'] / old_const['MG']['mh']
new_data = new_const['HG']['ml'] / new_const['MG']['mh']
data1 = np.nanmean((new_data - old_data), axis=1)
data2 = np.nanmean((new_data - old_data)/old_data*100, axis=1)
data = [data1, data2]
key ='Slopes ratio HG/MG'
plot_definition(data, g, key)
else:
for i, key in enumerate(old_const[g].keys()):
data1 = np.nanmean((new_const[g][key] - old_const[g][key]), axis=1)
data2 = np.nanmean((new_const[g][key] - old_const[g][key])/old_const[g][key]*100, axis=1)
data = [data1, data2]
plot_definition(data, g, key)
```
%% Cell type:code id: tags:
``` python
for gain in old_gain_data.keys():
plot_diff_consts(old_gain_data, gain_data, gain)
```
%% Cell type:code id: tags:
``` python
plot_diff_consts(old_gain_data, gain_data, 'Ratio HG/MG', ratio=True)
```
%% 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 = []
ratio_old = old_gain_data['HG']['ml'] / old_gain_data['MG']['mh']
for mod in modules:
table.append((mod,
f"{np.nanmean(ratio[mod]):0.1f} +- {np.nanstd(ratio[mod]):0.2f}",
f"{np.nanmean(ratio_old[mod]):0.1f} +- {np.nanstd(ratio_old[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 = []
all_HM_old = []
for mod in modules:
all_HM.extend(ratio[mod])
all_HM_old.extend(ratio_old[mod])
all_HM = np.array(all_HM)
all_HM_old = np.array(all_HM_old)
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_HM_old):0.1f} +- {np.nanstd(all_HM_old):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",
"HG/MG Ratio",
"Previous HG/MG Ratio",
"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()
```

Legal Notice | Contacts