source: sasview/DataLoader/data_info.py @ 4a7ad5f

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.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change on this file since 4a7ad5f was a7a5886, checked in by Gervaise Alina <gervyh@…>, 14 years ago

working pylint
  • Property mode set to 100644
File size: 28.3 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
12#This software was developed by the University of Tennessee as part of the
13#Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
14#project funded by the US National Science Foundation.
15#If you use DANSE applications to do scientific research that leads to
16#publication, we ask that you acknowledge the use of the software with the
17#following sentence:
18#This work benefited from DANSE software developed under NSF award DMR-0520547.
19#copyright 2008, University of Tennessee
20
21
22#TODO: Keep track of data manipulation in the 'process' data structure.
23#TODO: This module should be independent of plottables. We should write
24#        an adapter class for plottables when needed.
25
26#from sans.guitools.plottables import Data1D as plottable_1D
27from data_util.uncertainty import Uncertainty
28import numpy
29import math
30
31class plottable_1D:
32    """
33    Data1D is a place holder for 1D plottables.
34    """
35    # The presence of these should be mutually
36    # exclusive with the presence of Qdev (dx)
37    x = None
38    y = None
39    dx = None
40    dy = None
41    ## Slit smearing length
42    dxl = None
43    ## Slit smearing width
44    dxw = None
45   
46    # Units
47    _xaxis = ''
48    _xunit = ''
49    _yaxis = ''
50    _yunit = ''
51   
52    def __init__(self, x, y, dx=None, dy=None, dxl=None, dxw=None):
53        self.x = numpy.asarray(x)
54        self.y = numpy.asarray(y)
55        if dx is not None:
56             self.dx = numpy.asarray(dx)
57        if dy is not None: 
58            self.dy = numpy.asarray(dy)
59        if dxl is not None: 
60            self.dxl = numpy.asarray(dxl)
61        if dxw is not None: 
62            self.dxw = numpy.asarray(dxw)
63
64    def xaxis(self, label, unit):
65        """
66        set the x axis label and unit
67        """
68        self._xaxis = label
69        self._xunit = unit
70       
71    def yaxis(self, label, unit):
72        """
73        set the y axis label and unit
74        """
75        self._yaxis = label
76        self._yunit = unit
77
78class plottable_2D:
79    """
80    Data2D is a place holder for 2D plottables.
81    """
82    xmin = None
83    xmax = None
84    ymin = None
85    ymax = None
86    data = None
87    qx_data     = None
88    qy_data     = None
89    q_data      = None
90    err_data    = None
91    dqx_data    = None
92    dqy_data    = None
93    mask        = None
94   
95    # Units
96    _xaxis = ''
97    _xunit = ''
98    _yaxis = ''
99    _yunit = ''
100    _zaxis = ''
101    _zunit = ''
102   
103    def __init__(self, data=None, err_data=None, qx_data=None,
104                  qy_data=None, q_data=None, mask=None,
105                   dqx_data=None, dqy_data=None):
106        self.data = numpy.asarray(data)
107        self.qx_data = numpy.asarray(qx_data)
108        self.qy_data = numpy.asarray(qy_data)
109        self.q_data = numpy.asarray(q_data)
110        self.mask = numpy.asarray(mask)
111        self.err_data = numpy.asarray(err_data)
112        if dqx_data is not None: self.dqx_data = numpy.asarray(dqx_data) 
113        if dqy_data is not None: self.dqy_data = numpy.asarray(dqy_data) 
114               
115    def xaxis(self, label, unit):
116        """
117        set the x axis label and unit
118        """
119        self._xaxis = label
120        self._xunit = unit
121       
122    def yaxis(self, label, unit):
123        """
124        set the y axis label and unit
125        """
126        self._yaxis = label
127        self._yunit = unit
128           
129    def zaxis(self, label, unit):
130        """
131        set the z axis label and unit
132        """
133        self._zaxis = label
134        self._zunit = unit
135
136           
137class Vector:
138    """
139    Vector class to hold multi-dimensional objects
140    """
141    ## x component
142    x = None
143    ## y component
144    y = None
145    ## z component
146    z = None
147   
148    def __init__(self, x=None, y=None, z=None):
149        """
150        Initialization. Components that are not
151        set a set to None by default.
152       
153        :param x: x component
154        :param y: y component
155        :param z: z component
156       
157        """
158        self.x = x
159        self.y = y
160        self.z = z
161       
162    def __str__(self):
163        msg = "x = %s\ty = %s\tz = %s" % (str(self.x), str(self.y), str(self.z))
164        return msg
165       
166
167class Detector:
168    """
169    Class to hold detector information
170    """
171    ## Name of the instrument [string]
172    name = None
173    ## Sample to detector distance [float] [mm]
174    distance = None
175    distance_unit = 'mm'
176    ## Offset of this detector position in X, Y,
177    #(and Z if necessary) [Vector] [mm]
178    offset = None
179    offset_unit = 'm'
180    ## Orientation (rotation) of this detector in roll,
181    # pitch, and yaw [Vector] [degrees]
182    orientation = None
183    orientation_unit = 'degree'
184    ## Center of the beam on the detector in X and Y
185    #(and Z if necessary) [Vector] [mm]
186    beam_center = None
187    beam_center_unit = 'mm'
188    ## Pixel size in X, Y, (and Z if necessary) [Vector] [mm]
189    pixel_size = None
190    pixel_size_unit = 'mm'
191    ## Slit length of the instrument for this detector.[float] [mm]
192    slit_length = None
193    slit_length_unit = 'mm'
194   
195    def __init__(self):
196        """
197       
198        Initialize class attribute that are objects...
199       
200        """
201        self.offset      = Vector()
202        self.orientation = Vector()
203        self.beam_center = Vector()
204        self.pixel_size  = Vector()
205       
206   
207    def __str__(self):
208        _str  = "Detector:\n"
209        _str += "   Name:         %s\n" % self.name
210        _str += "   Distance:     %s [%s]\n" % \
211            (str(self.distance), str(self.distance_unit))
212        _str += "   Offset:       %s [%s]\n" % \
213            (str(self.offset), str(self.offset_unit))
214        _str += "   Orientation:  %s [%s]\n" % \
215            (str(self.orientation), str(self.orientation_unit))
216        _str += "   Beam center:  %s [%s]\n" % \
217            (str(self.beam_center), str(self.beam_center_unit))
218        _str += "   Pixel size:   %s [%s]\n" % \
219            (str(self.pixel_size), str(self.pixel_size_unit))
220        _str += "   Slit length:  %s [%s]\n" % \
221            (str(self.slit_length), str(self.slit_length_unit))
222        return _str
223
224class Aperture:
225    ## Name
226    name = None
227    ## Type
228    type = None
229    ## Size name
230    size_name = None
231    ## Aperture size [Vector]
232    size = None
233    size_unit = 'mm'
234    ## Aperture distance [float]
235    distance = None
236    distance_unit = 'mm'
237   
238    def __init__(self):
239        self.size = Vector()
240   
241class Collimation:
242    """
243    Class to hold collimation information
244    """
245    ## Name
246    name = None
247    ## Length [float] [mm]
248    length = None
249    length_unit = 'mm'
250    ## Aperture
251    aperture = None
252   
253    def __init__(self):
254        self.aperture = []
255   
256    def __str__(self):
257        _str = "Collimation:\n"
258        _str += "   Length:       %s [%s]\n" % \
259            (str(self.length), str(self.length_unit))
260        for item in self.aperture:
261            _str += "   Aperture size:%s [%s]\n" % \
262                (str(item.size), str(item.size_unit))
263            _str += "   Aperture_dist:%s [%s]\n" % \
264                (str(item.distance), str(item.distance_unit))
265        return _str
266
267class Source:
268    """
269    Class to hold source information
270    """ 
271    ## Name
272    name = None
273    ## Radiation type [string]
274    radiation = None
275    ## Beam size name
276    beam_size_name = None
277    ## Beam size [Vector] [mm]
278    beam_size = None
279    beam_size_unit = 'mm'
280    ## Beam shape [string]
281    beam_shape = None
282    ## Wavelength [float] [Angstrom]
283    wavelength = None
284    wavelength_unit = 'A'
285    ## Minimum wavelength [float] [Angstrom]
286    wavelength_min = None
287    wavelength_min_unit = 'nm'
288    ## Maximum wavelength [float] [Angstrom]
289    wavelength_max = None
290    wavelength_max_unit = 'nm'
291    ## Wavelength spread [float] [Angstrom]
292    wavelength_spread = None
293    wavelength_spread_unit = 'percent'
294   
295    def __init__(self):
296        self.beam_size = Vector()
297       
298   
299    def __str__(self):
300        _str  = "Source:\n"
301        _str += "   Radiation:    %s\n" % str(self.radiation)
302        _str += "   Shape:        %s\n" % str(self.beam_shape)
303        _str += "   Wavelength:   %s [%s]\n" % \
304            (str(self.wavelength), str(self.wavelength_unit))
305        _str += "   Waveln_min:   %s [%s]\n" % \
306            (str(self.wavelength_min), str(self.wavelength_min_unit))
307        _str += "   Waveln_max:   %s [%s]\n" % \
308            (str(self.wavelength_max), str(self.wavelength_max_unit))
309        _str += "   Waveln_spread:%s [%s]\n" % \
310            (str(self.wavelength_spread), str(self.wavelength_spread_unit))
311        _str += "   Beam_size:    %s [%s]\n" % \
312            (str(self.beam_size), str(self.beam_size_unit))
313        return _str
314   
315   
316"""
317Definitions of radiation types
318"""
319NEUTRON  = 'neutron'
320XRAY     = 'x-ray'
321MUON     = 'muon'
322ELECTRON = 'electron'
323   
324class Sample:
325    """
326    Class to hold the sample description
327    """
328    ## Short name for sample
329    name = ''
330    ## ID
331    ID = ''
332    ## Thickness [float] [mm]
333    thickness = None
334    thickness_unit = 'mm'
335    ## Transmission [float] [fraction]
336    transmission = None
337    ## Temperature [float] [C]
338    temperature = None
339    temperature_unit = 'C'
340    ## Position [Vector] [mm]
341    position = None
342    position_unit = 'mm'
343    ## Orientation [Vector] [degrees]
344    orientation = None
345    orientation_unit = 'degree'
346    ## Details
347    details = None
348   
349    def __init__(self):
350        self.position    = Vector()
351        self.orientation = Vector()
352        self.details     = []
353   
354    def __str__(self):
355        _str  = "Sample:\n"
356        _str += "   ID:           %s\n" % str(self.ID)
357        _str += "   Transmission: %s\n" % str(self.transmission)
358        _str += "   Thickness:    %s [%s]\n" % \
359            (str(self.thickness), str(self.thickness_unit))
360        _str += "   Temperature:  %s [%s]\n" % \
361            (str(self.temperature), str(self.temperature_unit))
362        _str += "   Position:     %s [%s]\n" % \
363            (str(self.position), str(self.position_unit))
364        _str += "   Orientation:  %s [%s]\n" % \
365            (str(self.orientation), str(self.orientation_unit))
366       
367        _str += "   Details:\n"
368        for item in self.details:
369            _str += "      %s\n" % item
370           
371        return _str
372 
373class Process:
374    """
375    Class that holds information about the processes
376    performed on the data.
377    """
378    name = ''
379    date = ''
380    description = ''
381    term = None
382    notes = None
383   
384    def __init__(self):
385        self.term = []
386        self.notes = []
387   
388    def __str__(self):
389        _str  = "Process:\n"
390        _str += "   Name:         %s\n" % self.name
391        _str += "   Date:         %s\n" % self.date
392        _str += "   Description:  %s\n" % self.description
393        for item in self.term:
394            _str += "   Term:         %s\n" % item
395        for item in self.notes:
396            _str += "   Note:         %s\n" % item
397        return _str
398   
399 
400class DataInfo:
401    """
402    Class to hold the data read from a file.
403    It includes four blocks of data for the
404    instrument description, the sample description,
405    the data itself and any other meta data.
406    """
407    ## Title
408    title      = ''
409    ## Run number
410    run        = None
411    ## Run name
412    run_name   = None
413    ## File name
414    filename   = ''
415    ## Notes
416    notes      = None
417    ## Processes (Action on the data)
418    process    = None
419    ## Instrument name
420    instrument = ''
421    ## Detector information
422    detector   = None
423    ## Sample information
424    sample     = None
425    ## Source information
426    source     = None
427    ## Collimation information
428    collimation = None
429    ## Additional meta-data
430    meta_data  = None
431    ## Loading errors
432    errors = None
433           
434    def __init__(self):
435        """
436        Initialization
437        """
438        ## Title
439        self.title      = ''
440        ## Run number
441        self.run        = []
442        self.run_name   = {}
443        ## File name
444        self.filename   = ''
445        ## Notes
446        self.notes      = []
447        ## Processes (Action on the data)
448        self.process    = []
449        ## Instrument name
450        self.instrument = ''
451        ## Detector information
452        self.detector   = []
453        ## Sample information
454        self.sample     = Sample()
455        ## Source information
456        self.source     = Source()
457        ## Collimation information
458        self.collimation = []
459        ## Additional meta-data
460        self.meta_data  = {}
461        ## Loading errors
462        self.errors = []       
463       
464    def __str__(self):
465        """
466        Nice printout
467        """
468        _str =  "File:            %s\n" % self.filename
469        _str += "Title:           %s\n" % self.title
470        _str += "Run:             %s\n" % str(self.run)
471        _str += "Instrument:      %s\n" % str(self.instrument)
472        _str += "%s\n" % str(self.sample)
473        _str += "%s\n" % str(self.source)
474        for item in self.detector:
475            _str += "%s\n" % str(item)
476        for item in self.collimation:
477            _str += "%s\n" % str(item)
478        for item in self.process:
479            _str += "%s\n" % str(item)
480        for item in self.notes:
481            _str += "%s\n" % str(item)
482
483        return _str
484           
485    # Private method to perform operation. Not implemented for DataInfo,
486    # but should be implemented for each data class inherited from DataInfo
487    # that holds actual data (ex.: Data1D)
488    def _perform_operation(self, other, operation): 
489        """
490        Private method to perform operation. Not implemented for DataInfo,
491        but should be implemented for each data class inherited from DataInfo
492        that holds actual data (ex.: Data1D)
493        """
494        return NotImplemented
495
496    def __add__(self, other):
497        """
498        Add two data sets
499       
500        :param other: data set to add to the current one
501        :return: new data set
502        :raise ValueError: raised when two data sets are incompatible
503        """
504        def operation(a, b):
505            return a + b
506        return self._perform_operation(other, operation)
507       
508    def __radd__(self, other):
509        """
510        Add two data sets
511       
512        :param other: data set to add to the current one
513       
514        :return: new data set
515       
516        :raise ValueError: raised when two data sets are incompatible
517       
518        """
519        def operation(a, b):
520            return b + a
521        return self._perform_operation(other, operation)
522       
523    def __sub__(self, other):
524        """
525        Subtract two data sets
526       
527        :param other: data set to subtract from the current one
528       
529        :return: new data set
530       
531        :raise ValueError: raised when two data sets are incompatible
532       
533        """
534        def operation(a, b):
535            return a - b
536        return self._perform_operation(other, operation)
537       
538    def __rsub__(self, other):
539        """
540        Subtract two data sets
541       
542        :param other: data set to subtract from the current one
543       
544        :return: new data set
545       
546        :raise ValueError: raised when two data sets are incompatible
547       
548        """
549        def operation(a, b):
550            return b - a
551        return self._perform_operation(other, operation)
552       
553    def __mul__(self, other):
554        """
555        Multiply two data sets
556       
557        :param other: data set to subtract from the current one
558       
559        :return: new data set
560       
561        :raise ValueError: raised when two data sets are incompatible
562       
563        """
564        def operation(a, b):
565            return a * b
566        return self._perform_operation(other, operation)
567       
568    def __rmul__(self, other):
569        """
570        Multiply two data sets
571       
572        :param other: data set to subtract from the current one
573       
574        :return: new data set
575       
576        :raise ValueError: raised when two data sets are incompatible
577        """
578        def operation(a, b):
579            return b * a
580        return self._perform_operation(other, operation)
581       
582    def __div__(self, other):
583        """
584        Divided a data set by another
585       
586        :param other: data set that the current one is divided by
587       
588        :return: new data set
589       
590        :raise ValueError: raised when two data sets are incompatible
591       
592        """
593        def operation(a, b):
594            return a/b
595        return self._perform_operation(other, operation)
596       
597    def __rdiv__(self, other):
598        """
599        Divided a data set by another
600       
601        :param other: data set that the current one is divided by
602       
603        :return: new data set
604       
605        :raise ValueError: raised when two data sets are incompatible
606       
607        """
608        def operation(a, b):
609            return b/a
610        return self._perform_operation(other, operation)           
611           
612class Data1D(plottable_1D, DataInfo):
613    """
614    1D data class
615    """
616    x_unit = '1/A'
617    y_unit = '1/cm'
618   
619    def __init__(self, x, y, dx=None, dy=None):
620        DataInfo.__init__(self)
621        plottable_1D.__init__(self, x, y, dx, dy)
622       
623       
624    def __str__(self):
625        """
626        Nice printout
627        """
628        _str =  "%s\n" % DataInfo.__str__(self)
629   
630        _str += "Data:\n"
631        _str += "   Type:         %s\n" % self.__class__.__name__
632        _str += "   X-axis:       %s\t[%s]\n" % (self._xaxis, self._xunit)
633        _str += "   Y-axis:       %s\t[%s]\n" % (self._yaxis, self._yunit)
634        _str += "   Length:       %g\n" % len(self.x)
635
636        return _str
637
638    def is_slit_smeared(self):
639        """
640        Check whether the data has slit smearing information
641       
642        :return: True is slit smearing info is present, False otherwise
643       
644        """
645        def _check(v):           
646            if (v.__class__ == list or v.__class__ == numpy.ndarray) \
647                and len(v) > 0 and min(v) > 0:
648                return True
649           
650            return False
651       
652        return _check(self.dxl) or _check(self.dxw)
653       
654    def clone_without_data(self, length=0, clone=None):
655        """
656        Clone the current object, without copying the data (which
657        will be filled out by a subsequent operation).
658        The data arrays will be initialized to zero.
659       
660        :param length: length of the data array to be initialized
661        :param clone: if provided, the data will be copied to clone
662        """
663        from copy import deepcopy
664       
665        if clone is None or not issubclass(clone.__class__, Data1D):
666            x  = numpy.zeros(length) 
667            dx = numpy.zeros(length) 
668            y  = numpy.zeros(length) 
669            dy = numpy.zeros(length) 
670            clone = Data1D(x, y, dx=dx, dy=dy)
671       
672        clone.title       = self.title
673        clone.run         = self.run
674        clone.filename    = self.filename
675        clone.notes       = deepcopy(self.notes) 
676        clone.process     = deepcopy(self.process) 
677        clone.detector    = deepcopy(self.detector) 
678        clone.sample      = deepcopy(self.sample) 
679        clone.source      = deepcopy(self.source) 
680        clone.collimation = deepcopy(self.collimation) 
681        clone.meta_data   = deepcopy(self.meta_data) 
682        clone.errors      = deepcopy(self.errors) 
683       
684        return clone
685
686    def _validity_check(self, other):
687        """
688        Checks that the data lengths are compatible.
689        Checks that the x vectors are compatible.
690        Returns errors vectors equal to original
691        errors vectors if they were present or vectors
692        of zeros when none was found.
693       
694        :param other: other data set for operation
695       
696        :return: dy for self, dy for other [numpy arrays]
697       
698        :raise ValueError: when lengths are not compatible
699       
700        """
701        dy_other = None
702        if isinstance(other, Data1D):
703            # Check that data lengths are the same
704            if len(self.x) != len(other.x) or \
705                len(self.y) != len(other.y):
706                msg =  "Unable to perform operation: data length are not equal"
707                raise ValueError, msg
708           
709            # Here we could also extrapolate between data points
710            for i in range(len(self.x)):
711                if self.x[i] != other.x[i]:
712                    msg = "Incompatible data sets: x-values do not match"
713                    raise ValueError, msg
714           
715            # Check that the other data set has errors, otherwise
716            # create zero vector
717            dy_other = other.dy
718            if other.dy == None or (len(other.dy) != len(other.y)):
719                dy_other = numpy.zeros(len(other.y))
720           
721        # Check that we have errors, otherwise create zero vector
722        dy = self.dy
723        if self.dy == None or (len(self.dy) != len(self.y)):
724            dy = numpy.zeros(len(self.y))           
725           
726        return dy, dy_other
727
728    def _perform_operation(self, other, operation):
729        """
730        """
731        # First, check the data compatibility
732        dy, dy_other = self._validity_check(other)
733        result = self.clone_without_data(len(self.x))
734       
735        for i in range(len(self.x)):
736            result.x[i] = self.x[i]
737            if self.dx is not None and len(self.x) == len(self.dx):
738                result.dx[i] = self.dx[i]
739           
740            a = Uncertainty(self.y[i], dy[i]**2)
741            if isinstance(other, Data1D):
742                b = Uncertainty(other.y[i], dy_other[i]**2)
743            else:
744                b = other
745           
746            output = operation(a, b)
747            result.y[i] = output.x
748            result.dy[i] = math.sqrt(math.fabs(output.variance))
749        return result
750       
751class Data2D(plottable_2D, DataInfo):
752    """
753    2D data class
754    """
755    ## Units for Q-values
756    Q_unit = '1/A'
757   
758    ## Units for I(Q) values
759    I_unit = '1/cm'
760   
761    ## Vector of Q-values at the center of each bin in x
762    x_bins = None
763   
764    ## Vector of Q-values at the center of each bin in y
765    y_bins = None
766   
767   
768    def __init__(self, data=None, err_data=None, qx_data=None,
769                  qy_data=None, q_data=None, mask=None, 
770                  dqx_data=None, dqy_data=None):
771        self.y_bins = []
772        self.x_bins = []
773        DataInfo.__init__(self)
774        plottable_2D.__init__(self, data, err_data, qx_data,
775                              qy_data, q_data,mask, dqx_data, dqy_data)
776        if len(self.detector) > 0:
777            raise RuntimeError, "Data2D: Detector bank already filled at init"
778
779    def __str__(self):
780        _str =  "%s\n" % DataInfo.__str__(self)
781       
782        _str += "Data:\n"
783        _str += "   Type:         %s\n" % self.__class__.__name__
784        _str += "   X- & Y-axis:  %s\t[%s]\n" % (self._yaxis, self._yunit)
785        _str += "   Z-axis:       %s\t[%s]\n" % (self._zaxis, self._zunit)
786        #leny = 0
787        #if len(self.data) > 0:
788        #    leny = len(self.data)
789        _str += "   Length:       %g \n" % (len(self.data))
790       
791        return _str
792 
793    def clone_without_data(self, length=0, clone=None):
794        """
795        Clone the current object, without copying the data (which
796        will be filled out by a subsequent operation).
797        The data arrays will be initialized to zero.
798       
799        :param length: length of the data array to be initialized
800        :param clone: if provided, the data will be copied to clone
801        """
802        from copy import deepcopy
803       
804        if clone is None or not issubclass(clone.__class__, Data2D): 
805            data     = numpy.zeros(length) 
806            err_data = numpy.zeros(length) 
807            qx_data = numpy.zeros(length)
808            qy_data = numpy.zeros(length)
809            q_data = numpy.zeros(length)
810            mask = numpy.zeros(length)
811            dqx_data = None
812            dqy_data = None
813            clone = Data2D(data, err_data, qx_data, qy_data,
814                            q_data,mask, dqx_data=dqx_data, dqy_data=dqy_data)
815
816        clone.title       = self.title
817        clone.run         = self.run
818        clone.filename    = self.filename
819        clone.notes       = deepcopy(self.notes) 
820        clone.process     = deepcopy(self.process) 
821        clone.detector    = deepcopy(self.detector) 
822        clone.sample      = deepcopy(self.sample) 
823        clone.source      = deepcopy(self.source) 
824        clone.collimation = deepcopy(self.collimation) 
825        clone.meta_data   = deepcopy(self.meta_data) 
826        clone.errors      = deepcopy(self.errors) 
827       
828        return clone
829 
830 
831    def _validity_check(self, other):
832        """
833        Checks that the data lengths are compatible.
834        Checks that the x vectors are compatible.
835        Returns errors vectors equal to original
836        errors vectors if they were present or vectors
837        of zeros when none was found.
838       
839        :param other: other data set for operation
840       
841        :return: dy for self, dy for other [numpy arrays]
842       
843        :raise ValueError: when lengths are not compatible
844       
845        """
846        err_other = None
847        if isinstance(other, Data2D):
848            # Check that data lengths are the same
849            if numpy.size(self.data) != numpy.size(other.data):
850                msg = "Unable to perform operation: data length are not equal"
851                raise ValueError, msg
852               
853            # Check that the scales match
854            #TODO: matching scales?     
855           
856            # Check that the other data set has errors, otherwise
857            # create zero vector
858            #TODO: test this
859            err_other = other.err_data
860            if other.err_data == None or \
861                (numpy.size(other.err_data) != numpy.size(other.data)):
862                err_other = numpy.zeros([numpy.size(other.data, 0),
863                                          numpy.size(other.data, 1)])
864           
865        # Check that we have errors, otherwise create zero vector
866        err = self.err_data
867        if self.err_data == None or \
868            (numpy.size(self.err_data) != numpy.size(self.data)):
869            err = numpy.zeros([numpy.size(self.data, 0),
870                                numpy.size(self.data, 1)])
871           
872        return err, err_other
873 
874 
875    def _perform_operation(self, other, operation):
876        """
877        Perform 2D operations between data sets
878       
879        :param other: other data set
880        :param operation: function defining the operation
881       
882        """
883        # First, check the data compatibility
884        dy, dy_other = self._validity_check(other)
885   
886        result = self.clone_without_data([numpy.size(self.data, 0), 
887                                          numpy.size(self.data, 1)])
888       
889        for i in range(numpy.size(self.data, 0)):
890            for j in range(numpy.size(self.data, 1)):
891                result.data[i][j] = self.data[i][j]
892                if self.err_data is not None and \
893                    numpy.size(self.data) == numpy.size(self.err_data):
894                    result.err_data[i][j] = self.err_data[i][j]
895               
896                a = Uncertainty(self.data[i][j], dy[i][j]**2)
897                if isinstance(other, Data2D):
898                    b = Uncertainty(other.data[i][j], dy_other[i][j]**2)
899                else:
900                    b = other
901               
902                output = operation(a, b)
903                result.data[i][j] = output.x
904                result.err_data[i][j] = math.sqrt(math.fabs(output.variance))
905        return result
906   
907 
908 
Note: See TracBrowser for help on using the repository browser.