Ignore:
Timestamp:
Jun 16, 2010 2:10:47 PM (12 years ago)
Author:
srkline
Message:

SA_Includes_v410 : now include Smear_2D

PeakGauss_2D, Sphere_2D : included threaded resolution smearing calculations for testing

DataSetHandling? : Included a quick and dirty batch converter for XML->6-col. See the top
of the file for the command to run

GaussUtils? : re-define the ResSemear_2D_AAOStruct. Relocated q-value and phi calculations from
RawWindowHook? to this file so they would be available for reduction and analysis

Smear_2D : now has a generic (non-threaded) smearing routine. Threading must be done in
individual functions since FUNCREF can't be passed to threads (plus a few other issues)

PlotUtils_2D : updated loader for new QxQy? columns. Fixes to Wrapper_2D to enable smeared fits

RawWindowHook? : removed Q-value calculation functions and moved these to GaussUtils?

WriteQIS : now writes out 8-columns for QxQy? data, defining the resolution
function in terms of directions parallel and perpendicular to Q. TEMPORARILY in the data
file an error in intensity is generated that is SQRT(I), being careful to
replace any NaN, inf, or zero with an average error value

MultiScatter_MonteCarlo_2D : 4-processor aware

NCNR_Utils : 2D resolution calculation is now in terms of parallel and perpendicular
rather than x and y. Gravity is included in the y-component

File:
1 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Common/Packages/PlotManager/PlotUtils2D_v40.ipf

    r665 r708  
    44 
    55// 
    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) 
     6// TO DO: 
     7// 
     8// - more intelligent beam stop masking 
     9// - constraints 
     10// - interactive masking 
     11// 
     12// 
     13// 
     14 
     15 
     16// 
     17// changed June 2010 to read in the resolution information (now 8! columns) 
    1818// -- subject to change -- 
    1919// 
     
    2525        Variable numCols = V_flag 
    2626 
    27         String w0,w1,w2,n0,n1,n2 
    28         String w3,w4,w5,w6,n3,n4,n5,n6 
     27        String w0,w1,w2,w3,w4,w5,w6,w7 
     28        String n0,n1,n2,n3,n4,n5,n6,n7 
    2929                 
    3030        // put the names of the three loaded waves into local names 
     
    3636        n5 = StringFromList(5, S_waveNames ,";" ) 
    3737        n6 = StringFromList(6, S_waveNames ,";" ) 
     38        n7 = StringFromList(7, S_waveNames ,";" ) 
    3839         
    3940        //remove the semicolon AND period from files from the VAX 
     
    4142        w1 = CleanupName((S_fileName + "_qy"),0) 
    4243        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) 
     44        w3 = CleanupName((S_fileName + "_iErr"),0) 
     45        w4 = CleanupName((S_fileName + "_qz"),0) 
     46        w5 = CleanupName((S_fileName + "_sQpl"),0) 
     47        w6 = CleanupName((S_fileName + "_sQpp"),0) 
     48        w7 = CleanupName((S_fileName + "_fs"),0) 
    4749 
    4850        String baseStr=w1[0,strlen(w1)-4] 
     
    5153                        if(V_flag==2)   //user selected No, don't load the data 
    5254                                SetDataFolder root: 
    53                                 KillWaves $n0,$n1,$n2,$n3,$n4,$n5,$n6           // kill the default waveX that were loaded 
     55                                KillWaves $n0,$n1,$n2,$n3,$n4,$n5,$n6,$n7               // kill the default waveX that were loaded 
    5456                                return          //quits the macro 
    5557                        endif 
     
    8486        Duplicate/O $("root:"+n5), $w5 
    8587        Duplicate/O $("root:"+n6), $w6 
     88        Duplicate/O $("root:"+n7), $w7 
    8689         
    8790        Variable/G gIsLogScale = 0 
     
    9497         
    9598        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 
    10299 
    103100        //clean up               
    104101        SetDataFolder root: 
    105         KillWaves/Z $n0,$n1,$n2,$n3,$n4,$n5,$n6 
     102        KillWaves/Z $n0,$n1,$n2,$n3,$n4,$n5,$n6,$n7 
    106103EndMacro 
    107104 
     
    547544        String DF="root:"+folderStr+":"  
    548545         
    549         Struct ResSmearAAOStruct fs 
    550 //      WAVE resW = $(DF+folderStr+"_res")               
    551 //      WAVE fs.resW =  resW 
    552546        WAVE inten=$(DF+folderStr+"_i") 
    553         WAVE Qx=$(DF+folderStr+"_qx") 
    554         WAVE Qy=$(DF+folderStr+"_qy") 
    555 //      Wave fs.coefW = cw 
    556 //      Wave fs.yW = yw 
    557 //      Wave fs.xW = xw 
    558          
    559         // generate my own error wave for I(qx,qy) 
    560         Duplicate/O inten sw 
    561         sw = sqrt(sw)           //assumes Poisson statistics for each cell (counter) 
    562 //      sw = 0.05*sw            // uniform 5% error? tends to favor the low intensity too strongly 
     547        WAVE sw=$(DF+folderStr+"_iErr") 
     548        WAVE qx=$(DF+folderStr+"_qx") 
     549        WAVE qy=$(DF+folderStr+"_qy") 
     550        WAVE qz=$(DF+folderStr+"_qz") 
     551        WAVE sQpl=$(DF+folderStr+"_sQpl") 
     552        WAVE sQpp=$(DF+folderStr+"_sQpp") 
     553        WAVE shad=$(DF+folderStr+"_fs") 
     554 
     555//just a dummy - I shouldn't need this 
     556        Duplicate/O qx resultW 
     557        resultW=0 
     558         
     559        STRUCT ResSmear_2D_AAOStruct s 
     560        WAVE s.coefW = cw        
     561        WAVE s.zw = resultW      
     562        WAVE s.xw[0] = qx 
     563        WAVE s.xw[1] = qy 
     564        WAVE s.qz = qz 
     565        WAVE s.sQpl = sQpl 
     566        WAVE s.sQpp = sQpp 
     567        WAVE s.fs = shad 
     568         
    563569 
    564570        // generate my own mask wave - as a matrix first, then redimension to N*N vector 
     
    571577        Endif 
    572578         
    573  
    574          
    575 //      Duplicate/O yw $(DF+"FitYw") 
    576 //      WAVE fitYw = $(DF+"FitYw") 
    577 //      fitYw = NaN 
     579         
     580        Duplicate/O inten inten_masked 
     581        inten_masked = (mask[p][q] == 0) ? NaN : inten[p][q] 
     582         
    578583 
    579584        //for now, use res will always be 0 for 2D functions     
     
    582587                useRes=1 
    583588        endif 
    584          
    585         // do not construct constraints for any of the coefficients that are being held 
    586         // -- this will generate an "unknown error" from the curve fitting 
    587         Make/O/T/N=0 constr 
     589 
     590        // can't use constraints in this way for multivariate fits. See the curve fitting help file 
     591        // and "Contraint Matrix and Vector" 
    588592        if(useConstr) 
    589                 String constraintExpression 
    590                 Variable i, nPnts=DimSize(lolim, 0),nextRow=0 
    591                 for (i=0; i < nPnts; i += 1) 
    592                         if (strlen(lolim[i]) > 0 && hold[i] == 0) 
    593                                 InsertPoints nextRow, 1, constr 
    594                                 sprintf constraintExpression, "K%d > %s", i, lolim[i] 
    595                                 constr[nextRow] = constraintExpression 
    596                                 nextRow += 1 
    597                         endif 
    598                         if (strlen(hilim[i]) > 0 && hold[i] == 0) 
    599                                 InsertPoints nextRow, 1, constr 
    600                                 sprintf constraintExpression, "K%d < %s", i, hilim[i] 
    601                                 constr[nextRow] = constraintExpression 
    602                                 nextRow += 1 
    603                         endif 
    604                 endfor 
    605         endif 
     593                Print "Constraints not yet implemented" 
     594                useConstr = 0 
     595        endif    
     596         
     597//      // do not construct constraints for any of the coefficients that are being held 
     598//      // -- this will generate an "unknown error" from the curve fitting 
     599//      Make/O/T/N=0 constr 
     600//      if(useConstr) 
     601//              String constraintExpression 
     602//              Variable i, nPnts=DimSize(lolim, 0),nextRow=0 
     603//              for (i=0; i < nPnts; i += 1) 
     604//                      if (strlen(lolim[i]) > 0 && hold[i] == 0) 
     605//                              InsertPoints nextRow, 1, constr 
     606//                              sprintf constraintExpression, "K%d > %s", i, lolim[i] 
     607//                              constr[nextRow] = constraintExpression 
     608//                              nextRow += 1 
     609//                      endif 
     610//                      if (strlen(hilim[i]) > 0 && hold[i] == 0) 
     611//                              InsertPoints nextRow, 1, constr 
     612//                              sprintf constraintExpression, "K%d < %s", i, hilim[i] 
     613//                              constr[nextRow] = constraintExpression 
     614//                              nextRow += 1 
     615//                      endif 
     616//              endfor 
     617//      endif 
    606618 
    607619///// NO CURSORS for 2D waves 
     
    631643// dispatch the fit 
    632644        //      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  
     645        Variable t1=StopMSTimer(-2) 
     646 
     647 
     648// /NTHR=1 means just one thread for the fit (since the function itself is threaded) 
    633649 
    634650        do 
    635651                if(useRes && useEps && useCursors && useConstr)         //do it all 
    636                         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 
     652                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /C=constr /STRC=s 
    637653                        break 
    638654                endif 
    639655                 
    640656                if(useRes && useEps && useCursors)              //no constr 
    641                         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 
     657                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /STRC=s 
    642658                        break 
    643659                endif 
    644660                 
    645661                if(useRes && useEps && useConstr)               //no crsr 
    646                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr 
     662                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /C=constr /STRC=s 
    647663                        break 
    648664                endif 
    649665                 
    650666                if(useRes && useCursors && useConstr)           //no eps 
    651                         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 
     667                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /C=constr /STRC=s 
    652668                        break 
    653669                endif 
    654670                 
    655671                if(useRes && useCursors)                //no eps, no constr 
    656                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D 
     672                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /STRC=s 
    657673                        break 
    658674                endif 
    659675                 
    660676                if(useRes && useEps)            //no crsr, no constr 
    661                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps 
     677                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /STRC=s 
    662678                        break 
    663679                endif 
    664680         
    665681                if(useRes && useConstr)         //no crsr, no eps 
    666                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr 
     682                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /C=constr /STRC=s 
    667683                        break 
    668684                endif 
    669685                 
    670686                if(useRes)              //just res 
    671                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D 
     687                        Print "useRes only" 
     688                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /STRC=s 
     689//                      FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /W=sw /I=1 /STRC=s 
     690//                      FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 
     691//                      FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten_masked /X={Qx,Qy} /W=sw /I=1 /STRC=s 
    672692                        break 
    673693                endif 
     
    675695/////   same as above, but all without useRes (no /STRC flag) 
    676696                if(useEps && useCursors && useConstr)           //do it all 
    677                         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 
     697                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /C=constr 
    678698                        break 
    679699                endif 
    680700                 
    681701                if(useEps && useCursors)                //no constr 
    682                         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 
     702                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps 
    683703                        break 
    684704                endif 
     
    686706                 
    687707                if(useEps && useConstr)         //no crsr 
    688                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps /C=constr 
     708                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps /C=constr 
    689709                        break 
    690710                endif 
    691711                 
    692712                if(useCursors && useConstr)             //no eps 
    693                         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 
     713                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /C=constr 
    694714                        break 
    695715                endif 
    696716                 
    697717                if(useCursors)          //no eps, no constr 
    698                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 /D 
     718                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten[pcsr(A),pcsr(B)] /X={Qx,Qy} /M=mask /W=sw /I=1 
    699719                        break 
    700720                endif 
    701721                 
    702722                if(useEps)              //no crsr, no constr 
    703                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /E=eps 
     723                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /E=eps 
    704724                        break 
    705725                endif 
    706726         
    707727                if(useConstr)           //no crsr, no eps 
    708                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D /C=constr 
     728                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /C=constr 
    709729                        break 
    710730                endif 
    711731                 
    712732                //just a plain vanilla fit 
    713                         FuncFit/H=getHStr(hold) /NTHR=0 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 /D 
     733                print "unsmeared vanilla" 
     734                        FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /M=mask /W=sw /I=1 
     735//                      FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten /X={Qx,Qy} /W=sw /I=1 
     736//                      FuncFit/H=getHStr(hold) /NTHR=1 $funcStr cw, inten_masked /X={Qx,Qy} /W=sw /I=1 
    714737         
    715738        while(0) 
    716739         
     740        Print "elapsed 2D fit time  = ",(StopMSTimer(-2) - t1)/1e6," s = ",(StopMSTimer(-2) - t1)/1e6/60," min" 
     741 
    717742        // append the fit 
    718743        // need to manage duplicate copies 
     
    754779         
    755780        if(yesReport) 
    756                 String parStr=GetWavesDataFolder(cw,1)+ WaveList("*param*"+suffix, "", "TEXT:1," )              //this is *hopefully* one wave 
     781                String parStr = getFunctionParams(funcStr) 
     782//              String parStr=GetWavesDataFolder(cw,1)+ WaveList("*param*"+suffix, "", "TEXT:1," )              //old way, but doesn't work in 2D folders 
    757783                String topGraph= TopGizmoWindow()       //this is the topmost Gizmo (XOP) window 
    758784                 
Note: See TracChangeset for help on using the changeset viewer.