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

Last change on this file was 1088, checked in by srkline, 5 years ago

corrected bad logic in reduction when reducing to 2D. Nothing was wrong, just a meaningless error dialog was presented.

updated the beamstop sizes in SASCALC for the NGB-10m SANS instrument. Now the 4 beamstops are 1", 1.5", 2" and 3".

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