diff --git a/cal_tools/cal_tools/agipdlib.py b/cal_tools/cal_tools/agipdlib.py index cb4a1ca9891671d372494883b90252ec26492fdb..54afabdc2d78b32c0808dd4336817878fd838507 100644 --- a/cal_tools/cal_tools/agipdlib.py +++ b/cal_tools/cal_tools/agipdlib.py @@ -242,7 +242,8 @@ class AgipdCorrections: self.corr_bools = corr_bools else: bools = list(set(corr_bools) - set(tot_corr_bools)) - raise Exception(f'Correction Booleans: {bools} are not available!') # noqa + raise Exception('Correction Booleans: ' + f'{bools} are not available!') # Flags allowing for pulse capacitor constant retrieval. self.pc_bools = [self.corr_bools.get("rel_gain"), @@ -853,15 +854,15 @@ class AgipdCorrections: # save INDEX contents (first, count) in CORR files outfile.create_dataset(idx_base + "{}/first".format(do), - fidxv.shape, - dtype=fidxv.dtype, - data=fidxv, - fletcher32=True) + fidxv.shape, + dtype=fidxv.dtype, + data=fidxv, + fletcher32=True) outfile.create_dataset(idx_base + "{}/count".format(do), - cntsv.shape, - dtype=cntsv.dtype, - data=cntsv, - fletcher32=True) + cntsv.shape, + dtype=cntsv.dtype, + data=cntsv, + fletcher32=True) def retrieve_constant_and_time(self, const_dict, device, cal_db_interface, creation_time): @@ -899,7 +900,7 @@ class AgipdCorrections: print_once=0) return cons_data, when - def init_constants(self, cons_data, module_idx): + def init_constants(self, cons_data, when, module_idx): """ For CI derived gain, a mean multiplication factor of 4.48 compared to medium gain is used, as no reliable CI data for all memory cells @@ -938,6 +939,8 @@ class AgipdCorrections: rel_low gain = _rel_medium gain * 4.48 :param cons_data: A dictionary for each retrieved constant value. + :param when: A dictionary for the creation time + of each retrieved constant. :param module_idx: A module_idx index :return: """ @@ -954,114 +957,129 @@ class AgipdCorrections: bpixels = cons_data["BadPixelsDark"].astype(np.uint32) if self.corr_bools.get("xray_corr"): - bpixels |= cons_data["BadPixelsFF"].astype(np.uint32)[..., :bpixels.shape[2], None] # noqa - slopesFF = cons_data["SlopesFF"] - # This could be used for backward compatibility - # for very old SlopesFF constants - if len(slopesFF.shape) == 4: - slopesFF = slopesFF[..., 0] - # This is for backward compatability for old FF constants - # (128, 512, mem_cells) - if slopesFF.shape[-1] == 2: - xray_cor = np.squeeze(slopesFF[...,0]) - xray_cor_med = np.nanmedian(xray_cor) - xray_cor[np.isnan(xray_cor)]= xray_cor_med - xray_cor[(xray_cor<0.8) | (xray_cor>1.2)] = xray_cor_med - xray_cor = np.dstack([xray_cor]*self.max_cells) + if when["BadPixelsFF"]: + bpixels |= cons_data["BadPixelsFF"].astype(np.uint32)[..., + :bpixels.shape[2], # noqa + None] + + if when["SlopesFF"]: # Checking if constant was retrieved + + slopesFF = cons_data["SlopesFF"] + # This could be used for backward compatibility + # for very old SlopesFF constants + if len(slopesFF.shape) == 4: + slopesFF = slopesFF[..., 0] + # This is for backward compatability for old FF constants + # (128, 512, mem_cells) + if slopesFF.shape[-1] == 2: + xray_cor = np.squeeze(slopesFF[...,0]) + xray_cor_med = np.nanmedian(xray_cor) + xray_cor[np.isnan(xray_cor)]= xray_cor_med + xray_cor[(xray_cor<0.8) | (xray_cor>1.2)] = xray_cor_med + xray_cor = np.dstack([xray_cor]*self.max_cells) + else: + # Memory cell resolved xray_cor correction + xray_cor = slopesFF # (128, 512, mem_cells) + if xray_cor.shape[-1] < self.max_cells: + # In case of having new constant with less memory cells, + # due to lack of enough FF data or during development. + # xray_cor should be expanded by last memory cell. + xray_cor = np.dstack(xray_cor, + np.dstack([xray_cor[..., -1]] + * (self.max_cells - xray_cor.shape[-1]))) # noqa + # This is already done for old constants, + # but new constant is absolute and we need to have + # global ADU output for the moment + xray_cor /= self.ff_gain else: - # Memory cell resolved xray_cor correction - xray_cor = slopesFF # (128, 512, mem_cells) - if xray_cor.shape[-1] < self.max_cells: - # In case of having new constant with less memory cells, - # due to lack of enough FF data or during development. - # xray_cor should be expanded by last memory cell. - xray_cor = np.dstack(xray_cor, - np.dstack([xray_cor[..., -1]] - * (self.max_cells - xray_cor.shape[-1]))) # noqa - # This is already done for old constants, - # but new constant is absolute and we need to have - # global ADU output for the moment - xray_cor /= self.ff_gain + xray_cor = np.ones((128, 512, self.max_cells), np.float32) self.xray_cor[module_idx][...] = xray_cor.transpose()[...] # add additional bad pixel information if any(self.pc_bools): - bppc = np.moveaxis(cons_data["BadPixelsPC"].astype(np.uint32), - 0, 2) - bpixels |= bppc[..., :bpixels.shape[2], None] - - slopesPC = cons_data["SlopesPC"].astype(np.float32) - - # This will handle some historical data in a different format - # constant dimension injected first - if slopesPC.shape[0] == 10 or slopesPC.shape[0] == 11: - slopesPC = np.moveaxis(slopesPC, 0, 3) - slopesPC = np.moveaxis(slopesPC, 0, 2) + if when["BadPixelsPC"]: + bppc = np.moveaxis(cons_data["BadPixelsPC"].astype(np.uint32), + 0, 2) + bpixels |= bppc[..., :bpixels.shape[2], None] # calculate relative gain from the constants rel_gain = np.ones((128, 512, self.max_cells, 3), np.float32) - # high gain slope from pulse capacitor data - pc_high_m = slopesPC[..., :self.max_cells, 0] - pc_high_l = slopesPC[..., :self.max_cells, 1] - # medium gain slope from pulse capacitor data - pc_med_m = slopesPC[..., :self.max_cells, 3] - pc_med_l = slopesPC[..., :self.max_cells, 4] - - # calculate median for slopes - pc_high_med = np.nanmedian(pc_high_m, axis=(0,1)) - pc_med_med = np.nanmedian(pc_med_m, axis=(0,1)) - # calculate median for intercepts: - pc_high_l_med = np.nanmedian(pc_high_l, axis=(0,1)) - pc_med_l_med = np.nanmedian(pc_med_l, axis=(0,1)) - - # sanitize PC data - # (it should be done already on the level of constants) - # In the following loop, - # replace `nan`s across memory cells with - # the median value calculated previously. - # Then, values outside of the valid range (0.8 and 1.2) - # are fixed to the median value. - # This is applied for high and medium gain stages - for i in range(self.max_cells): - pc_high_m[np.isnan(pc_high_m[..., i])] = pc_high_med[i] - pc_med_m[np.isnan(pc_med_m[..., i])] = pc_med_med[i] - - pc_high_l[np.isnan(pc_high_l[..., i])] = pc_high_l_med[i] - pc_med_l[np.isnan(pc_med_l[..., i])] = pc_med_l_med[i] - - pc_high_m[(pc_high_m[..., i] < 0.8 * pc_high_med[i]) | - (pc_high_m[..., i] > 1.2 * pc_high_med[i])] = pc_high_med[i] # noqa - pc_med_m[(pc_med_m[..., i] < 0.8 * pc_med_med[i]) | - (pc_med_m[..., i] > 1.2 * pc_med_med[i])] = pc_med_med[i] # noqa - - pc_high_l[(pc_high_l[..., i] < 0.8 * pc_high_l_med[i]) | - (pc_high_l[..., i] > 1.2 * pc_high_l_med[i])] = pc_high_l_med[i] # noqa - pc_med_l[(pc_med_l[..., i] < 0.8 * pc_med_l_med[i]) | - (pc_med_l[..., i] > 1.2 * pc_med_l_med[i])] = pc_med_l_med[i] # noqa - - # ration between HG and MG per pixel per mem cell used - # for rel gain calculation - frac_high_med_pix = pc_high_m / pc_med_m - # avarage ration between HG and MG as a function of - # mem cell (needed for bls_stripes) - # TODO: Per pixel would be more optimal correction - frac_high_med = pc_high_med / pc_med_med - # calculate additional medium-gain offset - md_additional_offset = pc_high_l - pc_med_l * pc_high_m / pc_med_m - - # Calculate relative gain. If FF constants are available, - # use them for high gain - # if not rel_gain is calculated using PC data only - # if self.corr_bools.get("xray_corr"): - # rel_gain[..., :self.max_cells, 0] /= xray_corr - - # PC data should be 'calibrated with X-ray data, - # if it is not done, it is better to use 1 instead of bias - # the results with PC arteffacts. - # rel_gain[..., 0] = 1./(pc_high_m / pc_high_ave) - rel_gain[..., 1] = rel_gain[..., 0] * frac_high_med_pix - rel_gain[..., 2] = rel_gain[..., 1] * 4.48 + + if when["SlopesPC"]: + slopesPC = cons_data["SlopesPC"].astype(np.float32) + + # This will handle some historical data in a different format + # constant dimension injected first + if slopesPC.shape[0] == 10 or slopesPC.shape[0] == 11: + slopesPC = np.moveaxis(slopesPC, 0, 3) + slopesPC = np.moveaxis(slopesPC, 0, 2) + + # high gain slope from pulse capacitor data + pc_high_m = slopesPC[..., :self.max_cells, 0] + pc_high_l = slopesPC[..., :self.max_cells, 1] + # medium gain slope from pulse capacitor data + pc_med_m = slopesPC[..., :self.max_cells, 3] + pc_med_l = slopesPC[..., :self.max_cells, 4] + + # calculate median for slopes + pc_high_med = np.nanmedian(pc_high_m, axis=(0,1)) + pc_med_med = np.nanmedian(pc_med_m, axis=(0,1)) + # calculate median for intercepts: + pc_high_l_med = np.nanmedian(pc_high_l, axis=(0,1)) + pc_med_l_med = np.nanmedian(pc_med_l, axis=(0,1)) + + # sanitize PC data + # (it should be done already on the level of constants) + # In the following loop, + # replace `nan`s across memory cells with + # the median value calculated previously. + # Then, values outside of the valid range (0.8 and 1.2) + # are fixed to the median value. + # This is applied for high and medium gain stages + for i in range(self.max_cells): + pc_high_m[np.isnan(pc_high_m[..., i])] = pc_high_med[i] + pc_med_m[np.isnan(pc_med_m[..., i])] = pc_med_med[i] + + pc_high_l[np.isnan(pc_high_l[..., i])] = pc_high_l_med[i] + pc_med_l[np.isnan(pc_med_l[..., i])] = pc_med_l_med[i] + + pc_high_m[(pc_high_m[..., i] < 0.8 * pc_high_med[i]) | + (pc_high_m[..., i] > 1.2 * pc_high_med[i])] = pc_high_med[i] # noqa + pc_med_m[(pc_med_m[..., i] < 0.8 * pc_med_med[i]) | + (pc_med_m[..., i] > 1.2 * pc_med_med[i])] = pc_med_med[i] # noqa + + pc_high_l[(pc_high_l[..., i] < 0.8 * pc_high_l_med[i]) | + (pc_high_l[..., i] > 1.2 * pc_high_l_med[i])] = pc_high_l_med[i] # noqa + pc_med_l[(pc_med_l[..., i] < 0.8 * pc_med_l_med[i]) | + (pc_med_l[..., i] > 1.2 * pc_med_l_med[i])] = pc_med_l_med[i] # noqa + + # ration between HG and MG per pixel per mem cell used + # for rel gain calculation + frac_high_med_pix = pc_high_m / pc_med_m + # average ratio between HG and MG as a function of + # mem cell (needed for bls_stripes) + # TODO: Per pixel would be more optimal correction + frac_high_med = pc_high_med / pc_med_med + # calculate additional medium-gain offset + md_additional_offset = pc_high_l - pc_med_l * pc_high_m / pc_med_m # noqa + + # Calculate relative gain. If FF constants are available, + # use them for high gain + # if not rel_gain is calculated using PC data only + # if self.corr_bools.get("xray_corr"): + # rel_gain[..., :self.max_cells, 0] /= xray_corr + + # PC data should be 'calibrated with X-ray data, + # if it is not done, it is better to use 1 instead of bias + # the results with PC arteffacts. + # rel_gain[..., 0] = 1./(pc_high_m / pc_high_ave) + rel_gain[..., 1] = rel_gain[..., 0] * frac_high_med_pix + rel_gain[..., 2] = rel_gain[..., 1] * 4.48 + else: + # Intialize with fake calculated parameters of Ones + md_additional_offset = rel_gain + frac_high_med = np.ones((self.max_cells,), np.float32) self.md_additional_offset[module_idx][...] = md_additional_offset.transpose()[...] # noqa self.rel_gain[module_idx][...] = rel_gain[...].transpose() @@ -1099,7 +1117,8 @@ class AgipdCorrections: else: # Create empty constant using the list elements cons_data[cname] = getattr(np, mdata["file-path"][0])(mdata["file-path"][1]) # noqa - self.init_constants(cons_data, module_idx) + + self.init_constants(cons_data, when, module_idx) return when @@ -1174,7 +1193,7 @@ class AgipdCorrections: cal_db_interface, creation_time) - self.init_constants(cons_data, module_idx) + self.init_constants(cons_data, when, module_idx) return when @@ -1182,7 +1201,7 @@ class AgipdCorrections: """ Allocate memory for correction constants - :param modules: Module indises + :param modules: Module indices :param constant_shape: Shape of expected constants (gain, cells, x, y) """ for module_idx in modules: