Skip to content
Snippets Groups Projects
Commit 5f464dc1 authored by David Hammer's avatar David Hammer
Browse files

Adding prototype gain mode detection from slow data

parent 47f79f00
No related branches found
No related tags found
1 merge request!438Feat/agipd add fixed gain mode to dark notebook
Hide whitespace changes
Inline Side-by-side
notebooks/AGIPD/Characterize_AGIPD_Gain_Darks_NBC.ipynb
+ 28
− 6
View file @ 5f464dc1
%% Cell type:markdown id: tags:
# AGIPD Characterize Dark Images #
Author: S. Hauf, Version: 0.1
The following code analyzes 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. Data for the detector's three gain stages needs to be present, separated into separate runs.
The evaluated calibration constants are stored locally and injected in the calibration data base.
%% Cell type:code id: tags:
``` python
in_folder = "/gpfs/exfel/d/raw/CALLAB/202031/p900113" # path to input data, required
out_folder = "/gpfs/exfel/data/scratch/hammerd/agipd-fixed-gain" # path to output to, required
sequences = [0] # sequence files to evaluate.
modules = [-1] # list of modules to evaluate, RANGE ALLOWED
run_high = 9985 # run number in which high gain data was recorded, required
run_med = 9984 # run number in which medium gain data was recorded, required
run_low = 9983 # run number in which low gain data was recorded, required
operation_mode = 'ADAPTIVE_GAIN' # Detector operation mode, optional
karabo_id = "HED_DET_AGIPD500K2G" # karabo karabo_id
karabo_da = ['-1'] # a list of data aggregators names, Default [-1] for selecting all data aggregators
receiver_id = "{}CH0" # inset for receiver devices
path_template = 'RAW-R{:04d}-{}-S{:05d}.h5' # the template to use to access data
h5path = '/INSTRUMENT/{}/DET/{}:xtdf/image' # path in the HDF5 file to images
h5path_idx = '/INDEX/{}/DET/{}:xtdf/image' # path in the HDF5 file to images
h5path_ctrl = '/CONTROL/{}/MDL/FPGA_COMP' # path to control information
karabo_id_control = "HED_EXP_AGIPD500K2G" # karabo-id for control device '
karabo_da_control = "AGIPD500K2G00" # karabo DA for control infromation
use_dir_creation_date = True # use dir creation date as data production reference date
cal_db_interface = "tcp://max-exfl016:8020" # the database interface to use
cal_db_timeout = 3000000 # timeout on caldb requests"
local_output = True # output constants locally
db_output = False # output constants to database
mem_cells = 0 # number of memory cells used, set to 0 to automatically infer
bias_voltage = 0 # detector bias voltage
gain_setting = 0.1 # the gain setting, use 0.1 to try to auto-determine
acq_rate = 0. # the detector acquisition rate, use 0 to try to auto-determine
interlaced = False # assume interlaced data format, for data prior to Dec. 2017
rawversion = 2 # RAW file format version
thresholds_offset_sigma = 3. # offset sigma thresholds for offset deduced bad pixels
thresholds_offset_hard = [0, 0] # For setting the same threshold offset for the 3 gains. Left for backcompatability. Default [0, 0] to take the following parameters.
thresholds_offset_hard_hg = [3000, 7000] # High-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_offset_hard_mg = [6000, 10000] # Medium-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_offset_hard_lg = [6000, 10000] # Low-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_noise_sigma = 5. # noise sigma thresholds for offset deduced bad pixels
thresholds_noise_hard = [0, 0] # For setting the same threshold noise for the 3 gains. Left for backcompatability. Default [0, 0] to take the following parameters.
thresholds_noise_hard_hg = [4, 20] # High-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_noise_hard_mg = [4, 20] # Medium-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_noise_hard_lg = [4, 20] # Low-gain thresholds in absolute ADU terms for offset deduced bad pixels
thresholds_gain_sigma = 5. # Gain separation sigma threshold
high_res_badpix_3d = False # set this to True if you need high-resolution 3d bad pixel plots. ~7mins extra time for 64 memory cells
```
%% Cell type:code id: tags:
``` python
# imports and things that do not usually need to be changed
from datetime import datetime
import dateutil.parser
import os
from collections import OrderedDict
from typing import List, Tuple
import h5py
import matplotlib
import numpy as np
import tabulate
from cal_tools.enums import BadPixels
matplotlib.use('agg')
import matplotlib.pyplot as plt
from IPython.display import Latex, Markdown, display
%matplotlib inline
import multiprocessing
from cal_tools.agipdlib import (get_acq_rate, get_bias_voltage,
get_gain_setting, get_num_cells)
from cal_tools.enums import AgipdGainMode
from cal_tools.plotting import (create_constant_overview, plot_badpix_3d,
show_overview, show_processed_modules)
from cal_tools.tools import (get_dir_creation_date, get_from_db,
get_pdu_from_db, get_random_db_interface,
get_report, map_gain_stages,
run_prop_seq_from_path, save_const_to_h5,
send_to_db)
from iCalibrationDB import Conditions, Constants, Detectors
```
%% Cell type:code id: tags:
``` python
def get_gain_mode(h5fn):
h5path_gainmode = f'{h5path_ctrl.replace("/CONTROL/", "/RUN/", 1)}/gainModeIndex/value'
with h5py.File(h5fn, 'r') as fd:
if h5path_gainmode in fd:
return AgipdGainMode(fd[h5path_gainmode][0])
return AgipdGainMode.ADAPTIVE_GAIN
```
%% Cell type:code id: tags:
``` python
# insert control device if format string (does nothing otherwise)
h5path_ctrl = h5path_ctrl.format(karabo_id_control)
gains = np.arange(3)
max_cells = mem_cells
offset_runs = OrderedDict()
offset_runs["high"] = run_high
offset_runs["med"] = run_med
offset_runs["low"] = run_low
creation_time=None
if use_dir_creation_date:
creation_time = get_dir_creation_date(in_folder, run_high)
print(f"Using {creation_time} as creation time of constant.")
run, prop, seq = run_prop_seq_from_path(in_folder)
cal_db_interface = get_random_db_interface(cal_db_interface)
print(f'Calibration database interface: {cal_db_interface}')
instrument = karabo_id.split("_")[0]
if instrument == "SPB":
dinstance = "AGIPD1M1"
nmods = 16
elif instrument == "MID":
dinstance = "AGIPD1M2"
nmods = 16
elif instrument == "HED":
dinstance = "AGIPD500K"
nmods = 8
control_names = [f'{in_folder}/r{r:04d}/RAW-R{r:04d}-{karabo_da_control}-S00000.h5'
for r in (run_high, run_med, run_low)]
run_gain_modes = [get_gain_mode(fn) for fn in control_names]
if all(gm == AgipdGainMode.ADAPTIVE_GAIN for gm in run_gain_modes):
fixed_gain_mode = False
elif run_gain_modes == [AgipdGainMode.FIXED_HIGH_GAIN, AgipdGainMode.FIXED_MEDIUM_GAIN, AgipdGainMode.FIXED_LOW_GAIN]:
fixed_gain_mode = True
else:
print(f'Something is clearly wrong; slow data indicates gain modes {run_gain_modes}')
print(f"Detector in use is {karabo_id}")
print(f"Instrument {instrument}")
print(f"Detector instance {dinstance}")
print(f"Operation mode is {operation_mode}")
print(f"Operation mode is {'fixed' if fixed_gain_mode else 'adaptive'} gain mode")
```
%% Cell type:code id: tags:
``` python
runs = [run_high, run_med, run_low]
# insert control device if format string (does nothing otherwise)
h5path_ctrl = h5path_ctrl.format(karabo_id_control)
if gain_setting == 0.1:
if creation_time.replace(tzinfo=None) < dateutil.parser.parse('2020-01-31'):
print("Set gain-setting to None for runs taken before 2020-01-31")
gain_setting = None
else:
try:
# extract gain setting and validate that all runs have the same setting
gsettings = []
for r in runs:
control_fname = '{}/r{:04d}/RAW-R{:04d}-{}-S00000.h5'.format(in_folder, r, r,
karabo_da_control)
gsettings.append(get_gain_setting(control_fname, h5path_ctrl))
if not all(g == gsettings[0] for g in gsettings):
raise ValueError(f"Different gain settings for the 3 input runs {gsettings}")
gain_setting = gsettings[0]
except Exception as e:
print(f'Error while reading gain setting from: \n{control_fname}')
print(f'Error: {e}')
if "component not found" in str(e):
print("Gain setting is not found in the control information")
print("Data will not be processed")
sequences = []
```
%% Cell type:code id: tags:
``` python
if karabo_da[0] == '-1':
if modules[0] == -1:
modules = list(range(nmods))
karabo_da = ["AGIPD{:02d}".format(i) for i in modules]
else:
modules = [int(x[-2:]) for x in karabo_da]
h5path = h5path.format(karabo_id, receiver_id)
h5path_idx = h5path_idx.format(karabo_id, receiver_id)
if bias_voltage == 0:
# Read the bias voltage from files, if recorded.
# If not available, make use of the historical voltage the detector is running at
bias_voltage = get_bias_voltage(control_names[0], karabo_id_control)
bias_voltage = bias_voltage if bias_voltage is not None else 300
print("Parameters are:")
print(f"Proposal: {prop}")
print(f"Memory cells: {mem_cells}/{max_cells}")
print("Runs: {}".format([ v for v in offset_runs.values()]))
print(f"Sequences: {sequences}")
print(f"Interlaced mode: {interlaced}")
print(f"Using DB: {db_output}")
print(f"Input: {in_folder}")
print(f"Output: {out_folder}")
print(f"Bias voltage: {bias_voltage}V")
print(f"Gain setting: {gain_setting}")
```
%% Cell type:markdown id: tags:
The following lines will create a queue of files which will the be executed module-parallel. Distiguishing between different gains.
%% Cell type:code id: tags:
``` python
# set everything up filewise
os.makedirs(out_folder, exist_ok=True)
gmf = map_gain_stages(in_folder, offset_runs, path_template, karabo_da, sequences)
gain_mapped_files, total_sequences, total_file_size = gmf
print(f"Will process a total of {total_sequences} files.")
```
%% Cell type:code id: tags:
``` python
def index_to_qm(i: int):
return f"Q{i//4+1}M{i%4+1}"
```
%% Cell type:markdown id: tags:
## Calculate Offsets, Noise and Thresholds ##
The calculation is performed per-pixel and per-memory-cell. Offsets are simply the median value for a set of dark data taken at a given gain, noise the standard deviation, and gain-bit values the medians of the gain array.
%% Cell type:code id: tags:
``` python
if thresholds_offset_hard == [0, 0]:
thresholds_offset_hard = [thresholds_offset_hard_hg, thresholds_offset_hard_mg, thresholds_offset_hard_lg]
else:
thresholds_offset_hard = [thresholds_offset_hard] * 3
if thresholds_noise_hard == [0, 0]:
thresholds_noise_hard = [thresholds_noise_hard_hg, thresholds_noise_hard_mg, thresholds_noise_hard_lg]
else:
thresholds_noise_hard = [thresholds_noise_hard] * 3
```
%% Cell type:code id: tags:
``` python
def characterize_module(fast_data_filename: str, channel: int, gg: int) -> Tuple[np.array, np.array, np.array, np.array, int, np.array, int, float]:
if max_cells == 0:
num_cells = get_num_cells(fast_data_filename, karabo_id, channel)
else:
num_cells = max_cells
print(f"Using {num_cells} memory cells")
if acq_rate == 0.:
slow_paths = control_names[gg], karabo_id_control
fast_paths = fast_data_filename, karabo_id, channel
local_acq_rate = get_acq_rate(fast_paths, slow_paths)
else:
local_acq_rate = acq_rate
local_thresholds_offset_hard = thresholds_offset_hard[gg]
local_thresholds_noise_hard = thresholds_noise_hard[gg]
h5path_f = h5path.format(channel)
h5path_idx_f = h5path_idx.format(channel)
with h5py.File(fast_data_filename, "r", driver="core") as infile:
if rawversion == 2:
count = np.squeeze(infile[f"{h5path_idx_f}/count"])
first = np.squeeze(infile[f"{h5path_idx_f}/first"])
last_index = int(first[count != 0][-1]+count[count != 0][-1])
first_index = int(first[count != 0][0])
else:
status = np.squeeze(infile[f"{h5path_idx_f}/status"])
if np.count_nonzero(status != 0) == 0:
return
last = np.squeeze(infile[f"{h5path_idx_f}/last"])
first = np.squeeze(infile[f"{h5path_idx_f}/first"])
last_index = int(last[status != 0][-1]) + 1
first_index = int(first[status != 0][0])
im = np.array(infile[f"{h5path_f}/data"][first_index:last_index,...])
cellIds = np.squeeze(infile[f"{h5path_f}/cellId"][first_index:last_index,...])
if interlaced:
ga = im[1::2, 0, ...]
im = im[0::2, 0, ...].astype(np.float32)
cellIds = cellIds[::2]
else:
ga = im[:, 1, ...]
im = im[:, 0, ...].astype(np.float32)
im = np.rollaxis(im, 2)
im = np.rollaxis(im, 2, 1)
ga = np.rollaxis(ga, 2)
ga = np.rollaxis(ga, 2, 1)
offset = np.zeros((im.shape[0], im.shape[1], num_cells))
gains = np.zeros((im.shape[0], im.shape[1], num_cells))
noise = np.zeros((im.shape[0], im.shape[1], num_cells))
gains_std = np.zeros((im.shape[0], im.shape[1], num_cells))
for cc in np.unique(cellIds[cellIds < num_cells]):
cellidx = cellIds == cc
offset[...,cc] = np.median(im[..., cellidx], axis=2)
noise[...,cc] = np.std(im[..., cellidx], axis=2)
gains[...,cc] = np.median(ga[..., cellidx], axis=2)
gains_std[...,cc] = np.std(ga[..., cellidx], axis=2)
# bad pixels
bp = np.zeros(offset.shape, np.uint32)
# offset related bad pixels
offset_mn = np.nanmedian(offset, axis=(0,1))
offset_std = np.nanstd(offset, axis=(0,1))
bp[(offset < offset_mn-thresholds_offset_sigma*offset_std) |
(offset > offset_mn+thresholds_offset_sigma*offset_std)] |= BadPixels.OFFSET_OUT_OF_THRESHOLD.value
bp[(offset < local_thresholds_offset_hard[0]) | (
offset > local_thresholds_offset_hard[1])] |= BadPixels.OFFSET_OUT_OF_THRESHOLD.value
bp[~np.isfinite(offset)] |= BadPixels.OFFSET_NOISE_EVAL_ERROR.value
# noise related bad pixels
noise_mn = np.nanmedian(noise, axis=(0,1))
noise_std = np.nanstd(noise, axis=(0,1))
bp[(noise < noise_mn-thresholds_noise_sigma*noise_std) |
(noise > noise_mn+thresholds_noise_sigma*noise_std)] |= BadPixels.NOISE_OUT_OF_THRESHOLD.value
bp[(noise < local_thresholds_noise_hard[0]) | (noise > local_thresholds_noise_hard[1])] |= BadPixels.NOISE_OUT_OF_THRESHOLD.value
bp[~np.isfinite(noise)] |= BadPixels.OFFSET_NOISE_EVAL_ERROR.value
return offset, noise, gains, gains_std, gg, bp, num_cells, local_acq_rate
offset_g = OrderedDict()
noise_g = OrderedDict()
gain_g = OrderedDict()
gainstd_g = OrderedDict()
badpix_g = OrderedDict()
start = datetime.now()
all_cells = []
all_acq_rate = []
inp = []
for gg, (gain, mapped_files) in enumerate(gain_mapped_files.items()):
dones = []
for i in modules:
qm = index_to_qm(i)
if qm in mapped_files and not mapped_files[qm].empty():
fname_in = mapped_files[qm].get()
print("Process file: ", fname_in)
dones.append(mapped_files[qm].empty())
else:
continue
inp.append((fname_in, i, gg))
with multiprocessing.Pool() as pool:
results = pool.starmap(characterize_module, inp)
for offset, noise, gains, gains_std, gg, bp, thiscell, thisacq in results:
all_cells.append(thiscell)
all_acq_rate.append(thisacq)
for i in modules:
qm = index_to_qm(i)
if qm not in offset_g:
offset_g[qm] = np.zeros((offset.shape[0], offset.shape[1], offset.shape[2], 3))
noise_g[qm] = np.zeros_like(offset_g[qm])
gain_g[qm] = np.zeros_like(offset_g[qm])
gainstd_g[qm] = np.zeros_like(offset_g[qm])
badpix_g[qm] = np.zeros_like(offset_g[qm], np.uint32)
offset_g[qm][...,gg] = offset
noise_g[qm][...,gg] = noise
gain_g[qm][...,gg] = gains
gainstd_g[qm][..., gg] = gains_std
badpix_g[qm][...,gg] = bp
duration = (datetime.now() - start).total_seconds()
max_cells = np.max(all_cells)
print(f"Using {max_cells} memory cells")
acq_rate = np.max(all_acq_rate)
print(f"Using {acq_rate} MHz acquisition rate")
```
%% Cell type:code id: tags:
``` python
# Add a badpixel due to bad gain separation
for g in range(2):
# Bad pixels during bad gain separation.
# Fraction of pixels in the module with separation lower than "thresholds_gain_sigma".
bad_sep = (gain_g[qm][..., g+1] - gain_g[qm][..., g]) / np.sqrt(gainstd_g[qm][..., g+1]**2 + gainstd_g[qm][..., g]**2)
badpix_g[qm][...,g+1][(bad_sep)<thresholds_gain_sigma]|= BadPixels.GAIN_THRESHOLDING_ERROR.value
```
%% Cell type:markdown id: tags:
The thresholds for gain switching are then defined as the mean value between in individual gain bit levels. Note that these thresholds need to be refined with charge induced thresholds, as the two are not the same.
%% Cell type:code id: tags:
``` python
thresholds_g = {}
for qm in gain_g.keys():
thresholds_g[qm] = np.zeros((gain_g[qm].shape[0], gain_g[qm].shape[1], gain_g[qm].shape[2], 5))
thresholds_g[qm][...,0] = (gain_g[qm][...,1]+gain_g[qm][...,0])/2
thresholds_g[qm][...,1] = (gain_g[qm][...,2]+gain_g[qm][...,1])/2
for i in range(3):
thresholds_g[qm][...,2+i] = gain_g[qm][...,i]
```
%% Cell type:code id: tags:
``` python
res = OrderedDict()
for i in modules:
qm = index_to_qm(i)
res[qm] = {'Offset': offset_g[qm],
'Noise': noise_g[qm],
'ThresholdsDark': thresholds_g[qm],
'BadPixelsDark': badpix_g[qm]
}
```
%% Cell type:code id: tags:
``` python
# Read report path and create file location tuple to add with the injection
proposal = list(filter(None, in_folder.strip('/').split('/')))[-2]
file_loc = 'proposal:{} runs:{} {} {}'.format(proposal, run_low, run_med, run_high)
report = get_report(out_folder)
```
%% Cell type:code id: tags:
``` python
# TODO: add db_module when received from myMDC
# Create the modules dict of karabo_das and PDUs
qm_dict = OrderedDict()
for i, k_da in zip(modules, karabo_da):
qm = index_to_qm(i)
qm_dict[qm] = {"karabo_da": k_da,
"db_module": ""}
```
%% Cell type:code id: tags:
``` python
# Retrieve existing constants for comparison
clist = ["Offset", "Noise", "ThresholdsDark", "BadPixelsDark"]
old_const = {}
old_mdata = {}
detinst = getattr(Detectors, dinstance)
print('Retrieve pre-existing constants for comparison.')
for qm in res:
qm_db = qm_dict[qm]
karabo_da = qm_db["karabo_da"]
for const in res[qm]:
dconst = getattr(Constants.AGIPD, const)()
dconst.data = res[qm][const]
# Setting conditions
condition = Conditions.Dark.AGIPD(memory_cells=max_cells,
bias_voltage=bias_voltage,
acquisition_rate=acq_rate,
gain_setting=gain_setting)
# This should be used in case of running notebook
# by a different method other than myMDC which already
# sends CalCat info.
# TODO: Set db_module to "" by default in the first cell
if not qm_db["db_module"]:
qm_db["db_module"] = get_pdu_from_db(karabo_id, karabo_da, dconst,
condition, cal_db_interface,
snapshot_at=creation_time)[0]
data, mdata = get_from_db(karabo_id, karabo_da,
dconst, condition,
None, cal_db_interface,
creation_time=creation_time,
verbosity=2, timeout=cal_db_timeout)
old_const[const] = data
if mdata is not None and data is not None:
time = mdata.calibration_constant_version.begin_at
old_mdata[const] = time.isoformat()
os.makedirs('{}/old/'.format(out_folder), exist_ok=True)
save_const_to_h5(qm_db["db_module"], karabo_id,
dconst, condition, data,
file_loc, report, creation_time,
f'{out_folder}/old/')
else:
old_mdata[const] = "Not found"
```
%% Cell type:code id: tags:
``` python
md = None
for qm in res:
karabo_da = qm_dict[qm]["karabo_da"]
db_module = qm_dict[qm]["db_module"]
for const in res[qm]:
dconst = getattr(Constants.AGIPD, const)()
dconst.data = res[qm][const]
# set the operating condition
condition = Conditions.Dark.AGIPD(memory_cells=max_cells,
bias_voltage=bias_voltage,
acquisition_rate=acq_rate,
gain_setting=gain_setting)
if db_output:
md = send_to_db(db_module, karabo_id, dconst, condition, file_loc,
report, cal_db_interface, creation_time=creation_time,
timeout=cal_db_timeout)
if local_output:
md = save_const_to_h5(db_module, karabo_id, dconst, condition, dconst.data,
file_loc, report, creation_time, out_folder)
print(f"Calibration constant {const} is stored locally.\n")
print("Constants parameter conditions are:\n")
print(f"• memory_cells: {max_cells}\n• bias_voltage: {bias_voltage}\n"
f"• acquisition_rate: {acq_rate}\n• gain_setting: {gain_setting}\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
mnames=[]
for i in modules:
qm = index_to_qm(i)
mnames.append(qm)
display(Markdown(f'## Position of the module {qm} and its ASICs##'))
show_processed_modules(dinstance, constants=None, mnames=mnames, mode="position")
```
%% Cell type:markdown id: tags:
## Single-Cell Overviews ##
Single cell overviews allow to identify potential effects on all memory cells, e.g. on sensor level. Additionally, they should serve as a first sanity check on expected behaviour, e.g. if structuring on the ASIC level is visible in the offsets, but otherwise no immediate artifacts are visible.
%% Cell type:markdown id: tags:
### High Gain ###
%% Cell type:code id: tags:
``` python
cell = 3
gain = 0
show_overview(res, cell, gain, infix="{}-{}-{}".format(*offset_runs.values()))
```
%% Cell type:markdown id: tags:
### Medium Gain ###
%% Cell type:code id: tags:
``` python
cell = 3
gain = 1
show_overview(res, cell, gain, infix="{}-{}-{}".format(*offset_runs.values()))
```
%% Cell type:markdown id: tags:
### Low Gain ###
%% Cell type:code id: tags:
``` python
cell = 3
gain = 2
show_overview(res, cell, gain, infix="{}-{}-{}".format(*offset_runs.values()))
```
%% Cell type:code id: tags:
``` python
cols = {BadPixels.NOISE_OUT_OF_THRESHOLD.value: (BadPixels.NOISE_OUT_OF_THRESHOLD.name, '#FF000080'),
BadPixels.OFFSET_NOISE_EVAL_ERROR.value: (BadPixels.OFFSET_NOISE_EVAL_ERROR.name, '#0000FF80'),
BadPixels.OFFSET_OUT_OF_THRESHOLD.value: (BadPixels.OFFSET_OUT_OF_THRESHOLD.name, '#00FF0080'),
BadPixels.GAIN_THRESHOLDING_ERROR.value: (BadPixels.GAIN_THRESHOLDING_ERROR.name, '#FF40FF40'),
BadPixels.OFFSET_OUT_OF_THRESHOLD.value | BadPixels.NOISE_OUT_OF_THRESHOLD.value: ('OFFSET_OUT_OF_THRESHOLD + NOISE_OUT_OF_THRESHOLD', '#DD00DD80'),
BadPixels.OFFSET_OUT_OF_THRESHOLD.value | BadPixels.NOISE_OUT_OF_THRESHOLD.value |
BadPixels.GAIN_THRESHOLDING_ERROR.value: ('MIXED', '#BFDF009F')}
if high_res_badpix_3d:
display(Markdown("""
## Global Bad Pixel Behaviour ##
The following plots show the results of bad pixel evaluation for all evaluated memory cells.
Cells are stacked in the Z-dimension, while pixels values in x/y are rebinned with a factor of 2.
This excludes single bad pixels present only in disconnected pixels.
Hence, any bad pixels spanning at least 4 pixels in the x/y-plane, or across at least two memory cells are indicated.
Colors encode the bad pixel type, or mixed type.
"""))
gnames = ['High Gain', 'Medium Gain', 'Low Gain']
for gain in range(3):
display(Markdown(f'### {gnames[gain]} ###'))
for mod, data in badpix_g.items():
plot_badpix_3d(data[...,gain], cols, title=mod, rebin_fac=1)
plt.show()
```
%% Cell type:markdown id: tags:
## Aggregate values, and per Cell behaviour ##
The following tables and plots give an overview of statistical aggregates for each constant, as well as per cell behavior.
%% Cell type:code id: tags:
``` python
create_constant_overview(offset_g, "Offset (ADU)", max_cells, 4000, 8000,
badpixels=[badpix_g, np.nan])
```
%% Cell type:code id: tags:
``` python
create_constant_overview(noise_g, "Noise (ADU)", max_cells, 0, 100,
badpixels=[badpix_g, np.nan])
```
%% Cell type:code id: tags:
``` python
# Plot only three gain threshold maps.
bp_thresh = OrderedDict()
for mod, con in badpix_g.items():
bp_thresh[mod] = np.zeros((con.shape[0], con.shape[1], con.shape[2], 5), dtype=con.dtype)
bp_thresh[mod][...,:2] = con[...,:2]
bp_thresh[mod][...,2:] = con
create_constant_overview(thresholds_g, "Threshold (ADU)", max_cells, 4000, 10000, 5,
badpixels=[bp_thresh, np.nan],
gmap=['HG-MG Threshold', 'MG-LG Threshold', 'High gain', 'Medium gain', 'low gain'],
marker=['d','d','','','']
)
```
%% Cell type:code id: tags:
``` python
bad_pixel_aggregate_g = OrderedDict()
for m, d in badpix_g.items():
bad_pixel_aggregate_g[m] = d.astype(np.bool).astype(np.float)
create_constant_overview(bad_pixel_aggregate_g, "Bad pixel fraction", max_cells, 0, 0.10, 3)
```
%% Cell type:markdown id: tags:
## Summary tables ##
The following tables show summary information for the evaluated module. Values for currently evaluated constants are compared with values for pre-existing constants retrieved from the calibration database.
%% Cell type:code id: tags:
``` python
table = []
gain_names = ['High', 'Medium', 'Low']
bits = [BadPixels.NOISE_OUT_OF_THRESHOLD, BadPixels.OFFSET_OUT_OF_THRESHOLD, BadPixels.OFFSET_NOISE_EVAL_ERROR, BadPixels.GAIN_THRESHOLDING_ERROR]
for qm in badpix_g.keys():
for gain in range(3):
l_data = []
l_data_old = []
data = np.copy(badpix_g[qm][:,:,:,gain])
datau32 = data.astype(np.uint32)
l_data.append(len(datau32[datau32>0].flatten()))
for bit in bits:
l_data.append(np.count_nonzero(badpix_g[qm][:,:,:,gain] & bit.value))
if old_const['BadPixelsDark'] is not None:
dataold = np.copy(old_const['BadPixelsDark'][:, :, :, gain])
datau32old = dataold.astype(np.uint32)
l_data_old.append(len(datau32old[datau32old>0].flatten()))
for bit in bits:
l_data_old.append(np.count_nonzero(old_const['BadPixelsDark'][:, :, :, gain] & bit.value))
l_data_name = ['All bad pixels', 'NOISE_OUT_OF_THRESHOLD',
'OFFSET_OUT_OF_THRESHOLD', 'OFFSET_NOISE_EVAL_ERROR', 'GAIN_THRESHOLDING_ERROR']
l_threshold = ['', f'{thresholds_noise_sigma}' f'{thresholds_noise_hard[gain]}',
f'{thresholds_offset_sigma}' f'{thresholds_offset_hard[gain]}',
'', f'{thresholds_gain_sigma}']
for i in range(len(l_data)):
line = [f'{l_data_name[i]}, {gain_names[gain]} gain', l_threshold[i], l_data[i]]
if old_const['BadPixelsDark'] is not None:
line += [l_data_old[i]]
else:
line += ['-']
table.append(line)
table.append(['', '', '', ''])
display(Markdown('''
### Number of bad pixels ###
One pixel can be bad for different reasons, therefore, the sum of all types of bad pixels can be more than the number of all bad pixels.
'''))
if len(table)>0:
md = display(Latex(tabulate.tabulate(table, tablefmt='latex',
headers=["Pixel type", "Threshold",
"New constant", "Old constant"])))
```
%% Cell type:code id: tags:
``` python
header = ['Parameter',
"New constant", "Old constant ",
"New constant", "Old constant ",
"New constant", "Old constant ",
"New constant", "Old constant "]
for const in ['Offset', 'Noise', 'ThresholdsDark']:
if const != 'ThresholdsDark':
table = [['','High gain', 'High gain', 'Medium gain', 'Medium gain', 'Low gain', 'Low gain']]
else:
table = [['','HG-MG threshold', 'HG-MG threshold', 'MG-LG threshold', 'MG-LG threshold']]
for qm in res.keys():
data = np.copy(res[qm][const])
if const == 'ThresholdsDark':
data[...,0][res[qm]['BadPixelsDark'][...,0]>0] = np.nan
data[...,1][res[qm]['BadPixelsDark'][...,1]>0] = np.nan
else:
data[res[qm]['BadPixelsDark']>0] = np.nan
if old_const[const] is not None and old_const['BadPixelsDark'] is not None:
dataold = np.copy(old_const[const])
if const == 'ThresholdsDark':
dataold[...,0][old_const['BadPixelsDark'][...,0]>0] = np.nan
dataold[...,1][old_const['BadPixelsDark'][...,1]>0] = np.nan
else:
dataold[old_const['BadPixelsDark']>0] = np.nan
f_list = [np.nanmedian, np.nanmean, np.nanstd, np.nanmin, np.nanmax]
n_list = ['Median', 'Mean', 'Std', 'Min', 'Max']
for i, f in enumerate(f_list):
line = [n_list[i]]
for gain in range(3):
# Compare only 3 threshold gain-maps
if gain == 2 and const == 'ThresholdsDark':
continue
line.append('{:6.1f}'.format(f(data[...,gain])))
if old_const[const] is not None and old_const['BadPixelsDark'] is not None:
line.append('{:6.1f}'.format(f(dataold[...,gain])))
else:
line.append('-')
table.append(line)
display(Markdown('### {} [ADU], good pixels only ###'.format(const)))
md = display(Latex(tabulate.tabulate(table, tablefmt='latex', headers=header)))
```
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment

Legal Notice | Contacts