Skip to content
Snippets Groups Projects
Commit 701662e0 authored by Cyril Danilevski's avatar Cyril Danilevski :scooter:
Browse files

[AGIPD][Summary] Add Summary Figures for mini-half

Closes calibration_workshop#219

See merge request detectors/pycalibration!393
parents a10c44e8 8b35a368
No related branches found
No related tags found
1 merge request!393[AGIPD][Summary] Add Summary Figures for mini-half
%% Cell type:code id: tags:
``` python
#Author: K. Ahmed, M. Karnevsky, Version: 0.1
#The following is a summary for the processing of dark images and calibration constants production.
out_folder = "/gpfs/exfel/data/scratch/ahmedk/test/miniHalfAGIPD" # path to output to, required
karabo_id = "DETLAB_DET_AGIPD500K2G" # detector instance
gain_names = ['High gain', 'Medium gain', 'Low gain'] # a list of gain names to be used in plotting
threshold_names = ['HG-MG threshold', 'MG_LG threshold'] # a list of gain names to be used in plotting
```
%% Cell type:code id: tags:
``` python
from collections import OrderedDict
import copy
from datetime import datetime
import os
import warnings
warnings.filterwarnings('ignore')
import glob
import h5py
from IPython.display import display, Markdown, Latex
import numpy as np
import matplotlib
matplotlib.use("agg")
import matplotlib.gridspec as gridspec
import matplotlib.patches as patches
import matplotlib.pyplot as plt
%matplotlib inline
import tabulate
from cal_tools.ana_tools import get_range
from cal_tools.plotting import show_processed_modules
import extra_geom
from iCalibrationDB import Detectors
from XFELDetAna.plotting.heatmap import heatmapPlot
from XFELDetAna.plotting.simpleplot import simplePlot
```
%% Cell type:code id: tags:
``` python
if "AGIPD" in karabo_id:
if "SPB" in karabo_id:
dinstance = "AGIPD1M1"
nmods = 16
elif "MID" in karabo_id:
dinstance = "AGIPD1M2"
nmods = 16
elif "HED" in karabo_id or "DETLAB" in karabo_id:
elif "HED" in karabo_id:
dinstance = "AGIPD500K"
nmods = 8
display(Markdown("""
# Summary of AGIPD dark characterization #
The following report shows a set of dark images taken with the AGIPD detector to deduce detector offsets, noise, bad-pixel maps and thresholding. All four types of constants are evaluated per-pixel and per-memory cell.
**The offset** ($O$) is defined as the median ($M$) of the dark signal ($Ds$) over trains ($t$) for a given pixel ($x,y$) and memory cell ($c$).
**The noise** $N$ is the standard deviation $\sigma$ of the dark signal.
$$ O_{x,y,c} = M(Ds)_{t} ,\,\,\,\,\,\, N_{x,y,c} = \sigma(Ds)_{t}$$
**The bad pixel** mask is encoded as a bit mask.
**"OFFSET_OUT_OF_THRESHOLD":**
Offset outside of bounds:
$$M(O)_{x,y} - \sigma(O)_{x,y} * \mathrm{thresholds\_offset\_sigma} < O < M(O)_{x,y} + \sigma(O)_{x,y} * \mathrm{thresholds\_offset\_sigma} $$
or offset outside of hard limits
$$ \mathrm{thresholds\_offset\_hard}_\mathrm{low} < O < \mathrm{thresholds\_offset\_hard}_\mathrm{high} $$
**"NOISE_OUT_OF_THRESHOLD":**
Noise outside of bounds:
$$M(N)_{x,y} - \sigma(N)_{x,y} * \mathrm{thresholds\_noise\_sigma} < N < M(N)_{x,y} + \sigma(N)_{x,y} * \mathrm{thresholds\_noise\_sigma} $$
or noise outside of hard limits
$$\mathrm{thresholds\_noise\_hard}_\mathrm{low} < N < \mathrm{thresholds\_noise\_hard}_\mathrm{high} $$
**"OFFSET_NOISE_EVAL_ERROR":**
Offset and Noise both not $nan$ values
Values: $\mathrm{thresholds\_offset\_sigma}$, $\mathrm{thresholds\_offset\_hard}$, $\mathrm{thresholds\_noise\_sigma}$, $\mathrm{thresholds\_noise\_hard}$ are given as parameters.
"**\"GAIN_THRESHOLDING_ERROR\":**
Bad gain separated pixels with sigma separation less than gain_separation_sigma_threshold
$$ sigma\_separation = \\frac{\mathrm{gain\_offset} - \mathrm{previous\_gain\_offset}}{\sqrt{\mathrm{gain\_offset_{std}}^\mathrm{2} + \mathrm{previuos\_gain\_offset_{std}}^\mathrm{2}}}$$
$$ Bad\_separation = sigma\_separation < \mathrm{gain\_separation\_sigma\_threshold} $$
"""))
elif "LPD" in karabo_id:
dinstance = "LPD1M1"
nmods = 16
display(Markdown("""
# Summary of LPD dark characterization #
The following report shows a set of dark images taken with the LPD detector to deduce detector offsets, noise, bad-pixel maps. All three types of constants are evaluated per-pixel and per-memory cell.
**The offset** ($O$) is defined as the median ($M$) of the dark signal ($Ds$) over trains ($t$) for a given pixel ($x,y$) and memory cell ($c$).
**The noise** $N$ is the standard deviation $\sigma$ of the dark signal.
$$ O_{x,y,c} = M(Ds)_{t} ,\,\,\,\,\,\, N_{x,y,c} = \sigma(Ds)_{t}$$
**The bad pixel** mask is encoded as a bit mask.
**"OFFSET_OUT_OF_THRESHOLD":**
Offset outside of bounds:
$$M(O)_{x,y} - \sigma(O)_{x,y} * \mathrm{thresholds\_offset\_sigma} < O < M(O)_{x,y} + \sigma(O)_{x,y} * \mathrm{thresholds\_offset\_sigma} $$
or offset outside of hard limits
$$ \mathrm{thresholds\_offset\_hard}_\mathrm{low} < O < \mathrm{thresholds\_offset\_hard}_\mathrm{high} $$
**"NOISE_OUT_OF_THRESHOLD":**
Noise outside of bounds:
$$M(N)_{x,y} - \sigma(N)_{x,y} * \mathrm{thresholds\_noise\_sigma} < N < M(N)_{x,y} + \sigma(N)_{x,y} * \mathrm{thresholds\_noise\_sigma} $$
or noise outside of hard limits
$$\mathrm{thresholds\_noise\_hard}_\mathrm{low} < N < \mathrm{thresholds\_noise\_hard}_\mathrm{high} $$
**"OFFSET_NOISE_EVAL_ERROR":**
Offset and Noise both not $nan$ values
"Values: $\\mathrm{thresholds\\_offset\\_sigma}$, $\\mathrm{thresholds\\_offset\\_hard}$, $\\mathrm{thresholds\\_noise\\_sigma}$, $\\mathrm{thresholds\\_noise\\_hard}$ are given as parameters.\n",
"""))
elif "DSSC" in karabo_id:
dinstance = "DSSC1M1"
nmods = 16
display(Markdown("""
# Summary of DSSC dark characterization #
"""))
```
%% Cell type:markdown id: tags:
Preparing newly injected and previous constants from produced local folder in out_folder.
%% Cell type:code id: tags:
``` python
# TODO: After changes in the Cal DB interface read files from cal repository
# Load constants from local files
data = OrderedDict()
old_cons = OrderedDict()
mod_names = []
# Loop over modules
for i in range(nmods):
qm = f"Q{i//4+1}M{i%4+1}"
# loop over constants
detinst = getattr(Detectors, dinstance)
for const in ['Offset', 'Noise', 'ThresholdsDark', 'BadPixelsDark']:
det = getattr(detinst, qm)
if det is None:
continue
det_name = det.device_name
fpath = '{}/const_{}_{}.h5'.format(out_folder, const, det_name)
oldfpath = '{}/old/const_{}_{}.h5'.format(out_folder, const, det_name)
if not os.path.isfile(fpath):
continue
with h5py.File(fpath, 'r') as f:
if qm not in data:
mod_names.append(qm)
data[qm] = OrderedDict()
data[qm][const] = f["data"][()]
if not os.path.isfile(oldfpath):
continue
with h5py.File(oldfpath, 'r') as oldf:
if qm not in old_cons:
old_cons[qm] = OrderedDict()
old_cons[qm][const] = oldf["data"][()]
```
%% Cell type:code id: tags:
``` python
cons_shape = {}
# extracting constant shape.
for qm, constant in data.items():
for cname, cons in constant.items():
shape = data[qm][cname].shape
if cname not in cons_shape:
cons_shape[cname] = shape
break
constants = {}
prev_const = {}
for cname, sh in cons_shape.items():
constants[cname]= np.zeros((len(mod_names),) + sh[:])
prev_const[cname]= np.zeros((len(mod_names),) + sh[:])
for i in range(len(mod_names)):
for cname, cval in constants.items():
cval[i] = data[mod_names[i]][cname]
if mod_names[i] in old_cons.keys():
prev_const[cname][i] = old_cons[mod_names[i]][cname]
else:
print(f"No previous {cname} found for {mod_names[i]}")
```
%% Cell type:code id: tags:
``` python
display(Markdown('## Processed modules ##'))
show_processed_modules(dinstance, constants, mod_names, mode="processed")
```
%% Cell type:markdown id: tags:
## Summary figures across Modules ##
The following plots give an overview of calibration constants averaged across pixels and memory cells. A bad pixel mask is applied.
%% Cell type:code id: tags:
``` python
if "LPD" in dinstance:
geom = extra_geom.LPD_1MGeometry.from_quad_positions(quad_pos=[(11.4, 299),
(-11.5, 8),
(254.5, -16),
(278.5, 275)])
pixels_y = 256
pixels_x = 256
elif "AGIPD" in dinstance:
elif dinstance in ('AGIPD1M1', 'AGIPD1M2'):
geom = extra_geom.AGIPD_1MGeometry.from_quad_positions(quad_pos=[(-525, 625),
(-550, -10),
(520, -160),
(542.5, 475)])
pixels_y = 128
pixels_x = 512
elif dinstance == "AGIPD500K":
geom = extra_geom.AGIPD_500K2GGeometry.from_origin()
pixels_y = 128
pixels_x = 512
elif "DSSC" in dinstance:
pixels_y = 512
pixels_x = 128
quadpos = [(-130, 5), (-130, -125), (5, -125), (5, 5)]
extrageom_pth = os.path.dirname(extra_geom.__file__)
geom = extra_geom.DSSC_1MGeometry.from_h5_file_and_quad_positions(
f"{extrageom_pth}/tests/dssc_geo_june19.h5", positions=quadpos)
Mod_data=OrderedDict()
gainstages = 1
for const_name, const in constants.items():
#TODO: add a Summary figure across modules
if dinstance == "AGIPD500K":
display(Markdown(f'WARNING: No summary figures are available for AGIPD mini-half at the moment.'))
break
if const_name == 'BadPixelsDark':
continue
# Check if constant gain available in constant e.g. AGIPD, LPD
if len(const.shape) == 5:
gainstages = 3
else:
gainstages = 1
for dname in ['{}', 'd-{}', 'dpct-{}']:
Mod_data[dname.format(const_name)] = OrderedDict()
display(Markdown(f'##### {const_name}'))
print_once = True
for gain in range(gainstages):
if const_name == 'ThresholdsDark':
if gain > 1:
continue
glabel = threshold_names
else:
glabel = gain_names
for i in range(nmods):
qm = f"Q{i//4+1}M{i%4+1}"
if qm in mod_names:
m_idx = mod_names.index(qm)
# Check if constant shape of 5 indices e.g. AGIPD, LPD
if len(const.shape) == 5:
values = np.nanmean(const[m_idx, :, :, :, gain], axis=2)
dval = np.nanmean(prev_const[const_name][m_idx, :, :, :, gain], axis=2)
else:
values = np.nanmean(const[m_idx, :, :, :], axis=2)
dval = np.nanmean(prev_const[const_name][m_idx, :, :, :], axis=2)
values[values == 0] = np.nan
dval[dval == 0] = np.nan
dval = values - dval
dval_pct = dval/values * 100
values = np.moveaxis(values, 0, -1).reshape(1, values.shape[1], values.shape[0])
dval = np.moveaxis(dval, 0, -1).reshape(1, dval.shape[1], dval.shape[0])
dval_pct = np.moveaxis(dval_pct, 0, -1).reshape(1, dval_pct.shape[1], dval_pct.shape[0])
else:
# if module not available fill arrays with nan
values = np.zeros((1, pixels_x, pixels_y),dtype=np.float64)
values[values == 0] = np.nan
dval = values
dval_pct = dval
for k, v in {'{}': values, 'd-{}': dval , 'dpct-{}': dval_pct}.items():
try:
Mod_data[k.format(const_name)][gain_names[gain]] = \
np.concatenate((Mod_data[k.format(const_name)][gain_names[gain]],
v), axis=0)
except:
Mod_data[k.format(const_name)][gain_names[gain]] = v
if np.count_nonzero(dval) == 0 and print_once:
display(Markdown(f'New and previous {const_name} are the same, hence there is no difference.'))
print_once = False
# Plotting constant overall modules.
display(Markdown(f'###### {glabel[gain]} ######'))
gs = gridspec.GridSpec(2, 2)
fig = plt.figure(figsize=(24, 32))
axis = OrderedDict()
axis = {"ax0": {"cname": "{}" ,"gs": gs[0, :], "shrink": 0.9, "pad": 0.05, "label": "ADUs", "title": '{}'},
"ax1": {"cname": "d-{}","gs": gs[1, 0], "shrink": 0.6, "pad": 0.1, "label": "ADUs", "title": 'Difference with previous {}'},
"ax2": {"cname": "dpct-{}", "gs": gs[1, 1], "shrink": 0.6, "pad": 0.1, "label": "%", "title": 'Difference with previous {} %'}}
for ax, axv in axis.items():
# Add the min and max plot values for each axis.
vmin, vmax = get_range(Mod_data[axv["cname"].format(const_name)][gain_names[gain]], 2)
ax = fig.add_subplot(axv["gs"])
geom.plot_data_fast(Mod_data[axv["cname"].format(const_name)][gain_names[gain]],
vmin=vmin, vmax=vmax, ax=ax,
colorbar={'shrink': axv["shrink"],
'pad': axv["pad"]
}
)
colorbar = ax.images[0].colorbar
colorbar.set_label(axv["label"])
ax.set_title(axv["title"].format(f"{const_name} {glabel[gain]}"), fontsize=15)
ax.set_xlabel('Columns', fontsize=15)
ax.set_ylabel('Rows', fontsize=15)
plt.show()
```
%% Cell type:code id: tags:
``` python
# Loop over modules and constants
for const_name, const in constants.items():
display(Markdown(f'### Summary across Modules - {const_name}'))
for gain in range(gainstages):
if const_name == 'ThresholdsDark':
if gain == 2:
continue
glabel = threshold_names
else:
glabel = gain_names
if len(const.shape) == 5:
data = np.copy(const[:, :, :, :, gain])
else:
data = np.copy(const[:, :, :, :])
if const_name != 'BadPixelsDark':
if "BadPixelsDark" in constants.keys():
label = f'{const_name} value [ADU], good pixels only'
if len(const.shape) == 5:
data[constants['BadPixelsDark'][:, :, :, :, gain] > 0] = np.nan
else:
data[constants['BadPixelsDark'][:, :, :, :] > 0] = np.nan
else:
label = f'{const_name} value [ADU], good and bad pixels'
datamean = np.nanmean(data, axis=(1, 2))
fig = plt.figure(figsize=(15, 6), tight_layout={
'pad': 0.2, 'w_pad': 1.3, 'h_pad': 1.3})
ax = fig.add_subplot(121)
else:
label = 'Fraction of bad pixels'
data[data > 0] = 1.0
datamean = np.nanmean(data, axis=(1, 2))
datamean[datamean == 1.0] = np.nan
fig = plt.figure(figsize=(15, 6), tight_layout={
'pad': 0.2, 'w_pad': 1.3, 'h_pad': 1.3})
ax = fig.add_subplot(111)
d = []
for im, mod in enumerate(datamean):
d.append({'x': np.arange(mod.shape[0]),
'y': mod,
'drawstyle': 'steps-pre',
'label': mod_names[im],
})
_ = simplePlot(d, figsize=(10, 10), xrange=(-12, 510),
x_label='Memory Cell ID',
y_label=label,
use_axis=ax,
title='{}'.format(glabel[gain]),
title_position=[0.5, 1.18],
legend='outside-top-ncol6-frame', legend_size='18%',
legend_pad=0.00)
if const_name != 'BadPixelsDark':
ax = fig.add_subplot(122)
if "BadPixelsDark" in constants.keys():
label = f'$\sigma$ {const_name} [ADU], good pixels only'
else:
label = f'$\sigma$ {const_name} [ADU], good and bad pixels'
d = []
for im, mod in enumerate(np.nanstd(data, axis=(1, 2))):
d.append({'x': np.arange(mod.shape[0]),
'y': mod,
'drawstyle': 'steps-pre',
'label': mod_names[im],
})
_ = simplePlot(d, figsize=(10, 10), xrange=(-12, 510),
x_label='Memory Cell ID',
y_label=label,
use_axis=ax,
title='{} $\sigma$'.format(glabel[gain]),
title_position=[0.5, 1.18],
legend='outside-top-ncol6-frame', legend_size='18%',
legend_pad=0.00)
plt.show()
```
%% Cell type:markdown id: tags:
## Summary tables across Modules ##
Tables show values averaged across all pixels and memory cells of a given detector module.
%% Cell type:code id: tags:
``` python
if u'$' in tabulate.LATEX_ESCAPE_RULES:
del(tabulate.LATEX_ESCAPE_RULES[u'$'])
if u'\\' in tabulate.LATEX_ESCAPE_RULES:
del(tabulate.LATEX_ESCAPE_RULES[u'\\'])
```
%% Cell type:code id: tags:
``` python
head = ['Module', 'High gain', 'Medium gain', 'Low gain']
head_th = ['Module', 'HG_MG threshold', 'MG_LG threshold']
for const_name, const in constants.items():
table = []
for i_mod, mod in enumerate(mod_names):
t_line = [mod]
for gain in range(gainstages):
if const_name == 'ThresholdsDark':
if gain == 2:
continue
header = head_th
else:
header = head
if len(const.shape) == 5:
data = np.copy(const[i_mod, :, :, :, gain])
else:
data = np.copy(const[i_mod, :, :, :])
if const_name == 'BadPixelsDark':
data[data > 0] = 1.0
datasum = np.nansum(data)
datamean = np.nanmean(data)
if datamean == 1.0:
datamean = np.nan
datasum = np.nan
t_line.append(f'{datasum:6.0f} ({datamean:6.3f}) ')
label = '## Number(fraction) of bad pixels'
else:
if "BadPixelsDark" in constants.keys():
data[constants['BadPixelsDark']
[i_mod, :, :, :, gain] > 0] = np.nan
label = f'### Average {const_name} [ADU], good pixels only'
else:
label = f'### Average {const_name} [ADU], good and bad pixels'
t_line.append(f'{np.nanmean(data):6.1f} $\\pm$ {np.nanstd(data):6.1f}')
label = f'## Average {const_name} [ADU], good pixels only'
table.append(t_line)
display(Markdown(label))
md = display(Latex(tabulate.tabulate(
table, tablefmt='latex', headers=header)))
```
%% Cell type:code id: tags:
``` python
```
%% Cell type:code id: tags:
``` python
```
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment