source: sasview/sansdataloader/src/sans/dataloader/data_info.py @ 015eb2b

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 015eb2b was e38e335, checked in by Jae Cho <jhjcho@…>, 13 years ago

added instrument (name) for writes xml
  • Property mode set to 100644
File size: 28.6 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 append_empty_process(self):
465        """
466        """
467        self.process.append(Process())
468       
469    def add_notes(self, message=""):
470        """
471        Add notes to datainfo
472        """
473        self.notes.append(message)
474       
475    def __str__(self):
476        """
477        Nice printout
478        """
479        _str =  "File:            %s\n" % self.filename
480        _str += "Title:           %s\n" % self.title
481        _str += "Run:             %s\n" % str(self.run)
482        _str += "Instrument:      %s\n" % str(self.instrument)
483        _str += "%s\n" % str(self.sample)
484        _str += "%s\n" % str(self.source)
485        for item in self.detector:
486            _str += "%s\n" % str(item)
487        for item in self.collimation:
488            _str += "%s\n" % str(item)
489        for item in self.process:
490            _str += "%s\n" % str(item)
491        for item in self.notes:
492            _str += "%s\n" % str(item)
493
494        return _str
495           
496    # Private method to perform operation. Not implemented for DataInfo,
497    # but should be implemented for each data class inherited from DataInfo
498    # that holds actual data (ex.: Data1D)
499    def _perform_operation(self, other, operation): 
500        """
501        Private method to perform operation. Not implemented for DataInfo,
502        but should be implemented for each data class inherited from DataInfo
503        that holds actual data (ex.: Data1D)
504        """
505        return NotImplemented
506
507    def __add__(self, other):
508        """
509        Add two data sets
510       
511        :param other: data set to add to the current one
512        :return: new data set
513        :raise ValueError: raised when two data sets are incompatible
514        """
515        def operation(a, b):
516            return a + b
517        return self._perform_operation(other, operation)
518       
519    def __radd__(self, other):
520        """
521        Add two data sets
522       
523        :param other: data set to add to the current one
524       
525        :return: new data set
526       
527        :raise ValueError: raised when two data sets are incompatible
528       
529        """
530        def operation(a, b):
531            return b + a
532        return self._perform_operation(other, operation)
533       
534    def __sub__(self, other):
535        """
536        Subtract two data sets
537       
538        :param other: data set to subtract from the current one
539       
540        :return: new data set
541       
542        :raise ValueError: raised when two data sets are incompatible
543       
544        """
545        def operation(a, b):
546            return a - b
547        return self._perform_operation(other, operation)
548       
549    def __rsub__(self, other):
550        """
551        Subtract two data sets
552       
553        :param other: data set to subtract from the current one
554       
555        :return: new data set
556       
557        :raise ValueError: raised when two data sets are incompatible
558       
559        """
560        def operation(a, b):
561            return b - a
562        return self._perform_operation(other, operation)
563       
564    def __mul__(self, other):
565        """
566        Multiply two data sets
567       
568        :param other: data set to subtract from the current one
569       
570        :return: new data set
571       
572        :raise ValueError: raised when two data sets are incompatible
573       
574        """
575        def operation(a, b):
576            return a * b
577        return self._perform_operation(other, operation)
578       
579    def __rmul__(self, other):
580        """
581        Multiply two data sets
582       
583        :param other: data set to subtract from the current one
584       
585        :return: new data set
586       
587        :raise ValueError: raised when two data sets are incompatible
588        """
589        def operation(a, b):
590            return b * a
591        return self._perform_operation(other, operation)
592       
593    def __div__(self, other):
594        """
595        Divided a data set by another
596       
597        :param other: data set that the current one is divided by
598       
599        :return: new data set
600       
601        :raise ValueError: raised when two data sets are incompatible
602       
603        """
604        def operation(a, b):
605            return a/b
606        return self._perform_operation(other, operation)
607       
608    def __rdiv__(self, other):
609        """
610        Divided a data set by another
611       
612        :param other: data set that the current one is divided by
613       
614        :return: new data set
615       
616        :raise ValueError: raised when two data sets are incompatible
617       
618        """
619        def operation(a, b):
620            return b/a
621        return self._perform_operation(other, operation)           
622           
623class Data1D(plottable_1D, DataInfo):
624    """
625    1D data class
626    """
627    x_unit = '1/A'
628    y_unit = '1/cm'
629   
630    def __init__(self, x, y, dx=None, dy=None):
631        DataInfo.__init__(self)
632        plottable_1D.__init__(self, x, y, dx, dy)
633       
634       
635    def __str__(self):
636        """
637        Nice printout
638        """
639        _str =  "%s\n" % DataInfo.__str__(self)
640   
641        _str += "Data:\n"
642        _str += "   Type:         %s\n" % self.__class__.__name__
643        _str += "   X-axis:       %s\t[%s]\n" % (self._xaxis, self._xunit)
644        _str += "   Y-axis:       %s\t[%s]\n" % (self._yaxis, self._yunit)
645        _str += "   Length:       %g\n" % len(self.x)
646
647        return _str
648
649    def is_slit_smeared(self):
650        """
651        Check whether the data has slit smearing information
652       
653        :return: True is slit smearing info is present, False otherwise
654       
655        """
656        def _check(v):           
657            if (v.__class__ == list or v.__class__ == numpy.ndarray) \
658                and len(v) > 0 and min(v) > 0:
659                return True
660           
661            return False
662       
663        return _check(self.dxl) or _check(self.dxw)
664       
665    def clone_without_data(self, length=0, clone=None):
666        """
667        Clone the current object, without copying the data (which
668        will be filled out by a subsequent operation).
669        The data arrays will be initialized to zero.
670       
671        :param length: length of the data array to be initialized
672        :param clone: if provided, the data will be copied to clone
673        """
674        from copy import deepcopy
675       
676        if clone is None or not issubclass(clone.__class__, Data1D):
677            x  = numpy.zeros(length) 
678            dx = numpy.zeros(length) 
679            y  = numpy.zeros(length) 
680            dy = numpy.zeros(length) 
681            clone = Data1D(x, y, dx=dx, dy=dy)
682       
683        clone.title       = self.title
684        clone.run         = self.run
685        clone.filename    = self.filename
686        clone.instrument   = self.instrument
687        clone.notes       = deepcopy(self.notes) 
688        clone.process     = deepcopy(self.process) 
689        clone.detector    = deepcopy(self.detector) 
690        clone.sample      = deepcopy(self.sample) 
691        clone.source      = deepcopy(self.source) 
692        clone.collimation = deepcopy(self.collimation) 
693        clone.meta_data   = deepcopy(self.meta_data) 
694        clone.errors      = deepcopy(self.errors) 
695       
696        return clone
697
698    def _validity_check(self, other):
699        """
700        Checks that the data lengths are compatible.
701        Checks that the x vectors are compatible.
702        Returns errors vectors equal to original
703        errors vectors if they were present or vectors
704        of zeros when none was found.
705       
706        :param other: other data set for operation
707       
708        :return: dy for self, dy for other [numpy arrays]
709       
710        :raise ValueError: when lengths are not compatible
711       
712        """
713        dy_other = None
714        if isinstance(other, Data1D):
715            # Check that data lengths are the same
716            if len(self.x) != len(other.x) or \
717                len(self.y) != len(other.y):
718                msg =  "Unable to perform operation: data length are not equal"
719                raise ValueError, msg
720           
721            # Here we could also extrapolate between data points
722            for i in range(len(self.x)):
723                if self.x[i] != other.x[i]:
724                    msg = "Incompatible data sets: x-values do not match"
725                    raise ValueError, msg
726           
727            # Check that the other data set has errors, otherwise
728            # create zero vector
729            dy_other = other.dy
730            if other.dy == None or (len(other.dy) != len(other.y)):
731                dy_other = numpy.zeros(len(other.y))
732           
733        # Check that we have errors, otherwise create zero vector
734        dy = self.dy
735        if self.dy == None or (len(self.dy) != len(self.y)):
736            dy = numpy.zeros(len(self.y))           
737           
738        return dy, dy_other
739
740    def _perform_operation(self, other, operation):
741        """
742        """
743        # First, check the data compatibility
744        dy, dy_other = self._validity_check(other)
745        result = self.clone_without_data(len(self.x))
746       
747        for i in range(len(self.x)):
748            result.x[i] = self.x[i]
749            if self.dx is not None and len(self.x) == len(self.dx):
750                result.dx[i] = self.dx[i]
751           
752            a = Uncertainty(self.y[i], dy[i]**2)
753            if isinstance(other, Data1D):
754                b = Uncertainty(other.y[i], dy_other[i]**2)
755            else:
756                b = other
757           
758            output = operation(a, b)
759            result.y[i] = output.x
760            result.dy[i] = math.sqrt(math.fabs(output.variance))
761        return result
762       
763class Data2D(plottable_2D, DataInfo):
764    """
765    2D data class
766    """
767    ## Units for Q-values
768    Q_unit = '1/A'
769   
770    ## Units for I(Q) values
771    I_unit = '1/cm'
772   
773    ## Vector of Q-values at the center of each bin in x
774    x_bins = None
775   
776    ## Vector of Q-values at the center of each bin in y
777    y_bins = None
778   
779   
780    def __init__(self, data=None, err_data=None, qx_data=None,
781                  qy_data=None, q_data=None, mask=None, 
782                  dqx_data=None, dqy_data=None):
783        self.y_bins = []
784        self.x_bins = []
785        DataInfo.__init__(self)
786        plottable_2D.__init__(self, data, err_data, qx_data,
787                              qy_data, q_data,mask, dqx_data, dqy_data)
788        if len(self.detector) > 0:
789            raise RuntimeError, "Data2D: Detector bank already filled at init"
790
791    def __str__(self):
792        _str =  "%s\n" % DataInfo.__str__(self)
793       
794        _str += "Data:\n"
795        _str += "   Type:         %s\n" % self.__class__.__name__
796        _str += "   X- & Y-axis:  %s\t[%s]\n" % (self._yaxis, self._yunit)
797        _str += "   Z-axis:       %s\t[%s]\n" % (self._zaxis, self._zunit)
798        #leny = 0
799        #if len(self.data) > 0:
800        #    leny = len(self.data)
801        _str += "   Length:       %g \n" % (len(self.data))
802       
803        return _str
804 
805    def clone_without_data(self, length=0, clone=None):
806        """
807        Clone the current object, without copying the data (which
808        will be filled out by a subsequent operation).
809        The data arrays will be initialized to zero.
810       
811        :param length: length of the data array to be initialized
812        :param clone: if provided, the data will be copied to clone
813        """
814        from copy import deepcopy
815       
816        if clone is None or not issubclass(clone.__class__, Data2D): 
817            data     = numpy.zeros(length) 
818            err_data = numpy.zeros(length) 
819            qx_data = numpy.zeros(length)
820            qy_data = numpy.zeros(length)
821            q_data = numpy.zeros(length)
822            mask = numpy.zeros(length)
823            dqx_data = None
824            dqy_data = None
825            clone = Data2D(data, err_data, qx_data, qy_data,
826                            q_data,mask, dqx_data=dqx_data, dqy_data=dqy_data)
827
828        clone.title       = self.title
829        clone.run         = self.run
830        clone.filename    = self.filename
831        clone.instrument  = self.instrument
832        clone.notes       = deepcopy(self.notes) 
833        clone.process     = deepcopy(self.process) 
834        clone.detector    = deepcopy(self.detector) 
835        clone.sample      = deepcopy(self.sample) 
836        clone.source      = deepcopy(self.source) 
837        clone.collimation = deepcopy(self.collimation) 
838        clone.meta_data   = deepcopy(self.meta_data) 
839        clone.errors      = deepcopy(self.errors) 
840       
841        return clone
842 
843 
844    def _validity_check(self, other):
845        """
846        Checks that the data lengths are compatible.
847        Checks that the x vectors are compatible.
848        Returns errors vectors equal to original
849        errors vectors if they were present or vectors
850        of zeros when none was found.
851       
852        :param other: other data set for operation
853       
854        :return: dy for self, dy for other [numpy arrays]
855       
856        :raise ValueError: when lengths are not compatible
857       
858        """
859        err_other = None
860        if isinstance(other, Data2D):
861            # Check that data lengths are the same
862            if numpy.size(self.data) != numpy.size(other.data):
863                msg = "Unable to perform operation: data length are not equal"
864                raise ValueError, msg
865               
866            # Check that the scales match
867            #TODO: matching scales?     
868           
869            # Check that the other data set has errors, otherwise
870            # create zero vector
871            #TODO: test this
872            err_other = other.err_data
873            if other.err_data == None or \
874                (numpy.size(other.err_data) != numpy.size(other.data)):
875                err_other = numpy.zeros([numpy.size(other.data, 0),
876                                          numpy.size(other.data, 1)])
877           
878        # Check that we have errors, otherwise create zero vector
879        err = self.err_data
880        if self.err_data == None or \
881            (numpy.size(self.err_data) != numpy.size(self.data)):
882            err = numpy.zeros([numpy.size(self.data, 0),
883                                numpy.size(self.data, 1)])
884           
885        return err, err_other
886 
887 
888    def _perform_operation(self, other, operation):
889        """
890        Perform 2D operations between data sets
891       
892        :param other: other data set
893        :param operation: function defining the operation
894       
895        """
896        # First, check the data compatibility
897        dy, dy_other = self._validity_check(other)
898   
899        result = self.clone_without_data([numpy.size(self.data, 0), 
900                                          numpy.size(self.data, 1)])
901       
902        for i in range(numpy.size(self.data, 0)):
903            for j in range(numpy.size(self.data, 1)):
904                result.data[i][j] = self.data[i][j]
905                if self.err_data is not None and \
906                    numpy.size(self.data) == numpy.size(self.err_data):
907                    result.err_data[i][j] = self.err_data[i][j]
908               
909                a = Uncertainty(self.data[i][j], dy[i][j]**2)
910                if isinstance(other, Data2D):
911                    b = Uncertainty(other.data[i][j], dy_other[i][j]**2)
912                else:
913                    b = other
914               
915                output = operation(a, b)
916                result.data[i][j] = output.x
917                result.err_data[i][j] = math.sqrt(math.fabs(output.variance))
918        return result
919   
920 
921 
Note: See TracBrowser for help on using the repository browser.