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

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.2.2ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change on this file since 9d786e5 was 9d786e5, checked in by lewis, 7 years ago

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