Changeset 6053bfd in sasview


Ignore:
Timestamp:
Jul 8, 2016 7:08:15 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:
cdd8910
Parents:
f4622db
Message:

Implement extract parameters

Location:
src/sas
Files:
2 edited

Legend:

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

    r62c0c64 r6053bfd  
    1313class CorfuncCalculator(object): 
    1414 
    15     # Helper class 
    16     class _Struct: 
    17         def __init__(self, **entries): 
    18             self.__dict__.update(entries) 
    19  
    2015    class _Interpolator(object): 
    2116        """ 
     
    147142        return transform 
    148143 
     144    def extract_parameters(self, transformed_data): 
     145        """ 
     146        Extract the interesting measurements from a correlation function 
     147        :param transformed_data: Fourier transformation of the 
     148            extrapolated data 
     149        """ 
     150        # Calculate indexes of maxima and minima 
     151        x = transformed_data.x 
     152        y = transformed_data.y 
     153        maxs = argrelextrema(y, np.greater)[0] 
     154        mins = argrelextrema(y, np.less)[0] 
     155 
     156        # If there are no maxima, return None 
     157        if len(maxs) == 0: 
     158            return None 
     159 
     160        GammaMin = y[mins[0]]  # The value at the first minimum 
     161 
     162        ddy = (y[:-2]+y[2:]-2*y[1:-1])/(x[2:]-x[:-2])**2  # 2nd derivative of y 
     163        dy = (y[2:]-y[:-2])/(x[2:]-x[:-2])  # 1st derivative of y 
     164        # Find where the second derivative goes to zero 
     165        zeros = argrelextrema(np.abs(ddy), np.less)[0] 
     166        # locate the first inflection point 
     167        linear_point = zeros[0] 
     168        linear_point = int(mins[0]/10) 
     169 
     170        # Try to calculate slope around linear_point using 80 data points 
     171        lower = linear_point - 40 
     172        upper = linear_point + 40 
     173 
     174        # If too few data points to the left, use linear_point*2 data points 
     175        if lower < 0: 
     176            lower = 0 
     177            upper = linear_point * 2 
     178        # If too few to right, use 2*(dy.size - linear_point) data points 
     179        elif upper > len(dy): 
     180            upper = len(dy) 
     181            width = len(dy) - linear_point 
     182            lower = 2*linear_point - dy.size 
     183 
     184        m = np.mean(dy[lower:upper])  # Linear slope 
     185        b = y[1:-1][linear_point]-m*x[1:-1][linear_point]  # Linear intercept 
     186 
     187        Lc = (GammaMin-b)/m  # Hard block thickness 
     188 
     189        # Find the data points where the graph is linear to within 1% 
     190        mask = np.where(np.abs((y-(m*x+b))/y) < 0.01)[0] 
     191        if len(mask) == 0:  # Return garbage for bad fits 
     192            return garbage 
     193        dtr = x[mask[0]]  # Beginning of Linear Section 
     194        d0 = x[mask[-1]]  # End of Linear Section 
     195        GammaMax = y[mask[-1]] 
     196        A = -GammaMin/GammaMax  # Normalized depth of minimum 
     197 
     198        params = { 
     199            'max': x[maxs[0]], 
     200            'dtr': dtr, 
     201            'Lc': Lc, 
     202            'd0': d0, 
     203            'A': A, 
     204            'Lc/max': Lc/x[maxs[0]] 
     205        } 
     206 
     207        for key, val in params.iteritems(): 
     208            params[key] = self._round_sig_figs(val, 6) 
     209 
     210        return params 
     211 
     212 
    149213    def _porod(self, q, K, sigma, bg): 
    150214        """Equation for the Porod region of the data""" 
     
    191255 
    192256        return s2 
     257 
     258    def _round_sig_figs(self, x, sigfigs): 
     259        """ 
     260        Round a number to a given number of significant figures. 
     261 
     262        :param x: The value to round 
     263        :param sigfigs: How many significant figures to round to 
     264        :return rounded_str: x rounded to the given number of significant 
     265            figures, as a string 
     266        """ 
     267        # Index of first significant digit 
     268        significant_digit = -int(np.floor(np.log10(np.abs(x)))) 
     269        # Number of digits required for correct number of sig figs 
     270        digits = significant_digit + (sigfigs - 1) 
     271        rounded = np.round(x, decimals=digits) 
     272        rounded_str = "{1:.{0}f}".format(sigfigs -1  + significant_digit, 
     273            rounded) 
     274        return rounded_str 
  • src/sas/sasgui/perspectives/corfunc/corfunc_panel.py

    rf4622db r6053bfd  
    1414from sas.sascalc.corfunc.corfunc_calculator import CorfuncCalculator 
    1515 
     16OUTPUT_STRINGS = { 
     17    'max': "Long Period (A): ", 
     18    'Lc': "Average Hard Block Thickness (A): ", 
     19    'dtr': "Average Interface Thickness (A): ", 
     20    'd0': "Average Core Thickness: ", 
     21    'A': "PolyDispersity: ", 
     22    'Lc/max': "Filling Fraction: " 
     23} 
     24 
    1625if sys.platform.count("win32") > 0: 
    1726    _STATICBOX_WIDTH = 350 
     
    4150        self._data = data # The data to be analysed (corrected fr background) 
    4251        self._extrapolated_data = None # The extrapolated data set 
     52        self._transformed_data = None # Fourier trans. of the extrapolated data 
    4353        self._calculator = CorfuncCalculator() 
    4454        self._data_name_box = None # Text box to show name of file 
     
    4858        self._qmax2_input = None 
    4959        self._transform_btn = None 
    50         self._compute_btn = None 
     60        self._extract_btn = None 
    5161        self.qmin = 0 
    5262        self.qmax = (0, 0) 
    5363        self.background = 0 
    54         # Dictionary for saving IDs of text boxes used to display output data 
    55         self._output_ids = None 
     64        # Dictionary for saving refs to text boxes used to display output data 
     65        self._output_boxes = None 
    5666        self.state = None 
    5767        self._do_layout() 
     
    119129        self._enable_inputs() 
    120130        self._transform_btn.Disable() 
     131        self._extract_btn.Disable() 
    121132        self._data_name_box.SetValue(str(data.title)) 
    122133        self._data = data 
     
    174185        transformed_data = self._calculator.compute_transform( 
    175186            self._extrapolated_data, self.background) 
     187        self._transformed_data = transformed_data 
    176188        from sas.sasgui.perspectives.corfunc.corfunc import TRANSFORM_LABEL 
    177189        import numpy as np 
     
    179191        plot_y = transformed_data.y[np.where(transformed_data.x <= 200)] 
    180192        self._manager.show_data(Data1D(plot_x, plot_y), TRANSFORM_LABEL) 
     193        self._extract_btn.Enable() 
     194 
     195    def extract_parameters(self, event=None): 
     196        params = None 
     197        try: 
     198            params = self._calculator.extract_parameters(self._transformed_data) 
     199        except: 
     200            params = None 
     201        if params is None: 
     202            msg = "Error extracting parameters." 
     203            wx.PostEvent(self._manager.parent, 
     204                StatusEvent(status=msg, info="error")) 
     205            return 
     206        for key in OUTPUT_STRINGS.keys(): 
     207            value = params[key] 
     208            self._output_boxes[key].SetValue(value) 
     209 
    181210 
    182211 
     
    284313                msg = "qmin must be less than qmax" 
    285314                qmin_valid = False 
    286             elif background < 0 or background > self._data.y.max(): 
    287                 msg = "background must be positive and less than highest I" 
     315            elif background > self._data.y.max(): 
     316                msg = "background must be less than highest I" 
    288317                background_valid = False 
    289318        if not qmin_valid: 
     
    420449        output_sizer = wx.GridBagSizer(5, 5) 
    421450 
    422         label_strings = [ 
    423             "Long Period (A): ", 
    424             "Average Hard Block Thickness (A): ", 
    425             "Average Interface Thickness (A): ", 
    426             "Average Core Thickness: ", 
    427             "PolyDispersity: ", 
    428             "Filling Fraction: " 
    429         ] 
    430         self._output_ids = dict() 
    431         for i in range(len(label_strings)): 
     451        self._output_boxes = dict() 
     452        i = 0 
     453        for key, value in OUTPUT_STRINGS.iteritems(): 
    432454            # Create a label and a text box for each poperty 
    433             label = wx.StaticText(self, -1, label_strings[i]) 
    434             output_box = OutputTextCtrl(self, wx.NewId(), size=(50, 20), 
     455            label = wx.StaticText(self, -1, value) 
     456            output_box = OutputTextCtrl(self, wx.NewId(), 
    435457                value="-", style=wx.ALIGN_CENTER_HORIZONTAL) 
    436458            # Save the ID of each of the text boxes for accessing after the 
    437459            # output data has been calculated 
    438             self._output_ids[label_strings[i]] = output_box.GetId() 
     460            self._output_boxes[key] = output_box 
    439461            output_sizer.Add(label, (i, 0), (1, 1), wx.LEFT | wx.EXPAND, 15) 
    440462            output_sizer.Add(output_box, (i, 2), (1, 1), 
    441463                wx.RIGHT | wx.EXPAND, 15) 
     464            i += 1 
    442465 
    443466        outputbox_sizer.Add(output_sizer, wx.TOP, 0) 
     
    454477        extrapolate_btn = wx.Button(self, wx.NewId(), "Extrapolate") 
    455478        self._transform_btn = wx.Button(self, wx.NewId(), "Transform") 
    456         self._compute_btn = wx.Button(self, wx.NewId(), "Compute Measuments") 
     479        self._extract_btn = wx.Button(self, wx.NewId(), "Compute Measuments") 
    457480 
    458481        self._transform_btn.Disable() 
    459         self._compute_btn.Disable() 
     482        self._extract_btn.Disable() 
    460483 
    461484        extrapolate_btn.Bind(wx.EVT_BUTTON, self.compute_extrapolation) 
    462485        self._transform_btn.Bind(wx.EVT_BUTTON, self.compute_transform) 
     486        self._extract_btn.Bind(wx.EVT_BUTTON, self.extract_parameters) 
    463487 
    464488        controls_sizer.Add(extrapolate_btn, wx.CENTER | wx.EXPAND) 
    465489        controls_sizer.Add(self._transform_btn, wx.CENTER | wx.EXPAND) 
    466         controls_sizer.Add(self._compute_btn, wx.CENTER | wx.EXPAND) 
     490        controls_sizer.Add(self._extract_btn, wx.CENTER | wx.EXPAND) 
    467491 
    468492        controlbox_sizer.Add(controls_sizer, wx.TOP | wx.EXPAND, 0) 
Note: See TracChangeset for help on using the changeset viewer.