Changeset 911dbe4 in sasview for src/sas


Ignore:
Timestamp:
Jul 11, 2016 4:36:55 AM (8 years ago)
Author:
lewis
Branches:
master, ESS_GUI, ESS_GUI_Docs, ESS_GUI_batch_fitting, ESS_GUI_bumps_abstraction, ESS_GUI_iss1116, ESS_GUI_iss879, ESS_GUI_iss959, ESS_GUI_opencl, ESS_GUI_ordering, ESS_GUI_sync_sascalc, costrafo411, magnetic_scatt, release-4.1.1, release-4.1.2, release-4.2.2, ticket-1009, ticket-1094-headless, ticket-1242-2d-resolution, ticket-1243, ticket-1249, ticket885, unittest-saveload
Children:
7a219e3e
Parents:
f2bbabf
git-author:
Lewis O'Driscoll <lewis.o'driscoll@…> (07/07/16 11:30:07)
git-committer:
Lewis O'Driscoll <lewis.o'driscoll@…> (07/11/16 04:36:55)
Message:

Add transform functionality

Location:
src/sas
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • src/sas/sascalc/corfunc/corfunc_calculator.py

    r77e9ac6 r911dbe4  
    5555 
    5656 
    57     def __init__(self, data, lowerq, upperq, background=0, scale=1): 
     57    def __init__(self, data=None, lowerq=None, upperq=None, scale=1): 
    5858        """ 
    5959        Initialize the class. 
    6060 
    6161        :param data: Data of the type DataLoader.Data1D 
    62         :param background: Background value. Will be subtracted from the data 
    63             before processing 
     62        :param lowerq: The Q value to use as the boundary for 
     63            Guinier extrapolation 
     64        :param upperq: A tuple of the form (lower, upper). 
     65            Values between lower and upper will be used for Porod extrapolation 
    6466        :param scale: Scaling factor for I(q) 
    6567        """ 
    66         print "Before: {}".format(data.y[0]) 
    67         self._data = self._get_data(data, background, scale) 
    68         print "After: {}".format(self._data.y[0]) 
     68        self._data = None 
     69        self.set_data(data, scale) 
    6970        self.lowerq = lowerq 
    7071        self.upperq = upperq 
     72        self.background = 0 
     73 
     74    def set_data(self, data, scale=1): 
     75        """ 
     76        Prepares the data for analysis 
     77 
     78        :return: new_data = data * scale - background 
     79        """ 
     80        if data is None: 
     81            return 
     82        # Only process data of the class Data1D 
     83        if not issubclass(data.__class__, Data1D): 
     84            raise ValueError, "Data must be of the type DataLoader.Data1D" 
     85 
     86        # Prepare the data 
     87        new_data = Data1D(x=data.x, y=data.y) 
     88        new_data *= scale 
     89 
     90        # Ensure the errors are set correctly 
     91        if new_data.dy is None or len(new_data.x) != len(new_data.dy) or \ 
     92            (min(new_data.dy) == 0 and max(new_data.dy) == 0): 
     93            new_data.dy = np.ones(len(new_data.x)) 
     94 
     95        self._data = new_data 
     96 
     97    def compute_background(self, upperq=None): 
     98        """ 
     99        Compute the background level from the Porod region of the data 
     100        """ 
     101        if self._data is None: return 0 
     102        elif upperq is None and self.upperq is not None: upperq = self.upperq 
     103        elif upperq == 0: return 0 
     104        q = self._data.x 
     105        mask = np.logical_and(q > upperq[0], q < upperq[1]) 
     106        _, _, bg = self._fit_porod(q[mask], self._data.y[mask]) 
     107 
     108        return bg 
    71109 
    72110    def compute_extrapolation(self): 
     
    82120        return extrapolation 
    83121 
    84     def compute_transform(self, extrapolation): 
     122    def compute_transform(self, extrapolation, background=None): 
    85123        """ 
    86124        Transform an extrapolated scattering curve into a correlation function. 
     
    88126        qs = extrapolation.x 
    89127        iqs = extrapolation.y 
     128        q = self._data.x 
     129        if background is None: background = self.background 
    90130 
    91         gamma = dct(iqs*qs**2) 
     131        gamma = dct((iqs-background)*qs**2) 
    92132        gamma = gamma / gamma.max() 
    93133        xs = np.pi*np.arange(len(qs),dtype=np.float32)/(q[1]-q[0])/len(qs) 
     
    97137        return transform 
    98138 
    99     def _porod(self, q, K, sigma): 
     139    def _porod(self, q, K, sigma, bg): 
    100140        """Equation for the Porod region of the data""" 
    101         return (K*q**(-4))*np.exp(-q**2*sigma**2) 
     141        return bg + (K*q**(-4))*np.exp(-q**2*sigma**2) 
    102142 
    103143    def _fit_guinier(self, q, iq): 
     
    106146        return lstsq(A, np.log(iq)) 
    107147 
    108     def _get_data(self, data, background, scale): 
    109         """ 
    110         Prepares the data for analysis 
    111  
    112         :return: new_data = data * scale - background 
    113         """ 
    114         # Only process data of the class Data1D 
    115         if not issubclass(data.__class__, Data1D): 
    116             raise ValueError, "Data must be of the type DataLoader.Data1D" 
    117  
    118         # Prepare the data 
    119         new_data = (scale * data) 
    120         new_data.y -= background 
    121  
    122         # Check the vector lengths are equal 
    123         assert len(new_data.x) == len(new_data.y) 
    124  
    125         # Ensure the errors are set correctly 
    126         if new_data.dy is None or len(new_data.x) != len(new_data.dy) or \ 
    127             (min(new_data.dy) == 0 and max(new_data.dy) == 0): 
    128             new_data.dy = np.ones(len(new_data.x)) 
    129  
    130         return new_data 
     148    def _fit_porod(self, q, iq): 
     149        """Fit the Porod region of the curve""" 
     150        fitp = curve_fit(lambda q, k, sig, bg: self._porod(q, k, sig, bg)*q**2, 
     151                         q, iq*q**2)[0] 
     152        k, sigma, bg = fitp 
     153        return k, sigma, bg 
    131154 
    132155    def _fit_data(self, q, iq): 
     
    137160        # Returns an array where the 1st and 2nd elements are the values of k 
    138161        # and sigma for the best-fit Porod function 
    139         fitp = curve_fit(lambda q, k, sig: self._porod(q, k, sig)*q**2, 
    140                          q[mask], iq[mask]*q[mask]**2)[0] 
     162        k, sigma, bg = self._fit_porod(q[mask], iq[mask]) 
     163        self.background = bg 
    141164 
    142165        # Smooths between the best-fit porod function and the data to produce a 
     
    144167        data = interp1d(q, iq) 
    145168        s1 = self._Interpolator(data, 
    146             lambda x: self._porod(x, fitp[0], fitp[1]), self.upperq[0], q[-1]) 
     169            lambda x: self._porod(x, k, sigma, bg), self.upperq[0], q[-1]) 
    147170 
    148171        mask = np.logical_and(q < self.lowerq, 0 < q) 
  • src/sas/sasgui/perspectives/corfunc/corfunc.py

    r3b8efec r911dbe4  
    2222IQ_DATA_LABEL = r"$I_{obs}(q)$" 
    2323IQ_EXTRAPOLATED_DATA_LABEL = r"$I_{extrap}(q)$" 
     24 
     25GROUP_ID_TRANSFORM = r"$\Gamma(x)$" 
     26TRANSFORM_LABEL = r"$\Gamma(x)$" 
    2427 
    2528 
     
    129132        new_plot = Data1D(data.x, data.y, dy=data.dy) 
    130133 
     134        group_id = "" 
    131135        if label == IQ_DATA_LABEL or label == IQ_EXTRAPOLATED_DATA_LABEL: 
    132136            new_plot.xaxis("\\rm{Q}", 'A^{-1}') 
    133137            new_plot.yaxis("\\rm{Intensity} ", "cm^{-1}") 
    134  
     138            group_id = GROUP_ID_IQ_DATA 
     139        elif label == TRANSFORM_LABEL: 
     140            new_plot.xaxis("{x}", 'A') 
     141            new_plot.yaxis("{\\Gamma}", '') 
     142            group_id = GROUP_ID_TRANSFORM 
    135143        new_plot.symbol = GUIFRAME_ID.CURVE_SYMBOL_NUM 
    136144        new_plot.id = label 
    137145        new_plot.name = label 
     146        new_plot.group_id = group_id 
    138147        new_plot.interactive = True 
    139         new_plot.group_id = GROUP_ID_IQ_DATA 
    140         new_plot.title = "I(q)" 
     148        new_plot.title = group_id.replace('$', '').replace('\\', '') 
    141149        # Show data on a linear scale 
    142150        new_plot.xtransform = 'x' 
    143151        new_plot.ytransform = 'y' 
    144152        wx.PostEvent(self.parent, 
    145                      NewPlotEvent(plot=new_plot, title="I(q)", reset=reset)) 
     153                     NewPlotEvent(plot=new_plot, title=new_plot.title, 
     154                        reset=reset)) 
  • src/sas/sasgui/perspectives/corfunc/corfunc_panel.py

    rf2bbabf r911dbe4  
    3838        self.SetWindowVariant(variant=FONT_VARIANT) 
    3939        self._manager = manager 
    40         self._data = data # The data to be analysed 
     40        # The data with no correction for background values 
     41        self._data = data # The data to be analysed (corrected fr background) 
    4142        self._extrapolated_data = None # The extrapolated data set 
     43        self._calculator = CorfuncCalculator() 
    4244        self._data_name_box = None # Text box to show name of file 
    4345        self._background_input = None 
     
    4547        self._qmax1_input = None 
    4648        self._qmax2_input = None 
     49        self._transform_btn = None 
     50        self._compute_btn = None 
    4751        self.qmin = 0 
    4852        self.qmax = (0, 0) 
     
    107111        :param data: The data that has been loaded 
    108112        """ 
     113        self._transform_btn.Disable() 
    109114        if data is None: 
    110115            return 
    111116        self._data_name_box.SetValue(str(data.title)) 
    112117        self._data = data 
     118        self._calculator.set_data(data) 
    113119        if self._manager is not None: 
    114120            from sas.sasgui.perspectives.corfunc.corfunc import IQ_DATA_LABEL 
    115             self._manager.show_data(data, IQ_DATA_LABEL, reset=True) 
     121            self._manager.show_data(self._data, IQ_DATA_LABEL, reset=True) 
    116122        if set_qrange: 
    117123            lower = data.x[-1]*0.05 
     
    120126            self.set_qmin(lower) 
    121127            self.set_qmax((upper1, upper2)) 
    122             self.set_background(0.0) 
     128            self.set_background(self._calculator.compute_background(self.qmax)) 
    123129 
    124130    def get_data(self): 
     
    134140            wx.PostEvent(self.parent.parent, StatusEvent(status=msg)) 
    135141            return 
    136         calculator = CorfuncCalculator(self._data, self.qmin, self.qmax, 
    137             background=self.background) 
    138         self._extrapolated_data = calculator.compute_extrapolation() 
     142        self._calculator.set_data(self._data) 
     143        self._calculator.lowerq = self.qmin 
     144        self._calculator.upperq = self.qmax 
     145        self._extrapolated_data = self._calculator.compute_extrapolation() 
    139146        # TODO: Find way to set xlim and ylim so full range of data can be 
    140147        # plotted 
     
    148155            IQ_EXTRAPOLATED_DATA_LABEL 
    149156        self._manager.show_data(to_plot, IQ_EXTRAPOLATED_DATA_LABEL) 
     157        self._transform_btn.Enable() 
     158 
     159    def compute_transform(self, event=None): 
     160        transformed_data = self._calculator.compute_transform( 
     161            self._extrapolated_data, self.background) 
     162        from sas.sasgui.perspectives.corfunc.corfunc import TRANSFORM_LABEL 
     163        import numpy as np 
     164        plot_x = transformed_data.x[np.where(transformed_data.x <= 200)] 
     165        plot_y = transformed_data.y[np.where(transformed_data.x <= 200)] 
     166        self._manager.show_data(Data1D(plot_x, plot_y), TRANSFORM_LABEL) 
    150167 
    151168 
     
    183200        self._background_input.SetValue(str(bg)) 
    184201 
     202    def _compute_background(self, event=None): 
     203        self.set_background(self._calculator.compute_background(self.qmax)) 
    185204 
    186205    def _on_enter_input(self, event=None): 
     
    200219        group_id = GROUP_ID_IQ_DATA 
    201220        if event is not None: 
     221            if self._manager is not None and\ 
     222                event.GetEventObject() == self._background_input: 
     223                from sas.sasgui.perspectives.corfunc.corfunc\ 
     224                    import IQ_DATA_LABEL 
     225                self._manager.show_data(self._data, IQ_DATA_LABEL, reset=True) 
    202226            wx.PostEvent(self._manager.parent, PlotQrangeEvent( 
    203227                ctrl=[self._qmin_input, self._qmax1_input, self._qmax2_input], 
     
    309333            wx.RIGHT | wx.EXPAND, 5) 
    310334 
     335        background_button = wx.Button(self, wx.NewId(), "Calculate", 
     336            size=(75, 25)) 
     337        background_button.Bind(wx.EVT_BUTTON, self._compute_background) 
     338        q_sizer.Add(background_button, (1, 2), (1, 1), wx.RIGHT | wx.EXPAND, 5) 
     339 
    311340        qrange_label = wx.StaticText(self, -1, "Q Range:", size=(50,20)) 
    312341        q_sizer.Add(qrange_label, (2,0), (1,1), wx.LEFT | wx.EXPAND, 5) 
     
    395424 
    396425        extrapolate_btn = wx.Button(self, wx.NewId(), "Extrapolate") 
    397         transform_btn = wx.Button(self, wx.NewId(), "Transform") 
    398         compute_btn = wx.Button(self, wx.NewId(), "Compute Measuments") 
     426        self._transform_btn = wx.Button(self, wx.NewId(), "Transform") 
     427        self._compute_btn = wx.Button(self, wx.NewId(), "Compute Measuments") 
     428 
     429        self._transform_btn.Disable() 
     430        self._compute_btn.Disable() 
    399431 
    400432        extrapolate_btn.Bind(wx.EVT_BUTTON, self.compute_extrapolation) 
     433        self._transform_btn.Bind(wx.EVT_BUTTON, self.compute_transform) 
    401434 
    402435        controls_sizer.Add(extrapolate_btn, wx.CENTER | wx.EXPAND) 
    403         controls_sizer.Add(transform_btn, wx.CENTER | wx.EXPAND) 
    404         controls_sizer.Add(compute_btn, wx.CENTER | wx.EXPAND) 
     436        controls_sizer.Add(self._transform_btn, wx.CENTER | wx.EXPAND) 
     437        controls_sizer.Add(self._compute_btn, wx.CENTER | wx.EXPAND) 
    405438 
    406439        controlbox_sizer.Add(controls_sizer, wx.TOP | wx.EXPAND, 0) 
Note: See TracChangeset for help on using the changeset viewer.