source: sans/Analysis/branches/ajj_23APR07/IGOR_Package_Files/Put in User Procedures/SANS_Models_v3.00/OblateForm.ipf @ 166

Last change on this file since 166 was 166, checked in by srkline, 15 years ago

Modified all procedure files to add to the keyword=value strings to identify the function, coefficients, and suffix once the model is plotted (and the objects will exist)

a one-liner in the Plot and PlotSmeared? macros.

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