source: sans/Dev/trunk/NCNR_User_Procedures/SANS/Analysis/Models/NewModels_2008/Fractal_PolySphere_v40.ipf @ 396

Last change on this file since 396 was 396, checked in by ajj, 14 years ago

Adding models:

  • Fractal model with polydisperse spheres as base unit
  • Core + NShells with n=1,2,3,4
  • PolyCore? + NShells with n=1,2,3,4
File size: 6.5 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4// plots scattering from a mass fractal object
5// uses the model of Teixeria
6//
7// REFERENCE: J. Appl. Cryst. vol 21, p781-785
8//  Uses eq.1, 4, and 16
9//
10// - basic subunit is a polydisperse (Schulz) sphere - SRK July 2008
11//
12// Macro for fractal parameters added JGB 2004
13
14Proc PlotFractalPolySphere(num,qmin,qmax)                                               
15        Variable num=128,qmin=0.001,qmax=0.5
16        Prompt num "Enter number of data points for model: "
17        Prompt qmin "Enter minimum q-value (A^-1) for model: "
18        Prompt qmax "Enter maximum q-value (A^-1) for model: "
19       
20        Make/O/D/n=(num) xwave_fraPolySph,ywave_fraPolySph                                     
21        xwave_fraPolySph = alog(log(qmin) + x*((log(qmax)-log(qmin))/num))                                     
22        Make/O/D coef_fraPolySph = {0.05,5,0.1,2,100,2e-6,6.35e-6,0}                                           
23        make/o/t parameters_fraPolySph = {"Volume Fraction (scale)","Block Radius (A)","block polydispersity (0,1)","fractal dimension","correlation length (A)","SLD block (A-2)","SLD solvent (A-2)","bkgd (cm^-1 sr^-1)"}           
24        Edit parameters_fraPolySph,coef_fraPolySph             
25       
26        Variable/G root:g_fraPolySph
27        g_fraPolySph := FractalPolySphere(coef_fraPolySph,ywave_fraPolySph,xwave_fraPolySph)                   
28        Display ywave_fraPolySph vs xwave_fraPolySph                                                   
29        ModifyGraph log=1,marker=29,msize=2,mode=4                     
30        Label bottom "q (A\\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("FractalPolySphere","coef_fraPolySph","fraPolySph")
35End
36
37// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
38Proc PlotSmearedFractalPolySphere(str)                                                         
39        String str
40        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
41       
42        // if any of the resolution waves are missing => abort
43        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
44                Abort
45        endif
46       
47        SetDataFolder $("root:"+str)
48       
49        // Setup parameter table for model function
50        Make/O/D smear_coef_fraPolySph =  {0.05,5,0.1,2,100,2e-6,6.35e-6,0}
51        make/o/t smear_parameters_fraPolySph = {"Volume Fraction (scale)","Block Radius (A)","block polydispersity (0,1)","fractal dimension","correlation length (A)","SLD block (A-2)","SLD solvent (A-2)","bkgd (cm^-1 sr^-1)"}
52        Edit smear_parameters_fraPolySph,smear_coef_fraPolySph                                  //display parameters in a table
53       
54        // output smeared intensity wave, dimensions are identical to experimental QSIG values
55        // make extra copy of experimental q-values for easy plotting
56        Duplicate/O $(str+"_q") smeared_fraPolySph,smeared_qvals                                //
57        SetScale d,0,0,"1/cm",smeared_fraPolySph                                                //
58                                       
59        Variable/G gs_fraPolySph=0
60        gs_fraPolySph := fSmearedFractalPolySphere(smear_coef_fraPolySph,smeared_fraPolySph,smeared_qvals)      //this wrapper fills the STRUCT
61       
62        Display smeared_fraPolySph vs smeared_qvals             
63        ModifyGraph log=1,marker=29,msize=2,mode=4
64        Label bottom "q (A\\S-1\\M)"
65        Label left "I(q) (cm\\S-1\\M)"
66        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
67       
68        SetDataFolder root:
69        AddModelToStrings("SmearedFractalPolySphere","smear_coef_fraPolySph","fraPolySph")
70End
71
72
73
74//calculates the physical parameters related to the
75//model parameters. See the reference at the top of the
76//file for details
77//
78// this macro is currently only appicable to the monodisperse case and must be appropriately
79// modified before use
80//
81//Macro NumberDensity_FractalPolySphere()
82//     
83//      Variable nden,phi,r0,Df,corr,s0,vpoly,i0,rg
84//     
85//      if(Exists("coef_fraPolySph")!=1)
86//              abort "You need to plot the unsmeared model first to create the coefficient table"
87//      Endif
88//     
89//      phi = coef_fraPolySph[0]   // volume fraction of building blocks
90//      r0 = coef_fraPolySph[1]    // building block radius
91//      Df = coef_fraPolySph[2]    // fractal dimension
92//      corr = coef_fraPolySph[3]  // fractal correlation length (of cluster)
93//     
94//      Print "mean building block radius (A) = ",r0
95//      Print "volume fraction = ",phi
96//     
97// //  average particle volume         
98//      vpoly = 4*Pi/3*r0^3
99//      nden = phi/vpoly                //nden in 1/A^3
100//      i0 = 1.0e8*phi*vpoly*(coef_fraPolySph[4]-coef_fraPolySph[5])^2  // 1/cm/sr
101//      rg = corr*( Df*(Df+1)/2 )^0.5
102//      s0 = exp(gammln(Df+1))*(corr/r0)^Df
103//      Print "number density (A^-3) = ",nden
104//      Print "Guinier radius (A) = ",rg
105//      Print "Aggregation number G = ",s0
106//      Print "Forward cross section of building blocks (cm-1 sr-1) I(0) = ",i0
107//      Print "Forward cross section of clusters (cm-1 sr-1) I(0) = ",i0*s0
108//End
109
110//AAO version, uses XOP if available
111// simply calls the original single point calculation with
112// a wave assignment (this will behave nicely if given point ranges)
113Function FractalPolySphere(cw,yw,xw) : FitFunc
114        Wave cw,yw,xw
115       
116#if exists("FractalPolySphereX")
117        yw = FractalPolySphereX(cw,xw)
118#else
119        yw = fFractalPolySphere(cw,xw)
120#endif
121        return(0)
122End
123
124//fractal scattering function
125Function fFractalPolySphere(w,x) :FitFunc
126        wave w
127        variable x
128       
129        variable r0,Df,corr,phi,sldp,sldm,bkg
130        variable pq,sq,ans,pd
131        phi=w[0]   // volume fraction of building block spheres...
132        r0=w[1]    //  radius of building block
133        pd = w[2]               // polydispersity of sphere
134        Df=w[3]     //  fractal dimension
135        corr=w[4] //  correlation length of fractal-like aggregates
136        sldp = w[5] // SLD of building block
137        sldm = w[6] // SLD of matrix or solution
138        bkg=w[7]  //  flat background
139       
140        //calculate P(q) for the spherical subunits, units cm-1 sr-1
141//      pq = 1.0e8*phi*4/3*pi*r0^3*(sldp-sldm)^2*(3*(sin(x*r0) - x*r0*cos(x*r0))/(x*r0)^3)^2
142        Make/O/D/N=6 tmp_SchSph
143        tmp_SchSph[0] = phi
144        tmp_SchSph[1] = r0
145        tmp_SchSph[2] = pd
146        tmp_SchSph[3] = sldp
147        tmp_SchSph[4] = sldm
148        tmp_SchSph[5] = 0
149
150#if exists("SchulzSpheresX")
151        pq = SchulzSpheresX(tmp_SchSph,x)
152#else
153        pq = fSchulzSpheres(tmp_SchSph,x)
154#endif
155
156        //calculate S(q)
157        sq = Df*exp(gammln(Df-1))*sin((Df-1)*atan(x*corr))
158        sq /= (x*r0)^Df * (1 + 1/(x*corr)^2)^((Df-1)/2)
159        sq += 1
160        //combine and return
161        ans = pq*sq + bkg
162
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 fSmearedFractalPolySphere(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 = SmearedFractalPolySphere(fs)
187       
188        return (0)
189End
190
191//the smeared model calculation
192Function SmearedFractalPolySphere(s) :FitFunc
193        Struct ResSmearAAOStruct &s
194
195//      the name of your unsmeared model (AAO) is the first argument
196        Smear_Model_20(FractalPolySphere,s.coefW,s.xW,s.yW,s.resW)
197
198        return(0)
199End
Note: See TracBrowser for help on using the repository browser.