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

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

Typo in dialog in every model...

File size: 6.3 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 an ellipsoid of rotation with uniform scattering length density
9//
10// 06 NOV 98 SRK
11////////////////////////////////////////////////
12
13Proc PlotEllipsoidForm(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 (^-1) for model: "
17        Prompt qmax "Enter maximum q-value (^-1) for model: "
18       
19        Make/O/D/n=(num) xwave_eor,ywave_eor
20        xwave_eor =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
21        Make/O/D coef_eor = {1.,20.,400,3.0e-6,0.01}
22        make/o/t parameters_eor = {"scale","R a (rotation axis) (A)","R b (A)","Contrast (A^-2)","incoh. bkg (cm^-1)"}
23        Edit parameters_eor,coef_eor
24        Variable/G root:g_eor
25        g_eor := EllipsoidForm(coef_eor,ywave_eor,xwave_eor)
26//      ywave_eor := EllipsoidForm(coef_eor,xwave_eor)
27        Display ywave_eor vs xwave_eor
28        ModifyGraph log=1,marker=29,msize=2,mode=4
29        Label bottom "q (\\S-1\\M)"
30        Label left "Intensity (cm\\S-1\\M)"
31        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
32End
33///////////////////////////////////////////////////////////
34
35// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
36Proc PlotSmearedEllipsoidForm(str)     
37        String str
38        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
39       
40        // if no gQvals wave, data must not have been loaded => abort
41        if(ResolutionWavesMissingDF(str))
42                Abort
43        endif
44       
45        SetDataFolder $("root:"+str)
46       
47        // Setup parameter table for model function
48        Make/O/D smear_coef_eor = {1.,20.,400,3.0e-6,0.01}
49        make/o/t smear_parameters_eor = {"scale","R a (rotation axis) (A)","R b (A)","Contrast (A^-2)","incoh. bkg (cm^-1)"}
50        Edit smear_parameters_eor,smear_coef_eor
51       
52        // output smeared intensity wave, dimensions are identical to experimental QSIG values
53        // make extra copy of experimental q-values for easy plotting
54        Duplicate/O $(str+"_q") smeared_eor,smeared_qvals
55        SetScale d,0,0,"1/cm",smeared_eor       
56
57        Variable/G root:gs_eor=0
58        gs_eor := fSmearedEllipsoidForm(smear_coef_eor,smeared_eor,smeared_qvals)
59       
60        Display smeared_eor vs smeared_qvals
61        ModifyGraph log=1,marker=29,msize=2,mode=4
62        Label bottom "q (\\S-1\\M)"
63        Label left "Intensity (cm\\S-1\\M)"
64        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
65       
66        SetDataFolder root:
67End
68
69//AAO version
70Function EllipsoidForm(cw,yw,xw) : FitFunc
71        Wave cw,yw,xw
72       
73#if exists("EllipsoidFormX")
74        yw = EllipsoidFormX(cw,xw)
75#else
76        yw = fEllipsoidForm(cw,xw)
77#endif
78        return(0)
79End
80
81///////////////////////////////////////////////////////////////
82// unsmeared model calculation
83///////////////////////////
84Function fEllipsoidForm(w,x) : FitFunc
85        Wave w
86        Variable x
87
88//The input variables are (and output)
89        //[0] scale
90        //[1] Axis of rotation
91        //[2] two equal radii
92        //[3] contrast (A^-2)
93        //[4] background (cm^-1)
94        Variable scale, ra,vra,delrho,bkg
95        scale = w[0]
96        vra = w[1]
97        ra = w[2]
98        delrho = w[3]
99        bkg = w[4]
100        //if vra < ra, OBLATE
101        //if vra > ra, PROLATE
102//
103// the OUTPUT form factor is <f^2>/Vell [cm-1]
104//
105
106// local variables
107        Variable nord,ii,va,vb,contr,vell,nden,summ,yyy,zi,qq,halfheight
108        Variable answer
109        String weightStr,zStr
110       
111        weightStr = "gauss76wt"
112        zStr = "gauss76z"
113
114       
115//      if wt,z waves don't exist, create them
116// 20 Gauss points is not enough for Ellipsoid calculation
117       
118        if (WaveExists($weightStr) == 0) // wave reference is not valid,
119                Make/D/N=76 $weightStr,$zStr
120                Wave w76 = $weightStr
121                Wave z76 = $zStr                // wave references to pass
122                Make76GaussPoints(w76,z76)     
123        //                  printf "w[0],z[0] = %g %g\r", w76[0],z76[0]
124        else
125                if(exists(weightStr) > 1)
126                         Abort "wave name is already in use"    // execute if condition is false
127                endif
128                Wave w76 = $weightStr
129                Wave z76 = $zStr                // Not sure why this has to be "declared" twice
130        //          printf "w[0],z[0] = %g %g\r", w76[0],z76[0]
131        endif
132
133
134// set up the integration
135        // end points and weights
136        nord = 76
137        va = 0
138        vb = 1
139
140// evaluate at Gauss points
141        // remember to index from 0,size-1
142
143        qq = x          //current x point is the q-value for evaluation
144      summ = 0.0                // initialize integral
145      ii=0
146      do
147                // Using 76 Gauss points
148                zi = ( z76[ii]*(vb-va) + vb + va )/2.0         
149                yyy = w76[ii] * eor(qq, ra,vra, zi)
150                summ += yyy
151
152                ii+=1
153        while (ii<nord)                         // end of loop over quadrature points
154//   
155// calculate value of integral to return
156
157      answer = (vb-va)/2.0*summ
158     
159// Multiply by contrast^2
160        answer *= delrho*delrho
161//normalize by Ellipsoid volume
162//NOTE that for this (Fournet) definition of the integral, one must MULTIPLY by Vell
163        vell=4*Pi/3*ra*ra*vra
164        answer *= vell
165//convert to [cm-1]
166        answer *= 1.0e8
167//Scale
168        answer *= scale
169// add in the background
170        answer += bkg
171
172        Return (answer)
173       
174End             //End of function EllipsoidForm()
175
176///////////////////////////////////////////////////////////////
177Function eor(qq,ra,vra,theta)
178        Variable qq,ra,vra,theta
179       
180// qq is the q-value for the calculation (1/A)
181// ra are the like radius of the Ellipsoid (A)
182// vra is the unlike semiaxis of the Ellipsoid = L/2 (A) (= the rotation axis!)
183// theta is the dummy variable of the integration
184
185   //Local variables
186        Variable retval,arg,t1, nu
187       
188        nu = vra/ra
189        arg = qq*ra*sqrt(1+theta^2*(nu^2-1))
190       
191        retval = 9*((sin(arg)-arg*cos(arg))/arg^3)^2
192       
193    return retval
194   
195End     //Function eor()
196
197// this is all there is to the smeared calculation!
198Function SmearedEllipsoidForm(s) :FitFunc
199        Struct ResSmearAAOStruct &s
200       
201        //the name of your unsmeared model is the first argument
202        s.yW = Smear_Model_20(EllipsoidForm,s.coefW,s.xW,s.resW)
203
204        return(0)
205End
206
207//wrapper to calculate the smeared model as an AAO-Struct
208// fills the struct and calls the ususal function with the STRUCT parameter
209//
210// used only for the dependency, not for fitting
211//
212Function fSmearedEllipsoidForm(coefW,yW,xW)
213        Wave coefW,yW,xW
214       
215        String str = getWavesDataFolder(yW,0)
216        String DF="root:"+str+":"
217       
218        WAVE resW = $(DF+str+"_res")
219       
220        STRUCT ResSmearAAOStruct fs
221        WAVE fs.coefW = coefW   
222        WAVE fs.yW = yW
223        WAVE fs.xW = xW
224        WAVE fs.resW = resW
225       
226        Variable err
227        err = SmearedEllipsoidForm(fs)
228       
229        return (0)
230End
Note: See TracBrowser for help on using the repository browser.