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

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

Changes to utility functions to properly list the raw data from the 10m SANS (SA4)

Added a utility function to patch BS xPos in 10m data files so that they will be recognized as transmission files. ICE is currently not writing out the BS xPos to the file. This is done in batch mode, similar to what was done for the 8m SANS.

Some changes to the simulation code to permit easier scripting.

Additions to SASCALC for the 10m SANS. Dimensions from John. Still need attenuator table and flux.

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