diff --git a/src/toolbox_scs/detectors/viking.py b/src/toolbox_scs/detectors/viking.py index 6c6c50aa4834e85622d9e614a4fda497232e0b46..dc795e50a818ae75aeb450aa203e14b6f15fd0cb 100644 --- a/src/toolbox_scs/detectors/viking.py +++ b/src/toolbox_scs/detectors/viking.py @@ -18,13 +18,14 @@ class Viking: # dimension for integration INTEGRATE_DIM = 'newt_y' USE_DARK = False - + dark_image = None + # polynomial degree for background subtraction POLY_DEG = 5 FIELDS = ['newton'] - ENERGY_CALIB = [9.8233e-7, 0.0240, 514.4795] + ENERGY_CALIB = [0, 1, 0] def set_params(self, **params): for key, value in params.items(): @@ -53,13 +54,11 @@ class Viking: data.attrs[k] = v return data - def load_dark(self, runNB=None, proposal=None): + def load_dark(self, runNB=None): if runNB is None: self.USE_DARK = False return - if proposal is None: - proposal = self.PROPOSAL - data = self.from_run(runNB, proposal, add_attrs=False) + data = self.from_run(runNB, add_attrs=False) self.dark_image = data['newton'].mean(dim='trainId') self.dark_image.attrs['runNB'] = runNB self.USE_DARK = True @@ -67,8 +66,8 @@ class Viking: def integrate(self, data): imgs = data['newton'] if self.USE_DARK: - imgs -= self.dark_image - data['spectrum'] = imgs.sum(dim=self.INTEGRATE_DIM) + imgs = imgs - self.dark_image + data['spectrum'] = imgs.mean(dim=self.INTEGRATE_DIM) return data def get_camera_gain(self, run): @@ -94,8 +93,9 @@ class Viking: """ Conversion factor from camera digital counts to photoelectrons per count. The values can be found in the camera datasheet - but they have been slightly corrected for High Sensitivity - mode after analysis of runs 1204, 1207 and 1208, proposal 2937. + (Andor Newton) but they have been slightly corrected for High + Sensitivity mode after analysis of runs 1204, 1207 and 1208, + proposal 2937. Parameters: ------ @@ -127,7 +127,7 @@ class Viking: 'startY': 'imageSpecifications.startY.value', 'endY': 'imageSpecifications.endY.value', 'temperature': 'CoolerActual.temperature.value', - 'high_sensitivity': 'HighCapacity.value', + 'high_capacity': 'HighCapacity.value', 'exposure_s': 'exposureTime.value' } ret = {} @@ -175,10 +175,10 @@ class Viking: bl = spectra[:, mask] fit = np.polyfit(x_bl, bl.T, self.POLY_DEG) if len(spectra.shape) == 1: - return spectra - np.poly1d(fit)(x) + return spectra - np.polyval(fit, x) final_bl = np.empty(spectra.shape) for t in range(spectra.shape[0]): - final_bl[t] = np.poly1d(fit[:, t])(x) + final_bl[t] = np.polyval(fit[:, t], x) data['spectrum_nobg'] = spectra - final_bl return spectra - final_bl @@ -197,23 +197,27 @@ class Viking: ds = xr.Dataset() ds['It'] = spectrum ds['It_std'] = std - ds['It_std_err'] = std_err + ds['It_stderr'] = std_err ds['I0'] = spectrum_ref ds['I0_std'] = std - ds['I0_std_err'] = std_err - ds['absorption'] = spectrum_ref / spectrum - ds['absorption_std'] = np.abs(ds['absorption']) * np.sqrt( + ds['I0_stderr'] = std_err + absorption = spectrum_ref / spectrum + # assume that It and I0 are independent variables: + absorption_std = np.abs(absorption) * np.sqrt( std_ref**2 / spectrum_ref**2 + std**2 / spectrum**2) - ds['absorption_stderr'] = np.abs(ds['absorption']) * np.sqrt( + absorption_stderr = np.abs(absorption) * np.sqrt( (std_err_ref / spectrum_ref)**2 + (std_err / spectrum)**2) - ds['absorptionCoef'] = np.log(ds['absorption']) / thickness - ds['absorptionCoef_std'] = ds['absorption_std'] / (thickness * np.abs(ds['absorption'])) - ds['absorptionCoef_stderr'] = ds['absorption_stderr'] / (thickness * np.abs(ds['absorption'])) + ds['absorptionCoef'] = np.log(absorption) / thickness + ds['absorptionCoef_std'] = absorption_std / (thickness * np.abs(absorption)) + ds['absorptionCoef_stderr'] = absorption_stderr / (thickness * np.abs(absorption)) + ds.attrs['n_It'] = data[key].sizes['trainId'] + ds.attrs['n_I0'] = data_ref[key].sizes['trainId'] if plot: plot_viking_xas(ds, plot_errors, xas_ylim) return ds + def plot_viking_xas(xas, plot_errors=True, xas_ylim=(-1, 3)): fig, ax = plt.subplots(3, 1, figsize=(7,7), sharex=True) ax[0].plot(xas.newt_x, xas['I0']) @@ -231,12 +235,12 @@ def plot_viking_xas(xas, plot_errors=True, xas_ylim=(-1, 3)): if plot_errors: ax[0].fill_between(xas.newt_x, - xas['I0'] - 1.96*xas['I0_std_err'], - xas['I0'] + 1.96*xas['I0_std_err'], + xas['I0'] - 1.96*xas['I0_stderr'], + xas['I0'] + 1.96*xas['I0_stderr'], alpha=0.2) ax[1].fill_between(xas.newt_x, - xas['It'] - 1.96*xas['It_std_err'], - xas['It'] + 1.96*xas['It_std_err'], + xas['It'] - 1.96*xas['It_stderr'], + xas['It'] + 1.96*xas['It_stderr'], alpha=0.2) ax[2].fill_between(xas.newt_x, xas['absorptionCoef'] - 1.96*xas['absorptionCoef_stderr'],