cansas_reader_HDF5.py @ 25dd9c9

Last change on this file since 25dd9c9 was c94280c, checked in by krzywon, 8 years ago
Further down PEP 8.
Property mode set to `100644`
File size: 27.6 KB

Line
1	"""
2	CanSAS 2D data reader for reading HDF5 formatted CanSAS files.
3	"""
4
5	import h5py
6	import numpy as np
7	import re
8	import os
9	import sys
10
11	from sas.sascalc.dataloader.data_info import plottable_1D, plottable_2D,\
12	Data1D, Data2D, DataInfo, Process, Aperture, Collimation, \
13	TransmissionSpectrum, Detector
14	from sas.sascalc.dataloader.data_info import combine_data_info_with_plottable
15
16
17	class Reader():
18	"""
19	A class for reading in CanSAS v2.0 data files. The existing iteration opens
20	Mantid generated HDF5 formatted files with file extension .h5/.H5. Any
21	number of data sets may be present within the file and any dimensionality
22	of data may be used. Currently 1D and 2D SAS data sets are supported, but
23	future implementations will include 1D and 2D SESANS data.
24
25	Any number of SASdata sets may be present in a SASentry and the data within
26	can be either 1D I(Q) or 2D I(Qx, Qy).
27
28	Also supports reading NXcanSAS formatted HDF5 files
29
30	:Dependencies:
31	The CanSAS HDF5 reader requires h5py => v2.5.0 or later.
32	"""
33
34	# CanSAS version
35	cansas_version = 2.0
36	# Logged warnings or messages
37	logging = None
38	# List of errors for the current data set
39	errors = None
40	# Raw file contents to be processed
41	raw_data = None
42	# Data info currently being read in
43	current_datainfo = None
44	# SASdata set currently being read in
45	current_dataset = None
46	# List of plottable1D objects that should be linked to the current_datainfo
47	data1d = None
48	# List of plottable2D objects that should be linked to the current_datainfo
49	data2d = None
50	# Data type name
51	type_name = "CanSAS 2.0"
52	# Wildcards
53	type = ["CanSAS 2.0 HDF5 Files (.h5)\|.h5"]
54	# List of allowed extensions
55	ext = ['.h5', '.H5']
56	# Flag to bypass extension check
57	allow_all = True
58	# List of files to return
59	output = None
60
61	def read(self, filename):
62	"""
63	This is the general read method that all SasView data_loaders must have.
64
65	:param filename: A path for an HDF5 formatted CanSAS 2D data file.
66	:return: List of Data1D/2D objects and/or a list of errors.
67	"""
68	# Reinitialize when loading a new data file to reset all class variables
69	self.reset_class_variables()
70	# Check that the file exists
71	if os.path.isfile(filename):
72	basename = os.path.basename(filename)
73	_, extension = os.path.splitext(basename)
74	# If the file type is not allowed, return empty list
75	if extension in self.ext or self.allow_all:
76	# Load the data file
77	self.raw_data = h5py.File(filename, 'r')
78	# Read in all child elements of top level SASroot
79	self.read_children(self.raw_data, [])
80	# Add the last data set to the list of outputs
81	self.add_data_set()
82	# Close the data file
83	self.raw_data.close()
84	# Return data set(s)
85	return self.output
86
87	def reset_class_variables(self):
88	"""
89	Create the reader object and define initial states for class variables
90	"""
91	self.current_datainfo = None
92	self.current_dataset = None
93	self.data1d = []
94	self.data2d = []
95	self.raw_data = None
96	self.errors = set()
97	self.logging = []
98	self.output = []
99	self.parent_class = u''
100	self.detector = Detector()
101	self.collimation = Collimation()
102	self.aperture = Aperture()
103	self.process = Process()
104	self.trans_spectrum = TransmissionSpectrum()
105
106	def read_children(self, data, parent_list):
107	"""
108	A recursive method for stepping through the hierarchical data file.
109
110	:param data: h5py Group object of any kind
111	:param parent: h5py Group parent name
112	"""
113
114	# Loop through each element of the parent and process accordingly
115	for key in data.keys():
116	# Get all information for the current key
117	value = data.get(key)
118	if value.attrs.get(u'canSAS_class') is not None:
119	class_name = value.attrs.get(u'canSAS_class')
120	else:
121	class_name = value.attrs.get(u'NX_class')
122	if class_name is not None:
123	class_prog = re.compile(class_name)
124	else:
125	class_prog = re.compile(value.name)
126
127	if isinstance(value, h5py.Group):
128	self.parent_class = class_name
129	parent_list.append(key)
130	# If a new sasentry, store the current data sets and create
131	# a fresh Data1D/2D object
132	if class_prog.match(u'SASentry'):
133	self.add_data_set(key)
134	elif class_prog.match(u'SASdata'):
135	self._initialize_new_data_set(parent_list)
136	# Recursion step to access data within the group
137	self.read_children(value, parent_list)
138	self.add_intermediate()
139	parent_list.remove(key)
140
141	elif isinstance(value, h5py.Dataset):
142	# If this is a dataset, store the data appropriately
143	data_set = data[key][:]
144	unit = self._get_unit(value)
145
146	# I and Q Data
147	if key == u'I':
148	if isinstance(self.current_dataset, plottable_2D):
149	self.current_dataset.data = data_set
150	self.current_dataset.zaxis("Intensity", unit)
151	else:
152	self.current_dataset.y = data_set.flatten()
153	self.current_dataset.yaxis("Intensity", unit)
154	continue
155	elif key == u'Idev':
156	if isinstance(self.current_dataset, plottable_2D):
157	self.current_dataset.err_data = data_set.flatten()
158	else:
159	self.current_dataset.dy = data_set.flatten()
160	continue
161	elif key == u'Q':
162	self.current_dataset.xaxis("Q", unit)
163	if isinstance(self.current_dataset, plottable_2D):
164	self.current_dataset.q = data_set.flatten()
165	else:
166	self.current_dataset.x = data_set.flatten()
167	continue
168	elif key == u'Qdev':
169	self.current_dataset.dx = data_set.flatten()
170	continue
171	elif key == u'dQw':
172	self.current_dataset.dxw = data_set.flatten()
173	continue
174	elif key == u'dQl':
175	self.current_dataset.dxl = data_set.flatten()
176	continue
177	elif key == u'Qy':
178	self.current_dataset.yaxis("Q_y", unit)
179	self.current_dataset.qy_data = data_set.flatten()
180	continue
181	elif key == u'Qydev':
182	self.current_dataset.dqy_data = data_set.flatten()
183	continue
184	elif key == u'Qx':
185	self.current_dataset.xaxis("Q_x", unit)
186	self.current_dataset.qx_data = data_set.flatten()
187	continue
188	elif key == u'Qxdev':
189	self.current_dataset.dqx_data = data_set.flatten()
190	continue
191	elif key == u'Mask':
192	self.current_dataset.mask = data_set.flatten()
193	continue
194	# Transmission Spectrum
195	elif (key == u'T'
196	and self.parent_class == u'SAStransmission_spectrum'):
197	self.trans_spectrum.transmission = data_set.flatten()
198	continue
199	elif (key == u'Tdev'
200	and self.parent_class == u'SAStransmission_spectrum'):
201	self.trans_spectrum.transmission_deviation = \
202	data_set.flatten()
203	continue
204	elif (key == u'lambda'
205	and self.parent_class == u'SAStransmission_spectrum'):
206	self.trans_spectrum.wavelength = data_set.flatten()
207	continue
208
209	for data_point in data_set:
210	# Top Level Meta Data
211	if key == u'definition':
212	self.current_datainfo.meta_data['reader'] = data_point
213	elif key == u'run':
214	self.current_datainfo.run.append(data_point)
215	try:
216	run_name = value.attrs['name']
217	run_dict = {data_point: run_name}
218	self.current_datainfo.run_name = run_dict
219	except:
220	pass
221	elif key == u'title':
222	self.current_datainfo.title = data_point
223	elif key == u'SASnote':
224	self.current_datainfo.notes.append(data_point)
225
226	# Sample Information
227	# CanSAS 2.0 format
228	elif key == u'Title' and self.parent_class == u'SASsample':
229	self.current_datainfo.sample.name = data_point
230	# NXcanSAS format
231	elif key == u'name' and self.parent_class == u'SASsample':
232	self.current_datainfo.sample.name = data_point
233	# NXcanSAS format
234	elif key == u'ID' and self.parent_class == u'SASsample':
235	self.current_datainfo.sample.name = data_point
236	elif (key == u'thickness'
237	and self.parent_class == u'SASsample'):
238	self.current_datainfo.sample.thickness = data_point
239	elif (key == u'temperature'
240	and self.parent_class == u'SASsample'):
241	self.current_datainfo.sample.temperature = data_point
242	elif (key == u'transmission'
243	and self.parent_class == u'SASsample'):
244	self.current_datainfo.sample.transmission = data_point
245	elif (key == u'x_position'
246	and self.parent_class == u'SASsample'):
247	self.current_datainfo.sample.position.x = data_point
248	elif (key == u'y_position'
249	and self.parent_class == u'SASsample'):
250	self.current_datainfo.sample.position.y = data_point
251	elif key == u'pitch' and self.parent_class == u'SASsample':
252	self.current_datainfo.sample.orientation.x = data_point
253	elif key == u'yaw' and self.parent_class == u'SASsample':
254	self.current_datainfo.sample.orientation.y = data_point
255	elif key == u'roll' and self.parent_class == u'SASsample':
256	self.current_datainfo.sample.orientation.z = data_point
257	elif (key == u'details'
258	and self.parent_class == u'SASsample'):
259	self.current_datainfo.sample.details.append(data_point)
260
261	# Instrumental Information
262	elif (key == u'name'
263	and self.parent_class == u'SASinstrument'):
264	self.current_datainfo.instrument = data_point
265	elif key == u'name' and self.parent_class == u'SASdetector':
266	self.detector.name = data_point
267	elif key == u'SDD' and self.parent_class == u'SASdetector':
268	self.detector.distance = float(data_point)
269	self.detector.distance_unit = unit
270	elif (key == u'slit_length'
271	and self.parent_class == u'SASdetector'):
272	self.detector.slit_length = float(data_point)
273	self.detector.slit_length_unit = unit
274	elif (key == u'x_position'
275	and self.parent_class == u'SASdetector'):
276	self.detector.offset.x = float(data_point)
277	self.detector.offset_unit = unit
278	elif (key == u'y_position'
279	and self.parent_class == u'SASdetector'):
280	self.detector.offset.y = float(data_point)
281	self.detector.offset_unit = unit
282	elif (key == u'pitch'
283	and self.parent_class == u'SASdetector'):
284	self.detector.orientation.x = float(data_point)
285	self.detector.orientation_unit = unit
286	elif key == u'roll' and self.parent_class == u'SASdetector':
287	self.detector.orientation.z = float(data_point)
288	self.detector.orientation_unit = unit
289	elif key == u'yaw' and self.parent_class == u'SASdetector':
290	self.detector.orientation.y = float(data_point)
291	self.detector.orientation_unit = unit
292	elif (key == u'beam_center_x'
293	and self.parent_class == u'SASdetector'):
294	self.detector.beam_center.x = float(data_point)
295	self.detector.beam_center_unit = unit
296	elif (key == u'beam_center_y'
297	and self.parent_class == u'SASdetector'):
298	self.detector.beam_center.y = float(data_point)
299	self.detector.beam_center_unit = unit
300	elif (key == u'x_pixel_size'
301	and self.parent_class == u'SASdetector'):
302	self.detector.pixel_size.x = float(data_point)
303	self.detector.pixel_size_unit = unit
304	elif (key == u'y_pixel_size'
305	and self.parent_class == u'SASdetector'):
306	self.detector.pixel_size.y = float(data_point)
307	self.detector.pixel_size_unit = unit
308	elif (key == u'distance'
309	and self.parent_class == u'SAScollimation'):
310	self.collimation.length = data_point
311	self.collimation.length_unit = unit
312	elif (key == u'name'
313	and self.parent_class == u'SAScollimation'):
314	self.collimation.name = data_point
315	elif (key == u'shape'
316	and self.parent_class == u'SASaperture'):
317	self.aperture.shape = data_point
318	elif (key == u'x_gap'
319	and self.parent_class == u'SASaperture'):
320	self.aperture.size.x = data_point
321	elif (key == u'y_gap'
322	and self.parent_class == u'SASaperture'):
323	self.aperture.size.y = data_point
324
325	# Process Information
326	elif (key == u'Title'
327	and self.parent_class == u'SASprocess'): # CanSAS 2.0
328	self.process.name = data_point
329	elif (key == u'name'
330	and self.parent_class == u'SASprocess'): # NXcanSAS
331	self.process.name = data_point
332	elif (key == u'description'
333	and self.parent_class == u'SASprocess'):
334	self.process.description = data_point
335	elif key == u'date' and self.parent_class == u'SASprocess':
336	self.process.date = data_point
337	elif key == u'term' and self.parent_class == u'SASprocess':
338	self.process.term = data_point
339	elif self.parent_class == u'SASprocess':
340	self.process.notes.append(data_point)
341
342	# Source
343	elif (key == u'wavelength'
344	and self.parent_class == u'SASdata'):
345	self.current_datainfo.source.wavelength = data_point
346	self.current_datainfo.source.wavelength_unit = unit
347	elif (key == u'incident_wavelength'
348	and self.parent_class == 'SASsource'):
349	self.current_datainfo.source.wavelength = data_point
350	self.current_datainfo.source.wavelength_unit = unit
351	elif (key == u'wavelength_max'
352	and self.parent_class == u'SASsource'):
353	self.current_datainfo.source.wavelength_max = data_point
354	self.current_datainfo.source.wavelength_max_unit = unit
355	elif (key == u'wavelength_min'
356	and self.parent_class == u'SASsource'):
357	self.current_datainfo.source.wavelength_min = data_point
358	self.current_datainfo.source.wavelength_min_unit = unit
359	elif (key == u'incident_wavelength_spread'
360	and self.parent_class == u'SASsource'):
361	self.current_datainfo.source.wavelength_spread = \
362	data_point
363	self.current_datainfo.source.wavelength_spread_unit = \
364	unit
365	elif (key == u'beam_size_x'
366	and self.parent_class == u'SASsource'):
367	self.current_datainfo.source.beam_size.x = data_point
368	self.current_datainfo.source.beam_size_unit = unit
369	elif (key == u'beam_size_y'
370	and self.parent_class == u'SASsource'):
371	self.current_datainfo.source.beam_size.y = data_point
372	self.current_datainfo.source.beam_size_unit = unit
373	elif (key == u'beam_shape'
374	and self.parent_class == u'SASsource'):
375	self.current_datainfo.source.beam_shape = data_point
376	elif (key == u'radiation'
377	and self.parent_class == u'SASsource'):
378	self.current_datainfo.source.radiation = data_point
379	elif (key == u'transmission'
380	and self.parent_class == u'SASdata'):
381	self.current_datainfo.sample.transmission = data_point
382
383	# Everything else goes in meta_data
384	else:
385	new_key = self._create_unique_key(
386	self.current_datainfo.meta_data, key)
387	self.current_datainfo.meta_data[new_key] = data_point
388
389	else:
390	# I don't know if this reachable code
391	self.errors.add("ShouldNeverHappenException")
392
393	def add_intermediate(self):
394	"""
395	This method stores any intermediate objects within the final data set
396	after fully reading the set.
397
398	:param parent: The NXclass name for the h5py Group object that just
399	finished being processed
400	"""
401
402	if self.parent_class == u'SASprocess':
403	self.current_datainfo.process.append(self.process)
404	self.process = Process()
405	elif self.parent_class == u'SASdetector':
406	self.current_datainfo.detector.append(self.detector)
407	self.detector = Detector()
408	elif self.parent_class == u'SAStransmission_spectrum':
409	self.current_datainfo.trans_spectrum.append(self.trans_spectrum)
410	self.trans_spectrum = TransmissionSpectrum()
411	elif self.parent_class == u'SAScollimation':
412	self.current_datainfo.collimation.append(self.collimation)
413	self.collimation = Collimation()
414	elif self.parent_class == u'SASaperture':
415	self.collimation.aperture.append(self.aperture)
416	self.aperture = Aperture()
417	elif self.parent_class == u'SASdata':
418	if isinstance(self.current_dataset, plottable_2D):
419	self.data2d.append(self.current_dataset)
420	elif isinstance(self.current_dataset, plottable_1D):
421	self.data1d.append(self.current_dataset)
422
423	def final_data_cleanup(self):
424	"""
425	Does some final cleanup and formatting on self.current_datainfo and
426	all data1D and data2D objects and then combines the data and info into
427	Data1D and Data2D objects
428	"""
429
430	# Type cast data arrays to float64
431	if len(self.current_datainfo.trans_spectrum) > 0:
432	spectrum_list = []
433	for spectrum in self.current_datainfo.trans_spectrum:
434	spectrum.transmission = np.delete(spectrum.transmission, [0])
435	spectrum.transmission = spectrum.transmission.astype(np.float64)
436	spectrum.transmission_deviation = np.delete(
437	spectrum.transmission_deviation, [0])
438	spectrum.transmission_deviation = \
439	spectrum.transmission_deviation.astype(np.float64)
440	spectrum.wavelength = np.delete(spectrum.wavelength, [0])
441	spectrum.wavelength = spectrum.wavelength.astype(np.float64)
442	if len(spectrum.transmission) > 0:
443	spectrum_list.append(spectrum)
444	self.current_datainfo.trans_spectrum = spectrum_list
445
446	# Append errors to dataset and reset class errors
447	self.current_datainfo.errors = self.errors
448	self.errors.clear()
449
450	# Combine all plottables with datainfo and append each to output
451	# Type cast data arrays to float64 and find min/max as appropriate
452	for dataset in self.data2d:
453	dataset.data = dataset.data.astype(np.float64)
454	dataset.err_data = dataset.err_data.astype(np.float64)
455	if dataset.qx_data is not None:
456	dataset.xmin = np.min(dataset.qx_data)
457	dataset.xmax = np.max(dataset.qx_data)
458	dataset.qx_data = dataset.qx_data.astype(np.float64)
459	if dataset.dqx_data is not None:
460	dataset.dqx_data = dataset.dqx_data.astype(np.float64)
461	if dataset.qy_data is not None:
462	dataset.ymin = np.min(dataset.qy_data)
463	dataset.ymax = np.max(dataset.qy_data)
464	dataset.qy_data = dataset.qy_data.astype(np.float64)
465	if dataset.dqy_data is not None:
466	dataset.dqy_data = dataset.dqy_data.astype(np.float64)
467	if dataset.q_data is not None:
468	dataset.q_data = dataset.q_data.astype(np.float64)
469	zeros = np.ones(dataset.data.size, dtype=bool)
470	try:
471	for i in range(0, dataset.mask.size - 1):
472	zeros[i] = dataset.mask[i]
473	except:
474	self.errors.add(sys.exc_value)
475	dataset.mask = zeros
476	# Calculate the actual Q matrix
477	try:
478	if dataset.q_data.size <= 1:
479	dataset.q_data = np.sqrt(dataset.qx_data
480	* dataset.qx_data
481	+ dataset.qy_data
482	* dataset.qy_data)
483	except:
484	dataset.q_data = None
485
486	if dataset.data.ndim == 2:
487	(n_rows, n_cols) = dataset.data.shape
488	dataset.y_bins = dataset.qy_data[0::n_cols]
489	dataset.x_bins = dataset.qx_data[:n_cols]
490	dataset.data = dataset.data.flatten()
491
492	final_dataset = combine_data_info_with_plottable(
493	dataset, self.current_datainfo)
494	self.output.append(final_dataset)
495
496	for dataset in self.data1d:
497	if dataset.x is not None:
498	dataset.x = dataset.x.astype(np.float64)
499	dataset.xmin = np.min(dataset.x)
500	dataset.xmax = np.max(dataset.x)
501	if dataset.y is not None:
502	dataset.y = dataset.y.astype(np.float64)
503	dataset.ymin = np.min(dataset.y)
504	dataset.ymax = np.max(dataset.y)
505	if dataset.dx is not None:
506	dataset.dx = dataset.dx.astype(np.float64)
507	if dataset.dxl is not None:
508	dataset.dxl = dataset.dxl.astype(np.float64)
509	if dataset.dxw is not None:
510	dataset.dxw = dataset.dxw.astype(np.float64)
511	if dataset.dy is not None:
512	dataset.dy = dataset.dy.astype(np.float64)
513	final_dataset = combine_data_info_with_plottable(
514	dataset, self.current_datainfo)
515	self.output.append(final_dataset)
516
517	def add_data_set(self, key=""):
518	"""
519	Adds the current_dataset to the list of outputs after preforming final
520	processing on the data and then calls a private method to generate a
521	new data set.
522
523	:param key: NeXus group name for current tree level
524	"""
525
526	if self.current_datainfo and self.current_dataset:
527	self.final_data_cleanup()
528	self.data1d = []
529	self.data2d = []
530	self.current_datainfo = DataInfo()
531
532
533	def _initialize_new_data_set(self, parent_list=None):
534	"""
535	A private class method to generate a new 1D or 2D data object based on
536	the type of data within the set. Outside methods should call
537	add_data_set() to be sure any existing data is stored properly.
538
539	:param parent_list: List of names of parent elements
540	"""
541
542	if parent_list is None:
543	parent_list = []
544	if self._find_intermediate(parent_list, "Qx"):
545	self.current_dataset = plottable_2D()
546	else:
547	x = np.array(0)
548	y = np.array(0)
549	self.current_dataset = plottable_1D(x, y)
550	self.current_datainfo.filename = self.raw_data.filename
551
552	def _find_intermediate(self, parent_list, basename=""):
553	"""
554	A private class used to find an entry by either using a direct key or
555	knowing the approximate basename.
556
557	:param parent_list: List of parents nodes in the HDF5 file
558	:param basename: Approximate name of an entry to search for
559	:return:
560	"""
561
562	entry = False
563	key_prog = re.compile(basename)
564	top = self.raw_data
565	for parent in parent_list:
566	top = top.get(parent)
567	for key in top.keys():
568	if key_prog.match(key):
569	entry = True
570	break
571	return entry
572
573	def _create_unique_key(self, dictionary, name, numb=0):
574	"""
575	Create a unique key value for any dictionary to prevent overwriting
576	Recurses until a unique key value is found.
577
578	:param dictionary: A dictionary with any number of entries
579	:param name: The index of the item to be added to dictionary
580	:param numb: The number to be appended to the name, starts at 0
581	:return: The new name for the dictionary entry
582	"""
583	if dictionary.get(name) is not None:
584	numb += 1
585	name = name.split("_")[0]
586	name += "_{0}".format(numb)
587	name = self._create_unique_key(dictionary, name, numb)
588	return name
589
590	def _get_unit(self, value):
591	"""
592	Find the unit for a particular value within the h5py dictionary
593
594	:param value: attribute dictionary for a particular value set
595	:return: unit for the value passed to the method
596	"""
597	unit = value.attrs.get(u'units')
598	if unit is None:
599	unit = value.attrs.get(u'unit')
600	# Convert the unit formats
601	if unit == "1/A":
602	unit = "A^{-1}"
603	elif unit == "1/cm":
604	unit = "cm^{-1}"
605	return unit

Note: See TracBrowser for help on using the repository browser.

SasView

source: sasview/src/sas/sascalc/dataloader/readers/cansas_reader_HDF5.py @ 25dd9c9

Download in other formats: