source: sans/Dev/trunk/NCNR_User_Procedures/Analysis/VSANS/V_BroadPeak_WB.ipf @ 1095

Last change on this file since 1095 was 1095, checked in by srkline, 4 years ago

Added two model functions for white beam smearing.

Many other small changes for processing of the back detector, shuffling of VSANS menu items, and consistent naming of V_ procedures.

File size: 6.4 KB
Line 
1#pragma TextEncoding = "MacRoman"
2#pragma rtGlobals=3             // Use modern global access method and strict wave access.
3
4
5////////////////////////////////////////////////////
6//
7// an empirical model containing power law scattering + a broad peak
8//
9// B. Hammouda OCT 2008
10//
11//
12// updated for use with latest macros SRK Nov 2008
13//
14// Updated 2018 to White Beam Smearing
15//
16//
17////////////////////////////////////////////////////
18
19//
20Proc PlotBroadPeakWB(num,qmin,qmax)
21        Variable num=200, qmin=0.001, qmax=0.7
22        Prompt num "Enter number of data points for model: "
23        Prompt qmin "Enter minimum q-value (^-1) for model: "
24        Prompt qmax "Enter maximum q-value (^-1) for model: "
25//
26        Make/O/D/n=(num) xwave_BroadPeakWB, ywave_BroadPeakWB
27        xwave_BroadPeakWB =  alog(log(qmin) + x*((log(qmax)-log(qmin))/num))
28        Make/O/D coef_BroadPeakWB = {1e-5, 3, 10, 50.0, 0.1,2,0.1}             
29        make/o/t parameters_BroadPeakWB = {"Porod Scale", "Porod Exponent","Lorentzian Scale","Lor Screening Length [A]","Qzero [1/A]","Lorentzian Exponent","Bgd [1/cm]"}      //CH#2
30        Edit parameters_BroadPeakWB, coef_BroadPeakWB
31       
32        Variable/G root:g_BroadPeakWB
33        g_BroadPeakWB := BroadPeakWB(coef_BroadPeakWB, ywave_BroadPeakWB, xwave_BroadPeakWB)
34        Display ywave_BroadPeakWB vs xwave_BroadPeakWB
35        ModifyGraph marker=29, msize=2, mode=4
36        ModifyGraph log=1,grid=1,mirror=2
37        Label bottom "q (\\S-1\\M) "
38        Label left "I(q) (cm\\S-1\\M)"
39        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
40       
41        AddModelToStrings("BroadPeakWB","coef_BroadPeakWB","parameters_BroadPeakWB","BroadPeakWB")
42//
43End
44
45
46//
47//no input parameters are necessary, it MUST use the experimental q-values
48// from the experimental data read in from an AVE/QSIG data file
49////////////////////////////////////////////////////
50// - sets up a dependency to a wrapper, not the actual SmearedModelFunction
51Proc PlotSmearedBroadPeakWB(str)                                                               
52        String str
53        Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4)
54       
55        // if any of the resolution waves are missing => abort
56        if(ResolutionWavesMissingDF(str))               //updated to NOT use global strings (in GaussUtils)
57                Abort
58        endif
59       
60        SetDataFolder $("root:"+str)
61       
62        // Setup parameter table for model function
63        Make/O/D smear_coef_BroadPeakWB = {1e-5, 3, 10, 50.0, 0.1,2,0.1}               
64        make/o/t smear_parameters_BroadPeakWB = {"Porod Scale", "Porod Exponent","Lorentzian Scale","Lor Screening Length [A]","Qzero [1/A]","Lorentzian Exponent","Bgd [1/cm]"}
65        Edit smear_parameters_BroadPeakWB,smear_coef_BroadPeakWB                                        //display parameters in a table
66       
67        // output smeared intensity wave, dimensions are identical to experimental QSIG values
68        // make extra copy of experimental q-values for easy plotting
69        Duplicate/O $(str+"_q") smeared_BroadPeakWB,smeared_qvals
70        SetScale d,0,0,"1/cm",smeared_BroadPeakWB
71                                       
72        Variable/G gs_BroadPeakWB=0
73        gs_BroadPeakWB := fSmearedBroadPeakWB(smear_coef_BroadPeakWB,smeared_BroadPeakWB,smeared_qvals) //this wrapper fills the STRUCT
74       
75        Display smeared_BroadPeakWB 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)
80       
81        SetDataFolder root:
82        AddModelToStrings("SmearedBroadPeakWB","smear_coef_BroadPeakWB","smear_parameters_BroadPeakWB","BroadPeakWB")
83End
84
85
86//
87//AAO version, uses XOP if available
88// simply calls the original single point calculation with
89// a wave assignment (this will behave nicely if given point ranges)
90Function BroadPeakWB(cw,yw,xw) : FitFunc
91        Wave cw,yw,xw
92       
93#if exists("BroadPeakX")
94//      yw = BroadPeakX(cw,xw)
95        yw = V_fBroadPeakWB(cw,xw)
96#else
97//      yw = fBroadPeakWB(cw,xw)
98        yw = 1
99#endif
100        return(0)
101End
102
103//
104// unsmeared model calculation
105//
106Function V_fBroadPeakWB(w,x) : FitFunc
107        Wave w
108        Variable x
109       
110        // variables are:                                                       
111        //[0] Porod term scaling
112        //[1] Porod exponent
113        //[2] Lorentzian term scaling
114        //[3] Lorentzian screening length [A]
115        //[4] peak location [1/A]
116        //[5] Lorentzian exponent
117        //[6] background
118       
119        Variable aa,nn,cc,LL,Qzero,mm,bgd,inten,lolim,uplim
120
121        // define limits based on lo/mean, hi/mean of the wavelength distribution
122        // using the empirical definition, "middle" of the peaks
123        loLim = 3.37/5.3
124        upLim = 8.37/5.3
125       
126        inten = V_IntegrBroadPeakWB_mid(w,loLim,upLim,x)
127
128// why do I need this? Is this because this is defined as the mean of the distribution
129//  and is needed to normalize the integral? verify this on paper.     
130        inten *= 5.3
131
132// normalize the integral       
133        inten /= 19933          // "middle"  of peaks
134
135// additional normalization???
136        inten /= 1.05           //
137        Return (inten)
138       
139End
140
141// the trick here is that declaring the last qVal wave as a variable
142// since this is implicitly called N times in the wave assignment of the answer wave
143Function V_IntegrBroadPeakWB_mid(cw,loLim,upLim,qVal)
144        Wave cw
145        Variable loLim,upLim
146        Variable qVal
147       
148        Variable/G root:qq = qval
149        Variable ans
150       
151//      ans = Integrate1D(V_intgrnd_top,lolim,uplim,2,0,cw)             //adaptive quadrature
152        ans = Integrate1D(V_integrand_BroadPeakWB,lolim,uplim,1,0,cw)           // Romberg integration
153       
154        return ans
155end
156
157Function V_integrand_BroadPeakWB(cw,dum)
158        Wave cw
159        Variable dum            // the dummy of the integration
160
161        Variable val
162        NVAR qq = root:qq               //the q-value of the integration, not part of cw, so pass global
163//      SVAR funcStr = root:gFunctionString
164//      FUNCREF SANSModel_proto func = $funcStr
165
166        val = V_WhiteBeamDist_mid(dum*5.3)*BroadPeakX(cw,qq/dum)
167       
168        return (val)
169End
170
171//CH#4 
172///////////////////////////////////////////////////////////////
173// smeared model calculation
174//
175// you don't need to do anything with this function, as long as
176// your BroadPeak works correctly, you get the resolution-smeared
177// version for free.
178//
179// this is all there is to the smeared model calculation!
180Function SmearedBroadPeakWB(s) : FitFunc
181        Struct ResSmearAAOStruct &s
182
183//      the name of your unsmeared model (AAO) is the first argument
184        Smear_Model_76(BroadPeakWB,s.coefW,s.xW,s.yW,s.resW)
185
186        return(0)
187End
188
189
190///////////////////////////////////////////////////////////////
191
192
193// nothing to change here
194//
195//wrapper to calculate the smeared model as an AAO-Struct
196// fills the struct and calls the ususal function with the STRUCT parameter
197//
198// used only for the dependency, not for fitting
199//
200Function fSmearedBroadPeakWB(coefW,yW,xW)
201        Wave coefW,yW,xW
202       
203        String str = getWavesDataFolder(yW,0)
204        String DF="root:"+str+":"
205       
206        WAVE resW = $(DF+str+"_res")
207       
208        STRUCT ResSmearAAOStruct fs
209        WAVE fs.coefW = coefW   
210        WAVE fs.yW = yW
211        WAVE fs.xW = xW
212        WAVE fs.resW = resW
213       
214        Variable err
215        err = SmearedBroadPeakWB(fs)
216       
217        return (0)
218End
Note: See TracBrowser for help on using the repository browser.