Coverage for larch/xrd/xrd.py: 13%

256 statements  

« prev     ^ index     » next       coverage.py v7.6.0, created at 2024-10-16 21:04 +0000

1''' 

2This module defines a device-independent XRD class. 

3 

4Authors/Modifications: 

5---------------------- 

6* Margaret Koker, koker@cars.uchicago.edu 

7* modeled after MCA class 

8''' 

9 

10########################################################################## 

11# IMPORT PYTHON PACKAGES 

12 

13import os 

14import numpy as np 

15import larch 

16 

17from .xrd_tools import (d_from_q,d_from_twth,twth_from_d,twth_from_q, 

18 q_from_d,q_from_twth,E_from_lambda,lambda_from_E) 

19from .xrd_pyFAI import integrate_xrd,calc_cake 

20from .xrd_bgr import xrd_background 

21from .xrd_fitting import peakfinder,peaklocater,peakfilter,peakfitter 

22from larch.io import tifffile 

23from larch.utils import get_cwd 

24 

25########################################################################## 

26# CLASSES 

27 

28class xrd1d(larch.Group): 

29 ''' 

30 1D XRD data class 

31 

32 --> all length units in m (unless otherwise noted) 

33 

34 Attributes: 

35 ------------ 

36 * self.filename = 'CeO2_Allende.xy' # file containing x-y data 

37 * self.label = 'Data: CeO2_Allende' # name of data set 

38 

39 # Data parameters 

40 * self.I = None or array # intensity (units: counts) 

41 * self.q = None or array # inverse lattice spacing (units: 1/A) 

42 * self.d = None or array # lattice spacing (units: A) 

43 * self.twth = None or array # 2-theta (units: degrees) 

44 

45 # Calibration/experiemental parameters 

46 * self.wavelength = 0.688801071778 # incident x-ray wavelength (units: A) 

47 * self.energy = 18.000 # incident x-ray energy (units: keV) 

48 * self.distance = 0.49212785 # distance from sample to detector 

49 * self.poni = [0.2286, 0.2547] # center/point of normal incidence 

50 * self.rotation = [0.010, -0.00857, 0.0] # detector rotations (units: radians) 

51 * self.pixelsize = [0.0004, 0.0004] # size of pixels 

52 * self.splinefile = None # spline file for detector 

53 * self.polarization = None # polarization of detector 

54 * self.normalization = 1.0 # normalization factor for detector 

55 

56 # Data fitting parameters 

57 * self.uvw = [0.313, -0.109, 0.019] # instrumental broadening parameters 

58 * self.D = None # particle size broadening (units: A) 

59 * self.pki = [8, 254, 3664] # list of peak indecides 

60 

61 * self.imin = None # range for trimmed data 

62 * self.imax = None # range for trimmed data 

63 * self.bkgd = None # fit background for data 

64 

65 * self.matches = None # list of amcsd matches from database 

66 

67 * self.xrd2d = None # 2D data 

68 * self.cake = None # 2D cake 

69 

70 mkak 2017.03.15 

71 ''' 

72 def __init__(self, file=None, label=None, x=None, xtype=None, I=None, 

73 wavelength=None, energy=None): 

74 

75 self.filename = file 

76 self.label = label 

77 self.wavelength = wavelength 

78 self.energy = energy 

79 if energy is None and wavelength is None: 

80 self.energy = 19.0 

81 self.wavelength = lambda_from_E(self.energy) 

82 if self.energy is None: 

83 self.energy = E_from_lambda(self.wavelength) 

84 if self.wavelength is None: 

85 self.wavelength = lambda_from_E(self.energy) 

86 

87 if file is not None: 

88 self.xrd_from_file(file) 

89 else: 

90 ## Default values 

91 self.distance = None 

92 self.poni = None 

93 self.rotation = None 

94 self.pixelsize = None 

95 self.splinefile = None 

96 self.polarization = None 

97 self.normalization = None 

98 

99 if I is not None and x is not None: 

100 self.xrd_from_2d([x,I], xtype) 

101 self.bkgd = np.zeros(np.shape(self.I)) 

102 else: 

103 self.q = None 

104 self.twth = None 

105 self.d = None 

106 self.I = None 

107 self.bkgd = None 

108 

109 

110 ## Analysis parameters - set defaults 

111 self.uvw = None 

112 self.D = None 

113 self.pki = [] 

114 

115 self.imin = None 

116 self.imax = None 

117 

118 self.matches = None 

119 

120 self.xrd2d = None 

121 self.cake = None 

122 

123 larch.Group.__init__(self) 

124 

125 

126 def xrd_from_2d(self,xy, xtype, verbose=True): 

127 self.set_xy_data(xy, xtype) 

128 

129 def xrd_from_file(self, filename, verbose=False): 

130 

131 try: 

132 from ..xrmmap.asciifiles import read1DXRDFile 

133 head, dat = read1DXRDFile(filename) 

134 if verbose: 

135 print('Opening xrd data file: %s' % os.path.split(filename)[-1]) 

136 if len(head) < 4: 

137 print('WARNING: Using default energy for data. None given in file.') 

138 except: 

139 print('incorrect xy file format: %s' % os.path.split(filename)[-1]) 

140 return 

141 

142 if self.label is None: self.label = os.path.split(filename)[-1] 

143 

144 ## header info 

145 for line in head: 

146 import re 

147 line = re.sub(',','',line) 

148 

149 if 'SplineFile' in line: 

150 self.splinefile = line.split()[-1] 

151 if 'PixelSize' in line: 

152 self.pixelsize = [float(line.split()[-3]),float(line.split()[-2])] 

153 if 'PONI' in line: 

154 self.poni = [float(line.split()[2]),float(line.split()[3])] 

155 if 'Distance Sample to Detector' in line: 

156 self.distance = float(line.split()[-2]) 

157 if 'Rotations' in line: 

158 self.rotation = [float(line.split()[2]),float(line.split()[3]),float(line.split()[4])] 

159 

160 if 'Wavelength' in line: 

161 try: 

162 self.wavelength = float(line.split()[-1])*1e10 

163 except: 

164 self.wavelength = float(line.split()[-2])*1e10 

165 self.energy = E_from_lambda(self.wavelength) 

166 if 'Polarization' in line: 

167 if line.split()[-1] != 'None': self.polarization = float(line.split()[-1]) 

168 if 'Normalization' in line: 

169 try: 

170 value = float(line.split()[-1]) 

171 except: 

172 value = 1.0 

173 self.normalization = value 

174 

175 if 'q_' in line or '2th_' in line: 

176 xtype = line.split()[1] 

177 

178 ## data 

179 self.set_xy_data(dat, xtype) 

180 

181 def set_xy_data(self, xy, xtype): 

182 if xy is not None: 

183 xy = np.array(xy) 

184 if xy.shape[0] > xy.shape[1]: 

185 x,y = np.split(xy,2,axis=1) 

186 else: 

187 x,y = np.split(xy, 2, axis=0) 

188 self.q, self.twth, self.d = calculate_xvalues(x, xtype, self.wavelength) 

189 self.I = np.array(y).squeeze() 

190 

191 self.imin,self.imax = 0,len(self.q) 

192 self.bkgd = np.zeros(np.shape(self.I)) 

193 

194 def set_wavelength(self, wavelength): 

195 self.wavelength = wavelength 

196 self.energy = E_from_lambda(self.wavelength) 

197 self.q, self.twth, self.d = calculate_xvalues(self.q, 'q', self.wavelength) 

198 

199 def plot(self,reset=False,bkgd=False): 

200 

201 if reset: self.imin,self.imax = 0,len(self.I) 

202 if bkgd: 

203 return [self.q[self.imin:self.imax], 

204 self.twth[self.imin:self.imax], 

205 self.d[self.imin:self.imax], 

206 self.I[self.imin:self.imax]-self.bkgd, 

207 self.bkgd] 

208 return [self.q[self.imin:self.imax], 

209 self.twth[self.imin:self.imax], 

210 self.d[self.imin:self.imax], 

211 self.I[self.imin:self.imax], 

212 self.bkgd] 

213 

214 def reset_bkgd(self): 

215 self.bkgd = np.zeros(np.shape(self.I)) 

216 

217 def slct_xaxis(self, xtype='', xi=None): 

218 

219 if xtype.startswith('q') or xi == 0: 

220 x = self.q 

221 elif xtype.startswith('2th') or xi == 1: 

222 x = self.twth 

223 elif xtype.startswith('d') or xi == 2: 

224 x = self.d 

225 else: 

226 print('The provided x-axis label (%s or &i) not correct.' % (xtype,xi)) 

227 return 

228 

229 return x 

230 

231 def set_trim(self,xmin,xmax,xtype='',xi=None): 

232 

233 x = self.slct_xaxis(xtype=xtype,xi=xi) 

234 

235 self.imin,self.imax = 0,len(x)-1 

236 if xmin > np.min(x): 

237 self.imin = (np.abs(x-xmin)).argmin() 

238 if xmax < np.max(x): 

239 self.imax = (np.abs(x-xmax)).argmin() 

240 

241 def all_data(self,reset=False,bkgd=False): 

242 

243 if reset: self.imin,self.imax = 0,len(self.I) 

244 if len(self.I[self.imin:self.imax]) != len(self.bkgd): 

245 self.bkgd = np.zeros(len(self.I[self.imin:self.imax])) 

246 if bkgd: 

247 return [self.q[self.imin:self.imax], 

248 self.twth[self.imin:self.imax], 

249 self.d[self.imin:self.imax], 

250 self.I[self.imin:self.imax]-self.bkgd, 

251 self.bkgd] 

252 return [self.q[self.imin:self.imax], 

253 self.twth[self.imin:self.imax], 

254 self.d[self.imin:self.imax], 

255 self.I[self.imin:self.imax], 

256 self.bkgd] 

257 

258 def fit_background(self, **kwargs): 

259 x = self.q[self.imin:self.imax], 

260 y = self.I[self.imin:self.imax] 

261 bkgd = xrd_background(x, y, **kwargs) 

262 self.bkgd = np.zeros(len(y)) 

263 self.bkgd[:len(bkgd)] = bkgd 

264 

265 def find_peaks(self,bkgd=False,threshold=None,**kwargs): 

266 

267 all_data = np.array(self.all_data(bkgd=bkgd)) 

268 

269 

270 self.pki = peakfinder(all_data[3],**kwargs) 

271 if threshold is not None: 

272 self.pki = peakfilter(threshold,self.pki,all_data[3]) 

273 

274 pk_data = np.zeros((5,len(self.pki))) 

275 for i,pki in enumerate(self.pki): pk_data[:,i] = all_data[:,pki] 

276 

277# def refine_peaks(self,trim=False,bkgd=False): 

278# 

279# q,twth,d,I = self.trim_all(trim,bkgd) 

280# 

281# pktwth,pkfwhm,self.Ipks = peakfitter(self.pki,twth,I,fittype='double') 

282# #self.peaks = zip(pkfwhm,pkI) 

283# 

284# self.qpks,self.twthpks,self.dpks = calculate_xvalues(pktwth,'2th',self.wavelength) 

285 

286# def fit_pattern(self): 

287# ## This isn't written yet 

288# ## mkak 2017.04.03 

289# fit = np.zeros(len(self.I)) 

290# for i,j in enumerate(self.pki): 

291# a = None 

292 

293 

294 

295def read_xrd_data(filepath): 

296 

297 if not os.path.exists(filepath): 

298 return 

299 

300 try: 

301 data = np.array(tifffile.imread(filepath)) 

302 except: # TypeError: 

303 try: 

304 from ..xrmmap import read_xrd_netcdf 

305 data = np.array(read_xrd_netcdf(filepath)) 

306 except: 

307 try: 

308 data = xrd1d(file=filepath).I 

309 except: 

310 return 

311 return data 

312 

313 

314class XRD(larch.Group): 

315 ''' 

316 X-Ray Diffraction (XRD) class 

317 

318 Attributes: 

319 ----------- 

320 * self.name = 'xrd' # Name of the object 

321 * self.xpix = 2048 # Number of x pixels 

322 * self.ypix = 2048 # Number of y pixels 

323 * self.data2D = None # 2D XRD data 

324 * self.data1D = None # 1D XRD data 

325 

326 Notes: 

327 ------ 

328 

329 mkak 2016.08.20 

330 ''' 

331 

332 def __init__(self, data2D=None, xpixels=2048, ypixels=2048, 

333 data1D=None, nwedge=0, title=None, 

334 steps=5001, name='xrd', filename=None, 

335 calfile=None, energy=None, wavelength=None, 

336 npixels=None, _larch=None, **kws): 

337 

338 self.name = name 

339 self.xpix = xpixels 

340 self.ypix = ypixels 

341 self.data2D = data2D 

342 self.nwedge = nwedge 

343 self.steps = steps 

344 self.data1D = data1D 

345 self.data2D = data2D 

346 self.cake = None 

347 

348 self.calfile = calfile 

349 

350 if energy is None and wavelength is not None: 

351 self.wavelegth = wavelength 

352 self.energy = E_from_lambda(wavelength) 

353 elif energy is not None and wavelength is None: 

354 self.energy = energy 

355 self.wavelength = lambda_from_E(self.energy) 

356 else: 

357 self.energy = energy 

358 self.wavelength = wavelength 

359 

360 self.filename = filename 

361 self.title = title 

362 self.npixels = npixels 

363 

364 larch.Group.__init__(self) 

365 

366 def __repr__(self): 

367 if self.data2D is not None: 

368 form = "<2DXRD %s: pixels = %d, %d>" 

369 return form % (self.name, self.xpix, self.ypix) 

370 elif self.data1D is not None: 

371 form = "<1DXRD %s: channels = %d>" 

372 return form % (self.name, self.steps) 

373 else: 

374 form = "<no 1D or 2D XRD pattern given>" 

375 return form 

376 

377 def add_environ(self, desc='', val='', addr=''): 

378 '''add an Environment setting''' 

379 if len(desc) > 0 and len(val) > 0: 

380 self.environ.append(Environment(desc=desc, val=val, addr=addr)) 

381 

382 

383 def calc_1D(self,save=False,calccake=True,calc1d=True,verbose=False): 

384 

385 kwargs = {'steps':self.steps} 

386 

387 if save: 

388 file = self.save_1D() 

389 kwargs.update({'file':file}) 

390 

391 if os.path.exists(self.calfile): 

392 if len(self.data2D) > 0: 

393 self.data1D = integrate_xrd(self.data2D,self.calfile, 

394 verbose=verbose,**kwargs) 

395 self.cake = calc_cake(self.data2D,self.calfile, unit='q') 

396 else: 

397 if verbose: 

398 print('No 2D XRD data provided.') 

399 else: 

400 if verbose: 

401 print('Could not locate file %s' % self.calfile) 

402 

403 def save_1D(self,file=None): 

404 

405 if file is None: 

406 counter = 1 

407 while os.path.exists('%s/%s_%03d.xy' % (get_cwd(),self.title,counter)): 

408 counter += 1 

409 file = '%s/%s_%03d.xy' % (get_cwd(),self.title,counter) 

410 

411 return file 

412 

413 

414 def save_2D(self,file=None,verbose=False): 

415 

416 if file is None: 

417 counter = 1 

418 while os.path.exists('%s/%s_%03d.tiff' % (get_cwd(),self.title,counter)): 

419 counter += 1 

420 file = '%s/%s_%03d.tiff' % (get_cwd(),self.title,counter) 

421 

422 tifffile.imsave(file,self.data2D) 

423 

424 

425 

426########################################################################## 

427# FUNCTIONS 

428 

429def calculate_xvalues(x, xtype, wavelength): 

430 ''' 

431 projects given x-axis onto q-, 2theta-, and d-axes 

432 

433 x : list or array (expected units: 1/A, deg, or A) 

434 xtype : options 'q', '2th', or 'd' 

435 wavelength : incident x-ray wavelength (units: A) 

436 

437 q, twth, d : returned with same dimensions as x (units: 1/A, deg, A) 

438 ''' 

439 

440 x = np.array(x).squeeze() 

441 if xtype.startswith('q'): 

442 q = x 

443 d = d_from_q(q) 

444 if wavelength is not None: 

445 twth = twth_from_q(q,wavelength) 

446 else: 

447 twth = np.zeros(len(q)) 

448 

449 elif xtype.startswith('2th'): 

450 twth = x 

451 if wavelength is not None: 

452 q = q_from_twth(twth,wavelength) 

453 d = d_from_twth(twth,wavelength) 

454 else: 

455 q = np.zeros(len(twth)) 

456 d = np.zeros(len(twth)) 

457 

458 elif xtype.startswith('d'): 

459 

460 d = x 

461 q = q_from_d(d) 

462 if wavelength is not None: 

463 twth = twth_from_d(d,wavelength) 

464 else: 

465 twth = np.zeros(len(d)) 

466 

467 else: 

468 print('The provided x-axis label (%s) not correct. Check data.' % xtype) 

469 return None,None,None 

470 

471 

472 return q,twth,d 

473 

474 

475def create_xrd(data2D=None, xpixels=2048, ypixels=2048, 

476 data1D=None, nwedge=0, steps=5001, 

477 name='xrd', _larch=None, **kws): 

478 

479 ''' 

480 create an XRD object 

481 

482 Parameters: 

483 ------------ 

484 data2D: 2D diffraction patterns 

485 data1D: 1D diffraction patterns 

486 xpixels: number of x pixels 

487 ypixels: number of y pixels 

488 

489 Returns: 

490 ---------- 

491 an XRD object 

492 

493 ''' 

494 return XRD(data2D=data2D, data1D=data1D, xpixels=xpixels, ypixels=ypixels, 

495 name=name, **kws) 

496 

497 

498def create_xrd1d(file, _larch=None, **kws): 

499 

500 ''' 

501 create an XRD object 

502 

503 Parameters: 

504 ------------ 

505 data2D: 2D diffraction patterns 

506 data1D: 1D diffraction patterns 

507 xpixels: number of x pixels 

508 ypixels: number of y pixels 

509 

510 Returns: 

511 ---------- 

512 an XRD object 

513 

514 ''' 

515 return xrd1d(file=file, **kws)