cansas_reader_HDF5.py @ f7d720f

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalcmagnetic_scattrelease-4.2.2ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload

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

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

SasView

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

Download in other formats: