source: sasview/src/sas/sascalc/dataloader/data_info.py @ a24eacf

ESS_GUIESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalc
Last change on this file since a24eacf was b8080e1, checked in by Piotr Rozyczko <rozyczko@…>, 6 years ago

cherry picking sascalc changes from master SASVIEW-996
minor unit test fixes

  • Property mode set to 100644
File size: 39.8 KB
Line 
1"""
2    Module that contains classes to hold information read from
3    reduced data files.
4
5    A good description of the data members can be found in
6    the CanSAS 1D XML data format:
7
8    http://www.smallangles.net/wgwiki/index.php/cansas1d_documentation
9"""
10#####################################################################
11#This software was developed by the University of Tennessee as part of the
12#Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
13#project funded by the US National Science Foundation.
14#See the license text in license.txt
15#copyright 2008, University of Tennessee
16######################################################################
17
18from __future__ import print_function
19
20#TODO: Keep track of data manipulation in the 'process' data structure.
21#TODO: This module should be independent of plottables. We should write
22#        an adapter class for plottables when needed.
23
24#from sas.guitools.plottables import Data1D as plottable_1D
25from sas.sascalc.data_util.uncertainty import Uncertainty
26import numpy as np
27import math
28
29class plottable_1D(object):
30    """
31    Data1D is a place holder for 1D plottables.
32    """
33    # The presence of these should be mutually
34    # exclusive with the presence of Qdev (dx)
35    x = None
36    y = None
37    dx = None
38    dy = None
39    ## Slit smearing length
40    dxl = None
41    ## Slit smearing width
42    dxw = None
43    ## SESANS specific params (wavelengths for spin echo length calculation)
44    lam = None
45    dlam = None
46
47    # Units
48    _xaxis = ''
49    _xunit = ''
50    _yaxis = ''
51    _yunit = ''
52
53    def __init__(self, x, y, dx=None, dy=None, dxl=None, dxw=None, lam=None, dlam=None):
54        self.x = np.asarray(x)
55        self.y = np.asarray(y)
56        if dx is not None:
57            self.dx = np.asarray(dx)
58        if dy is not None:
59            self.dy = np.asarray(dy)
60        if dxl is not None:
61            self.dxl = np.asarray(dxl)
62        if dxw is not None:
63            self.dxw = np.asarray(dxw)
64        if lam is not None:
65            self.lam = np.asarray(lam)
66        if dlam is not None:
67            self.dlam = np.asarray(dlam)
68
69    def xaxis(self, label, unit):
70        """
71        set the x axis label and unit
72        """
73        self._xaxis = label
74        self._xunit = unit
75
76    def yaxis(self, label, unit):
77        """
78        set the y axis label and unit
79        """
80        self._yaxis = label
81        self._yunit = unit
82
83
84class plottable_2D(object):
85    """
86    Data2D is a place holder for 2D plottables.
87    """
88    xmin = None
89    xmax = None
90    ymin = None
91    ymax = None
92    data = None
93    qx_data = None
94    qy_data = None
95    q_data = None
96    err_data = None
97    dqx_data = None
98    dqy_data = None
99    mask = None
100
101    # Units
102    _xaxis = ''
103    _xunit = ''
104    _yaxis = ''
105    _yunit = ''
106    _zaxis = ''
107    _zunit = ''
108
109    def __init__(self, data=None, err_data=None, qx_data=None,
110                 qy_data=None, q_data=None, mask=None,
111                 dqx_data=None, dqy_data=None):
112        self.data = np.asarray(data)
113        self.qx_data = np.asarray(qx_data)
114        self.qy_data = np.asarray(qy_data)
115        self.q_data = np.asarray(q_data)
116        self.mask = np.asarray(mask)
117        self.err_data = np.asarray(err_data)
118        if dqx_data is not None:
119            self.dqx_data = np.asarray(dqx_data)
120        if dqy_data is not None:
121            self.dqy_data = np.asarray(dqy_data)
122
123    def xaxis(self, label, unit):
124        """
125        set the x axis label and unit
126        """
127        self._xaxis = label
128        self._xunit = unit
129
130    def yaxis(self, label, unit):
131        """
132        set the y axis label and unit
133        """
134        self._yaxis = label
135        self._yunit = unit
136
137    def zaxis(self, label, unit):
138        """
139        set the z axis label and unit
140        """
141        self._zaxis = label
142        self._zunit = unit
143
144
145class Vector(object):
146    """
147    Vector class to hold multi-dimensional objects
148    """
149    ## x component
150    x = None
151    ## y component
152    y = None
153    ## z component
154    z = None
155
156    def __init__(self, x=None, y=None, z=None):
157        """
158        Initialization. Components that are not
159        set a set to None by default.
160
161        :param x: x component
162        :param y: y component
163        :param z: z component
164        """
165        self.x = x
166        self.y = y
167        self.z = z
168
169    def __str__(self):
170        msg = "x = %s\ty = %s\tz = %s" % (str(self.x), str(self.y), str(self.z))
171        return msg
172
173
174class Detector(object):
175    """
176    Class to hold detector information
177    """
178    ## Name of the instrument [string]
179    name = None
180    ## Sample to detector distance [float] [mm]
181    distance = None
182    distance_unit = 'mm'
183    ## Offset of this detector position in X, Y,
184    #(and Z if necessary) [Vector] [mm]
185    offset = None
186    offset_unit = 'm'
187    ## Orientation (rotation) of this detector in roll,
188    # pitch, and yaw [Vector] [degrees]
189    orientation = None
190    orientation_unit = 'degree'
191    ## Center of the beam on the detector in X and Y
192    #(and Z if necessary) [Vector] [mm]
193    beam_center = None
194    beam_center_unit = 'mm'
195    ## Pixel size in X, Y, (and Z if necessary) [Vector] [mm]
196    pixel_size = None
197    pixel_size_unit = 'mm'
198    ## Slit length of the instrument for this detector.[float] [mm]
199    slit_length = None
200    slit_length_unit = 'mm'
201
202    def __init__(self):
203        """
204        Initialize class attribute that are objects...
205        """
206        self.offset = Vector()
207        self.orientation = Vector()
208        self.beam_center = Vector()
209        self.pixel_size = Vector()
210
211    def __str__(self):
212        _str = "Detector:\n"
213        _str += "   Name:         %s\n" % self.name
214        _str += "   Distance:     %s [%s]\n" % \
215            (str(self.distance), str(self.distance_unit))
216        _str += "   Offset:       %s [%s]\n" % \
217            (str(self.offset), str(self.offset_unit))
218        _str += "   Orientation:  %s [%s]\n" % \
219            (str(self.orientation), str(self.orientation_unit))
220        _str += "   Beam center:  %s [%s]\n" % \
221            (str(self.beam_center), str(self.beam_center_unit))
222        _str += "   Pixel size:   %s [%s]\n" % \
223            (str(self.pixel_size), str(self.pixel_size_unit))
224        _str += "   Slit length:  %s [%s]\n" % \
225            (str(self.slit_length), str(self.slit_length_unit))
226        return _str
227
228
229class Aperture(object):
230    ## Name
231    name = None
232    ## Type
233    type = None
234    ## Size name
235    size_name = None
236    ## Aperture size [Vector]
237    size = None
238    size_unit = 'mm'
239    ## Aperture distance [float]
240    distance = None
241    distance_unit = 'mm'
242
243    def __init__(self):
244        self.size = Vector()
245
246
247class Collimation(object):
248    """
249    Class to hold collimation information
250    """
251    ## Name
252    name = None
253    ## Length [float] [mm]
254    length = None
255    length_unit = 'mm'
256    ## Aperture
257    aperture = None
258
259    def __init__(self):
260        self.aperture = []
261
262    def __str__(self):
263        _str = "Collimation:\n"
264        _str += "   Length:       %s [%s]\n" % \
265            (str(self.length), str(self.length_unit))
266        for item in self.aperture:
267            _str += "   Aperture size:%s [%s]\n" % \
268                (str(item.size), str(item.size_unit))
269            _str += "   Aperture_dist:%s [%s]\n" % \
270                (str(item.distance), str(item.distance_unit))
271        return _str
272
273
274class Source(object):
275    """
276    Class to hold source information
277    """
278    ## Name
279    name = None
280    ## Radiation type [string]
281    radiation = None
282    ## Beam size name
283    beam_size_name = None
284    ## Beam size [Vector] [mm]
285    beam_size = None
286    beam_size_unit = 'mm'
287    ## Beam shape [string]
288    beam_shape = None
289    ## Wavelength [float] [Angstrom]
290    wavelength = None
291    wavelength_unit = 'A'
292    ## Minimum wavelength [float] [Angstrom]
293    wavelength_min = None
294    wavelength_min_unit = 'nm'
295    ## Maximum wavelength [float] [Angstrom]
296    wavelength_max = None
297    wavelength_max_unit = 'nm'
298    ## Wavelength spread [float] [Angstrom]
299    wavelength_spread = None
300    wavelength_spread_unit = 'percent'
301
302    def __init__(self):
303        self.beam_size = Vector()
304
305    def __str__(self):
306        _str = "Source:\n"
307        _str += "   Radiation:    %s\n" % str(self.radiation)
308        _str += "   Shape:        %s\n" % str(self.beam_shape)
309        _str += "   Wavelength:   %s [%s]\n" % \
310            (str(self.wavelength), str(self.wavelength_unit))
311        _str += "   Waveln_min:   %s [%s]\n" % \
312            (str(self.wavelength_min), str(self.wavelength_min_unit))
313        _str += "   Waveln_max:   %s [%s]\n" % \
314            (str(self.wavelength_max), str(self.wavelength_max_unit))
315        _str += "   Waveln_spread:%s [%s]\n" % \
316            (str(self.wavelength_spread), str(self.wavelength_spread_unit))
317        _str += "   Beam_size:    %s [%s]\n" % \
318            (str(self.beam_size), str(self.beam_size_unit))
319        return _str
320
321
322"""
323Definitions of radiation types
324"""
325NEUTRON = 'neutron'
326XRAY = 'x-ray'
327MUON = 'muon'
328ELECTRON = 'electron'
329
330
331class Sample(object):
332    """
333    Class to hold the sample description
334    """
335    ## Short name for sample
336    name = ''
337    ## ID
338    ID = ''
339    ## Thickness [float] [mm]
340    thickness = None
341    thickness_unit = 'mm'
342    ## Transmission [float] [fraction]
343    transmission = None
344    ## Temperature [float] [No Default]
345    temperature = None
346    temperature_unit = None
347    ## Position [Vector] [mm]
348    position = None
349    position_unit = 'mm'
350    ## Orientation [Vector] [degrees]
351    orientation = None
352    orientation_unit = 'degree'
353    ## Details
354    details = None
355    ## SESANS zacceptance
356    zacceptance = (0,"")
357    yacceptance = (0,"")
358
359    def __init__(self):
360        self.position = Vector()
361        self.orientation = Vector()
362        self.details = []
363
364    def __str__(self):
365        _str = "Sample:\n"
366        _str += "   ID:           %s\n" % str(self.ID)
367        _str += "   Transmission: %s\n" % str(self.transmission)
368        _str += "   Thickness:    %s [%s]\n" % \
369            (str(self.thickness), str(self.thickness_unit))
370        _str += "   Temperature:  %s [%s]\n" % \
371            (str(self.temperature), str(self.temperature_unit))
372        _str += "   Position:     %s [%s]\n" % \
373            (str(self.position), str(self.position_unit))
374        _str += "   Orientation:  %s [%s]\n" % \
375            (str(self.orientation), str(self.orientation_unit))
376
377        _str += "   Details:\n"
378        for item in self.details:
379            _str += "      %s\n" % item
380
381        return _str
382
383
384class Process(object):
385    """
386    Class that holds information about the processes
387    performed on the data.
388    """
389    name = ''
390    date = ''
391    description = ''
392    term = None
393    notes = None
394
395    def __init__(self):
396        self.term = []
397        self.notes = []
398
399    def is_empty(self):
400        """
401            Return True if the object is empty
402        """
403        return len(self.name) == 0 and len(self.date) == 0 and len(self.description) == 0 \
404            and len(self.term) == 0 and len(self.notes) == 0
405
406    def single_line_desc(self):
407        """
408            Return a single line string representing the process
409        """
410        return "%s %s %s" % (self.name, self.date, self.description)
411
412    def __str__(self):
413        _str = "Process:\n"
414        _str += "   Name:         %s\n" % self.name
415        _str += "   Date:         %s\n" % self.date
416        _str += "   Description:  %s\n" % self.description
417        for item in self.term:
418            _str += "   Term:         %s\n" % item
419        for item in self.notes:
420            _str += "   Note:         %s\n" % item
421        return _str
422
423
424class TransmissionSpectrum(object):
425    """
426    Class that holds information about transmission spectrum
427    for white beams and spallation sources.
428    """
429    name = ''
430    timestamp = ''
431    ## Wavelength (float) [A]
432    wavelength = None
433    wavelength_unit = 'A'
434    ## Transmission (float) [unit less]
435    transmission = None
436    transmission_unit = ''
437    ## Transmission Deviation (float) [unit less]
438    transmission_deviation = None
439    transmission_deviation_unit = ''
440
441    def __init__(self):
442        self.wavelength = []
443        self.transmission = []
444        self.transmission_deviation = []
445
446    def __str__(self):
447        _str = "Transmission Spectrum:\n"
448        _str += "   Name:             \t{0}\n".format(self.name)
449        _str += "   Timestamp:        \t{0}\n".format(self.timestamp)
450        _str += "   Wavelength unit:  \t{0}\n".format(self.wavelength_unit)
451        _str += "   Transmission unit:\t{0}\n".format(self.transmission_unit)
452        _str += "   Trans. Dev. unit:  \t{0}\n".format(\
453                                            self.transmission_deviation_unit)
454        length_list = [len(self.wavelength), len(self.transmission), \
455                len(self.transmission_deviation)]
456        _str += "   Number of Pts:    \t{0}\n".format(max(length_list))
457        return _str
458
459
460class DataInfo(object):
461    """
462    Class to hold the data read from a file.
463    It includes four blocks of data for the
464    instrument description, the sample description,
465    the data itself and any other meta data.
466    """
467    ## Title
468    title = ''
469    ## Run number
470    run = None
471    ## Run name
472    run_name = None
473    ## File name
474    filename = ''
475    ## Notes
476    notes = None
477    ## Processes (Action on the data)
478    process = None
479    ## Instrument name
480    instrument = ''
481    ## Detector information
482    detector = None
483    ## Sample information
484    sample = None
485    ## Source information
486    source = None
487    ## Collimation information
488    collimation = None
489    ## Transmission Spectrum INfo
490    trans_spectrum = None
491    ## Additional meta-data
492    meta_data = None
493    ## Loading errors
494    errors = None
495    ## SESANS data check
496    isSesans = None
497
498
499    def __init__(self):
500        """
501        Initialization
502        """
503        ## Title
504        self.title = ''
505        ## Run number
506        self.run = []
507        self.run_name = {}
508        ## File name
509        self.filename = ''
510        ## Notes
511        self.notes = []
512        ## Processes (Action on the data)
513        self.process = []
514        ## Instrument name
515        self.instrument = ''
516        ## Detector information
517        self.detector = []
518        ## Sample information
519        self.sample = Sample()
520        ## Source information
521        self.source = Source()
522        ## Collimation information
523        self.collimation = []
524        ## Transmission Spectrum
525        self.trans_spectrum = []
526        ## Additional meta-data
527        self.meta_data = {}
528        ## Loading errors
529        self.errors = []
530        ## SESANS data check
531        self.isSesans = False
532
533    def append_empty_process(self):
534        """
535        """
536        self.process.append(Process())
537
538    def add_notes(self, message=""):
539        """
540        Add notes to datainfo
541        """
542        self.notes.append(message)
543
544    def __str__(self):
545        """
546        Nice printout
547        """
548        _str = "File:            %s\n" % self.filename
549        _str += "Title:           %s\n" % self.title
550        _str += "Run:             %s\n" % str(self.run)
551        _str += "SESANS:          %s\n" % str(self.isSesans)
552        _str += "Instrument:      %s\n" % str(self.instrument)
553        _str += "%s\n" % str(self.sample)
554        _str += "%s\n" % str(self.source)
555        for item in self.detector:
556            _str += "%s\n" % str(item)
557        for item in self.collimation:
558            _str += "%s\n" % str(item)
559        for item in self.process:
560            _str += "%s\n" % str(item)
561        for item in self.notes:
562            _str += "%s\n" % str(item)
563        for item in self.trans_spectrum:
564            _str += "%s\n" % str(item)
565        return _str
566
567    # Private method to perform operation. Not implemented for DataInfo,
568    # but should be implemented for each data class inherited from DataInfo
569    # that holds actual data (ex.: Data1D)
570    def _perform_operation(self, other, operation):
571        """
572        Private method to perform operation. Not implemented for DataInfo,
573        but should be implemented for each data class inherited from DataInfo
574        that holds actual data (ex.: Data1D)
575        """
576        return NotImplemented
577
578    def _perform_union(self, other):
579        """
580        Private method to perform union operation. Not implemented for DataInfo,
581        but should be implemented for each data class inherited from DataInfo
582        that holds actual data (ex.: Data1D)
583        """
584        return NotImplemented
585
586    def __add__(self, other):
587        """
588        Add two data sets
589
590        :param other: data set to add to the current one
591        :return: new data set
592        :raise ValueError: raised when two data sets are incompatible
593        """
594        def operation(a, b):
595            return a + b
596        return self._perform_operation(other, operation)
597
598    def __radd__(self, other):
599        """
600        Add two data sets
601
602        :param other: data set to add to the current one
603        :return: new data set
604        :raise ValueError: raised when two data sets are incompatible
605        """
606        def operation(a, b):
607            return b + a
608        return self._perform_operation(other, operation)
609
610    def __sub__(self, other):
611        """
612        Subtract two data sets
613
614        :param other: data set to subtract from the current one
615        :return: new data set
616        :raise ValueError: raised when two data sets are incompatible
617        """
618        def operation(a, b):
619            return a - b
620        return self._perform_operation(other, operation)
621
622    def __rsub__(self, other):
623        """
624        Subtract two data sets
625
626        :param other: data set to subtract from the current one
627        :return: new data set
628        :raise ValueError: raised when two data sets are incompatible
629        """
630        def operation(a, b):
631            return b - a
632        return self._perform_operation(other, operation)
633
634    def __mul__(self, other):
635        """
636        Multiply two data sets
637
638        :param other: data set to subtract from the current one
639        :return: new data set
640        :raise ValueError: raised when two data sets are incompatible
641        """
642        def operation(a, b):
643            return a * b
644        return self._perform_operation(other, operation)
645
646    def __rmul__(self, other):
647        """
648        Multiply two data sets
649
650        :param other: data set to subtract from the current one
651        :return: new data set
652        :raise ValueError: raised when two data sets are incompatible
653        """
654        def operation(a, b):
655            return b * a
656        return self._perform_operation(other, operation)
657
658    def __div__(self, other):
659        """
660        Divided a data set by another
661
662        :param other: data set that the current one is divided by
663        :return: new data set
664        :raise ValueError: raised when two data sets are incompatible
665        """
666        def operation(a, b):
667            return a/b
668        return self._perform_operation(other, operation)
669
670    def __rdiv__(self, other):
671        """
672        Divided a data set by another
673
674        :param other: data set that the current one is divided by
675        :return: new data set
676        :raise ValueError: raised when two data sets are incompatible
677        """
678        def operation(a, b):
679            return b/a
680        return self._perform_operation(other, operation)
681
682    def __or__(self, other):
683        """
684        Union a data set with another
685
686        :param other: data set to be unified
687        :return: new data set
688        :raise ValueError: raised when two data sets are incompatible
689        """
690        return self._perform_union(other)
691
692    def __ror__(self, other):
693        """
694        Union a data set with another
695
696        :param other: data set to be unified
697        :return: new data set
698        :raise ValueError: raised when two data sets are incompatible
699        """
700        return self._perform_union(other)
701
702class Data1D(plottable_1D, DataInfo):
703    """
704    1D data class
705    """
706    def __init__(self, x=None, y=None, dx=None, dy=None, lam=None, dlam=None, isSesans=None):
707        DataInfo.__init__(self)
708        plottable_1D.__init__(self, x, y, dx, dy,None, None, lam, dlam)
709        self.isSesans = isSesans
710        try:
711            if self.isSesans: # the data is SESANS
712                self.x_unit = 'A'
713                self.y_unit = 'pol'
714            elif not self.isSesans: # the data is SANS
715                self.x_unit = '1/A'
716                self.y_unit = '1/cm'
717        except: # the data is not recognized/supported, and the user is notified
718            raise TypeError('data not recognized, check documentation for supported 1D data formats')
719
720    def __str__(self):
721        """
722        Nice printout
723        """
724        _str = "%s\n" % DataInfo.__str__(self)
725        _str += "Data:\n"
726        _str += "   Type:         %s\n" % self.__class__.__name__
727        _str += "   X-axis:       %s\t[%s]\n" % (self._xaxis, self._xunit)
728        _str += "   Y-axis:       %s\t[%s]\n" % (self._yaxis, self._yunit)
729        _str += "   Length:       %g\n" % len(self.x)
730        return _str
731
732    def is_slit_smeared(self):
733        """
734        Check whether the data has slit smearing information
735        :return: True is slit smearing info is present, False otherwise
736        """
737        def _check(v):
738            if (v.__class__ == list or v.__class__ == np.ndarray) \
739                and len(v) > 0 and min(v) > 0:
740                return True
741            return False
742        return _check(self.dxl) or _check(self.dxw)
743
744    def clone_without_data(self, length=0, clone=None):
745        """
746        Clone the current object, without copying the data (which
747        will be filled out by a subsequent operation).
748        The data arrays will be initialized to zero.
749
750        :param length: length of the data array to be initialized
751        :param clone: if provided, the data will be copied to clone
752        """
753        from copy import deepcopy
754
755        if clone is None or not issubclass(clone.__class__, Data1D):
756            x = np.zeros(length)
757            dx = np.zeros(length)
758            y = np.zeros(length)
759            dy = np.zeros(length)
760            lam = np.zeros(length)
761            dlam = np.zeros(length)
762            clone = Data1D(x, y, lam=lam, dx=dx, dy=dy, dlam=dlam)
763
764        clone.title = self.title
765        clone.run = self.run
766        clone.filename = self.filename
767        clone.instrument = self.instrument
768        clone.notes = deepcopy(self.notes)
769        clone.process = deepcopy(self.process)
770        clone.detector = deepcopy(self.detector)
771        clone.sample = deepcopy(self.sample)
772        clone.source = deepcopy(self.source)
773        clone.collimation = deepcopy(self.collimation)
774        clone.trans_spectrum = deepcopy(self.trans_spectrum)
775        clone.meta_data = deepcopy(self.meta_data)
776        clone.errors = deepcopy(self.errors)
777
778        return clone
779
780    def _validity_check(self, other):
781        """
782        Checks that the data lengths are compatible.
783        Checks that the x vectors are compatible.
784        Returns errors vectors equal to original
785        errors vectors if they were present or vectors
786        of zeros when none was found.
787
788        :param other: other data set for operation
789        :return: dy for self, dy for other [numpy arrays]
790        :raise ValueError: when lengths are not compatible
791        """
792        dy_other = None
793        if isinstance(other, Data1D):
794            # Check that data lengths are the same
795            if len(self.x) != len(other.x) or \
796                len(self.y) != len(other.y):
797                msg = "Unable to perform operation: data length are not equal"
798                raise ValueError(msg)
799            # Here we could also extrapolate between data points
800            TOLERANCE = 0.01
801            for i in range(len(self.x)):
802                if math.fabs((self.x[i] - other.x[i])/self.x[i]) > TOLERANCE:
803                    msg = "Incompatible data sets: x-values do not match"
804                    raise ValueError(msg)
805
806            # Check that the other data set has errors, otherwise
807            # create zero vector
808            dy_other = other.dy
809            if other.dy is None or (len(other.dy) != len(other.y)):
810                dy_other = np.zeros(len(other.y))
811
812        # Check that we have errors, otherwise create zero vector
813        dy = self.dy
814        if self.dy is None or (len(self.dy) != len(self.y)):
815            dy = np.zeros(len(self.y))
816
817        return dy, dy_other
818
819    def _perform_operation(self, other, operation):
820        """
821        """
822        # First, check the data compatibility
823        dy, dy_other = self._validity_check(other)
824        result = self.clone_without_data(len(self.x))
825        if self.dxw is None:
826            result.dxw = None
827        else:
828            result.dxw = np.zeros(len(self.x))
829        if self.dxl is None:
830            result.dxl = None
831        else:
832            result.dxl = np.zeros(len(self.x))
833
834        for i in range(len(self.x)):
835            result.x[i] = self.x[i]
836            if self.dx is not None and len(self.x) == len(self.dx):
837                result.dx[i] = self.dx[i]
838            if self.dxw is not None and len(self.x) == len(self.dxw):
839                result.dxw[i] = self.dxw[i]
840            if self.dxl is not None and len(self.x) == len(self.dxl):
841                result.dxl[i] = self.dxl[i]
842
843            a = Uncertainty(self.y[i], dy[i]**2)
844            if isinstance(other, Data1D):
845                b = Uncertainty(other.y[i], dy_other[i]**2)
846                if other.dx is not None:
847                    result.dx[i] *= self.dx[i]
848                    result.dx[i] += (other.dx[i]**2)
849                    result.dx[i] /= 2
850                    result.dx[i] = math.sqrt(result.dx[i])
851                if result.dxl is not None and other.dxl is not None:
852                    result.dxl[i] *= self.dxl[i]
853                    result.dxl[i] += (other.dxl[i]**2)
854                    result.dxl[i] /= 2
855                    result.dxl[i] = math.sqrt(result.dxl[i])
856            else:
857                b = other
858
859            output = operation(a, b)
860            result.y[i] = output.x
861            result.dy[i] = math.sqrt(math.fabs(output.variance))
862        return result
863
864    def _validity_check_union(self, other):
865        """
866        Checks that the data lengths are compatible.
867        Checks that the x vectors are compatible.
868        Returns errors vectors equal to original
869        errors vectors if they were present or vectors
870        of zeros when none was found.
871
872        :param other: other data set for operation
873        :return: bool
874        :raise ValueError: when data types are not compatible
875        """
876        if not isinstance(other, Data1D):
877            msg = "Unable to perform operation: different types of data set"
878            raise ValueError(msg)
879        return True
880
881    def _perform_union(self, other):
882        """
883        """
884        # First, check the data compatibility
885        self._validity_check_union(other)
886        result = self.clone_without_data(len(self.x) + len(other.x))
887        if self.dy is None or other.dy is None:
888            result.dy = None
889        else:
890            result.dy = np.zeros(len(self.x) + len(other.x))
891        if self.dx is None or other.dx is None:
892            result.dx = None
893        else:
894            result.dx = np.zeros(len(self.x) + len(other.x))
895        if self.dxw is None or other.dxw is None:
896            result.dxw = None
897        else:
898            result.dxw = np.zeros(len(self.x) + len(other.x))
899        if self.dxl is None or other.dxl is None:
900            result.dxl = None
901        else:
902            result.dxl = np.zeros(len(self.x) + len(other.x))
903
904        result.x = np.append(self.x, other.x)
905        #argsorting
906        ind = np.argsort(result.x)
907        result.x = result.x[ind]
908        result.y = np.append(self.y, other.y)
909        result.y = result.y[ind]
910        if result.dy is not None:
911            result.dy = np.append(self.dy, other.dy)
912            result.dy = result.dy[ind]
913        if result.dx is not None:
914            result.dx = np.append(self.dx, other.dx)
915            result.dx = result.dx[ind]
916        if result.dxw is not None:
917            result.dxw = np.append(self.dxw, other.dxw)
918            result.dxw = result.dxw[ind]
919        if result.dxl is not None:
920            result.dxl = np.append(self.dxl, other.dxl)
921            result.dxl = result.dxl[ind]
922        return result
923
924
925class Data2D(plottable_2D, DataInfo):
926    """
927    2D data class
928    """
929    ## Units for Q-values
930    Q_unit = '1/A'
931    ## Units for I(Q) values
932    I_unit = '1/cm'
933    ## Vector of Q-values at the center of each bin in x
934    x_bins = None
935    ## Vector of Q-values at the center of each bin in y
936    y_bins = None
937    ## No 2D SESANS data as of yet. Always set it to False
938    isSesans = False
939
940    def __init__(self, data=None, err_data=None, qx_data=None,
941                 qy_data=None, q_data=None, mask=None,
942                 dqx_data=None, dqy_data=None):
943        DataInfo.__init__(self)
944        plottable_2D.__init__(self, data, err_data, qx_data,
945                              qy_data, q_data, mask, dqx_data, dqy_data)
946        self.y_bins = []
947        self.x_bins = []
948
949        if len(self.detector) > 0:
950            raise RuntimeError("Data2D: Detector bank already filled at init")
951
952    def __str__(self):
953        _str = "%s\n" % DataInfo.__str__(self)
954        _str += "Data:\n"
955        _str += "   Type:         %s\n" % self.__class__.__name__
956        _str += "   X- & Y-axis:  %s\t[%s]\n" % (self._yaxis, self._yunit)
957        _str += "   Z-axis:       %s\t[%s]\n" % (self._zaxis, self._zunit)
958        _str += "   Length:       %g \n" % (len(self.data))
959        _str += "   Shape:        (%d, %d)\n" % (len(self.y_bins), len(self.x_bins))
960        return _str
961
962    def clone_without_data(self, length=0, clone=None):
963        """
964        Clone the current object, without copying the data (which
965        will be filled out by a subsequent operation).
966        The data arrays will be initialized to zero.
967
968        :param length: length of the data array to be initialized
969        :param clone: if provided, the data will be copied to clone
970        """
971        from copy import deepcopy
972
973        if clone is None or not issubclass(clone.__class__, Data2D):
974            data = np.zeros(length)
975            err_data = np.zeros(length)
976            qx_data = np.zeros(length)
977            qy_data = np.zeros(length)
978            q_data = np.zeros(length)
979            mask = np.zeros(length)
980            dqx_data = None
981            dqy_data = None
982            clone = Data2D(data=data, err_data=err_data,
983                           qx_data=qx_data, qy_data=qy_data,
984                           q_data=q_data, mask=mask)
985
986        clone.title = self.title
987        clone.run = self.run
988        clone.filename = self.filename
989        clone.instrument = self.instrument
990        clone.notes = deepcopy(self.notes)
991        clone.process = deepcopy(self.process)
992        clone.detector = deepcopy(self.detector)
993        clone.sample = deepcopy(self.sample)
994        clone.source = deepcopy(self.source)
995        clone.collimation = deepcopy(self.collimation)
996        clone.trans_spectrum = deepcopy(self.trans_spectrum)
997        clone.meta_data = deepcopy(self.meta_data)
998        clone.errors = deepcopy(self.errors)
999
1000        return clone
1001
1002    def _validity_check(self, other):
1003        """
1004        Checks that the data lengths are compatible.
1005        Checks that the x vectors are compatible.
1006        Returns errors vectors equal to original
1007        errors vectors if they were present or vectors
1008        of zeros when none was found.
1009
1010        :param other: other data set for operation
1011        :return: dy for self, dy for other [numpy arrays]
1012        :raise ValueError: when lengths are not compatible
1013        """
1014        err_other = None
1015        TOLERANCE = 0.01
1016        if isinstance(other, Data2D):
1017            # Check that data lengths are the same
1018            if len(self.data) != len(other.data) or \
1019                len(self.qx_data) != len(other.qx_data) or \
1020                len(self.qy_data) != len(other.qy_data):
1021                msg = "Unable to perform operation: data length are not equal"
1022                raise ValueError(msg)
1023            for ind in range(len(self.data)):
1024                if math.fabs((self.qx_data[ind] - other.qx_data[ind])/self.qx_data[ind]) > TOLERANCE:
1025                    msg = "Incompatible data sets: qx-values do not match: %s %s" % (self.qx_data[ind], other.qx_data[ind])
1026                    raise ValueError(msg)
1027                if math.fabs((self.qy_data[ind] - other.qy_data[ind])/self.qy_data[ind]) > TOLERANCE:
1028                    msg = "Incompatible data sets: qy-values do not match: %s %s" % (self.qy_data[ind], other.qy_data[ind])
1029                    raise ValueError(msg)
1030
1031            # Check that the scales match
1032            err_other = other.err_data
1033            if other.err_data is None or \
1034                (len(other.err_data) != len(other.data)):
1035                err_other = np.zeros(len(other.data))
1036
1037        # Check that we have errors, otherwise create zero vector
1038        err = self.err_data
1039        if self.err_data is None or \
1040            (len(self.err_data) != len(self.data)):
1041            err = np.zeros(len(other.data))
1042        return err, err_other
1043
1044    def _perform_operation(self, other, operation):
1045        """
1046        Perform 2D operations between data sets
1047
1048        :param other: other data set
1049        :param operation: function defining the operation
1050        """
1051        # First, check the data compatibility
1052        dy, dy_other = self._validity_check(other)
1053        result = self.clone_without_data(np.size(self.data))
1054        if self.dqx_data is None or self.dqy_data is None:
1055            result.dqx_data = None
1056            result.dqy_data = None
1057        else:
1058            result.dqx_data = np.zeros(len(self.data))
1059            result.dqy_data = np.zeros(len(self.data))
1060        for i in range(np.size(self.data)):
1061            result.data[i] = self.data[i]
1062            if self.err_data is not None and \
1063                            np.size(self.data) == np.size(self.err_data):
1064                result.err_data[i] = self.err_data[i]
1065            if self.dqx_data is not None:
1066                result.dqx_data[i] = self.dqx_data[i]
1067            if self.dqy_data is not None:
1068                result.dqy_data[i] = self.dqy_data[i]
1069            result.qx_data[i] = self.qx_data[i]
1070            result.qy_data[i] = self.qy_data[i]
1071            result.q_data[i] = self.q_data[i]
1072            result.mask[i] = self.mask[i]
1073
1074            a = Uncertainty(self.data[i], dy[i]**2)
1075            if isinstance(other, Data2D):
1076                b = Uncertainty(other.data[i], dy_other[i]**2)
1077                if other.dqx_data is not None and \
1078                        result.dqx_data is not None:
1079                    result.dqx_data[i] *= self.dqx_data[i]
1080                    result.dqx_data[i] += (other.dqx_data[i]**2)
1081                    result.dqx_data[i] /= 2
1082                    result.dqx_data[i] = math.sqrt(result.dqx_data[i])
1083                if other.dqy_data is not None and \
1084                        result.dqy_data is not None:
1085                    result.dqy_data[i] *= self.dqy_data[i]
1086                    result.dqy_data[i] += (other.dqy_data[i]**2)
1087                    result.dqy_data[i] /= 2
1088                    result.dqy_data[i] = math.sqrt(result.dqy_data[i])
1089            else:
1090                b = other
1091            output = operation(a, b)
1092            result.data[i] = output.x
1093            result.err_data[i] = math.sqrt(math.fabs(output.variance))
1094        return result
1095
1096    def _validity_check_union(self, other):
1097        """
1098        Checks that the data lengths are compatible.
1099        Checks that the x vectors are compatible.
1100        Returns errors vectors equal to original
1101        errors vectors if they were present or vectors
1102        of zeros when none was found.
1103
1104        :param other: other data set for operation
1105        :return: bool
1106        :raise ValueError: when data types are not compatible
1107        """
1108        if not isinstance(other, Data2D):
1109            msg = "Unable to perform operation: different types of data set"
1110            raise ValueError(msg)
1111        return True
1112
1113    def _perform_union(self, other):
1114        """
1115        Perform 2D operations between data sets
1116
1117        :param other: other data set
1118        :param operation: function defining the operation
1119        """
1120        # First, check the data compatibility
1121        self._validity_check_union(other)
1122        result = self.clone_without_data(np.size(self.data) + \
1123                                         np.size(other.data))
1124        result.xmin = self.xmin
1125        result.xmax = self.xmax
1126        result.ymin = self.ymin
1127        result.ymax = self.ymax
1128        if self.dqx_data is None or self.dqy_data is None or \
1129                other.dqx_data is None or other.dqy_data is None:
1130            result.dqx_data = None
1131            result.dqy_data = None
1132        else:
1133            result.dqx_data = np.zeros(len(self.data) + \
1134                                       np.size(other.data))
1135            result.dqy_data = np.zeros(len(self.data) + \
1136                                       np.size(other.data))
1137
1138        result.data = np.append(self.data, other.data)
1139        result.qx_data = np.append(self.qx_data, other.qx_data)
1140        result.qy_data = np.append(self.qy_data, other.qy_data)
1141        result.q_data = np.append(self.q_data, other.q_data)
1142        result.mask = np.append(self.mask, other.mask)
1143        if result.err_data is not None:
1144            result.err_data = np.append(self.err_data, other.err_data)
1145        if self.dqx_data is not None:
1146            result.dqx_data = np.append(self.dqx_data, other.dqx_data)
1147        if self.dqy_data is not None:
1148            result.dqy_data = np.append(self.dqy_data, other.dqy_data)
1149
1150        return result
1151
1152
1153def combine_data_info_with_plottable(data, datainfo):
1154    """
1155    A function that combines the DataInfo data in self.current_datainto with a plottable_1D or 2D data object.
1156
1157    :param data: A plottable_1D or plottable_2D data object
1158    :return: A fully specified Data1D or Data2D object
1159    """
1160
1161    final_dataset = None
1162    if isinstance(data, plottable_1D):
1163        final_dataset = Data1D(data.x, data.y, isSesans=datainfo.isSesans)
1164        final_dataset.dx = data.dx
1165        final_dataset.dy = data.dy
1166        final_dataset.dxl = data.dxl
1167        final_dataset.dxw = data.dxw
1168        final_dataset.x_unit = data._xunit
1169        final_dataset.y_unit = data._yunit
1170        final_dataset.xaxis(data._xaxis, data._xunit)
1171        final_dataset.yaxis(data._yaxis, data._yunit)
1172    elif isinstance(data, plottable_2D):
1173        final_dataset = Data2D(data.data, data.err_data, data.qx_data, data.qy_data, data.q_data,
1174                               data.mask, data.dqx_data, data.dqy_data)
1175        final_dataset.xaxis(data._xaxis, data._xunit)
1176        final_dataset.yaxis(data._yaxis, data._yunit)
1177        final_dataset.zaxis(data._zaxis, data._zunit)
1178        if len(data.data.shape) == 2:
1179            n_rows, n_cols = data.data.shape
1180            final_dataset.y_bins = data.qy_data[0::int(n_cols)]
1181            final_dataset.x_bins = data.qx_data[:int(n_cols)]
1182    else:
1183        return_string = "Should Never Happen: _combine_data_info_with_plottable input is not a plottable1d or " + \
1184                        "plottable2d data object"
1185        return return_string
1186
1187    if hasattr(data, "xmax"):
1188        final_dataset.xmax = data.xmax
1189    if hasattr(data, "ymax"):
1190        final_dataset.ymax = data.ymax
1191    if hasattr(data, "xmin"):
1192        final_dataset.xmin = data.xmin
1193    if hasattr(data, "ymin"):
1194        final_dataset.ymin = data.ymin
1195    final_dataset.isSesans = datainfo.isSesans
1196    final_dataset.title = datainfo.title
1197    final_dataset.run = datainfo.run
1198    final_dataset.run_name = datainfo.run_name
1199    final_dataset.filename = datainfo.filename
1200    final_dataset.notes = datainfo.notes
1201    final_dataset.process = datainfo.process
1202    final_dataset.instrument = datainfo.instrument
1203    final_dataset.detector = datainfo.detector
1204    final_dataset.sample = datainfo.sample
1205    final_dataset.source = datainfo.source
1206    final_dataset.collimation = datainfo.collimation
1207    final_dataset.trans_spectrum = datainfo.trans_spectrum
1208    final_dataset.meta_data = datainfo.meta_data
1209    final_dataset.errors = datainfo.errors
1210    return final_dataset
Note: See TracBrowser for help on using the repository browser.