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

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 5c5e7fd 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
5import h5py
6import numpy as np
7import re
8import os
9import sys
10
11from sas.sascalc.dataloader.data_info import plottable_1D, plottable_2D,\
12    Data1D, Data2D, DataInfo, Process, Aperture, Collimation, \
13    TransmissionSpectrum, Detector
14from sas.sascalc.dataloader.data_info import combine_data_info_with_plottable
15from sas.sascalc.dataloader.loader_exceptions import FileContentsException, DefaultReaderException
16from sas.sascalc.dataloader.file_reader_base_class import FileReader
17
18def decode(s):
19    return s.decode() if isinstance(s, bytes) else s
20
21def h5attr(node, key, default=None):
22    return decode(node.attrs.get(key, default))
23
24class 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.