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

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

Updated the reduction Read/Write? and corresponding factility stubs to accomodate detector dead time written to the VAX header. It is currently not written to the header, but may be with NICE (hopefully).

With the move of NG3 SANS to CGB(upper), the NG3 designation in the account name [NGxSANSxx] has been replaced with CGB. Folders on charlotte and radio button on SASCALC are labeled "NGB30" to be more obvious which instrument it is.

FFT routines have minor typos cleaned up.

File size: 79.8 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//
2328Function beamDiameter(direction)
2329        String direction
2330
2331        NVAR lens = root:Packages:NIST:SAS:gUsingLenses
2332        if(lens)
2333                return sourceApertureDiam()
2334        endif
2335       
2336    Variable l1 = sourceToSampleDist()
2337    Variable l2 //= sampleAperToDetDist()
2338    Variable d1,d2,bh,bv,bm,umbra,a1,a2
2339   
2340    NVAR L2diff = root:Packages:NIST:SAS:L2diff
2341    NVAR lambda = root:Packages:NIST:SAS:gLambda
2342    NVAR lambda_width = root:Packages:NIST:SAS:gDeltaLambda
2343    NVAR bs_factor = root:Packages:NIST:SAS:bs_factor
2344   
2345    l2 = sampleToDetectorDist() + L2diff
2346    a1 = sourceApertureDiam()
2347    a2 = sampleApertureDiam()
2348   
2349    d1 = a1*l2/l1
2350    d2 = a2*(l1+l2)/l1
2351    bh = d1+d2          //beam size in horizontal direction
2352    umbra = abs(d1-d2)
2353    //vertical spreading due to gravity
2354    bv = bh + 1.25e-8*(l1+l2)*l2*lambda*lambda*lambda_width
2355    bm = (bs_factor*bh > bv) ? bs_factor*bh : bv //use the larger of horiz*safety or vertical
2356   
2357    strswitch(direction)        // string switch
2358        case "vertical":                // execute if case matches expression
2359                return(bv)
2360                break                                           // exit from switch
2361        case "horizontal":              // execute if case matches expression
2362                return(bh)
2363                break
2364        case "maximum":         // execute if case matches expression
2365                return(bm)
2366                break
2367        default:                                                        // optional default expression executed
2368                return(bm)                                              // when no case matches
2369    endswitch
2370End
2371
2372//on NG3 and NG7, allowable sizes are 1,2,3,4" diameter
2373//will return values larger than 4.0*2.54 if a larger beam is needed
2374//
2375// - in an approximate way, account for lenses
2376Function beamstopDiam()
2377
2378        NVAR yesLens = root:Packages:NIST:SAS:gUsingLenses
2379        Variable bm=0
2380        Variable bs=0.0
2381   
2382        if(yesLens)
2383                //bm = sourceApertureDiam()             //ideal result, not needed
2384                bs = 1                                                          //force the diameter to 1"
2385        else
2386                bm = beamDiameter("maximum")
2387                do
2388                bs += 1
2389           while ( (bs*2.54 < bm) || (bs > 30.0))                       //30 = ridiculous limit to avoid inf loop
2390        endif
2391
2392        //update the wave
2393        WAVE rw=root:Packages:NIST:SAS:realsRead
2394        rw[21] = bs*25.4                //store the BS diameter in mm
2395       
2396    return (bs*2.54)            //return diameter in cm, not inches for txt
2397End
2398
2399//returns the projected diameter of the beamstop at the anode plane.
2400// most noticeable at short SDD
2401//if flag == 0 use conservative estimate = largest diameter (for SASCALC, default)
2402//if flag != 0 use point aperture = average diameter (for resolution calculation)
2403Function beamstopDiamProjection(flag)
2404        Variable flag
2405       
2406        NVAR L2diff = root:Packages:NIST:SAS:L2diff
2407        Variable a2 = sampleApertureDiam()
2408        Variable bs = beamstopDiam()
2409        Variable l2, LB, BS_P
2410   
2411        l2 = sampleToDetectorDist() + L2diff
2412        LB = 20.1 + 1.61*BS                     //distance in cm from beamstop to anode plane (empirical)
2413        if(flag==0)
2414                BS_P = bs + (bs+a2)*lb/(l2-lb)          //diameter of shadow from parallax
2415        else
2416                BS_P = bs + bs*lb/(l2-lb)               //diameter of shadow, point A2
2417        endif
2418        return (bs_p)           //return projected diameter in cm
2419End
2420
2421// 19MAR07 - using correction from John for an estimate of the shadow of the beamstop
2422// at the detection plane. This is a noticeable effect at short SDD, where the projected
2423// diameter of the beamstop is much larger than the physical diameter.
2424Function qMin()
2425
2426    Variable l2s = sampleToDetectorDist()       //distance from sample to detector in cm
2427//    Variable bs = beamstopDiam()              //beamstop diameter in cm
2428    Variable bs_p = beamstopDiamProjection(0)           //projected beamstop diameter in cm
2429    NVAR lambda = root:Packages:NIST:SAS:gLambda
2430    NVAR d_det = root:Packages:NIST:SAS:d_det           //cm
2431    NVAR a_pixel = root:Packages:NIST:SAS:a_pixel       //cm
2432
2433    return( (pi/lambda)*(bs_p + d_det + a_pixel)/l2s )          //use bs_p rather than bs
2434   // return( (pi/lambda)*(bs + d_det + a_pixel)/l2s )          //use bs (incorrect)
2435End
2436
2437Function qMaxVert()
2438
2439    Variable theta
2440    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2441    NVAR lambda = root:Packages:NIST:SAS:gLambda
2442        NVAR det_width = root:Packages:NIST:SAS:det_width
2443       
2444    theta = atan( (det_width/2.0)/l2s )
2445   
2446    return ( 4.0*pi/lambda * sin(theta/2.0) )
2447end
2448
2449Function qMaxCorner()
2450
2451    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2452    Variable radial
2453    NVAR lambda = root:Packages:NIST:SAS:gLambda
2454        NVAR det_off = root:Packages:NIST:SAS:gOffset
2455        NVAR det_width = root:Packages:NIST:SAS:det_width
2456
2457    radial=sqrt( (0.5*det_width)*(0.5*det_width)+(0.5*det_width+det_off)*(0.5*det_width+det_off) )
2458   
2459    return ( 4*pi/lambda*sin(0.5*atan(radial/l2s)) )
2460End
2461
2462Function qMaxHoriz()
2463
2464    Variable theta
2465    Variable l2s = sampleToDetectorDist()       //distance from sample to detector
2466    NVAR lambda = root:Packages:NIST:SAS:gLambda
2467        NVAR det_off = root:Packages:NIST:SAS:gOffset
2468        NVAR det_width = root:Packages:NIST:SAS:det_width
2469
2470    theta = atan( ((det_width/2.0)+det_off)/l2s )       //from the instance variables
2471   
2472    return ( 4.0*pi/lambda * sin(theta/2.0) )
2473End
2474
2475// calculate sigma for the resolution function at either limit of q-range
2476Function deltaQ(atQ)
2477        Variable atQ
2478       
2479    Variable k02,lp,l1,l2,sig_02,sigQ2,a1,a2
2480    NVAR l2Diff = root:Packages:NIST:SAS:L2diff
2481        NVAR lambda = root:Packages:NIST:SAS:gLambda
2482        NVAR lambda_width = root:Packages:NIST:SAS:gDeltaLambda
2483        NVAR d_det = root:Packages:NIST:SAS:d_det
2484        NVAR del_r = root:Packages:NIST:SAS:del_r
2485       
2486       
2487    l1 = sourceToSampleDist()
2488    l2 = sampleToDetectorDist() + L2diff
2489    a1 = sourceApertureDiam()
2490    a2 = sampleApertureDiam()
2491   
2492    k02 = (6.2832/lambda)*(6.2832/lambda)
2493    lp = 1/(1/l1 + 1/l2)
2494   
2495    sig_02 = (0.25*a1/l1)*(0.25*a1/l1)
2496    sig_02 += (0.25*a2/lp)*(0.25*a2/lp)
2497    sig_02 += (d_det/(2.355*l2))*(d_det/(2.355*l2))
2498    sig_02 += (del_r/l2)*(del_r/l2)/12
2499    sig_02 *= k02
2500   
2501    sigQ2 = sig_02 + (atQ*lambda_width)*(atQ*lambda_width)/6
2502
2503    return(100*sqrt(sigQ2)/atQ)
2504End
2505
2506// updated with new flux numbers from John Barker
2507// NG3 - Feb 2009
2508// NG7 - July 2009
2509//
2510// guide loss has been changed to 0.95 rather than the old value of 0.95
2511//
2512// other values are changed in the initialization routines
2513//
2514Function beamIntensity()
2515
2516    Variable alpha,f,t,t4,t5,t6,as,solid_angle,l1,d2_phi
2517    Variable a1,a2,retVal
2518    SetDataFolder root:Packages:NIST:SAS
2519    NVAR l_gap=l_gap,guide_width =guide_width,ng = gNg
2520    NVAR lambda_t=lambda_t,b=b,c=c
2521    NVAR lambda=gLambda,t1=t1,t2=t2,t3=t3,phi_0=phi_0
2522    NVAR lambda_width=gDeltaLambda
2523    NVAR guide_loss=gGuide_loss
2524   
2525    l1 = sourceToSampleDist()
2526    a1 = sourceApertureDiam()
2527    a2 = sampleApertureDiam()
2528   
2529   
2530    alpha = (a1+a2)/(2*l1)      //angular divergence of beam
2531    f = l_gap*alpha/(2*guide_width)
2532    t4 = (1-f)*(1-f)
2533    t5 = exp(ng*ln(guide_loss)) // trans losses of guides in pre-sample flight
2534    t6 = 1 - lambda*(b-(ng/8)*(b-c))            //experimental correction factor
2535    t = t1*t2*t3*t4*t5*t6
2536   
2537    as = pi/4*a2*a2             //area of sample in the beam
2538    d2_phi = phi_0/(2*pi)
2539    d2_phi *= exp(4*ln(lambda_t/lambda))
2540    d2_phi *= exp(-1*(lambda_t*lambda_t/lambda/lambda))
2541
2542    solid_angle = pi/4* (a1/l1)*(a1/l1)
2543
2544    retVal = as * d2_phi * lambda_width * solid_angle * t
2545    SetDataFolder root:
2546    return (retVal)
2547end
2548
2549Function figureOfMerit()
2550
2551        Variable bi = beamIntensity()
2552        NVAR lambda = root:Packages:NIST:SAS:gLambda
2553       
2554    return (lambda*lambda*bi)
2555End
2556
2557//estimate the number of pixels in the beam, and enforce the maximum countrate per pixel (idmax)
2558Function attenuatorTransmission()
2559
2560    Variable num_pixels,i_pix           //i_pix = id in John's notation
2561    Variable bDiam = beamDiameter("horizontal") //!! note that prev calculations used bh (horizontal)
2562    Variable atten,a2
2563    SetDataFolder root:Packages:NIST:SAS
2564    NVAR a_pixel=a_pixel,idmax=idmax
2565   
2566    a2 = sampleApertureDiam()
2567   
2568    num_pixels = pi/4*(0.5*(a2+bDiam))*(0.5*(a2+bDiam))/a_pixel/a_pixel
2569    i_pix = ( beamIntensity() )/num_pixels
2570   
2571    atten = (i_pix < idmax) ? 1.0 : idmax/i_pix
2572    SetDataFolder root:
2573    return(atten)
2574End
2575
2576// TODO_10m: is this going to be close enough? Are the steps the same?
2577Function attenuatorNumber()
2578
2579    Variable atten = attenuatorTransmission()
2580    Variable af,nf,numAtten
2581    SetDataFolder root:Packages:NIST:SAS
2582    NVAR lambda=gLambda
2583   
2584    af = 0.498 + 0.0792*lambda - 1.66e-3*lambda*lambda
2585    nf = -ln(atten)/af          //floating point
2586   
2587    numAtten = trunc(nf) + 1                    //in c, (int)nf
2588    //correct for larger step thickness at n > 6
2589    if(numAtten > 6)
2590        numAtten = 7 + trunc( (numAtten-6)/2 )          //in c, numAtten = 7 + (int)( (numAtten-6)/2 )
2591    endif
2592   
2593    SetDatafolder root:
2594    return (numAtten)
2595End
Note: See TracBrowser for help on using the repository browser.