-                      r240a2d2
+                      rf8aa738
 import logging
 import sys
-import sas.models.sas_extension.smearer as smearer
-from sas.sascalc.data_util.smearing_2d import Smearer2D
+def smear_selection(data1D, model = None):
+from sasmodels.resolution import Slit1D, Pinhole1D
+from sasmodels.resolution2d import Pinhole2D
+def smear_selection(data, model = None):
     """
     Creates the right type of smearer according
+    Creates the right type of smearer according
     to the data.
     The canSAS format has a rule that either
     slit smearing data OR resolution smearing data
     is available.
+    is available.
     For the present purpose, we choose the one that
     has none-zero data. If both slit and resolution
     smearing arrays are filled with good data
+    smearing arrays are filled with good data
     (which should not happen), then we choose the
     resolution smearing data.
     :param data1D: Data1D object
+    resolution smearing data.
+    :param data: Data1D object
     :param model: sas.model instance
     """
     # Sanity check. If we are not dealing with a SAS Data1D
     # object, just return None
     if  data1D.__class__.__name__ not in ['Data1D', 'Theory1D']:
         if data1D == None:
+    if  data.__class__.__name__ not in ['Data1D', 'Theory1D']:
+        if data == None:
             return None
         elif data1D.dqx_data == None or data1D.dqy_data == None:
+        elif data.dqx_data == None or data.dqy_data == None:
             return None
         return Smearer2D(data1D)
     if  not hasattr(data1D, "dx") and not hasattr(data1D, "dxl")\
          and not hasattr(data1D, "dxw"):
+        return Pinhole2D(data)
+    if  not hasattr(data, "dx") and not hasattr(data, "dxl")\
+         and not hasattr(data, "dxw"):
         return None
     # Look for resolution smearing data
     _found_resolution = False
     if data1D.dx is not None and len(data1D.dx) == len(data1D.x):
+    if data.dx is not None and len(data.dx) == len(data.x):
         # Check that we have non-zero data
         if data1D.dx[0] > 0.0:
+        if data.dx[0] > 0.0:
             _found_resolution = True
             #print "_found_resolution",_found_resolution
 …
     # If we found resolution smearing data, return a QSmearer
     if _found_resolution == True:
+        return QSmearer(data1D, model)
+        #return pinhole_smear(data1D, model)
+         return pinhole_smear(data, model)
     # Look for slit smearing data
     _found_slit = False
     if data1D.dxl is not None and len(data1D.dxl) == len(data1D.x) \
         and data1D.dxw is not None and len(data1D.dxw) == len(data1D.x):
+    if data.dxl is not None and len(data.dxl) == len(data.x) \
+        and data.dxw is not None and len(data.dxw) == len(data.x):
         # Check that we have non-zero data
         if data1D.dxl[0] > 0.0 or data1D.dxw[0] > 0.0:
+        if data.dxl[0] > 0.0 or data.dxw[0] > 0.0:
             _found_slit = True
         # Sanity check: all data should be the same as a function of Q
         for item in data1D.dxl:
             if data1D.dxl[0] != item:
+        for item in data.dxl:
+            if data.dxl[0] != item:
                 _found_resolution = False
                 break
         for item in data1D.dxw:
             if data1D.dxw[0] != item:
+        for item in data.dxw:
+            if data.dxw[0] != item:
                 _found_resolution = False
                 break
     # If we found slit smearing data, return a slit smearer
     if _found_slit == True:
+        #return SlitSmearer(data1D, model)
+        return slit_smear(data1D, model)
+        return slit_smear(data, model)
     return None
-class _BaseSmearer(object):
-    """
-        Base class for smearers
-    """
-    def __init__(self):
-        self.nbins = 0
-        self.nbins_low = 0
-        self.nbins_high = 0
-        self._weights = None
-        ## Internal flag to keep track of C++ smearer initialization
-        self._init_complete = False
-        self._smearer = None
-        self.model = None
-        self.min = None
-        self.max = None
-        self.qvalues = []
-    def __deepcopy__(self, memo=None):
-        """
-        Return a valid copy of self.
-        Avoid copying the _smearer C object and force a matrix recompute
-        when the copy is used.
-        """
-        result = _BaseSmearer()
-        result.nbins = self.nbins
-        return result
-    def _compute_matrix(self):
-        """
-            Place holder for matrix computation
-        """
-        return NotImplemented
-    def get_unsmeared_range(self, q_min=None, q_max=None):
-        """
-            Place holder for method returning unsmeared range
-        """
-        return NotImplemented
-    def get_bin_range(self, q_min=None, q_max=None):
-        """
-        :param q_min: minimum q-value to smear
-        :param q_max: maximum q-value to smear
-        """
-        # If this is the first time we call for smearing,
-        # initialize the C++ smearer object first
-        if not self._init_complete:
-            self._initialize_smearer()
-        if q_min == None:
-            q_min = self.min
-        if q_max == None:
-            q_max = self.max
-        _qmin_unsmeared, _qmax_unsmeared = self.get_unsmeared_range(q_min,
-                                                                     q_max)
-        _first_bin = None
-        _last_bin  = None
-        #step = (self.max - self.min) / (self.nbins - 1.0)
-        # Find the first and last bin number in all extrapolated and real data
-        try:
-            for i in range(self.nbins):
-                q_i = smearer.get_q(self._smearer, i)
-                if (q_i >= _qmin_unsmeared) and (q_i <= _qmax_unsmeared):
-                    # Identify first and last bin
-                    if _first_bin is None:
-                        _first_bin = i
-                    else:
-                        _last_bin  = i
-        except:
-            msg = "_BaseSmearer.get_bin_range: "
-            msg += " error getting range\n  %s" % sys.exc_value
-            raise RuntimeError, msg
-        #  Find the first and last bin number only in the real data
-        _first_bin, _last_bin = self._get_unextrapolated_bin( \
-                                        _first_bin, _last_bin)
-        return _first_bin, _last_bin
-    def __call__(self, iq_in, first_bin = 0, last_bin = None):
-        """
-        Perform smearing
-        """
-        # If this is the first time we call for smearing,
-        # initialize the C++ smearer object first
-        if not self._init_complete:
-            self._initialize_smearer()
-        if last_bin is None or last_bin >= len(iq_in):
-            last_bin = len(iq_in) - 1
-        # Check that the first bin is positive
-        if first_bin < 0:
-            first_bin = 0
-        # With a model given, compute I for the extrapolated points and append
-        # to the iq_in
-        iq_in_temp = iq_in
-        if self.model != None:
-            temp_first, temp_last = self._get_extrapolated_bin( \
-                                                        first_bin, last_bin)
-            if self.nbins_low > 0:
-                iq_in_low = self.model.evalDistribution( \
-                                    numpy.fabs(self.qvalues[0:self.nbins_low]))
-            iq_in_high = self.model.evalDistribution( \
-                                            self.qvalues[(len(self.qvalues) - \
-                                            self.nbins_high - 1):])
-            # Todo: find out who is sending iq[last_poin] = 0.
-            if iq_in[len(iq_in) - 1] == 0:
-                iq_in[len(iq_in) - 1] = iq_in_high[0]
-            # Append the extrapolated points to the data points
-            if self.nbins_low > 0:
-                iq_in_temp = numpy.append(iq_in_low, iq_in)
-            if self.nbins_high > 0:
-                iq_in_temp = numpy.append(iq_in_temp, iq_in_high[1:])
-        else:
-            temp_first = first_bin
-            temp_last = last_bin
-            #iq_in_temp = iq_in
-        # Sanity check
-        if len(iq_in_temp) != self.nbins:
-            msg = "Invalid I(q) vector: inconsistent array "
-            msg += " length %d != %s" % (len(iq_in_temp), str(self.nbins))
-            raise RuntimeError, msg
-        # Storage for smeared I(q)
-        iq_out = numpy.zeros(self.nbins)
-        smear_output = smearer.smear(self._smearer, iq_in_temp, iq_out,
-                                      #0, self.nbins - 1)
-                                      temp_first, temp_last)
-                                      #first_bin, last_bin)
-        if smear_output < 0:
-            msg = "_BaseSmearer: could not smear, code = %g" % smear_output
-            raise RuntimeError, msg
-        temp_first = first_bin + self.nbins_low
-        temp_last = self.nbins - self.nbins_high
-        out = iq_out[temp_first: temp_last]
-        return out
-    def _initialize_smearer(self):
-        """
-            Place holder for initializing data smearer
-        """
-        return NotImplemented
-    def _get_unextrapolated_bin(self, first_bin = 0, last_bin = 0):
-        """
-        Get unextrapolated first bin and the last bin
-        : param first_bin: extrapolated first_bin
-        : param last_bin: extrapolated last_bin
-        : return fist_bin, last_bin: unextrapolated first and last bin
-        """
-        #  For first bin
-        if first_bin <= self.nbins_low:
-            first_bin = 0
-        else:
-            first_bin = first_bin - self.nbins_low
-        # For last bin
-        if last_bin >= (self.nbins - self.nbins_high):
-            last_bin  = self.nbins - (self.nbins_high + self.nbins_low + 1)
-        elif last_bin >= self.nbins_low:
-            last_bin = last_bin - self.nbins_low
-        else:
-            last_bin = 0
-        return first_bin, last_bin
-    def _get_extrapolated_bin(self, first_bin = 0, last_bin = 0):
-        """
-        Get extrapolated first bin and the last bin
-        : param first_bin: unextrapolated first_bin
-        : param last_bin: unextrapolated last_bin
-        : return first_bin, last_bin: extrapolated first and last bin
-        """
-        #  For the first bin
-        # In the case that needs low extrapolation data
-        first_bin = 0
-        # For last bin
-        if last_bin >= self.nbins - (self.nbins_high + self.nbins_low + 1):
-            # In the case that needs higher q extrapolation data
-            last_bin = self.nbins - 1
-        else:
-            # In the case that doesn't need higher q extrapolation data
-            last_bin += self.nbins_low
-        return first_bin, last_bin
-class _SlitSmearer(_BaseSmearer):
-    """
-    Slit smearing for I(q) array
-    """
-    def __init__(self, nbins=None, width=None, height=None, min=None, max=None):
-        """
-        Initialization
-        :param iq: I(q) array [cm-1]
-        :param width: slit width [A-1]
-        :param height: slit height [A-1]
-        :param min: Q_min [A-1]
-        :param max: Q_max [A-1]
-        """
-        _BaseSmearer.__init__(self)
-        ## Slit width in Q units
-        self.width  = width
-        ## Slit height in Q units
-        self.height = height
-        ## Q_min (Min Q-value for I(q))
-        self.min    = min
-        ## Q_max (Max Q_value for I(q))
-        self.max    = max
-        ## Number of Q bins
-        self.nbins  = nbins
-        ## Number of points used in the smearing computation
-        self.npts   = 3000
-        ## Smearing matrix
-        self._weights = None
-        self.qvalues  = None
-    def _initialize_smearer(self):
-        """
-        Initialize the C++ smearer object.
-        This method HAS to be called before smearing
-        """
-        #self._smearer = smearer.new_slit_smearer(self.width,
-        # self.height, self.min, self.max, self.nbins)
-        self._smearer = smearer.new_slit_smearer_with_q(self.width,
-                                                    self.height, self.qvalues)
-        self._init_complete = True
-    def get_unsmeared_range(self, q_min, q_max):
-        """
-        Determine the range needed in unsmeared-Q to cover
-        the smeared Q range
-        """
-        # Range used for input to smearing
-        _qmin_unsmeared = q_min
-        _qmax_unsmeared = q_max
-        try:
-            _qmin_unsmeared = self.min
-            _qmax_unsmeared = self.max
-        except:
-            logging.error("_SlitSmearer.get_bin_range: %s" % sys.exc_value)
-        return _qmin_unsmeared, _qmax_unsmeared
-class SlitSmearer(_SlitSmearer):
-    """
-    Adaptor for slit smearing class and SAS data
-    """
-    def __init__(self, data1D, model = None):
-        """
-        Assumption: equally spaced bins of increasing q-values.
-        :param data1D: data used to set the smearing parameters
-        """
-        # Initialization from parent class
-        super(SlitSmearer, self).__init__()
-        ## Slit width
-        self.width = 0
-        self.nbins_low = 0
-        self.nbins_high = 0
-        self.model = model
-        if data1D.dxw is not None and len(data1D.dxw) == len(data1D.x):
-            self.width = data1D.dxw[0]
-            # Sanity check
-            for value in data1D.dxw:
-                if value != self.width:
-                    msg = "Slit smearing parameters must "
-                    msg += " be the same for all data"
-                    raise RuntimeError, msg
-        ## Slit height
-        self.height = 0
-        if data1D.dxl is not None and len(data1D.dxl) == len(data1D.x):
-            self.height = data1D.dxl[0]
-            # Sanity check
-            for value in data1D.dxl:
-                if value != self.height:
-                    msg = "Slit smearing parameters must be"
-                    msg += " the same for all data"
-                    raise RuntimeError, msg
-        # If a model is given, get the q extrapolation
-        if self.model == None:
-            data1d_x = data1D.x
-        else:
-            # Take larger sigma
-            if self.height > self.width:
-                # The denominator (2.0) covers all the possible w^2 + h^2 range
-                sigma_in = data1D.dxl / 2.0
-            elif self.width > 0:
-                sigma_in = data1D.dxw / 2.0
-            else:
-                sigma_in = []
-            self.nbins_low, self.nbins_high, _, data1d_x = \
-                                get_qextrapolate(sigma_in, data1D.x)
-        ## Number of Q bins
-        self.nbins = len(data1d_x)
-        ## Minimum Q
-        self.min = min(data1d_x)
-        ## Maximum
-        self.max = max(data1d_x)
-        ## Q-values
-        self.qvalues = data1d_x
-class _QSmearer(_BaseSmearer):
-    """
-    Perform Gaussian Q smearing
-    """
-    def __init__(self, nbins=None, width=None, min=None, max=None):
-        """
-        Initialization
-        :param nbins: number of Q bins
-        :param width: array standard deviation in Q [A-1]
-        :param min: Q_min [A-1]
-        :param max: Q_max [A-1]
-        """
-        _BaseSmearer.__init__(self)
-        ## Standard deviation in Q [A-1]
-        self.width = width
-        ## Q_min (Min Q-value for I(q))
-        self.min = min
-        ## Q_max (Max Q_value for I(q))
-        self.max = max
-        ## Number of Q bins
-        self.nbins = nbins
-        ## Smearing matrix
-        self._weights = None
-        self.qvalues  = None
-    def _initialize_smearer(self):
-        """
-        Initialize the C++ smearer object.
-        This method HAS to be called before smearing
-        """
-        #self._smearer = smearer.new_q_smearer(numpy.asarray(self.width),
-        # self.min, self.max, self.nbins)
-        self._smearer = smearer.new_q_smearer_with_q(numpy.asarray(self.width),
-                                                      self.qvalues)
-        self._init_complete = True
-    def get_unsmeared_range(self, q_min, q_max):
-        """
-        Determine the range needed in unsmeared-Q to cover
-        the smeared Q range
-        Take 3 sigmas as the offset between smeared and unsmeared space
-        """
-        # Range used for input to smearing
-        _qmin_unsmeared = q_min
-        _qmax_unsmeared = q_max
-        try:
-            offset = 3.0 * max(self.width)
-            _qmin_unsmeared = self.min#max([self.min, q_min - offset])
-            _qmax_unsmeared = self.max#min([self.max, q_max + offset])
-        except:
-            logging.error("_QSmearer.get_bin_range: %s" % sys.exc_value)
-        return _qmin_unsmeared, _qmax_unsmeared
-class QSmearer(_QSmearer):
-    """
-    Adaptor for Gaussian Q smearing class and SAS data
-    """
-    def __init__(self, data1D, model = None):
-        """
-        Assumption: equally spaced bins of increasing q-values.
-        :param data1D: data used to set the smearing parameters
-        """
-        # Initialization from parent class
-        super(QSmearer, self).__init__()
-        data1d_x = []
-        self.nbins_low = 0
-        self.nbins_high = 0
-        self.model = model
-        ## Resolution
-        #self.width = numpy.zeros(len(data1D.x))
-        if data1D.dx is not None and len(data1D.dx) == len(data1D.x):
-            self.width = data1D.dx
-        if self.model == None:
-            data1d_x = data1D.x
-        else:
-            self.nbins_low, self.nbins_high, self.width, data1d_x = \
-                                get_qextrapolate(self.width, data1D.x)
-        ## Number of Q bins
-        self.nbins = len(data1d_x)
-        ## Minimum Q
-        self.min = min(data1d_x)
-        ## Maximum
-        self.max = max(data1d_x)
-        ## Q-values
-        self.qvalues = data1d_x
-def get_qextrapolate(width, data_x):
-    """
-    Make fake data_x points extrapolated outside of the data_x points
-    :param width: array of std of q resolution
-    :param Data1D.x: Data1D.x array
-    :return new_width, data_x_ext: extrapolated width array and x array
-    :assumption1: data_x is ordered from lower q to higher q
-    :assumption2: len(data) = len(width)
-    :assumption3: the distance between the data points is more compact than the size of width
-    :Todo1: Make sure that the assumptions are correct for Data1D
-    :Todo2: This fixes the edge problem in Qsmearer but still needs to make smearer interface
-    """
-    # Length of the width
-    length = len(width)
-    width_low = math.fabs(width[0])
-    width_high = math.fabs(width[length -1])
-    nbins_low = 0.0
-    nbins_high = 0.0
-    # Compare width(dQ) to the data bin size and take smaller one as the bin
-    # size of the extrapolation; this will correct some weird behavior
-    # at the edge: This method was out (commented)
-    # because it becomes very expansive when
-    # bin size is very small comparing to the width.
-    # Now on, we will just give the bin size of the extrapolated points
-    # based on the width.
-    # Find bin sizes
-    #bin_size_low = math.fabs(data_x[1] - data_x[0])
-    #bin_size_high = math.fabs(data_x[length - 1] - data_x[length - 2])
-    # Let's set the bin size 1/3 of the width(sigma), it is good as long as
-    # the scattering is monotonous.
-    #if width_low < (bin_size_low):
-    bin_size_low = width_low / 10.0
-    #if width_high < (bin_size_high):
-    bin_size_high = width_high / 10.0
-    # Number of q points required below the 1st data point in order to extend
-    # them 3 times of the width (std)
-    if width_low > 0.0:
-        nbins_low = math.ceil(3.0 * width_low / bin_size_low)
-    # Number of q points required above the last data point
-    if width_high > 0.0:
-        nbins_high = math.ceil(3.0 * width_high / bin_size_high)
-    # Make null q points
-    extra_low = numpy.zeros(nbins_low)
-    extra_high = numpy.zeros(nbins_high)
-    # Give extrapolated values
-    ind = 0
-    qvalue = data_x[0] - bin_size_low
-    #if qvalue > 0:
-    while(ind < nbins_low):
-        extra_low[nbins_low - (ind + 1)] = qvalue
-        qvalue -= bin_size_low
-        ind += 1
-        #if qvalue <= 0:
-        #    break
-    # Redefine nbins_low
-    nbins_low = ind
-    # Reset ind for another extrapolation
-    ind = 0
-    qvalue = data_x[length -1] + bin_size_high
-    while(ind < nbins_high):
-        extra_high[ind] = qvalue
-        qvalue += bin_size_high
-        ind += 1
-    # Make a new qx array
-    if nbins_low > 0:
-        data_x_ext = numpy.append(extra_low, data_x)
-    else:
-        data_x_ext = data_x
-    data_x_ext = numpy.append(data_x_ext, extra_high)
-    # Redefine extra_low and high based on corrected nbins
-    # And note that it is not necessary for extra_width to be a non-zero
-    if nbins_low > 0:
-        extra_low = numpy.zeros(nbins_low)
-    extra_high = numpy.zeros(nbins_high)
-    # Make new width array
-    new_width = numpy.append(extra_low, width)
-    new_width = numpy.append(new_width, extra_high)
-    # nbins corrections due to the negative q value
-    nbins_low = nbins_low - len(data_x_ext[data_x_ext <= 0])
-    return nbins_low, nbins_high, \
-             new_width[data_x_ext > 0], data_x_ext[data_x_ext > 0]
-from .resolution import Slit1D, Pinhole1D
 class PySmear(object):
     """
 …
         """
         Apply the resolution function to the data.
         Note that this is called with iq_in matching data.x, but with
         iq_in[first_bin:last_bin] set to theory values for these bins,
         and the remainder left undefined.  The first_bin, last_bin values
         should be those returned from get_bin_range.
         The returned value is of the same length as iq_in, with the range
         first_bin:last_bin set to the resolution smeared values.
 …
         """
         For a given q_min, q_max, find the corresponding indices in the data.
         Returns first, last.
         Note that these are indexes into q from the data, not the q_calc
         needed by the resolution function.  Note also that these are the

Note: See TracChangeset for help on using the changeset viewer.

SasView

Changeset f8aa738 in sasview for src/sas/sascalc/data_util/qsmearing.py

Legend:

src/sas/sascalc/data_util/qsmearing.py

Download in other formats: