source: sasview/sansmodels/src/sans/models/ReflectivityModel.py @ 20f00bed

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 20f00bed was 339ce67, checked in by Jae Cho <jhjcho@…>, 14 years ago

added some models and tests
  • Property mode set to 100644
File size: 12.9 KB
Line 
1   
2from sans.models.BaseComponent import BaseComponent
3from sans.models.ReflModel import ReflModel
4import copy
5from math import floor
6from scipy.special import erf
7func_list = {'Erf':0, 'Linear':1, 'RParabolic':2, \
8                     'LParabola':3, 'RCubic':4, 'LCubic':5}
9max_nshells = 10
10class ReflectivityModel(BaseComponent):
11    """
12    This multi-model is based on Parratt formalism and provides the capability
13    of changing the number of layers between 0 and 10.
14    """
15    def __init__(self, multfactor=1):
16        BaseComponent.__init__(self)
17        """
18        :param multfactor: number of layers in the model, assumes 0<= n_shells <=10.
19        """
20
21        ## Setting  model name model description
22        self.description=""
23        model = ReflModel()
24        self.model = model
25        self.name = "ReflectivityModel"
26        self.description=model.description
27        self.n_layers = multfactor
28        ## Define parameters
29        self.params = {}
30
31        ## Parameter details [units, min, max]
32        self.details = {}
33       
34        # non-fittable parameters
35        self.non_fittable = model.non_fittable
36       
37        # list of function in order of the function number
38        self.fun_list = self._get_func_list()
39        ## dispersion
40        self._set_dispersion()
41        ## Define parameters
42        self._set_params()
43       
44        ## Parameter details [units, min, max]
45        self._set_details()
46       
47        #list of parameter that can be fitted
48        self._set_fixed_params() 
49        self.model.params['n_layers'] = self.n_layers
50       
51        ## functional multiplicity info of the model
52        # [int(maximum no. of functionality),"str(Titl),
53        # [str(name of function0),...], [str(x-asix name of sld),...]]
54        self.multiplicity_info = [max_nshells,"No. of Layers:",[],['Depth']]
55
56   
57    def _clone(self, obj):
58        """
59        Internal utility function to copy the internal
60        data members to a fresh copy.
61        """
62        obj.params     = copy.deepcopy(self.params)
63        obj.non_fittable     = copy.deepcopy(self.non_fittable)
64        obj.description     = copy.deepcopy(self.description)
65        obj.details    = copy.deepcopy(self.details)
66        obj.dispersion = copy.deepcopy(self.dispersion)
67        obj.model  = self.model.clone()
68
69        return obj
70   
71   
72    def _set_dispersion(self):
73        """
74        model dispersions
75        """ 
76        ##set dispersion from model
77        self.dispersion = {}
78                   
79
80    def _set_params(self):
81        """
82        Concatenate the parameters of the model to create
83        this model parameters
84        """
85        # rearrange the parameters for the given # of shells
86        for name , value in self.model.params.iteritems():
87            n = 0
88            pos = len(name.split('_'))-1
89            if name.split('_')[0] == 'sldIM':
90                continue
91            elif name.split('_')[0] == 'func':
92                n= -1
93                while n<self.n_layers:
94                    n += 1
95                    if name.split('_')[pos] == 'inter%s' % str(n):
96                        self.params[name]=value
97                        continue
98                #continue
99            elif name.split('_')[pos][0:5] == 'inter':
100                n= -1
101                while n<self.n_layers:
102                    n += 1
103                    if name.split('_')[pos] == 'inter%s' % str(n):
104                        self.params[name]= value
105                        continue
106            elif name.split('_')[pos][0:4] == 'flat':
107                while n<self.n_layers:
108                    n += 1
109                    if name.split('_')[pos] == 'flat%s' % str(n):
110                        self.params[name]= value
111                        continue
112            elif name == 'n_layers':
113                continue
114            else:
115                self.params[name]= value
116               
117        self.model.params['n_layers'] = self.n_layers   
118
119        # set constrained values for the original model params
120        self._set_xtra_model_param()       
121
122    def _set_details(self):
123        """
124        Concatenate details of the original model to create
125        this model details
126        """
127        for name ,detail in self.model.details.iteritems():
128            if name in self.params.iterkeys():
129                self.details[name]= detail
130           
131   
132    def _set_xtra_model_param(self):
133        """
134        Set params of original model that are hidden from this model
135        """
136        # look for the model parameters that are not in param list
137        for key in self.model.params.iterkeys():
138            if key not in self.params.keys():
139                if  key.split('_')[0] == 'thick':
140                    self.model.setParam(key, 0)
141                    continue
142                if  key.split('_')[0] == 'func':
143                    self.model.setParam(key, 0)
144                    continue
145                for nshell in range(self.n_layers,max_nshells):
146                    if key.split('_')[1] == 'flat%s' % str(nshell+1):
147                        try:
148                            if key.split('_')[0] == 'sld':
149                                value = self.model.params['sld_medium']
150                            elif key.split('_')[0] == 'sldIM':
151                                value = self.model.params['sldIM_medium']
152                            self.model.setParam(key, value)
153                        except: pass
154   
155    def _get_func_list(self):
156        """
157        Get the list of functions in each layer (shell)
158        """
159        #func_list = {}
160        return func_list
161       
162    def getProfile(self):
163        """
164        Get SLD profile
165       
166        : return: (z, beta) where z is a list of depth of the transition points
167                beta is a list of the corresponding SLD values
168        """
169        # max_pts for each layers
170        n_sub = 21
171        z = []
172        beta = []
173        sub_range = floor(n_sub/2.0)
174        z.append(0)
175        beta.append(self.params['sld_sub0']) 
176       
177        z0 = 0
178        # for layers from the top
179        for n in range(1,self.n_layers+2):
180            i = n
181
182            for j in range(0,2):
183                for n_s in range(-sub_range,sub_range+1):
184                    if j==1:
185                        if i==self.n_layers+1:
186                            break
187                        # shift half sub thickness for the first point
188                        z0 += dz/2.0
189                        z.append(z0)
190                        #z0 -= dz/2.0
191                        z0 += self.params['thick_flat%s'% str(i)]
192                       
193                        sld_i = self.params['sld_flat%s'% str(i)]
194                        beta.append(self.params['sld_flat%s'% str(i)])
195                    else:
196                       
197                        dz = self.params['thick_inter%s'% str(i-1)]/n_sub
198                       
199                        if n_s == -sub_range:
200                            # shift half sub thickness for the first point
201                            z0 -= dz/2.0
202                        #exec "dz = self.params['thick_inter[%s-1]'% str(i)]/9"
203                        #print "%d = %g \n"% (i,self.params['thick_inter3'])
204                        z0 += dz
205
206                        if i == 1:
207                            sld_l = self.params['sld_sub0']
208                        else:
209                            sld_l = self.params['sld_flat%s'% str(i-1)]
210                        if i == self.n_layers+1:
211                            sld_r = self.params['sld_medium']
212                        else:
213                            sld_r = self.params['sld_flat%s'% str(i)]
214                        func_idx = self.params['func_inter%s'% str(i-1)]
215                        func = self._get_func(n_s, n_sub, func_idx)
216                        if sld_r>sld_l:
217                            sld_i = (sld_r-sld_l)*func+sld_l
218                        elif sld_r<sld_l:
219                            sld_i = (sld_l-sld_r)*(1-func)+sld_r
220                        else:
221                             sld_i = sld_r
222                    z.append(z0)
223                    beta.append(sld_i)
224                    if j==1: break
225        # put substrate and superstrate profile
226        # shift half sub thickness for the first point
227        z0 += dz/2.0
228        z.append(z0)
229        beta.append(self.params['sld_medium']) 
230        z_ext = z0/6.0
231       
232        # put the extra points for the substrate
233        # and superstrate
234        z.append(z0+z_ext)
235        beta.append(self.params['sld_medium']) 
236        z.insert(0,-z_ext)
237        beta.insert(0,self.params['sld_sub0']) 
238        z = [z0 - x for x in z]
239        z.reverse()
240        beta.reverse() 
241        return z, beta
242   
243    def _get_func(self, index, n_sub, func_idx):
244        """
245        Get the function asked to buil sld profile
246        : param index: index of sub_layer
247        : param n_sub: total number of sub_layer
248        : param func_idx: an integer to identify a function
249       
250        : return out: the output from the function, float
251        """
252        # cal bin_size
253        bin_size = 1.0/n_sub
254        # erf
255        if func_idx == 0:
256            out = erf(index/(n_sub/5.0))/2.0 + 0.5
257            return out
258        else:
259            index += 0.5
260        # linear
261        if func_idx == 1:
262            out = ((index + floor(n_sub/2.0))*bin_size)
263        # r_parabolic
264        elif func_idx == 2:
265            out = ((index + floor(n_sub/2.0))*bin_size)* \
266                    ((index + floor(n_sub/2.0))*bin_size)
267        # l_parabolic
268        elif func_idx == 3:
269            out = 1.0-(((index + floor(n_sub/2.0))*bin_size) - 1.0) *\
270                    (((index + floor(n_sub/2.0))*bin_size) - 1.0)
271        # r_cubic
272        elif func_idx == 4:
273            out = ((index + floor(n_sub/2.0))*bin_size)* \
274                    ((index + floor(n_sub/2.0))*bin_size)* \
275                    ((index + floor(n_sub/2.0))*bin_size)
276        # l_cubic
277        elif func_idx == 5:
278            out = 1.0+(((index + floor(n_sub/2.0)))*bin_size - 1.0) *\
279                    (((index + floor(n_sub/2.0)))*bin_size - 1.0) *\
280                    (((index + floor(n_sub/2.0)))*bin_size - 1.0)
281        # return output
282        return out
283   
284    def setParam(self, name, value):
285        """
286        Set the value of a model parameter
287   
288        : param name: name of the parameter
289        : param value: value of the parameter
290        """
291        # set param to new model
292        self._setParamHelper( name, value)
293       
294        ## setParam to model
295        if name=='sld_medium':
296            # the sld_*** model.params not in params must set to value of sld_solv
297            for key in self.model.params.iterkeys():
298                if key not in self.params.keys()and key.split('_')[0] == 'sld':
299                        self.model.setParam(key, value)
300           
301        self.model.setParam( name, value)
302
303    def _setParamHelper(self, name, value):
304        """
305        Helper function to setParam
306        """
307
308        # Look for standard parameter
309        for item in self.params.keys():
310            if item.lower()==name.lower():
311                self.params[item] = value
312                return
313       
314        raise ValueError, "Model does not contain parameter %s" % name
315             
316   
317    def _set_fixed_params(self):
318        """
319        Fill the self.fixed list with the model fixed list
320        """
321        pass         
322
323    def run(self, x = 0.0):
324        """
325        Evaluate the model
326       
327        :param x: input q, or [q,phi]
328       
329        :return: scattering function P(q)
330       
331        """
332
333        return self.model.run(x)
334
335    def runXY(self, x = 0.0):
336        """
337        Evaluate the model
338       
339        : param x: input q-value (float or [float, float] as [qx, qy])
340        : return: scattering function value
341        """ 
342
343        return self.model.runXY(x)
344   
345    ## Now (May27,10) directly uses the model eval function
346    ## instead of the for-loop in Base Component.
347    def evalDistribution(self, x = []):
348        """
349        Evaluate the model in cartesian coordinates
350       
351        : param x: input q[], or [qx[], qy[]]
352        : return: scattering function P(q[])
353        """
354        # set effective radius and scaling factor before run
355        return self.model.evalDistribution(x)
356    def calculate_ER(self):
357        """
358        """
359        return self.model.calculate_ER()
360    def set_dispersion(self, parameter, dispersion):
361        """
362        Set the dispersion object for a model parameter
363       
364        : param parameter: name of the parameter [string]
365        :dispersion: dispersion object of type DispersionModel
366        """
367        pass
Note: See TracBrowser for help on using the repository browser.