source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/USANS/U_CALC.ipf @ 926

Last change on this file since 926 was 926, checked in by srkline, 9 years ago

Added utilities to RealSpaceModeling? to allow previously saved matrices to be added together. currently they are strictly added - that is overlapping voxels are added together, which may not be the optimal behavior, but this can be changed as needed.

Also added utility to allow arbitrary rotation of a 3D object. No interpolation is done, so the rotation is lossy. Rotation is applied in xyz order, which can make a difference. Beware.

Documentation to follow (eventually).

File size: 59.0 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion=6.1
3
4
5// USANS version of SASCALC
6
7// to simulate the intensity from a USANS experiment for planning
8// see John's instrument paper:
9//
10// J. Appl. Cryst. (2005) 38 1004-1011.
11//
12// SRK JUL 2009
13
14
15// ideas:
16//
17// - more presets
18//
19// - NEED -
20// - printable output that makes sense to an instrument scientist at least
21//   for instrument setup
22//
23// - make sure that the # points <-> increment relation meshes with ICP
24
25//
26// X plot as countrate, not absolute scale
27// X 3e-5 cutoff
28// X ? don't plot lowest angle range (but needs to be in the count time)
29// - need empty beam and empty cell count rate vs. aperture (Cd vs. Gd?)
30//
31
32
33// need to add in empty and background corrections to see "reduced" data
34// or at least compare to what the empty cell would give in the same count time
35//
36// -? model the direct beam?? currently the "red" region from -1 to 0.6 is almost entirely
37// the primary beam, so it's a bit artificial (qmin is really ~ 3e-5)
38//
39
40//
41// Need T_wide, T_rock, I peak for proper absolute scaling, but I don't know if it's really important
42// to be able to simlulate to this extent
43//
44
45
46//#include "MultScatter_MonteCarlo_2D"
47
48
49
50// Bring the UCALC panel to the front
51// ALWAYS initializes the folders and variables
52// then draws the panel if necessary
53Proc Show_UCALC()
54
55        DoWindow/F UCALC
56        if(V_Flag==0)
57                Init_UCALC()
58                UCALC_Panel()
59                CalcTotalCountTime()
60                ControlUpdate/W=UCALC U_popup0          //force a pop of the function list
61        Endif
62       
63End
64
65
66//set up the global values for the angles, # points, and count time
67Proc Init_UCALC()
68       
69        //
70        NewDataFolder/O root:Simulation                         // required to calculate the RandomDeviate
71        NewDataFolder/O root:Packages:NIST:SAS                          // required to calculate the RandomDeviate
72
73//
74        NewDataFolder/O root:Packages:NIST:USANS:SIM            //for the fake raw data
75        NewDataFolder/O/S root:Packages:NIST:USANS:Globals:U_Sim                //for constants, panel, etc
76
77        Variable/G gAngLow1 = -1
78        Variable/G gAngHigh1 = 0.6
79        Variable/G gNumPts1 = 33
80        Variable/G gCtTime1 = 25
81        Variable/G gIncr1 = 0.05       
82       
83        Variable/G gAngLow2 = 0.7
84        Variable/G gAngHigh2 = 1.9
85        Variable/G gNumPts2 = 13
86        Variable/G gCtTime2 = 100
87        Variable/G gIncr2 = 0.1
88       
89        Variable/G gAngLow3 = 2
90        Variable/G gAngHigh3 = 4.8
91        Variable/G gNumPts3 = 15
92        Variable/G gCtTime3 = 300
93        Variable/G gIncr3 = 0.2
94       
95        Variable/G gAngLow4 = 5
96        Variable/G gAngHigh4 = 9.5
97        Variable/G gNumPts4 = 10
98        Variable/G gCtTime4 = 600
99        Variable/G gIncr4 = 0.5
100       
101        Variable/G gAngLow5 = 10
102        Variable/G gAngHigh5 = 19
103        Variable/G gNumPts5 = 10
104        Variable/G gCtTime5 = 1200
105        Variable/G gIncr5 = 1
106       
107        Variable/G gAngLow6 = 20
108        Variable/G gAngHigh6 = 48
109        Variable/G gNumPts6 = 15
110        Variable/G gCtTime6 = 2000
111        Variable/G gIncr6 = 2   
112       
113        Variable/G gAngLow7 = 50
114        Variable/G gAngHigh7 = 95
115        Variable/G gNumPts7 = 10
116        Variable/G gCtTime7 = 3000
117        Variable/G gIncr7 = 5   
118       
119        // results, setup values
120        String/G gFuncStr=""
121        String/G gTotTimeStr=""
122        Variable/G gAnalyzerOmega = 7.1e-7              //solid angle of the analyzer, in steradians
123//      Variable/G gBeamCurrent=25000           // (?? not used) beam current Ed*I      (n/s) for 5/8" diam = 25000 n/s
124        Variable/G gThick=0.1           //sample thickness (cm)
125        Variable/G gSamTrans=0.8
126        Variable/G g_1D_DoABS = 0               //=1 for abs scale, 0 for just counts
127        Variable/G g_1D_PlotCR = 1              //=1 to plot countrate
128        Variable/G g_1D_AddNoise = 1            // add in appropriate noise to simulation
129       
130// a box for the results
131        // ??   maybe not useful to report
132        Variable/G g_1DEstDetCR  = 0            // estimated detector count rate
133        Variable/G g_1DTotCts = 0               
134        Variable/G g_1DFracScatt= 0     // ??
135        Variable/G g_1DEstTrans = 0     // ? can I calculate this?
136        Variable/G g_MultScattFraction = 0
137
138
139// not on panel yet
140        Variable/G g_Empirical_EMP = 0          // use an emperical model for empty cell subtraction
141        Variable/G g_EmptyLevel = 0.7
142        Variable/G g_BkgLevel = 0.6
143
144        SetDataFolder root:
145       
146End
147
148//// make the display panel a graph with a control bar just as in SASCALC
149//// so that the subwindow syntax doesn't break all of the other functionality
150////
151//Window UCALC_Panel() : Graph
152//      PauseUpdate; Silent 1           // building window...
153//      Display /W=(55,44,670,850) /K=1
154//      ModifyGraph cbRGB=(36929,50412,31845)
155//      DoWindow/C UCALC
156//      DoWindow/T UCALC,"USANS Simulation"
157//      ControlBar 320
158//     
159//      GroupBox group0,pos={5,0},size={577,159},title="Instrument Setup"
160//      GroupBox group1,pos={5,165},size={240,147},title="Sample Setup"
161//      GroupBox group2,pos={327,165},size={259,147},title="Results"
162//     
163//      PopupMenu popup0,pos={17,18},size={165,20},title="Sample Aperture Diam (in)"
164//      PopupMenu popup0,mode=3,popvalue="0.625",value="0.25;0.50;0.625;0.75;1.0;1.75;2.0;"
165//      PopupMenu popup2,pos={220,18},size={165,20},title="Presets"
166//      PopupMenu popup2,mode=3,popvalue="Long Count",value="Short Count;Medium Count;Long Count;"
167//      PopupMenu popup2,proc=UCALC_PresetPopup
168//
169//      SetDataFolder root:Packages:NIST:USANS:Globals:U_Sim
170//     
171//      Variable top=44,pt=0,inc=18
172//      SetVariable setvar1a,pos={12,top},size={100,15},title="theta min",value= gAngLow1
173//      SetVariable setvar1b,pos={119,top},size={100,15},title="theta max",value= gAngHigh1
174//      SetVariable setvar1c,pos={227,top},size={100,15},title="increm",value= gIncr1
175//      SetVariable setvar1d,pos={335,top},size={100,15},title="# points",value= gNumPts1
176//      SetVariable setvar1e,pos={443,top},size={100,15},title="count (s)",value= gCtTime1
177//      SetVariable setvar1a,labelBack=(65535,32768,32768)
178//     
179//      pt += inc
180//      SetVariable setvar2a,pos={12,top+pt},size={100,15},title="theta min",value= gAngLow2
181//      SetVariable setvar2b,pos={119,top+pt},size={100,15},title="theta max",value= gAngHigh2
182//      SetVariable setvar2c,pos={227,top+pt},size={100,15},title="increm",value= gIncr2
183//      SetVariable setvar2d,pos={335,top+pt},size={100,15},title="# points",value= gNumPts2
184//      SetVariable setvar2e,pos={443,top+pt},size={100,15},title="count (s)",value= gCtTime2
185//      SetVariable setvar2a labelBack=(65535,65533,32768)
186//     
187//      pt += inc
188//      SetVariable setvar3a,pos={12,top+pt},size={100,15},title="theta min",value= gAngLow3
189//      SetVariable setvar3b,pos={119,top+pt},size={100,15},title="theta max",value= gAngHigh3
190//      SetVariable setvar3c,pos={227,top+pt},size={100,15},title="increm",value= gIncr3
191//      SetVariable setvar3d,pos={335,top+pt},size={100,15},title="# points",value= gNumPts3
192//      SetVariable setvar3e,pos={443,top+pt},size={100,15},title="count (s)",value= gCtTime3
193//      SetVariable setvar3a labelBack=(32769,65535,32768)
194//     
195//      pt += inc
196//      SetVariable setvar4a,pos={12,top+pt},size={100,15},title="theta min",value= gAngLow4
197//      SetVariable setvar4b,pos={119,top+pt},size={100,15},title="theta max",value= gAngHigh4
198//      SetVariable setvar4c,pos={227,top+pt},size={100,15},title="increm",value= gIncr4
199//      SetVariable setvar4d,pos={335,top+pt},size={100,15},title="# points",value= gNumPts4
200//      SetVariable setvar4e,pos={443,top+pt},size={100,15},title="count (s)",value= gCtTime4
201//      SetVariable setvar4a labelBack=(32768,65535,65535)
202//     
203//      pt += inc
204//      SetVariable setvar5a,pos={12,top+pt},size={100,15},title="theta min",value= gAngLow5
205//      SetVariable setvar5b,pos={119,top+pt},size={100,15},title="theta max",value= gAngHigh5
206//      SetVariable setvar5c,pos={227,top+pt},size={100,15},title="increm",value= gIncr5
207//      SetVariable setvar5d,pos={335,top+pt},size={100,15},title="# points",value= gNumPts5
208//      SetVariable setvar5e,pos={443,top+pt},size={100,15},title="count (s)",value= gCtTime5
209//      SetVariable setvar5a labelBack=(32768,54615,65535)
210//     
211//      pt += inc
212//      SetVariable setvar6a,pos={12,top+pt},size={100,15},title="theta min",value= gAngLow6
213//      SetVariable setvar6b,pos={119,top+pt},size={100,15},title="theta max",value= gAngHigh6
214//      SetVariable setvar6c,pos={227,top+pt},size={100,15},title="increm",value= gIncr6
215//      SetVariable setvar6d,pos={335,top+pt},size={100,15},title="# points",value= gNumPts6
216//      SetVariable setvar6e,pos={443,top+pt},size={100,15},title="count (s)",value= gCtTime6
217//      SetVariable setvar6a labelBack=(44253,29492,58982)
218//
219//// the action procedures and limits/increments
220//      SetVariable setvar1a proc=ThetaMinSetVarProc            //,limits={-2,0,0.1}
221//      SetVariable setvar2a proc=ThetaMinSetVarProc
222//      SetVariable setvar3a proc=ThetaMinSetVarProc
223//      SetVariable setvar4a proc=ThetaMinSetVarProc
224//      SetVariable setvar5a proc=ThetaMinSetVarProc
225//      SetVariable setvar6a proc=ThetaMinSetVarProc
226//
227////
228//      SetVariable setvar1b proc=ThetaMaxSetVarProc            //,limits={0.4,1,0.1}
229//      SetVariable setvar2b proc=ThetaMaxSetVarProc
230//      SetVariable setvar3b proc=ThetaMaxSetVarProc
231//      SetVariable setvar4b proc=ThetaMaxSetVarProc
232//      SetVariable setvar5b proc=ThetaMaxSetVarProc
233//      SetVariable setvar6b proc=ThetaMaxSetVarProc
234////
235//      SetVariable setvar1c proc=IncrSetVarProc,limits={0.01,0.1,0.01}
236//      SetVariable setvar2c proc=IncrSetVarProc,limits={0.02,0.2,0.02}
237//      SetVariable setvar3c proc=IncrSetVarProc,limits={0.05,0.4,0.05}
238//      SetVariable setvar4c proc=IncrSetVarProc,limits={0.1,1,0.1}
239//      SetVariable setvar5c proc=IncrSetVarProc,limits={0.5,5,1}
240//      SetVariable setvar6c proc=IncrSetVarProc,limits={1,10,2}
241////
242//      SetVariable setvar1d proc=NumPtsSetVarProc,limits={2,50,1}
243//      SetVariable setvar2d proc=NumPtsSetVarProc,limits={2,50,1}
244//      SetVariable setvar3d proc=NumPtsSetVarProc,limits={2,50,1}
245//      SetVariable setvar4d proc=NumPtsSetVarProc,limits={2,50,1}
246//      SetVariable setvar5d proc=NumPtsSetVarProc,limits={2,50,1}
247//      SetVariable setvar6d proc=NumPtsSetVarProc,limits={2,50,1}
248////
249//      SetVariable setvar1e proc=CtTimeSetVarProc,limits={-1,50000,1}
250//      SetVariable setvar2e proc=CtTimeSetVarProc,limits={-1,50000,10}
251//      SetVariable setvar3e proc=CtTimeSetVarProc,limits={-1,50000,10}
252//      SetVariable setvar4e proc=CtTimeSetVarProc,limits={-1,50000,30}
253//      SetVariable setvar5e proc=CtTimeSetVarProc,limits={-1,50000,100}
254//      SetVariable setvar6e proc=CtTimeSetVarProc,limits={-1,50000,100}
255//     
256//      Button button0,pos={255,180},size={60,20},fColor=(65535,65535,0),proc=U_SimPlotButtonProc,title="Plot"
257//      Button button1,pos={260,286},size={50,20},proc=U_SaveButtonProc,title="Save"
258//
259////checkbox for "easy" mode
260//      CheckBox check0 title="Simple mode?",pos={400,19},proc=EnterModeCheckProc,value=1
261//      ThetaEditMode(2)                //checked on startup
262//     
263////    instrument setup
264//      SetVariable U_setvar0_1,pos={20,211},size={160,15},title="Thickness (cm)"
265//      SetVariable U_setvar0_1,limits={0,inf,0.1},value= root:Packages:NIST:USANS:Globals:U_Sim:gThick
266//      SetVariable U_setvar0_3,pos={20,235},size={160,15},title="Sample Transmission"
267//      SetVariable U_setvar0_3,limits={0,1,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:gSamTrans
268//      PopupMenu U_popup0,pos={20,185},size={165,20},proc=Sim_USANS_ModelPopMenuProc,title="Model"
269//      PopupMenu U_popup0,mode=1,value= #"U_FunctionPopupList()"
270//      SetVariable setvar0,pos={20,259},size={120,15},title="Empty Level"
271//      SetVariable setvar0,limits={0,10,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:g_EmptyLevel
272//      SetVariable setvar0_1,pos={20,284},size={120,15},title="Bkg Level"
273//      SetVariable setvar0_1,limits={0,10,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:g_BkgLevel
274//     
275//      CheckBox check0_4 title="Show EMP?",pos={160,260},proc=ShowEMPCheckProc,value=0
276//     
277//      CheckBox check0_2,pos={253,239},size={60,14},title="CountRate?",variable= root:Packages:NIST:USANS:Globals:U_Sim:g_1D_PlotCR
278//      CheckBox check0_3,pos={262,264},size={60,14},title="Noise?",variable= root:Packages:NIST:USANS:Globals:U_Sim:g_1D_AddNoise
279//     
280//// a box for the results
281//      SetVariable totalTime,pos={338,185},size={150,15},title="Count time (h:m)",value= gTotTimeStr
282////    ValDisplay valdisp0,pos={338,210},size={220,13},title="Total detector counts"
283////    ValDisplay valdisp0,limits={0,0,0},barmisc={0,1000},value= root:Packages:NIST:USANS:Globals:U_Sim:g_1DTotCts
284//      ValDisplay valdisp0_2,pos={338,234},size={220,13},title="Fraction of beam scattered"
285//      ValDisplay valdisp0_2,limits={0,0,0},barmisc={0,1000},value= root:Packages:NIST:USANS:Globals:U_Sim:g_1DFracScatt
286//      ValDisplay valdisp0_3,pos={338,259},size={220,13},title="Estimated transmission"
287//      ValDisplay valdisp0_3,limits={0,0,0},barmisc={0,1000},value=root:Packages:NIST:USANS:Globals:U_Sim:g_1DEstTrans
288//
289//     
290//      SetDataFolder root:
291//
292//EndMacro
293
294// ??make the USANS simulation results into a graph and a separate panel.
295// the control bar at the top makes the whole thing too large, and control bars are limited
296// to 500 pix wide, which is really tight
297// left of 1055 is good for Mac, 955 is better for Win
298//
299Window UCALC_Panel()
300        PauseUpdate; Silent 1           // building window...
301        String platform=UpperStr(IgorInfo(2))
302        Variable pos=strsearch(platform,"WINDOWS",0)
303        if(pos >= 0)            //windows
304                Display /W=(55,44,855,450) /K=1
305        else            //mac
306                Display /W=(55,44,1055,544) /K=1
307        endif
308       
309        ModifyGraph cbRGB=(36929,50412,31845)
310        DoWindow/C UCALC
311        DoWindow/T UCALC,"USANS Simulation"
312        ControlBar/L 500
313       
314        GroupBox group0,pos={5,1},size={493,177},title="Instrument Setup"
315        GroupBox group1,pos={5,183},size={240,110},title="Sample Setup"
316        GroupBox group2,pos={5,310},size={240,170},title="Results"
317       
318        PopupMenu popup0,pos={17,19},size={165,20},title="Sample Aperture Diam (in)"
319        PopupMenu popup0,mode=3,popvalue="0.625",value="0.25;0.50;0.625;0.75;1.0;1.75;2.0;"
320        PopupMenu popup0,proc=UCALC_SampleAperturePopup
321        PopupMenu popup2,pos={220,19},size={165,20},title="Presets"
322        PopupMenu popup2,mode=3,popvalue="Long Count",value="Short Count;Medium Count;Long Count;"
323        PopupMenu popup2,proc=UCALC_PresetPopup
324
325        SetDataFolder root:Packages:NIST:USANS:Globals:U_Sim
326       
327        Variable top=46,pt=0,inc=18,left=0//left=533
328        SetVariable setvar1a,pos={left+17,top},size={90,15},title="theta min",value= gAngLow1
329        SetVariable setvar1b,pos={left+113,top},size={89,15},title="theta max",value= gAngHigh1
330        SetVariable setvar1c,pos={left+209,top},size={89,15},title="increm",value= gIncr1
331        SetVariable setvar1d,pos={left+299,top},size={100,15},title="# points",value= gNumPts1
332        SetVariable setvar1e,pos={left+399,top},size={93,15},title="count (s)",value= gCtTime1
333//      SetVariable setvar1a,labelBack=(65535,32768,32768)              //old rainbow
334        SetVariable setvar1a,labelBack=(49858,65535,65535)
335       
336        pt += inc
337        SetVariable setvar2a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow2
338        SetVariable setvar2b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh2
339        SetVariable setvar2c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr2
340        SetVariable setvar2d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts2
341        SetVariable setvar2e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime2
342//      SetVariable setvar2a labelBack=(65535,65533,32768)              //old rainbow
343        SetVariable setvar2a labelBack=(21074,8995,21074)
344       
345        pt += inc
346        SetVariable setvar3a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow3
347        SetVariable setvar3b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh3
348        SetVariable setvar3c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr3
349        SetVariable setvar3d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts3
350        SetVariable setvar3e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime3
351//      SetVariable setvar3a labelBack=(32769,65535,32768)              //old rainbow
352        SetVariable setvar3a labelBack=(0,60652,60652)
353       
354        pt += inc
355        SetVariable setvar4a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow4
356        SetVariable setvar4b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh4
357        SetVariable setvar4c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr4
358        SetVariable setvar4d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts4
359        SetVariable setvar4e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime4
360//      SetVariable setvar4a labelBack=(32768,65535,65535)              //old rainbow
361        SetVariable setvar4a labelBack=(0,51400,0)
362       
363        pt += inc
364        SetVariable setvar5a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow5
365        SetVariable setvar5b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh5
366        SetVariable setvar5c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr5
367        SetVariable setvar5d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts5
368        SetVariable setvar5e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime5
369//      SetVariable setvar5a labelBack=(32768,54615,65535)              //old rainbow
370        SetVariable setvar5a labelBack=(59367,49344,0)
371       
372        pt += inc
373        SetVariable setvar6a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow6
374        SetVariable setvar6b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh6
375        SetVariable setvar6c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr6
376        SetVariable setvar6d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts6
377        SetVariable setvar6e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime6
378//      SetVariable setvar6a labelBack=(44253,29492,58982)              //old rainbow
379        SetVariable setvar6a labelBack=(54998,0,0)
380
381        pt += inc
382        SetVariable setvar7a,pos={left+17,top+pt},size={90,15},title="theta min",value= gAngLow7
383        SetVariable setvar7b,pos={left+113,top+pt},size={89,15},title="theta max",value= gAngHigh7
384        SetVariable setvar7c,pos={left+209,top+pt},size={89,15},title="increm",value= gIncr7
385        SetVariable setvar7d,pos={left+299,top+pt},size={100,15},title="# points",value= gNumPts7
386        SetVariable setvar7e,pos={left+399,top+pt},size={93,15},title="count (s)",value= gCtTime7
387//      SetVariable setvar7a labelBack=(44253,29492,58982)              //old rainbow
388        SetVariable setvar7a labelBack=(39321,21845,51657)
389
390// the action procedures and limits/increments
391        SetVariable setvar1a proc=ThetaMinSetVarProc            //,limits={-2,0,0.1}
392        SetVariable setvar2a proc=ThetaMinSetVarProc
393        SetVariable setvar3a proc=ThetaMinSetVarProc
394        SetVariable setvar4a proc=ThetaMinSetVarProc
395        SetVariable setvar5a proc=ThetaMinSetVarProc
396        SetVariable setvar6a proc=ThetaMinSetVarProc
397        SetVariable setvar7a proc=ThetaMinSetVarProc
398
399//
400        SetVariable setvar1b proc=ThetaMaxSetVarProc            //,limits={0.4,1,0.1}
401        SetVariable setvar2b proc=ThetaMaxSetVarProc
402        SetVariable setvar3b proc=ThetaMaxSetVarProc
403        SetVariable setvar4b proc=ThetaMaxSetVarProc
404        SetVariable setvar5b proc=ThetaMaxSetVarProc
405        SetVariable setvar6b proc=ThetaMaxSetVarProc
406        SetVariable setvar7b proc=ThetaMaxSetVarProc
407//
408        SetVariable setvar1c proc=IncrSetVarProc,limits={0.01,0.1,0.01}
409        SetVariable setvar2c proc=IncrSetVarProc,limits={0.02,0.2,0.02}
410        SetVariable setvar3c proc=IncrSetVarProc,limits={0.05,0.4,0.05}
411        SetVariable setvar4c proc=IncrSetVarProc,limits={0.1,1,0.1}
412        SetVariable setvar5c proc=IncrSetVarProc,limits={0.5,5,1}
413        SetVariable setvar6c proc=IncrSetVarProc,limits={1,10,2}
414        SetVariable setvar7c proc=IncrSetVarProc,limits={1,20,2}
415//
416        SetVariable setvar1d proc=NumPtsSetVarProc,limits={2,50,1}
417        SetVariable setvar2d proc=NumPtsSetVarProc,limits={2,50,1}
418        SetVariable setvar3d proc=NumPtsSetVarProc,limits={2,50,1}
419        SetVariable setvar4d proc=NumPtsSetVarProc,limits={2,50,1}
420        SetVariable setvar5d proc=NumPtsSetVarProc,limits={2,50,1}
421        SetVariable setvar6d proc=NumPtsSetVarProc,limits={2,50,1}
422        SetVariable setvar7d proc=NumPtsSetVarProc,limits={2,50,1}
423//
424        SetVariable setvar1e proc=CtTimeSetVarProc,limits={-1,50000,1}
425        SetVariable setvar2e proc=CtTimeSetVarProc,limits={-1,50000,10}
426        SetVariable setvar3e proc=CtTimeSetVarProc,limits={-1,50000,10}
427        SetVariable setvar4e proc=CtTimeSetVarProc,limits={-1,50000,30}
428        SetVariable setvar5e proc=CtTimeSetVarProc,limits={-1,50000,100}
429        SetVariable setvar6e proc=CtTimeSetVarProc,limits={-1,50000,100}
430        SetVariable setvar7e proc=CtTimeSetVarProc,limits={-1,50000,100}
431       
432        Button button0,pos={left+280,200},size={130,20},fColor=(65535,65535,0),proc=U_SimPlotButtonProc,title="Simulate USANS"
433        CheckBox check0_2,pos={left+280,250},size={60,14},title="CountRate?",variable= root:Packages:NIST:USANS:Globals:U_Sim:g_1D_PlotCR
434        CheckBox check0_3,pos={left+280,270},size={60,14},title="Noise?",variable= root:Packages:NIST:USANS:Globals:U_Sim:g_1D_AddNoise
435        CheckBox check0_4 title="Show EMP?",pos={left+280,290},proc=ShowEMPCheckProc,value=0
436
437
438//checkbox for "easy" mode
439        CheckBox check0 title="Simple mode?",pos={left+400,19},proc=EnterModeCheckProc,value=1
440        ThetaEditMode(2)                //checked on startup
441       
442//      sample setup
443        SetVariable U_setvar0_1,pos={left+20,231},size={160,15},title="Thickness (cm)"
444        SetVariable U_setvar0_1,limits={0,inf,0.1},value= root:Packages:NIST:USANS:Globals:U_Sim:gThick
445        SetVariable U_setvar0_3,pos={left+20,255},size={160,15},title="Sample Transmission"
446        SetVariable U_setvar0_3,limits={0,1,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:gSamTrans
447        PopupMenu U_popup0,pos={left+20,205},size={165,20},proc=Sim_USANS_ModelPopMenuProc,title="Model"
448        PopupMenu U_popup0,mode=1,value= #"U_FunctionPopupList()"
449       
450        //unused as of yet, not sure I'll ever need these
451//      SetVariable setvar0,pos={left+20,279},size={120,15},title="Empty Level",disable=2
452//      SetVariable setvar0,limits={0,10,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:g_EmptyLevel
453//      SetVariable setvar0_1,pos={left+20,304},size={120,15},title="Bkg Level",disable=2
454//      SetVariable setvar0_1,limits={0,10,0.01},value= root:Packages:NIST:USANS:Globals:U_Sim:g_BkgLevel
455       
456       
457       
458// a box for the results
459        SetVariable totalTime,pos={left+20,335},size={150,15},title="Count time (h:m)",value= gTotTimeStr
460        SetVariable totalTime,limits={-inf,inf,0},noedit=1
461        SetVariable valdisp0_2,pos={left+20,365},size={200,15},title="Multiple Coherent Scattering"
462        SetVariable valdisp0_2,limits={-inf,inf,0},format="%8f",noedit=1,value=g_MultScattFraction
463//      ValDisplay valdisp0_3,pos={left+20,420},size={220,13},title="Estimated transmission"
464//      ValDisplay valdisp0_3,limits={0,0,0},barmisc={0,1000},value=root:Packages:NIST:USANS:Globals:U_Sim:g_1DEstTrans
465//      ValDisplay valdisp0_3,disable=1
466        Button button1,pos={left+20,390},size={150,20},proc=U_SavePanelProc,title="Save PNG"
467        Button button2,pos={left+20,420},size={150,20},proc=U_ConfigTextProc,title="Config Text"
468        Button button3,pos={left+20,450},size={150,20},proc=U_SaveButtonProc,title="Save Simulated Data"
469
470// help, done buttons
471        Button U_helpButton,pos={300,440},size={25,20},proc=showUCALCHelp,title="?"
472        Button U_helpButton,help={"Show help file for simulation of USANS Data"}
473        Button U_DoneButton,pos={350,440},size={50,20},proc=UCALCDoneButton,title="Done"
474        Button U_DoneButton,help={"This button will close the panel"}
475               
476        SetDataFolder root:
477
478EndMacro
479
480Proc UCALCDoneButton(ctrlName): ButtonControl
481        String ctrlName
482        DoWindow/K UCALC
483end
484
485Proc showUCALCHelp(ctrlName): ButtonControl
486        String ctrlName
487        DisplayHelpTopic/K=1/Z "UCALC"
488        if(V_flag !=0)
489                DoAlert 0,"The USANS Simulation Help file could not be found"
490        endif
491end
492
493// changing theta min - hold incr and #, result is new theta max
494Function ThetaMinSetVarProc(sva) : SetVariableControl
495        STRUCT WMSetVariableAction &sva
496
497        switch( sva.eventCode )
498                case 1: // mouse up
499                case 2: // Enter key
500                case 3: // Live update
501                        Variable ThetaMin = sva.dval
502                        String ns=CtrlNumber(sva.ctrlName)              //control number as a string
503                        NVAR incr = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gIncr"+ns)
504                        NVAR NumPts = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gNumPts"+ns)
505                        NVAR thetaMax = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngHigh"+ns)
506                       
507                        thetaMax = ThetaMin + incr*NumPts
508                       
509                        break
510        endswitch
511
512        return 0
513End
514
515
516// changing theta max - hold min and incr, result is new # of points
517// then need to recalculate the total counting time
518Function ThetaMaxSetVarProc(sva) : SetVariableControl
519        STRUCT WMSetVariableAction &sva
520
521        switch( sva.eventCode )
522                case 1: // mouse up
523                case 2: // Enter key
524                case 3: // Live update
525                        Variable ThetaMax = sva.dval
526                        String ns=CtrlNumber(sva.ctrlName)              //control number as a string
527                        NVAR thetaMin = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngLow"+ns)
528                        NVAR NumPts = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gNumPts"+ns)
529                        NVAR Incr = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gIncr"+ns)
530                       
531                        NumPts = trunc( (thetaMax - ThetaMin) / incr ) + 1                      //+1 to count both ends of the interval
532               
533                        CalcTotalCountTime()
534                       
535                        break
536        endswitch
537
538        return 0
539End
540
541// changing increment - hold min and #, result is new max
542Function IncrSetVarProc(sva) : SetVariableControl
543        STRUCT WMSetVariableAction &sva
544
545        switch( sva.eventCode )
546                case 1: // mouse up
547                case 2: // Enter key
548                case 3: // Live update
549                        Variable incr = sva.dval
550                        String ns=CtrlNumber(sva.ctrlName)              //control number as a string
551                        NVAR thetaMin = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngLow"+ns)
552                        NVAR NumPts = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gNumPts"+ns)
553                        NVAR thetaMax = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngHigh"+ns)
554                       
555                        thetaMax = ThetaMin + incr*NumPts
556               
557                        break
558        endswitch
559
560        return 0
561End
562
563// changing #pts - hold min and incr, result is new max
564// then recalculate the total count time
565Function NumPtsSetVarProc(sva) : SetVariableControl
566        STRUCT WMSetVariableAction &sva
567
568        switch( sva.eventCode )
569                case 1: // mouse up
570                case 2: // Enter key
571                case 3: // Live update
572                        Variable NumPts = sva.dval
573                        String ns=CtrlNumber(sva.ctrlName)              //control number as a string
574                        NVAR thetaMin = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngLow"+ns)
575                        NVAR incr = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gIncr"+ns)
576                        NVAR thetaMax = $("root:Packages:NIST:USANS:Globals:U_Sim:"+"gAngHigh"+ns)
577
578// old way, like ICP                   
579//                      thetaMax = ThetaMin + incr*NumPts
580
581// new way - to spread the points out over the specified angle range
582                        incr = (thetaMax - thetaMin) /( numpts -1 )                     //-1 since end points are both counted
583                       
584                        CalcTotalCountTime()
585                       
586                        break
587        endswitch
588
589        return 0
590End
591
592// changing count time -
593// then recalculate the total count time
594Function CtTimeSetVarProc(sva) : SetVariableControl
595        STRUCT WMSetVariableAction &sva
596
597        switch( sva.eventCode )
598                case 1: // mouse up
599                case 2: // Enter key
600                case 3: // Live update
601                        CalcTotalCountTime()
602                        break
603        endswitch
604
605        return 0
606End
607
608// hard-wired for 7 controls
609// return value is in seconds
610// global display string is set with hrs:min
611Function CalcTotalCountTime()
612
613        Variable ii,num,totTime=0
614        String pathStr="root:Packages:NIST:USANS:Globals:U_Sim:"
615        num=7
616       
617        for(ii=1;ii<=num;ii+=1)
618                NVAR ctTime = $(pathStr+"gCtTime"+num2str(ii))
619                NVAR numPts = $(pathStr+"gNumPts"+num2str(ii))
620                if(ctTime>0)
621                        totTime += ctTime*numPts
622                endif
623        endfor
624        Variable hrs,mins
625        hrs = trunc(totTime/3600)
626        mins = trunc(mod(totTime,3600)/60)
627//      Printf "Counting time (hr:min) = %d:%d\r",hrs,mins
628       
629        SVAR str = $(pathStr+"gTotTimeStr")
630        sprintf str,"%d:%d",hrs,mins
631       
632        return(totTime)
633End
634
635//returns the control number from the name string
636// all are setvarNa
637Function/S CtrlNumber(str)
638        String str
639       
640        return(str[6])
641End
642
643// changes edit mode of the theta min/max boxes and increment for simplified setup
644// val = 2 = disable
645// val = 0 = edit enabled
646Function ThetaEditMode(val)
647        Variable val
648       
649        SetVariable setvar1a,win=UCALC,disable=val
650        SetVariable setvar1b,win=UCALC,disable=val
651        SetVariable setvar1c,win=UCALC,disable=val
652       
653        SetVariable setvar2a,win=UCALC,disable=val
654        SetVariable setvar2b,win=UCALC,disable=val
655        SetVariable setvar2c,win=UCALC,disable=val
656       
657        SetVariable setvar3a,win=UCALC,disable=val
658        SetVariable setvar3b,win=UCALC,disable=val
659        SetVariable setvar3c,win=UCALC,disable=val
660       
661        SetVariable setvar4a,win=UCALC,disable=val
662        SetVariable setvar4b,win=UCALC,disable=val
663        SetVariable setvar4c,win=UCALC,disable=val
664       
665        SetVariable setvar5a,win=UCALC,disable=val
666        SetVariable setvar5b,win=UCALC,disable=val
667        SetVariable setvar5c,win=UCALC,disable=val
668       
669        SetVariable setvar6a,win=UCALC,disable=val
670        SetVariable setvar6b,win=UCALC,disable=val
671        SetVariable setvar6c,win=UCALC,disable=val
672       
673        SetVariable setvar7a,win=UCALC,disable=val
674        SetVariable setvar7b,win=UCALC,disable=val
675        SetVariable setvar7c,win=UCALC,disable=val
676        return(0)
677End
678
679
680
681Function EnterModeCheckProc(cba) : CheckBoxControl
682        STRUCT WMCheckboxAction &cba
683
684        switch( cba.eventCode )
685                case 2: // mouse up
686                        Variable checked = cba.checked
687                       
688                        if(checked)
689                                ThetaEditMode(2)
690                        else
691                                ThetaEditMode(0)
692                        endif
693                       
694                        break
695        endswitch
696
697        return 0
698End
699
700
701Function ShowEMPCheckProc(cba) : CheckBoxControl
702        STRUCT WMCheckboxAction &cba
703
704        switch( cba.eventCode )
705                case 2: // mouse up
706                        Variable checked = cba.checked
707
708                        String list,item,popStr,qval,CR
709                        Variable OK=1
710
711                        if(exists("root:Packages:NIST:USANS:Globals:Q_2p0")==0)                 //
712                                MakeUSANSEmptyWaves()
713                        endif
714                       
715                       
716                        if(checked)
717                                // put it on the graph
718                                SetDataFolder root:Packages:NIST:USANS:Globals
719                                ControlInfo/W=UCALC popup0
720                                popStr = S_Value
721                                strswitch(popStr)       // string switch
722                                        case "0.25":            // execute if case matches expression
723                                                qval = "Q_0p25"
724                                                CR = "CR_0p25"
725                                                break           
726                                        case "0.50":   
727                                                qval = "Q_0p50"
728                                                CR = "CR_0p50"
729                                                break   
730                                        case "0.625":           
731                                                qval = "Q_0p625"
732                                                CR = "CR_0p625"
733                                                break   
734                                        case "1.75":           
735                                                qval = "Q_1p75"
736                                                CR = "CR_1p75"
737                                                break   
738                                        case "2.0":             
739                                                qval = "Q_2p0"
740                                                CR = "CR_2p0"
741                                                break   
742                                        default:                                                        // optional default expression executed
743                                                OK=0
744                                endswitch
745                               
746                                if(OK)
747                                        AppendToGraph/W=UCALC $CR vs $qval
748                                        ModifyGraph marker=19,mode($CR)=4,msize($CR)=3,rgb($CR)=(0,0,0)
749                                endif
750                               
751                        else
752                                //take it off of the graph
753                                SetDataFolder root:Packages:NIST:USANS:Globals
754                                list=WaveList("CR*", ";", "WIN:UCALC")
755                                item=StringFromList(0, list ,";")               //should be one item
756                                if(strlen(item) != 0)
757                                        RemoveFromGraph/W=UCALC $item
758                                endif
759                        endif
760                       
761                        break
762        endswitch
763
764        SetDataFolder root:
765        return 0
766End
767
768
769
770// based on the angle ranges above (with non-zero count times)
771// plot where the data points would be
772//
773Function U_SimPlotButtonProc(ba) : ButtonControl
774        STRUCT WMButtonAction &ba
775
776       
777        switch( ba.eventCode )
778                case 2: // mouse up
779                        // click code here
780                       
781                        CalcUSANS()
782                       
783                        break
784        endswitch
785
786        return 0
787End
788
789
790
791// Fills a fake USANS folder with the concatenated ranges, and converts to Q
792// root:Packages:NIST:USANS:SIM
793//
794// does the work of calculating and smearing the simluated USANS data set
795Function CalcUSANS()   
796
797        Variable num,ii,firstSet=0
798        String pathStr="root:Packages:NIST:USANS:Globals:U_Sim:gCtTime"
799        String fromType="SWAP",toType="SIM"
800        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
801       
802        num = 7                 //# of angular ranges
803       
804        // only try to plot ranges with non-zero count times
805        firstSet=0
806        for(ii=1;ii<=num;ii+=1)
807                NVAR dum = $(pathStr+num2str(ii))
808                if(dum>0)                                               
809//                      print "CtTime = ",dum
810                        firstSet += 1
811                        LoadSimulatedAngleRange(ii,"SWAP")      //overwrite what's in the SWAP folder
812                       
813                        // did not do this step
814                        //Convert2Countrate("SWAP",1)
815                        //
816                       
817                        if(firstSet==1) //first time, overwrite
818                                NewDataWaves("SWAP","SIM")
819                                //plus my two waves
820                                Duplicate/O $(USANSFolder+":"+fromType+":countingTime"),$(USANSFolder+":"+toType+":countingTime")
821                                Duplicate/O $(USANSFolder+":"+fromType+":SetNumber"),$(USANSFolder+":"+toType+":SetNumber")
822
823                        else            //append to waves in "SIM"
824                                AppendDataWaves("SWAP","SIM")
825                                // and my two waves
826                                ConcatenateData( (USANSFolder+":"+toType+":countingTime"),(USANSFolder+":"+fromType+":countingTime") )
827                                ConcatenateData( (USANSFolder+":"+toType+":SetNumber"),(USANSFolder+":"+fromType+":SetNumber") )
828
829                        endif
830                       
831                endif
832        endfor
833
834
835        //sort after all loaded - not by angle, but by Q
836// Get rid of the negative angles - the smearing integration does not like these!
837// (may add them back in later, but probably not)
838        UDoAngleSort("SIM")
839       
840        ConvertAngle2Qvals("SIM",0)
841       
842       
843       
844        //fill the data with something
845        WAVE qvals = root:Packages:NIST:USANS:SIM:qvals
846        WAVE DetCts = root:Packages:NIST:USANS:SIM:DetCts
847       
848
849        //find the Trans Cts for T_Wide
850//      FindTWideCts("EMP")
851
852//generate a "fake" 1d data folder/set named "Sim_USANS" that can be used for smearing
853// DOES NOT calculate a matrix, but instead fills 4-5-6 columns with -dQv
854// so that the trapezoid rule is used
855//     
856
857        FakeUSANSDataFolder(qvals,DetCts,0.117,"Sim_USANS")
858
859        // now calculate the smearing... instead of the counts
860        SVAR funcStr = root:Packages:NIST:USANS:Globals:U_Sim:gFuncStr          //set by the popup     
861       
862        Wave inten = root:Sim_USANS:Sim_USANS_i
863        Wave sigave = root:Sim_USANS:Sim_USANS_s
864        Wave countingTime = root:Sim_USANS:countingTime         //counting time per point in the set
865
866        Duplicate/O inten root:Sim_USANS:Smeared_inten          // a place for the smeared result for probabilities
867        Wave Smeared_inten = root:Sim_USANS:Smeared_inten
868
869        String coefStr=""
870        Variable sig_sas=0,wavelength = 2.4
871        Variable Imon
872       
873
874        NVAR omega = root:Packages:NIST:USANS:Globals:U_Sim:gAnalyzerOmega
875       
876//      if(exists(funcStr) != 0)
877        if(cmpstr(funcStr,"default") != 0)
878                FUNCREF SANSModelAAO_proto func=$("fSmeared"+funcStr)                   //a wrapper for the structure version
879                FUNCREF SANSModelAAO_proto funcUnsmeared=$(funcStr)             //unsmeared
880                coefStr = MC_getFunctionCoef(funcStr)
881               
882                if(!MC_CheckFunctionAndCoef(funcStr,coefStr))
883                        Abort "Function and coefficients do not match. You must plot the unsmeared function before simulation."
884                endif
885               
886                // do the smearing calculation
887                func($coefStr,Smeared_inten,qvals)
888
889
890                NVAR thick = root:Packages:NIST:USANS:Globals:U_Sim:gThick
891                NVAR trans = root:Packages:NIST:USANS:Globals:U_Sim:gSamTrans
892                NVAR SimDetCts = root:Packages:NIST:USANS:Globals:U_Sim:g_1DTotCts                      //summed counts (simulated)
893                NVAR estDetCR = root:Packages:NIST:USANS:Globals:U_Sim:g_1DEstDetCR                     // estimated detector count rate
894                NVAR fracScat = root:Packages:NIST:USANS:Globals:U_Sim:g_1DFracScatt            // fraction of beam captured on detector
895                NVAR estTrans = root:Packages:NIST:USANS:Globals:U_Sim:g_1DEstTrans             // estimated transmission of sample
896//              NVAR SimCountTime = root:Packages:NIST:USANS:Globals:U_Sim:gCntTime             //counting time used for simulation
897               
898                Imon = GetUSANSBeamIntensity()                          //based on the aperture size, select the beam intensity
899                Print "imon=",imon
900               
901                // calculate the scattering cross section simply to be able to estimate the transmission
902                // unfortunately, this calculation is useless in the USANS range. So although it's calculated, the
903                // results are never reported.
904                // -- the main issue is that the integration range is optimized for SANS, and is not useful for USANS
905                // -there are only a few points in the USANS range, and extending it lower caused issues before.
906                CalculateRandomDeviate(funcUnsmeared,$coefStr,wavelength,"root:Packages:NIST:SAS:ran_dev",sig_sas)
907
908                if(sig_sas > 100)
909                        DoAlert 0,"SAS cross section > 100. Estimates of multiple scattering are unreliable. Choosing a model with a well-defined Rg may help"
910                endif                   
911//              if(sig_sas > 100)
912//                      sprintf abortStr,"sig_sas = %g. Please check that the model coefficients have a zero background, or the low q is well-behaved.",sig_sas
913//              endif
914                // calculate the multiple scattering fraction for display (10/2009)
915                NVAR mScat = root:Packages:NIST:USANS:Globals:U_Sim:g_MultScattFraction
916                Variable nMax=10,tau
917                mScat=0
918                tau = thick*sig_sas
919                // this sums the normalized scattering P', so the result is the fraction of multiply coherently scattered
920                // neutrons out of those that were scattered
921                for(ii=2;ii<nMax;ii+=1)
922                        mScat += tau^(ii)/factorial(ii)
923        //              print tau^(ii)/factorial(ii)
924                endfor
925                estTrans = exp(-1*thick*sig_sas)                //thickness and sigma both in units of cm
926                mscat *= (estTrans)/(1-estTrans)
927               
928                if(mScat > 0.1)
929                        SetVariable valdisp0_2 win=UCALC,labelBack=(65535,32768,32768)
930                else
931                        SetVariable valdisp0_2 win=UCALC,labelBack=0
932                endif
933
934
935                Print "Sig_sas = ",sig_sas
936                Print "Estimated Transmission = ",estTrans
937               
938                Duplicate/O qvals prob_i
939                                       
940                prob_i = trans*thick*omega*Smeared_inten                        //probability of a neutron in q-bin(i)
941               
942//              Variable P_on = sum(prob_i,-inf,inf)
943//              Print "P_on = ",P_on
944                fracScat = 1-estTrans
945               
946                inten = (Imon*countingTime)*prob_i
947               
948                // do I round to an integer?
949                inten = round(inten)
950
951//              SimDetCts = sum(inten,-inf,inf)
952//              estDetCR = SimDetCts/SimCountTime
953               
954               
955                NVAR doABS = root:Packages:NIST:USANS:Globals:U_Sim:g_1D_DoABS
956                NVAR plotCR = root:Packages:NIST:USANS:Globals:U_Sim:g_1D_PlotCR
957                NVAR addNoise = root:Packages:NIST:USANS:Globals:U_Sim:g_1D_AddNoise
958                                       
959                sigave = sqrt(inten)            // assuming that N is large
960               
961                // add in random error in aveint based on the sigave
962                if(addNoise)
963                        inten += gnoise(sigave)
964                       
965                        //round to an integer again
966                        inten = round(inten)           
967                endif
968
969                // convert to absolute scale? Maybe not needed
970                // does nothing yet - need Ipeak, Twide
971//              if(doABS)
972//                      Variable kappa = thick*omega*trans*iMon*ctTime
973//                      inten /= kappa
974//                      inten /= kappa
975//              endif
976
977                // plot as countrate - maybe easier to visualize, and all of the data overlaps
978                if(plotCR)
979                        inten /= countingTime
980                        sigave /= countingTime
981                endif
982               
983                GraphSIM()
984
985        else
986                //no function plotted, no simulation can be done
987                DoAlert 0,"No function is selected or plotted, so no simulation is done. The default power law function is used."
988
989                inten = U_Power_Law_Model(1e-6,3,0,qvals)
990//              inten = U_SphereForm(1,9000,6e-6,0,qvals)               
991       
992                GraphSIM()
993
994        endif
995       
996
997end
998
999
1000//sort the data in the "type"folder, based on angle
1001//carry along all associated waves
1002//
1003// ---a duplicate of DoAngleSort(), by modified to
1004// include counting time and setNumber
1005//
1006// also trims the beginning of each data set so that it does not include any negative or zero angles
1007//
1008Function UDoAngleSort(type)
1009        String type
1010       
1011        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
1012       
1013        Wave angle = $(USANSFolder+":"+Type+":Angle")
1014        Wave detCts = $(USANSFolder+":"+Type+":DetCts")
1015        Wave ErrdetCts = $(USANSFolder+":"+Type+":ErrDetCts")
1016        Wave MonCts = $(USANSFolder+":"+Type+":MonCts")
1017        Wave TransCts = $(USANSFolder+":"+Type+":TransCts")
1018        Wave countingTime = $(USANSFolder+":"+Type+":countingTime")
1019        Wave SetNumber = $(USANSFolder+":"+Type+":SetNumber")
1020       
1021        Sort Angle DetCts,ErrDetCts,MonCts,TransCts,Angle,countingTime,SetNumber
1022       
1023        Variable ii,num,numBad,ang,val
1024        num=numpnts(angle)
1025        ii=0
1026        numBad=0
1027        val = 0         //cutoff value
1028        do
1029                ang = angle[ii]
1030                if(ang <= val)
1031                        numBad += 1
1032                else            //keep the points
1033                        Angle[ii-numBad] = ang
1034                        DetCts[ii-numBad] = DetCts[ii]
1035                        ErrDetCts[ii-numBad] = ErrDetCts[ii]
1036                        MonCts[ii-numBad] = MonCts[ii]
1037                        TransCts[ii-numBad] = TransCts[ii]
1038                        countingTime[ii-numBad] = countingTime[ii]
1039                        SetNumber[ii-numBad] = SetNumber[ii]
1040                endif
1041                ii += 1
1042        while(ii<num)
1043        //trim the end of the waves
1044        DeletePoints num-numBad, numBad, DetCts,ErrDetCts,MonCts,TransCts,Angle,countingTime,SetNumber
1045       
1046       
1047        return(0)
1048End
1049
1050// a simple default function
1051//
1052Function U_SphereForm(scale,radius,delrho,bkg,x)                               
1053        Variable scale,radius,delrho,bkg
1054        Variable x
1055       
1056        // variables are:                                                       
1057        //[0] scale
1058        //[1] radius (A)
1059        //[2] delrho (A-2)
1060        //[3] background (cm-1)
1061
1062       
1063        // calculates scale * f^2/Vol where f=Vol*3*delrho*((sin(qr)-qrcos(qr))/qr^3
1064        // and is rescaled to give [=] cm^-1
1065       
1066        Variable bes,f,vol,f2
1067        //
1068        //handle q==0 separately
1069        If(x==0)
1070                f = 4/3*pi*radius^3*delrho*delrho*scale*1e8 + bkg
1071                return(f)
1072        Endif
1073       
1074        bes = 3*(sin(x*radius)-x*radius*cos(x*radius))/x^3/radius^3
1075        vol = 4*pi/3*radius^3
1076        f = vol*bes*delrho              // [=] A
1077        // normalize to single particle volume, convert to 1/cm
1078        f2 = f * f / vol * 1.0e8                // [=] 1/cm
1079       
1080        return (scale*f2+bkg)   // Scale, then add in the background
1081       
1082End
1083
1084// better default function
1085Function U_Power_Law_Model(A,m,bgd,x) : FitFunc
1086        Variable A, m,bgd,x
1087//       Input (fitting) variables are:
1088        //[0] Coefficient
1089        //[1] (-) Power
1090        //[2] incoherent background
1091       
1092//      local variables
1093        Variable inten, qval
1094//      x is the q-value for the calculation
1095        qval = x
1096//      do the calculation and return the function value
1097       
1098        inten = A*qval^-m + bgd
1099        Return (inten)
1100End
1101
1102
1103// mimics LoadBT5File
1104// creates two other waves to identify the set and the counting time for that set
1105//
1106Function LoadSimulatedAngleRange(set,type)
1107        Variable set
1108        String type
1109
1110        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
1111       
1112        String s=num2str(set)
1113        NVAR angLow = $("root:Packages:NIST:USANS:Globals:U_Sim:gAngLow"+s)
1114        NVAR angHigh = $("root:Packages:NIST:USANS:Globals:U_Sim:gAngHigh"+s)
1115        NVAR numPts = $("root:Packages:NIST:USANS:Globals:U_Sim:gNumPts"+s)
1116        NVAR ctTime = $("root:Packages:NIST:USANS:Globals:U_Sim:gCtTime"+s)
1117        NVAR incr = $("root:Packages:NIST:USANS:Globals:U_Sim:gIncr"+s)
1118       
1119        Variable ii, err=0
1120       
1121        // generate q-points based on angular range from panel
1122       
1123        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":Angle")
1124        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":DetCts")
1125        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":ErrDetCts")
1126        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":MonCts")
1127        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":TransCts")
1128        //
1129        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":CountingTime")
1130        Make/O/D/N=(numPts) $(USANSFolder+":"+type+":SetNumber")
1131       
1132        Wave Angle = $(USANSFolder+":"+type+":Angle")
1133        Wave DetCts = $(USANSFolder+":"+type+":DetCts")
1134        Wave ErrDetCts = $(USANSFolder+":"+type+":ErrDetCts")
1135        Wave MonCts = $(USANSFolder+":"+type+":MonCts")
1136        Wave TransCts = $(USANSFolder+":"+type+":TransCts")
1137        Wave countingTime = $(USANSFolder+":"+type+":countingTime")
1138        Wave SetNumber = $(USANSFolder+":"+type+":SetNumber")
1139       
1140        countingTime = ctTime
1141        SetNumber = set
1142       
1143        for(ii=0;ii<numPts;ii+=1)
1144                Angle[ii] = angLow + ii*incr
1145                DetCts[ii] = set
1146        endfor
1147       
1148        //set the wave note for the DetCts
1149        String str=""
1150        str = "FILE:Sim "+s+";"
1151        str += "TIMEPT:"+num2str(ctTime)+";"
1152        str += "PEAKANG:0;"             //no value yet
1153        str += "STATUS:;"               //no value yet
1154        str += "LABEL:SimSet "+s+";"
1155        str += "PEAKVAL:;"              //no value yet
1156        str += "TWIDE:0;"               //no value yet
1157        Note DetCts,str
1158       
1159        return err                      // Zero signifies no error.     
1160End
1161
1162// add SIM data to the graph if it exists and is not already on the graph
1163//
1164//
1165// ** currently, I have changed the graph to not display the angle at the top. I have not yet
1166// found a way to keep the scaling of the two axes in-sync when the empty cell data is added to the
1167// graph (versus Q)
1168//
1169Function GraphSIM()
1170
1171//      SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
1172//      SetDataFolder $(USANSFolder+":SIM")
1173
1174        SetDataFolder root:Sim_USANS
1175        //is it already on the graph?
1176        String list=""
1177        list = Wavelist("Sim_USANS*",";","WIN:UCALC")
1178       
1179        if(strlen(list)!=0)
1180//              Print "SIM already on graph"
1181                SetDataFolder root:
1182                return(0)
1183        Endif
1184       
1185        //append the data if it exists
1186        If(waveExists($"Sim_USANS_i")==1)
1187                DoWindow/F UCALC
1188
1189// lines for the no-noise result vs angle
1190//              AppendToGraph/T Smeared_inten vs angle
1191                AppendToGraph Smeared_inten vs Sim_USANS_q
1192                ModifyGraph rgb(Smeared_inten)=(1,12815,52428)
1193                ModifyGraph mode(Smeared_inten)=4,marker(Smeared_inten)=19,msize(Smeared_inten)=2
1194                ModifyGraph tickUnit=1
1195
1196// colored points for the simulation with noise on top
1197                AppendToGraph Sim_USANS_i vs Sim_USANS_q
1198                ModifyGraph mode(Sim_USANS_i)=3,marker(Sim_USANS_i)=19,msize(Sim_USANS_i)=4
1199                //don't reverse the dbZ21 or the highest angle will have "invisible" light blue error bars
1200                ModifyGraph zColor(Sim_USANS_i)={setNumber,1,7,dBZ21,0}                         //better for 7 colors
1201//              ModifyGraph zColor(Sim_USANS_i)={setNumber,1,7,Rainbow,0}                       //force the colors from 1->7
1202                ModifyGraph useMrkStrokeRGB(Sim_USANS_i)=1
1203                ErrorBars/T=0 Sim_USANS_i Y,wave=(Sim_USANS_s,Sim_USANS_s)
1204               
1205                ModifyGraph log=1
1206                ModifyGraph mirror(left)=1
1207                // if no top axis, then
1208//              ModifyGraph mirror(bottom)=1
1209                ModifyGraph grid=2
1210                ModifyGraph standoff=0
1211               
1212                // to make sure that the scales are the same (but fails on zoom)
1213                NewFreeAxis/O/T top_angle
1214                ModifyFreeAxis/Z top_angle,master= bottom,hook= TransformAngleAxisHook
1215                ModifyGraph lblPos(top_angle)=50,freePos(top_angle)=0
1216                ModifyGraph log(top_angle)=1
1217               
1218                SetDrawEnv linefgc= (39321,1,1),dash= 3,linethick= 3
1219                SetDrawEnv xcoord= bottom
1220                DrawLine 3e-05,0.01,3e-05,0.99
1221               
1222                Label top_angle "Angle"
1223                Label bottom "Q (1/A)"
1224                Label left "Counts or Count Rate"
1225               
1226                Legend
1227        endif
1228       
1229        SetDataFolder root:
1230End
1231
1232Function TransformAngleAxisHook(s)
1233        STRUCT WMAxisHookStruct &s
1234
1235        s.max= s.max/5.55e-5
1236        s.min= s.min/5.55e-5
1237       
1238        return 0
1239End
1240
1241//fakes a folder with loaded 1-d usans data, no calculation of the matrix
1242Function        FakeUSANSDataFolder(qval,aveint,dqv,dataFolder)
1243        WAVE qval,aveint
1244        Variable dqv
1245        String dataFolder
1246       
1247
1248        String baseStr=dataFolder
1249        if(DataFolderExists("root:"+baseStr))
1250                SetDataFolder $("root:"+baseStr)
1251        else
1252                NewDataFolder/S $("root:"+baseStr)
1253        endif
1254
1255        ////overwrite the existing data, if it exists
1256        Duplicate/O qval, $(baseStr+"_q")
1257        Duplicate/O aveint, $(baseStr+"_i")
1258       
1259        Duplicate/O qval, $(baseStr+"_s")
1260        Wave sigave = $(baseStr+"_s")
1261
1262        sigave = sqrt(aveint)
1263
1264        // make a resolution matrix for SANS data
1265        Variable np=numpnts(qval)
1266        Make/D/O/N=(np,4) $(baseStr+"_res")
1267        Wave res=$(baseStr+"_res")
1268       
1269        res[][0] = -dQv         //sigQ
1270        res[][1] = -dQv         //qBar
1271        res[][2] = -dQv         //fShad
1272        res[][3] = qval[p]              //Qvalues
1273       
1274        // extra waves of set number and counting time for the simulation
1275        WAVE ctW = root:Packages:NIST:USANS:SIM:countingTime
1276        WAVE setW = root:Packages:NIST:USANS:SIM:setNumber
1277        WAVE ang = root:Packages:NIST:USANS:SIM:Angle
1278        Duplicate/O ctW countingTime
1279        Duplicate/O setW setNumber
1280        Duplicate/O ang Angle
1281
1282
1283        //clean up             
1284        SetDataFolder root:
1285       
1286End
1287
1288
1289Function/S U_FunctionPopupList()
1290        String list,tmp
1291        list = User_FunctionPopupList()
1292       
1293        //simplify the display, forcing smeared calculations behind the scenes
1294        tmp = FunctionList("Smear*",";","NPARAMS:1")    //smeared dependency calculations
1295        list = RemoveFromList(tmp, list ,";")
1296
1297
1298        if(strlen(list)==0)
1299                list = "No functions plotted"
1300        endif
1301       
1302        list = SortList(list)
1303       
1304        list = "default;"+ list
1305        return(list)
1306End     
1307
1308Function Sim_USANS_ModelPopMenuProc(pa) : PopupMenuControl
1309        STRUCT WMPopupAction &pa
1310
1311        switch( pa.eventCode )
1312                case 2: // mouse up
1313                        Variable popNum = pa.popNum
1314                        String popStr = pa.popStr
1315                        SVAR gStr = root:Packages:NIST:USANS:Globals:U_Sim:gFuncStr
1316                        gStr = popStr
1317                       
1318                        break
1319        endswitch
1320
1321        return 0
1322End         
1323
1324// if the sample aperture is changed, AND the empty data is displayed, change to the proper data
1325Function UCALC_SampleAperturePopup(pa) : PopupMenuControl
1326        STRUCT WMPopupAction &pa
1327
1328        switch( pa.eventCode )
1329                case 2: // mouse up
1330                        Variable popNum = pa.popNum
1331                        String popStr = pa.popStr
1332
1333                        ControlInfo/W=UCALC check0_4
1334                        if(V_Value==1)          //currently checked, need to update
1335                               
1336                                STRUCT WMCheckboxAction cba
1337                                cba.checked=0           //"un-check"
1338                                cba.eventCode=2
1339                                ShowEMPCheckProc(cba)
1340                               
1341                                cba.checked=1           //"re-check"
1342                                ShowEMPCheckProc(cba)
1343                       
1344                        endif
1345                        break
1346                       
1347        endswitch
1348       
1349        return 0
1350End 
1351
1352
1353Function UCALC_PresetPopup(pa) : PopupMenuControl
1354        STRUCT WMPopupAction &pa
1355
1356        switch( pa.eventCode )
1357                case 2: // mouse up
1358                        Variable popNum = pa.popNum
1359                        String popStr = pa.popStr
1360
1361                        SetDataFolder root:Packages:NIST:USANS:Globals:U_Sim
1362                       
1363                        NVAR gAngLow1 = gAngLow1
1364                        NVAR gAngHigh1 = gAngHigh1
1365                        NVAR gNumPts1 = gNumPts1
1366                        NVAR gCtTime1 = gCtTime1
1367                        NVAR gIncr1 = gIncr1
1368                       
1369                        NVAR gAngLow2 = gAngLow2
1370                        NVAR gAngHigh2 = gAngHigh2
1371                        NVAR gNumPts2 = gNumPts2
1372                        NVAR gCtTime2 = gCtTime2
1373                        NVAR gIncr2 = gIncr2
1374                       
1375                        NVAR gAngLow3 = gAngLow3
1376                        NVAR gAngHigh3 = gAngHigh3
1377                        NVAR gNumPts3 = gNumPts3
1378                        NVAR gCtTime3 = gCtTime3
1379                        NVAR gIncr3 = gIncr3
1380                       
1381                        NVAR gAngLow4 = gAngLow4
1382                        NVAR gAngHigh4 = gAngHigh4
1383                        NVAR gNumPts4 = gNumPts4
1384                        NVAR gCtTime4 = gCtTime4
1385                        NVAR gIncr4 = gIncr4
1386                       
1387                        NVAR gAngLow5 = gAngLow5
1388                        NVAR gAngHigh5 = gAngHigh5
1389                        NVAR gNumPts5 = gNumPts5
1390                        NVAR gCtTime5 = gCtTime5
1391                        NVAR gIncr5 = gIncr5
1392                       
1393                        NVAR gAngLow6 = gAngLow6
1394                        NVAR gAngHigh6 = gAngHigh6
1395                        NVAR gNumPts6 = gNumPts6
1396                        NVAR gCtTime6 = gCtTime6
1397                        NVAR gIncr6 = gIncr6
1398                       
1399                        NVAR gAngLow7 = gAngLow7
1400                        NVAR gAngHigh7 = gAngHigh7
1401                        NVAR gNumPts7 = gNumPts7
1402                        NVAR gCtTime7 = gCtTime7
1403                        NVAR gIncr7 = gIncr7
1404                       
1405                        strswitch(popStr)       // string switch
1406                                case "Short Count":             // execute if case matches expression
1407                                       
1408                                        gAngLow1 = -1
1409                                        gAngHigh1 = 0.6
1410                                        gNumPts1 = 33
1411                                        gCtTime1 = 10
1412                                        gIncr1 = 0.05   
1413                                       
1414                                        gAngLow2 = 0.7
1415                                        gAngHigh2 = 1.9
1416                                        gNumPts2 = 6
1417                                        gCtTime2 = 60
1418                                        gIncr2 = 0.2   
1419                                       
1420                                        gAngLow3 = 2
1421                                        gAngHigh3 = 4.8
1422                                        gNumPts3 = 7
1423                                        gCtTime3 = 120
1424                                        gIncr3 = 0.4   
1425                                       
1426                                        gAngLow4 = 5
1427                                        gAngHigh4 = 9.5
1428                                        gNumPts4 = 5
1429                                        gCtTime4 = 180
1430                                        gIncr4 = 0.9   
1431                                       
1432                                        gAngLow5 = 10
1433                                        gAngHigh5 = 19
1434                                        gNumPts5 = 5
1435                                        gCtTime5 = 240
1436                                        gIncr5 = 1.8
1437                                       
1438                                        gAngLow6 = 20
1439                                        gAngHigh6 = 48
1440                                        gNumPts6 = 15
1441                                        gCtTime6 = 0
1442                                        gIncr6 = 2
1443                                       
1444                                        gAngLow7 = 50
1445                                        gAngHigh7 = 95
1446                                        gNumPts7 = 10
1447                                        gCtTime7 = 0
1448                                        gIncr7 = 5
1449                                       
1450                                        break                                           
1451                                case "Medium Count":   
1452                       
1453                                        gAngLow1 = -1
1454                                        gAngHigh1 = 0.6
1455                                        gNumPts1 = 33
1456                                        gCtTime1 = 10
1457                                        gIncr1 = 0.05   
1458                                       
1459                                        gAngLow2 = 0.7
1460                                        gAngHigh2 = 1.9
1461                                        gNumPts2 = 13
1462                                        gCtTime2 = 60
1463                                        gIncr2 = 0.1   
1464                                       
1465                                        gAngLow3 = 2
1466                                        gAngHigh3 = 4.8
1467                                        gNumPts3 = 15
1468                                        gCtTime3 = 120
1469                                        gIncr3 = 0.2   
1470                                       
1471                                        gAngLow4 = 5
1472                                        gAngHigh4 = 9.5
1473                                        gNumPts4 = 10
1474                                        gCtTime4 = 300
1475                                        gIncr4 = 0.5   
1476                                       
1477                                        gAngLow5 = 10
1478                                        gAngHigh5 = 19
1479                                        gNumPts5 = 10
1480                                        gCtTime5 = 600
1481                                        gIncr5 = 1
1482                                       
1483                                        gAngLow6 = 20
1484                                        gAngHigh6 = 48
1485                                        gNumPts6 = 15
1486                                        gCtTime6 = 1200
1487                                        gIncr6 = 2     
1488                                       
1489                                        gAngLow7 = 50
1490                                        gAngHigh7 = 95
1491                                        gNumPts7 = 10
1492                                        gCtTime7 = 0
1493                                        gIncr7 = 5
1494                                                               
1495                                        break
1496                                case "Long Count":     
1497                                       
1498                                        gAngLow1 = -1
1499                                        gAngHigh1 = 0.6
1500                                        gNumPts1 = 33
1501                                        gCtTime1 = 25
1502                                        gIncr1 = 0.05   
1503                                       
1504                                        gAngLow2 = 0.7
1505                                        gAngHigh2 = 1.9
1506                                        gNumPts2 = 13
1507                                        gCtTime2 = 100
1508                                        gIncr2 = 0.1   
1509                                       
1510                                        gAngLow3 = 2
1511                                        gAngHigh3 = 4.8
1512                                        gNumPts3 = 15
1513                                        gCtTime3 = 300
1514                                        gIncr3 = 0.2   
1515                                       
1516                                        gAngLow4 = 5
1517                                        gAngHigh4 = 9.5
1518                                        gNumPts4 = 10
1519                                        gCtTime4 = 600
1520                                        gIncr4 = 0.5   
1521                                       
1522                                        gAngLow5 = 10
1523                                        gAngHigh5 = 19
1524                                        gNumPts5 = 10
1525                                        gCtTime5 = 1200
1526                                        gIncr5 = 1
1527                                       
1528                                        gAngLow6 = 20
1529                                        gAngHigh6 = 48
1530                                        gNumPts6 = 15
1531                                        gCtTime6 = 2000
1532                                        gIncr6 = 2
1533                                       
1534                                        gAngLow7 = 50
1535                                        gAngHigh7 = 95
1536                                        gNumPts7 = 10
1537                                        gCtTime7 = 3000
1538                                        gIncr7 = 5
1539                                       
1540                                        break
1541                                default:                       
1542                        endswitch
1543                       
1544                        break
1545        endswitch
1546
1547        //update the count time
1548        CalcTotalCountTime()
1549
1550        SetDataFolder root:
1551       
1552        return 0
1553End     
1554
1555
1556// return the beam intensity based on the sample aperture diameter
1557//
1558// based on the equation in John's instrument paper
1559//
1560Function GetUSANSBeamIntensity()       
1561
1562        String popStr
1563        Variable flux,diam,rad
1564       
1565        ControlInfo/W=UCALC popup0
1566        popStr = S_Value
1567        diam=str2num(popStr)            //in inches
1568        rad = diam/2*25.4                       //radius in mm
1569       
1570        flux = 662*rad*rad-39.9*rad*rad*rad+0.70*rad*rad*rad*rad
1571       
1572        return(flux)
1573End
1574
1575// based on the angle ranges above (with non-zero count times)
1576// save fake data points into a fake BT5 data file
1577//
1578Function U_SaveButtonProc(ba) : ButtonControl
1579        STRUCT WMButtonAction &ba
1580
1581       
1582        switch( ba.eventCode )
1583                case 2: // mouse up
1584                        // click code here
1585                       
1586                        Variable num,ii,baseNumber=101,firstSet=0
1587                        String pathStr="root:Packages:NIST:USANS:Globals:U_Sim:gCtTime"
1588                        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
1589                        String baseName="SIMUL"
1590                       
1591                        num = 7                 //# of angular ranges
1592                       
1593                        // only try to plot ranges with non-zero count times
1594                        firstSet=0
1595                        for(ii=1;ii<=num;ii+=1)
1596                                NVAR dum = $(pathStr+num2str(ii))
1597                                if(dum>0)                                               
1598                                        firstSet += 1
1599                                        if(firstSet==1) //first time, ask for base name
1600                                                Prompt baseName,"Enter a base name for the files"
1601                                                Prompt baseNumber,"Enter the starting index"
1602                                                DoPrompt "Enter Information for Save",baseName,baseNumber
1603                                               
1604                                                if(V_Flag==1)           //user canceled
1605                                                        return(1)
1606                                                endif
1607                                        endif
1608                                               
1609                                        SaveFakeUSANS(baseName,baseNumber-1,ii)
1610                                       
1611                                endif
1612                        endfor
1613                       
1614                        // no save the fake "COR" file that is the combination and rescaling of the data sets
1615                        // - the combination is in root:Sim_USANS:
1616                        NVAR g_1D_PlotCR = root:Packages:NIST:USANS:Globals:U_Sim:g_1D_PlotCR
1617
1618                        if(g_1D_PlotCR) // right now, save COR only if the CR is plotted
1619                                SaveFakeCOR(baseName)                   
1620                        endif
1621                       
1622                        break
1623        endswitch
1624
1625        return 0
1626End
1627
1628
1629Function SaveFakeCOR(nameStr)
1630        String nameStr
1631       
1632        NVAR gThick = root:Packages:NIST:USANS:Globals:U_Sim:gThick
1633        NVAR gSamTrans = root:Packages:NIST:USANS:Globals:U_Sim:gSamTrans
1634        NVAR gAnalyzerOmega = root:Packages:NIST:USANS:Globals:U_Sim:gAnalyzerOmega
1635        NVAR g_MultScattFraction = root:Packages:NIST:USANS:Globals:U_Sim:g_MultScattFraction
1636        SVAR funcStr = root:Packages:NIST:USANS:Globals:U_Sim:gFuncStr
1637       
1638        String folder = "root:Sim_USANS:",fileStr=""
1639        String termStr="\r\n"           //VAX uses only <CR> as terminator, but only CRLF seems to FTP correctly to VAX
1640        String destStr="",formatStr = "%15.6g %15.6g %15.6g %15.6g %15.6g %15.6g"+termStr
1641
1642        Variable refNum,integer,realval
1643        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
1644        WAVE qvals =$(folder + "Sim_USANS_q")
1645        WAVE inten=$(folder + "Sim_USANS_i")
1646        WAVE sig=$(folder + "Sim_USANS_s")
1647       
1648        //check each wave
1649        If(!(WaveExists(qvals)))
1650                Abort "qvals DNExist in SaveFakeCOR()"
1651        Endif
1652        If(!(WaveExists(inten)))
1653                Abort "inten DNExist in SaveFakeCOR()"
1654        Endif
1655        If(!(WaveExists(sig)))
1656                Abort "sig DNExist in SaveFakeCOR()"
1657        Endif
1658       
1659        //dummy wave for the resolution
1660        Duplicate/O qvals,dumWave
1661        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
1662        NVAR DQv=$(USANSFolder+":Globals:MainPanel:gDQv")
1663        dumWave = - DQv
1664       
1665        //temp waves for the rescaling
1666        Duplicate/O qvals,tq,ti,te
1667        ti=inten
1668        te=sig
1669        tq=qvals
1670
1671        //rescale
1672        ti *= 1.0e6 / GetUSANSBeamIntensity()           //converts to "SAM" scaled data
1673        ti *= 1/(gSamTrans*gThick*gAnalyzerOmega*1.0e6)         //converts to "COR" scaling (approx, since I don't know Twide)
1674        //rescale error appropriately too
1675        te *= 1.0e6 / GetUSANSBeamIntensity()           //converts to "SAM" scaled data
1676        te *= 1/(gSamTrans*gThick*gAnalyzerOmega*1.0e6)         //converts to "COR" scaling (approx, since I don't know Twide)
1677
1678        // add fake error to every few data points to mark it as simulated data
1679        te[2,numpnts(te)-1;7] *= 10
1680       
1681        // now trim off the q-values less than 3e-5
1682        if(tq[0] < 3e-5)
1683                do
1684                        DeletePoints 0, 1, ti,tq,te,dumWave
1685                while(tq[0] < 3e-5)
1686        endif
1687       
1688        String samStr="",empStr="",dateStr="",samLabelStr="",paramStr="",empLevStr="",bkgLevStr="",pkStr=""
1689        samStr = "SIM FILES: "+StringByKey("FILE",note(inten),":",";")
1690        empStr = "EMP FILES: "+"none"
1691        empLevStr = "EMP LEVEL: " + "none"
1692        bkgLevStr = "BKG LEVEL: " + "not used"
1693        paramStr = "Ds = "+num2str(gthick)+" cm ; "
1694        paramStr += "Twide = "+"unknown"+" ; "
1695        paramStr += "Trock = "+num2str(gSamTrans)       
1696        pkStr += "Multiple coherent scattering fraction: "+num2str(g_MultScattFraction)
1697       
1698        //these strings are always the same
1699        samLabelStr ="LABEL: "+"Simulated USANS Data from: "+funcStr
1700        dateStr = date()+" at  "+time()// + " "+Secs2Time(DateTime,2)
1701       
1702       
1703        fileStr=nameStr+".cor"
1704        //actually open the file
1705        PathInfo savePathName
1706        if(V_flag==1)
1707                Open/P=savePathName refNum as fileStr
1708        else
1709                Open refNum as fileStr
1710        endif
1711       
1712        fprintf refnum,"%s"+termStr,samStr
1713        fprintf refnum,"%s"+termStr,dateStr
1714        fprintf refnum,"%s"+termStr,samLabelStr
1715        fprintf refnum,"%s"+termStr,empStr
1716        fprintf refnum,"%s"+termStr,paramStr
1717        fprintf refnum,"%s"+termStr,pkStr
1718        fprintf refnum,"%s"+termStr,empLevStr + " ; "+bkglevStr
1719       
1720        //
1721        wfprintf refnum, formatStr, tq,ti,te,dumWave,dumWave,dumWave
1722       
1723        Close refnum
1724       
1725        Killwaves/Z dumWave,ti,te,tq
1726       
1727       
1728       
1729        return(0)
1730End
1731
1732
1733//
1734// duplicates the BT5 file format
1735//
1736Function SaveFakeUSANS(nameStr,num,set)
1737        String nameStr
1738        Variable num,set
1739               
1740        String folder = "root:Packages:NIST:USANS:SIM:"
1741        String pathStr="root:Packages:NIST:USANS:Globals:U_Sim:"
1742        String termStr="\r\n"           //VAX uses only <CR> as terminator, but only CRLF seems to FTP correctly to VAX
1743       
1744//      WAVE DetCts = $(folder+"DetCts")                //these are only dummy values
1745        WAVE DetCts = root:Sim_USANS:Sim_USANS_i
1746        WAVE angle = $(folder+"Angle")
1747        WAVE SetNumber = $(folder+"SetNumber")
1748        WAVE countingTime = $(folder+"countingTime")
1749       
1750        NVAR ang1 = $(pathStr+"gAngLow"+num2str(set))
1751        NVAR ang2 = $(pathStr+"gAngHigh"+num2str(set))
1752        NVAR incr = $(pathStr+"gIncr"+num2str(set))
1753        SVAR ext = root:Packages:NIST:USANS:Globals:MainPanel:gUExt
1754
1755       
1756        Variable refNum,ii,wavePts,numPts,first,last,ctTime,monCt,transDet
1757        String str,fileStr,dateStr
1758       
1759        wavePts = numpnts(angle)
1760        for(ii=0;ii<wavePts;ii+=1)
1761                if(setNumber[ii] == set)
1762                        first = ii
1763                        break
1764                endif
1765        endfor
1766               
1767        for(ii=wavePts-1;ii>=0;ii-=1)
1768                if(setNumber[ii] == set)
1769                        last = ii
1770                        break
1771                endif
1772        endfor
1773       
1774//      Print "First, last = ",first,last
1775       
1776        // set up some of the strings needed
1777        fileStr=nameStr+num2str(num+set)+ext
1778        dateStr = date()
1779        ctTime = countingTime[first]
1780        numPts = last-first+1
1781       
1782        MonCt = ctTime*GetUSANSBeamIntensity()
1783       
1784        transDet = 999          //bogus and constant value
1785               
1786        //actually open the file
1787        PathInfo savePathName
1788        if(V_flag==1)
1789                Open/P=savePathName refNum as fileStr
1790        else
1791                Open refNum as fileStr
1792        endif
1793       
1794        sprintf str,"'%s' '%s' 'I'        %d    1  'TIME'   %d  'RAW'",fileStr,dateStr,ctTime,numpts
1795        fprintf refnum,"%s"+termStr,str
1796        fprintf refnum,"%s"+termStr,"  Filename         Date            Scan       Mon    Prf  Base   #pts  Type"
1797       
1798        fprintf refnum,"Simulated USANS data"+termStr
1799
1800        fprintf refnum,"%s"+termStr,"   0    0    0    0   0  0  0     0.0000     0.00000 0.00000   0.00000    2"
1801        fprintf refnum,"%s"+termStr," Collimation      Mosaic    Wavelength   T-Start   Incr.   H-field #Det    "
1802        fprintf refnum,"%s"+termStr,"  1       0.00000    0.00000    0.00000"
1803        fprintf refnum,"  2       %9.5f    %9.5f    %9.5f"+termStr,ang1,incr,ang2
1804        fprintf refnum,"%s"+termStr,"  3      10.00000    0.00000   10.00000"
1805        fprintf refnum,"%s"+termStr,"  4      10.00000    0.00000   10.00000"
1806        fprintf refnum,"%s"+termStr,"  5       0.00000    0.00000    0.00000"
1807        fprintf refnum,"%s"+termStr,"  6       0.00000    0.00000    0.00000"
1808        fprintf refnum,"%s"+termStr," Mot:    Start       Step      End"
1809        fprintf refnum,"%s"+termStr,"     A2       MIN      MONITOR      COUNTS      EXTRA  "
1810
1811        //loop over the waves, picking out the desired set
1812        //write 2 lines each time
1813        for(ii=first;ii<=last;ii+=1)
1814                sprintf str,"      %6.3f    %6.2f    %d     %d     %d",angle[ii],ctTime/60,MonCt,DetCts[ii]*ctTime,transDet
1815                fprintf refnum,"%s"+termStr,str
1816
1817                sprintf str,"%d,%d,0,%d,0,0,0,0",MonCt,DetCts[ii]*ctTime,transDet
1818                fprintf refnum,"%s"+termStr,str
1819        endfor
1820       
1821        Close refnum
1822
1823        return(0)
1824end
1825
1826
1827// print out the USANS configuration in some reasonable format
1828Function/S USANSConfigurationText()
1829
1830        String str="",temp
1831
1832        SetDataFolder root:Packages:NIST:USANS:Globals:U_Sim
1833       
1834        // results, setup values
1835        SVAR gTotTimeStr=gTotTimeStr
1836       
1837        Variable ii,num
1838        String pathStr="root:Packages:NIST:USANS:Globals:U_Sim:"
1839        num=7
1840       
1841        str += "USANS Instrument Configuration:\r\r"
1842        str += "Theta Min  Theta Max   Increment   # Points   Count Time\r"
1843       
1844       
1845        for(ii=1;ii<=num;ii+=1)
1846                NVAR ctTime = $(pathStr+"gCtTime"+num2str(ii))
1847                if(ctTime>0)
1848                        NVAR angLow = $(pathStr+"gAngLow"+num2str(ii))
1849                        NVAR angHigh = $(pathStr+"gAngHigh"+num2str(ii))
1850                        NVAR incr = $(pathStr+"gIncr"+num2str(ii))
1851                        NVAR numPts = $(pathStr+"gNumPts"+num2str(ii))
1852                       
1853                        sprintf temp,"%9.3f  %9.3f  %9.3f  %9d  %9d\r",angLow,angHigh,incr,numPts,ctTime
1854                        str += temp
1855                endif
1856        endfor
1857       
1858       
1859        sprintf temp,"\r\rTotal Counting Time (HR:MIN) = %s\r",gTotTimeStr
1860        str += temp
1861
1862       
1863   setDataFolder root:
1864   return str                   
1865End
1866
1867Function DisplayUCALCText()
1868
1869        if(WinType("USANS_Configuration")==0)
1870                NewNotebook/F=0/K=1/N=USANS_Configuration /W=(480,44,880,369)
1871        endif
1872        //replace the text
1873        Notebook USANS_Configuration selection={startOfFile, endOfFile}
1874        Notebook USANS_Configuration font="Monaco",fSize=10,text=USANSConfigurationText()
1875        return(0)
1876end
1877
1878
1879//
1880Function U_ConfigTextProc(ba) : ButtonControl
1881        STRUCT WMButtonAction &ba
1882
1883       
1884        switch( ba.eventCode )
1885                case 2: // mouse up
1886                        // click code here
1887                       
1888                                DisplayUCALCText()                                     
1889                        break
1890        endswitch
1891
1892        return 0
1893End
1894
1895// this will save a graphic of the whole panel that then needs to be opened and printed
1896// must be a PNG @ screen resolution
1897//
1898Function U_SavePanelProc(ba) : ButtonControl
1899        STRUCT WMButtonAction &ba
1900
1901       
1902        switch( ba.eventCode )
1903                case 2: // mouse up
1904                        // click code here
1905                       
1906                        SavePICT/P=home/E=-5/B=72/SNAP=1       
1907                       
1908                        // can I reload and print?
1909                        // how will the users know where this went and what to do with it?
1910                        //                     
1911                        break
1912        endswitch
1913
1914        return 0
1915End
Note: See TracBrowser for help on using the repository browser.