# source:sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/NewModels_2008/Two_Power_Law_v40.ipf@445

Last change on this file since 445 was 445, checked in by srkline, 14 years ago

Adding a pile of new model functions. Some are from Boualem, some old code from Danilo, some from me.

As far as I can tell, these replicate the original references (not the empirical models).

Help file and XOPs have net been done yet.

File size: 4.7 KB
Line
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4//
5// empirical model to fit two power law regions. Model is a "v". No attempt is made to
6// smooth the transition at the crossover q-value. The two power law slopes are the important
7// results, not fudging some crossover function.
8//
9// JUN 2008 SRK
10
11//////////////////////////////////
12Proc PlotTwoPowerLaw(num,qmin,qmax)
13        Variable num=512, qmin=.001, qmax=.2
14        Prompt num "Enter number of data points for model: "
15        Prompt qmin "Enter minimum q-value (^1) for model: "
16         Prompt qmax "Enter maximum q-value (^1) for model: "
17//
18        Make/O/D/n=(num) xwave_TwoPowerLaw, ywave_TwoPowerLaw
19        xwave_TwoPowerLaw =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
20        Make/O/D coef_TwoPowerLaw = {1e-6, 4, 1, 0.01, 0}
21        make/o/t parameters_TwoPowerLaw = {"Coefficient, A ", "(-)Low Q Power","(-) high Q Power","Crossover Qc (A-1)","Incoherent Bgd (cm-1)"}
22        Edit parameters_TwoPowerLaw, coef_TwoPowerLaw
23        Variable/G root:g_TwoPowerLaw
24        g_TwoPowerLaw  := TwoPowerLaw(coef_TwoPowerLaw, ywave_TwoPowerLaw, xwave_TwoPowerLaw)
25//      ywave_TwoPowerLaw  := TwoPowerLaw(coef_TwoPowerLaw, xwave_TwoPowerLaw)
26        Display ywave_TwoPowerLaw vs xwave_TwoPowerLaw
27        ModifyGraph marker=29, msize=2, mode=4
28        ModifyGraph log(left)=1
29        ModifyGraph log(bottom)=1
30        Label bottom "q (\\S-1\\M) "
31        Label left "Power-Law (cm\\S-1\\M)"
32        AutoPositionWindow/M=1/R=\$(WinName(0,1)) \$WinName(0,2)
33
35//
36End
37
38////////////////////////////////////////////////////
39// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
40Proc PlotSmearedTwoPowerLaw(str)
41        String str
42        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
43
44        // if any of the resolution waves are missing => abort
45        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
46                Abort
47        endif
48
49        SetDataFolder \$("root:"+str)
50
51        // Setup parameter table for model function
52        Make/O/D smear_coef_TwoPowerLaw = {1e-6, 4, 1, 0.01, 0}
53        make/o/t smear_parameters_TwoPowerLaw = {"Coefficient, A ", "(-)Low Q Power","(-) high Q Power","Crossover Qc (A-1)","Incoherent Bgd (cm-1)"}
54        Edit smear_parameters_TwoPowerLaw,smear_coef_TwoPowerLaw                                        //display parameters in a table
55
56        // output smeared intensity wave, dimensions are identical to experimental QSIG values
57        // make extra copy of experimental q-values for easy plotting
58        Duplicate/O \$(str+"_q") smeared_TwoPowerLaw,smeared_qvals                               //
59        SetScale d,0,0,"1/cm",smeared_TwoPowerLaw                                                       //
60
61
62        Variable/G gs_TwoPowerLaw=0
63        gs_TwoPowerLaw := fSmearedTwoPowerLaw(smear_coef_TwoPowerLaw,smeared_TwoPowerLaw,smeared_qvals) //this wrapper fills the STRUCT
64
65        Display smeared_TwoPowerLaw vs smeared_qvals                                                                    //
66        ModifyGraph log=1,marker=29,msize=2,mode=4
67        Label bottom "q (\\S-1\\M)"
68        Label left "TwoPowerLaw (cm\\S-1\\M)"
69        AutoPositionWindow/M=1/R=\$(WinName(0,1)) \$WinName(0,2)
70
71        SetDataFolder root:
73//
74End     // end macro PlotSmearedTwoPowerLaw
75
76//AAO version
77Function TwoPowerLaw(cw,yw,xw) : FitFunc
78        Wave cw,yw,xw
79
80#if exists("TwoPowerLawX")
81        yw = TwoPowerLawX(cw,xw)
82#else
83        yw = fTwoPowerLaw(cw,xw)
84#endif
85        return(0)
86End
87
88Function fTwoPowerLaw(w,x) : FitFunc
89        Wave w
90        Variable x
91//       Input (fitting) variables are:
92        //[0] Coefficient
93        //[1] (-) Power @ low Q
94        //[2] (-) Power @ high Q
95        //[3] crossover Q-value
96        //[4] incoherent background
97//      give them nice names
98        Variable A, m1,m2,qc,bgd,scale
99        A = w[0]
100        m1 = w[1]
101        m2 = w[2]
102        qc = w[3]
103        bgd = w[4]
104
105//      local variables
106        Variable inten, qval
107//      x is the q-value for the calculation
108        qval = x
109//      do the calculation and return the function value
110
111        if(qval<=qc)
112                inten = A*qval^-m1
113        else
114                scale = A*qc^-m1 / qc^-m2
115                inten = scale*qval^-m2
116        endif
117
118        inten += bgd
119        Return (inten)
120End
121/////////////////////////////////////////////////////////////////////////////////
122
123
124// this is all there is to the smeared calculation!
125Function SmearedTwoPowerLaw(s) :FitFunc
126        Struct ResSmearAAOStruct &s
127
128////the name of your unsmeared model is the first argument
129        Smear_Model_20(TwoPowerLaw,s.coefW,s.xW,s.yW,s.resW)
130
131        return(0)
132End
133
134//wrapper to calculate the smeared model as an AAO-Struct
135// fills the struct and calls the ususal function with the STRUCT parameter
136//
137// used only for the dependency, not for fitting
138//
139Function fSmearedTwoPowerLaw(coefW,yW,xW)
140        Wave coefW,yW,xW
141
142        String str = getWavesDataFolder(yW,0)
143        String DF="root:"+str+":"
144
145        WAVE resW = \$(DF+str+"_res")
146
147        STRUCT ResSmearAAOStruct fs
148        WAVE fs.coefW = coefW
149        WAVE fs.yW = yW
150        WAVE fs.xW = xW
151        WAVE fs.resW = resW
152
153        Variable err
154        err = SmearedTwoPowerLaw(fs)
155
156        return (0)
157End
Note: See TracBrowser for help on using the repository browser.