source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/SASCALC.ipf @ 891

Last change on this file since 891 was 891, checked in by srkline, 10 years ago

Changes so that the simulation in SASCALC will work well together - it makes sure that the analysis package is loaded if the simulation checkbox is selected.

File size: 74.2 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma version=1.0
3#pragma IgorVersion=6.1
4#pragma ModuleName=SASCALC
5
6// SASCALC.ipf
7//
8// 04 OCT 2006 SRK
9// 30 OCT 2006 SRK - corrected beamDiameter size in attenuator calculation (Bh vs Bm)
10// 11 DEC 2006 SRK - added 2.5 cm A1 option for NG3@7guide config
11// 09 MAR 2007 SRK - now appends text of each frozen configuration for printing
12//                                        - colorized frozen traces so that they aren't all red (unfrozen is black)
13// 19 MAR 2007 SRK - corrections added for projected BS diameter at anode plane
14// 11 APR 2007 SRK - default aperture offset of 5 cm added to match VAX implementation
15// nn AUG 2007 SRK - added defulat sample aperture size
16//                                               added option of lenses, approximated beamDiam=sourceDiam, BSdiam=1"
17//                                               Lens flux, trans is not corrected for lens/prism transmission
18//                                               Lenses can still be inserted in incorrect cases, and are not automatically taken out
19// 27 JAN 2009 SRK - Changed behavior of Lens checkbox. Now, it SETS parameters as needed for proper
20//                                               configuration. 17.2 can be typed in for lens/prism on NG3. Invalid conditions
21//                                               will automatically uncheck the box
22// 03 MAR 2012 SRK - Added new source apertures in 6th and 7th boxes of NG3
23//
24// 05 OCT 2012 SRK - (not visible) but added the skeleton bits for the 10m SANS instrument. Don't have the
25//                                                      details of the distances, etc, but I'll fill that in as needed
26//                                              -- to "un-hide" the 10m SANS, uncomment the CheckBox control in the panel (see the constant declared below)
27// 03 JAN 2013 SRK -- settled on "A" for new 10m SANS and "B" for the NG3 when it's relocated. Both instruments share
28//                                                      NGB, so this will have to do. files will be SA4 and SA5 respectively, if the VAX naming continues
29//
30//
31//
32// calculate what q-values you get based on the instruments settings
33// that you would typically input to SASCALC
34// calculation is for (80A radius spheres) * (beam stop shadow factor)
35// or a Debye function (Rg=50A) * (beam stop shadow factor)
36// - NOT true intensity, not counts, just a display
37//
38// TODO_10m:
39//
40// -- For the 10m SANS...
41//
42//              x- how many guides are there
43//              x- no lenses? focusing pinholes?
44//              x- aperture sizes?
45//              - flux values
46//              x- huber/sample chamber is inline, so get rid of this for 10m (gray out?)
47//              x- does this hook into anywhere in the main program? Attenuators, etc.
48//                      where knowledge of the 10m SANS is necessary??
49//              x- setting the instrument "number" is done here as "10". This will have a trickle-down effect for the places that use
50//              the global -- root:Packages:NIST:SAS:instrument (search for it everywhere...)
51//              x- repair the naming/numbering scheme to look for the string "NGx" and switch on that, not a number, which is unreliable and confusing.
52//
53// use str2hex(str) (my function) to convert "A" to 10
54//
55// Optional:
56// - freeze configurations with a user defined tag
57// - different model functions (+ change simple parameters)
58// - resolution smeared models
59// - "simulation" of data and error bars given a model and a total number of detector counts
60// - streamline code (globals needed in panel vs. wave needed for calculation)
61//
62// - there is a lot of re-calculation of things (a consequence of the fake-OO) that could be streamlined
63//
64// Done:
65// - include resolution effects (includes BS effect, smeared model)
66// - (data = 1) then multiply by a typical form factor
67// - masked two pixels around the edge, as default
68// - conversion from # guides to SSD from sascalc
69// - show multiple configurations at once
70// - interactive graphics on panel
71// - full capabilities of SASCALC
72// - correct beamstop diameter
73// - get the sample/huber position properly incorporated (ssd and sdd)
74// - get SDD value to update when switching NG7->NG3 and illegal value
75// - disallow 6 guides at NG3 (post a warning)
76//
77//
78
79Constant show10mSANS = 1
80
81
82
83Proc SASCALC()
84        DoWindow/F SASCALC
85        if(V_flag==0)
86                S_initialize_space()
87                initNG3()               //start life as NG3
88                Sascalc_Panel()
89                ReCalculateInten(1)             //will use defaults
90        Endif
91
92End
93
94Proc S_initialize_space()
95        NewDataFolder/O root:Packages
96        NewDataFolder/O root:Packages:NIST
97        NewDataFolder/O root:Packages:NIST:SAS
98       
99        Make/O/D/N=23 root:Packages:NIST:SAS:integersRead
100        Make/O/D/N=52 root:Packages:NIST:SAS:realsRead
101        Make/O/T/N=11 root:Packages:NIST:SAS:textRead
102        // data
103        Make/O/D/N=(128,128) root:Packages:NIST:SAS:data,root:Packages:NIST:SAS:linear_data
104        Make/O/D/N=2 root:Packages:NIST:SAS:aveint,root:Packages:NIST:SAS:qval,root:Packages:NIST:SAS:sigave
105        root:Packages:NIST:SAS:data = 1
106        root:Packages:NIST:SAS:linear_data = 1
107
108       
109        // other variables
110        // -(hard coded right now - look for NVAR declarations)
111        Variable/G root:Packages:NIST:gBinWidth=1               //uses global preference
112        Variable/G root:Packages:NIST:SAS:gisLogScale=0
113        String/G root:Packages:NIST:SAS:FileList = "SASCALC"
114       
115        // for the panel
116//      Variable/G root:Packages:NIST:SAS:gInst=3               //or 7 for NG7
117        String/G root:Packages:NIST:SAS:gInstStr="NG3"          //or "NG7" or "NGA"=10m or "NGB"=NG3 moved
118        Variable/G root:Packages:NIST:SAS:gNg=0
119        Variable/G root:Packages:NIST:SAS:gTable=2              //2=chamber, 1=table
120        Variable/G root:Packages:NIST:SAS:gDetDist=1000         //sample chamber to detector in cm
121        Variable/G root:Packages:NIST:SAS:gSSD=1632             //!!SSD in cm fo 0 guides (derived from Ng)
122        Variable/G root:Packages:NIST:SAS:gOffset=0
123        Variable/G root:Packages:NIST:SAS:gSamAp=1.27           //samAp diameter in cm
124        Variable/G root:Packages:NIST:SAS:gLambda=6
125        Variable/G root:Packages:NIST:SAS:gDeltaLambda=0.125            //default value
126        String/G root:Packages:NIST:SAS:gSourceApString = "1.43 cm;2.54 cm;3.81 cm;"
127        String/G root:Packages:NIST:SAS:gDeltaLambdaStr = "0.109;0.125;0.236;"          //ng3 defaults
128        String/G root:Packages:NIST:SAS:gApPopStr = "1/16\";1/8\";3/16\";1/4\";5/16\";3/8\";7/16\";1/2\";9/16\";5/8\";11/16\";3/4\";other;"
129        Variable/G root:Packages:NIST:SAS:gSamApOther = 10              //non-standard aperture diameter, in mm
130        Variable/G root:Packages:NIST:SAS:gUsingLenses = 0              //0=no lenses, 1=lenses(or prisms)
131        Variable/G root:Packages:NIST:SAS:gModelOffsetFactor = 1
132       
133        // for the MC simulation
134        Variable/G root:Packages:NIST:SAS:doSimulation  =0              // == 1 if 1D simulated data, 0 if other from the checkbox
135        Variable/G root:Packages:NIST:SAS:gRanDateTime=datetime
136        Variable/G root:Packages:NIST:SAS:gImon = 10000
137        Variable/G root:Packages:NIST:SAS:gThick = 0.1
138        Variable/G root:Packages:NIST:SAS:gSig_incoh = 0.1
139        String/G root:Packages:NIST:SAS:gFuncStr = ""
140        Variable/G root:Packages:NIST:SAS:gR2 = 2.54/2 
141        Variable/G root:Packages:NIST:SAS:gSamTrans=0.8                 //for 1D, default value
142        Variable/G root:Packages:NIST:SAS:gCntTime = 300
143        Variable/G root:Packages:NIST:SAS:gDoMonteCarlo = 0
144        Variable/G root:Packages:NIST:SAS:gUse_MC_XOP = 1                               //set to zero to use Igor code
145        Variable/G root:Packages:NIST:SAS:gBeamStopIn = 1                       //set to zero for beamstop out (transmission)
146        Variable/G root:Packages:NIST:SAS:gRawCounts = 0
147        Variable/G root:Packages:NIST:SAS:gSaveIndex = 100
148        String/G root:Packages:NIST:SAS:gSavePrefix = "SIMUL"
149        Variable/G root:Packages:NIST:SAS:gAutoSaveIndex = 100                  //a way to set the index for automated saves
150        String/G root:Packages:NIST:SAS:gAutoSaveLabel = ""                             //a way to set the "sample" label for automated saves
151        Make/O/D/N=10 root:Packages:NIST:SAS:results = 0
152        Make/O/T/N=10 root:Packages:NIST:SAS:results_desc = {"total X-section (1/cm)","SAS X-section (1/cm)","number that scatter","number that reach detector","avg # times scattered","fraction single coherent","fraction multiple coherent","fraction multiple scattered","fraction transmitted","detector counts w/beamstop"}
153
154        Variable/G root:Packages:NIST:SAS:g_1DTotCts = 0                        //summed counts (simulated)
155        Variable/G root:Packages:NIST:SAS:g_1DEstDetCR = 0              // estimated detector count rate
156        Variable/G root:Packages:NIST:SAS:g_1DFracScatt = 0             // fraction of beam captured on detector
157        Variable/G root:Packages:NIST:SAS:g_1DEstTrans = 0              // estimated transmission of sample (OUTPUT)
158        Variable/G root:Packages:NIST:SAS:g_1D_DoABS = 1
159        Variable/G root:Packages:NIST:SAS:g_1D_AddNoise = 1
160        Variable/G root:Packages:NIST:SAS:g_MultScattFraction=0
161        Variable/G root:Packages:NIST:SAS:g_detectorEff=0.75                    //average value for most wavelengths
162        Variable/G root:Packages:NIST:SAS:g_actSimTime = 0                              //for the save
163        Variable/G root:Packages:NIST:SAS:g_SimTimeWarn = 10                    //manually set to a very large value for scripted operation
164        Variable/G root:Packages:NIST:SAS:g_estimateOnly = 0                    // == 1 for just a time estimate, == 0 (default) to do the simulation
165        Variable/G root:Packages:NIST:SAS:g_estimatedMCTime             // estimated MC sim time
166       
167        PlotEC_Empirical(100,0.001,0.7)         // plot an empirical empty cell, only used for the simulation, not needed for SASCALC
168//      KillTopGraphAndTable(1)                                 // then kill the graph and table that the macro generates (but the function is in Wrapper, which is often not loaded
169        String topGraph= WinName(0,1)   //this is the topmost graph     
170        String topTable= WinName(0,2)   //this is the topmost table
171        KillWindow $topGraph
172        KillWindow $topTable
173       
174       
175        //tick labels for SDD slider
176        //userTicks={tvWave,tlblWave }
177        Make/O/D/N=5 root:Packages:NIST:SAS:tickSDDNG3,root:Packages:NIST:SAS:tickSDDNG7
178        Make/O/T/N=5 root:Packages:NIST:SAS:lblSDDNG3,root:Packages:NIST:SAS:lblSDDNG7
179        Make/O/D/N=5 root:Packages:NIST:SAS:tickSDDNGA
180        Make/O/T/N=5 root:Packages:NIST:SAS:lblSDDNGA
181        root:Packages:NIST:SAS:tickSDDNG3 = {133,400,700,1000,1317}
182        root:Packages:NIST:SAS:lblSDDNG3 = {"133","400","700","1000","1317"}
183        root:Packages:NIST:SAS:tickSDDNG7 = {100,450,800,1150,1530}
184        root:Packages:NIST:SAS:lblSDDNG7 = {"100","450","800","1150","1530"}
185        root:Packages:NIST:SAS:tickSDDNGA = {106,200,300,400,525}
186        root:Packages:NIST:SAS:lblSDDNGA = {"106","200","300","400","525"}
187//      root:Packages:NIST:SAS:tickSDDNGA = {106,200,300,400,500,525}
188//      root:Packages:NIST:SAS:lblSDDNGA = {"106","200","300","400","500","525"}
189               
190        //for the fake dependency
191        Variable/G root:Packages:NIST:SAS:gTouched=0
192        Variable/G root:Packages:NIST:SAS:gCalculate=0
193        //for plotting
194        Variable/G root:Packages:NIST:SAS:gFreezeCount=1                //start the count at 1 to keep Jeff happy
195        Variable/G root:Packages:NIST:SAS:gDoTraceOffset=0              // (1==Yes, offset 2^n), 0==turn off the offset
196       
197        // fill a fake data header
198        // fill w/default values
199        S_fillDefaultHeader(root:Packages:NIST:SAS:integersRead,root:Packages:NIST:SAS:realsRead,root:Packages:NIST:SAS:textRead)
200       
201End
202
203Function initNG3()
204
205        SetDataFolder root:Packages:NIST:SAS
206       
207        String/G gSelectedInstrument="checkNG3"
208        String/G gInstStr = "NG3"
209       
210        Variable/G s12 = 54.8
211        Variable/G d_det = 0.508
212        Variable/G a_pixel = 0.508
213        Variable/G del_r = 0.5
214        Variable/G det_width = 64.0
215        Variable/G lambda_t = 5.50
216        Variable/G l2r_lower = 132.3
217        Variable/G l2r_upper =  1317
218        Variable/G lambda_lower = 2.5
219        Variable/G lambda_upper = 20.0
220        Variable/G d_upper = 25.0
221        Variable/G bs_factor = 1.05
222        Variable/G t1 = 0.63
223        Variable/G t2 = 1.0
224        Variable/G t3 = 0.75
225        Variable/G l_gap = 100.0
226        Variable/G guide_width = 6.0
227        Variable/G idmax = 100.0
228//      //old values, from 3/2002
229//      Variable/G phi_0 = 2.95e13
230//      Variable/G b = 0.023
231//      Variable/G c = 0.023
232
233        //new values, from 11/2009 --- BeamFluxReport_2009.ifn
234        Variable/G phi_0 = 2.42e13
235        Variable/G b = 0.0
236        Variable/G c = -0.0243
237        Variable/G gGuide_loss = 0.924
238       
239//      //fwhm values (new variables) (+3, 0, -3, calibrated 2009)
240//      Variable/G fwhm_narrow = 0.109
241//      Variable/G fwhm_mid = 0.125
242//      Variable/G fwhm_wide = 0.236
243//     
244//      //source apertures (cm)
245//      Variable/G a1_0_0 = 1.43
246//      Variable/G a1_0_1 = 2.54
247//      Variable/G a1_0_2 = 3.81
248//      Variable/G a1_7_0 = 2.5 // after the polarizer         
249//      Variable/G a1_7_1 = 5.0
250//      Variable/G a1_7_1 = 0.95                // new Mar 2012 9.5 mm diameter source aperture
251//      Variable/G a1_def = 5.00
252       
253        //default configuration values
254//      ng = 0
255//      a1 = 3.81
256//      pos_table = 2
257//      l2r = 1310.0
258//      a2 = 1.27
259//      det_off = 0.0
260//      lambda = 6.0
261//      lambda_width = 0.15
262        Variable/G      l2diff = 5.0
263//     
264
265
266        UpdateSASHeader()
267       
268        SetDataFolder root:
269end
270
271Function initNG7()
272
273        SetDataFolder root:Packages:NIST:SAS
274
275        String/G gSelectedInstrument="checkNG7"
276        String/G gInstStr = "NG7"
277       
278        Variable/G s12 = 54.8
279        Variable/G d_det = 0.508
280        Variable/G a_pixel = 0.508
281        Variable/G del_r = 0.5
282        Variable/G det_width = 64.0
283        Variable/G lambda_t = 5.50
284        Variable/G l2r_lower = 100
285        Variable/G l2r_upper =  1531
286        Variable/G lambda_lower = 4.0
287        Variable/G lambda_upper = 20.0
288        Variable/G d_upper = 25.0
289        Variable/G bs_factor = 1.05
290        Variable/G t1 = 0.63
291        Variable/G t2 = 0.7
292        Variable/G t3 = 0.75
293        Variable/G l_gap = 188.0
294        Variable/G guide_width = 5.0
295        Variable/G idmax = 100.0
296       
297//      //old values, from 3/2002
298//      Variable/G phi_0 = 2.3e13
299//      Variable/G b = 0.028
300//      Variable/G c = 0.028
301
302        //new values, from 11/2009     
303        Variable/G phi_0 = 2.55e13
304        Variable/G b = 0.0395
305        Variable/G c = 0.0442
306        Variable/G gGuide_loss = 0.974
307       
308//      //fwhm values (new variables)
309//      Variable/G fwhm_narrow = 0.09
310//      Variable/G fwhm_mid = 0.115             //2009
311//      Variable/G fwhm_wide = 0.22
312//     
313//      //source apertures (cm)
314//      Variable/G a1_0_0 = 1.43
315//      Variable/G a1_0_1 = 2.22
316//      Variable/G a1_0_2 = 3.81
317//      Variable/G a1_7_0 = 5.0         // don't apply to NG7
318//      Variable/G a1_7_1 = 5.0
319//      Variable/G a1_def = 5.00
320       
321        //default configuration values
322//      ng = 0
323//      a1 = 2.22
324//      pos_table = 2
325//      l2r = 1530.0
326//      a2 = 1.27
327//      det_off = 0.0
328//      lambda = 6.0
329//      lambda_width = 0.11
330        Variable/G      l2diff = 5.0
331//     
332
333        UpdateSASHeader()
334        SetDataFolder root:
335end
336
337/// this is the (incomplete) definition of the 10m SANS instrument
338// on NG-B, which will be referred to as NGA here to keep the NG(number) notation
339// reserving NGB for the moved NG3 instrument
340// which may be simpler here to keep functions from breaking...
341//
342// Updated 24 JAN 2013 with current flux numbers from John
343//
344Function initNGA()
345
346        SetDataFolder root:Packages:NIST:SAS
347
348        String/G gSelectedInstrument="checkNGA"
349        String/G gInstStr = "NGA"
350
351        Variable/G s12 = 0                      //**            no difference between sample and huber position
352        Variable/G d_det = 0.508
353        Variable/G a_pixel = 0.508
354        Variable/G del_r = 0.5
355        Variable/G det_width = 64.0
356        Variable/G lambda_t = 5.50                              //good, Jan 2013
357        Variable/G l2r_lower = 106
358        Variable/G l2r_upper =  525
359        Variable/G lambda_lower = 3.0
360        Variable/G lambda_upper = 30.0
361        Variable/G d_upper = 25.0                       
362        Variable/G bs_factor = 1.05             
363        Variable/G t1 = 0.63                            //**good, Jan 2013
364        Variable/G t2 = 1.0                                     //**good, Jan 2013
365        Variable/G t3 = 0.75                            //**good, Jan 2013
366        Variable/G l_gap = 165                          //**good, Jan 2013
367        Variable/G guide_width = 5.0            //**good, Jan 2013
368        Variable/G idmax = 100.0                        //** max count rate per pixel
369//      current values -- Jan 2013
370        Variable/G phi_0 = 2.5e13
371        Variable/G b = 0.03
372        Variable/G c = 0.03
373
374        Variable/G gGuide_loss = 0.95
375       
376//      //fwhm values (new variables) (+3, 0, -3, calibrated 2009)
377//      Variable/G fwhm_narrow = 0.10
378//      Variable/G fwhm_mid = 0.20
379//      Variable/G fwhm_wide = 0.30
380//     
381//      //source apertures (cm)
382//      Variable/G a1_0_0 = 1.43
383//      Variable/G a1_0_1 = 2.54
384//      Variable/G a1_0_2 = 3.81
385//      Variable/G a1_7_0 = 2.5 // after the polarizer         
386//      Variable/G a1_7_1 = 5.0
387//      Variable/G a1_7_1 = 0.95                // new Mar 2012 9.5 mm diameter source aperture
388//      Variable/G a1_def = 5.00
389//     
390        //default configuration values
391//      ng = 0
392//      a1 = 3.81
393//      pos_table = 2
394//      l2r = 1310.0
395//      a2 = 1.27
396//      det_off = 0.0
397//      lambda = 6.0
398//      lambda_width = 0.15
399        Variable/G      l2diff = 5.0
400//     
401        UpdateSASHeader()
402        SetDataFolder root:
403end
404
405// as the instrument changes, fill in appropriate values into the header waves of the SAS folder
406// so that any simulation data (2D VAX) that is written out will be as correct as possible
407// -- this should be taken care of in SimulationVAXHeader(), so this is really redundant
408Function UpdateSASHeader()
409
410        Wave/T tw=root:Packages:NIST:SAS:textRead
411        // fill in the instrument
412        SVAR gInstStr = root:Packages:NIST:SAS:gInstStr
413        tw[3] = "["+gInstStr+"SANS99]"
414
415End
416
417Function S_fillDefaultHeader(iW,rW,tW)
418        Wave iW,rW
419        Wave/T tW
420
421        // text wave
422        // fill in the instrument
423        SVAR gInstStr = root:Packages:NIST:SAS:gInstStr
424        tw[3] = "["+gInstStr+"SANS99]"
425       
426        // integer wave
427        // don't need anything
428       
429        // real wave
430        rw = 0
431        rw[16] = 65             // beamcenter X (pixels)
432        rw[17] = 64             // beamcenter Y
433       
434        rw[10]  = 5.08                  //detector resolution (5mm) and calibration constants (linearity)
435        rw[11] = 10000
436        rw[12] = 0
437        rw[13] = 5.08
438        rw[14] = 10000
439        rw[15] = 0
440       
441        rw[20] = 65             // det size in cm
442        rw[18] = 6              // SDD in meters (=L2)
443       
444        rw[26] = 6              //lambda in Angstroms
445        rw[4] = 1               //transmission
446       
447        rw[21] = 76.2                   //BS diameter in mm
448        rw[23] = 50                     //A1 diameter in mm
449        rw[24] = 12.7                   //A2 diameter in mm
450        rw[25] = 7.02                   //L1 distance in meters (derived from number of guides)
451        rw[27] = 0.11                   //DL/L wavelength spread
452       
453        return(0)
454End
455
456//
457Window SASCALC_Panel()
458
459        PauseUpdate; Silent 1           // building window...
460
461// if I make the graph a subwindow in a panel, it breaks the "append 1d" from the wrapper       
462//      NewPanel/W=(5,44,463,570)/K=1 as "SASCALC"
463//      DoWindow/C SASCALC
464//      ModifyPanel cbRGB=(49151,53155,65535)
465//
466//     
467//      String fldrSav0= GetDataFolder(1)
468//      SetDataFolder root:Packages:NIST:SAS:
469//     
470//      Display/HOST=#/W=(5,200,463,570) aveint vs qval
471//      ModifyGraph mode=3
472//      ModifyGraph marker=19
473//      ModifyGraph rgb=(0,0,0)
474//      Modifygraph log=1
475//      Modifygraph grid=1
476//      Modifygraph mirror=2
477//      ModifyGraph msize(aveint)=2
478//      ErrorBars/T=0 aveint Y,wave=(sigave,sigave)
479//      Label bottom, "Q (1/A)"
480//      Label left, "Relative Intensity"
481//      legend
482//     
483//      RenameWindow #,G_aveint
484//      SetActiveSubwindow ##
485
486//// end panel commands
487
488/// draw as a graph
489        String fldrSav0= GetDataFolder(1)
490        SetDataFolder root:Packages:NIST:SAS:
491       
492        Display/W=(5,44,463,570)/K=1  aveint vs qval as "SASCALC"
493        DoWindow/C SASCALC
494        ModifyGraph cbRGB=(49151,53155,65535)
495        ModifyGraph mode=3
496        ModifyGraph marker=19
497        ModifyGraph rgb=(0,0,0)
498        Modifygraph log=1
499        Modifygraph grid=1
500        Modifygraph mirror=2
501        ModifyGraph msize(aveint)=2
502        ErrorBars/T=0 aveint Y,wave=(sigave,sigave)
503        Label bottom, "Q (1/A)"
504        Label left, "Relative Intensity"
505        legend
506       
507        ControlBar/T 200
508
509        SetDataFolder fldrSav0
510       
511        Slider SC_Slider,pos={11,46},size={150,45},proc=GuideSliderProc,live=0
512        Slider SC_Slider,limits={0,8,1},variable= root:Packages:NIST:SAS:gNg,vert= 0//,thumbColor= (1,16019,65535)
513        Slider SC_Slider_1,pos={234,45},size={150,45},proc=DetDistSliderProc,live=0
514        Slider SC_Slider_1,tkLblRot= 90,userTicks={root:Packages:NIST:SAS:tickSDDNG3,root:Packages:NIST:SAS:lblSDDNG3 }
515        Slider SC_Slider_1,limits={133,1317,1},variable= root:Packages:NIST:SAS:gDetDist,vert= 0//,thumbColor= (1,16019,65535)
516        Slider SC_Slider_2,pos={394,21},size={47,65},proc=OffsetSliderProc,live=0,ticks=4
517        Slider SC_Slider_2,limits={0,25,1},variable= root:Packages:NIST:SAS:gOffset//,thumbColor= (1,16019,65535)
518        CheckBox checkNG3,pos={20,19},size={36,14},proc=SelectInstrumentCheckProc,title="NG3"
519        CheckBox checkNG3,value=1,mode=1
520        CheckBox checkNG7,pos={66,19},size={36,14},proc=SelectInstrumentCheckProc,title="NG7"
521        CheckBox checkNG7,value=0,mode=1
522
523        CheckBox checkChamber,pos={172,48},size={57,14},proc=TableCheckProc,title="Chamber"
524        CheckBox checkChamber,value=1,mode=1
525        CheckBox checkHuber,pos={172,27},size={44,14},proc=TableCheckProc,title="Huber"
526        CheckBox checkHuber,value=0,mode=1
527//      -- hide/unhide the 10m SANS
528        if(show10mSANS)
529                CheckBox checkNGA,pos={110,19},size={40,14},proc=SelectInstrumentCheckProc,title="NGA"
530                CheckBox checkNGA,value=0,mode=1
531        endif
532//             
533        PopupMenu popup0,pos={6,94},size={76,20},proc=SourceAperturePopMenuProc
534        PopupMenu popup0,mode=1,popvalue="3.81 cm",value= root:Packages:NIST:SAS:gSourceApString
535        PopupMenu popup0_1,pos={172,72},size={49,20},proc=SampleAperturePopMenuProc
536        PopupMenu popup0_1,mode=8,popvalue="1/2\"",value= root:Packages:NIST:SAS:gApPopStr
537        SetVariable setvar0,pos={301,107},size={130,15},title="Det Dist (cm)",proc=SDDSetVarProc
538        SetVariable setvar0,limits={133,1317,1},value=root:Packages:NIST:SAS:gDetDist
539        SetVariable setvar0_1,pos={321,129},size={110,15},title="Offset (cm)",proc=OffsetSetVarProc
540        SetVariable setvar0_1,limits={0,25,1},value= root:Packages:NIST:SAS:gOffset
541        SetVariable setvar0_2,pos={6,130},size={90,15},title="Lambda",proc=LambdaSetVarProc
542        SetVariable setvar0_2,limits={4,20,0.1},value= root:Packages:NIST:SAS:gLambda
543        PopupMenu popup0_2,pos={108,127},size={55,20},proc=DeltaLambdaPopMenuProc
544        PopupMenu popup0_2,mode=2,value= root:Packages:NIST:SAS:gDeltaLambdaStr
545        Button FreezeButton title="Freeze",size={60,20},pos={180,166}
546        Button FreezeButton proc=FreezeButtonProc
547        Button ClearButton title="Clear",size={60,20},pos={250,166}
548        Button ClearButton proc=S_ClearButtonProc
549//      GroupBox group0,pos={6,1},size={108,36},title="Instrument"
550        GroupBox group0,pos={6,1},size={160,36},title="Instrument"
551        SetDataFolder fldrSav0
552       
553        SetVariable setvar0_3,pos={140,94},size={110,15},title="Diam (mm)",disable=1
554        SetVariable setvar0_3,limits={0,100,0.1},value= root:Packages:NIST:SAS:gSamApOther,proc=SampleApOtherSetVarProc
555       
556        CheckBox checkLens,pos={6,155},size={44,14},proc=LensCheckProc,title="Lenses?"
557        CheckBox checkLens,value=root:Packages:NIST:SAS:gUsingLenses
558       
559        CheckBox checkSim,pos={6,175},size={44,14},proc=SimCheckProc,title="Simulation?"
560        CheckBox checkSim,value=0
561       
562        CheckBox check0_1,pos={80,155},size={90,14},title="Offset Traces?",variable= root:Packages:NIST:SAS:gDoTraceOffset
563
564        // help, done buttons
565        Button SC_helpButton,pos={340,166},size={25,20},proc=showSASCALCHelp,title="?"
566        Button SC_helpButton,help={"Show help file for simulation of SANS Data"}
567        Button SC_DoneButton,pos={380,166},size={50,20},proc=SASCALCDoneButton,title="Done"
568        Button SC_DoneButton,help={"This button will close the panel"}
569
570        // set up a fake dependency to trigger recalculation
571        //root:Packages:NIST:SAS:gCalculate := ReCalculateInten(root:Packages:NIST:SAS:gTouched)
572EndMacro
573
574//clean up
575Proc SASCALCDoneButton(ctrlName): ButtonControl
576        String ctrlName
577        DoWindow/K SASCALC
578        DoWindow/K Trial_Configuration
579        DoWindow/K Saved_Configurations
580        DoWindow/K MC_SASCALC
581        DoWindow/K Sim_1D_Panel
582end
583
584//
585Proc showSASCALCHelp(ctrlName): ButtonControl
586        String ctrlName
587        DisplayHelpTopic/K=1/Z "SASCALC"
588        if(V_flag !=0)
589                DoAlert 0,"The SANS Simulation Help file could not be found"
590        endif
591end
592
593// based on the instrument selected...
594// set the SDD range
595// set the source aperture popup (based on NGx and number of guides)
596// set the wavelength spread popup
597//
598Function UpdateControls()
599        //poll the controls on the panel, and change needed values
600        Variable isNG3,Ng,mode
601
602        ControlInfo/W=SASCALC SC_slider
603        Ng=V_value
604        NVAR gNg = root:Packages:NIST:SAS:gNg
605        gNg = Ng                //update the global from the panel
606       
607        SVAR A1str= root:Packages:NIST:SAS:gSourceApString// = "1.43 cm;2.54 cm;3.81 cm;"
608        SVAR dlStr = root:Packages:NIST:SAS:gDeltaLambdaStr
609        SVAR selInstr = root:Packages:NIST:SAS:gInstStr
610        NVAR lambda = root:Packages:NIST:SAS:gLambda
611
612        strswitch(selInstr)     // string switch
613                case "NG3":                     //
614                        switch(ng)     
615                                case 0:
616                                        ControlInfo/W=SASCALC popup0
617                                        mode=V_value
618                                        A1str="1.43 cm;2.54 cm;3.81 cm;"
619                                        break
620        //                      case 6:                                                                 //Mar 2012 - 6 guides is OK now, 5 cm aperture present in 7th box
621        //                              A1str = "! 6 Guides invalid;"
622        //                              mode=1
623        //                              break
624                                case 7:                                                 // switched order in 2009 to keep 5 cm as default, 2.5 cm for polarizer
625                                        A1Str = "5.00 cm;2.50 cm;0.95 cm;"              // Mar 2012 - 3rd aperture added
626                                        mode = 1
627                                        break
628                                default:
629                                        A1str = "5 cm;"
630                                        mode=1
631                        endswitch
632                        //wavelength spread
633                        dlStr = "0.109;0.125;0.236;"            //updated calibration 2009
634                        Slider SC_Slider win=SASCALC,limits={0,8,1},ticks=1
635                        //detector limits
636                        Slider SC_Slider_1,win=SASCALC,limits={133,1317,1},userTicks={root:Packages:NIST:SAS:tickSDDNG3,root:Packages:NIST:SAS:lblSDDNG3 }
637                        SetVariable setvar0,win=SASCALC,limits={133,1317,1}
638                        NVAR detDist=root:Packages:NIST:SAS:gDetDist
639                        if(detDist < 133 )
640                                detDist = 133
641                        elseif (detDist > 1317 )
642                                detDist = 1317
643                        endif
644
645                        // wavelength limits
646                        SetVariable setvar0_2 win=SASCALC,limits={4,20,0.1}
647                        if(lambda < 4)
648                                lambda = 4
649                        endif
650                        if(lambda > 20)
651                                lambda = 20
652                        endif
653
654                        // un-disable the sample table position
655                        CheckBox checkChamber,win=SASCALC,disable=0
656                        CheckBox checkHuber,win=SASCALC,disable=0
657                                               
658                        break
659                                                               
660                case "NG7":                     //
661                        switch(ng)     
662                                case 0:
663                                        ControlInfo/W=SASCALC popup0
664                                        mode=V_value
665                                        A1str="1.43 cm;2.22 cm;3.81 cm;"
666                                        break
667                                default:
668                                        A1str = "5.08 cm;"
669                                        mode=1
670                        endswitch
671                       
672                        dlStr = "0.09;0.115;0.22;"
673                        Slider SC_Slider win=SASCALC,limits={0,8,1},ticks=1
674                        Slider SC_Slider_1,win=SASCALC,limits={100,1531,1},userTicks={root:Packages:NIST:SAS:tickSDDNG7,root:Packages:NIST:SAS:lblSDDNG7 }
675                        SetVariable setvar0,win=SASCALC,limits={100,1531,1}
676                       
677                        if(detDist < 100 )
678                                detDist = 100
679                        elseif (detDist > 1531 )
680                                detDist = 1531
681                        endif
682                        // wavelength limits
683                        SetVariable setvar0_2 win=SASCALC,limits={4,20,0.1}
684                        if(lambda < 4)
685                                lambda = 4
686                        endif
687                        if(lambda > 20)
688                                lambda = 20
689                        endif
690                                               
691                        // un-disable the sample table position
692                        CheckBox checkChamber,win=SASCALC,disable=0
693                        CheckBox checkHuber,win=SASCALC,disable=0
694                       
695                        break
696                       
697                case "NGA":             // 10m SANS
698                        if(ng>2)
699                                ng=2
700                                gNg = Ng                //update the global
701                        endif
702                        switch(ng)     
703                                case 0:
704                                        ControlInfo/W=SASCALC popup0
705                                        mode=V_value
706                                        A1str="1.3 cm;2.5 cm;3.8 cm;"
707                                        break
708                                default:
709                                        A1str = "5.0 cm;"
710                                        mode=1
711                        endswitch
712                       
713                        dlStr = "0.10;0.132;0.154;0.25;"
714                        Slider SC_Slider win=SASCALC,limits={0,2,1},ticks=2                     //number of guides different on 10m SANS, 3 ticks
715                        Slider SC_Slider_1,win=SASCALC,limits={106,525,1},userTicks={root:Packages:NIST:SAS:tickSDDNGA,root:Packages:NIST:SAS:lblSDDNGA }
716                        SetVariable setvar0,win=SASCALC,limits={106,525,1}
717
718                        NVAR detDist=root:Packages:NIST:SAS:gDetDist
719                        if(detDist < 106 )
720                                detDist = 106
721                        elseif (detDist > 525 )
722                                detDist = 525
723                        endif
724                       
725                        // wavelength limits
726                        SetVariable setvar0_2 win=SASCALC,limits={3,30,0.1}
727                        if(lambda < 3)
728                                lambda = 3
729                        endif
730                        if(lambda > 30)
731                                lambda = 30
732                        endif
733                                       
734                        // disable the sample table position, as long as the sample positions are identical for the 10 m
735                        CheckBox checkChamber,win=SASCALC,disable=2
736                        CheckBox checkHuber,win=SASCALC,disable=2
737
738                        break
739                       
740                default:                                                        // optional default expression executed
741                                DoAlert 0,"No matching instrument!  UpdateControls"                                     // when no case matches
742        endswitch
743       
744        // force update of the sliders
745        Slider SC_Slider,win=SASCALC,variable=root:Packages:NIST:SAS:gNg                                //guides
746        Slider SC_Slider_1,win=SASCALC,variable=root:Packages:NIST:SAS:gDetDist         //detector distance
747        Slider SC_Slider_2,win=SASCALC,variable=root:Packages:NIST:SAS:gOffset          //detector offset
748
749
750        ControlUpdate popup0
751        PopupMenu popup0,win=SASCALC,mode=mode          //source Ap
752        ControlInfo/W=SASCALC popup0
753        SourceAperturePopMenuProc("",0,S_Value)                 //send popNum==0 so recalculation won't be done
754       
755        ControlUpdate popup0_2          // delta lambda pop
756        ControlInfo/W=SASCALC popup0_2
757        DeltaLambdaPopMenuProc("",0,S_Value)            //sets the global and the wave
758End
759
760//changing the number of guides changes the SSD
761// the source aperture popup may need to be updated
762//
763Function GuideSliderProc(ctrlName,sliderValue,event)
764        String ctrlName
765        Variable sliderValue
766        Variable event  // bit field: bit 0: value set, 1: mouse down, 2: mouse up, 3: mouse moved
767       
768        Variable recalc=0
769       
770        if(event %& 0x1)        // bit 0, value set
771                if(cmpstr(ctrlName,"") != 0)            //here by direct action, so do LensCheck and recalculate
772                        recalc=1
773                        LensCheckProc("",2)             //make sure lenses are deselected
774                endif
775                sourceToSampleDist()            //updates the SSD global and wave
776                //change the sourceAp popup, SDD range, etc
777                UpdateControls()                                // the Ng global is actually set inside this function
778                ReCalculateInten(recalc)
779        endif
780        return 0
781End
782
783// changing the detector position changes the SDD
784//
785Function DetDistSliderProc(ctrlName,sliderValue,event)
786        String ctrlName
787        Variable sliderValue
788        Variable event  // bit field: bit 0: value set, 1: mouse down, 2: mouse up, 3: mouse moved
789
790        Variable recalc=0
791        if(event %& 0x1)        // bit 0, value set
792                if(cmpstr(ctrlName,"") != 0)
793                        recalc=1
794                        LensCheckProc("",2)             //make sure lenses are only selected for valid configurations
795                endif
796                sampleToDetectorDist()          //changes the SDD and wave (DetDist is the global)
797                ReCalculateInten(recalc)
798        endif
799
800        return 0
801End
802
803// change the offset
804// - changes the beamcenter (x,y) position too
805Function OffsetSliderProc(ctrlName,sliderValue,event)
806        String ctrlName
807        Variable sliderValue
808        Variable event  // bit field: bit 0: value set, 1: mouse down, 2: mouse up, 3: mouse moved
809
810        if(event %& 0x1)        // bit 0, value set
811                detectorOffset()
812                ReCalculateInten(1)
813        endif
814
815        return 0
816End
817
818// changing the instrument -
819// re-initialize the global values
820// update the controls to appropriate limits/choices
821//
822// updated for 10m SANS
823Function SelectInstrumentCheckProc(ctrlName,checked) : CheckBoxControl
824        String ctrlName
825        Variable checked
826
827
828        strswitch(ctrlName)     // string switch
829                case "checkNG3":                        //
830                        checkBox checkNG3,win=SASCALC, value=1
831                        checkBox checkNG7,win=SASCALC, value=0
832                        checkBox checkNGA,win=SASCALC, value=0
833                        initNG3()
834                        break                                           
835                case "checkNG7":                        //
836                        checkBox checkNG3,win=SASCALC, value=0
837                        checkBox checkNG7,win=SASCALC, value=1
838                        checkBox checkNGA,win=SASCALC, value=0
839                        initNG7()
840                        break
841                case "checkNGA":                // 10m SANS
842                        checkBox checkNG3,win=SASCALC, value=0
843                        checkBox checkNG7,win=SASCALC, value=0
844                        checkBox checkNGA,win=SASCALC, value=1
845                        initNGA()
846                        break
847                default:                                                        // optional default expression executed
848                                                                        // when no case matches
849        endswitch
850
851
852        LensCheckProc("",2)             //check if lenses are still valid (they won't be)
853        UpdateControls()
854        ReCalculateInten(1)
855End
856
857//changing the table position
858// changes the SSD and SDD
859Function TableCheckProc(ctrlName,checked) : CheckBoxControl
860        String ctrlName
861        Variable checked
862
863        NVAR table=root:Packages:NIST:SAS:gTable
864        if(cmpstr(ctrlName,"checkHuber")==0)
865                checkBox checkHuber,win=SASCALC, value=1
866                checkBox checkChamber,win=SASCALC, value=0
867                table=1         //in Huber position
868        else
869                checkBox checkHuber,win=SASCALC, value=0
870                checkBox checkChamber,win=SASCALC, value=1
871                table = 2               //in Sample chamber
872        endif
873        sampleToDetectorDist()
874        sourceToSampleDist()            //update
875        ReCalculateInten(1)
876End
877
878
879//lenses (or prisms) in/out changes resolution
880//
881// passing in a checked == 2 will do a check of the configuration without
882// affecting the box state
883//
884// if an invalid configuration is detected, the box is unchecked
885// and the lens flag is set to zero (=out).
886//
887// When necessary controls are "popped", a ctrlName="" is passed to signal the control
888// to NOT recalculate the intensity, and to NOT (recursively) call LensCheckProc again
889//
890// currently, the 17.2 A for lens/prism @ ng3 must be typed in
891//
892Function LensCheckProc(ctrlName,checked) : CheckBoxControl
893        String ctrlName
894        Variable checked
895
896        // don't let the box get checked if the conditions are wrong
897        // lambda != 8.09,8.4,17.2
898        // Ng != 0
899        NVAR lens = root:Packages:NIST:SAS:gUsingLenses
900        NVAR ng = root:Packages:NIST:SAS:gNg
901        NVAR lam = root:Packages:NIST:SAS:gLambda
902        NVAR dist = root:Packages:NIST:SAS:gDetDist
903//      NVAR instrument = root:Packages:NIST:SAS:instrument
904        Wave rw=root:Packages:NIST:SAS:realsRead
905       
906        SVAR selInstr = root:Packages:NIST:SAS:gInstStr
907               
908        // directly uncheck the box, just set the flag and get out
909        if(checked == 0)
910                lens = 0
911                CheckBox checkLens,win=SASCALC,value=0
912                rw[28]=0                //flag for lenses out
913                ReCalculateInten(1)
914                return(0)
915        endif
916       
917        // check the box, enforce the proper conditions
918        if(checked == 1)
919                lens = checked 
920                strswitch(selInstr)     // string switch
921                        case "NG3":
922                                dist = 1317
923                                DetDistSliderProc("",1317,1)
924                               
925                                lam = 8.4
926                                LambdaSetVarProc("",8.4,"8.4","")
927       
928                                ng=0
929                                Slider SC_Slider,win=SASCALC,variable=root:Packages:NIST:SAS:gNg
930                                GuideSliderProc("",0,1)
931
932                                PopupMenu popup0,win=SASCALC,mode=1,popvalue="1.43 cm"          //first item in source aperture menu
933                               
934                                PopupMenu popup0_2,win=SASCALC,mode=2           //deltaLambda
935                                ControlInfo/W=SASCALC popup0_2
936                                DeltaLambdaPopMenuProc("",0,S_value)                    //zero as 2nd param skips recalculation
937                                rw[28]=1                //flag for lenses in (not the true number, but OK)
938               
939                                break
940                        case "NG7":
941                                dist = 1531
942                                DetDistSliderProc("",1531,1)
943                               
944                                lam = 8.09
945                                LambdaSetVarProc("",8.09,"8.09","")
946                               
947                                ng=0
948                                Slider SC_Slider,win=SASCALC,variable=root:Packages:NIST:SAS:gNg
949                                GuideSliderProc("",0,1)
950                               
951                                PopupMenu popup0,win=SASCALC,mode=1,popvalue="1.43 cm"          //first item
952                               
953                                PopupMenu popup0_2,win=SASCALC,mode=2           //deltaLambda
954                                ControlInfo/W=SASCALC popup0_2
955                                DeltaLambdaPopMenuProc("",0,S_value)                    //zero as 2nd param skips recalculation
956                                rw[28]=1                //flag for lenses in (not the true number, but OK)
957                               
958                                break
959                        case "NGA":
960                                // 10m SANS - force no lenses for now
961                                // TODO:  -- put in CORRECT VALUES -- THESE ARE FICTIONAL
962                                lens = 0                //no lenses
963                                CheckBox checkLens,win=SASCALC,value=0
964                                rw[28]=0                //flag for lenses out
965
966                                break
967                        default:
968                                DoAlert 0,"No matching instrument! LCP"
969                endswitch
970        endif
971       
972        // this is my internal check to see if conditions are still valid
973        // I get to this point by calling this myself and passing the fictional value of (2) for the
974        // checkbox value...
975        //
976       
977        // I'll uncheck as needed
978        Variable lensNotAllowed=0
979        if(checked == 2)
980
981                // source aperture must be 1.43 cm diameter
982                // number of guides must be zero
983                Variable a1 = sourceApertureDiam()
984                if(a1 != 1.43  || Ng !=0)
985                        lensNotAllowed=1
986                endif
987       
988                // instrument specific distance requirements
989                if(cmpstr(selInstr,"NG3") == 0 && dist != 1317)
990                        lensNotAllowed=1
991                endif
992       
993                if(cmpstr(selInstr,"NG7") == 0 && dist != 1531)
994                        lensNotAllowed=1
995                endif
996               
997                // instrument specific wavelength requirements
998                if(cmpstr(selInstr,"NG3") == 0 && !(lam == 8.4 || lam == 17.2) )
999                        lensNotAllowed=1
1000                endif
1001               
1002                if(cmpstr(selInstr,"NG7") == 0 && lam != 8.09 )
1003                        lensNotAllowed=1
1004                endif
1005
1006        // right now, if 10m instrument, no lenses allowed             
1007                if(cmpstr(selInstr,"NGA") == 0)
1008                        lensNotAllowed=1
1009                endif
1010               
1011        endif
1012       
1013        if(lensNotAllowed)
1014                lens = 0
1015                CheckBox checkLens,win=SASCALC,value=0
1016                rw[28]=0                //flag for lenses out
1017                return(0)
1018        endif
1019
1020        ReCalculateInten(1)
1021
1022
1023        return(1)               //return value not used
1024End
1025
1026//
1027//simulation controls as a control bar that toggles on/off to the right
1028//
1029// depending on the state of the 2D flag, open the 1d or 2d control panel
1030// AND - if doing simulation, check first that the Analysis package is loaded - otheriwse lots of errors
1031Function SimCheckProc(CB_Struct) : CheckBoxControl
1032        STRUCT WMCheckboxAction &CB_Struct
1033
1034
1035        if(CB_Struct.checked)
1036                if(exists("SANS_ANA_VERSION")!=2)               //if the analysis package is not loaded
1037                        // be sure that the SANS Analysis is loaded and compiles
1038                        NCNR_AnalysisLoader("Load NCNR Analysis Macros")
1039                endif
1040               
1041                NVAR do2D = root:Packages:NIST:SAS:gDoMonteCarlo
1042               
1043                if((CB_Struct.eventMod & 2^1) != 0 )            //if the shift key is down (bit 1)- go to 2D mode
1044                        do2D = 1
1045                endif
1046               
1047                if(do2D)
1048                        DoWindow/F MC_SASCALC
1049                        if(V_flag==0)
1050                                Execute "MC_SASCALC()"          //sets the variable
1051                                AutoPositionWindow/M=1/R=SASCALC MC_SASCALC
1052                        endif
1053                else
1054                        //draw 1D panel
1055                        DoWindow/F Sim_1D_Panel
1056                        if(V_flag==0)
1057                                Execute "Sim_1D_Panel()"
1058                                AutoPositionWindow/M=1/R=SASCALC Sim_1D_Panel
1059                        endif
1060                endif
1061                NVAR doSim = root:Packages:NIST:SAS:doSimulation
1062                doSim=1
1063        else
1064                //get rid of the controls (just try both)
1065                DoWindow MC_SASCALC
1066                if(V_flag!=0)
1067                        NVAR do2D = root:Packages:NIST:SAS:gDoMonteCarlo
1068                        do2D = 0
1069                        KillWindow MC_SASCALC
1070                endif
1071                DoWindow Sim_1D_Panel
1072                if(V_flag!=0)
1073                        KillWindow Sim_1D_Panel
1074                endif
1075               
1076                NVAR doSim = root:Packages:NIST:SAS:doSimulation
1077                doSim=0
1078        endif
1079
1080        return(0)
1081End
1082
1083//
1084// change the source aperture
1085//
1086Function SourceAperturePopMenuProc(ctrlName,popNum,popStr) : PopupMenuControl
1087        String ctrlName
1088        Variable popNum
1089        String popStr
1090
1091        Variable a1 = sourceApertureDiam()              //sets the new value in the wave
1092       
1093        // if LensCheckProc is calling this, don't call it again.
1094        // if this is called by a direct pop, then check
1095        if(cmpstr(ctrlName,"") != 0)
1096                LensCheckProc("",2)             //make sure lenses are only selected for valid configurations
1097        endif
1098       
1099        ReCalculateInten(popnum)                //skip the recalculation if I pass in a zero
1100End
1101
1102// set sample aperture
1103//
1104Function SampleAperturePopMenuProc(ctrlName,popNum,popStr) : PopupMenuControl
1105        String ctrlName
1106        Variable popNum
1107        String popStr
1108
1109        sampleApertureDiam()
1110        ReCalculateInten(1)
1111End
1112
1113// set sample to detector distance
1114//
1115Function SDDSetVarProc(ctrlName,varNum,varStr,varName) : SetVariableControl
1116        String ctrlName
1117        Variable varNum
1118        String varStr
1119        String varName
1120       
1121        sampleToDetectorDist()
1122        if(cmpstr(ctrlName,"") != 0)
1123                LensCheckProc("",2)             //make sure lenses are only selected for valid configurations
1124        endif
1125        ReCalculateInten(1)
1126End
1127
1128// set offset
1129// -- also changes the beamcenter (x,y)
1130//
1131Function OffsetSetVarProc(ctrlName,varNum,varStr,varName) : SetVariableControl
1132        String ctrlName
1133        Variable varNum
1134        String varStr
1135        String varName
1136
1137        detectorOffset()                //sets the offset in the wave and also the new (x,y) beamcenter
1138       
1139        ReCalculateInten(1)
1140        return(0)
1141End
1142
1143// change the wavelength
1144Function LambdaSetVarProc(ctrlName,varNum,varStr,varName) : SetVariableControl
1145        String ctrlName
1146        Variable varNum
1147        String varStr
1148        String varName
1149
1150        WAVE rw=root:Packages:NIST:SAS:realsRead
1151        Variable recalc=0
1152       
1153        rw[26] = str2num(varStr)
1154        if(cmpstr(ctrlName,"") != 0)
1155                recalc=1
1156                LensCheckProc("",2)             //make sure lenses are only selected for valid configurations
1157        endif
1158        ReCalculateInten(recalc)
1159        return(0)
1160End
1161
1162
1163// change the wavelength spread
1164Function DeltaLambdaPopMenuProc(ctrlName,popNum,popStr) : PopupMenuControl
1165        String ctrlName
1166        Variable popNum
1167        String popStr
1168
1169        NVAR dl=root:Packages:NIST:SAS:gDeltaLambda
1170        dl=str2num(popStr)
1171        WAVE rw=root:Packages:NIST:SAS:realsRead
1172        rw[27] = dl
1173        ReCalculateInten(popnum)                //skip the calculation if I pass in  zero
1174        return(0)
1175End
1176
1177
1178// calculate the intensity
1179// - either do MC or the straight calculation.
1180//
1181// *****currently the smeared calculation is turned off, by not prepending "fSmeared" to the FUNCREFs)
1182// --- the random deviate can't be calculated for the smeared model since I don't know the resolution
1183// function over an infinite q-range, just the detector. Maybe the interpolation is OK, but I
1184// don't really have a good way of testing this. Also, the resolution calculation explicitly multiplies
1185// by fShad, and this wrecks the random deviate calculation. The 2D looks great, but the probabilities
1186// are all wrong. fShad is only appropriate post-simulation.
1187//
1188//
1189// (--NO--) ALWAYS DOES THE RESOLUTION SMEARED CALCULATION
1190//  even though the unsmeared model is plotted. this is more realistic to
1191//  present to users that are planning base on what they will see in an experiment.
1192//
1193// some bits of the calculation are in a root:Simulation folder that are needed for the resolution smeared calculation
1194// all other bits (as possible) are in the SAS folder (a 2D work folder)
1195//
1196// passing in one does the calculation, "normal" or MC, depending on the global. Normal calculation is the default
1197// passing in zero from a control skips the calculation
1198//
1199Function ReCalculateInten(doIt)
1200        Variable doIt
1201       
1202        if(doIt==0)                     
1203                return(0)
1204        endif
1205       
1206        // update the wave with the beamstop diameter here, since I don't know what
1207        // combinations of parameters will change the BS - but anytime the curve is
1208        // recalculated, or the text displayed, the right BS must be present
1209        beamstopDiam()
1210       
1211        // generate the resolution waves first, so they are present for a smearing calculation
1212        // average the "fake" 2d data now to generate the smearing information
1213        S_CircularAverageTo1D("SAS")
1214        WAVE aveint=root:Packages:NIST:SAS:aveint
1215        WAVE qval=root:Packages:NIST:SAS:qval
1216        WAVE sigave=root:Packages:NIST:SAS:sigave
1217        WAVE SigmaQ=root:Packages:NIST:SAS:sigmaQ
1218        WAVE qbar=root:Packages:NIST:SAS:qbar
1219        WAVE fSubS=root:Packages:NIST:SAS:fSubS
1220       
1221        //generate a "fake" 1d data folder/set named "Simulation"
1222        Fake1DDataFolder(qval,aveint,sigave,sigmaQ,qbar,fSubs,"Simulation")
1223       
1224        // do the simulation here, or not
1225        Variable r1,xCtr,yCtr,sdd,pixSize,wavelength
1226        String coefStr,abortStr,str
1227
1228        // now the cases are: simulation (0|1), 1D (default) or 2D (hidden)
1229        NVAR doMonteCarlo = root:Packages:NIST:SAS:gDoMonteCarlo                // == 1 if 2D MonteCarlo set by hidden flag
1230        NVAR doSimulation = root:Packages:NIST:SAS:doSimulation         // == 1 if 1D simulated data, 0 if other from the checkbox
1231        NVAR g_estimateOnly = root:Packages:NIST:SAS:g_estimateOnly             // == 1 for just a time estimate, == 0 (default) to do the simulation
1232        NVAR g_estimatedMCTime = root:Packages:NIST:SAS:g_estimatedMCTime               // estimated MC sim time
1233        SVAR funcStr = root:Packages:NIST:SAS:gFuncStr          //set by the popup
1234
1235        if(doSimulation == 1)
1236                if(doMonteCarlo == 1)
1237                        //2D simulation (in MultiScatter_MonteCarlo_2D.ipf)
1238                       
1239                        g_estimatedMCTime = Simulate_2D_MC(funcStr,aveint,qval,sigave,sigmaq,qbar,fsubs,g_estimateOnly)
1240                       
1241                        //end 2D simulation
1242                else
1243                        //1D simulation
1244                       
1245                        if(exists(funcStr) != 0)
1246                               
1247                                Simulate_1D(funcStr,aveint,qval,sigave,sigmaq,qbar,fsubs)
1248                               
1249                        else
1250                                //no function plotted, no simulation can be done
1251                                DoAlert 0,"No function is selected or plotted, so no simulation is done. The default Debye function is used."
1252               
1253                                aveint = S_Debye(1000,100,0.0,qval)
1254                                aveint *= fSubS         // multiply either estimate by beamstop shadowing
1255                                sigave = 0              //reset for model calculation
1256                        endif
1257                       
1258                endif // end 1D simulation
1259        else
1260                //no simulation
1261               
1262                aveint = S_Debye(1000,100,0.0,qval)
1263                aveint *= fSubS         // multiply either estimate by beamstop shadowing
1264                sigave = 0              //reset for model calculation
1265               
1266                //end no simulation     
1267        endif
1268       
1269
1270        //display the configuration text in a separate notebook
1271        DisplayConfigurationText()
1272       
1273        return(0)
1274End
1275
1276//freezes the current configuration
1277// -1- duplicates the trace on the graph
1278// -2- copies the configuration text to a second notebook window for printing
1279Function FreezeButtonProc(ctrlName) : ButtonControl
1280        String ctrlName
1281
1282        String str=""
1283        NVAR ct=root:Packages:NIST:SAS:gFreezeCount
1284       
1285       
1286        SetDataFolder root:Packages:NIST:SAS
1287       
1288        Duplicate/O aveint,$("aveint_"+num2str(ct))
1289        Duplicate/O qval,$("qval_"+num2str(ct))
1290        Duplicate/O sigave,$("sigave_"+num2str(ct))
1291        Appendtograph $("aveint_"+num2str(ct)) vs $("qval_"+num2str(ct))
1292        ModifyGraph mode($("aveint_"+num2str(ct)))=3
1293        ModifyGraph marker($("aveint_"+num2str(ct)))=19
1294        ModifyGraph msize($("aveint_"+num2str(ct)))=2
1295        ErrorBars/T=0 $("aveint_"+num2str(ct)) Y,wave=($("sigave_"+num2str(ct)),$("sigave_"+num2str(ct)))
1296       
1297        switch(mod(ct,10))      // 10 different colors - black is the unfrozen color
1298                case 0:
1299                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(65535,16385,16385)
1300                        break
1301                case 1:
1302                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(2,39321,1)
1303                        break
1304                case 2:
1305                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(0,0,65535)
1306                        break
1307                case 3:
1308                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(39321,1,31457)
1309                        break
1310                case 4:
1311                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(48059,48059,48059)
1312                        break
1313                case 5:
1314                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(65535,32768,32768)
1315                        break
1316                case 6:
1317                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(0,65535,0)
1318                        break
1319                case 7:
1320                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(16385,65535,65535)
1321                        break
1322                case 8:
1323                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(65535,32768,58981)
1324                        break
1325                case 9:
1326                        ModifyGraph rgb($("aveint_"+num2str(ct)))=(36873,14755,58982)
1327                        break
1328        endswitch
1329       
1330        NVAR doTraceOffset = root:Packages:NIST:SAS:gDoTraceOffset
1331        NVAR offset = root:Packages:NIST:SAS:gModelOffsetFactor
1332        if(doTraceOffset)
1333                offset = 2^ct
1334                //multiply by current offset (>=1)
1335                Wave inten = $("aveint_"+num2str(ct))
1336                inten *= offset
1337                //      Print "new offset = ",offset
1338        endif
1339       
1340        ct +=1
1341        SetDataFolder root:
1342       
1343        // create or append the configuration to a new notebook
1344        if(WinType("Saved_Configurations")==0)
1345                NewNotebook/F=1/N=Saved_Configurations /W=(480,400,880,725)
1346                DoWindow/F SASCALC              //return focus to SASCALC window
1347        endif
1348        //append the text
1349        sprintf str,"\rConfiguration #%d\r",ct-1
1350        Notebook Saved_Configurations showRuler=0,defaultTab=20,selection={endOfFile, endOfFile}
1351        Notebook Saved_Configurations font="Monaco",fSize=10,fstyle=1,text=str          //bold
1352        Notebook Saved_Configurations font="Monaco",fSize=10,fstyle=0,text=SetConfigurationText()
1353        return(0)
1354End
1355
1356//clears the frozen traces on the graph, asks if you want to clear the saved text
1357//
1358Function S_ClearButtonProc(ctrlName) : ButtonControl
1359        String ctrlName
1360
1361        NVAR ct=root:Packages:NIST:SAS:gFreezeCount
1362        Variable ii
1363        Setdatafolder root:Packages:NIST:SAS
1364        for(ii=ct-1;ii>=1;ii-=1)
1365        //remove all traces, replace aveint
1366        // kill all waves _ct
1367                RemoveFromGraph $("aveint_"+num2str(ii))
1368                Killwaves/Z $("aveint_"+num2str(ii)),$("qval_"+num2str(ii))
1369        endfor
1370        ct=1
1371        setdatafolder root:
1372       
1373        DoAlert 1,"Do you also want to clear the \"Saved Configurations\" text?"
1374        if(V_flag == 1)         // yes
1375                DoWindow/K Saved_Configurations
1376        endif
1377       
1378        //reset offset value
1379        NVAR offset = root:Packages:NIST:SAS:gModelOffsetFactor
1380        offset = 1
1381        ReCalculateInten(1)
1382        return(0)
1383End
1384
1385
1386///////////////////////////////////////////////////////////
1387// 19MAR07 uses correction for beamstop diameter projection to get shadow factor correct
1388//
1389Function S_CircularAverageTo1D(type)
1390        String type
1391       
1392        Variable isCircular = 1
1393       
1394        //type is the data type to do the averaging on, and will be set as the current folder
1395        //get the current displayed data (so the correct folder is used)
1396        String destPath = "root:Packages:NIST:"+type
1397       
1398        //
1399        Variable xcenter,ycenter,x0,y0,sx,sx3,sy,sy3,dtsize,dtdist,dr,ddr
1400        Variable lambda,trans
1401        WAVE reals = $(destPath + ":RealsRead")
1402//      WAVE/T textread = $(destPath + ":TextRead")
1403//      String fileStr = textread[3]
1404       
1405        // center of detector, for non-linear corrections
1406        Variable pixelsX=128,pixelsY=128
1407       
1408        xcenter = pixelsX/2 + 0.5               // == 64.5 for 128x128 Ordela
1409        ycenter = pixelsY/2 + 0.5               // == 64.5 for 128x128 Ordela
1410       
1411        // beam center, in pixels
1412        x0 = reals[16]
1413        y0 = reals[17]
1414        //detector calibration constants
1415        sx = reals[10]          //mm/pixel (x)
1416        sx3 = reals[11]         //nonlinear coeff
1417        sy = reals[13]          //mm/pixel (y)
1418        sy3 = reals[14]         //nonlinear coeff
1419       
1420        dtsize = 10*reals[20]           //det size in mm
1421        dtdist = 1000*reals[18]         // det distance in mm
1422       
1423        NVAR binWidth=root:Packages:NIST:gBinWidth
1424//      Variable binWidth = 1
1425       
1426        dr = binWidth           // ***********annulus width set by user, default is one***********
1427        ddr = dr*sx             //step size, in mm (this value should be passed to the resolution calculation, not dr 18NOV03)
1428               
1429        Variable rcentr,large_num,small_num,dtdis2,nq,xoffst,dxbm,dybm,ii
1430        Variable phi_rad,dphi_rad,phi_x,phi_y
1431        Variable forward,mirror
1432       
1433        String side = "both"
1434//      side = StringByKey("SIDE",keyListStr,"=",";")
1435//      Print "side = ",side
1436       
1437        /// data wave is data in the current folder which was set at the top of the function
1438        WAVE data=$(destPath + ":data")
1439
1440// fake mask that uses all of the detector
1441        Make/D/O/N=(pixelsX,pixelsY) $(destPath + ":mask")
1442        Wave mask = $(destPath + ":mask")
1443        mask = 0
1444        //three pixels all around
1445        mask[0,2][] = 1
1446        mask[125,127][] = 1
1447        mask[][0,2] = 1
1448        mask[][125,127] = 1
1449        //
1450        //pixels within rcentr of beam center are broken into 9 parts (units of mm)
1451        rcentr = 100            //original
1452//      rcentr = 0
1453        // values for error if unable to estimate value
1454        //large_num = 1e10
1455        large_num = 1           //1e10 value (typically sig of last data point) plots poorly, arb set to 1
1456        small_num = 1e-10
1457       
1458        // output wave are expected to exist (?) initialized to zero, what length?
1459        // 200 points on VAX --- use 300 here, or more if SAXS data is used with 1024x1024 detector (1000 pts seems good)
1460        Variable defWavePts=500
1461        Make/O/D/N=(defWavePts) $(destPath + ":qval"),$(destPath + ":aveint")
1462        Make/O/D/N=(defWavePts) $(destPath + ":ncells"),$(destPath + ":dsq"),$(destPath + ":sigave")
1463        Make/O/D/N=(defWavePts) $(destPath + ":SigmaQ"),$(destPath + ":fSubS"),$(destPath + ":QBar")
1464
1465        WAVE qval = $(destPath + ":qval")
1466        WAVE aveint = $(destPath + ":aveint")
1467        WAVE ncells = $(destPath + ":ncells")
1468        WAVE dsq = $(destPath + ":dsq")
1469        WAVE sigave = $(destPath + ":sigave")
1470        WAVE qbar = $(destPath + ":QBar")
1471        WAVE sigmaq = $(destPath + ":SigmaQ")
1472        WAVE fsubs = $(destPath + ":fSubS")
1473       
1474        qval = 0
1475        aveint = 0
1476        ncells = 0
1477        dsq = 0
1478        sigave = 0
1479        qbar = 0
1480        sigmaq = 0
1481        fsubs = 0
1482
1483        dtdis2 = dtdist^2
1484        nq = 1
1485        xoffst=0
1486        //distance of beam center from detector center
1487        dxbm = S_FX(x0,sx3,xcenter,sx)
1488        dybm = S_FY(y0,sy3,ycenter,sy)
1489               
1490        //BEGIN AVERAGE **********
1491        Variable xi,dxi,dx,jj,data_pixel,yj,dyj,dy,mask_val=0.1
1492        Variable dr2,nd,fd,nd2,ll,kk,dxx,dyy,ir,dphi_p
1493       
1494        // IGOR arrays are indexed from [0][0], FORTAN from (1,1) (and the detector too)
1495        // loop index corresponds to FORTRAN (old code)
1496        // and the IGOR array indices must be adjusted (-1) to the correct address
1497        ii=1
1498        do
1499                xi = ii
1500                dxi = S_FX(xi,sx3,xcenter,sx)
1501                dx = dxi-dxbm           //dx and dy are in mm
1502               
1503                jj = 1
1504                do
1505                        data_pixel = data[ii-1][jj-1]           //assign to local variable
1506                        yj = jj
1507                        dyj = S_FY(yj,sy3,ycenter,sy)
1508                        dy = dyj - dybm
1509                        if(!(mask[ii-1][jj-1]))                 //masked pixels = 1, skip if masked (this way works...)
1510                                dr2 = (dx^2 + dy^2)^(0.5)               //distance from beam center NOTE dr2 used here - dr used above
1511                                if(dr2>rcentr)          //keep pixel whole
1512                                        nd = 1
1513                                        fd = 1
1514                                else                            //break pixel into 9 equal parts
1515                                        nd = 3
1516                                        fd = 2
1517                                endif
1518                                nd2 = nd^2
1519                                ll = 1          //"el-el" loop index
1520                                do
1521                                        dxx = dx + (ll - fd)*sx/3
1522                                        kk = 1
1523                                        do
1524                                                dyy = dy + (kk - fd)*sy/3
1525                                                if(isCircular)
1526                                                        //circular average, use all pixels
1527                                                        //(increment)
1528                                                        nq = S_IncrementPixel(data_pixel,ddr,dxx,dyy,aveint,dsq,ncells,nq,nd2)
1529                                                else
1530                                                        //a sector average - determine azimuthal angle
1531                                                        dphi_p = S_dphi_pixel(dxx,dyy,phi_x,phi_y)
1532                                                        if(dphi_p < dphi_rad)
1533                                                                forward = 1                     //within forward sector
1534                                                        else
1535                                                                forward = 0
1536                                                        Endif
1537                                                        if((Pi - dphi_p) < dphi_rad)
1538                                                                mirror = 1              //within mirror sector
1539                                                        else
1540                                                                mirror = 0
1541                                                        Endif
1542                                                        //check if pixel lies within allowed sector(s)
1543                                                        if(cmpstr(side,"both")==0)              //both sectors
1544                                                                if ( mirror || forward)
1545                                                                        //increment
1546                                                                        nq = S_IncrementPixel(data_pixel,ddr,dxx,dyy,aveint,dsq,ncells,nq,nd2)
1547                                                                Endif
1548                                                        else
1549                                                                if(cmpstr(side,"right")==0)             //forward sector only
1550                                                                        if(forward)
1551                                                                                //increment
1552                                                                                nq = S_IncrementPixel(data_pixel,ddr,dxx,dyy,aveint,dsq,ncells,nq,nd2)
1553                                                                        Endif
1554                                                                else                    //mirror sector only
1555                                                                        if(mirror)
1556                                                                                //increment
1557                                                                                nq = S_IncrementPixel(data_pixel,ddr,dxx,dyy,aveint,dsq,ncells,nq,nd2)
1558                                                                        Endif
1559                                                                Endif
1560                                                        Endif           //allowable sectors
1561                                                Endif   //circular or sector check
1562                                                kk+=1
1563                                        while(kk<=nd)
1564                                        ll += 1
1565                                while(ll<=nd)
1566                        Endif           //masked pixel check
1567                        jj += 1
1568                while (jj<=pixelsY)
1569                ii += 1
1570        while(ii<=pixelsX)              //end of the averaging
1571               
1572        //compute q-values and errors
1573        Variable ntotal,rr,theta,avesq,aveisq,var
1574       
1575        lambda = reals[26]
1576        ntotal = 0
1577        kk = 1
1578        do
1579                rr = (2*kk-1)*ddr/2
1580                theta = 0.5*atan(rr/dtdist)
1581                qval[kk-1] = (4*Pi/lambda)*sin(theta)
1582                if(ncells[kk-1] == 0)
1583                        //no pixels in annuli, data unknown
1584                        aveint[kk-1] = 0
1585                        sigave[kk-1] = large_num
1586                else
1587                        if(ncells[kk-1] <= 1)
1588                                //need more than one pixel to determine error
1589                                aveint[kk-1] = aveint[kk-1]/ncells[kk-1]
1590                                sigave[kk-1] = large_num
1591                        else
1592                                //assume that the intensity in each pixel in annuli is normally
1593                                // distributed about mean...
1594                                aveint[kk-1] = aveint[kk-1]/ncells[kk-1]
1595                                avesq = aveint[kk-1]^2
1596                                aveisq = dsq[kk-1]/ncells[kk-1]
1597                                var = aveisq-avesq
1598                                if(var<=0)
1599                                        sigave[kk-1] = small_num
1600                                else
1601                                        sigave[kk-1] = sqrt(var/(ncells[kk-1] - 1))
1602                                endif
1603                        endif
1604                endif
1605                ntotal += ncells[kk-1]
1606                kk+=1
1607        while(kk<=nq)
1608       
1609        //Print "NQ = ",nq
1610        // data waves were defined as 300 points (=defWavePts), but now have less than that (nq) points
1611        // use DeletePoints to remove junk from end of waves
1612        //WaveStats would be a more foolproof implementation, to get the # points in the wave
1613        Variable startElement,numElements
1614        startElement = nq
1615        numElements = defWavePts - startElement
1616        DeletePoints startElement,numElements, qval,aveint,ncells,dsq,sigave
1617       
1618        //////////////end of VAX sector_ave()
1619               
1620
1621
1622// ***************************************************************
1623//
1624// Do the extra 3 columns of resolution calculations starting here.
1625//
1626// ***************************************************************
1627
1628        Variable L2 = reals[18]
1629        Variable BS = reals[21]         //this the diameter is stored in mm
1630//  SRK - why was I using the projected diameter of the beam stop?? I added a step at the beginning of every recalculation
1631// of the intensity to get the right beamstop diameter into RealsRead...
1632//      Variable BS = beamstopDiamProjection(1) * 10            //calculated projection in cm *10 = mm
1633        Variable S1 = reals[23]
1634        Variable S2 = reals[24]
1635        Variable L1 = reals[25]
1636        lambda = reals[26]
1637        Variable lambdaWidth = reals[27]
1638
1639        Variable DDet, apOff
1640        //typical value for NG3 and NG7 - distance between sample aperture and sample in (cm)
1641        apOff=5.0
1642        // hard wire value for Ordela detectors
1643        DDet = 0.5              // resolution in cm
1644        //      String detStr=textRead[9]
1645        //      DDet = DetectorPixelResolution(fileStr,detStr)          //needs detector type and beamline
1646
1647        //Go from 0 to nq doing the calc for all three values at
1648        //every Q value
1649
1650        ii=0
1651        Variable ret1,ret2,ret3
1652        do
1653                S_getResolution(qval[ii],lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,ddr,ret1,ret2,ret3)
1654                sigmaq[ii] = ret1       //res_wave[0]
1655                qbar[ii] = ret2         //res_wave[1]
1656                fsubs[ii] = ret3                //res_wave[2]
1657                ii+=1
1658        while(ii<nq)
1659        DeletePoints startElement,numElements, sigmaq, qbar, fsubs
1660
1661        fsubs += 1e-8           //keep the values from being too small
1662
1663// End of resolution calculations
1664// ***************************************************************
1665       
1666        //get rid of the default mask, if one was created (it is in the current folder)
1667        //don't just kill "mask" since it might be pointing to the one in the MSK folder
1668        Killwaves/Z $(destPath+":mask")
1669       
1670        //return to root folder (redundant)
1671        SetDataFolder root:
1672       
1673        Return 0
1674End
1675
1676//returns nq, new number of q-values
1677//arrays aveint,dsq,ncells are also changed by this function
1678//
1679Function S_IncrementPixel(dataPixel,ddr,dxx,dyy,aveint,dsq,ncells,nq,nd2)
1680        Variable dataPixel,ddr,dxx,dyy
1681        Wave aveint,dsq,ncells
1682        Variable nq,nd2
1683       
1684        Variable ir
1685       
1686        ir = trunc(sqrt(dxx*dxx+dyy*dyy)/ddr)+1
1687        if (ir>nq)
1688                nq = ir         //resets maximum number of q-values
1689        endif
1690        aveint[ir-1] += dataPixel/nd2           //ir-1 must be used, since ir is physical
1691        dsq[ir-1] += dataPixel*dataPixel/nd2
1692        ncells[ir-1] += 1/nd2
1693       
1694        Return nq
1695End
1696
1697//function determines azimuthal angle dphi that a vector connecting
1698//center of detector to pixel makes with respect to vector
1699//at chosen azimuthal angle phi -> [cos(phi),sin(phi)] = [phi_x,phi_y]
1700//dphi is always positive, varying from 0 to Pi
1701//
1702Function S_dphi_pixel(dxx,dyy,phi_x,phi_y)
1703        Variable dxx,dyy,phi_x,phi_y
1704       
1705        Variable val,rr,dot_prod
1706       
1707        rr = sqrt(dxx^2 + dyy^2)
1708        dot_prod = (dxx*phi_x + dyy*phi_y)/rr
1709        //? correct for roundoff error? - is this necessary in IGOR, w/ double precision?
1710        if(dot_prod > 1)
1711                dot_prod =1
1712        Endif
1713        if(dot_prod < -1)
1714                dot_prod = -1
1715        Endif
1716       
1717        val = acos(dot_prod)
1718       
1719        return val
1720
1721End
1722
1723//calculates the x distance from the center of the detector, w/nonlinear corrections
1724//
1725Function S_FX(xx,sx3,xcenter,sx)               
1726        Variable xx,sx3,xcenter,sx
1727       
1728        Variable retval
1729       
1730        retval = sx3*tan((xx-xcenter)*sx/sx3)
1731        Return retval
1732End
1733
1734//calculates the y distance from the center of the detector, w/nonlinear corrections
1735//
1736Function S_FY(yy,sy3,ycenter,sy)               
1737        Variable yy,sy3,ycenter,sy
1738       
1739        Variable retval
1740       
1741        retval = sy3*tan((yy-ycenter)*sy/sy3)
1742        Return retval
1743End
1744
1745//**********************
1746// Resolution calculation - used by the averaging routines
1747// to calculate the resolution function at each q-value
1748// - the return value is not used
1749//
1750// equivalent to John's routine on the VAX Q_SIGMA_AVE.FOR
1751// Incorporates eqn. 3-15 from J. Appl. Cryst. (1995) v. 28 p105-114
1752//
1753Function/S S_getResolution(inQ,lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,del_r,SigmaQ,QBar,fSubS)
1754        Variable inQ, lambda, lambdaWidth, DDet, apOff, S1, S2, L1, L2, BS, del_r
1755        Variable &fSubS, &QBar, &SigmaQ         //these are the output quantities at the input Q value
1756       
1757        //lots of calculation variables
1758        Variable a2, q_small, lp, v_lambda, v_b, v_d, vz, yg, v_g
1759        Variable r0, delta, inc_gamma, fr, fv, rmd, v_r1, rm, v_r
1760
1761        //Constants
1762        //Variable del_r = .1
1763        Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron
1764        Variable g = 981.0                              //gravity acceleration [cm/s^2]
1765
1766        String results
1767        results ="Failure"
1768       
1769        NVAR usingLenses = root:Packages:NIST:SAS:gUsingLenses
1770
1771        //rename for working variables,  these must be gotten from global
1772        //variables
1773
1774//      Variable wLam, wLW, wL1, wL2, wS1, wS2
1775//      Variable wBS, wDDet, wApOff
1776//      wLam = lambda
1777//      wLW = lambdaWidth
1778//      wDDet = DDet
1779//      wApOff = apOff
1780        S1 *= 0.5*0.1                   //convert to radius and [cm]
1781        S2 *= 0.5*0.1
1782
1783        L1 *= 100.0                     // [cm]
1784        L1 -= apOff                             //correct the distance
1785
1786        L2 *= 100.0
1787        L2 += apOff
1788
1789        BS *= 0.5*0.1                   //convert to radius and [cm]
1790        del_r *= 0.1                            //width of annulus, convert mm to [cm]
1791       
1792        //Start resolution calculation
1793        a2 = S1*L2/L1 + S2*(L1+L2)/L1
1794        q_small = 2.0*Pi*(BS-a2)*(1.0-lambdaWidth)/(lambda*L2)
1795        lp = 1.0/( 1.0/L1 + 1.0/L2)
1796
1797        v_lambda = lambdaWidth^2/6.0
1798       
1799        if(usingLenses==1)                      //SRK 2007
1800                v_b = 0.25*(S1*L2/L1)^2 +0.25*(2/3)*(lambdaWidth/lambda)^2*(S2*L2/lp)^2         //correction to 2nd term
1801        else
1802                v_b = 0.25*(S1*L2/L1)^2 +0.25*(S2*L2/lp)^2              //original form
1803        endif
1804       
1805        v_d = (DDet/2.3548)^2 + del_r^2/12.0
1806        vz = vz_1 / lambda
1807        yg = 0.5*g*L2*(L1+L2)/vz^2
1808        v_g = 2.0*(2.0*yg^2*v_lambda)                                   //factor of 2 correction, B. Hammouda, 2007
1809
1810        r0 = L2*tan(2.0*asin(lambda*inQ/(4.0*Pi) ))
1811        delta = 0.5*(BS - r0)^2/v_d
1812
1813        if (r0 < BS)
1814                inc_gamma=exp(gammln(1.5))*(1-gammp(1.5,delta))
1815        else
1816                inc_gamma=exp(gammln(1.5))*(1+gammp(1.5,delta))
1817        endif
1818
1819        fSubS = 0.5*(1.0+erf( (r0-BS)/sqrt(2.0*v_d) ) )
1820        if (fSubS <= 0.0)
1821                fSubS = 1.e-10
1822        endif
1823        fr = 1.0 + sqrt(v_d)*exp(-1.0*delta) /(r0*fSubS*sqrt(2.0*Pi))
1824        fv = inc_gamma/(fSubS*sqrt(Pi)) - r0^2*(fr-1.0)^2/v_d
1825
1826        rmd = fr*r0
1827        v_r1 = v_b + fv*v_d +v_g
1828
1829        rm = rmd + 0.5*v_r1/rmd
1830        v_r = v_r1 - 0.5*(v_r1/rmd)^2
1831        if (v_r < 0.0)
1832                v_r = 0.0
1833        endif
1834        QBar = (4.0*Pi/lambda)*sin(0.5*atan(rm/L2))
1835        SigmaQ = QBar*sqrt(v_r/rmd^2 +v_lambda)
1836       
1837        results = "success"
1838        Return results
1839End
1840
1841Function S_Debye(scale,rg,bkg,x)
1842        Variable scale,rg,bkg
1843        Variable x
1844       
1845        // variables are:
1846        //[0] scale factor
1847        //[1] radius of gyration [A]
1848        //[2] background        [cm-1]
1849       
1850        // calculates (scale*debye)+bkg
1851        Variable Pq,qr2
1852       
1853        qr2=(x*rg)^2
1854        Pq = 2*(exp(-(qr2))-1+qr2)/qr2^2
1855       
1856        //scale
1857        Pq *= scale
1858        // then add in the background
1859        return (Pq+bkg)
1860End
1861
1862
1863Function S_SphereForm(scale,radius,delrho,bkg,x)                               
1864        Variable scale,radius,delrho,bkg
1865        Variable x
1866       
1867        // variables are:                                                       
1868        //[0] scale
1869        //[1] radius (A)
1870        //[2] delrho (A-2)
1871        //[3] background (cm-1)
1872       
1873//      Variable scale,radius,delrho,bkg                               
1874//      scale = w[0]
1875//      radius = w[1]
1876//      delrho = w[2]
1877//      bkg = w[3]
1878       
1879       
1880        // calculates scale * f^2/Vol where f=Vol*3*delrho*((sin(qr)-qrcos(qr))/qr^3
1881        // and is rescaled to give [=] cm^-1
1882       
1883        Variable bes,f,vol,f2
1884        //
1885        //handle q==0 separately
1886        If(x==0)
1887                f = 4/3*pi*radius^3*delrho*delrho*scale*1e8 + bkg
1888                return(f)
1889        Endif
1890       
1891        bes = 3*(sin(x*radius)-x*radius*cos(x*radius))/x^3/radius^3
1892        vol = 4*pi/3*radius^3
1893        f = vol*bes*delrho              // [=] A
1894        // normalize to single particle volume, convert to 1/cm
1895        f2 = f * f / vol * 1.0e8                // [=] 1/cm
1896       
1897        return (scale*f2+bkg)   // Scale, then add in the background
1898       
1899End
1900
1901Function/S SetConfigurationText()
1902
1903        String str="",temp
1904
1905        SetDataFolder root:Packages:NIST:SAS
1906       
1907        NVAR numberOfGuides=gNg
1908        NVAR gTable=gTable              //2=chamber, 1=table
1909        NVAR wavelength=gLambda
1910        NVAR lambdaWidth=gDeltaLambda
1911//      NVAR instrument = instrument
1912        NVAR L2diff = L2diff
1913   NVAR lens = root:Packages:NIST:SAS:gUsingLenses
1914        SVAR/Z aStr = root:Packages:NIST:gAngstStr
1915        SVAR selInstr = root:Packages:NIST:SAS:gInstStr
1916
1917       
1918        sprintf temp,"Source Aperture Diameter =\t\t%6.2f cm\r",sourceApertureDiam()
1919        str += temp
1920        sprintf temp,"Source to Sample =\t\t\t\t%6.0f cm\r",sourceToSampleDist()
1921        str += temp
1922        sprintf temp,"Sample Aperture to Detector =\t%6.0f cm\r",sampleToDetectorDist() + L2diff
1923        str += temp
1924        sprintf temp,"Beam diameter =\t\t\t\t\t%6.2f cm\r",beamDiameter("maximum")
1925        str += temp
1926        sprintf temp,"Beamstop diameter =\t\t\t\t%6.2f inches\r",beamstopDiam()/2.54
1927        str += temp
1928        sprintf temp,"Minimum Q-value =\t\t\t\t%6.4f 1/%s (sigQ/Q = %3.1f %%)\r",qMin(),aStr,deltaQ(qMin())
1929        str += temp
1930        sprintf temp,"Maximum Horizontal Q-value =\t%6.4f 1/%s\r",qMaxHoriz(),aStr
1931        str += temp
1932        sprintf temp,"Maximum Vertical Q-value =\t\t%6.4f 1/%s\r",qMaxVert(),aStr
1933        str += temp
1934        sprintf temp,"Maximum Q-value =\t\t\t\t%6.4f 1/%s (sigQ/Q = %3.1f %%)\r",qMaxCorner(),aStr,deltaQ(qMaxCorner())
1935        str += temp
1936        sprintf temp,"Beam Intensity =\t\t\t\t%.0f counts/s\r",beamIntensity()
1937        str += temp
1938        sprintf temp,"Figure of Merit =\t\t\t\t%3.3g %s^2/s\r",figureOfMerit(),aStr
1939        str += temp
1940        sprintf temp,"Attenuator transmission =\t\t%3.3g = Atten # %d\r",attenuatorTransmission(),attenuatorNumber()
1941        str += temp
1942//     
1943//      // add text of the user-edited values
1944//      //
1945        sprintf temp,"***************** %s *** %s *****************\r",selInstr,selInstr
1946        str += temp
1947        sprintf temp,"Sample Aperture Diameter =\t\t\t\t%.2f cm\r",sampleApertureDiam()
1948        str += temp
1949        sprintf temp,"Number of Guides =\t\t\t\t\t\t%d \r", numberOfGuides
1950        str += temp
1951        sprintf temp,"Sample Chamber to Detector =\t\t\t%.1f cm\r", chamberToDetectorDist()
1952        str += temp
1953        if(gTable==1)
1954                sprintf temp,"Sample Position is \t\t\t\t\t\tHuber\r"
1955        else
1956                sprintf temp,"Sample Position is \t\t\t\t\t\tChamber\r"
1957        endif
1958        str += temp
1959        sprintf temp,"Detector Offset =\t\t\t\t\t\t%.1f cm\r", detectorOffset()
1960        str += temp
1961        sprintf temp,"Neutron Wavelength =\t\t\t\t\t%.2f %s\r", wavelength,aStr
1962        str += temp
1963        sprintf temp,"Wavelength Spread, FWHM =\t\t\t\t%.3f\r", lambdaWidth
1964        str += temp
1965        sprintf temp,"Sample Aperture to Sample Position =\t%.2f cm\r", L2Diff
1966        str += temp
1967        if(lens==1)
1968                sprintf temp,"Lenses are IN\r"
1969        else
1970                sprintf temp,"Lenses are OUT\r"
1971        endif
1972        str += temp
1973       
1974   setDataFolder root:
1975   return str                   
1976End
1977
1978Function DisplayConfigurationText()
1979
1980        if(WinType("Trial_Configuration")==0)
1981                NewNotebook/F=0/K=1/N=Trial_Configuration /W=(480,44,880,369)
1982        endif
1983        //replace the text
1984        Notebook Trial_Configuration selection={startOfFile, endOfFile}
1985        Notebook Trial_Configuration font="Monaco",fSize=10,text=SetConfigurationText()
1986        return(0)
1987end
1988
1989//parses the control for A1 diam
1990// updates the wave
1991Function sourceApertureDiam()
1992        ControlInfo/W=SASCALC popup0
1993        Variable diam
1994        sscanf S_Value, "%g cm", diam
1995        WAVE rw=root:Packages:NIST:SAS:realsRead
1996        rw[23] = diam*10                        //sample aperture diameter in mm
1997        return(diam)
1998end
1999
2000// change the sample aperture to a non-standard value
2001Function SampleApOtherSetVarProc(ctrlName,varNum,varStr,varName) : SetVariableControl
2002        String ctrlName
2003        Variable varNum
2004        String varStr
2005        String varName
2006               
2007        WAVE rw=root:Packages:NIST:SAS:realsRead
2008        rw[24] = varNum                 //sample aperture diameter in mm
2009        ReCalculateInten(1)
2010        return(0)
2011End
2012
2013//parses the control for A2 diam
2014// updates the wave and global
2015// returns a2 in cm
2016Function sampleApertureDiam()
2017        ControlInfo/W=SASCALC popup0_1
2018       
2019//      Print "In sampleApertureDiam()"
2020        //set the global
2021        NVAR a2=root:Packages:NIST:SAS:gSamAp
2022       
2023        if(cmpstr(S_Value,"other") == 0)                // "other" selected
2024                //enable the setvar, diameter in mm!
2025                SetVariable setvar0_3,win=SASCALC, disable=0
2026                // read its value (a global)
2027                NVAR a2other = root:Packages:NIST:SAS:gSamApOther
2028                a2=a2other/10                           //a2 in cm
2029        else
2030                SetVariable setvar0_3,win=SASCALC, disable=1
2031                //1st item is 1/16", popup steps by 1/16"
2032                a2 = 2.54/16.0 * (V_Value)                      //convert to cm         
2033        endif
2034        WAVE rw=root:Packages:NIST:SAS:realsRead
2035        rw[24] = a2*10                  //sample aperture diameter in mm
2036       
2037        return(a2)
2038end
2039
2040//1=huber 2=chamber
2041Function tableposition()
2042        ControlInfo/W=SASCALC checkHuber
2043        if(V_Value)
2044                return(1)               //huber
2045        else
2046                return(2)               //chamber
2047        endif
2048End
2049
2050//compute SSD and update both the global and the wave
2051//
2052// for the 10m SANS, table position is inline, so S12 = 0
2053//
2054//
2055Function sourceToSampleDist()
2056
2057        NVAR NG=root:Packages:NIST:SAS:gNg
2058        NVAR S12 = root:Packages:NIST:SAS:S12
2059        NVAR L2Diff = root:Packages:NIST:SAS:L2Diff
2060        NVAR SSD = root:Packages:NIST:SAS:gSSD
2061//      NVAR instrument = root:Packages:NIST:SAS:instrument
2062        SVAR selInstr = root:Packages:NIST:SAS:gInstStr
2063
2064        Variable NGA_gap = 61.9         // extra distance between a1 and beginning of guide 1 on NGA
2065       
2066        strswitch(selInstr)     // string switch
2067                case "NG3":
2068                case "NG7":
2069                        // NG3 and NG7 are both the same
2070                        SSD = 1632 - 155*NG - s12*(2-tableposition()) - L2Diff
2071                        break
2072                case "NGA":
2073                        // 10m SANS handled differently
2074                        // s12 == 0 by definition
2075                        // -- 16JAN13 - these are now correct values
2076                        if(ng==0)
2077                                SSD = 513 - L2diff
2078                        else                   
2079                                SSD = 513 - NGA_gap - 150*NG - s12*(2-tableposition()) - L2Diff
2080                        endif
2081                        break
2082                default:
2083                        DoAlert 0,"No matching instrument! sourceToSampleDist"
2084        endswitch
2085
2086       
2087        WAVE rw=root:Packages:NIST:SAS:realsRead
2088        rw[25] = SSD/100                // in meters
2089        return(SSD)
2090End
2091
2092
2093// not part of SASCALC, but can be used to convert the SSD to number of guides
2094//
2095// SSD in meters
2096//
2097Function numGuides(SSD)
2098        variable SSD
2099       
2100        Variable Ng
2101       
2102        //NVAR instrument = root:Packages:NIST:SAS:instrument
2103        SVAR selInstr = root:Packages:NIST:SAS:gInstStr
2104
2105        Variable NGA_gap = 61.9         // extra distance between a1 and beginning of guide 1 on NGA
2106
2107        strswitch(selInstr)     // string switch
2108                case "NG3":
2109                case "NG7":
2110                        // NG3 and NG7 are both the same
2111                        Ng = SSD*100 + 5 - 1632
2112                        Ng /= -155
2113       
2114                        break
2115                case "NGA":
2116                        // 10m SANS handled differently
2117                        // -- 16JAN13 - these are now correct values
2118                        Ng = 513 - NGA_gap - SSD*100 -5
2119                        Ng /= 150
2120       
2121                        if(ng < 0)
2122                                ng=0
2123                        endif
2124                       
2125                        break
2126                default:
2127                       
2128        endswitch
2129//      Print "Ng = ",Ng
2130
2131        Ng = abs(round(Ng))
2132        return(Ng)
2133End
2134
2135
2136//returns the offset value
2137// slider and setVar are linked to the same global
2138// updates the wave and changes the beamcenter (x,y) in the wave
2139//
2140// does not account for the fall due to gravity in y-direction
2141Function detectorOffset()
2142       
2143        WAVE rw=root:Packages:NIST:SAS:RealsRead
2144        NVAR val = root:Packages:NIST:SAS:gOffset
2145        rw[19] = val            // already in cm
2146        //move the beamcenter, make it an integer value for the MC simulation
2147        rw[16] = 64  + round(2*rw[19]) + 0.5            //approximate beam X is 64 w/no offset, 114 w/25 cm offset
2148        rw[17] = 64     + 0.5 //typical value
2149
2150        return(val)
2151end
2152
2153//returns the detector distance (slider and setVar are linked by the global)
2154//
2155Function chamberToDetectorDist()
2156
2157        NVAR val = root:Packages:NIST:SAS:gDetDist
2158        return(val)
2159End
2160
2161//sets the SDD (slider and setVar are linked by the global and is the detector position
2162//  relative to the chamber)
2163// updates the wave
2164Function sampleToDetectorDist()
2165
2166        NVAR detDist = root:Packages:NIST:SAS:gDetDist
2167        NVAR S12 = root:Packages:NIST:SAS:S12
2168        WAVE rw=root:Packages:NIST:SAS:RealsRead       
2169        Variable SDD   
2170       
2171        SDD = detDist + s12*(2-tableposition())
2172        rw[18] = SDD/100                // convert to meters for header
2173        return(SDD)
2174End
2175
2176//direction = one of "vertical;horizontal;maximum;"
2177// all of this is bypassed if the lenses are in
2178//
2179Function beamDiameter(direction)
2180        String direction
2181
2182        NVAR lens = root:Packages:NIST:SAS:gUsingLenses
2183        if(lens)
2184                return sourceApertureDiam()
2185        endif
2186       
2187    Variable l1 = sourceToSampleDist()
2188    Variable l2 //= sampleAperToDetDist()
2189    Variable d1,d2,bh,bv,bm,umbra,a1,a2
2190   
2191    NVAR L2diff = root:Packages:NIST:SAS:L2diff
2192    NVAR lambda = root:Packages:NIST:SAS:gLambda
2193    NVAR lambda_width = root:Packages:NIST:SAS:gDeltaLambda
2194    NVAR bs_factor = root:Packages:NIST:SAS:bs_factor
2195   
2196    l2 = sampleToDetectorDist() + L2diff
2197    a1 = sourceApertureDiam()
2198    a2 = sampleApertureDiam()
2199   
2200    d1 = a1*l2/l1
2201    d2 = a2*(l1+l2)/l1
2202    bh = d1+d2          //beam size in horizontal direction
2203    umbra = abs(d1-d2)
2204    //vertical spreading due to gravity
2205    bv = bh + 1.25e-8*(l1+l2)*l2*lambda*lambda*lambda_width
2206    bm = (bs_factor*bh > bv) ? bs_factor*bh : bv //use the larger of horiz*safety or vertical
2207   
2208    strswitch(direction)        // string switch
2209        case "vertical":                // execute if case matches expression
2210                return(bv)
2211                break                                           // exit from switch
2212        case "horizontal":              // execute if case matches expression
2213                return(bh)
2214                break
2215        case "maximum":         // execute if case matches expression
2216                return(bm)
2217                break
2218        default:                                                        // optional default expression executed
2219                return(bm)                                              // when no case matches
2220    endswitch
2221End
2222
2223//on NG3 and NG7, allowable sizes are 1,2,3,4" diameter
2224//will return values larger than 4.0*2.54 if a larger beam is needed
2225//
2226// - in an approximate way, account for lenses
2227Function beamstopDiam()
2228
2229        NVAR yesLens = root:Packages:NIST:SAS:gUsingLenses
2230        Variable bm=0
2231        Variable bs=0.0
2232   
2233        if(yesLens)
2234                //bm = sourceApertureDiam()             //ideal result, not needed
2235                bs = 1                                                          //force the diameter to 1"
2236        else
2237                bm = beamDiameter("maximum")
2238                do
2239                bs += 1
2240           while ( (bs*2.54 < bm) || (bs > 30.0))                       //30 = ridiculous limit to avoid inf loop
2241        endif
2242
2243        //update the wave
2244        WAVE rw=root:Packages:NIST:SAS:realsRead
2245        rw[21] = bs*25.4                //store the BS diameter in mm
2246       
2247    return (bs*2.54)            //return diameter in cm, not inches for txt
2248End
2249
2250//returns the projected diameter of the beamstop at the anode plane.
2251// most noticeable at short SDD
2252//if flag == 0 use conservative estimate = largest diameter (for SASCALC, default)
2253//if flag != 0 use point aperture = average diameter (for resolution calculation)
2254Function beamstopDiamProjection(flag)
2255        Variable flag
2256       
2257        NVAR L2diff = root:Packages:NIST:SAS:L2diff
2258        Variable a2 = sampleApertureDiam()
2259        Variable bs = beamstopDiam()
2260        Variable l2, LB, BS_P
2261   
2262        l2 = sampleToDetectorDist() + L2diff
2263        LB = 20.1 + 1.61*BS                     //distance in cm from beamstop to anode plane (empirical)
2264        if(flag==0)
2265                BS_P = bs + (bs+a2)*lb/(l2-lb)          //diameter of shadow from parallax
2266        else
2267                BS_P = bs + bs*lb/(l2-lb)               //diameter of shadow, point A2
2268        endif
2269        return (bs_p)           //return projected diameter in cm
2270End
2271
2272// 19MAR07 - using correction from John for an estimate of the shadow of the beamstop
2273// at the detection plane. This is a noticeable effect at short SDD, where the projected
2274// diameter of the beamstop is much larger than the physical diameter.
2275Function qMin()
2276
2277    Variable l2s = sampleToDetectorDist()       //distance from sample to detector in cm
2278//    Variable bs = beamstopDiam()              //beamstop diameter in cm
2279    Variable bs_p = beamstopDiamProjection(0)           //projected beamstop diameter in cm
2280    NVAR lambda = root:Packages:NIST:SAS:gLambda
2281    NVAR d_det = root:Packages:NIST:SAS:d_det           //cm
2282    NVAR a_pixel = root:Packages:NIST:SAS:a_pixel       //cm
2283
2284    return( (pi/lambda)*(bs_p + d_det + a_pixel)/l2s )          //use bs_p rather than bs
2285   // return( (pi/lambda)*(bs + d_det + a_pixel)/l2s )          //use bs (incorrect)
2286End
2287
2288Function qMaxVert()
2289
2290    Variable theta
2291    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2292    NVAR lambda = root:Packages:NIST:SAS:gLambda
2293        NVAR det_width = root:Packages:NIST:SAS:det_width
2294       
2295    theta = atan( (det_width/2.0)/l2s )
2296   
2297    return ( 4.0*pi/lambda * sin(theta/2.0) )
2298end
2299
2300Function qMaxCorner()
2301
2302    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2303    Variable radial
2304    NVAR lambda = root:Packages:NIST:SAS:gLambda
2305        NVAR det_off = root:Packages:NIST:SAS:gOffset
2306        NVAR det_width = root:Packages:NIST:SAS:det_width
2307
2308    radial=sqrt( (0.5*det_width)*(0.5*det_width)+(0.5*det_width+det_off)*(0.5*det_width+det_off) )
2309   
2310    return ( 4*pi/lambda*sin(0.5*atan(radial/l2s)) )
2311End
2312
2313Function qMaxHoriz()
2314
2315    Variable theta
2316    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2317    NVAR lambda = root:Packages:NIST:SAS:gLambda
2318        NVAR det_off = root:Packages:NIST:SAS:gOffset
2319        NVAR det_width = root:Packages:NIST:SAS:det_width
2320
2321    theta = atan( ((det_width/2.0)+det_off)/l2s )       //from the instance variables
2322   
2323    return ( 4.0*pi/lambda * sin(theta/2.0) )
2324End
2325
2326// calculate sigma for the resolution function at either limit of q-range
2327Function deltaQ(atQ)
2328        Variable atQ
2329       
2330    Variable k02,lp,l1,l2,sig_02,sigQ2,a1,a2
2331    NVAR l2Diff = root:Packages:NIST:SAS:L2diff
2332        NVAR lambda = root:Packages:NIST:SAS:gLambda
2333        NVAR lambda_width = root:Packages:NIST:SAS:gDeltaLambda
2334        NVAR d_det = root:Packages:NIST:SAS:d_det
2335        NVAR del_r = root:Packages:NIST:SAS:del_r
2336       
2337       
2338    l1 = sourceToSampleDist()
2339    l2 = sampleToDetectorDist() + L2diff
2340    a1 = sourceApertureDiam()
2341    a2 = sampleApertureDiam()
2342   
2343    k02 = (6.2832/lambda)*(6.2832/lambda)
2344    lp = 1/(1/l1 + 1/l2)
2345   
2346    sig_02 = (0.25*a1/l1)*(0.25*a1/l1)
2347    sig_02 += (0.25*a2/lp)*(0.25*a2/lp)
2348    sig_02 += (d_det/(2.355*l2))*(d_det/(2.355*l2))
2349    sig_02 += (del_r/l2)*(del_r/l2)/12
2350    sig_02 *= k02
2351   
2352    sigQ2 = sig_02 + (atQ*lambda_width)*(atQ*lambda_width)/6
2353
2354    return(100*sqrt(sigQ2)/atQ)
2355End
2356
2357// updated with new flux numbers from John Barker
2358// NG3 - Feb 2009
2359// NG7 - July 2009
2360//
2361// guide loss has been changed to 0.95 rather than the old value of 0.95
2362//
2363// other values are changed in the initialization routines
2364//
2365Function beamIntensity()
2366
2367    Variable alpha,f,t,t4,t5,t6,as,solid_angle,l1,d2_phi
2368    Variable a1,a2,retVal
2369    SetDataFolder root:Packages:NIST:SAS
2370    NVAR l_gap=l_gap,guide_width =guide_width,ng = gNg
2371    NVAR lambda_t=lambda_t,b=b,c=c
2372    NVAR lambda=gLambda,t1=t1,t2=t2,t3=t3,phi_0=phi_0
2373    NVAR lambda_width=gDeltaLambda
2374    NVAR guide_loss=gGuide_loss
2375   
2376    l1 = sourceToSampleDist()
2377    a1 = sourceApertureDiam()
2378    a2 = sampleApertureDiam()
2379   
2380   
2381    alpha = (a1+a2)/(2*l1)      //angular divergence of beam
2382    f = l_gap*alpha/(2*guide_width)
2383    t4 = (1-f)*(1-f)
2384    t5 = exp(ng*ln(guide_loss)) // trans losses of guides in pre-sample flight
2385    t6 = 1 - lambda*(b-(ng/8)*(b-c))            //experimental correction factor
2386    t = t1*t2*t3*t4*t5*t6
2387   
2388    as = pi/4*a2*a2             //area of sample in the beam
2389    d2_phi = phi_0/(2*pi)
2390    d2_phi *= exp(4*ln(lambda_t/lambda))
2391    d2_phi *= exp(-1*(lambda_t*lambda_t/lambda/lambda))
2392
2393    solid_angle = pi/4* (a1/l1)*(a1/l1)
2394
2395    retVal = as * d2_phi * lambda_width * solid_angle * t
2396    SetDataFolder root:
2397    return (retVal)
2398end
2399
2400Function figureOfMerit()
2401
2402        Variable bi = beamIntensity()
2403        NVAR lambda = root:Packages:NIST:SAS:gLambda
2404       
2405    return (lambda*lambda*bi)
2406End
2407
2408//estimate the number of pixels in the beam, and enforce the maximum countrate per pixel (idmax)
2409Function attenuatorTransmission()
2410
2411    Variable num_pixels,i_pix           //i_pix = id in John's notation
2412    Variable bDiam = beamDiameter("horizontal") //!! note that prev calculations used bh (horizontal)
2413    Variable atten,a2
2414    SetDataFolder root:Packages:NIST:SAS
2415    NVAR a_pixel=a_pixel,idmax=idmax
2416   
2417    a2 = sampleApertureDiam()
2418   
2419    num_pixels = pi/4*(0.5*(a2+bDiam))*(0.5*(a2+bDiam))/a_pixel/a_pixel
2420    i_pix = ( beamIntensity() )/num_pixels
2421   
2422    atten = (i_pix < idmax) ? 1.0 : idmax/i_pix
2423    SetDataFolder root:
2424    return(atten)
2425End
2426
2427// TODO_10m: is this going to be close enough? Are the steps the same?
2428Function attenuatorNumber()
2429
2430    Variable atten = attenuatorTransmission()
2431    Variable af,nf,numAtten
2432    SetDataFolder root:Packages:NIST:SAS
2433    NVAR lambda=gLambda
2434   
2435    af = 0.498 + 0.0792*lambda - 1.66e-3*lambda*lambda
2436    nf = -ln(atten)/af          //floating point
2437   
2438    numAtten = trunc(nf) + 1                    //in c, (int)nf
2439    //correct for larger step thickness at n > 6
2440    if(numAtten > 6)
2441        numAtten = 7 + trunc( (numAtten-6)/2 )          //in c, numAtten = 7 + (int)( (numAtten-6)/2 )
2442    endif
2443   
2444    SetDatafolder root:
2445    return (numAtten)
2446End
Note: See TracBrowser for help on using the repository browser.