source: sans/Dev/trunk/NCNR_User_Procedures/Common/Packages/PlotManager/PlotUtils2D_v40.ipf @ 509

Last change on this file since 509 was 509, checked in by srkline, 14 years ago

(1) fix in ProtocolAsPanel? that correctly calculates Kappa for data taken with the Cerca detectors. In the Kappa calculation, the data is manually scaled, bypassing the normal "ADD" step where the 4x detector data is correctly renormalized. Added a flag to divide by 4 if ILL detector type found.

(2) lots of bits to try to accomodate 2D resolution smearing. Changes to the QxQy? output to have everything needed: Qz, sigX, sigY, and shadowing. So presumably, all of the information is in the reduced data, where it should be.

  • added a 2D resolution calculation based on David Mildner's notes.
  • added 2D quadrature calculation for smearing, and new structure definitions

File size: 26.4 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma version=4.00
3#pragma IgorVersion=6.0
4
5//
6// WARNING: there are a lot of assumptions that went into this code at the first pass
7// that will need to be made generic before a release to the wild
8//
9//
10// - some manipulations assume a square matrix
11// - naming schemes are assumed (to be documented)
12//
13//
14
15
16//
17// changed May 2009 to read in the resolution information (now 7 columns)
18// -- subject to change --
19//
20Proc LoadQxQy()
21
22        LoadWave/G/D/W/A
23        String fileName = S_fileName
24        String path = S_Path
25        Variable numCols = V_flag
26
27        String w0,w1,w2,n0,n1,n2
28        String w3,w4,w5,w6,n3,n4,n5,n6
29               
30        // put the names of the three loaded waves into local names
31        n0 = StringFromList(0, S_waveNames ,";" )
32        n1 = StringFromList(1, S_waveNames ,";" )
33        n2 = StringFromList(2, S_waveNames ,";" )
34        n3 = StringFromList(3, S_waveNames ,";" )
35        n4 = StringFromList(4, S_waveNames ,";" )
36        n5 = StringFromList(5, S_waveNames ,";" )
37        n6 = StringFromList(6, S_waveNames ,";" )
38       
39        //remove the semicolon AND period from files from the VAX
40        w0 = CleanupName((S_fileName + "_qx"),0)
41        w1 = CleanupName((S_fileName + "_qy"),0)
42        w2 = CleanupName((S_fileName + "_i"),0)
43        w3 = CleanupName((S_fileName + "_qz"),0)
44        w4 = CleanupName((S_fileName + "_sigQx"),0)
45        w5 = CleanupName((S_fileName + "_sigQy"),0)
46        w6 = CleanupName((S_fileName + "_fs"),0)
47
48        String baseStr=w1[0,strlen(w1)-4]
49        if(DataFolderExists("root:"+baseStr))
50                        DoAlert 1,"The file "+S_filename+" has already been loaded. Do you want to load the new data file, overwriting the data in memory?"
51                        if(V_flag==2)   //user selected No, don't load the data
52                                SetDataFolder root:
53                                KillWaves $n0,$n1,$n2,$n3,$n4,$n5,$n6           // kill the default waveX that were loaded
54                                return          //quits the macro
55                        endif
56                        SetDataFolder $("root:"+baseStr)
57        else
58                NewDataFolder/S $("root:"+baseStr)
59        endif
60       
61        //read in the 18 lines of header (18th line starts w/ ASCII... 19th line is blank)
62        Make/O/T/N=18 header
63        Variable refnum,ii
64        string tmpStr=""
65        Open/R refNum  as (path+filename)
66        ii=0
67        do
68                tmpStr = ""
69                FReadLine refNum, tmpStr
70                header[ii] = tmpStr
71                ii+=1
72        while(ii < 18)         
73        Close refnum
74       
75        // ? parse to get nice variable names? put these in a structure just for fun?
76       
77       
78        ////overwrite the existing data, if it exists
79        Duplicate/O $("root:"+n0), $w0
80        Duplicate/O $("root:"+n1), $w1
81        Duplicate/O $("root:"+n2), $w2
82        Duplicate/O $("root:"+n3), $w3
83        Duplicate/O $("root:"+n4), $w4
84        Duplicate/O $("root:"+n5), $w5
85        Duplicate/O $("root:"+n6), $w6
86       
87        Variable/G gIsLogScale = 0
88       
89        Variable num=numpnts($w0)
90        // assume that the Q-grid is "uniform enough" for DISPLAY ONLY
91        // use the 3 original waves for all of the fitting...
92        ConvertQxQy2Mat($w0,$w1,$w2,baseStr+"_mat")
93        Duplicate/O $(baseStr+"_mat"),$(baseStr+"_lin")                 //keep a linear-scaled version of the data
94       
95        PlotQxQy(baseStr)               //this sets the data folder back to root:!!
96       
97        //don't create the triplet - not really of much use
98//      Make/D/O/N=(num,3) $(baseStr+"_tri")
99//      $(baseStr+"_tri")[][0] = $w0[p]         // Qx
100//      $(baseStr+"_tri")[][1] = $w1[p]         // Qy
101//      $(baseStr+"_tri")[][2] = $w2[p]         // Intensity
102
103        //clean up             
104        SetDataFolder root:
105        KillWaves/Z $n0,$n1,$n2,$n3,$n4,$n5,$n6
106EndMacro
107
108//does not seem to need to be flipped at all from the standard QxQy output
109//
110// plots the data, as surface0, YellowHot color table
111//
112Proc PlotQxQy(str)
113        String str
114        Prompt str,"Pick the data folder containing 2D data",popup,getAList(4)
115
116        PauseUpdate; Silent 1   // Building Gizmo 6 window...
117
118        // Do nothing if the Gizmo XOP is not available.
119        if(exists("NewGizmo")!=4)
120                DoAlert 0, "Gizmo XOP must be installed"
121                return
122        endif
123
124        String data = "root:"+str+":"+str+"_mat"
125       
126        NewGizmo/N=$str/T=str /W=(169,44,495,359)
127        ModifyGizmo startRecMacro
128        AppendToGizmo Surface=$data,name=surface0
129
130        // for constant color (data is a blue color, lighter on the bottom)
131//      ModifyGizmo ModifyObject=surface0 property={ surfaceColorType,1}
132//      ModifyGizmo ModifyObject=surface0 property={ srcMode,0}
133//      ModifyGizmo ModifyObject=surface0 property={ frontColor,0.749996,0.811093,1,1}
134//      ModifyGizmo ModifyObject=surface0 property={ backColor,0.250019,0.433326,1,1}
135        //e constant color
136       
137        //for color table       
138        ModifyGizmo ModifyObject=surface0 property={ srcMode,0}
139//      ModifyGizmo ModifyObject=surface0 property={ surfaceCTab,Blue}
140        ModifyGizmo ModifyObject=surface0 property={ surfaceCTab,ColdWarm}
141        ModifyGizmo ModifyObject=surface0 property={ SurfaceCTABScaling,4}
142        //e color table
143       
144        AppendToGizmo Axes=boxAxes,name=axes0
145        ModifyGizmo ModifyObject=axes0,property={0,axisRange,-1,-1,-1,1,-1,-1}
146        ModifyGizmo ModifyObject=axes0,property={1,axisRange,-1,-1,-1,-1,1,-1}
147        ModifyGizmo ModifyObject=axes0,property={2,axisRange,-1,-1,-1,-1,-1,1}
148        ModifyGizmo ModifyObject=axes0,property={3,axisRange,-1,1,-1,-1,1,1}
149        ModifyGizmo ModifyObject=axes0,property={4,axisRange,1,1,-1,1,1,1}
150        ModifyGizmo ModifyObject=axes0,property={5,axisRange,1,-1,-1,1,-1,1}
151        ModifyGizmo ModifyObject=axes0,property={6,axisRange,-1,-1,1,-1,1,1}
152        ModifyGizmo ModifyObject=axes0,property={7,axisRange,1,-1,1,1,1,1}
153        ModifyGizmo ModifyObject=axes0,property={8,axisRange,1,-1,-1,1,1,-1}
154        ModifyGizmo ModifyObject=axes0,property={9,axisRange,-1,1,-1,1,1,-1}
155        ModifyGizmo ModifyObject=axes0,property={10,axisRange,-1,1,1,1,1,1}
156        ModifyGizmo ModifyObject=axes0,property={11,axisRange,-1,-1,1,1,-1,1}
157        ModifyGizmo ModifyObject=axes0,property={-1,axisScalingMode,1}
158        ModifyGizmo ModifyObject=axes0,property={-1,axisColor,0,0,0,1}
159        ModifyGizmo ModifyObject=axes0,property={0,ticks,3}
160        ModifyGizmo ModifyObject=axes0,property={1,ticks,3}
161        ModifyGizmo ModifyObject=axes0,property={2,ticks,3}
162        ModifyGizmo ModifyObject=axes0,property={0,fontScaleFactor,1.5}
163        ModifyGizmo ModifyObject=axes0,property={1,fontScaleFactor,1.5}
164        ModifyGizmo ModifyObject=axes0,property={2,fontScaleFactor,1.5}
165        ModifyGizmo setDisplayList=0, opName=clearColor0, operation=clearColor, data={1,0.917,0.75,1}
166        ModifyGizmo setDisplayList=1, object=surface0
167        ModifyGizmo setDisplayList=2, object=axes0
168        ModifyGizmo SETQUATERNION={0.521287,-0.097088,-0.138769,0.836408}
169        ModifyGizmo autoscaling=1
170        ModifyGizmo currentGroupObject=""
171        ModifyGizmo compile
172
173        ModifyGizmo NamedHook={GizmoContours,WMGizmoContoursNamedHook}
174        ModifyGizmo endRecMacro
175       
176// don't bother with the flat image plot right now
177       
178//      Display $yw vs $xw
179//      modifygraph log=0
180//      ModifyGraph mode=3,marker=16,zColor($yw)={$zw,*,*,YellowHot,0}
181//      ModifyGraph standoff=0
182//      ModifyGraph width={Aspect,1}
183//      ModifyGraph lowTrip=0.001
184
185End
186
187Proc FakeQxQy(minQx,maxQx,minQy,maxQy,numPix,dataFolder)
188        Variable minQx=-0.1,maxQx=0.1,minQy=-0.1,maxQy=0.1
189        Variable numPix=64
190        String dataFolder="fake2DData"
191
192        String baseStr=dataFolder
193       
194
195        if(DataFolderExists("root:"+baseStr))
196                        DoAlert 1,"Fake data "+baseStr+" has already been created. Do you want to overwrite this fake data set?"
197                        if(V_flag==2)   //user selected No, don't load the data
198                                SetDataFolder root:
199                                return          //quits the macro
200                        endif
201                        SetDataFolder $("root:"+baseStr)
202        else
203                NewDataFolder/S $("root:"+baseStr)
204        endif
205
206        Make/O/D/N=(numPix*numPix) $(baseStr+"_qx"),$(baseStr+"_qy"),$(baseStr+"_i")
207       
208        $(baseStr+"_i") = 1
209       
210        Make/O/D/N=(numPix) tmpX,tmpY
211       
212        variable delX,delY
213        delX = (maxQx - minQx)/numPix
214        delY = (maxQy - minQy)/numPix
215        tmpX = minQx + x*delX
216        tmpY = minQy + x*delY
217       
218        // get rid of Q=0 values in the waves (may be a bad actor in model calculations)
219        //
220        tmpX = (tmpX == 0) ? 1e-6 : tmpX
221        tmpY = (tmpY == 0) ? 1e-6 : tmpY
222       
223        // X wave varies more rapidly
224        // Y wave is blocks of values
225        Variable ii=0
226        do
227                $(baseStr+"_qy")[ii*numPix,(ii+1)*numPix] = tmpY[ii]
228                ii+=1
229        while(ii<numPix)
230       
231        $(baseStr+"_qx") = tmpX[mod(p,numPix)]
232       
233        Variable/G gIsLogScale = 0
234       
235        // assume that the Q-grid is "uniform enough" for DISPLAY ONLY
236        // use the 3 original waves for all of the fitting...
237        ConvertQxQy2Mat($(baseStr+"_qx"),$(baseStr+"_qy"),$(baseStr+"_i"),baseStr+"_mat")
238        Duplicate/O $(baseStr+"_mat"),$(baseStr+"_lin")                 //keep a linear-scaled version of the data
239       
240        PlotQxQy(baseStr)               //this sets the data folder back to root:!!
241
242        //clean up             
243        SetDataFolder root:
244       
245EndMacro
246
247// this assumes that:
248// --QxQy data was written out in the format specified by the Igor macros, that is the x varies most rapidly
249// --the matrix is square!
250//
251// probably some other stuff...
252//
253Function UpdateQxQy2Mat(Qx,Qy,inten,linMat,mat)
254        Wave Qx,Qy,inten,linMat,mat
255       
256        Variable xrows=DimSize(mat, 0 )                 //assumes square matrix!!
257       
258        String folderStr=GetWavesDataFolder(Qx,1)
259        NVAR gIsLogScale=$(folderStr+"gIsLogScale")
260       
261        linMat = inten[q*xrows+p]
262       
263        if(gIsLogScale)
264                mat = log(linMat)
265        else
266                mat = linMat
267        endif
268       
269        return(0)
270End
271
272// this assumes that:
273// --QxQy data was written out in the format specified by the Igor macros, that is the x varies most rapidly
274// --the matrix is square!
275//
276// probably some other stuff...
277//
278Function ConvertQxQy2Mat(Qx,Qy,inten,matStr)
279        Wave Qx,Qy,inten
280        String matStr
281       
282        String folderStr=GetWavesDataFolder(Qx,1)
283       
284        Variable num=sqrt(numpnts(Qx))          //assumes square matrix, Qx = num x num points long
285        Make/O/D/N=(num,num) $(folderStr + matStr)
286        Wave mat=$matStr
287       
288        WaveStats/Q Qx
289        SetScale/I x, V_min, V_max, "", mat
290        WaveStats/Q Qy
291        SetScale/I y, V_min, V_max, "", mat
292       
293        Variable xrows=num
294       
295        mat = inten[q*xrows+p]
296       
297        return(0)
298End
299
300
301//str is the full path to the surface to append
302Proc AppendSurfaceToGizmo(str)
303        String str
304       
305        PauseUpdate; Silent 1   // Building Gizmo 6 window...
306
307        AppendToGizmo/Z Surface=$str,name=surface1
308       
309        // for a constant color (model is darker on top, lighter on the bottom)
310        //need these two lines, plus a color
311        ModifyGizmo ModifyObject=surface1 property={ surfaceColorType,1}
312        ModifyGizmo ModifyObject=surface1 property={ srcMode,0}
313        //green
314        ModifyGizmo ModifyObject=surface1 property={ frontColor,0.528923,0.882353,0.321584,1}
315        ModifyGizmo ModifyObject=surface1 property={ backColor,0.300221,0.6,1.5259e-05,1}
316        //red
317//      ModifyGizmo ModifyObject=surface1 property={ frontColor,1,0.749996,0.749996,1}
318//      ModifyGizmo ModifyObject=surface1 property={ backColor,1,0.250019,0.250019,1}
319
320        ModifyGizmo ModifyObject=surface1 property={ SurfaceCTABScaling,4}
321       
322        // for a color table (maybe not a good choice for vizualization if data uses a color table)
323//      ModifyGizmo ModifyObject=surface1 property={ srcMode,0}
324//      ModifyGizmo ModifyObject=surface1 property={ surfaceCTab,Red}
325       
326//      ModifyGizmo setDisplayList=0, object=surface0
327//      ModifyGizmo setDisplayList=1, object=axes0
328        ModifyGizmo setDisplayList=3, object=surface1
329        ModifyGizmo SETQUATERNION={0.565517,-0.103105,-0.139134,0.806350}
330//      ModifyGizmo autoscaling=1
331//      ModifyGizmo currentGroupObject=""
332        ModifyGizmo compile
333
334//      ModifyGizmo endRecMacro
335End
336
337// would be nice, but I can't get this to work...
338//
339//Macro AdjustColorTables()
340//
341//      ModifyGizmo ModifyObject=surface1 property={ srcMode,0}
342//      ModifyGizmo ModifyObject=surface1 property={ surfaceCTab,Rainbow}
343//      ModifyGizmo setDisplayList=3, object=surface1
344////    ModifyGizmo SETQUATERNION={0.565517,-0.103105,-0.139134,0.806350}
345//      ModifyGizmo autoscaling=1
346//      ModifyGizmo currentGroupObject=""
347//      ModifyGizmo compile
348//     
349//end
350
351Function LogToggle2DButtonProc(ba) : ButtonControl
352        STRUCT WMButtonAction &ba
353
354        switch( ba.eventCode )
355                case 2: // mouse up
356                        // click code here
357                        String winStr=TopGizmoWindow()
358                        if(strlen(winstr)==0)
359                                Print "no gizmo window"
360                                break
361                        endif
362                        // the winStr is also the data folder
363                        // toggle everything in the data folder
364                        ToggleFolderScaling(winStr)             
365                        break
366        endswitch
367
368        return 0
369End
370
371Function Plot2DButtonProc(ba) : ButtonControl
372        STRUCT WMButtonAction &ba
373
374        switch( ba.eventCode )
375                case 2: // mouse up
376                        // click code here
377                        Execute "PlotQxQy()"
378                        break
379        endswitch
380
381        return 0
382End
383
384Function Append2DModelButtonProc(ba) : ButtonControl
385        STRUCT WMButtonAction &ba
386
387        switch( ba.eventCode )
388                case 2: // mouse up
389                        // click code here
390                        //str is the full path to the surface to append
391                        String DF=TopGizmoWindow()
392                        String str,funcStr,suffix
393                        ControlInfo/W=WrapperPanel popup_1
394                        funcStr = S_Value
395                        suffix = getModelSuffix(getFunctionCoef(funcStr))
396                        str = "root:"+DF+":"+suffix+"_mat"
397                       
398                        Execute "AppendSurfaceToGizmo(\""+str+"\")"
399//                      Print str
400                        break
401        endswitch
402
403        return 0
404End
405
406// toggle the scaling of every matrix in the folder
407Function ToggleFolderScaling(DF)
408        String DF
409
410        // look for waves DF+"_mat" (the data)
411        // and models prefix + "_mat"
412        SetDataFolder $("root:"+DF)
413        //check the global to see the state of the data
414        NVAR gIsLogScale = gIsLogScale
415       
416        String matrixList=WaveList("*_mat",";",""),item
417        Variable ii=0,num=itemsinlist(matrixlist,";"),len
418       
419        for(ii=0;ii<num;ii+=1)
420                        item = StringFromList(ii,matrixList,";")
421                        len = strlen(item)
422                        Wave w = $item
423                        Wave linW = $(item[0,len-5]+"_lin")
424                        if(gisLogScale)
425                                //make linear
426                                w=linW
427                        else   
428                                //make log scale
429                                w=log(linW)
430                                w = w[p][q] == -inf ? NaN : w[p][q]             //remove the -inf for display
431                        endif
432        endfor 
433       
434        // toggle the global
435        gIsLogScale = !gIsLogScale
436       
437        SetDataFolder root:
438        return(0)
439End
440
441
442Function Load2DDataButtonProc(ba) : ButtonControl
443        STRUCT WMButtonAction &ba
444
445        switch( ba.eventCode )
446                case 2: // mouse up
447                        // click code here
448                        Execute "LoadQxQy()"
449                        break
450        endswitch
451
452        return 0
453End
454
455
456// dispatch the fit
457Function Do2DFitButtonProc(ba) : ButtonControl
458        STRUCT WMButtonAction &ba
459
460        String folderStr,funcStr,coefStr
461        Variable useCursors,useEps,useConstr
462       
463        switch( ba.eventCode )
464                case 2: // mouse up
465                        ControlInfo/W=WrapperPanel popup_0
466                        folderStr=S_Value
467                       
468                        ControlInfo/W=WrapperPanel popup_1
469                        funcStr=S_Value
470                       
471                        ControlInfo/W=WrapperPanel popup_2
472                        coefStr=S_Value
473                       
474                        ControlInfo/W=WrapperPanel check_0
475                        useCursors=V_Value
476                        //
477                        // NO CURSORS for 2D waves - force to zero
478                        //
479                        useCursors = 0
480                        //
481                        ControlInfo/W=WrapperPanel check_1
482                        useEps=V_Value
483                        ControlInfo/W=WrapperPanel check_2
484                        useConstr=V_Value
485                       
486                        if(!CheckFunctionAndCoef(funcStr,coefStr))
487                                DoAlert 0,"The coefficients and function type do not match. Please correct the selections in the popup menus."
488                                break
489                        endif
490                       
491                        FitWrapper2D(folderStr,funcStr,coefStr,useCursors,useEps,useConstr)
492                       
493                        //      DoUpdate (does not work!)
494                        //?? why do I need to force an update ??
495                        if(!exists("root:"+folderStr+":"+coefStr))
496                                Wave w=$coefStr
497                        else
498                                Wave w=$("root:"+folderStr+":"+coefStr) //smeared coefs in data folder
499                        endif
500                        w[0] += 1e-6
501                        w[0] -= 1e-6
502       
503                        break
504        endswitch
505
506        return 0
507End
508
509// wrapper to do the desired fit, pretty much a clone of FitWrapper, the 1D version
510//
511// folderStr is the data folder for the desired data set
512//
513// Currently the limitations are:
514// - I have no error waves for the intensity
515// - There is no smeared model
516// - Cursors can't be used
517// - the report *probably* won't work
518Function FitWrapper2D(folderStr,funcStr,coefStr,useCursors,useEps,useConstr)
519        String folderStr,funcStr,coefStr
520        Variable useCursors,useEps,useConstr
521
522        String suffix=getModelSuffix(coefStr)
523       
524        SetDataFolder $("root:"+folderStr)
525        if(!exists(coefStr))
526                // must be unsmeared model, work in the root folder
527                SetDataFolder root:                             
528                if(!exists(coefStr))            //this should be fine if the coef filter is working, but check anyhow
529                        DoAlert 0,"the coefficient and data sets do not match"
530                        return 0
531                endif
532        endif
533               
534        WAVE cw=$(coefStr)     
535        Wave hold=$("Hold_"+suffix)
536        Wave/T lolim=$("LoLim_"+suffix)
537        Wave/T hilim=$("HiLim_"+suffix)
538        Wave eps=$("epsilon_"+suffix)
539       
540// fill a struct instance whether I need one or not
541        String DF="root:"+folderStr+":"
542       
543        Struct ResSmearAAOStruct fs
544//      WAVE resW = $(DF+folderStr+"_res")             
545//      WAVE fs.resW =  resW
546        WAVE inten=$(DF+folderStr+"_i")
547        WAVE Qx=$(DF+folderStr+"_qx")
548        WAVE Qy=$(DF+folderStr+"_qy")
549//      Wave fs.coefW = cw
550//      Wave fs.yW = yw
551//      Wave fs.xW = xw
552       
553        // generate my own error wave for I(qx,qy)
554        Duplicate/O inten sw
555        sw = sqrt(sw)           //assumes Poisson statistics for each cell (counter)
556//      sw = 0.05*sw            // uniform 5% error? tends to favor the low intensity too strongly
557
558        // generate my own mask wave - as a matrix first, then redimension to N*N vector
559        // default mask is two pixels all the way around, (0 is excluded, 1 is included)
560        WAVE DataMat=$(DF+folderStr+"_lin")
561        if(exists(DF+"mask") == 0)
562                Duplicate/O dataMat mask
563                Variable bsRadius=8             //pixels?
564                MakeBSMask(mask,bsRadius)
565        Endif
566       
567
568       
569//      Duplicate/O yw $(DF+"FitYw")
570//      WAVE fitYw = $(DF+"FitYw")
571//      fitYw = NaN
572
573        //for now, use res will always be 0 for 2D functions   
574        Variable useRes=0
575        if(stringmatch(funcStr, "Smear*"))              // if it's a smeared function, need a struct
576                useRes=1
577        endif
578       
579        // do not construct constraints for any of the coefficients that are being held
580        // -- this will generate an "unknown error" from the curve fitting
581        Make/O/T/N=0 constr
582        if(useConstr)
583                String constraintExpression
584                Variable i, nPnts=DimSize(lolim, 0),nextRow=0
585                for (i=0; i < nPnts; i += 1)
586                        if (strlen(lolim[i]) > 0 && hold[i] == 0)
587                                InsertPoints nextRow, 1, constr
588                                sprintf constraintExpression, "K%d > %s", i, lolim[i]
589                                constr[nextRow] = constraintExpression
590                                nextRow += 1
591                        endif
592                        if (strlen(hilim[i]) > 0 && hold[i] == 0)
593                                InsertPoints nextRow, 1, constr
594                                sprintf constraintExpression, "K%d < %s", i, hilim[i]
595                                constr[nextRow] = constraintExpression
596                                nextRow += 1
597                        endif
598                endfor
599        endif
600
601///// NO CURSORS for 2D waves
602        //if useCursors, and the data is USANS, need to feed a (reassigned) trimmed matrix to the fit
603//      Variable pt1,pt2,newN
604//      if(useCursors && (dimsize(resW,1) > 4) )
605//              if(pcsr(A) > pcsr(B))
606//                      pt1 = pcsr(B)
607//                      pt2 = pcsr(A)
608//              else
609//                      pt1 = pcsr(A)
610//                      pt2 = pcsr(B)
611//              endif
612//              newN = pt2 - pt1 + 1            // +1 includes both cursors in the fit
613//              Make/O/D/N=(newN,newN) $(DF+"crsrResW")
614//              WAVE crsrResW = $(DF+"crsrResW")
615//              crsrResW = resW[p+pt1][q+pt1]
616//              //assign to the struct
617//              WAVE fs.resW =  crsrResW               
618//      endif
619
620// create these variables so that FuncFit will set them on exit
621        Variable/G V_FitError=0                         //0=no err, 1=error,(2^1+2^0)=3=singular matrix
622        Variable/G V_FitQuitReason=0            //0=ok,1=maxiter,2=user stop,3=no chisq decrease
623       
624// don't use the auto-destination with no flag, it doesn't appear to work correctly
625// dispatch the fit
626        //      FuncFit/H="11110111111"/NTHR=0 Cylinder2D_D :cyl2d_c_txt:coef_Cyl2D_D  :cyl2d_c_txt:cyl2d_c_txt_i /X={:cyl2d_c_txt:cyl2d_c_txt_qy,:cyl2d_c_txt:cyl2d_c_txt_qx} /W=:cyl2d_c_txt:sw /I=1 /M=:cyl2d_c_txt:mask /D
627
628        do
629                if(useRes && useEps && useCursors && useConstr)         //do it all
630                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr
631                        break
632                endif
633               
634                if(useRes && useEps && useCursors)              //no constr
635                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps
636                        break
637                endif
638               
639                if(useRes && useEps && useConstr)               //no crsr
640                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr
641                        break
642                endif
643               
644                if(useRes && useCursors && useConstr)           //no eps
645                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr
646                        break
647                endif
648               
649                if(useRes && useCursors)                //no eps, no constr
650                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D
651                        break
652                endif
653               
654                if(useRes && useEps)            //no crsr, no constr
655                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps
656                        break
657                endif
658       
659                if(useRes && useConstr)         //no crsr, no eps
660                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr
661                        break
662                endif
663               
664                if(useRes)              //just res
665                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D
666                        break
667                endif
668               
669/////   same as above, but all without useRes (no /STRC flag)
670                if(useEps && useCursors && useConstr)           //do it all
671                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr
672                        break
673                endif
674               
675                if(useEps && useCursors)                //no constr
676                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps
677                        break
678                endif
679               
680               
681                if(useEps && useConstr)         //no crsr
682                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr
683                        break
684                endif
685               
686                if(useCursors && useConstr)             //no eps
687                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr
688                        break
689                endif
690               
691                if(useCursors)          //no eps, no constr
692                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D
693                        break
694                endif
695               
696                if(useEps)              //no crsr, no constr
697                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps
698                        break
699                endif
700       
701                if(useConstr)           //no crsr, no eps
702                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr
703                        break
704                endif
705               
706                //just a plain vanilla fit
707                        FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D
708       
709        while(0)
710       
711        // append the fit
712        // need to manage duplicate copies
713        // Don't plot the full curve if cursors were used (set fitYw to NaN on entry...)
714//      String traces=TraceNameList("", ";", 1 )                //"" as first parameter == look on the target graph
715//      if(strsearch(traces,"FitYw",0) == -1)
716//              AppendToGraph FitYw vs xw
717//      else
718//              RemoveFromGraph FitYw
719//              AppendToGraph FitYw vs xw
720//      endif
721//      ModifyGraph lsize(FitYw)=2,rgb(FitYw)=(0,0,0)
722       
723//      DoUpdate                //force update of table and graph with fitted values (why doesn't this work? - the table still does not update)
724       
725        // report the results (to the panel?)
726        print "V_chisq = ",V_chisq
727        print cw
728        WAVE w_sigma
729        print w_sigma
730        String resultStr=""
731               
732        //now re-write the results
733        sprintf resultStr,"Chi^2 = %g  Sqrt(X^2/N) = %g",V_chisq,sqrt(V_chisq/V_Npnts)
734        resultStr = PadString(resultStr,63,0x20)
735        GroupBox grpBox_2 title=resultStr
736        ControlUpdate/W=WrapperPanel grpBox_2
737        sprintf resultStr,"FitErr = %s : FitQuit = %s",W_ErrorMessage(V_FitError),W_QuitMessage(V_FitQuitReason)
738        resultStr = PadString(resultStr,63,0x20)
739        GroupBox grpBox_3 title=resultStr
740        ControlUpdate/W=WrapperPanel grpBox_3
741       
742        Variable yesSave=0,yesReport=0
743        ControlInfo/W=WrapperPanel check_4
744        yesReport = V_Value
745        ControlInfo/W=WrapperPanel check_5
746        yesSave = V_Value
747       
748       
749        if(yesReport)
750                String parStr=GetWavesDataFolder(cw,1)+ WaveList("*param*"+suffix, "", "TEXT:1," )              //this is *hopefully* one wave
751                String topGraph= TopGizmoWindow()       //this is the topmost Gizmo (XOP) window
752               
753                DoUpdate                //force an update of the graph before making a copy of it for the report
754       
755                W_GenerateReport(funcStr,folderStr,$parStr,cw,yesSave,V_chisq,W_sigma,V_npnts,V_FitError,V_FitQuitReason,V_startRow,V_endRow,topGraph)
756        endif
757       
758        SetDataFolder root:
759        return(0)
760End
761
762
763// right now, there are only unsmeared models, but I want experimental points
764// for the QxQy of the calculation - so get the folder string
765Function Plot2DFunctionButtonProc(ba) : ButtonControl
766        STRUCT WMButtonAction &ba
767
768        String folderStr,funcStr,coefStr,cmdStr=""
769       
770        Variable killWhat=0             //kill nothing as default
771
772        switch( ba.eventCode )
773                case 2: // mouse up
774                        ControlInfo/W=WrapperPanel popup_0
775                        folderStr=S_Value
776                       
777                        ControlInfo/W=WrapperPanel popup_1
778                        funcStr=S_Value
779                       
780                        // check for smeared or smeared function
781//                      if(stringmatch(funcStr, "Smear*" )==1)
782                                //it's a smeared model
783                                // check for the special case of RPA that has an extra parameter
784//                              if(strsearch(funcStr, "RPAForm", 0 ,0) == -1)
785                                        sprintf cmdStr, "Plot%s(\"%s\")",funcStr,folderStr              //not RPA
786//                              else
787//                                      sprintf cmdStr, "Plot%s(\"%s\",)",funcStr,folderStr             //yes RPA, leave a comma for input
788//                              endif
789//                      else
790//                              // it's not,                   
791//                              sprintf cmdStr, "Plot%s()",funcStr
792//                      endif
793                       
794                        //Print cmdStr
795                        Execute cmdStr
796                       
797                        //pop the function menu to set the proper coefficients
798                        DoWindow/F WrapperPanel
799                        STRUCT WMPopupAction pa
800                        pa.popStr = funcStr
801                        pa.eventcode = 2
802                        Function_PopMenuProc(pa)
803       
804                        KillWhat = 2            //kill just the table, leave the 2d visible for now
805                        KillTopGraphAndTable(killWhat)          // crude
806
807                        break
808        endswitch
809
810        return 0
811End
812
813// unused - all of this functionality has been added to the Wrapper Panel
814Window Plot_2D_Controls() : Panel
815        PauseUpdate; Silent 1           // building window...
816        NewPanel /W=(950,44,1250,244)
817        Button button0,pos={172,37},size={70,20},proc=LogToggle2DButtonProc,title="Log/Lin"
818        Button button1,pos={146,9},size={100,20},proc=Plot2DButtonProc,title="Plot 2D Data"
819        Button button2,pos={164,118},size={120,20},proc=Append2DModelButtonProc,title="Append Model"
820        PopupMenu popup_1,pos={9,84},size={132,20},proc=Function_PopMenuProc,title="Function"
821        PopupMenu popup_1,mode=1,popvalue="Cylinder2D",value= #"W_FunctionPopupList()"
822        Button button3,pos={11,9},size={100,20},proc=Load2DDataButtonProc,title="Load 2D Data"
823        Button button4,pos={164,161},size={120,20},proc=Do2DFitButtonProc,title="Do 2D Fit"
824        Button button5,pos={164,85},size={120,20},proc=Plot2DFunctionButtonProc,title="Plot 2D Function"
825EndMacro
826
827
828Function/S TopGizmoWindow()
829        Return(StringFromList(0,WinList("*",";","WIN:4096")))
830end
831
832////
833//
834// unused proc for testing
835Proc Plot2DVsPointNumber(str)
836        String str
837        Prompt str,"Pick the data folder containing 2D data",popup,getAList(4)
838
839        PauseUpdate; Silent 1           // building window...
840        String fldrSav0= GetDataFolder(1)
841        SetDataFolder $("root:"+str)
842        Display /W=(241,352,810,842) $(str+"_i")
843        SetDataFolder fldrSav0
844        ModifyGraph mode=3
845        ModifyGraph marker=19
846        ModifyGraph lSize=2
847        ModifyGraph rgb($(str+"_i"))=(0,0,0)
848        ModifyGraph msize=1
849        ModifyGraph grid=1
850        ModifyGraph log(left)=1
851        ModifyGraph mirror=2
852        if(exists("root:"+str+":sw"))
853                ErrorBars/T=0 $(str+"_i") Y,wave=($("root:"+str+":sw"),$("root:"+str+":sw"))
854        endif
855        SetDataFolder fldrSav0
856EndMacro
857
858// testing procedure, used when I had x,y swapped and was getting nonsensical results
859//Macro CalculateChiSquared(str)
860//      String str
861//      Prompt str,"Pick the data folder containing 2D data",popup,getAList(4)
862//
863//
864//      String fldrSav0= GetDataFolder(1)
865//      SetDataFolder $("root:"+str)
866//
867//      Duplicate/O $(str+"_i") chi
868//      chi = ((zwave_cyl2D_D - $(str+"_i"))/sw )^2
869//     
870//      chi = (mask == 1) ? chi : 0
871//     
872//      Print sum(chi,-inf,inf)
873//     
874//      SetDataFolder fldrSav0
875//EndMacro
876
877
878Function MakeBSMask(mask,rad)
879        Wave mask
880        Variable rad
881
882
883// find the center based on the wave scaling
884        Variable xCtr,yCtr,Qzero=0
885        xCtr = (Qzero - DimOffset(Mask, 0))/DimDelta(Mask,0)
886        yCtr = (Qzero - DimOffset(Mask, 1))/DimDelta(Mask,1)
887
888        Print xctr,yctr
889
890//      Variable center = sqrt(numpnts(mask))/2 -0.5
891
892        mask = (sqrt((p-xCtr)^2+(q-yCtr)^2) < rad) ? 0 : 1
893
894        Variable xDim, yDim
895        xDim = DimSize(mask,0)
896        yDim = DimSize(mask,1)
897        mask[][0] = 0
898        mask[][1] = 0
899        mask[][yDim-2] = 0
900        mask[][yDim-1] = 0
901        mask[0][] = 0
902        mask[1][] = 0
903        mask[xDim-2][] = 0
904        mask[xDim-1][] = 0
905       
906        Redimension/N=(xDim*yDim) mask          //now 1D
907       
908       
909End
Note: See TracBrowser for help on using the repository browser.