Skip to content
Snippets Groups Projects
Commit 92711a70 authored by Steffen Hauf's avatar Steffen Hauf
Browse files

Propagate backlog of cal_tools related changes from production system as of 06/2019

parent 6958d2b1
No related branches found
No related tags found
1 merge request!79Propagate backlog of cal_tools related changes from production system as of 06/2019
This diff is collapsed.
......@@ -5,24 +5,22 @@ class BadPixels(Enum):
""" The European XFEL Bad Pixel Encoding
"""
OFFSET_OUT_OF_THRESHOLD = 0b00000000000000000001 # bit 1
NOISE_OUT_OF_THRESHOLD = 0b00000000000000000010 # bit 2
OFFSET_NOISE_EVAL_ERROR = 0b00000000000000000100 # bit 3
NO_DARK_DATA = 0b00000000000000001000 # bit 4
CI_GAIN_OF_OF_THRESHOLD = 0b00000000000000010000 # bit 5
CI_LINEAR_DEVIATION = 0b00000000000000100000 # bit 6
CI_EVAL_ERROR = 0b00000000000001000000 # bit 7
FF_GAIN_EVAL_ERROR = 0b00000000000010000000 # bit 8
FF_GAIN_DEVIATION = 0b00000000000100000000 # bit 9
FF_NO_ENTRIES = 0b00000000001000000000 # bit 10
CI2_EVAL_ERROR = 0b00000000010000000000 # bit 11
VALUE_IS_NAN = 0b00000000100000000000 # bit 12
VALUE_OUT_OF_RANGE = 0b00000001000000000000 # bit 13
GAIN_THRESHOLDING_ERROR = 0b00000010000000000000 # bit 14
DATA_STD_IS_ZERO = 0b00000100000000000000 # bit 15
ASIC_STD_BELOW_NOISE = 0b00001000000000000000 # bit 16
INTERPOLATED = 0b00010000000000000000 # bit 17
NOISY_ADC = 0b00100000000000000000 # bit 18
OVERSCAN = 0b01000000000000000000 # bit 19
NON_SENSITIVE = 0b10000000000000000000 # bit 20
OFFSET_OUT_OF_THRESHOLD = 0b000000000000000001
NOISE_OUT_OF_THRESHOLD = 0b000000000000000010
OFFSET_NOISE_EVAL_ERROR = 0b000000000000000100
NO_DARK_DATA = 0b000000000000001000
CI_GAIN_OF_OF_THRESHOLD = 0b000000000000010000
CI_LINEAR_DEVIATION = 0b000000000000100000
CI_EVAL_ERROR = 0b000000000001000000
FF_GAIN_EVAL_ERROR = 0b000000000010000000
FF_GAIN_DEVIATION = 0b000000000100000000
FF_NO_ENTRIES = 0b000000001000000000
CI2_EVAL_ERROR = 0b000000010000000000
VALUE_IS_NAN = 0b000000100000000000
VALUE_OUT_OF_RANGE = 0b000001000000000000
GAIN_THRESHOLDING_ERROR = 0b000010000000000000
DATA_STD_IS_ZERO = 0b000100000000000000
ASIC_STD_BELOW_NOISE = 0b001000000000000000
INTERPOLATED = 0b010000000000000000
NOISY_ADC = 0b100000000000000000
......@@ -2,8 +2,8 @@ import copy
import h5py
import numpy as np
from cal_tools.tools import get_constant_from_db
from cal_tools.enums import BadPixels
from cal_tools.tools import get_constant_from_db, get_constant_from_db_and_time
from iCalibrationDB import Constants, Conditions, Detectors
......@@ -41,7 +41,7 @@ class LpdCorrections:
raw_fmt_version=2, chunk_size=512,
h5_data_path="INSTRUMENT/FXE_DET_LPD1M-1/DET/{}CH0:xtdf/",
h5_index_path="INDEX/FXE_DET_LPD1M-1/DET/{}CH0:xtdf/",
do_ff=True, correct_non_linear=False, karabo_data_mode=False):
do_ff=True, correct_non_linear=True, karabo_data_mode=False):
"""
Initialize an LpdCorrections Class
......@@ -85,24 +85,19 @@ class LpdCorrections:
self.bins_signal_high_range = bins_signal_high_range
self.cidx = 0
self.do_ff = do_ff
filter_modules = [11]
filter_modules = []
self.filter_cells = [0, 1] if channel in filter_modules else []
self.cnl = correct_non_linear
self.cnl = True # correct_non_linear
self.karabo_data_mode = karabo_data_mode
# emprically determined from APD datasets p900038, r155,r156
self.cnl_const = {'med':
{'lin': (0.10960777594700886,
859.4778021212404),
'repo': (0.11046365063379418,
751.5109159229758,
1436280.868413923,
9966.92701736398)},
'high': {'sigm': (1333.6836067255365,
8370.936473950766,
0.003912416547243528,
0.00013306099710934553),
'lin': (0.16284255864598732,
-29.830797125784763)}}
# emprically determined from APD datasets p900038, r155,r156
self.cnl_const = {
'high': {'A': -0.000815934, 'lam': 0.00811281, 'c': 1908.89,
'm': 0, 'b': 0}, # noqa
'med': {'A': 0.0999894, 'lam': -0.00137652, 'c': 3107.83,
'm': 3.89982e-06, 'b': -0.116811}, # noqa
'low': {'A': 0.0119132, 'lam': -0.0002, 'c': 36738.6,
'm': 2.00273e-08, 'b': 0.245537}} # noqa
def get_iteration_range(self):
"""Returns a range expression over which to iterate in chunks
......@@ -157,7 +152,7 @@ class LpdCorrections:
if not hasattr(self, "mask"):
self.mask = mask
else:
self.mask |= mask
self.mask |= mask[:self.mask.shape[0], ...]
if noise is not None:
self.noise = noise
if flatfield is not None:
......@@ -222,7 +217,6 @@ class LpdCorrections:
cells = np.squeeze(irange['image.cellId'])
length = np.squeeze(irange['image.length'])
# split gain and image info into separate arrays
im, gain = self.split_gain(im[:, 0, ...])
......@@ -269,35 +263,30 @@ class LpdCorrections:
im /= self.flatfield[None, :, :]
# hacky way of smoothening transition region between med and low
#im[gain == 2] -= im[gain == 2] * cfac[gain == 2]
# im[gain == 2] -= im[gain == 2] * cfac[gain == 2]
# perform non-linear corrections if requested
if self.cnl:
largs = self.cnl_const['high']['lin']
sargs = self.cnl_const['high']['sigm']
n, c, b, k = sargs
y = im[(gain==0) & (im > 1000)]
x = (-np.arctanh(b-y/n+1)+c*k)/k
dy = largs[0]*x+largs[1] - n*(1+b+np.tanh((x-c)*k))
im[(gain==0) & (im > 1000)] += dy
mhg = largs[0]
bhg = largs[1]
largs = self.cnl_const['med']['lin']
rargs = self.cnl_const['med']['repo']
c, k, j, n = rargs
j *= -1
mf = mhg/largs[0]
y = im[(gain==1)]
x = (np.sqrt(c**2*n**2+4*c*j+2*c*n*(k-y)+(k-y)**2)+c*n-k+y)/(2*c) - (bhg-largs[1])/largs[0]
dy = (mf*largs[0]-c)*x + largs[1] - k + j/(x-n)
im[(gain==1)] += (bhg-largs[1]) + dy
#im[(gain==1) & (im >= 3500)] *= mf
#im[(gain==1) & (im >= 3500)] += (bhg-largs[1])
def lin_exp_fun(x, m, b, A, lam, c):
return m * x + b + A * np.exp(lam * (x - c))
x = im[(gain == 0)]
cnl = self.cnl_const['high']
cf = lin_exp_fun(x, cnl['m'], cnl['b'], cnl['A'], cnl['lam'],
cnl['c'])
im[(gain == 0)] -= np.maximum(cf, -0.2) * x
x = im[(gain == 1)]
cnl = self.cnl_const['med']
cf = lin_exp_fun(x, cnl['m'], cnl['b'], cnl['A'], cnl['lam'],
cnl['c'])
im[(gain == 1)] -= np.minimum(cf, 0.2) * x
x = im[(gain == 2)]
cnl = self.cnl_const['low']
cf = lin_exp_fun(x, cnl['m'], cnl['b'], cnl['A'], cnl['lam'],
cnl['c'])
im[(gain == 2)] -= np.minimum(cf, 0.45) * x
# create bad pixels masks, here non-finite values
bidx = ~np.isfinite(im)
......@@ -378,11 +367,18 @@ class LpdCorrections:
status = np.squeeze(self.infile[lpd_base + "image/status"])
if np.count_nonzero(status != 0) == 0:
raise IOError("File {} has no valid counts".format(
self.infile))
self.infile))
last = np.squeeze(self.infile[lpd_base + "image/last"])
valid = status != 0
last_index = int(last[status != 0][-1])
first_index = int(last[status != 0][0])
idxtrains = np.squeeze(self.infile["/INDEX/trainId"])
medianTrain = np.nanmedian(idxtrains)
lowok = (idxtrains > medianTrain - 1e4)
highok = (idxtrains < medianTrain + 1e4)
valid &= lowok & highok
last_index = int(last[valid][-1])
first_index = int(last[valid][0])
else:
raise AttributeError(
"Not a known raw format version: {}".format(self.index_v))
......@@ -457,8 +453,32 @@ class LpdCorrections:
self.infile.visititems(visitor)
# sanitize indices
for do in ["image", ]:
uq, fidx, cnts = np.unique(alltrains[firange], return_index=True,
return_counts=True)
uq, fidxv, cntsv = np.unique(alltrains[firange - firange[0]],
return_index=True,
return_counts=True)
duq = (uq[1:] - uq[:-1]).astype(np.int64)
cfidxv = [fidxv[0], ]
ccntsv = [cntsv[0], ]
for i, du in enumerate(duq.tolist()):
if du > 1000:
du = 1
cntsv[i] = 0
cfidxv += [0] * (du - 1) + [fidxv[i + 1], ]
ccntsv += [0] * (du - 1) + [cntsv[i + 1], ]
mv = len(cfidxv)
fidx = np.zeros(len(cfidxv), fidxv.dtype)
fidx[self.valid[:mv]] = np.array(cfidxv)[self.valid[:mv]]
for i in range(len(fidx) - 1, 2, -1):
if fidx[i - 1] == 0 and fidx[i] != 0:
fidx[i - 1] = fidx[i]
cnts = np.zeros(len(cfidxv), cntsv.dtype)
cnts[self.valid[:mv]] = np.array(ccntsv)[self.valid[:mv]]
self.outfile.create_dataset(idx_base + "{}/first".format(do),
fidx.shape,
dtype=fidx.dtype,
......@@ -585,23 +605,24 @@ class LpdCorrections:
(cal_db_interface, creation_time, max_cells_db, capacitor,
bias_voltage, photon_energy, timeout) = dbparms
if "#" in cal_db_interface:
prot, serv, ran = cal_db_interface.split(":")
r1, r2 = ran.split("#")
cal_db_interface = ":".join(
[prot, serv, str(np.random.randint(int(r1), int(r2)))])
print("Connecting to interface: {}".format(cal_db_interface))
offset = get_constant_from_db(getattr(Detectors.LPD1M1, qm),
Constants.LPD.Offset(),
Conditions.Dark.LPD(
memory_cells=max_cells_db,
bias_voltage=bias_voltage,
capacitor=capacitor),
np.zeros((256, 256, max_cells_db, 3)),
cal_db_interface,
creation_time=creation_time,
timeout=timeout)
# if "#" in cal_db_interface:
# prot, serv, ran = cal_db_interface.split(":")
# r1, r2 = ran.split("#")
# cal_db_interface = ":".join(
# [prot, serv, str(np.random.randint(int(r1), int(r2)))])
# print("Connecting to interface: {}".format(cal_db_interface))
offset, when = get_constant_from_db_and_time(
getattr(Detectors.LPD1M1, qm),
Constants.LPD.Offset(),
Conditions.Dark.LPD(
memory_cells=max_cells_db,
bias_voltage=bias_voltage,
capacitor=capacitor),
np.zeros((256, 256, max_cells_db, 3)),
cal_db_interface,
creation_time=creation_time,
timeout=timeout)
noise = get_constant_from_db(getattr(Detectors.LPD1M1, qm),
Constants.LPD.Noise(),
......@@ -628,7 +649,7 @@ class LpdCorrections:
# done loading constants and returning
if only_dark:
self.initialize(offset, None, None, bpixels, noise, None)
return
return when
slopesCI = get_constant_from_db(getattr(Detectors.LPD1M1, qm),
Constants.LPD.SlopesCI(),
......@@ -653,7 +674,8 @@ class LpdCorrections:
pixels_x=256,
pixels_y=256,
beam_energy=None,
capacitor=capacitor), # noqa
capacitor=capacitor),
# noqa
np.ones((256, 256)),
cal_db_interface,
creation_time=creation_time,
......@@ -685,6 +707,7 @@ class LpdCorrections:
bpixels |= bpix[..., None]
self.initialize(offset, rel_gains, rel_gain_offset, bpixels, noise,
flat_fields)
return when
def initialize_from_file(self, filename, qm, with_dark=True):
""" Initialize calibration constants from a calibration file
......
This diff is collapsed.
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