# source:sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/NewModels_2008/Dumbbell_v40.ipf@510

Last change on this file since 510 was 510, checked in by srkline, 14 years ago

Simple change in all of the model function files to include the name of the parameter wave in the Keyword=list that is generated when a model is plotted. This is becoming an issue where the proper parameter wave can't be deduced from just the suffix, then there is nothing to put in the table.

I should have added this when I initially wrote the wrapper...

File size: 7.1 KB
Line
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4////////////////////////////////////////////////////
5//
6// calculates the scattering of a dumbbell, two spheres in contact, possibly
7// with overlap (a di-colloid). The length of the cylinder is zero, but is
8// programmed as 0.01 A to avoid numerical errors (I'm too lazy to go through
9// the math to find the proper limits)
10//
11//
12// a double integral is used, both using Gaussian quadrature
13// routines that are now included with GaussUtils
14//
16//
17//
18// REFERENCE:
19//              H. Kaya, J. Appl. Cryst. (2004) 37, 223-230.
20//              H. Kaya and N-R deSouza, J. Appl. Cryst. (2004) 37, 508-509. (addenda and errata)
21//
22////////////////////////////////////////////////////
23
24//this macro sets up all the necessary parameters and waves that are
25//needed to calculate the model function.
26//
27Proc PlotDumbbell(num,qmin,qmax)
28        Variable num=100, qmin=.001, qmax=.7
29        Prompt num "Enter number of data points for model: "
30        Prompt qmin "Enter minimum q-value (^1) for model: "
31        Prompt qmax "Enter maximum q-value (^1) for model: "
32        //
33        Make/O/D/n=(num) xwave_Dumb, ywave_Dumb
34        xwave_Dumb =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
35        Make/O/D coef_Dumb = {1,20,40,1e-6,6.3e-6,0}                    //CH#2
36        make/o/t parameters_Dumb = {"Scale Factor","cylinder radius rc (A)","end cap radius R >= rc (A)","SLD cylinder (A^-2)","SLD solvent (A^-2)","Incoherent Bgd (cm-1)"}    //CH#3
37        Edit parameters_Dumb, coef_Dumb
38
39        Variable/G root:g_Dumb
40        g_Dumb := Dumbbell(coef_Dumb, ywave_Dumb, xwave_Dumb)
41        Display ywave_Dumb vs xwave_Dumb
42        ModifyGraph marker=29, msize=2, mode=4
43        ModifyGraph log=1
44        Label bottom "q (A\\S-1\\M)"
45        Label left "I(q) (cm\\S-1\\M)"
46        AutoPositionWindow/M=1/R=\$(WinName(0,1)) \$WinName(0,2)
47
49//
50End
51
52
53// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
54Proc PlotSmearedDumbbell(str)
55        String str
56        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
57
58        // if any of the resolution waves are missing => abort
59        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
60                Abort
61        endif
62
63        SetDataFolder \$("root:"+str)
64
65        // Setup parameter table for model function
66        Make/O/D smear_coef_Dumb = {1,20,40,1e-6,6.3e-6,0}              //CH#4
67        make/o/t smear_parameters_Dumb = {"Scale Factor","cylinder radius rc (A)","end cap radius R >= rc (A)","SLD cylinder (A^-2)","SLD solvent (A^-2)","Incoherent Bgd (cm-1)"}
68        Edit smear_parameters_Dumb,smear_coef_Dumb                                      //display parameters in a table
69
70        // output smeared intensity wave, dimensions are identical to experimental QSIG values
71        // make extra copy of experimental q-values for easy plotting
72        Duplicate/O \$(str+"_q") smeared_Dumb,smeared_qvals                              //
73        SetScale d,0,0,"1/cm",smeared_Dumb                                                      //
74
75        Variable/G gs_Dumb=0
76        gs_Dumb := fSmearedDumbbell(smear_coef_Dumb,smeared_Dumb,smeared_qvals) //this wrapper fills the STRUCT
77
78        Display smeared_Dumb vs smeared_qvals                                                                   //
79        ModifyGraph log=1,marker=29,msize=2,mode=4
80        Label bottom "q (A\\S-1\\M)"
81        Label left "I(q) (cm\\S-1\\M)"
82        AutoPositionWindow/M=1/R=\$(WinName(0,1)) \$WinName(0,2)
83
84        SetDataFolder root:
86End
87
88
89
90//AAO version, uses XOP if available
91// simply calls the original single point calculation with
92// a wave assignment (this will behave nicely if given point ranges)
93Function Dumbbell(cw,yw,xw) : FitFunc
94        Wave cw,yw,xw
95
96#if exists("DumbbellX")
97        yw = DumbbellX(cw,xw)
98#else
99        yw = fDumbbell(cw,xw)
100#endif
101        return(0)
102End
103
104//
105// - a double integral - choose points wisely - 76 for both...
106//
107Function fDumbbell(w,x) : FitFunc
108        Wave w
109        Variable x
110//       Input (fitting) variables are:
111        //[0] scale factor
112        //[1] cylinder radius (little r)
113        //[2] cylinder length (big L)
114        //[3] end cap radius (big R)
115        //[4] sld cylinder (A^-2)
116        //[5] sld solvent
117        //[6] incoherent background (cm^-1)
118//      give them nice names
119        Variable scale,contr,bkg,inten,sldc,slds
121        scale = w[0]
123//      len = w[2]
125        sldc = w[3]
126        slds = w[4]
127        bkg = w[5]
128
130
131        Make/O/D/N=7 Dumb_tmp
132        Dumb_tmp[0] = w[0]
133        Dumb_tmp[1] = w[1]
134        Dumb_tmp[2] = 0.01              //length is some small number, essentially zero
135        Dumb_tmp[3] = w[2]
136        Dumb_tmp[4] = w[3]
137        Dumb_tmp[5] = w[4]
138        Dumb_tmp[6] = w[5]
139
140
141        contr = sldc-slds
142
143        Variable/G root:gDumTheta=0,root:gDumT=0
144
145        inten = IntegrateFn76(Dumb_Outer,0,pi/2,Dumb_tmp,x)
146
148        inten *= 1e8            //convert to cm^-1
149        inten *= contr*contr
150        inten *= scale
151        inten += bkg
152
153        Return (inten)
154End
155
156// outer integral
157// x is the q-value
158Function Dumb_Outer(w,x,dum)
159        Wave w
160        Variable x,dum
161
162        Variable retVal
163        Variable scale,contr,bkg,inten,sldc,slds
165        scale = w[0]
167        len = w[2]
169        sldc = w[4]
170        slds = w[5]
171        bkg = w[6]
172
174
175        NVAR dTheta = root:gDumTheta
176        NVAR dt = root:gDumT
177        dTheta = dum
179
180        Variable arg1,arg2
181        arg1 = x*len/2*cos(dum)
183
185
186        retVal *= retval*sin(dum)               // = |A(q)|^2*sin(theta)
187
188        return(retVal)
189End
190
191//returns the value of the integrand of the inner integral
192Function Dumb_Inner(w,x,dum)
193        Wave w
194        Variable x,dum
195
196        Variable retVal
197        Variable scale,contr,bkg,inten,sldc,slds
199        scale = w[0]
201        len = w[2]
203        sldc = w[4]
204        slds = w[5]
205        bkg = w[6]
207
208        NVAR dTheta = root:gDumTheta
209        NVAR dt = root:gDumT
210        dt = dum
211
212        retVal = Dumb(w,x,dt,dTheta)
213
215
216        return(retVal)
217End
218
219Function Dumb(w,x,tt,Theta)
220        Wave w
221        Variable x,tt,Theta
222
223        Variable val,arg1,arg2
224        Variable scale,contr,bkg,inten,sldc,slds
226        scale = w[0]
228        len = w[2]
230        sldc = w[4]
231        slds = w[5]
232        bkg = w[6]
233
235
238
239        val = cos(arg1)*(1-tt*tt)*Besselj(1,arg2)/arg2
240
241        return(val)
242end
243
244//wrapper to calculate the smeared model as an AAO-Struct
245// fills the struct and calls the ususal function with the STRUCT parameter
246//
247// used only for the dependency, not for fitting
248//
249Function fSmearedDumbbell(coefW,yW,xW)
250        Wave coefW,yW,xW
251
252        String str = getWavesDataFolder(yW,0)
253        String DF="root:"+str+":"
254
255        WAVE resW = \$(DF+str+"_res")
256
257        STRUCT ResSmearAAOStruct fs
258        WAVE fs.coefW = coefW
259        WAVE fs.yW = yW
260        WAVE fs.xW = xW
261        WAVE fs.resW = resW
262
263        Variable err
264        err = SmearedDumbbell(fs)
265
266        return (0)
267End
268
269// this is all there is to the smeared calculation!
270//
271// 20 points should be fine here. This function is not much different than cylinders, where 20 is sufficient
272Function SmearedDumbbell(s) :FitFunc
273        Struct ResSmearAAOStruct &s
274
275//      the name of your unsmeared model (AAO) is the first argument
276        Smear_Model_20(Dumbbell,s.coefW,s.xW,s.yW,s.resW)
277
278        return(0)
279End
Note: See TracBrowser for help on using the repository browser.