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

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

(1) - cursors can now be used to select a subrange of USANS data to fit. This is done by th fit wrapper, assigning a subrange of resW to the struct

(2) all of the smeared model functions are now in the latest form of Smear_Model_N() that is NOT a pointwise calculation anymore, since the USANS matrix smearing in inherently not so.

File size: 5.3 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 PlotEllipCylinderForm(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 := EllipCyl20(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 PlotSmearedEllipCylForm(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 := fSmearedEllipCylForm(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 EllipCyl20(cw,yw,xw) : FitFunc
83        Wave cw,yw,xw
84       
85#if exists("EllipCyl20X")
86        yw = EllipCyl20X(cw,xw)
87#else
88        yw = fEllipCyl20(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 fEllipCyl20(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 fSmearedEllipCylForm(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 = SmearedEllipCylForm(fs)
187       
188        return (0)
189End
190
191// this is all there is to the smeared calculation!
192Function SmearedEllipCylForm(s) :FitFunc
193        Struct ResSmearAAOStruct &s
194
195//      the name of your unsmeared model (AAO) is the first argument
196        Smear_Model_20(EllipCyl20,s.coefW,s.xW,s.yW,s.resW)
197
198        return(0)
199End
200       
Note: See TracBrowser for help on using the repository browser.