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

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

Changes to SASCALC and other locations in the code that identify the instrument from the account name that is stored in the file header. New designation of NGA for the 10-m SANS, and NGB is reserved for the NG3 30-m SANS when it is moved into the new guide hall. Changes to incorporate the 10m SANS are functional, but INCOMPLETE since many of the instrument details have not been filled in yet (they haven't been measured yet).

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