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

Last change on this file since 127 was 127, checked in by srkline, 16 years ago

LOTS of changes to the analysis ipf files:

-- see sphere.ipf for the simplest example of the changes --

  • #pragma Igor 6
  • #if directive to look for XOP
  • AAO unsmeared functions
  • STRUCT functions for smearing (also AAO)
  • new wrappers for dependencies to struct functions

(2006 models have NOT been completed yet, only the old models)

  • loading data files into data folders (PlotUtils?) + some streamlining of the loaders
  • Smear_Model_N is now AAO + some streamlining of the quadrature code

-- SHS and SW structure factor XOPs are crashing (need DP wave, I may have old XOP)
-- this breaks fitting of the smeared models until wrappers can be devised
-- all packages will be broken due to the new data folder structure
-- multiple instances of functions will now cause problems (MSA)
-- RPA model is a problem with its odd functional form (extra wave)

-- lots of other carnage to follow as the bugs and typos are shaken out

24 JUL 2007 SRK

File size: 5.5 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4////////////////////////////////////////////////////
5//
6// calculates the scattering from a rectangular solid
7// i.e. a parallelepiped with sides a < b < c
8//
9// - the user must make sure that the constraints are not violated
10// otherwise the calculation will not be correct
11//
12// From: Mittelbach and Porod, Acta Phys. Austriaca 14 (1961) 185-211.
13//                              equations (1), (13), and (14) (in German!)
14//
15// note that the equations listed in Feigin and Svergun appears
16// to be wrong - they use equation (12), which does not appear to
17// be a complete orientational average (?)
18//
19// a double integral is used, both using Gaussian quadrature
20// routines that are now included with GaussUtils
21// 20-pt quadrature appears to be enough, 76 pt is available
22// by changing the function calls
23//
24////////////////////////////////////////////////////
25
26//this macro sets up all the necessary parameters and waves that are
27//needed to calculate the model function.
28//
29Proc Plot_Parallelepiped(num,qmin,qmax)
30        Variable num=100, qmin=.001, qmax=.7
31        Prompt num "Enter number of data points for model: "
32        Prompt qmin "Enter minimum q-value (^1) for model: "
33        Prompt qmax "Enter maximum q-value (^1) for model: "
34        //
35        Make/O/D/n=(num) xwave_Parallelepiped, ywave_Parallelepiped
36        xwave_Parallelepiped =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
37        Make/O/D coef_Parallelepiped = {1,35,75,400,6e-6,0}                     //CH#2
38        make/o/t parameters_Parallelepiped = {"Scale Factor","Shortest Edge A ()","B ()","Longest Edge C ()","Contrast (^-2)","Incoherent Bgd (cm-1)"}      //CH#3
39        Edit parameters_Parallelepiped, coef_Parallelepiped
40        ywave_Parallelepiped  := Parallelepiped(coef_Parallelepiped, xwave_Parallelepiped)
41        Display ywave_Parallelepiped vs xwave_Parallelepiped
42        ModifyGraph marker=29, msize=2, mode=4
43        ModifyGraph log=1
44        Label bottom "q (\\S-1\\M) "
45        Label left "I(q) (cm\\S-1\\M)"
46        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
47//
48End
49
50//
51//this macro sets up all the necessary parameters and waves that are
52//needed to calculate the  smeared model function.
53//
54//no input parameters are necessary, it MUST use the experimental q-values
55// from the experimental data read in from an AVE/QSIG data file
56////////////////////////////////////////////////////
57Proc PlotSmeared_Parallelepiped()                                                               //Parallelepiped
58        // if no gQvals wave, data must not have been loaded => abort
59        if(ResolutionWavesMissing())
60                Abort
61        endif
62       
63       
64        // Setup parameter table for model function
65        Make/O/D smear_coef_Parallelepiped = {1,35,75,400,6e-6,0}               //CH#4
66        make/o/t smear_parameters_Parallelepiped = {"Scale Factor","Shortest Edge A ()","B ()","Longest Edge C ()","Contrast (^-2)","Incoherent Bgd (cm-1)"}
67        Edit smear_parameters_Parallelepiped,smear_coef_Parallelepiped                                  //display parameters in a table
68       
69        // output smeared intensity wave, dimensions are identical to experimental QSIG values
70        // make extra copy of experimental q-values for easy plotting
71        Duplicate/O $gQvals smeared_Parallelepiped,smeared_qvals                                //
72        SetScale d,0,0,"1/cm",smeared_Parallelepiped                                                    //
73
74        smeared_Parallelepiped := Parallelepiped_Smeared(smear_coef_Parallelepiped,$gQvals)             // SMEARED function name
75        Display smeared_Parallelepiped vs smeared_qvals                                                                 //
76        ModifyGraph log=1,marker=29,msize=2,mode=4
77        Label bottom "q (\\S-1\\M)"
78        Label left "I(q) (cm\\S-1\\M)"
79        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
80End     // end macro
81
82// calculates the form factor of a rectangular solid
83// - a double integral - choose points wisely
84//
85Function Parallelepiped(w,x) : FitFunc
86        Wave w
87        Variable x
88//       Input (fitting) variables are:
89        //[0] scale factor
90        //[1] Edge A (A)
91        //[2] Edge B (A)
92        //[3] Edge C (A)
93        //[4] contrast (A^-2)
94        //[5] incoherent background (cm^-1)
95//      give them nice names
96        Variable scale,aa,bb,cc,contr,bkg,inten,qq,ii,arg,mu
97        scale = w[0]
98        aa = w[1]
99        bb = w[2]
100        cc = w[3]
101        contr = w[4]
102        bkg = w[5]
103       
104//      mu = bb*x               //scale in terms of B
105//      aa = aa/bb
106//      cc = cc/bb
107       
108        inten = IntegrateFn20(PP_Outer,0,1,w,x)
109//      inten = IntegrateFn76(PP_Outer,0,1,w,x)
110       
111        inten *= aa*cc*bb               //multiply by volume
112        inten *= 1e8            //convert to cm^-1
113        inten *= contr*contr
114        inten *= scale
115        inten += bkg
116       
117        Return (inten)
118End
119
120// outer integral
121
122// x is the q-value - remember that "mu" in the notation = B*Q
123Function PP_Outer(w,x,dum)
124        Wave w
125        Variable x,dum
126       
127        Variable retVal,mu,aa,bb,cc,mudum,arg
128        aa = w[1]
129        bb = w[2]
130        cc = w[3]
131        mu = bb*x
132       
133        mudum = mu*sqrt(1-dum^2)
134        retval = IntegrateFn20(PP_inner,0,1,w,mudum)
135//      retval = IntegrateFn76(PP_inner,0,1,w,mudum)
136       
137        cc = cc/bb
138        arg = mu*cc*dum/2
139        if(arg==0)
140                retval *= 1
141        else
142                retval *= sin(arg)*sin(arg)/arg/arg
143        endif
144       
145        return(retVal)
146End
147
148//returns the integrand of the inner integral
149Function PP_Inner(w,mu,uu)
150        Wave w
151        Variable mu,uu
152       
153        Variable aa,bb,retVal,arg1,arg2,tmp1,tmp2
154       
155        //NVAR mu = root:gEvalQval              //already has been converted to S=2*pi*q
156        aa = w[1]
157        bb = w[2]
158        aa = aa/bb
159       
160        //Mu*(1-x^2)^(0.5)
161       
162        //handle arg=0 separately, as sin(t)/t -> 1 as t->0
163        arg1 = (mu/2)*cos(Pi*uu/2)
164        arg2 = (mu*aa/2)*sin(Pi*uu/2)
165        if(arg1==0)
166                tmp1 = 1
167        else
168                tmp1 = sin(arg1)*sin(arg1)/arg1/arg1
169        endif
170        if(arg2==0)
171                tmp2 = 1
172        else
173                tmp2 = sin(arg2)*sin(arg2)/arg2/arg2
174        endif
175        retval = tmp1*tmp2
176       
177        return(retVal)
178End
179
180
181// this is all there is to the smeared calculation!
182Function Parallelepiped_Smeared(w,x) :FitFunc
183        Wave w
184        Variable x
185       
186        Variable ans
187        SVAR sq = gSig_Q
188        SVAR qb = gQ_bar
189        SVAR sh = gShadow
190        SVAR gQ = gQVals
191       
192        //the name of your unsmeared model is the first argument
193        ans = Smear_Model_20(Parallelepiped,$sq,$qb,$sh,$gQ,w,x)
194
195        return(ans)
196End
Note: See TracBrowser for help on using the repository browser.