source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_TubeAdjustments.ipf @ 1079

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

cleanup of TODO items in code, no noteworthy changes

prepare a test release package for the January startup of VSANS, not a general release (since no changes to SANS or USANS)

File size: 31.0 KB
Line 
1#pragma rtGlobals=3             // Use modern global access method and strict wave access.
2
3//
4// functions for testing and then actually applying the nonlinear corrections to the
5// tube detectors. These routines are for a test bank of 8 tubes (vertical) that were
6// run at a subdivision of 1024. VSANS will be different in practice
7//
8// but the fundamental process is the same, and can be translated into proper functions as needed
9//
10//
11
12
13//
14//
15// x- need a way to generate the known, physical dimensions of the slots
16// Make/O/D/N=5 peak_spacing_mm_ctr
17// peak_spacing_mm_ctr = {-350,-190,0,190,350} (to be filled in with the correct measurements,
18//   possibly different for each panel)
19//
20// x- a 128 point wave "tube_pixel" (=p) is made in V_ArrayToTubes(), and is needed for the WM
21//   procedures to identify the peak positions.
22//
23// x- fit with either gauss or lor function to get non-integer pixel values for the peak locations
24//
25// x- do I fit each individually to "tweak" the located values, or fit all 5 at once with a
26//    custom fit function and guess some good starting values for peak height, location, etc.
27//
28//
29// x- find a way to display all of the results - in a way that can quickly identify any fits
30//    that may be incorrect
31//
32//
33
34
35
36// the main routines are:
37
38//(1)
39//
40// to get from an array to individual tubes
41// V_ArrayToTubes(detector)
42//
43// (not needed) to get from individual tubes to an array
44//      V_Tubes_to_Array()                     
45
46//(2)
47// then to locate all of the peak positions
48//      V_MakeTableForPeaks(numTube,numPeak)           
49//      V_Identify_AllPeaks()
50//              AutoFindPeaksCustom()           // if Identify_AllPeaks  doesn't work -try this, setting the "noise" to 1 and smoothing to 2
51
52// (3) Refine the fitted peak positions
53//
54
55//(4)
56// fit to find all of the quadratic coefficients
57//      MakeTableForFitCoefs(numTube,numCoef)
58//      PlotFit_AllPeaks()
59
60
61//(5)
62// then pick a display method
63//
64//      Make_mm_tubes()
65//      Append_Adjusted_mm()
66//
67//      MakeMatrix_PixCenter()
68//      FillMatrix_Pix_Center(ind)
69//
70//
71// -or- (note that the pack_image wave that is generated here is for display ONLY)
72// --since it is interpolated
73//
74// Interpolate_mm_tubes()
75
76
77// The function most used externally is:
78// V_TubePix_to_mm(c0,c1,c2,pix)
79//
80// which will return the real space distance (in mm?) for a given pixel
81// and the tube's coefficients. The distance is relative to the zero position on the
82// detector (which was defined when the coefficients were determined)
83
84
85
86//
87// (0) -- what I start with:
88// -- a table of the mm spacing of the slots (20 of them)
89// -- masked data from each of the (8) tubes
90// -- the table of slots may need to be corrected for parallax, depending on the geometry of the test
91// ** In the table of slots, pick a slot near the center, and SET that to ZERO. Then all of the other
92//   distances are relative to that zero point. This is a necessary reference point.
93//
94
95//
96// x- need a routine to set up the actual measurements of the slot positions
97//
98//
99//
100// x- the slot positioning is different for the L/R and T/B detectors
101//
102Proc V_SetupSlotDimensions()
103        Make/O/D/N=5 peak_spacing_mm_ctr_TB,peak_spacing_mm_ctr_LR
104        peak_spacing_mm_ctr_TB = {-159.54,-80.17,0,80.17,159.54}
105        peak_spacing_mm_ctr_LR = {-379.4,-189.7,0,189.7,380.2}
106        DoWindow/F Real_mm_Table
107        if(V_Flag == 0)
108                Edit/N=Real_mm_Table peak_spacing_mm_ctr_TB,peak_spacing_mm_ctr_LR
109        endif
110End
111
112
113
114//
115// (1) -- get the individual tubes into an array
116//
117//
118Proc V_Tubes_to_Array()
119        Make/O/D/N=(8,1127) pack
120        edit pack
121        display;appendimage pack
122        pack[0][] = tube1[q]
123        pack[1][] = tube2[q]
124        pack[2][] = tube3[q]
125        pack[3][] = tube4[q]
126        pack[4][] = tube5[q]
127        pack[5][] = tube6[q]
128        pack[6][] = tube7[q]
129        pack[7][] = tube8[q]
130        ModifyImage pack ctab= {*,*,ColdWarm,0}
131End
132
133// or the other way around
134// - get the array into individual tubes ready for fitting.
135//
136Proc V_ArrayToTubes(detStr)
137        String detStr
138//      Prompt wStr,"Select detector panel",popup,WaveList("data_*",";","")
139        Prompt detStr,"Select detector panel",popup,ksDetectorListAll
140       
141        String/G root:detUsed = detStr
142       
143        Variable ii,numTubes=48
144        String str="tube"
145       
146        Variable dim0,dim1
147       
148        detStr = "root:Packages:NIST:VSANS:RAW:entry:instrument:detector_"+detStr+":data"
149       
150        dim0 = DimSize($detStr,0)
151        dim1 = DimSize($detStr,1)
152
153       
154        Make/O/D/N=128 tube_pixel
155        tube_pixel = p
156       
157       
158        ii=0
159        do
160                Make/O/D/N=128 $(str+num2str(ii))
161               
162                if(dim0 == 128)
163                        $(str+num2str(ii)) = $(detStr)[p][ii]
164                else
165                        $(str+num2str(ii)) = $(detStr)[ii][p]
166                endif
167               
168                ii+=1
169        while(ii < numTubes)
170
171End
172
173
174// (2) -- for each of the tubes, find the x-position (in pixels) of each of the (20) peaks
175// -- automatically loads the Analysis Package "MultiPeakFit 2"
176//
177// automatically find the peaks (after including MultiPeakFit 2)
178//              AutomaticallyFindPeaks()
179//
180//-- or if having difficulty
181//              AutoFindPeaksCustom()           // try this, setting the "noise" to 1 and smoothing to 2
182//
183// -- if really having difficulty, you can do the "full" MultiPeak Fit
184//
185// -- If (hopefully) using the easy way, the results are in:
186// root:WA_PeakCentersY,root:WA_PeakCentersX
187//
188// -- so to see the results:
189//¥Edit/K=0  root:WA_PeakCentersY,root:WA_PeakCentersX
190//
191// -- then sort the results - they seem to be in no real order...
192//¥Sort WA_PeakCentersX WA_PeakCentersY,WA_PeakCentersX
193//
194Proc V_MakeTableForPeaks(numTube,numPeak)
195        Variable numTube=48,numPeak=5
196       
197        Make/O/D/N=(numPeak,numTube) PeakTableX,PeakTableY              //*2 to store x-location and peak height (y)
198       
199        DoWindow/F Peak_Pixel_Loc
200        if(V_flag == 0)
201                Edit/N=Peak_Pixel_Loc peakTableX
202        endif
203       
204        Execute/P "INSERTINCLUDE <Multi-peak fitting 2.0>"
205        DoWindow/K MultiPeak2StarterPanel
206
207//      DoAlert 0, "Load the Package: Analysis->MultiPeak Fitting->MultiPeak Fitting 2"
208End
209
210Proc V_Identify_AllPeaks()
211
212        Variable ii,numTubes=48
213        String str="tube"
214       
215        ii=0
216        do
217                V_Identify_Peaks(str+num2str(ii),ii)
218                ii+=1
219        while(ii < numTubes)
220
221End
222
223Proc V_Identify_Peaks(tubeStr,ind)
224        String tubeStr
225        Variable ind
226       
227        // must use a wave of pixels rather than "calculated" -- if calculated is used it only
228        // returns integer values for the peak locations
229       
230//      AutomaticallyFindPeaks() //-- this is a function that doesn't take any parameters - so
231// I need to pull this from the WM function to call the worker directly
232        Variable pBegin=0, pEnd= numpnts($(tubeStr))-1
233        Variable/C estimates= EstPeakNoiseAndSmfact($(tubeStr),pBegin, pEnd)
234        Variable noiselevel=real(estimates)
235        Variable smoothingFactor=imag(estimates)
236        Variable maxPeaks = 20
237        Variable minPeakPercent = 10
238       
239        AutoFindPeaksWorker($(tubeStr), $("tube_pixel"), pBegin, pEnd, maxPeaks, minPeakPercent, noiseLevel, smoothingFactor)
240// end WM function call
241
242        Sort WA_PeakCentersX WA_PeakCentersY,WA_PeakCentersX
243       
244        peakTableX[][ind] = WA_PeakCentersX[p]          // the peak position
245        peakTableY[][ind] = WA_PeakCentersY[p]          // the peak height
246       
247End
248
249
250
251
252// ADD
253// a step to refine the peak positioning - currently an integer value
254//  fit with a gauss or lorentzian
255
256// CurveFit/M=2/W=0/TBOX=(0x310) lor, tube47[29,53]/X=tube_pixel[29,53]/D
257
258//CurveFit/M=2/W=0 lor, tube47[29,53]/X=tube_pixel[29,53]/D
259//fit_tube47= W_coef[0]+W_coef[1]/((x-W_coef[2])^2+W_coef[3])
260//W_coef={-20.37,876.94,40.078,0.5201}
261//W_sigma={6.52,47.3,0.0241,0.0308}
262
263Proc V_MakeTableForRefinedFit(numTube,numPeak)
264        Variable numTube=48,numPeak=5
265       
266        Make/O/D/N=(numPeak,numTube) position_refined,position_refined_err              //
267       
268        DoWindow/F Refined_Positions
269        if(V_flag == 0)
270                Edit/N=Refined_Positions position_refined
271        endif
272       
273End
274
275Proc V_Refine_All_PeakPos()
276
277        Variable ii,numTubes=48
278       
279        ii=0
280        do
281                V_Refine_PeakPos(ii)
282                ii+=1
283        while(ii<numTubes)
284
285End
286
287
288//CurveFit/M=2/W=0 lor, tube47[29,53]/X=tube_pixel[29,53]/D
289//fit_tube47= W_coef[0]+W_coef[1]/((x-W_coef[2])^2+W_coef[3])
290
291Proc V_Refine_PeakPos(ind)
292        Variable ind
293       
294//
295// x- hard-wired for 5 peaks
296
297        Variable ii,lo,hi
298       
299       
300        ii=0
301        do
302       
303                if(ii==0)
304                // 1st peak
305                // define fitting range pixels (integer)
306                        lo = 0
307                else
308                        lo = trunc(0.5*(peakTableX[ii-1][ind] + peakTableX[ii][ind]))
309                endif
310               
311                if(ii==4)
312                        hi = numpnts(tube_pixel)-1
313                else
314                        hi = trunc(0.5*(peakTableX[ii][ind] + peakTableX[ii+1][ind]))
315                endif
316               
317                // do I need initial guesses?
318                CurveFit/M=0/W=2 lor, $("tube"+num2str(ind))[lo,hi]/X=tube_pixel[lo,hi]/D
319               
320                position_refined[ii][ind] = W_coef[2]
321                position_refined_err[ii][ind] = W_sigma[2]
322
323                ii += 1
324
325        while(ii < 5)
326       
327End
328
329
330
331
332
333// -- save a copy of the root:WA_PeakCentersY,root:WA_PeakCentersX values
334//    for later in case the fitting failed, then you can go back and re-do
335//
336// -- then plot:
337//
338//      Display peak_spacing_mm_ctr vs WA_PeakCentersX
339//
340//      Then do a "QuickFit" of the peak position to the data using a polynomial of order 3 (= quadratic)
341//
342// result is in W_coef, W_sigma
343//
344// -- an example of the "quickFit" command is below, so it can be programmed rather than the menu every time
345//¥Display peak_spacing_mm_ctr vs WA_PeakCentersX3
346//¥CurveFit/M=2/W=0/TBOX=(0x310) poly 3, peak_spacing_mm_ctr/X=WA_PeakCentersX3/D
347//  fit_peak_spacing_mm_ctr= poly(W_coef,x)
348//  W_coef={-571.42,1.1135,-4.2444e-05}
349//  V_chisq= 8.5841;V_npnts= 20;V_numNaNs= 0;V_numINFs= 0;
350//  V_startRow= 0;V_endRow= 19;
351//  W_sigma={0.595,0.00246,2.15e-06}
352//  Coefficient values ± one standard deviation
353//      K0      =-571.42 ± 0.595
354//      K1      =1.1135 ± 0.00246
355//      K2      =-4.2444e-05 ± 2.15e-06
356//
357//
358//
359// for (8) tubes, keep all of the fit coefficients
360//
361//¥make/O/D/N=(3,8) fit_coef
362//¥edit fit_coef
363//¥make/O/D/N=(3,8) fit_sigma
364//¥edit fit_sigma
365//
366// -- copy and paste in the W_coef and W_sigma values (or by a command)
367//
368
369
370Proc V_MakeTableForFitCoefs(numTube,numCoef)
371        Variable numTube=48,numCoef=3
372       
373        Make/O/D/N=(numTube,numCoef) TubeCoefTable,TubeSigmaTable               //
374       
375        DoWindow/F Quad_Coefficients
376        if(V_flag == 0)
377                Edit/N=Quad_Coefficients TubeCoefTable
378        endif
379
380        String detUsed = root:detUsed
381       
382        if(strsearch(detUsed,"L",0) >= 0 || strsearch(detUsed,"R",0) >= 0)
383                Duplicate/O     peak_spacing_mm_ctr_LR, peak_spacing_mm_ctr
384                DoAlert 0,"Using peak spacing for L/R"
385        else
386                Duplicate/O     peak_spacing_mm_ctr_TB, peak_spacing_mm_ctr
387                DoAlert 0,"Using peak spacing for T/B"
388        endif
389End
390
391Proc V_PlotFit_AllPeakPosition()
392
393        Variable ii,numTubes=48
394       
395        ii=0
396        do
397                V_PlotFit_PeakPosition(ii)
398                ii+=1
399        while(ii<numTubes)
400
401End
402
403Proc V_PlotFit_PeakPosition(ind)
404        Variable ind
405       
406        Duplicate/O WA_PeakCentersX, tmpX
407       
408//      tmpX = peakTableX[p][ind]
409        tmpX = position_refined[p][ind]
410//      Display peak_spacing_mm_ctr vs tmpX
411       
412//      CurveFit/M=2/W=0/TBOX=(0x310) poly 3, peak_spacing_mm_ctr/X=tmpX/D
413        CurveFit/M=0/W=2 poly 3, peak_spacing_mm_ctr/X=tmpX/D
414       
415        TubeCoefTable[ind][] = W_coef[q]
416        TubeSigmaTable[ind][] = W_sigma[q]
417       
418End
419
420
421
422
423//¥Duplicate tube1 tube1_mm
424//¥tube1_mm = V_TubePix_to_mm(fit_coef[0][0],fit_coef[1][0],fit_coef[2][0],p)
425
426
427////////
428// then there are various display options:
429
430// adjust the center (pixel) of the tube:
431// - measCtr is the pixel location of the DEFINED "zero" peak
432// nominal Ctr is the pixel number of this DEFINED zer0 position, nominally nPix/2-1
433//
434// ( be sure to pick better names, and use a loop...)
435//      adj_tube = raw_tube[p+(measCtr-nominalCtr)]
436//
437// then fill and display a new matrix. The center will be reasonably well aligned, and will
438// get worse towards the ends of the tubes
439// (this may be the "preferred" way of displaying the raw data if the centers are far off)
440// -- this also may be what I need to start with to fit the data to locate the beam center.
441
442
443// I can also display the fully corrected tubes, where the y-axis is now in real-space mm
444// rather than arbitrary pixels. The x-axis is still tube nubmer.
445// -- do this with the procedure"
446//   Append_Adjusted_mm()               // name may change...
447//
448// -- the point of the appending is that it allows each tube to be plotted on an image plot
449// with its own y-axis. Every one will be different and will be non-linear. These conditions
450// BOTH prevent using any of the "normal" image plotting or manipulation routines.
451// - the gist is below:
452//
453//      duplicate tube1_mm adjusted_mm_edge
454//      InsertPoints 0,1, adjusted_mm_edge
455//      // be sure to use the correct set of coefficients
456//      adjusted_mm_edge[0] = V_TubePix_to_mm(fit_coef[0][0],fit_coef[1][0],fit_coef[2][0],-1)
457//     
458//      Display;AppendImage adjusted_pack vs {*,adjusted_mm_edge}
459
460
461Proc V_Make_mm_tubes()
462
463        Variable ii,numTubes=8
464       
465        ii=1
466        do
467                Duplicate $("tube"+num2str(ii)) $("tube"+num2str(ii)+"_mm")
468                $("tube"+num2str(ii)+"_mm") = V_TubePix_to_mm(TubeCoefTable[ii-1][0],TubeCoefTable[ii-1][1],TubeCoefTable[ii-1][2],p)
469                ii+=1
470        while(ii<=numTubes)
471       
472End
473
474
475Proc V_Append_Adjusted_mm()
476
477// a blank base image
478        Duplicate/O pack junk
479        junk=0
480        SetScale/I y -600,600,"", junk          // -600,600 is tooo large and not general
481        Display;Appendimage junk
482
483        Variable ii
484       
485        ii=1
486        do
487                make/O/D/N=(1,1127) $("tube"+num2str(ii)+"_mm_mat")=0   
488       
489                $("tube"+num2str(ii)+"_mm_mat")[0][] = $("tube"+num2str(ii))
490                SetScale/I x (ii-1),ii,"", $("tube"+num2str(ii)+"_mm_mat")              //builds up the x-axis
491               
492                duplicate/O $("tube"+num2str(ii)+"_mm") $("edge"+num2str(ii)+"_mm")
493                InsertPoints 0,1, $("edge"+num2str(ii)+"_mm")           //needs to be one point longer
494        // be sure to use the correct set of coefficients
495                $("edge"+num2str(ii)+"_mm")[0] = V_TubePix_to_mm(TubeCoefTable[0][0],TubeCoefTable[0][1],TubeCoefTable[0][2],-1)
496       
497                AppendImage $("tube"+num2str(ii)+"_mm_mat") vs {*,$("edge"+num2str(ii)+"_mm")}
498                ModifyImage $("tube"+num2str(ii)+"_mm_mat") ctab= {*,*,ColdWarm,0}
499       
500                ii+=1
501        while(ii < 9)
502       
503end
504
505
506
507////////////////////////////////
508
509Proc V_MakeMatrix_PixCenter()
510        Duplicate/O pack pack_centered
511       
512        Variable ii,numTubes=8
513       
514        ii=1
515        do
516                V_FillMatrix_Pix_Center(ii)
517                ii+=1
518        while(ii<=numTubes)
519       
520        Display;AppendImage pack_centered
521        ModifyImage pack_centered ctab= {*,*,ColdWarm,0}
522
523end
524
525//
526// this fills a matrix with the tubes center aligned ONLY, with the y-axis in pixels
527//
528// adj_tube = raw_tube[p+(measCtr-nominalCtr)]
529// finds the center automatically from the tubeN_mm wave, where it crosses zero
530//
531// Tube #1 is the "base" ans others are shifted to match that one's "zero"
532//
533// FindRoots/P=W_coef           can also be used to find the roots (but which one?)
534//
535Function V_FillMatrix_Pix_Center(ind)
536        Variable ind
537       
538        Wave adj=root:pack_centered
539        Wave tube1_mm = $("root:tube1_mm")
540        Wave tube = $("root:tube"+num2str(ind))
541        wave tube_mm = $("root:tube"+num2str(ind)+"_mm")
542
543        Variable base,shift,ii,num,pt
544       
545        num=numpnts(tube)
546       
547        FindLevel tube1_mm 0
548        base = round(V_LevelX)
549       
550       
551        FindLevel tube_mm 0
552        shift = round(V_LevelX)
553       
554        for(ii=0;ii<num;ii+=1)
555                pt = ii + (shift-base)
556                if(pt >= 0 && pt < num)
557                        adj[ind-1][ii] = tube[pt]
558                endif
559        endfor
560       
561        return(0)
562End
563
564
565// this fills a matrix with the tubes center aligned ONLY, with the y-axis in mm
566// -- there seems to be little reason to do this --
567// -- either keep pixels and align centers
568// -- OR -- use mm and append each tube with its own y-axis
569//
570Function V_FillAdjusted(ind)
571        Variable ind
572       
573        Wave adj=root:adjusted_pack
574        Wave tube1_mm
575        Wave tube = $("root:tube"+num2str(ind))
576        wave tube_mm = $("root:tube"+num2str(ind)+"_mm")
577
578        Variable base,shift,ii,num,pt
579       
580        num=numpnts(tube1_mm)
581       
582        FindLevel tube1_mm 0
583        base = round(V_LevelX)
584       
585       
586        FindLevel tube_mm 0
587        shift = round(V_LevelX)
588       
589        for(ii=0;ii<num;ii+=1)
590                pt = ii + (shift-base)
591                if(pt >= 0 && pt < num)
592                        adj[ind-1][ii] = tube[pt]
593                endif
594        endfor
595       
596        return(0)
597End
598
599
600// c0,c1,c2,pix
601// c0-c2 are the fit coefficients
602// pix is the test pixel
603//
604// returns the distance in mm (relative to ctr pixel)
605// ctr is the center pixel, as defined when fitting to quadratic
606//
607Function V_TubePix_to_mm(c0,c1,c2,pix)
608        Variable c0,c1,c2,pix
609       
610        Variable dist
611        dist = c0 + c1*pix + c2*pix*pix
612       
613        return(dist)
614End
615
616////
617
618
619
620// set the (linear) range of the y-axis of the matrix to be the
621// range of the 1st tube. This is completely arbitrary
622//
623Proc V_Interpolate_mm_tubes()
624
625        Duplicate/O pack pack_image
626
627        Variable ii,numTubes=8
628        Variable p1,p2
629        p1 = tube1_mm[0]
630        p2 = tube1_mm[numpnts(tube1_mm)-1]
631       
632        SetScale/I y p1,p2,"", pack_image
633       
634        // then make a temporary 1D wave to help with the interpolation
635        Duplicate/O tube1_mm lin_mm,lin_val
636        SetScale/I x p1,p2,"", lin_mm
637        lin_mm = x                      //fill with the linear mm spacing
638        lin_val=0
639       
640        ii=1
641        do
642                lin_val = interp(lin_mm, $("tube"+num2str(ii)+"_mm"), $("tube"+num2str(ii)))
643                pack_image[ii-1][] = lin_val[q]
644               
645                ii+=1
646        while(ii<=numTubes)
647       
648        display;appendimage pack_image
649        ModifyImage pack_image ctab= {*,*,ColdWarm,0}
650       
651End
652
653
654
655
656
657
658
659
660
661
662
663
664
665// this doesn't work - the interpolation step doesn't do what I want
666// and the plot of the triplet with f(z) for color doesn't fill space like I want
667Proc V_AnotherExample()
668
669        Concatenate/O/NP {tube1_mm,tube2_mm,tube3_mm,tube4_mm,tube5_mm,tube6_mm,tube7_mm,tube8_mm},cat_mm
670        Concatenate/O/NP {tube1,tube2,tube3,tube4,tube5,tube6,tube7,tube8},cat_tubes
671        Duplicate/O cat_mm,cat_num
672        Variable num=1127
673        cat_num[0,num-1]=1
674        cat_num[num,2*num-1]=2
675        cat_num[2*num,3*num-1]=3
676        cat_num[3*num,4*num-1]=4
677        cat_num[4*num,5*num-1]=5
678        cat_num[5*num,6*num-1]=6
679        cat_num[6*num,7*num-1]=7
680        cat_num[7*num,8*num-1]=8
681
682        Display cat_mm vs cat_num
683        ModifyGraph mode=3,marker=9
684        ModifyGraph zColor(cat_mm)={cat_tubes,*,*,ColdWarm,0}
685
686        Concatenate/O {cat_num,cat_mm,cat_tubes}, tripletWave
687        ImageInterpolate Kriging tripletWave
688        AppendImage M_InterpolatedImage
689
690//      Make/O/N=20 xWave=enoise(4),yWave=enoise(5),zWave=enoise(6)  // Random points
691//      Display yWave vs xWave
692//      ModifyGraph mode=3,marker=19
693//      ModifyGraph zColor(yWave)={zWave,*,*,Rainbow,0}
694//
695//      Concatenate/O {xWave,yWave,zWave}, tripletWave
696//      ImageInterpolate/S={-5,0.1,5,-5,0.1,5} voronoi tripletWave
697//      AppendImage M_InterpolatedImage
698
699end
700
701// this desn't work either...
702// (same y-axis for the entire image, which is not the case for the tubes)
703//
704// from the WM help file:
705// Plotting a 2D Z Wave With 1D X and Y Center Data
706//
707Function V_MakeEdgesWave(centers, edgesWave)
708        Wave centers    // Input
709        Wave edgesWave  // Receives output
710       
711        Variable N=numpnts(centers)
712        Redimension/N=(N+1) edgesWave
713
714        edgesWave[0]=centers[0]-0.5*(centers[1]-centers[0])
715        edgesWave[N]=centers[N-1]+0.5*(centers[N-1]-centers[N-2])
716        edgesWave[1,N-1]=centers[p]-0.5*(centers[p]-centers[p-1])
717End
718
719//This function demonstrates the use of MakeEdgesWave:
720Function V_DemoPlotXYZAsImage()
721        Make/O mat={{0,1,2},{2,3,4},{3,4,5}}    // Matrix containing Z values
722        Make/O centersX = {1, 2.5, 5}           // X centers wave
723        Make/O centersY = {300, 400, 600}               // Y centers wave
724        Make/O edgesX; V_MakeEdgesWave(centersX, edgesX)        // Create X edges wave
725        Make/O edgesY; V_MakeEdgesWave(centersY, edgesY)        // Create Y edges wave
726        Display; AppendImage mat vs {edgesX,edgesY}
727End
728
729
730
731////////////////////////////
732//
733// Main entry - open the panel and go through
734// each of the steps for each of the detector panels
735//
736Proc V_TubeCoefPanel() : Panel
737        PauseUpdate; Silent 1           // building window...
738        NewPanel /W=(973,45,1156,535)/K=1
739        DoWindow/C V_TubeCoefPanel
740//      ShowTools/A
741
742        SetDrawLayer UserBack
743        SetDrawEnv fsize= 14,fstyle= 1
744        DrawText 5,58,"(1)"
745        SetDrawEnv fsize= 14,fstyle= 1
746        DrawText 5,108,"(2)"
747        SetDrawEnv fsize= 14,fstyle= 1
748        DrawText 5,158,"(3)"
749        SetDrawEnv fsize= 14,fstyle= 1
750        DrawText 5,208,"(4)"
751        SetDrawEnv fsize= 14,fstyle= 1
752        DrawText 5,258,"(5)"
753        SetDrawEnv fsize= 14,fstyle= 1
754        DrawText 5,308,"(6)"
755        SetDrawEnv fsize= 14,fstyle= 1
756        DrawText 5,358,"(7)"
757        SetDrawEnv fsize= 14,fstyle= 1
758        DrawText 5,408,"(8)"
759        SetDrawEnv fsize= 14,fstyle= 1
760        DrawText 5,458,"(9)"
761                       
762        Button button_0,pos={30.00,40.00},size={120.00,20.00},proc=V_Setup_MasksButton,title="Setup"
763        Button button_1,pos={30.00,90.00},size={120.00,20.00},proc=V_ArrayToTubesButton,title="Array to Tubes"
764        Button button_2,pos={30.00,140.00},size={120.00,20.00},proc=V_TableForPeaksButton,title="Table for Peaks"
765        Button button_3,pos={30.00,190.00},size={120.00,20.00},proc=V_IdentifyPeaksButton,title="Identify Peaks"
766        Button button_4,pos={30.00,240.00},size={120.00,20.00},proc=V_RefineTableButton,title="Refine Peak Table"
767        Button button_5,pos={30.00,290.00},size={120.00,20.00},proc=V_RefinePeaksButton,title="Refine Peaks"
768
769        Button button_6,pos={30.00,340.00},size={120.00,20.00},proc=V_QuadFitTableButton,title="Table for Quad"
770        Button button_7,pos={30.00,390.00},size={120.00,20.00},proc=V_QuadFitButton,title="Fit to Quad"
771        Button button_8,pos={30.00,440},size={120.00,20.00},proc=V_PeakPlotButton,title="Plot Peaks"
772       
773EndMacro
774
775
776// a simple display of the refined results
777//
778Function V_PeakPlotButton(ba) : ButtonControl
779        STRUCT WMButtonAction &ba
780
781        switch( ba.eventCode )
782                case 2: // mouse up
783                        // click code here
784                        Execute "V_OpenPeakResultsGraph()"
785                        break
786                case -1: // control being killed
787                        break
788        endswitch
789
790        return 0
791End
792
793// generate the waves and the table
794//
795Function V_TableForPeaksButton(ba) : ButtonControl
796        STRUCT WMButtonAction &ba
797
798        switch( ba.eventCode )
799                case 2: // mouse up
800                        // click code here
801                        Execute "V_MakeTableForPeaks()"
802                        break
803                case -1: // control being killed
804                        break
805        endswitch
806
807        return 0
808End
809
810// use the WM procedures to quickly identify the peak position (and height)
811// to be used in the refining fits
812//
813Function V_IdentifyPeaksButton(ba) : ButtonControl
814        STRUCT WMButtonAction &ba
815
816        switch( ba.eventCode )
817                case 2: // mouse up
818                        // click code here
819                        Execute "V_Identify_AllPeaks()"
820                        break
821                case -1: // control being killed
822                        break
823        endswitch
824
825        return 0
826End
827
828// generate the waves and the table
829//
830Function V_RefineTableButton(ba) : ButtonControl
831        STRUCT WMButtonAction &ba
832
833        switch( ba.eventCode )
834                case 2: // mouse up
835                        // click code here
836                        Execute "V_MakeTableForRefinedFit()"
837                        break
838                case -1: // control being killed
839                        break
840        endswitch
841
842        return 0
843End
844
845// using the initial peak locations from WM, refine the values
846// by fitting each individual peak
847//
848Function V_RefinePeaksButton(ba) : ButtonControl
849        STRUCT WMButtonAction &ba
850
851        switch( ba.eventCode )
852                case 2: // mouse up
853                        // click code here
854                        Execute "V_Refine_All_PeakPos()"
855                        break
856                case -1: // control being killed
857                        break
858        endswitch
859
860        return 0
861End
862
863// finally, with the peak positions, make waves and a table for the
864// quadratic coefficients
865//
866Function V_QuadFitTableButton(ba) : ButtonControl
867        STRUCT WMButtonAction &ba
868
869        switch( ba.eventCode )
870                case 2: // mouse up
871                        // click code here
872                        Execute "V_MakeTableForFitCoefs()"
873                        break
874                case -1: // control being killed
875                        break
876        endswitch
877
878        return 0
879End
880
881// fit all of the peak positions per tube vs. the actual mm locations to
882// obtain the quadratic coefficients
883//
884Function V_QuadFitButton(ba) : ButtonControl
885        STRUCT WMButtonAction &ba
886
887        switch( ba.eventCode )
888                case 2: // mouse up
889                        // click code here
890                        Execute "V_PlotFit_AllPeakPosition()"
891                        break
892                case -1: // control being killed
893                        break
894        endswitch
895
896        return 0
897End
898
899
900// fill the waves and table with the hard-wired slot positions (mm)
901//
902Function V_Setup_MasksButton(ba) : ButtonControl
903        STRUCT WMButtonAction &ba
904
905        switch( ba.eventCode )
906                case 2: // mouse up
907                        // click code here
908                        Execute "V_SetupSlotDimensions()"
909                        break
910                case -1: // control being killed
911                        break
912        endswitch
913
914        return 0
915End
916
917
918// convert the named detector array to 48 individual tube waves
919//
920Function V_ArrayToTubesButton(ba) : ButtonControl
921        STRUCT WMButtonAction &ba
922
923        switch( ba.eventCode )
924                case 2: // mouse up
925                        // click code here
926                        Execute "V_ArrayToTubes()"
927                        break
928                case -1: // control being killed
929                        break
930        endswitch
931
932        return 0
933End
934
935
936///////////////////////////
937//
938// unused - a simple line graph for each tube is much simpler
939//
940Window Gizmo_refinedPositions() : GizmoPlot
941        PauseUpdate; Silent 1           // building window...
942        // Building Gizmo 7 window...
943        NewGizmo/W=(232,448,747,908)
944        ModifyGizmo startRecMacro=700
945        ModifyGizmo scalingOption=63
946        AppendToGizmo Surface=root:position_refined,name=surface0
947        ModifyGizmo ModifyObject=surface0,objectType=surface,property={ srcMode,0}
948        ModifyGizmo ModifyObject=surface0,objectType=surface,property={ surfaceCTab,Rainbow}
949        AppendToGizmo Axes=boxAxes,name=axes0
950        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={0,axisRange,-1,-1,-1,1,-1,-1}
951        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={1,axisRange,-1,-1,-1,-1,1,-1}
952        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={2,axisRange,-1,-1,-1,-1,-1,1}
953        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={3,axisRange,-1,1,-1,-1,1,1}
954        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={4,axisRange,1,1,-1,1,1,1}
955        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={5,axisRange,1,-1,-1,1,-1,1}
956        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={6,axisRange,-1,-1,1,-1,1,1}
957        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={7,axisRange,1,-1,1,1,1,1}
958        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={8,axisRange,1,-1,-1,1,1,-1}
959        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={9,axisRange,-1,1,-1,1,1,-1}
960        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={10,axisRange,-1,1,1,1,1,1}
961        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={11,axisRange,-1,-1,1,1,-1,1}
962        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={-1,axisScalingMode,1}
963        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={-1,axisColor,0,0,0,1}
964        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={0,ticks,2}
965        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={1,ticks,2}
966        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={2,ticks,2}
967        ModifyGizmo modifyObject=axes0,objectType=Axes,property={-1,Clipped,0}
968        ModifyGizmo setDisplayList=0, object=surface0
969        ModifyGizmo setDisplayList=1, object=axes0
970        ModifyGizmo autoscaling=1
971        ModifyGizmo currentGroupObject=""
972        ModifyGizmo showInfo
973        ModifyGizmo infoWindow={651,303,1468,602}
974        ModifyGizmo endRecMacro
975        ModifyGizmo SETQUATERNION={0.573113,-0.115160,-0.275160,0.763255}
976EndMacro
977
978///////////////////////////
979//
980// unused - a simple line graph for each tube is much simpler
981//
982Window Gizmo_DetPanel() : GizmoPlot
983        PauseUpdate; Silent 1           // building window...
984        // Building Gizmo 7 window...
985        NewGizmo/W=(96,290,611,750)
986        ModifyGizmo startRecMacro=700
987        ModifyGizmo scalingOption=63
988        AppendToGizmo Surface=root:slices_L,name=surface0
989        ModifyGizmo ModifyObject=surface0,objectType=surface,property={ srcMode,0}
990        ModifyGizmo ModifyObject=surface0,objectType=surface,property={ surfaceCTab,ColdWarm}
991        AppendToGizmo Axes=boxAxes,name=axes0
992        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={0,axisRange,-1,-1,-1,1,-1,-1}
993        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={1,axisRange,-1,-1,-1,-1,1,-1}
994        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={2,axisRange,-1,-1,-1,-1,-1,1}
995        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={3,axisRange,-1,1,-1,-1,1,1}
996        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={4,axisRange,1,1,-1,1,1,1}
997        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={5,axisRange,1,-1,-1,1,-1,1}
998        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={6,axisRange,-1,-1,1,-1,1,1}
999        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={7,axisRange,1,-1,1,1,1,1}
1000        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={8,axisRange,1,-1,-1,1,1,-1}
1001        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={9,axisRange,-1,1,-1,1,1,-1}
1002        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={10,axisRange,-1,1,1,1,1,1}
1003        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={11,axisRange,-1,-1,1,1,-1,1}
1004        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={-1,axisScalingMode,1}
1005        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={-1,axisColor,0,0,0,1}
1006        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={0,ticks,3}
1007        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={1,ticks,3}
1008        ModifyGizmo ModifyObject=axes0,objectType=Axes,property={2,ticks,3}
1009        ModifyGizmo modifyObject=axes0,objectType=Axes,property={-1,Clipped,0}
1010        AppendToGizmo Surface=root:position_refined,name=surface1
1011        ModifyGizmo ModifyObject=surface1,objectType=surface,property={ fillMode,4}
1012        ModifyGizmo ModifyObject=surface1,objectType=surface,property={ srcMode,0}
1013        ModifyGizmo ModifyObject=surface1,objectType=surface,property={ surfaceCTab,Rainbow}
1014        ModifyGizmo setDisplayList=0, object=axes0
1015        ModifyGizmo setDisplayList=1, object=surface0
1016        ModifyGizmo autoscaling=1
1017        ModifyGizmo currentGroupObject=""
1018        ModifyGizmo showInfo
1019        ModifyGizmo infoWindow={550,23,1367,322}
1020        ModifyGizmo endRecMacro
1021        ModifyGizmo SETQUATERNION={0.499484,-0.278571,-0.448869,0.686609}
1022EndMacro
1023
1024
1025////////////////////////////////////
1026//
1027// An easy way to see the fit results to check if the peak locations all make sense.
1028//
1029Proc V_OpenPeakResultsGraph()
1030
1031        DoWindow/F V_PeakResultsGraph
1032        if(V_flag == 0)
1033                Make/O/D/N=5 tmpPeak,dummyLevel
1034                Make/O/D/N=128 tmpTube
1035               
1036                tmpPeak = position_refined[p][0]
1037                dummyLevel = WaveMax(tube0)
1038                tmpTube = tube0
1039               
1040                V_PeakResultsGraph()
1041        endif
1042
1043End
1044
1045
1046///////////////
1047//
1048// An easy way to see the fit results to check if the peak locations all make sense.
1049//
1050Window V_PeakResultsGraph() : Graph
1051        PauseUpdate; Silent 1           // building window...
1052        Display /W=(750,45,1161,376)/K=1 tmpTube vs tube_pixel
1053       
1054        ControlBar 50
1055       
1056       
1057        AppendToGraph dummyLevel vs tmpPeak
1058        ModifyGraph mode(dummyLevel)=3
1059        ModifyGraph marker(dummyLevel)=19
1060        ModifyGraph rgb(dummyLevel)=(1,16019,65535)
1061       
1062        SetVariable setvar0,pos={10.00,10.00},size={120.00,14.00},proc=V_TubePeakSetVarProc,title="Tube"
1063        SetVariable setvar0,limits={0,47,1},value= _NUM:0
1064       
1065        Label left "Counts"
1066        Label bottom "Pixel Number"
1067EndMacro
1068
1069//
1070// cycle through the tubes
1071//
1072Function V_TubePeakSetVarProc(sva) : SetVariableControl
1073        STRUCT WMSetVariableAction &sva
1074
1075        switch( sva.eventCode )
1076                case 1: // mouse up
1077                case 2: // Enter key
1078                case 3: // Live update
1079                        Variable dval = sva.dval
1080                        String sval = sva.sval
1081                       
1082                        Wave tmpPeak = tmpPeak
1083                        Wave dummyLevel = dummyLevel
1084                        Wave tmpTube = tmpTube
1085                       
1086                        Wave pos_ref = position_refined
1087                        Wave tube = $("tube"+num2str(dval))
1088                       
1089                        tmpPeak = pos_ref[p][dval]
1090                        dummyLevel = WaveMax(tube)
1091                        tmpTube = tube
1092               
1093                        break
1094                case -1: // control being killed
1095                        break
1096        endswitch
1097
1098        return 0
1099End
1100
1101////////////////////////////////////
1102//
1103// DONE
1104// x- document the "simple" save of the detector panels for import and subsequent fitting.
1105//   Documentation is done in the "main" VSANS documentation, and is largely not needed, since Phil
1106//   is doing the nonlinear calibration calculations, not me.
1107//
1108//
1109// takes the data from RAW, by default. This is OK, since even though whatever is in the calibration data
1110// of the file is used when loading into RAW, it is only used for the calculation of q. The actual data
1111// array is unchanged. Alternatively, the data could be pulled from the RawVSANS folder after a
1112// file catalog operation
1113//
1114Proc V_CopyDetectorsToRoot()
1115
1116        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_B:data data_B
1117
1118        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_ML:data data_ML
1119        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_MR:data data_MR
1120        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_MT:data data_MT
1121        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_MB:data data_MB
1122
1123        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_FL:data data_FL
1124        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_FR:data data_FR
1125        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_FT:data data_FT
1126        Duplicate/O root:Packages:NIST:VSANS:RAW:entry:instrument:detector_FB:data data_FB
1127       
1128End
1129
1130//
1131//
1132Proc V_SaveDetectorsITX()
1133// binary save makes each wave an individual file. Igor text groups them all into one file.
1134//      Save/C data_B,data_FB,data_FL,data_FR,data_FT,data_MB,data_ML,data_MR,data_MT
1135        Save/T/M="\r\n" data_B,data_FB,data_FL,data_FR,data_FT,data_MB,data_ML,data_MR,data_MT as "data_B++.itx"
1136
1137End
Note: See TracBrowser for help on using the repository browser.