diff --git a/cal_tools/cal_tools/agipdlib.py b/cal_tools/cal_tools/agipdlib.py index be8033e5a99e7ada109cbce79fc99a37f9c7ff0b..a098b07b83a7af795b8000c0aced7f4eca62c6bb 100644 --- a/cal_tools/cal_tools/agipdlib.py +++ b/cal_tools/cal_tools/agipdlib.py @@ -24,8 +24,7 @@ from cal_tools.cython import agipdalgs as calgs def get_num_cells(fname, loc, module): with h5py.File(fname, "r") as f: - cells = \ - f[f"INSTRUMENT/{loc}/DET/{module}CH0:xtdf/image/cellId"][()] + cells = f[f"INSTRUMENT/{loc}/DET/{module}CH0:xtdf/image/cellId"][()] if cells.shape[0] == 0: return None maxcell = np.max(cells) @@ -77,7 +76,7 @@ def get_acq_rate(fast_paths: Tuple[str, str, int], fast_data_file = Path(fast_data_file) if fast_data_file.is_file(): fast_data_path = f'INSTRUMENT/{karabo_id}/DET/{module}CH0:xtdf/image/pulseId' # noqa - with h5py.File(fast_data_file) as fin: + with h5py.File(fast_data_file, "r") as fin: if fast_data_path in fin: # pulses is of shape (NNNN, 1), of type uint8. # Squeeze out the data, and subtract the 3rd entry from the 2nd @@ -98,7 +97,7 @@ def get_gain_setting(fname: str, h5path_ctrl: str) -> int: gain-setting 1: setupr@dark=8, setupr@slopespc=40 gain-setting 0: setupr@dark=0, setupr@slopespc=32 - patternTypeIndex 1: High-gian + patternTypeIndex 1: High-gain patternTypeIndex 2: Medium-gain patternTypeIndex 3: Low-gain patternTypeIndex 4: SlopesPC @@ -134,7 +133,7 @@ def get_gain_mode(fname: str, h5path_ctrl: str) -> AgipdGainMode: h5path_run = h5path_ctrl.replace("CONTROL/", "RUN/", 1) h5path_gainmode = f'{h5path_run}/gainModeIndex/value' - with h5py.File(fname, 'r') as fd: + with h5py.File(fname, "r") as fd: if h5path_gainmode in fd: return AgipdGainMode(fd[h5path_gainmode][0]) return AgipdGainMode.ADAPTIVE_GAIN @@ -166,10 +165,15 @@ def get_bias_voltage(fname: str, karabo_id_control: str, class AgipdCorrections: - def __init__(self, max_cells, max_pulses, - h5_data_path="INSTRUMENT/SPB_DET_AGIPD1M-1/DET/{}CH0:xtdf/", - h5_index_path="INDEX/SPB_DET_AGIPD1M-1/DET/{}CH0:xtdf/", - corr_bools: Optional[dict] = None): + def __init__( + self, + max_cells, + max_pulses, + h5_data_path="INSTRUMENT/SPB_DET_AGIPD1M-1/DET/{}CH0:xtdf/", + h5_index_path="INDEX/SPB_DET_AGIPD1M-1/DET/{}CH0:xtdf/", + corr_bools: Optional[dict] = None, + gain_mode: AgipdGainMode = AgipdGainMode.ADAPTIVE_GAIN, + ): """ Initialize an AgipdCorrections Class @@ -198,9 +202,7 @@ class AgipdCorrections: const_yaml = metadata["retrieved-constants"] for mod in modules: - qm = f"Q{mod // 4 + 1}M{mod % 4 + 1}" - agipd_corr.initialize_from_yaml(karabo_da, - const_yaml, mod) + agipd_corr.initialize_from_yaml(karabo_da, const_yaml, mod) data_shape = (n_images_max, 512, 128) agipd_corr.allocate_images(data_shape, n_cores_files) @@ -226,6 +228,7 @@ class AgipdCorrections: self.pulses_lst = list(range(*max_pulses)) \ if not (len(max_pulses) == 1 and max_pulses[0] == 0) else max_pulses # noqa self.max_cells = max_cells + self.gain_mode = gain_mode # Correction parameters self.baseline_corr_noise_threshold = -1000 @@ -297,7 +300,7 @@ class AgipdCorrections: data_dict = self.shared_dict[i_proc] data_dict['moduleIdx'][0] = module_idx try: - f = h5py.File(file_name, 'r') + f = h5py.File(file_name, "r") group = f[agipd_base]["image"] (_, first_index, last_index, @@ -326,8 +329,7 @@ class AgipdCorrections: data_dict['rawgain'][:n_img] = raw_data[:, 1] data_dict['cellId'][:n_img] = allcells[firange] data_dict['pulseId'][:n_img] = allpulses[firange] - data_dict['trainId'][:n_img] = np.squeeze( - group['trainId'][:][firange]) + data_dict['trainId'][:n_img] = np.squeeze(group['trainId'][:][firange]) # noqa except Exception as e: print(f'Error during reading data from file {file_name}: {e}') print(f'Error traceback: {traceback.format_exc()}') @@ -361,14 +363,14 @@ class AgipdCorrections: return trains = data_dict['trainId'][:n_img] - with h5py.File(ofile_name, 'w') as outfile: + with h5py.File(ofile_name, "w") as outfile: # Copy any other data from the input file. # This includes indexes, so it's important that the corrected data # we write is aligned with the raw data. - with h5py.File(file_name, 'r') as infile: - self.copy_and_sanitize_non_cal_data(infile, outfile, - agipd_base, - idx_base, trains) + with h5py.File(file_name, "r") as infile: + self.copy_and_sanitize_non_cal_data( + infile, outfile, agipd_base, idx_base, trains + ) # All corrected data goes in a /INSTRUMENT/.../image group image_grp = outfile[data_path] @@ -495,35 +497,35 @@ class AgipdCorrections: gain = self.shared_dict[i_proc]['gain'][first:last] cellid = self.shared_dict[i_proc]['cellId'][first:last] - # first evaluate the gain into 0, 1, 2 --> high, medium, low - t0 = self.thresholds[module_idx][0] - t1 = self.thresholds[module_idx][1] + if self.gain_mode is AgipdGainMode.ADAPTIVE_GAIN: + # first evaluate the gain into 0, 1, 2 --> high, medium, low + t0 = self.thresholds[module_idx][0] + t1 = self.thresholds[module_idx][1] - # load raw_data and rgain to be used during gain_correction, - # if requested - if self.corr_bools.get('melt_snow'): - self.shared_dict[i_proc]['t0_rgain'][first:last] = \ - rawgain / t0[cellid, ...] - self.shared_dict[i_proc]['raw_data'][first:last] = np.copy(data) - - # Often most pixels are in high-gain, so it's more efficient to - # set the whole output block to zero than select the right pixels. - gain[:] = 0 - # exceeding first threshold means data is medium or low gain - gain[rawgain > t0[cellid, ...]] = 1 - # exceeding also second threshold means data is low gain - gain[rawgain > t1[cellid, ...]] = 2 + if self.corr_bools.get("melt_snow"): + # load raw_data and rgain to be used during gain_correction + self.shared_dict[i_proc]["t0_rgain"][first:last] = rawgain / t0[cellid, ...] # noqa + self.shared_dict[i_proc]["raw_data"][first:last] = np.copy(data) # noqa + + # Often most pixels are in high-gain, so it's more efficient to + # set the whole output block to zero than select the right pixels. + gain[:] = 0 + # exceeding first threshold means data is medium or low gain + gain[rawgain > t0[cellid, ...]] = 1 + # exceeding also second threshold means data is low gain + gain[rawgain > t1[cellid, ...]] = 2 + else: + # the enum values map 1, 2, 3 to (fixed) gain modes + gain[:] = self.gain_mode - 1 offsetb = self.offset[module_idx][:, cellid] # force into high or medium gain if requested - if self.corr_bools.get('force_mg_if_below'): - gain[(gain == 2) & ( - (data - offsetb[1]) < self.mg_hard_threshold)] = 1 + if self.corr_bools.get("force_mg_if_below"): + gain[(gain == 2) & ((data - offsetb[1]) < self.mg_hard_threshold)] = 1 # noqa - if self.corr_bools.get('force_hg_if_below'): - gain[(gain > 0) & ( - (data - offsetb[0]) < self.hg_hard_threshold)] = 0 + if self.corr_bools.get("force_hg_if_below"): + gain[(gain > 0) & ((data - offsetb[0]) < self.hg_hard_threshold)] = 0 # noqa # choose constants according to gain setting off = calgs.gain_choose(gain, offsetb) @@ -556,8 +558,9 @@ class AgipdCorrections: cellid = self.shared_dict[i_proc]['cellId'][first:last] # output is saved in sharedmem to pass for correct_agipd() # as this function takes about 3 seconds. - self.shared_dict[i_proc]['msk'][first:last] = \ - calgs.gain_choose_int(gain, self.mask[module_idx][:, cellid]) + self.shared_dict[i_proc]["msk"][first:last] = calgs.gain_choose_int( + gain, self.mask[module_idx][:, cellid] + ) if hasattr(self, "rel_gain"): # Get the correct rel_gain depending on cell-id @@ -629,14 +632,12 @@ class AgipdCorrections: # if baseline correction was not requested # msk and rel_corr will still be empty shared_mem arrays if not any(self.blc_bools): - msk = calgs.gain_choose_int(gain, - self.mask[module_idx][:, cellid]) + msk = calgs.gain_choose_int(gain, self.mask[module_idx][:, cellid]) # same for relative gain and then bad pixel mask if hasattr(self, "rel_gain"): # Get the correct rel_gain depending on cell-id - rel_corr = calgs.gain_choose(gain, - self.rel_gain[module_idx][:, cellid]) # noqa + rel_corr = calgs.gain_choose(gain, self.rel_gain[module_idx][:, cellid]) # noqa # Correct for relative gain if self.corr_bools.get("rel_gain") and hasattr(self, "rel_gain"): @@ -699,11 +700,11 @@ class AgipdCorrections: # Copy the data across into the existing shared-memory array mask[...] = msk[...] - def get_valid_image_idx(self, idx_base: str, infile: str, - index_v: Optional[int] = 2): - """ Return the indices of valid data - """ - if index_v == 2: + def get_valid_image_idx( + self, idx_base: str, infile: str, raw_format_version: int = 2 + ): + """Return the indices of valid data""" + if raw_format_version == 2: count = np.squeeze(infile[idx_base + "image/count"]) first = np.squeeze(infile[idx_base + "image/first"]) if np.count_nonzero(count != 0) == 0: @@ -728,13 +729,16 @@ class AgipdCorrections: # Creating an array of validated indices. # If all indices were validated this array will be the same, # as what is stored at /DET/image/trainId - valid_indices = np.concatenate([np.arange(validf[i], - validf[i]+validc[i]) - for i in range(validf.size)], - axis=0) + valid_indices = np.concatenate( + [ + np.arange(validf[i], validf[i] + validc[i]) + for i in range(validf.size) + ], + axis=0, + ) valid_indices = np.squeeze(valid_indices).astype(np.int32) - elif index_v == 1: + elif raw_format_version == 1: status = np.squeeze(infile[idx_base + "image/status"]) if np.count_nonzero(status != 0) == 0: raise IOError(f"File {infile} has no valid counts") @@ -752,10 +756,9 @@ class AgipdCorrections: valid_indices = None else: raise AttributeError( - f"Not a known raw format version: {index_v}") + f"Not a known raw format version: {raw_format_version}") - return (valid, first_index, last_index, idxtrains, - valid_indices) + return (valid, first_index, last_index, idxtrains, valid_indices) def apply_selected_pulses(self, i_proc: int) -> int: """Select sharedmem data indices to correct based on selected @@ -961,13 +964,11 @@ class AgipdCorrections: # Extract parameters through identifying # unique trains, index and numbers. - uq, fidxv, cntsv = np.unique(trains, return_index=True, - return_counts=True) + uq, fidxv, cntsv = np.unique(trains, return_index=True, return_counts=True) # noqa # Validate calculated CORR INDEX contents by checking # difference between trainId stored in RAW data and trains from - train_diff = np.isin(np.array(infile["/INDEX/trainId"]), uq, - invert=True) + train_diff = np.isin(np.array(infile["/INDEX/trainId"]), uq, invert=True) # noqa # Insert zeros for missing trains. # fidxv and cntsv should have same length as @@ -1094,7 +1095,8 @@ class AgipdCorrections: self.offset[module_idx][...] = cons_data["Offset"].transpose()[...] self.noise[module_idx][...] = cons_data["Noise"].transpose()[...] - self.thresholds[module_idx][...] = cons_data["ThresholdsDark"].transpose()[:3, ...] # noqa + if self.gain_mode is AgipdGainMode.ADAPTIVE_GAIN: + self.thresholds[module_idx][...] = cons_data["ThresholdsDark"].transpose()[:3,...] # noqa if self.corr_bools.get("low_medium_gap"): t0 = self.thresholds[module_idx][0] @@ -1238,20 +1240,18 @@ class AgipdCorrections: return - def initialize_from_yaml(self, karabo_da: str, - const_yaml: Dict[str, Any], - module_idx: int - ) -> Dict[str, Any]: + def initialize_from_yaml( + self, karabo_da: str, const_yaml: Dict[str, Any], module_idx: int + ) -> Dict[str, Any]: """Initialize calibration constants from a yaml file - :param karabo-da: karabo data aggerator - :param const_yaml: (Dict) from the "retrieved-constants" part of a yaml - file in pre-notebook, which consists of metadata of either the constant + :param karabo_da: a karabo data aggregator + :param const_yaml: from the "retrieved-constants" part of a yaml file + from pre-notebook, which consists of metadata of either the constant file path or the empty constant shape, and the creation-time of the retrieved constants - :param module_idx: Index of module. - :return when: Dictionary of retrieved constants with - their creation-time. + :param module_idx: Index of module + :return when: dict of retrieved constants with their creation-time """ # string of the device name. @@ -1261,10 +1261,6 @@ class AgipdCorrections: for cname, mdata in const_yaml[karabo_da]["constants"].items(): when[cname] = mdata["creation-time"] if when[cname]: - # This path is only used when testing new flat fields from - # file during development: it takes ages to test using all - # cells. Consequently, the shape needs to be fixed when less - # cells are used. with h5py.File(mdata["file-path"], "r") as cf: cons_data[cname] = np.copy(cf[f"{db_module}/{cname}/0/data"]) # noqa else: @@ -1335,15 +1331,21 @@ class AgipdCorrections: """ - const_dict = \ - assemble_constant_dict(self.corr_bools, self.pc_bools, - memory_cells, bias_voltage, gain_setting, - acquisition_rate, photon_energy, - beam_energy=None, only_dark=only_dark) - - cons_data, when = \ - self.retrieve_constant_and_time(karabo_id, karabo_da, const_dict, - cal_db_interface, creation_time) + const_dict = assemble_constant_dict( + self.corr_bools, + self.pc_bools, + memory_cells, + bias_voltage, + gain_setting, + acquisition_rate, + photon_energy, + beam_energy=None, + only_dark=only_dark, + ) + + cons_data, when = self.retrieve_constant_and_time( + karabo_id, karabo_da, const_dict, cal_db_interface, creation_time + ) self.init_constants(cons_data, when, module_idx) @@ -1357,23 +1359,17 @@ class AgipdCorrections: :param constant_shape: Shape of expected constants (gain, cells, x, y) """ for module_idx in modules: - self.offset[module_idx] = sharedmem.empty(constant_shape, - dtype='f4') - self.thresholds[module_idx] = sharedmem.empty(constant_shape, - dtype='f4') - self.noise[module_idx] = sharedmem.empty(constant_shape, - dtype='f4') - - self.md_additional_offset[module_idx] = sharedmem.empty( - constant_shape[1:], dtype='f4') - self.rel_gain[module_idx] = sharedmem.empty(constant_shape, - dtype='f4') - self.frac_high_med[module_idx] = sharedmem.empty(constant_shape[1], - dtype='f4') - - self.mask[module_idx] = sharedmem.empty(constant_shape, dtype='i4') - self.xray_cor[module_idx] = sharedmem.empty(constant_shape[1:], - dtype='f4') + self.offset[module_idx] = sharedmem.empty(constant_shape, dtype="f4") # noqa + if self.gain_mode is AgipdGainMode.ADAPTIVE_GAIN: + self.thresholds[module_idx] = sharedmem.empty(constant_shape, dtype="f4") # noqa + self.noise[module_idx] = sharedmem.empty(constant_shape, dtype="f4") # noqa + + self.md_additional_offset[module_idx] = sharedmem.empty(constant_shape[1:], dtype="f4") # noqa + self.rel_gain[module_idx] = sharedmem.empty(constant_shape, dtype="f4") # noqa + self.frac_high_med[module_idx] = sharedmem.empty(constant_shape[1], dtype="f4") # noqa + + self.mask[module_idx] = sharedmem.empty(constant_shape, dtype="i4") + self.xray_cor[module_idx] = sharedmem.empty(constant_shape[1:], dtype="f4") # noqa def allocate_images(self, shape, n_cores_files): """ @@ -1387,26 +1383,18 @@ class AgipdCorrections: self.shared_dict = [] for i in range(n_cores_files): self.shared_dict.append({}) - self.shared_dict[i]['cellId'] = sharedmem.empty(shape[0], - dtype='u2') - self.shared_dict[i]['pulseId'] = sharedmem.empty(shape[0], - dtype='u8') - self.shared_dict[i]['trainId'] = sharedmem.empty(shape[0], - dtype='u8') - self.shared_dict[i]['moduleIdx'] = sharedmem.empty(1, dtype='i4') - self.shared_dict[i]['nImg'] = sharedmem.empty(1, dtype='i4') - self.shared_dict[i]['mask'] = sharedmem.empty(shape, dtype='u4') - self.shared_dict[i]['data'] = sharedmem.empty(shape, dtype='f4') - self.shared_dict[i]['rawgain'] = sharedmem.empty(shape, - dtype='u2') - self.shared_dict[i]['gain'] = sharedmem.empty(shape, dtype='u1') - self.shared_dict[i]['blShift'] = sharedmem.empty(shape[0], - dtype='f4') + self.shared_dict[i]["cellId"] = sharedmem.empty(shape[0], dtype="u2") # noqa + self.shared_dict[i]["pulseId"] = sharedmem.empty(shape[0], dtype="u8") # noqa + self.shared_dict[i]["trainId"] = sharedmem.empty(shape[0], dtype="u8") # noqa + self.shared_dict[i]["moduleIdx"] = sharedmem.empty(1, dtype="i4") + self.shared_dict[i]["nImg"] = sharedmem.empty(1, dtype="i4") + self.shared_dict[i]["mask"] = sharedmem.empty(shape, dtype="u4") + self.shared_dict[i]["data"] = sharedmem.empty(shape, dtype="f4") + self.shared_dict[i]["rawgain"] = sharedmem.empty(shape, dtype="u2") + self.shared_dict[i]["gain"] = sharedmem.empty(shape, dtype="u1") + self.shared_dict[i]["blShift"] = sharedmem.empty(shape[0], dtype="f4") # noqa # Parameters shared between image-wise correction functions - self.shared_dict[i]['msk'] = sharedmem.empty(shape, dtype='i4') - self.shared_dict[i]['raw_data'] = sharedmem.empty(shape, - dtype='f4') - self.shared_dict[i]['rel_corr'] = sharedmem.empty(shape, - dtype='f4') - self.shared_dict[i]['t0_rgain'] = sharedmem.empty(shape, - dtype='u2') + self.shared_dict[i]["msk"] = sharedmem.empty(shape, dtype="i4") + self.shared_dict[i]["raw_data"] = sharedmem.empty(shape, dtype="f4") # noqa + self.shared_dict[i]["rel_corr"] = sharedmem.empty(shape, dtype="f4") # noqa + self.shared_dict[i]["t0_rgain"] = sharedmem.empty(shape, dtype="u2") # noqa diff --git a/cal_tools/cal_tools/agipdutils.py b/cal_tools/cal_tools/agipdutils.py index 9d58ccf83ba9cc692a043c64c562d78f1c1afb99..44cac4c87b1a0ffc8c8ba4628e38d1f68056672e 100644 --- a/cal_tools/cal_tools/agipdutils.py +++ b/cal_tools/cal_tools/agipdutils.py @@ -2,18 +2,27 @@ import copy from typing import Tuple import numpy as np -from cal_tools.enums import BadPixels, SnowResolution -from scipy.signal import cwt, ricker +from cal_tools.enums import AgipdGainMode, BadPixels, SnowResolution +from scipy.signal import cwt, find_peaks_cwt, ricker from sklearn.mixture import GaussianMixture from sklearn.preprocessing import StandardScaler -def assemble_constant_dict(corr_bools, pc_bools, memory_cells, bias_voltage, - gain_setting, acquisition_rate, - photon_energy, beam_energy=None, only_dark=False): +def assemble_constant_dict( + corr_bools, + pc_bools, + memory_cells, + bias_voltage, + gain_setting, + acquisition_rate, + photon_energy, + beam_energy=None, + only_dark=False, + gain_mode=AgipdGainMode.ADAPTIVE_GAIN, +): """ Assemble a dictionary with the iCalibrationDB constant names and - the operating conditions for retrieveing the required constants + the operating conditions for retrieving the required constants for correction. :param corr_bools: (Dict) A dict of booleans for applying @@ -23,22 +32,25 @@ def assemble_constant_dict(corr_bools, pc_bools, memory_cells, bias_voltage, :param bias_voltage: (Int) Bias Voltage :param gain_setting: (Float) Gain setting :param acquisition_rate: (Float) Acquisition rate - :param photon_energy: (Float) Photong energy + :param photon_energy: (Float) Photon energy :param beam_energy: (Float) Beam Energy - :param only_dark: (Bool) Indicating a retrieval for dark - constants only from db - :return: const_dict: (Dict) An assembeld dictionary that can be used + :param only_dark: (Bool) Indicating a retrieval for dark constants only from db + :param gain_mode: Operation mode of the detector (default to adaptive gain) + :return: const_dict: (Dict) An assembled dictionary that can be used to retrieve the required constants """ darkcond = [ - "Dark", { + "Dark", + { "memory_cells": memory_cells, "bias_voltage": bias_voltage, "acquisition_rate": acquisition_rate, "gain_setting": gain_setting, + "gain_mode": gain_mode, "pixels_x": 512, - "pixels_y": 128, } + "pixels_y": 128, + }, ] const_dict = { "Offset": ["zeros", (128, 512, memory_cells, 3), darkcond], @@ -47,28 +59,21 @@ def assemble_constant_dict(corr_bools, pc_bools, memory_cells, bias_voltage, "BadPixelsDark": ["zeros", (128, 512, memory_cells, 3), darkcond], } - if not (corr_bools.get('only_offset') or only_dark): - + if not (corr_bools.get("only_offset") or only_dark): if any(pc_bools): - const_dict["BadPixelsPC"] = \ - ["zeros", (memory_cells, 128, 512), darkcond] - const_dict["SlopesPC"] = \ - ["ones", (128, 512, memory_cells, 10), darkcond] + const_dict["BadPixelsPC"] = ["zeros", (memory_cells, 128, 512), darkcond] + const_dict["SlopesPC"] = ["ones", (128, 512, memory_cells, 10), darkcond] if corr_bools.get("xray_corr"): # Add illuminated conditions illumcond = [ - "Illuminated", { - "beam_energy": beam_energy, - "photon_energy": photon_energy - } + "Illuminated", + {"beam_energy": beam_energy, "photon_energy": photon_energy}, ] illumcond[1].update(darkcond[1]) - const_dict["BadPixelsFF"] = ["zeros", (128, 512, memory_cells), - illumcond] - const_dict["SlopesFF"] = ["ones", (128, 512, memory_cells, 2), - illumcond] + const_dict["BadPixelsFF"] = ["zeros", (128, 512, memory_cells), illumcond] + const_dict["SlopesFF"] = ["ones", (128, 512, memory_cells, 2), illumcond] return const_dict diff --git a/cal_tools/cal_tools/tools.py b/cal_tools/cal_tools/tools.py index 0ddf2c13264ee0eaddd8a26217b382376dcf80ef..7572e7d18911d92abfba522c31dcf8f453745461 100644 --- a/cal_tools/cal_tools/tools.py +++ b/cal_tools/cal_tools/tools.py @@ -539,8 +539,7 @@ def get_from_db(karabo_id: str, karabo_da: str, if ntries > 0: if load_data and meta_only: mdata_const = metadata.calibration_constant_version - fpath = Path(mdata_const.hdf5path, - mdata_const.filename) + fpath = Path(mdata_const.hdf5path, mdata_const.filename) with h5py.File(fpath, "r") as f: arr = f[f"{mdata_const.h5path}/data"][()] metadata.calibration_constant.data = arr diff --git a/notebooks/AGIPD/AGIPD_Correct_and_Verify.ipynb b/notebooks/AGIPD/AGIPD_Correct_and_Verify.ipynb index c4cbf1ef2759f6ade29cb6f4861bca470707e085..27de082393d224644b97335bc69b1c0ca92266a4 100644 --- a/notebooks/AGIPD/AGIPD_Correct_and_Verify.ipynb +++ b/notebooks/AGIPD/AGIPD_Correct_and_Verify.ipynb @@ -14,12 +14,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-21T11:30:06.730220Z", - "start_time": "2019-02-21T11:30:06.658286Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "in_folder = \"/gpfs/exfel/exp/HED/202031/p900174/raw\" # the folder to read data from, required\n", @@ -51,7 +46,7 @@ "gain_setting = 0.1 # the gain setting, use 0.1 to try to auto-determine\n", "photon_energy = 9.2 # photon energy in keV\n", "overwrite = True # set to True if existing data should be overwritten\n", - "max_pulses = [0, 500, 1] # range list [st, end, step] of maximum pulse indices within a train. 3 allowed maximum list input elements. \n", + "max_pulses = [0, 500, 1] # range list [st, end, step] of maximum pulse indices within a train. 3 allowed maximum list input elements.\n", "mem_cells_db = 0 # set to a value different than 0 to use this value for DB queries\n", "cell_id_preview = 1 # cell Id used for preview in single-shot plots\n", "\n", @@ -71,7 +66,7 @@ "xray_gain = False # do relative gain correction based on xray data\n", "blc_noise = False # if set, baseline correction via noise peak location is attempted\n", "blc_stripes = False # if set, baseline corrected via stripes\n", - "blc_hmatch = False # if set, base line correction via histogram matching is attempted \n", + "blc_hmatch = False # if set, base line correction via histogram matching is attempted\n", "match_asics = False # if set, inner ASIC borders are matched to the same signal level\n", "adjust_mg_baseline = False # adjust medium gain baseline to match highest high gain value\n", "zero_nans = False # set NaN values in corrected data to 0\n", @@ -104,21 +99,19 @@ "metadata": {}, "outputs": [], "source": [ - "import copy\n", - "import gc\n", "import itertools\n", "import math\n", + "import multiprocessing\n", "import re\n", "import traceback\n", "import warnings\n", "from datetime import timedelta\n", - "from multiprocessing import Pool\n", "from pathlib import Path\n", - "from time import perf_counter, sleep, time\n", + "from time import perf_counter\n", "\n", "import tabulate\n", "from dateutil import parser\n", - "from IPython.display import HTML, Latex, Markdown, display\n", + "from IPython.display import Latex, Markdown, display\n", "\n", "warnings.filterwarnings('ignore')\n", "import matplotlib\n", @@ -126,11 +119,8 @@ "import yaml\n", "from extra_data import RunDirectory, stack_detector_data\n", "from extra_geom import AGIPD_1MGeometry, AGIPD_500K2GGeometry\n", - "from iCalibrationDB import Detectors\n", "from matplotlib import cm as colormap\n", "from matplotlib.colors import LogNorm\n", - "from matplotlib.ticker import FormatStrFormatter, LinearLocator\n", - "from mpl_toolkits.mplot3d import Axes3D\n", "\n", "matplotlib.use(\"agg\")\n", "%matplotlib inline\n", @@ -141,22 +131,19 @@ "sns.set_context(\"paper\", font_scale=1.4)\n", "sns.set_style(\"ticks\")\n", "\n", + "import cal_tools\n", "import seaborn as sns\n", "from cal_tools.agipdlib import (\n", " AgipdCorrections,\n", " get_acq_rate,\n", + " get_gain_mode,\n", " get_gain_setting,\n", " get_num_cells,\n", ")\n", "from cal_tools.ana_tools import get_range\n", "from cal_tools.cython import agipdalgs as calgs\n", - "from cal_tools.enums import BadPixels\n", + "from cal_tools.enums import AgipdGainMode, BadPixels\n", "from cal_tools.step_timing import StepTimer\n", - "from cal_tools.tools import (\n", - " CalibrationMetadata,\n", - " get_dir_creation_date,\n", - " map_modules_from_folder,\n", - ")\n", "\n", "sns.set()\n", "sns.set_context(\"paper\", font_scale=1.4)\n", @@ -188,7 +175,7 @@ "source": [ "# Fill dictionaries comprising bools and arguments for correction and data analysis\n", "\n", - "# Here the herarichy and dependability for correction booleans are defined\n", + "# Here the hierarchy and dependability for correction booleans are defined\n", "corr_bools = {}\n", "\n", "# offset is at the bottom of AGIPD correction pyramid.\n", @@ -213,8 +200,17 @@ " corr_bools[\"common_mode\"] = common_mode\n", " corr_bools[\"melt_snow\"] = melt_snow\n", " corr_bools[\"mask_zero_std\"] = mask_zero_std\n", - " corr_bools[\"low_medium_gap\"] = low_medium_gap \n", - " " + " corr_bools[\"low_medium_gap\"] = low_medium_gap\n", + "\n", + "# Many corrections don't apply to fixed gain mode; will explicitly disable later if detected\n", + "disable_for_fixed_gain = [\n", + " \"adjust_mg_baseline\",\n", + " \"blc_set_min\",\n", + " \"force_hg_if_below\",\n", + " \"force_mg_if_below\",\n", + " \"low_medium_gap\",\n", + " \"melt_snow\"\n", + "]" ] }, { @@ -262,12 +258,8 @@ " karabo_da = [\"AGIPD{:02d}\".format(i) for i in modules]\n", "else:\n", " modules = [int(x[-2:]) for x in karabo_da]\n", - " \n", - "def mod_name(modno):\n", - " return f\"Q{modno // 4 + 1}M{modno % 4 + 1}\"\n", "\n", - "print(\"Process modules: \", ', '.join(\n", - " [mod_name(x) for x in modules]))\n", + "print(\"Process modules:\", ', '.join(cal_tools.tools.module_index_to_qm(x) for x in modules))\n", "print(f\"Detector in use is {karabo_id}\")\n", "print(f\"Instrument {instrument}\")\n", "print(f\"Detector instance {dinstance}\")" @@ -280,17 +272,17 @@ "outputs": [], "source": [ "# Display Information about the selected pulses indices for correction.\n", - "pulses_lst = list(range(*max_pulses)) if not (len(max_pulses)==1 and max_pulses[0]==0) else max_pulses \n", + "pulses_lst = list(range(*max_pulses)) if not (len(max_pulses)==1 and max_pulses[0]==0) else max_pulses\n", "\n", "try:\n", - " if len(pulses_lst) > 1: \n", + " if len(pulses_lst) > 1:\n", " print(\"A range of {} pulse indices is selected: from {} to {} with a step of {}\"\n", " .format(len(pulses_lst), pulses_lst[0] , pulses_lst[-1] + (pulses_lst[1] - pulses_lst[0]),\n", " pulses_lst[1] - pulses_lst[0]))\n", " else:\n", - " print(\"one pulse is selected: a pulse of idx {}\".format(pulses_lst[0]))\n", + " print(f\"one pulse is selected: a pulse of idx {pulses_lst[0]}\")\n", "except Exception as e:\n", - " raise ValueError('max_pulses input Error: {}'.format(e))" + " raise ValueError(f\"max_pulses input Error: {e}\")" ] }, { @@ -300,9 +292,9 @@ "outputs": [], "source": [ "# set everything up filewise\n", - "mmf = map_modules_from_folder(str(in_folder), run, path_template,\n", - " karabo_da, sequences)\n", - "mapped_files, mod_ids, total_sequences, sequences_qm, _ = mmf\n", + "mapped_files, _, total_sequences, _, _ = cal_tools.tools.map_modules_from_folder(\n", + " str(in_folder), run, path_template, karabo_da, sequences\n", + ")\n", "file_list = []\n", "\n", "# ToDo: Split table over pages\n", @@ -357,12 +349,11 @@ "# Evaluate creation time\n", "creation_time = None\n", "if use_dir_creation_date:\n", - " creation_time = get_dir_creation_date(str(in_folder), run)\n", + " creation_time = cal_tools.tools.get_dir_creation_date(str(in_folder), run)\n", " offset = parser.parse(creation_date_offset)\n", - " delta = timedelta(hours=offset.hour,\n", - " minutes=offset.minute, seconds=offset.second)\n", + " delta = timedelta(hours=offset.hour, minutes=offset.minute, seconds=offset.second)\n", " creation_time += delta\n", - " \n", + "\n", "# Evaluate gain setting\n", "if gain_setting == 0.1:\n", " if creation_time.replace(tzinfo=None) < parser.parse('2020-01-31'):\n", @@ -376,7 +367,9 @@ " print(e)\n", " print(\"Set gain setting to 0\")\n", " gain_setting = 0\n", - " " + "\n", + "# Evaluate gain mode (operation mode)\n", + "gain_mode = get_gain_mode(control_fn, h5path_ctrl)" ] }, { @@ -386,8 +379,26 @@ "outputs": [], "source": [ "print(f\"Using {creation_time} as creation time\")\n", - "print(f\"Operating conditions are:\\n• Bias voltage: {bias_voltage}\\n• Memory cells: {mem_cells_db}\\n\"\n", - " f\"• Acquisition rate: {acq_rate}\\n• Gain setting: {gain_setting}\\n• Photon Energy: {photon_energy}\\n\")" + "print(\"Operating conditions are:\")\n", + "print(f\"• Bias voltage: {bias_voltage}\")\n", + "print(f\"• Memory cells: {mem_cells_db}\")\n", + "print(f\"• Acquisition rate: {acq_rate}\")\n", + "print(f\"• Gain setting: {gain_setting}\")\n", + "print(f\"• Gain mode: {gain_mode.name}\")\n", + "print(f\"• Photon Energy: {photon_energy}\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "if gain_mode:\n", + " for to_disable in disable_for_fixed_gain:\n", + " if corr_bools.get(to_disable, False):\n", + " print(f\"Warning: {to_disable} correction was requested, but does not apply to fixed gain mode\")\n", + " corr_bools[to_disable] = False" ] }, { @@ -403,10 +414,14 @@ "metadata": {}, "outputs": [], "source": [ - "agipd_corr = AgipdCorrections(max_cells, max_pulses,\n", - " h5_data_path=h5path,\n", - " h5_index_path=h5path_idx,\n", - " corr_bools=corr_bools)\n", + "agipd_corr = AgipdCorrections(\n", + " max_cells,\n", + " max_pulses,\n", + " h5_data_path=h5path,\n", + " h5_index_path=h5path_idx,\n", + " corr_bools=corr_bools,\n", + " gain_mode=gain_mode\n", + ")\n", "\n", "agipd_corr.baseline_corr_noise_threshold = -blc_noise_threshold\n", "agipd_corr.hg_hard_threshold = hg_hard_threshold\n", @@ -420,6 +435,15 @@ "agipd_corr.ff_gain = ff_gain" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "module_index_to_karabo_da = {mod: da for (mod, da) in zip(modules, karabo_da)}" + ] + }, { "cell_type": "code", "execution_count": null, @@ -429,28 +453,37 @@ "# Retrieve calibration constants to RAM\n", "agipd_corr.allocate_constants(modules, (3, mem_cells_db, 512, 128))\n", "\n", - "metadata = CalibrationMetadata(out_folder)\n", + "metadata = cal_tools.tools.CalibrationMetadata(out_folder)\n", "# NOTE: this notebook will not overwrite calibration metadata file\n", "const_yaml = metadata.get(\"retrieved-constants\", {})\n", "\n", - "# retrive constants\n", "def retrieve_constants(mod):\n", " \"\"\"\n", " Retrieve calibration constants and load them to shared memory\n", - " \n", + "\n", " Metadata for constants is taken from yml file or retrieved from the DB\n", " \"\"\"\n", - " err = ''\n", - " # TODO: parallelize over modules\n", - " k_da = karabo_da[mod]\n", + " err = \"\"\n", + " k_da = module_index_to_karabo_da[mod]\n", " try:\n", " # check if there is a yaml file in out_folder that has the device constants.\n", " if k_da in const_yaml:\n", " when = agipd_corr.initialize_from_yaml(k_da, const_yaml, mod)\n", " else:\n", - " # TODO: should we save what is found here in metadata?\n", - " when = agipd_corr.initialize_from_db(karabo_id, k_da, cal_db_interface, creation_time, mem_cells_db, bias_voltage,\n", - " photon_energy, gain_setting, acq_rate, mod, False)\n", + " # TODO: replace with proper retrieval (as done in pre-correction)\n", + " when = agipd_corr.initialize_from_db(\n", + " karabo_id,\n", + " k_da,\n", + " cal_db_interface,\n", + " creation_time,\n", + " mem_cells_db,\n", + " bias_voltage,\n", + " photon_energy,\n", + " gain_setting,\n", + " acq_rate,\n", + " mod,\n", + " False,\n", + " )\n", " except Exception as e:\n", " err = f\"Error: {e}\\nError traceback: {traceback.format_exc()}\"\n", " when = None\n", @@ -458,7 +491,7 @@ "\n", "\n", "ts = perf_counter()\n", - "with Pool(processes=len(modules)) as pool:\n", + "with multiprocessing.Pool(processes=len(modules)) as pool:\n", " const_out = pool.map(retrieve_constants, modules)\n", "print(f\"Constants were loaded in {perf_counter()-ts:.01f}s\")" ] @@ -470,7 +503,7 @@ "outputs": [], "source": [ "# allocate memory for images and hists\n", - "n_images_max = max_cells*256\n", + "n_images_max = max_cells * 256\n", "data_shape = (n_images_max, 512, 128)\n", "agipd_corr.allocate_images(data_shape, n_cores_files)" ] @@ -497,7 +530,7 @@ "source": [ "def imagewise_chunks(img_counts):\n", " \"\"\"Break up the loaded data into chunks of up to chunk_size\n", - " \n", + "\n", " Yields (file data slot, start index, stop index)\n", " \"\"\"\n", " for i_proc, n_img in enumerate(img_counts):\n", @@ -518,17 +551,13 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "scrolled": true - }, + "metadata": {}, "outputs": [], "source": [ - "with Pool() as pool:\n", + "with multiprocessing.Pool() as pool:\n", " for file_batch in batches(file_list, n_cores_files):\n", " # TODO: Move some printed output to logging or similar\n", - " print(f\"Processing next {len(file_batch)} files:\")\n", - " for file_name in file_batch:\n", - " print(\" \", file_name)\n", + " print(f\"Processing next {len(file_batch)} files\")\n", " step_timer.start()\n", " img_counts = pool.starmap(agipd_corr.read_file, zip(range(len(file_batch)), file_batch,\n", " [not common_mode]*len(file_batch)))\n", @@ -549,7 +578,7 @@ " # Perform image-wise correction\n", " pool.starmap(agipd_corr.baseline_correction, imagewise_chunks(img_counts))\n", " step_timer.done_step(\"Base-line shift correction\")\n", - " \n", + "\n", " if common_mode:\n", " # Perform cross-file correction parallel over asics\n", " pool.starmap(agipd_corr.cm_correction, itertools.product(\n", @@ -587,9 +616,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "scrolled": true - }, + "metadata": {}, "outputs": [], "source": [ "# if the yml file contains \"retrieved-constants\", that means a leading\n", @@ -598,7 +625,7 @@ "timestamps = {}\n", "\n", "for i, (error, modno, when, k_da) in enumerate(const_out):\n", - " qm = mod_name(modno)\n", + " qm = cal_tools.tools.module_index_to_qm(modno)\n", " # expose errors while applying correction\n", " if error:\n", " print(\"Error: {}\".format(error) )\n", @@ -607,7 +634,7 @@ " if fst_print:\n", " print(\"Constants are retrieved with creation time: \")\n", " fst_print = False\n", - " \n", + "\n", " module_timestamps = {}\n", "\n", " # If correction is crashed\n", @@ -655,8 +682,7 @@ " Z = data.T\n", "\n", " # Plot the surface.\n", - " surf = ax.plot_surface(X, Y, Z, cmap=colormap.coolwarm,\n", - " linewidth=0, antialiased=False)\n", + " ax.plot_surface(X, Y, Z, cmap=colormap.coolwarm, linewidth=0, antialiased=False)\n", " ax.set_xlabel(x_axis)\n", " ax.set_ylabel(y_axis)\n", " ax.set_zlabel(\"Counts\")\n", @@ -683,7 +709,7 @@ "source": [ "def get_trains_data(run_folder, source, include, detector_id, tid=None, modules=16, fillvalue=np.nan):\n", " \"\"\"Load single train for all module\n", - " \n", + "\n", " :param run_folder: Path to folder with data\n", " :param source: Data source to be loaded\n", " :param include: Inset of file name to be considered\n", @@ -696,7 +722,7 @@ " tid, data = run_data.select(f'{detector_id}/DET/*', source).train_from_id(tid)\n", " else:\n", " tid, data = next(iter(run_data.select(f'{detector_id}/DET/*', source).trains(require_all=True)))\n", - " \n", + "\n", " return tid, stack_detector_data(train=data, data=source, fillvalue=fillvalue, modules=modules)" ] }, @@ -794,7 +820,7 @@ "print(f\"Gain statistics in %\")\n", "table = [[f'{gains[gains==0].size/gains.size*100:.02f}',\n", " f'{gains[gains==1].size/gains.size*100:.03f}',\n", - " f'{gains[gains==2].size/gains.size*100:.03f}']] \n", + " f'{gains[gains==2].size/gains.size*100:.03f}']]\n", "md = display(Latex(tabulate.tabulate(table, tablefmt='latex',\n", " headers=[\"High\", \"Medium\", \"Low\"])))" ] @@ -809,9 +835,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "scrolled": false - }, + "metadata": {}, "outputs": [], "source": [ "pulse_range = [np.min(pulseId[pulseId>=0]), np.max(pulseId[pulseId>=0])]\n", @@ -885,12 +909,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:33.226396Z", - "start_time": "2019-02-18T17:29:27.027758Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", @@ -926,12 +945,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:33.761015Z", - "start_time": "2019-02-18T17:29:33.227922Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", @@ -944,24 +958,19 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:35.903487Z", - "start_time": "2019-02-18T17:29:33.762568Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", "ax = fig.add_subplot(111)\n", "vmin, vmax = get_range(corrected[cell_id_preview], 5, -50)\n", "nbins = np.int((vmax + 50) / 2)\n", - "h = ax.hist(corrected[cell_id_preview].flatten(), \n", - " bins=nbins, range=(-50, vmax), \n", + "h = ax.hist(corrected[cell_id_preview].flatten(),\n", + " bins=nbins, range=(-50, vmax),\n", " histtype='stepfilled', log=True)\n", - "_ = plt.xlabel('[ADU]')\n", - "_ = plt.ylabel('Counts')\n", - "_ = ax.grid()" + "plt.xlabel('[ADU]')\n", + "plt.ylabel('Counts')\n", + "ax.grid()" ] }, { @@ -971,18 +980,13 @@ "outputs": [], "source": [ "display(Markdown('### Mean CORRECTED Preview ###\\n'))\n", - "display(Markdown(f'A mean across one train \\n'))" + "display(Markdown(f'A mean across one train\\n'))" ] }, { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:39.369686Z", - "start_time": "2019-02-18T17:29:35.905152Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", @@ -995,12 +999,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:49.217848Z", - "start_time": "2019-02-18T17:29:39.371232Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", @@ -1012,16 +1011,16 @@ "nbins = np.int((vmax + 100) / 5)\n", "h = ax.hist(corrected.flatten(), bins=nbins,\n", " range=(-100, vmax), histtype='step', log=True, label = 'All')\n", - "_ = ax.hist(corrected[gains == 0].flatten(), bins=nbins, range=(-100, vmax),\n", - " alpha=0.5, log=True, label='High gain', color='green')\n", - "_ = ax.hist(corrected[gains == 1].flatten(), bins=nbins, range=(-100, vmax),\n", - " alpha=0.5, log=True, label='Medium gain', color='red')\n", - "_ = ax.hist(corrected[gains == 2].flatten(), bins=nbins,\n", - " range=(-100, vmax), alpha=0.5, log=True, label='Low gain', color='yellow')\n", - "_ = ax.legend()\n", - "_ = ax.grid()\n", - "_ = plt.xlabel('[ADU]')\n", - "_ = plt.ylabel('Counts')" + "ax.hist(corrected[gains == 0].flatten(), bins=nbins, range=(-100, vmax),\n", + " alpha=0.5, log=True, label='High gain', color='green')\n", + "ax.hist(corrected[gains == 1].flatten(), bins=nbins, range=(-100, vmax),\n", + " alpha=0.5, log=True, label='Medium gain', color='red')\n", + "ax.hist(corrected[gains == 2].flatten(), bins=nbins, range=(-100, vmax),\n", + " alpha=0.5, log=True, label='Low gain', color='yellow')\n", + "ax.legend()\n", + "ax.grid()\n", + "plt.xlabel('[ADU]')\n", + "plt.ylabel('Counts')" ] }, { @@ -1037,12 +1036,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:49.641675Z", - "start_time": "2019-02-18T17:29:49.224167Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", @@ -1053,9 +1047,7 @@ }, { "cell_type": "markdown", - "metadata": { - "collapsed": true - }, + "metadata": {}, "source": [ "## Bad Pixels ##\n", "The mask contains dedicated entries for all pixels and memory cells as well as all three gains stages. Each mask entry is encoded in 32 bits as:" @@ -1064,12 +1056,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:49.651913Z", - "start_time": "2019-02-18T17:29:49.643556Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "table = []\n", @@ -1092,24 +1079,17 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:50.086169Z", - "start_time": "2019-02-18T17:29:49.653391Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", "ax = fig.add_subplot(111)\n", - "ax = geom.plot_data_fast(np.log2(mask[cell_id_preview]), ax=ax, vmin=0, vmax=32, cmap=\"jet\")" + "geom.plot_data_fast(np.log2(mask[cell_id_preview]), ax=ax, vmin=0, vmax=32, cmap=\"jet\")" ] }, { "cell_type": "markdown", - "metadata": { - "collapsed": true - }, + "metadata": {}, "source": [ "### Percentage of Bad Pixels across one train ###" ] @@ -1117,18 +1097,12 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:51.686562Z", - "start_time": "2019-02-18T17:29:50.088883Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", "ax = fig.add_subplot(111)\n", - "ax = geom.plot_data_fast(np.mean(mask>0, axis=0),\n", - " vmin=0, ax=ax, vmax=1, cmap=\"jet\")" + "geom.plot_data_fast(np.mean(mask>0, axis=0), vmin=0, ax=ax, vmax=1, cmap=\"jet\")" ] }, { @@ -1141,12 +1115,7 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-18T17:29:55.483270Z", - "start_time": "2019-02-18T17:29:53.664226Z" - } - }, + "metadata": {}, "outputs": [], "source": [ "fig = plt.figure(figsize=(20, 10))\n", diff --git a/notebooks/AGIPD/AGIPD_Retrieve_Constants_Precorrection.ipynb b/notebooks/AGIPD/AGIPD_Retrieve_Constants_Precorrection.ipynb index 8126a0a01a063d298eb2d950f61978d08f923d4f..b92a901a9c82422651601445255b6ab2df4b0d96 100644 --- a/notebooks/AGIPD/AGIPD_Retrieve_Constants_Precorrection.ipynb +++ b/notebooks/AGIPD/AGIPD_Retrieve_Constants_Precorrection.ipynb @@ -14,15 +14,9 @@ { "cell_type": "code", "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-21T11:30:06.730220Z", - "start_time": "2019-02-21T11:30:06.658286Z" - } - }, + "metadata": {}, "outputs": [], "source": [ - "cluster_profile = \"noDB\"\n", "in_folder = \"/gpfs/exfel/exp/SPB/202030/p900119/raw\" # the folder to read data from, required\n", "out_folder = \"/gpfs/exfel/data/scratch/ahmedk/test/AGIPD_\" # the folder to output to, required\n", "sequences = [-1] # sequences to correct, set to -1 for all, range allowed\n", @@ -57,7 +51,7 @@ "xray_gain = True # do relative gain correction based on xray data\n", "blc_noise = False # if set, baseline correction via noise peak location is attempted\n", "blc_stripes = False # if set, baseline corrected via stripes\n", - "blc_hmatch = False # if set, base line correction via histogram matching is attempted \n", + "blc_hmatch = False # if set, base line correction via histogram matching is attempted\n", "match_asics = False # if set, inner ASIC borders are matched to the same signal level\n", "adjust_mg_baseline = False # adjust medium gain baseline to match highest high gain value" ] @@ -69,8 +63,7 @@ "outputs": [], "source": [ "# Fill dictionaries comprising bools and arguments for correction and data analysis\n", - "\n", - "# Here the herarichy and dependability for correction booleans are defined \n", + "# Here the hierarichy and dependencies for correction booleans are defined \n", "corr_bools = {}\n", "\n", "# offset is at the bottom of AGIPD correction pyramid.\n", @@ -91,27 +84,18 @@ "metadata": {}, "outputs": [], "source": [ - "import sys\n", - "from collections import OrderedDict\n", - "from functools import partial\n", "from typing import List, Tuple\n", "\n", - "import h5py\n", "import matplotlib\n", "import numpy as np\n", "\n", "matplotlib.use(\"agg\")\n", - "import multiprocessing as mp\n", + "import multiprocessing\n", "from datetime import timedelta\n", "from pathlib import Path\n", "\n", "import matplotlib.pyplot as plt\n", - "from cal_tools.agipdlib import get_gain_setting\n", - "from cal_tools.tools import (\n", - " CalibrationMetadata,\n", - " get_dir_creation_date,\n", - " map_modules_from_folder,\n", - ")\n", + "from cal_tools import agipdlib, tools\n", "from dateutil import parser\n", "from iCalibrationDB import Conditions, Constants, Detectors" ] @@ -125,7 +109,7 @@ "# slopes_ff_from_files left as str for now\n", "in_folder = Path(in_folder)\n", "out_folder = Path(out_folder)\n", - "metadata = CalibrationMetadata(out_folder)" + "metadata = tools.CalibrationMetadata(out_folder)" ] }, { @@ -138,7 +122,7 @@ "\n", "creation_time = None\n", "if use_dir_creation_date:\n", - " creation_time = get_dir_creation_date(str(in_folder), run)\n", + " creation_time = tools.get_dir_creation_date(str(in_folder), run)\n", " offset = parser.parse(creation_date_offset)\n", " delta = timedelta(hours=offset.hour, minutes=offset.minute, seconds=offset.second)\n", " creation_time += delta\n", @@ -150,13 +134,7 @@ "print(f\"Outputting to {out_folder}\")\n", "out_folder.mkdir(parents=True, exist_ok=True)\n", "\n", - "import warnings\n", - "\n", - "warnings.filterwarnings('ignore')\n", - "\n", - "from cal_tools.agipdlib import SnowResolution\n", - "\n", - "melt_snow = False if corr_bools[\"only_offset\"] else SnowResolution.NONE" + "melt_snow = False if corr_bools[\"only_offset\"] else agipdlib.SnowResolution.NONE" ] }, { @@ -174,14 +152,18 @@ " gain_setting = None\n", " else:\n", " try:\n", - " gain_setting = get_gain_setting(str(control_fn), h5path_ctrl)\n", + " gain_setting = agipdlib.get_gain_setting(str(control_fn), h5path_ctrl)\n", " except Exception as e:\n", " print(f'ERROR: while reading gain setting from: \\n{control_fn}')\n", " print(e)\n", " print(\"Set gain setting to 0\")\n", " gain_setting = 0\n", "\n", + "# Evaluate gain mode (operation mode)\n", + "gain_mode = agipdlib.get_gain_mode(control_fn, h5path_ctrl)\n", + " \n", "print(f\"Gain setting: {gain_setting}\")\n", + "print(f\"Gain mode: {gain_mode.name}\")\n", "print(f\"Detector in use is {karabo_id}\")\n", "\n", "\n", @@ -210,24 +192,6 @@ " modules = [int(x[-2:]) for x in karabo_da]" ] }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "ExecuteTime": { - "end_time": "2019-02-21T11:30:07.974174Z", - "start_time": "2019-02-21T11:30:07.914832Z" - } - }, - "outputs": [], - "source": [ - "# set everything up filewise\n", - "print(f\"Checking the files before retrieving constants\")\n", - "mmf = map_modules_from_folder(str(in_folder), run, path_template, karabo_da, sequences)\n", - "\n", - "mapped_files, mod_ids, total_sequences, sequences_qm, _ = mmf" - ] - }, { "cell_type": "markdown", "metadata": {}, @@ -241,30 +205,12 @@ "metadata": {}, "outputs": [], "source": [ - "def retrieve_constants(karabo_id: str, bias_voltage: int, max_cells: float,\n", - " acq_rate: float, gain_setting: float, photon_energy: float,\n", - " only_dark: bool, nodb_with_dark: bool, \n", - " cal_db_interface: str, creation_time: str, \n", - " corr_bools: dict, pc_bools: List[bool],\n", - " inp: Tuple[str, str, str, int]\n", - " ) -> Tuple[str, str, float, float, str, dict]:\n", + "def retrieve_constants(\n", + " qm_files: List[Path], qm: str, karabo_da: str, idx: int\n", + ") -> Tuple[str, str, float, float, str, dict]:\n", " \"\"\"\n", - " Retreive constant for each module in parallel and produce a dictionary\n", - " with the creation-time and constant file path.\n", - " \n", - " :param karabo_id: (STR) Karabo ID\n", - " :param bias_voltage: (FLOAT) Bias Voltage\n", - " :param max_cells: (INT) Memory cells\n", - " :param acq_rate: (FLOAT) Acquisition Rate\n", - " :param gain_setting: (FLOAT) Gain setting\n", - " :param photon_energy: (FLOAT) Photon Energy\n", - " :param only_dark: (BOOL) only retrieve dark constants\n", - " :param nodb_with_dark: (BOOL) no constant retrieval even for dark\n", - " :param cal_db_interface: (STR) the database interface port\n", - " :param creation_time: (STR) raw data creation time\n", - " :param corr_bools: (DICT) A dictionary with bools for applying requested corrections\n", - " :param pc_bools: (LIST) list of bools to retrieve pulse capacitor constants\n", - " :param inp: (LIST) input for the parallel cluster of the partial function\n", + " Retrieve constants for a module.\n", + "\n", " :return:\n", " qm: module virtual name i.e. Q1M1.\n", " karabo_da: karabo data aggregator.\n", @@ -274,78 +220,104 @@ " mdata_dict: (DICT) dictionary with the metadata for the retrieved constants.\n", " \"\"\"\n", "\n", - " import sys\n", - " import traceback\n", - "\n", - " import numpy as np\n", - " from cal_tools.agipdlib import get_acq_rate, get_num_cells\n", - " from cal_tools.agipdutils import assemble_constant_dict\n", - " from cal_tools.tools import get_from_db\n", - " from iCalibrationDB import Conditions, Constants, Detectors\n", - "\n", " err = None\n", - "\n", - " qm_files, qm, karabo_da, idx = inp\n", - " # get number of memory cells from a sequence file with image data\n", - " for f in qm_files:\n", - " if not max_cells:\n", - " max_cells = get_num_cells(f, karabo_id, idx)\n", - " if max_cells is None:\n", - " if f != qm_files[-1]:\n", - " continue\n", - " else:\n", - " raise ValueError(f\"No raw images found for {qm} for all sequences\")\n", - " else:\n", - " cells = np.arange(max_cells)\n", - " # get out of the loop,\n", - " # if max_cells is successfully calculated. \n", + " if max_cells != 0:\n", + " # either use overriding notebook parameter\n", + " local_max_cells = max_cells\n", + " else:\n", + " # or look around in sequence files\n", + " for f in qm_files:\n", + " local_max_cells = agipdlib.get_num_cells(f, karabo_id, idx)\n", + " if local_max_cells is not None:\n", " break\n", + " # maybe we never found this in a sequence file...\n", + " if local_max_cells is None:\n", + " raise ValueError(f\"No raw images found for {qm} for all sequences\")\n", "\n", - " if acq_rate == 0.:\n", - " acq_rate = get_acq_rate((f, karabo_id, idx))\n", - "\n", - " # avoid creating retireving constant, if requested.\n", - " if not nodb_with_dark:\n", - " const_dict = assemble_constant_dict(corr_bools, pc_bools, max_cells, bias_voltage,\n", - " gain_setting, acq_rate, photon_energy,\n", - " beam_energy=None, only_dark=only_dark)\n", - "\n", - " # Retrieve multiple constants through an input dictionary\n", - " # to return a dict of useful metadata.\n", - " mdata_dict = dict()\n", - " mdata_dict['constants'] = dict()\n", - " mdata_dict['physical-detector-unit'] = None # initialization\n", - "\n", - " for cname, cval in const_dict.items():\n", - " # saving metadata in a dict\n", - " mdata_dict['constants'][cname] = dict()\n", + " if acq_rate == 0:\n", + " local_acq_rate = agipdlib.get_acq_rate(fast_paths=(f, karabo_id, idx))\n", + " else:\n", + " local_acq_rate = acq_rate\n", + "\n", + " # avoid retrieving constant, if requested.\n", + " if nodb_with_dark:\n", + " return\n", + "\n", + " const_dict = agipdlib.assemble_constant_dict(\n", + " corr_bools,\n", + " pc_bools,\n", + " local_max_cells,\n", + " bias_voltage,\n", + " gain_setting,\n", + " local_acq_rate,\n", + " photon_energy,\n", + " gain_mode=gain_mode,\n", + " beam_energy=None,\n", + " only_dark=only_dark,\n", + " )\n", + "\n", + " # Retrieve multiple constants through an input dictionary\n", + " # to return a dict of useful metadata.\n", + " mdata_dict = dict()\n", + " mdata_dict[\"constants\"] = dict()\n", + " mdata_dict[\"physical-detector-unit\"] = None # initialization\n", + "\n", + " for const_name, (const_init_fun, const_shape, (cond_type, cond_param)) in const_dict.items():\n", + " if gain_mode and const_name in (\"ThresholdsDark\",):\n", + " continue\n", + " \n", + " # saving metadata in a dict\n", + " const_mdata = dict()\n", + " mdata_dict[\"constants\"][const_name] = const_mdata\n", + "\n", + " if slopes_ff_from_files and const_name in [\"SlopesFF\", \"BadPixelsFF\"]:\n", + " const_mdata[\"file-path\"] = f\"{slopes_ff_from_files}/slopesff_bpmask_module_{qm}.h5\"\n", + " const_mdata[\"creation-time\"] = \"00:00:00\"\n", + " continue\n", + " \n", + " if gain_mode and const_name in (\"BadPixelsPC\", \"SlopesPC\", \"BadPixelsFF\", \"SlopesFF\"):\n", + " param_copy = cond_param.copy()\n", + " del param_copy[\"gain_mode\"]\n", + " condition = getattr(Conditions, cond_type).AGIPD(**param_copy)\n", + " else:\n", + " condition = getattr(Conditions, cond_type).AGIPD(**cond_param)\n", + "\n", + " _, mdata = tools.get_from_db(\n", + " karabo_id,\n", + " karabo_da,\n", + " getattr(Constants.AGIPD, const_name)(),\n", + " condition,\n", + " getattr(np, const_init_fun)(const_shape),\n", + " cal_db_interface,\n", + " creation_time,\n", + " meta_only=True,\n", + " verbosity=0,\n", + " )\n", + " mdata_const = mdata.calibration_constant_version\n", + " # check if constant was sucessfully retrieved.\n", + " if mdata.comm_db_success:\n", + " const_mdata[\"file-path\"] = (\n", + " f\"{mdata_const.hdf5path}\" f\"{mdata_const.filename}\"\n", + " )\n", + " const_mdata[\"creation-time\"] = f\"{mdata_const.begin_at}\"\n", + " mdata_dict[\"physical-detector-unit\"] = mdata_const.device_name\n", + " else:\n", + " const_mdata[\"file-path\"] = const_dict[const_name][:2]\n", + " const_mdata[\"creation-time\"] = None\n", "\n", - " if slopes_ff_from_files and cname in [\"SlopesFF\", \"BadPixelsFF\"]:\n", - " mdata_dict['constants'][cname][\"file-path\"] = f\"{slopes_ff_from_files}/slopesff_bpmask_module_{qm}.h5\"\n", - " mdata_dict['constants'][cname][\"creation-time\"] = \"00:00:00\"\n", - " else:\n", - " try:\n", - " condition = getattr(Conditions, cval[2][0]).AGIPD(**cval[2][1])\n", - " co, mdata = \\\n", - " get_from_db(karabo_id, karabo_da, getattr(Constants.AGIPD, cname)(),\n", - " condition, getattr(np, cval[0])(cval[1]),\n", - " cal_db_interface, creation_time, meta_only=True, verbosity=0)\n", - " mdata_const = mdata.calibration_constant_version\n", - " device_name = mdata.calibration_constant_version.device_name\n", - " # check if constant was sucessfully retrieved.\n", - " if mdata.comm_db_success:\n", - " mdata_dict['constants'][cname][\"file-path\"] = f\"{mdata_const.hdf5path}\" \\\n", - " f\"{mdata_const.filename}\"\n", - " mdata_dict['constants'][cname][\"creation-time\"] = f\"{mdata_const.begin_at}\"\n", - " mdata_dict['physical-detector-unit'] = mdata_const.device_name\n", - " else:\n", - " mdata_dict['constants'][cname][\"file-path\"] = const_dict[cname][:2]\n", - " mdata_dict['constants'][cname][\"creation-time\"] = None\n", - " except Exception as e:\n", - " err = f\"Error: {e}, Traceback: {traceback.format_exc()}\"\n", - " print(err)\n", - "\n", - " return qm, mdata_dict, karabo_da, acq_rate, max_cells, err\n", + " return qm, mdata_dict, karabo_da, acq_rate, local_max_cells, err" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# set everything up filewise\n", + "mapped_files, _, _, _, _ = tools.map_modules_from_folder(\n", + " str(in_folder), run, path_template, karabo_da, sequences\n", + ")\n", "\n", "pc_bools = [corr_bools.get(\"rel_gain\"),\n", " corr_bools.get(\"adjust_mg_baseline\"),\n", @@ -358,52 +330,63 @@ "only_dark = False\n", "nodb_with_dark = False\n", "if not nodb:\n", - " only_dark=(calfile != \"\")\n", + " only_dark = (calfile != \"\")\n", "if calfile != \"\" and not corr_bools[\"only_offset\"]:\n", " nodb_with_dark = nodb\n", "\n", "# A dict to connect virtual device\n", "# to actual device name.\n", - "for i, k_da in zip(modules, karabo_da):\n", - " qm = f\"Q{i//4+1}M{i%4+1}\"\n", + "for module_index, k_da in zip(modules, karabo_da):\n", + " qm = tools.module_index_to_qm(module_index)\n", " if qm in mapped_files and not mapped_files[qm].empty():\n", " device = getattr(getattr(Detectors, dinstance), qm)\n", - " qm_files = [str(mapped_files[qm].get()) for _ in range(mapped_files[qm].qsize())]\n", - "\n", + " # TODO: make map_modules_from_folder just return list(s)\n", + " qm_files = [Path(mapped_files[qm].get()) for _ in range(mapped_files[qm].qsize())]\n", " else:\n", - " print(f\"Skipping {qm}\")\n", " continue\n", "\n", - " inp.append((qm_files, qm, k_da, i))\n", - "\n", - "p = partial(retrieve_constants, karabo_id, bias_voltage, max_cells, \n", - " acq_rate, gain_setting, photon_energy, only_dark, nodb_with_dark, \n", - " cal_db_interface, creation_time, \n", - " corr_bools, pc_bools)\n", - "\n", - "with mp.Pool(processes=nmods) as pool:\n", - " results = pool.map(p, inp)\n", - "\n", + " inp.append((qm_files, qm, k_da, module_index))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "with multiprocessing.Pool(processes=nmods) as pool:\n", + " results = pool.starmap(retrieve_constants, inp)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ "mod_dev = dict()\n", "mdata_dict = dict()\n", - "for r in results:\n", - " if r:\n", - " qm, md_dict, karabo_da, acq_rate, max_cells, err = r\n", - " mod_dev[karabo_da] = {\"mod\": qm, \"err\": err}\n", - " if err:\n", - " print(f\"Error for module {qm}: {err}\")\n", - " mdata_dict[karabo_da] = md_dict\n", + "for qm, md_dict, karabo_da, acq_rate, max_cells, err in results:\n", + " mod_dev[karabo_da] = {\"mod\": qm, \"err\": err}\n", + " if err:\n", + " print(f\"Error for module {qm}: {err}\")\n", + " mdata_dict[karabo_da] = md_dict\n", "# check if it is requested not to retrieve any constants from the database\n", - "if not nodb_with_dark:\n", - " metadata.update({\"retrieved-constants\": mdata_dict})\n", - " \n", - " print(\"\\nRetrieved constants for modules: \",\n", - " f\"{[', '.join([f'Q{x//4+1}M{x%4+1}' for x in modules])]}\")\n", - " print(f\"Operating conditions are:\\n• Bias voltage: {bias_voltage}\\n• Memory cells: {max_cells}\\n\"\n", - " f\"• Acquisition rate: {acq_rate}\\n• Gain setting: {gain_setting}\\n• Photon Energy: {photon_energy}\\n\")\n", - " print(f\"Constant metadata is saved under \\\"retrieved-constants\\\" in calibration_metadata.yml\\n\")\n", + "if nodb_with_dark:\n", + " print(\"No constants were retrieved as calibrated files will be used.\")\n", "else:\n", - " print(\"No constants were retrieved as calibrated files will be used.\")" + " metadata.update({\"retrieved-constants\": mdata_dict})\n", + "\n", + " print(\"\\nRetrieved constants for modules:\",\n", + " ', '.join([tools.module_index_to_qm(x) for x in modules]))\n", + " print(f\"Operating conditions are:\")\n", + " print(f\"• Bias voltage: {bias_voltage}\")\n", + " print(f\"• Memory cells: {max_cells}\")\n", + " print(f\"• Acquisition rate: {acq_rate}\")\n", + " print(f\"• Gain mode: {gain_mode.name}\")\n", + " print(f\"• Gain setting: {gain_setting}\")\n", + " print(f\"• Photon Energy: {photon_energy}\")\n", + " print(\"Constant metadata is saved under \\\"retrieved-constants\\\" in calibration_metadata.yml\\n\")" ] }, { @@ -412,12 +395,11 @@ "metadata": {}, "outputs": [], "source": [ - "print(\"Constants are retrieved with creation time: \")\n", - "i = 0\n", + "print(\"Constants are retrieved with creation time:\")\n", "timestamps = {}\n", "\n", "for k_da, dinfo in mod_dev.items():\n", - " print(dinfo[\"mod\"], \":\")\n", + " print(f\"{dinfo['mod']}:\")\n", " module_timestamps = {}\n", " module_name = dinfo[\"mod\"]\n", " if k_da in mdata_dict:\n", @@ -425,10 +407,8 @@ " if hasattr(mdata[\"creation-time\"], 'strftime'):\n", " mdata[\"creation-time\"] = mdata[\"creation-time\"].strftime('%y-%m-%d %H:%M')\n", " print(f'{cname:.<12s}', mdata[\"creation-time\"])\n", - " # Store few time stamps if exists\n", - " # Add NA to keep array structure\n", " for cname in ['Offset', 'SlopesPC', 'SlopesFF']:\n", - " if not k_da in mdata_dict or dinfo[\"err\"]:\n", + " if k_da not in mdata_dict or dinfo[\"err\"]:\n", " module_timestamps[cname] = \"Err\"\n", " else:\n", " if cname in mdata_dict[k_da]:\n", @@ -440,13 +420,6 @@ " module_timestamps[cname] = \"NA\"\n", " timestamps[module_name] = module_timestamps\n", "\n", - " i += 1\n", - " if sequences:\n", - " seq_num = sequences[0]\n", - " else:\n", - " # if sequences[0] changed to None as it was -1\n", - " seq_num = 0\n", - "\n", "time_summary = metadata.setdefault(\"retrieved-constants\", {}).setdefault(\"time-summary\", {})\n", "time_summary[\"SAll\"] = timestamps\n", "\n", diff --git a/setup.py b/setup.py index 2e1236355e348595f1d11d4b53ef3909ea6db652..f4310aefeec6475f4a78c47f3bc2160a5cb5a87f 100644 --- a/setup.py +++ b/setup.py @@ -79,7 +79,7 @@ setup( }, ext_modules=extensions, install_requires=[ - "iCalibrationDB @ git+ssh://git@git.xfel.eu:10022/detectors/cal_db_interactive.git@2.0.1", # noqa + "iCalibrationDB @ git+ssh://git@git.xfel.eu:10022/detectors/cal_db_interactive.git@2.0.4", # noqa "nbparameterise @ git+ssh://git@git.xfel.eu:10022/detectors/nbparameterise.git@0.3", # noqa "XFELDetectorAnalysis @ git+ssh://git@git.xfel.eu:10022/karaboDevices/pyDetLib.git@2.5.6-2.10.0#subdirectory=lib", # noqa "Cython==0.29.21",