source: sans/Analysis/branches/ajj_23APR07/IGOR_Package_Files/Put in User Procedures/SANS_Models_v3.00/NewModels_2006/EllipticalCylinder.ipf @ 153

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

Changed Plot* and PlotSmeared?* naming schemes to be all consistent prefixes for the actual function name, so that the macros can be constructed from the function name, or vice versa.

also some tweaks to the wrapper to make sure that plot and append really work

File size: 5.4 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4///////////////////////////////////////////
5// this function is for the form factor of an cylinder with an ellipsoidal cross-section
6// and a uniform scattering length density
7//
8// 06 NOV 98 SRK
9//
10// re-written to not use MacOS XOP for calculation
11// now requires the "new" version of GaussUtils that includes the generic quadrature routines
12//
13// 09 SEP 03 SRK
14////////////////////////////////////////////////
15
16Proc PlotEllipticalCylinder(num,qmin,qmax)
17        Variable num=50,qmin=0.001,qmax=0.7
18        Prompt num "Enter number of data points for model: "
19        Prompt qmin "Enter minimum q-value (^-1) for model: "
20        Prompt qmax "Enter maximum q-value (^-1) for model: "
21       
22//      //constants needed for the integration if qtrap is used (in a separate procedure file!)
23//      Variable/G root:gNumPoints=200
24//      Variable/G root:gTol=1e-5
25//      Variable/G root:gMaxIter=20
26        //
27        Make/O/D/n=(num) xwave_ecf,ywave_ecf
28        xwave_ecf =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
29        Make/O/D coef_ecf = {1.,20.,1.5,400,3.0e-6,0.0}
30        make/o/t parameters_ecf = {"scale","minor radius (A)","nu = major/minor (-)","length ()","SLD diff (A^-2)","incoh. bkg (cm^-1)"}
31        Edit parameters_ecf,coef_ecf
32       
33        Variable/G root:g_ecf
34        g_ecf := EllipticalCylinder(coef_ecf,ywave_ecf,xwave_ecf)
35        Display ywave_ecf vs xwave_ecf
36        ModifyGraph log=1,marker=29,msize=2,mode=4
37        Label bottom "q (\\S-1\\M)"
38        Label left "Intensity (cm\\S-1\\M)"
39        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
40End
41///////////////////////////////////////////////////////////
42
43// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
44Proc PlotSmearedEllipticalCylinder(str)                                                         
45        String str
46        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
47       
48        // if any of the resolution waves are missing => abort
49        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
50                Abort
51        endif
52       
53        SetDataFolder $("root:"+str)
54       
55        // Setup parameter table for model function
56        Make/O/D smear_coef_ecf =  {1.,20.,1.5,400,3.0e-6,0.0}
57        make/o/t smear_parameters_ecf = {"scale","minor radius (A)","nu = major/minor (-)","length ()","SLD diff (A^-2)","incoh. bkg (cm^-1)"}
58        Edit smear_parameters_ecf,smear_coef_ecf
59       
60        // output smeared intensity wave, dimensions are identical to experimental QSIG values
61        // make extra copy of experimental q-values for easy plotting
62        Duplicate/O $(str+"_q") smeared_ecf,smeared_qvals
63        SetScale d,0,0,"1/cm",smeared_ecf       
64                                       
65        Variable/G gs_ecf=0
66        gs_ecf := fSmearedEllipticalCylinder(smear_coef_ecf,smeared_ecf,smeared_qvals)  //this wrapper fills the STRUCT
67       
68        Display smeared_ecf vs smeared_qvals
69        ModifyGraph log=1,marker=29,msize=2,mode=4
70        Label bottom "q (\\S-1\\M)"
71        Label left "Intensity (cm\\S-1\\M)"
72        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
73       
74        SetDataFolder root:
75End
76
77
78
79//AAO version, uses XOP if available
80// simply calls the original single point calculation with
81// a wave assignment (this will behave nicely if given point ranges)
82Function EllipticalCylinder(cw,yw,xw) : FitFunc
83        Wave cw,yw,xw
84       
85#if exists("EllipticalCylinderX")
86        yw = EllipticalCylinderX(cw,xw)
87#else
88        yw = fEllipticalCylinder(cw,xw)
89#endif
90        return(0)
91End
92
93// the main function that calculates the form factor of an elliptical cylinder
94// integrates EllipCyl_integrand, which itself is an integral function
95// 20 points of quadrature seems to be sufficient for both integrals
96//
97Function fEllipticalCylinder(w,x)       : FitFunc
98        Wave w
99        Variable x
100       
101        Variable inten,scale,rad,nu,len,contr,bkg,ii
102        scale = w[0]
103        rad = w[1]
104        nu = w[2]
105        len = w[3]
106        contr = w[4]
107        bkg = w[5]
108       
109        inten = IntegrateFn20(EllipCyl_Integrand,0,1,w,x)
110       
111        //multiply by volume
112        inten *= Pi*rad*rad*nu*len
113        inten *= 1e8    //convert to 1/cm
114        inten *= contr*contr
115        inten *= scale
116        inten += bkg
117       
118        return(inten)
119End
120
121//the outer integral
122Function EllipCyl_Integrand(w,x,dum)
123        Wave w
124        Variable x,dum
125       
126        Variable val,rad,arg,len
127        rad = w[1]
128        len = w[3]
129       
130        arg = rad*sqrt(1-dum^2)
131        duplicate/O w temp_w
132        Wave temp_w=temp_w
133        temp_w[1] = arg         //replace radius with transformed variable
134        val = (1/pi)*IntegrateFn20(Phi_EC,0,Pi,temp_w,x)
135       
136        // equivalent to the 20-pt quadrature
137//      val = (1/pi)*qtrap(Phi_EC,temp_w,x,0,Pi,1e-3,20)
138       
139        arg = x*len*dum/2
140        if(arg==0)
141                val *= 1
142        else
143                val *= sin(arg)*sin(arg)/arg/arg
144        endif
145        //Print "val=",val
146        return(val)
147End
148
149//the inner integral
150Function Phi_EC(w,x,dum)
151        Wave w
152        Variable x,dum
153       
154        Variable ans,arg,aa,nu
155        aa = w[1]               // = rad*sqrt(1-dum^2)
156        nu = w[2]
157        arg = x*aa*sqrt( (1+nu^2)/2 + (1-nu^2)/2*cos(dum) )
158        if(arg==0)
159                ans = (2*0.5)^2         // == 1
160        else
161                ans = 2*2*bessJ(1,arg)*bessJ(1,arg)/arg/arg
162        endif
163        return(ans)
164End
165
166//wrapper to calculate the smeared model as an AAO-Struct
167// fills the struct and calls the ususal function with the STRUCT parameter
168//
169// used only for the dependency, not for fitting
170//
171Function fSmearedEllipticalCylinder(coefW,yW,xW)
172        Wave coefW,yW,xW
173       
174        String str = getWavesDataFolder(yW,0)
175        String DF="root:"+str+":"
176       
177        WAVE resW = $(DF+str+"_res")
178       
179        STRUCT ResSmearAAOStruct fs
180        WAVE fs.coefW = coefW   
181        WAVE fs.yW = yW
182        WAVE fs.xW = xW
183        WAVE fs.resW = resW
184       
185        Variable err
186        err = SmearedEllipticalCylinder(fs)
187       
188        return (0)
189End
190
191// this is all there is to the smeared calculation!
192Function SmearedEllipticalCylinder(s) :FitFunc
193        Struct ResSmearAAOStruct &s
194
195//      the name of your unsmeared model (AAO) is the first argument
196        Smear_Model_20(EllipticalCylinder,s.coefW,s.xW,s.yW,s.resW)
197
198        return(0)
199End
200       
Note: See TracBrowser for help on using the repository browser.