Ignore:
Timestamp:
Apr 4, 2011 11:12:38 AM (12 years ago)
Author:
srkline
Message:

Lots of changes:
-2D resolution smearing
-2D error propagation

1) 2D resolution smearing has been corrected to use sigma (perp) correctly
rather than phi. This shows up in the quadrature loop in all of the 2D models
and in the Smear_2D "generic" function.

2) 1D resolutionn smearing is now corrected to account for the integration
range of +/- 3 sigma (99.73% of distribution). Now the result is divided by
0.9973 to rescale it to the correct value.

3) Smeared models are now AAO to improve speed and to allow easier use with
functions that are inherently AAO. No changes are needed, since the call is
behind the scenes, replacing Smear_Model_N() with Smear_Model_N_AAO().

4) in PlotUtils_2D added functions to re-bin the QxQy? data into a 1D format
BinQxQy_to_1D(). This also re-bins the errors in two ways, adding the per-pixel
errors in quadrature, or the deviation from the mean of the intensity. Manually
editing the intensity declaration allows 2D->1D binning of smeared models.

5) Per-pixel propagation of errors has been added through the entire data
reduction sequence. Data errors are generated on loading using Poisson
statistics (specifically tailored for accuracy at low counts) and then is
propagated through each manipulation of the data using standard error
propagation. The error matrix is linear_data_error. As a by-product, all of
the math operations on data are explicitly done on linear_data, to avoid
any potential mistakes of log/linear scaling. Results of this propagation
largely match J. Barker's /ERROR routines from the VAX, with some differences
at low data count values (as expected) and at higher count values near the
beam stop (somewhat unexpected). This per-pixel error is ouput in the QxQy_ASCII
data files. NO CHANGE has been made to the 1D data, which uses the deviation from
the mean as the error - since this is correct.

6) Added tables for the uncertainty in attenuator transmission (from JGB)

7) added two more REAL values to the VAX header to store information
necessary for error propagation. These are couting error that are part of
the transmission error and of the absolute scaling error. These add Read/Write?
functions in NCNR_DataReadWrite

The transmission error (one standard deviation) is written at byte 396 (4 bytes)

RealsRead?[41]

The Box Sum error (one standard deviation) is written at byte 400 (4 bytes)

RealsRead?[42]

File:
1 edited

Legend:

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

    r771 r794  
    1717//                                              2D model calculations (model_lin) 2D data files as well. 
    1818 
     19// Call the procedure BinQxQy_to_1D() to re-bin the QxQy data into I(q). good for testing to check the 2D error propagation 
     20// 
     21// 
     22 
    1923 
    2024// 
     
    2529Proc LoadQxQy() 
    2630 
     31        SetDataFolder root: 
     32         
    2733        LoadWave/G/D/W/A 
    2834        String fileName = S_fileName 
     
    10181024// be made too, if needed. 
    10191025// 
    1020 // X- need error on I(q) 
    1021 // -- need to set "proper" number of data points (delta and qMax?) 
     1026// X- need error on I(q) - now calculated both ways, from 2D error propagation, and from deviation from the mean 
     1027// X- need to set "proper" number of data points (delta and qMax?) 
    10221028// X- need to remove points at high Q end 
    10231029// 
     
    10271033// the results are in iBin_qxqy, qBin_qxqy, and eBin_qxqy, in the folder passed 
    10281034//  
     1035 
     1036Proc BinQxQy_to_1D(folderStr) 
     1037        String folderStr 
     1038        Prompt folderStr,"Pick the data folder containing 2D data",popup,getAList(4) 
     1039 
     1040        fDoBinning_QxQy2D(folderStr) 
     1041End 
     1042 
     1043 
     1044 
    10291045//Function fDoBinning_QxQy2D(inten,qx,qy,qz) 
     1046//      Wave inten,qx,qy,qz 
     1047 
    10301048Function fDoBinning_QxQy2D(folderStr) 
    10311049        String folderStr 
    10321050 
    1033 //      Wave inten,qx,qy,qz 
    10341051 
    10351052        SetDataFolder $("root:"+folderStr) 
    10361053         
    1037         WAVE inten = $(folderStr + "_i") 
     1054        WAVE inten = $("smeared_sf2D") 
     1055 
     1056//      WAVE inten = $(folderStr + "_i") 
     1057        WAVE iErr = $(folderStr + "_iErr") 
    10381058        WAVE qx = $(folderStr + "_qx") 
    10391059        WAVE qy = $(folderStr + "_qy") 
     
    10481068         
    10491069        yDim = XDim 
    1050         Make/O/D/N=(nq) iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy 
     1070        Make/O/D/N=(nq) iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy,eBin2D_qxqy 
    10511071        delQ = abs(sqrt(qx[2]^2+qy[2]^2+qz[2]^2) - sqrt(qx[1]^2+qy[1]^2+qz[1]^2))               //use bins of 1 pixel width  
    10521072        qBin_qxqy[] =  p*       delQ     
     
    10561076        iBin2_qxqy = 0 
    10571077        eBin_qxqy = 0 
     1078        eBin2D_qxqy = 0 
    10581079        nBin_qxqy = 0   //number of intensities added to each bin 
    10591080         
     
    10651086                        iBin_qxqy[binIndex] += val 
    10661087                        iBin2_qxqy[binIndex] += val*val 
     1088                        eBin2D_qxqy[binIndex] += iErr[ii]*iErr[ii] 
    10671089                        nBin_qxqy[binIndex] += 1 
    10681090                endif 
     
    10751097                        iBin_qxqy[ii] = 0 
    10761098                        eBin_qxqy[ii] = 1 
     1099                        eBin2D_qxqy[ii] = NaN 
    10771100                else 
    10781101                        if(nBin_qxqy[ii] <= 1) 
     
    10801103                                iBin_qxqy[ii] /= nBin_qxqy[ii] 
    10811104                                eBin_qxqy[ii] = 1 
     1105                                eBin2D_qxqy[ii] /= (nBin_qxqy[ii])^2 
    10821106                        else 
    10831107                                //assume that the intensity in each pixel in annuli is normally 
     
    10921116                                        eBin_qxqy[ii] = sqrt(var/(nBin_qxqy[ii] - 1)) 
    10931117                                endif 
     1118                                // and calculate as it is propagated pixel-by-pixel 
     1119                                eBin2D_qxqy[ii] /= (nBin_qxqy[ii])^2 
    10941120                        endif 
    10951121                endif 
    10961122        endfor 
     1123         
     1124        eBin2D_qxqy = sqrt(eBin2D_qxqy)         // as equation (3) of John's memo 
    10971125         
    10981126        // find the last non-zero point, working backwards 
     
    11031131         
    11041132//      print val, nBin_qxqy[val] 
    1105         DeletePoints val, nq-val, iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy 
     1133        DeletePoints val, nq-val, iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy,eBin2D_qxqy 
    11061134         
    11071135        SetDataFolder root: 
Note: See TracChangeset for help on using the changeset viewer.