source: sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/CoreShell_v40.ipf @ 633

Revision 633, 5.1 KB checked in by srkline, 4 years ago (diff)

Corrected models to explicitly return proper values for I(q=0). There are some models that just can't be fixed, and these typically return NaN. Some, however, are simply numerically unstable at extreme conditions. Beware.

Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion=6.1
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//
9// this function is for the form factor of a sphere with a core-shell structure
10//
11// 06 NOV 98 SRK
12////////////////////////////////////////////////
13
14Proc PlotCoreShellSphere(num,qmin,qmax)
15        Variable num=128,qmin=0.001,qmax=0.7
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_css,ywave_css
21        xwave_css =alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
22        make/o/d coef_css = {1.,60,10,1e-6,2e-6,3e-6,0.001}
23        make/o/t parameters_css = {"scale","core radius (A)","shell thickness (A)","Core SLD (A-2)","Shell SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"}
24        Edit parameters_css,coef_css
25        Variable/G root:g_css
26        g_css := CoreShellSphere(coef_css,ywave_css,xwave_css)
27//      ywave_css := CoreShellSphere(coef_css,xwave_css)
28        Display ywave_css vs xwave_css
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("CoreShellSphere","coef_css","parameters_css","css")
35End
36
37///////////////////////////////////////////////////////////
38// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
39Proc PlotSmearedCoreShellSphere(str)                                                           
40        String str
41        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
42       
43        // if any of the resolution waves are missing => abort
44        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
45                Abort
46        endif
47       
48        SetDataFolder $("root:"+str)
49       
50        // Setup parameter table for model function
51        make/o/d smear_coef_css = {1.,60,10,1e-6,2e-6,3e-6,0.001}
52        make/o/t smear_parameters_css = {"scale","core radius (A)","shell thickness (A)","Core SLD (A-2)","Shell SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"}
53        Edit smear_parameters_css,smear_coef_css
54       
55        // output smeared intensity wave, dimensions are identical to experimental QSIG values
56        // make extra copy of experimental q-values for easy plotting
57       
58        Duplicate/O $(str+"_q") smeared_css,smeared_qvals                               
59        SetScale d,0,0,"1/cm",smeared_css                                                       
60                                       
61        Variable/G gs_css=0
62        gs_css := fSmearedCoreShellSphere(smear_coef_css,smeared_css,smeared_qvals)     //this wrapper fills the STRUCT
63       
64        Display smeared_css vs smeared_qvals                                                                   
65        ModifyGraph log=1,marker=29,msize=2,mode=4
66        Label bottom "q (A\\S-1\\M)"
67        Label left "Intensity (cm\\S-1\\M)"
68        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
69       
70        SetDataFolder root:
71        AddModelToStrings("SmearedCoreShellSphere","smear_coef_css","smear_parameters_css","css")
72End
73
74
75//AAO version
76Function CoreShellSphere(cw,yw,xw) : FitFunc
77        Wave cw,yw,xw
78
79#if exists("CoreShellSphereX")
80        yw = CoreShellSphereX(cw,xw)
81#else
82        yw = fCoreShellSphere(cw,xw)
83#endif
84        return(0)
85End
86
87///////////////////////////////////////////////////////////////
88// unsmeared model calculation
89///////////////////////////
90Function fCoreShellSphere(w,x) : FitFunc
91        Wave w
92        Variable x
93       
94        // variables are:
95        //[0] scale factor
96        //[1] radius of core [A]
97        //[2] thickness of the shell    [A]
98        //[3] SLD of the core   [A-2]
99        //[4] SLD of the shell
100        //[5] SLD of the solvent
101        //[6] background        [cm-1]
102       
103        // All inputs are in ANGSTROMS
104        //OUTPUT is normalized by the particle volume, and converted to [cm-1]
105       
106       
107        Variable scale,rcore,thick,rhocore,rhoshel,rhosolv,bkg
108        scale = w[0]
109        rcore = w[1]
110        thick = w[2]
111        rhocore = w[3]
112        rhoshel = w[4]
113        rhosolv = w[5]
114        bkg = w[6]
115       
116        // calculates scale *( f^2 + bkg)
117        Variable bes,f,vol,qr,contr,f2
118       
119        // core first, then add in shell
120        qr=x*rcore
121        contr = rhocore-rhoshel
122        if(qr==0)
123                bes = 1
124        else
125                bes = 3*(sin(qr)-qr*cos(qr))/qr^3
126        endif
127        vol = 4*pi/3*rcore^3
128        f = vol*bes*contr
129        //now the shell
130        qr=x*(rcore+thick)
131        if(qr==0)
132                bes = 1
133        else
134                bes = 3*(sin(qr)-qr*cos(qr))/qr^3
135        endif
136        contr = rhoshel-rhosolv
137        vol = 4*pi/3*(rcore+thick)^3
138        f += vol*bes*contr
139       
140        // normalize to particle volume and rescale from [A-1] to [cm-1]
141        f2 = f*f/vol*1.0e8
142       
143        //scale if desired
144        f2 *= scale
145        // then add in the background
146        f2 += bkg
147       
148        return (f2)
149End
150
151// this is all there is to the smeared calculation!
152Function SmearedCoreShellSphere(s) :FitFunc
153        Struct ResSmearAAOStruct &s
154
155////the name of your unsmeared model is the first argument
156        Smear_Model_20(CoreShellSphere,s.coefW,s.xW,s.yW,s.resW)
157
158        return(0)
159End
160
161//wrapper to calculate the smeared model as an AAO-Struct
162// fills the struct and calls the ususal function with the STRUCT parameter
163//
164// used only for the dependency, not for fitting
165//
166Function fSmearedCoreShellSphere(coefW,yW,xW)
167        Wave coefW,yW,xW
168       
169        String str = getWavesDataFolder(yW,0)
170        String DF="root:"+str+":"
171       
172        WAVE resW = $(DF+str+"_res")
173       
174        STRUCT ResSmearAAOStruct fs
175        WAVE fs.coefW = coefW   
176        WAVE fs.yW = yW
177        WAVE fs.xW = xW
178        WAVE fs.resW = resW
179       
180        Variable err
181        err = SmearedCoreShellSphere(fs)
182       
183        return (0)
184End
Note: See TracBrowser for help on using the repository browser.