source: sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/ProlateCoreShell_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.8 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4////////////////////////////////////////////////
5// GaussUtils.proc and PlotUtils.proc MUST be included for the smearing calculation to compile
6// Adopting these into the experiment will insure that they are always present
7////////////////////////////////////////////////
8// this function is for the form factor of a prolate ellipsoid
9//
10// 06 NOV 98 SRK
11////////////////////////////////////////////////
12
13Proc PlotProlateForm(num,qmin,qmax)
14        Variable num=128,qmin=0.001,qmax=0.7
15        Prompt num "Enter number of data points for model: "
16        Prompt qmin "Enter minimum q-value (A^-1) for model: "
17        Prompt qmax "Enter maximum q-value (A^-1) for model: "
18       
19        Make/O/D/n=(num) xwave_pef,ywave_pef
20        xwave_pef =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
21        Make/O/D coef_pef = {1.,100,50,110,60,1e-6,2e-6,6.3e-6,0.001}
22        make/o/t parameters_pef = {"scale","major core radius (A)","minor core radius (A)","major shell radius (A)","minor shell radius (A)","SLD core (A-2)","SLD shell (A-2)","SLD solvent (A-2)","bkg (cm-1)"}
23        Edit parameters_pef,coef_pef
24        Variable/G root:g_pef
25        g_pef := ProlateForm(coef_pef,ywave_pef,xwave_pef)
26//      ywave_pef := ProlateForm(coef_pef,xwave_pef)
27        Display ywave_pef vs xwave_pef
28        ModifyGraph log=1,marker=29,msize=2,mode=4
29        Label bottom "q (A\\S-1\\M)"
30        Label left "Intensity (cm\\S-1\\M)"
31        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
32       
33        AddModelToStrings("ProlateForm","coef_pef","parameters_pef","pef")
34End
35
36///////////////////////////////////////////////////////////
37// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
38Proc PlotSmearedProlateForm(str)                                                               
39        String str
40        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
41       
42        // if any of the resolution waves are missing => abort
43        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
44                Abort
45        endif
46       
47        SetDataFolder $("root:"+str)
48       
49        // Setup parameter table for model function
50        Make/O/D smear_coef_pef = {1.,100,50,110,60,1e-6,2e-6,6.3e-6,0.001}
51        make/o/t smear_parameters_pef = {"scale","major core radius (A)","minor core radius (A)","major shell radius (A)","minor shell radius (A)","SLD core (A-2)","SLD shell (A-2)","SLD solvent (A-2)","bkg (cm-1)"}
52        Edit smear_parameters_pef,smear_coef_pef
53       
54        // output smeared intensity wave, dimensions are identical to experimental QSIG values
55        // make extra copy of experimental q-values for easy plotting
56        Duplicate/O $(str+"_q") smeared_pef,smeared_qvals                               //**** mod
57        SetScale d,0,0,"1/cm",smeared_pef                                               
58               
59        Variable/G gs_pef=0
60        gs_pef := fSmearedProlateForm(smear_coef_pef,smeared_pef,smeared_qvals) //this wrapper fills the STRUCT
61
62        Display smeared_pef vs smeared_qvals
63        ModifyGraph log=1,marker=29,msize=2,mode=4
64        Label bottom "q (A\\S-1\\M)"
65        Label left "Intensity (cm\\S-1\\M)"
66        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 
67       
68        SetDataFolder root:
69        AddModelToStrings("SmearedProlateForm","smear_coef_pef","smear_parameters_pef","pef")
70End
71
72//AAO version
73Function ProlateForm(cw,yw,xw) : FitFunc
74        Wave cw,yw,xw
75
76#if exists("ProlateFormX")
77        yw = ProlateFormX(cw,xw)
78#else
79        yw = fProlateForm(cw,xw)
80#endif
81        return(0)
82End
83
84///////////////////////////////////////////////////////////////
85// unsmeared model calculation
86///////////////////////////
87Function fProlateForm(w,x) : FitFunc
88        Wave w
89        Variable x
90
91//The input variables are (and output)
92        //[0] scale
93        //[1] crmaj, major radius of core       [A]
94        //[2] crmin, minor radius of core
95        //[3] trmaj, overall major radius
96        //[4] trmin, overall minor radius
97        //[5] sld core, [A^-2]
98        //[6] sld shell,
99        //[7] sld solvent
100        //[8] bkg [cm-1]
101        Variable scale,crmaj,crmin,trmaj,trmin,delpc,delps,bkg,sldc,slds,sld
102        scale = w[0]
103        crmaj = w[1]
104        crmin = w[2]
105        trmaj = w[3]
106        trmin = w[4]
107        sldc = w[5]
108        slds = w[6]
109        sld = w[7]
110        bkg = w[8]
111
112        delpc = sldc - slds                     //core - shell
113        delps = slds - sld                      //shell - solvent
114// local variables
115        Variable yyy,va,vb,ii,nord,zi,qq,summ,nfn,npro,answer,prolatevol
116        String weightStr,zStr
117       
118        weightStr = "gauss76wt"
119        zStr = "gauss76z"
120
121//      if wt,z waves don't exist, create them
122
123        if (WaveExists($weightStr) == 0) // wave reference is not valid,
124                Make/D/N=76 $weightStr,$zStr
125                Wave w76 = $weightStr
126                Wave z76 = $zStr                // wave references to pass
127                Make76GaussPoints(w76,z76)     
128        else
129                if(exists(weightStr) > 1)
130                         Abort "wave name is already in use"    // execute if condition is false
131                endif
132                Wave w76 = $weightStr
133                Wave z76 = $zStr
134        endif
135
136// set up the integration
137        // end points and weights
138        nord = 76
139        nfn = 2         //only <f^2> is calculated
140        npro = 1        // PROLATE ELLIPSOIDS
141        va =0
142        vb =1
143//move this zi(i) evaluation inside other nord loop, since I don't have an array
144//      i=0
145//      do
146//       zi[i] = ( z76[i]*(vb-va) + vb + va )/2.0
147 //       i +=1
148 //     while (i<nord)
149//
150// evaluate at Gauss points
151        // remember to index from 0,size-1
152
153        qq = x          //current x point is the q-value for evaluation
154        summ = 0.0
155        ii=0
156        do
157                //printf "top of nord loop, i = %g\r",i
158                if(nfn ==1) //then              // "f1" required for beta factor
159                        if(npro ==1) //then     // prolate
160                                zi = ( z76[ii]*(vb-va) + vb + va )/2.0 
161//           yyy = w76[ii]*gfn1(zi,crmaj,crmin,trmaj,trmin,delpc,delps,qq)
162                        Endif
163//
164                        if(npro ==0) //then     // oblate 
165                                zi = ( z76[ii]*(vb-va) + vb + va )/2.0
166//            yyy = w76[i]*gfn3(zi,crmaj,crmin,trmaj,trmin,delpc,delps,qq)
167                        Endif
168                Endif           //nfn = 1
169          //
170                if(nfn !=1) //then              //calculate"f2" = <f^2> = averaged form factor
171                        if(npro ==1) //then     //prolate
172                                zi = ( z76[ii]*(vb-va) + vb + va )/2.0
173                                yyy = w76[ii]*gfn2(zi,crmaj,crmin,trmaj,trmin,delpc,delps,qq)
174                                //printf "yyy = %g\r",yyy
175                        Endif
176//
177                        if(npro ==0) //then     //oblate
178                                zi = ( z76[ii]*(vb-va) + vb + va )/2.0
179//              yyy = w76[ii]*gfn4(zi,crmaj,crmin,trmaj,trmin,delpc,delps,qq)
180                        Endif
181                Endif           //nfn <>1
182         
183                summ = yyy + summ               // get running total of integral
184                ii+=1
185        while (ii<nord)                         // end of loop over quadrature points
186//   
187// calculate value of integral to return
188
189        answer = (vb-va)/2.0*summ
190       
191        //normailze by particle volume
192        prolatevol = 4*Pi/3*trmaj*trmin*trmin
193        answer /= prolatevol
194       
195        // rescale from 1/A to 1/cm
196        answer *= 1.0e8
197        //scale (arb)
198        answer *= scale
199        ////then add in background
200        answer += bkg
201
202        Return (answer)
203End     //prolate form factor
204
205//
206//     FUNCTION gfn2:    CONTAINS F(Q,A,B,mu)**2  AS GIVEN
207//                       BY (53) AND (56,57) IN CHEN AND
208//                       KOTLARCHYK REFERENCE
209//
210//     <PROLATE ELLIPSOIDS>
211//
212Function gfn2(xx,crmaj,crmin,trmaj,trmin,delpc,delps,qq)
213        Variable xx,crmaj,crmin,trmaj,trmin,delpc,delps,qq
214        // local variables
215        Variable aa,bb,u2,ut2,uq,ut,vc,vt,gfnc,gfnt,tgfn,gfn2,pi43,gfn
216
217        PI43=4.0/3.0*PI
218        aa = crmaj
219        bb = crmin
220        u2 = (aa*aa*xx*xx + bb*bb*(1.0-xx*xx))
221        ut2 = (trmaj*trmaj*xx*xx + trmin*trmin*(1.0-xx*xx))
222        uq = sqrt(u2)*qq
223        ut= sqrt(ut2)*qq
224        vc = PI43*aa*bb*bb
225        vt = PI43*trmaj*trmin*trmin
226        gfnc = 3.0*(sin(uq)/uq/uq - cos(uq)/uq)/uq*vc*delpc
227        gfnt = 3.0*(sin(ut)/ut/ut - cos(ut)/ut)/ut*vt*delps
228        gfn = gfnc+gfnt
229        gfn2 = gfn*gfn
230       
231        return gfn2
232End             //function gfn2 for prolate ellipsoids
233
234
235// this is all there is to the smeared calculation!
236Function SmearedProlateForm(s) :FitFunc
237        Struct ResSmearAAOStruct &s
238
239////the name of your unsmeared model is the first argument
240        Smear_Model_20(ProlateForm,s.coefW,s.xW,s.yW,s.resW)
241
242        return(0)
243End
244
245//wrapper to calculate the smeared model as an AAO-Struct
246// fills the struct and calls the ususal function with the STRUCT parameter
247//
248// used only for the dependency, not for fitting
249//
250Function fSmearedProlateForm(coefW,yW,xW)
251        Wave coefW,yW,xW
252       
253        String str = getWavesDataFolder(yW,0)
254        String DF="root:"+str+":"
255       
256        WAVE resW = $(DF+str+"_res")
257       
258        STRUCT ResSmearAAOStruct fs
259        WAVE fs.coefW = coefW   
260        WAVE fs.yW = yW
261        WAVE fs.xW = xW
262        WAVE fs.resW = resW
263       
264        Variable err
265        err = SmearedProlateForm(fs)
266       
267        return (0)
268End
Note: See TracBrowser for help on using the repository browser.