Mercurial > repos > kls286 > chap_test_20230328

comparison build/bdist.linux-x86_64/egg/CHAP/models/edd.py @ 0:cbbe42422d56 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/CHESSComputing/ChessAnalysisPipeline/tree/galaxy commit 1401a7e1ae007a6bda260d147f9b879e789b73e0-dirty
author kls286
date Tue, 28 Mar 2023 15:07:30 +0000
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:cbbe42422d56
1 from msnctools.general import create_mask
2 from msnctools.material import Material
3 from msnctools.scanparsers import SMBMCAScanParser as ScanParser
4 import numpy as np
5 from pathlib import PosixPath
6 from pydantic import (BaseModel,
7 confloat,
8 conint,
9 conlist,
10 constr,
11 FilePath,
12 validator)
13 from scipy.interpolate import interp1d
14 from typing import Optional
15
16
17 class MCACeriaCalibrationConfig(BaseModel):
18 '''Class representing metadata required to perform a Ceria calibration for an
19 MCA detector.
20
21 :ivar spec_file: Path to the SPEC file containing the CeO2 scan
22 :ivar scan_number: Number of the CeO2 scan in `spec_file`
23 :ivar scan_step_index: Index of the scan step to use for calibration,
24 optional. If not specified, the calibration routine will be performed on
25 the average of all MCA spectra for the scan.
26
27 :ivar flux_file: csv file containing station beam energy in eV (column 0)
28 and flux (column 1)
29
30 :ivar detector_name: name of the MCA to calibrate
31 :ivar num_bins: number of channels on the MCA to calibrate
32 :ivar max_energy_kev: maximum channel energy of the MCA in keV
33
34 :ivar hexrd_h5_material_file: path to a HEXRD materials.h5 file containing an
35 entry for the material properties.
36 :ivar hexrd_h5_material_name: Name of the material entry in
37 `hexrd_h5_material_file`, defaults to `'CeO2'`.
38 :ivar lattice_parameter_angstrom: lattice spacing in angstrom to use for
39 the cubic CeO2 crystal, defaults to `5.41153`.
40
41 :ivar tth_max: detector rotation about hutch x axis, defaults to `90`.
42 :ivar hkl_tth_tol: minimum resolvable difference in 2&theta between two
43 unique HKL peaks, defaults to `0.15`.
44
45 :ivar fit_include_bin_ranges: list of MCA channel index ranges whose data
46 will be included in the calibration routine
47 :ivar fit_hkls: list of unique HKL indices to fit peaks for in the
48 calibration routine
49
50 :ivar tth_initial_guess: initial guess for 2&theta
51 :ivar slope_initial_guess: initial guess for detector channel energy
52 correction linear slope, defaults to `1.0`.
53 :ivar intercept_initial_guess: initial guess for detector channel energy
54 correction y-intercept, defaults to `0.0`.
55
56 :ivar tth_calibrated: calibrated value for 2&theta, defaults to None
57 :ivar slope_calibrated: calibrated value for detector channel energy
58 correction linear slope, defaults to `None`
59 :ivar intercept_calibrated: calibrated value for detector channel energy
60 correction y-intercept, defaluts to None
61
62 :ivar max_iter: maximum number of iterations of the calibration routine,
63 defaults to `10`.
64 :ivar tune_tth_tol: stop iteratively tuning 2&theta when an iteration
65 produces a change in the tuned value of 2&theta that is smaller than this
66 value, defaults to `1e-8`.
67 '''
68
69 spec_file: FilePath
70 scan_number: conint(gt=0)
71 scan_step_index: Optional[conint(ge=0)]
72
73 flux_file: FilePath
74
75 detector_name: constr(strip_whitespace=True, min_length=1)
76 num_bins: conint(gt=0)
77 max_energy_kev: confloat(gt=0)
78
79 hexrd_h5_material_file: FilePath
80 hexrd_h5_material_name: constr(strip_whitespace=True, min_length=1) = 'CeO2'
81 lattice_parameter_angstrom: confloat(gt=0) = 5.41153
82
83 tth_max: confloat(gt=0, allow_inf_nan=False) = 90.0
84 hkl_tth_tol: confloat(gt=0, allow_inf_nan=False) = 0.15
85
86 fit_include_bin_ranges: conlist(min_items=1,
87 item_type=conlist(item_type=conint(ge=0),
88 min_items=2,
89 max_items=2))
90 fit_hkls: conlist(item_type=conint(ge=0), min_items=1)
91
92 tth_initial_guess: confloat(gt=0, le=tth_max, allow_inf_nan=False)
93 slope_initial_guess: float = 1.0
94 intercept_initial_guess: float = 0.0
95 tth_calibrated: Optional[confloat(gt=0, allow_inf_nan=False)]
96 slope_calibrated: Optional[confloat(allow_inf_nan=False)]
97 intercept_calibrated: Optional[confloat(allow_inf_nan=False)]
98
99 max_iter: conint(gt=0) = 10
100 tune_tth_tol: confloat(ge=0) = 1e-8
101
102 @validator('fit_include_bin_ranges', each_item=True)
103 def validate_include_bin_range(cls, value, values):
104 '''Ensure no bin ranges are outside the boundary of the detector'''
105
106 num_bins = values.get('num_bins')
107 value[1] = min(value[1], num_bins)
108 return(value)
109
110 def mca_data(self):
111 '''Get the 1D array of MCA data to use for calibration.
112
113 :return: MCA data
114 :rtype: np.ndarray
115 '''
116
117 scanparser = ScanParser(self.spec_file, self.scan_number)
118 if self.scan_step_index is None:
119 data = scanparser.get_all_detector_data(self.detector_name)
120 if scanparser.spec_scan_npts > 1:
121 data = np.average(data, axis=1)
122 else:
123 data = data[0]
124 else:
125 data = scanparser.get_detector_data(self.detector_name, self.scan_step_index)
126
127 return(np.array(data))
128
129 def mca_mask(self):
130 '''Get a boolean mask array to use on MCA data before fitting.
131
132 :return: boolean mask array
133 :rtype: numpy.ndarray
134 '''
135
136 mask = None
137 bin_indices = np.arange(self.num_bins)
138 for bin_range in self.fit_include_bin_ranges:
139 mask = create_mask(bin_indices,
140 bounds=bin_range,
141 exclude_bounds=False,
142 current_mask=mask)
143
144 return(mask)
145
146 def flux_correction_interpolation_function(self):
147 '''Get an interpolation function to correct MCA data for relative energy
148 flux of the incident beam.
149
150 :return: energy flux correction interpolation function
151 :rtype: scipy.interpolate._polyint._Interpolator1D
152 '''
153
154 flux = np.loadtxt(self.flux_file)
155 energies = flux[:,0]/1.e3
156 relative_intensities = flux[:,1]/np.max(flux[:,1])
157 interpolation_function = interp1d(energies, relative_intensities)
158 return(interpolation_function)
159
160 def material(self):
161 '''Get CeO2 as a `msnctools.materials.Material` object.
162
163 :return: CeO2 material
164 :rtype: msnctools.material.Material
165 '''
166
167 material = Material(material_name=self.hexrd_h5_material_name,
168 material_file=self.hexrd_h5_material_file,
169 lattice_parameters_angstroms=self.lattice_parameter_angstrom)
170 # The following kwargs will be needed if we allow the material to be
171 # built using xrayutilities (for now, we only allow hexrd to make the
172 # material):
173 # sgnum=225,
174 # atoms=['Ce4p', 'O2mdot'],
175 # pos=[(0.,0.,0.), (0.25,0.75,0.75)],
176 # enrgy=50000.) # Why do we need to specify an energy to get HKLs when using xrayutilities?
177 return(material)
178
179 def unique_ds(self):
180 '''Get a list of unique HKLs and their lattice spacings
181
182 :return: unique HKLs and their lattice spacings in angstroms
183 :rtype: np.ndarray, np.ndarray
184 '''
185
186 unique_hkls, unique_ds = self.material().get_unique_ds(tth_tol=self.hkl_tth_tol, tth_max=self.tth_max)
187
188 return(unique_hkls, unique_ds)
189
190 def fit_ds(self):
191 '''Get a list of HKLs and their lattice spacings that will be fit in the
192 calibration routine
193
194 :return: HKLs to fit and their lattice spacings in angstroms
195 :rtype: np.ndarray, np.ndarray
196 '''
197
198 unique_hkls, unique_ds = self.unique_ds()
199
200 fit_hkls = np.array([unique_hkls[i] for i in self.fit_hkls])
201 fit_ds = np.array([unique_ds[i] for i in self.fit_hkls])
202
203 return(fit_hkls, fit_ds)
204
205 def dict(self):
206 '''Return a representation of this configuration in a dictionary that is
207 suitable for dumping to a YAML file (one that converts all instances of
208 fields with type `PosixPath` to `str`).
209
210 :return: dictionary representation of the configuration.
211 :rtype: dict
212 '''
213
214 d = super().dict()
215 for k,v in d.items():
216 if isinstance(v, PosixPath):
217 d[k] = str(v)
218 return(d)