Ignore:
Timestamp:
Jul 3, 2018 2:30:41 PM (4 years ago)
Author:
srkline
Message:

Many changes:

Made the VCALC panel aware of all of the binning options
Corrected the behavior of the VCALC preset conditions
Adjusted how the Slit data is binned so that there are not duplicated q-values in the output

Made Absolute scaling aware of the back detector. Now the ABS String in the protocol has a second
set of scaling constants tagged with "_B" for the back detector. There is an added button
on the protocol panel to set the second set of constants. For the back detector, the read noise
is subtracted by reading it from the empty beam file (shifting over to the right by one box width)
All of the associated abs procedures are now aware of this.
More error checking needs to be added.

Back detector image is now shifted upon loading of the data. the default mask takes this into account
and masks out the padded (zero) regions.

in the protocol, DIV and MSK do not use grep any longer. it was just way too slow. Now it depends on

the file name having DIV or MASK respectively.



Raw data files can now be added together, in the usual way from the protocol panel.



File:
1 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_DetectorCorrections.ipf

    r1101 r1109  
    965965 
    966966 
    967 //////////// 
    968 // TODO: all of below is untested code 
    969 //   copied from SANS 
    970 // 
    971 // 
    972 // NOV 2017 
    973 // Currently, this is not called from any VSANS routines. it is only referenced 
    974 // from V_Add_raw_to_work(), which would add two VSANS raw data files together. This has 
    975 // not yet been implemented. I am only keeping this function around to be sure that  
    976 // if/when V_Add_raw_to_work() is implemented, all of the functionality of V_DetCorr() is 
    977 // properly duplicated. 
    978 // 
    979 // 
    980 // 
    981 //performs solid angle and non-linear detector corrections to raw data as it is "added" to a work folder 
    982 //function is called by Raw_to_work() and Add_raw_to_work() functions 
    983 //works on the actual data array, assumes that is is already on LINEAR scale 
    984 // 
    985 Function V_DetCorr(data,data_err,realsread,doEfficiency,doTrans) 
    986         Wave data,data_err,realsread 
    987         Variable doEfficiency,doTrans 
    988  
    989         DoAlert 0,"This has not yet been updated for VSANS" 
    990          
    991         Variable xcenter,ycenter,x0,y0,sx,sx3,sy,sy3,xx0,yy0 
    992         Variable ii,jj,dtdist,dtdis2 
    993         Variable xi,xd,yd,rad,ratio,domega,xy 
    994         Variable lambda,trans,trans_err,lat_err,tmp_err,lat_corr 
    995          
    996 //      Print "...doing jacobian and non-linear corrections" 
    997  
    998         NVAR pixelsX = root:myGlobals:gNPixelsX 
    999         NVAR pixelsY = root:myGlobals:gNPixelsY 
    1000          
    1001         //set up values to send to auxiliary trig functions 
    1002         xcenter = pixelsX/2 + 0.5               // == 64.5 for 128x128 Ordela 
    1003         ycenter = pixelsY/2 + 0.5               // == 64.5 for 128x128 Ordela 
    1004  
    1005         x0 = realsread[16] 
    1006         y0 = realsread[17] 
    1007         sx = realsread[10] 
    1008         sx3 = realsread[11] 
    1009         sy = realsread[13] 
    1010         sy3 = realsread[14] 
    1011          
    1012         dtdist = 1000*realsread[18]     //sdd in mm 
    1013         dtdis2 = dtdist^2 
    1014          
    1015         lambda = realsRead[26] 
    1016         trans = RealsRead[4] 
    1017         trans_err = RealsRead[41]               //new, March 2011 
    1018          
    1019  
    1020         //waves to contain repeated function calls 
    1021         Make/O/N=(pixelsX) fyy,xx,yy            //Assumes square detector !!! 
    1022         ii=0 
    1023         do 
    1024                 xi = ii 
    1025 //              fyy[ii] = dc_fy(ii+1,sy,sy3,ycenter) 
    1026 //              xx[ii] = dc_fxn(ii+1,sx,sx3,xcenter) 
    1027 //              yy[ii] = dc_fym(ii+1,sy,sy3,ycenter) 
    1028                 ii+=1 
    1029         while(ii<pixelsX) 
    1030          
    1031         Make/O/N=(pixelsX,pixelsY) SolidAngle           // testing only 
    1032          
    1033         ii=0 
    1034         do 
    1035                 xi = ii 
    1036 //              xd = dc_fx(ii+1,sx,sx3,xcenter)-xx0 
    1037                 jj=0 
    1038                 do 
    1039                         yd = fyy[jj]-yy0 
    1040                         //rad is the distance of pixel ij from the sample 
    1041                         //domega is the ratio of the solid angle of pixel ij versus center pixel 
    1042                         // product xy = 1 for a detector with a linear spatial response (modern Ordela) 
    1043                         // solid angle calculated, dW^3 >=1, so multiply data to raise measured values to correct values. 
    1044                         rad = sqrt(dtdis2 + xd^2 + yd^2) 
    1045                         domega = rad/dtdist 
    1046                         ratio = domega^3 
    1047                         xy = xx[ii]*yy[jj] 
    1048                          
    1049                         data[ii][jj] *= xy*ratio 
    1050                          
    1051                         solidAngle[ii][jj] = xy*ratio           //testing only   
    1052                         data_err[ii][jj] *= xy*ratio                    //error propagation assumes that SA and Jacobian are exact, so simply scale error 
    1053                          
    1054                          
    1055                         // correction factor for detector efficiency JBG memo det_eff_cor2.doc 3/20/07 
    1056                         // correction inserted 11/2007 SRK 
    1057                         // large angle detector efficiency is >= 1 and will "bump up" the measured value at the highest angles 
    1058                         // so divide here to get the correct answer (5/22/08 SRK) 
    1059                         if(doEfficiency) 
    1060 //                              data[ii][jj] /= DetEffCorr(lambda,dtdist,xd,yd) 
    1061 //                              data_err[ii][jj] /= DetEffCorr(lambda,dtdist,xd,yd) 
    1062 //                              solidAngle[ii][jj] /= DetEffCorr(lambda,dtdist,xd,yd)           //testing only 
    1063                         endif 
    1064                          
    1065                         // large angle transmission calculation is <= 1 and will "bump down" the measured value at the highest angles 
    1066                         // so divide here to get the correct answer 
    1067                         if(doTrans) 
    1068                          
    1069                                 if(trans<0.1 && ii==0 && jj==0) 
    1070                                         Print "***transmission is less than 0.1*** and is a significant correction" 
    1071                                 endif 
    1072                                  
    1073                                 if(trans==0) 
    1074                                         if(ii==0 && jj==0) 
    1075                                                 Print "***transmission is ZERO*** and has been reset to 1.0 for the averaging calculation" 
    1076                                         endif 
    1077                                         trans = 1 
    1078                                 endif 
    1079                                  
    1080                                 // pass in the transmission error, and the error in the correction is returned as the last parameter 
    1081  
    1082 //                              lat_corr = V_LargeAngleTransmissionCorr(trans,dtdist,xd,yd,trans_err,lat_err)           //moved from 1D avg SRK 11/2007 
    1083  
    1084                                 data[ii][jj] /= lat_corr                        //divide by the correction factor 
    1085                                 // 
    1086                                 // 
    1087                                 // 
    1088                                 // relative errors add in quadrature 
    1089                                 tmp_err = (data_err[ii][jj]/lat_corr)^2 + (lat_err/lat_corr)^2*data[ii][jj]*data[ii][jj]/lat_corr^2 
    1090                                 tmp_err = sqrt(tmp_err) 
    1091                                  
    1092                                 data_err[ii][jj] = tmp_err 
    1093                                  
    1094 //                              solidAngle[ii][jj] = lat_err 
    1095  
    1096                                  
    1097                                 //solidAngle[ii][jj] = LargeAngleTransmissionCorr(trans,dtdist,xd,yd)           //testing only 
    1098                         endif 
    1099                          
    1100                         jj+=1 
    1101                 while(jj<pixelsX) 
    1102                 ii+=1 
    1103         while(ii<pixelsX) 
    1104          
    1105         //clean up waves 
    1106          
    1107         Return(0) 
    1108 End 
    1109967 
    1110968 
     
    12591117//w_ is the "work" file 
    12601118//both are work files and should already be normalized to 10^8 monitor counts 
    1261 Function V_Absolute_Scale(type,w_trans,w_thick,s_trans,s_thick,s_izero,s_cross,kappa_err) 
    1262         String type 
     1119Function V_Absolute_Scale(type,absStr) 
     1120        String type,absStr 
     1121         
     1122         
    12631123        Variable w_trans,w_thick,s_trans,s_thick,s_izero,s_cross,kappa_err 
    1264  
    12651124 
    12661125        Variable defmon = 1e8,w_moncount,s1,s2,s3,s4 
     
    12831142         
    12841143        w_moncount = V_getBeamMonNormData(type) 
    1285          
    1286          
     1144                 
    12871145        if(w_moncount == 0) 
    12881146                //zero monitor counts will give divide by zero --- 
     
    12901148                Return(1)               //report error 
    12911149        Endif 
     1150 
     1151        w_trans = V_getSampleTransmission(type)         //sample transmission 
     1152        w_thick = V_getSampleThickness(type)            //sample thickness 
     1153        trans_err = V_getSampleTransError(type)  
     1154         
     1155         
     1156        //get the parames from the list 
     1157        s_trans = NumberByKey("TSTAND", absStr, "=", ";")       //parse the list of values 
     1158        s_thick = NumberByKey("DSTAND", absStr, "=", ";") 
     1159        s_izero = NumberByKey("IZERO", absStr, "=", ";") 
     1160        s_cross = NumberByKey("XSECT", absStr, "=", ";") 
     1161        kappa_err = NumberByKey("SDEV", absStr, "=", ";") 
     1162 
    12921163         
    12931164        //calculate scale factor 
     
    12981169        scale = s1*s2*s3*s4 
    12991170 
    1300         trans_err = V_getSampleTransError(type)  
    13011171         
    13021172        // kappa comes in as s_izero, so be sure to use 1/kappa_err 
     
    13041174        // and now loop through all of the detectors 
    13051175        //do the actual absolute scaling here, modifying the data in ABS 
    1306         for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
    1307                 detStr = StringFromList(ii, ksDetectorListAll, ";") 
     1176        for(ii=0;ii<ItemsInList(ksDetectorListNoB);ii+=1) 
     1177                detStr = StringFromList(ii, ksDetectorListNoB, ";") 
    13081178                Wave data = V_getDetectorDataW("ABS",detStr) 
    13091179                Wave data_err = V_getDetectorDataErrW("ABS",detStr) 
     
    13121182                data_err = sqrt(scale^2*data_err^2 + scale^2*data^2*(kappa_err^2/s_izero^2 +trans_err^2/w_trans^2)) 
    13131183        endfor 
     1184 
     1185        // do the back detector separately, if it is set to be used 
     1186        NVAR gIgnoreDetB = root:Packages:NIST:VSANS:Globals:gIgnoreDetB 
     1187        if(gIgnoreDetB == 0) 
     1188                detStr = "B" 
     1189                Wave data = V_getDetectorDataW("ABS",detStr) 
     1190                Wave data_err = V_getDetectorDataErrW("ABS",detStr) 
     1191                 
     1192                //get the parames from the list 
     1193                s_trans = NumberByKey("TSTAND_B", absStr, "=", ";")     //parse the list of values 
     1194                s_thick = NumberByKey("DSTAND_B", absStr, "=", ";") 
     1195                s_izero = NumberByKey("IZERO_B", absStr, "=", ";") 
     1196                s_cross = NumberByKey("XSECT_B", absStr, "=", ";") 
     1197                kappa_err = NumberByKey("SDEV_B", absStr, "=", ";") 
     1198 
     1199                //calculate scale factor 
     1200                s1 = defmon/w_moncount          // monitor count (s1 should be 1) 
     1201                s2 = s_thick/w_thick 
     1202                s3 = s_trans/w_trans 
     1203                s4 = s_cross/s_izero 
     1204                scale = s1*s2*s3*s4 
     1205                 
     1206                data *= scale 
     1207                data_err = sqrt(scale^2*data_err^2 + scale^2*data^2*(kappa_err^2/s_izero^2 +trans_err^2/w_trans^2)) 
     1208        endif 
    13141209         
    13151210        //********* 15APR02 
     
    14761371 
    14771372////////////////////////// 
     1373// detector corrections to stitch the back detector into one proper image 
     1374// 
     1375// 
     1376// 
     1377 
     1378 
     1379// 
     1380// to register the image on the back detector panel 
     1381// 
     1382// middle portion (552 pix in Y) is held fixed 
     1383// top portion of image is shifted right and down 
     1384// bottom portion of image is shifted right and up 
     1385// 
     1386// remainder of image is filled with Zero (NaN causes problems converting to WORK) 
     1387// 
     1388// currently, data is not added together and averaged, but it could be 
     1389// 
     1390Function V_ShiftBackDetImage(w,adjW) 
     1391        Wave w,adjW 
     1392         
     1393        adjW=0 
     1394         
     1395        Variable topX,bottomX 
     1396        Variable topY,bottomY 
     1397         
     1398        topX = 7 
     1399        topY = 105 
     1400         
     1401        bottomX = 7 
     1402        bottomY = 35 
     1403         
     1404        // middle 
     1405        adjW[][552,552+552] = w[p][q] 
     1406 
     1407        //top 
     1408        adjW[0+topX,679][552+552,1655-topY] = w[p-topX][q+topY] 
     1409         
     1410        //bottom 
     1411        adjW[0+bottomX,679][0+bottomY,551] = w[p-bottomX][q-bottomY] 
     1412         
     1413        return(0) 
     1414End 
     1415 
     1416 
Note: See TracChangeset for help on using the changeset viewer.