Ignore:
Timestamp:
Jan 19, 2017 7:53:30 AM (6 years ago)
Author:
srkline
Message:

changes to V_Correct.ipf to account for VSANS (9) detector banks.

File:
1 edited

Legend:

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

    r1019 r1020  
    1 #pragma rtGlobals=1             // Use modern global access method. 
    2 #pragma version=5.0 
    3 #pragma IgorVersion=6.1 
     1#pragma rtGlobals=3             // Use modern global access method and strict wave access. 
     2//#pragma version=5.0 
     3//#pragma IgorVersion=6.1 
     4 
     5// Updated for VSANS Jan2017 
     6// 
     7// largely duplication of the SANS subtractions and error propagation. 
     8// Changes: (1) the beam center mismatch is ALWAYS ignored. It is flagged, and alerted, but nothing is shifted 
     9//          (2) the condition of trans == 1 is not flagged, and there is no stopping for user input 
     10//  
     11// TODO -- verify the operation of all modes 
     12// -- decide if/how to implement/re-implement the trans == 1 check and dialog 
     13// -- decide if the beam center mismatch is ever to be re-implemented 
     14// -- check the monitor count calls and rescaled values (correct monitor? where is rescaling written?) 
     15// 
     16 
     17 
    418 
    519 
     
    2943 
    3044// 
    31 //unused test procedure for Correct() function 
     45// test procedure for Correct() function 
    3246//must be updated to include "mode" parameter before re-use 
    3347// 
     
    139153                                return(err) 
    140154                        Endif 
    141 // TODO                 err = V_CorrectMode_1() 
     155                        err = V_CorrectMode_1() 
    142156                        break 
    143157                case 2: 
     
    146160                                return(err) 
    147161                        Endif 
    148 // TODO                 err = V_CorrectMode_2() 
     162                        err = V_CorrectMode_2() 
    149163                        break 
    150164                case 3: 
     
    165179//                      endif 
    166180 
    167 // TODO                 err = V_CorrectMode_3() 
     181                        err = V_CorrectMode_3() 
    168182                        break 
    169183                case 4: 
     
    194208                                return(err) 
    195209                        Endif 
    196 // TODO                 err = V_CorrectMode_11() 
     210                        err = V_CorrectMode_11() 
    197211                        break 
    198212                case 12: 
     
    205219                                return(err) 
    206220                        Endif 
    207 // TODO                 err = V_CorrectMode_12() 
     221                        err = V_CorrectMode_12() 
    208222                        break 
    209223                case 13: 
     
    227241                                return(err) 
    228242                        Endif 
    229 // TODO                 err = V_CorrectMode_13() 
     243                        err = V_CorrectMode_13() 
    230244                        break 
    231245                case 14: 
     
    234248                                return(err) 
    235249                        Endif 
    236 // TODO                 err = V_CorrectMode_14() 
     250                        err = V_CorrectMode_14() 
    237251                        break 
    238252                default:        //something wrong 
     
    253267//                                      was freshly loaded. added final copy of cor result to cor:data and cor:linear_data 
    254268// 
    255 xFunction V_CorrectMode_1() 
    256          
    257         //create the necessary wave references 
    258         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    259         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    260         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    261         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    262         WAVE bgd_data=$"root:Packages:NIST:BGD:linear_data" 
    263         WAVE bgd_reals=$"root:Packages:NIST:BGD:realsread" 
    264         WAVE bgd_ints=$"root:Packages:NIST:BGD:integersread" 
    265         WAVE/T bgd_text=$"root:Packages:NIST:BGD:textread" 
    266         WAVE emp_data=$"root:Packages:NIST:EMP:linear_data" 
    267         WAVE emp_reals=$"root:Packages:NIST:EMP:realsread" 
    268         WAVE emp_ints=$"root:Packages:NIST:EMP:integersread" 
    269         WAVE/T emp_text=$"root:Packages:NIST:EMP:textread" 
    270         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    271         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    272          
    273         // needed to propagate error 
    274         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    275         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    276         WAVE bgd_err =$"root:Packages:NIST:BGD:linear_data_error" 
    277         WAVE emp_err =$"root:Packages:NIST:EMP:linear_data_error" 
    278         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    279          
    280         Variable sam_trans_err,emp_trans_err 
    281         sam_trans_err = sam_reals[41] 
    282         emp_trans_err = emp_reals[41] 
    283          
    284          
    285         //get sam and bgd attenuation factors 
    286         String fileStr="" 
    287         Variable lambda,attenNo,sam_AttenFactor,bgd_attenFactor,emp_AttenFactor 
     269Function V_CorrectMode_1() 
     270 
     271        //get SAM, BGD, EMP attenuation factor 
     272        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     273        Variable bgd_AttenFactor,bgd_atten_err 
     274        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     275        Variable ii 
     276        String detStr 
    288277        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
    289278        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
    290         Variable sam_atten_err,emp_atten_err,bgd_atten_err 
    291         fileStr = sam_text[3] 
    292         lambda = sam_reals[26] 
    293         attenNo = sam_reals[3] 
    294         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    295         fileStr = bgd_text[3] 
    296         lambda = bgd_reals[26] 
    297         attenNo = bgd_reals[3] 
    298         bgd_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,bgd_atten_err) 
    299         fileStr = emp_text[3] 
    300         lambda = emp_reals[26] 
    301         attenNo = emp_reals[3] 
    302         emp_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,emp_atten_err) 
    303          
    304         //get relative monitor counts (should all be 10^8, since normalized in add step) 
    305         tmonsam = sam_reals[0]          //monitor count in SAM 
    306         tsam = sam_reals[4]             //SAM transmission 
    307         csam = sam_reals[16]            //x center 
    308         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    309         tmonbgd = bgd_reals[0]          //monitor count in BGD 
    310         cbgd = bgd_reals[16] 
    311         rbgd = bgd_reals[17] 
    312         tmonemp = emp_reals[0]          //monitor count in EMP 
    313         temp = emp_reals[4]                     //trans emp 
    314         cemp = emp_reals[16]            //beamcenter of EMP 
    315         remp = emp_reals[17] 
    316          
    317         if(temp==0) 
    318                 DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
    319                 temp=1 
    320         Endif 
    321          
    322         NVAR pixelsX = root:myGlobals:gNPixelsX 
    323         NVAR pixelsY = root:myGlobals:gNPixelsY 
    324          
    325         //get the shifted data arrays, EMP and BGD, each relative to SAM 
    326         Make/D/O/N=(pixelsX,pixelsY) cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
    327         xshift = cbgd-csam 
    328         yshift = rbgd-rsam 
    329         if(abs(xshift) <= wcen) 
    330                 xshift = 0 
    331         Endif 
    332         if(abs(yshift) <= wcen) 
    333                 yshift = 0 
    334         Endif 
    335         GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)              //bgd_temp 
    336          
    337         xshift = cemp-csam 
    338         yshift = remp-rsam 
    339         if(abs(xshift) <= wcen) 
    340                 xshift = 0 
    341         Endif 
    342         if(abs(yshift) <= wcen) 
    343                 yshift = 0 
    344         Endif 
    345         GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)              //emp_temp 
    346  
    347         //do the subtraction 
    348         fsam=1 
    349         femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
    350         fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
    351         cor1 = fsam*sam_data/sam_attenFactor - fbgd*bgd_temp/bgd_attenFactor 
    352         cor1 -= (tsam/temp)*(femp*emp_temp/emp_attenFactor - fbgd*bgd_temp/bgd_attenFactor) 
    353         cor1 *= noadd_bgd*noadd_emp             //zero out the array mismatch values 
    354  
    355 // do the error propagation piecewise    
    356         Duplicate/O sam_err, tmp_a, tmp_b, tmp_c, tmp_d,c_val,d_val 
    357         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    358          
    359         tmp_b = (bgd_err/bgd_attenFactor)^2*(tsam/temp - 1)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2*(1-tsam/temp)^2            //sig b ^2 
    360  
    361         tmp_c = (sam_trans_err/temp)^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
    362         tmp_c += (tsam/temp^2)^2*emp_trans_err^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
    363          
    364         tmp_d = (tsam/(temp*emp_attenFactor))^2*(emp_err)^2 + (tsam*emp_data/(temp*emp_attenFactor^2))^2*(emp_atten_err)^2 
    365  
    366         cor_err = sqrt(tmp_a + tmp_b + tmp_c + tmp_d) 
    367          
    368         //we're done, get out w/no error 
    369         //set the COR data and linear_data to the result 
    370         cor_data = cor1 
    371         cor_data_display = cor1 
    372          
    373         //update COR header 
    374         cor_text[1] = date() + " " + time()             //date + time stamp 
    375  
    376         KillWaves/Z cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
    377         Killwaves/Z tmp_a,tmp_b,tmp_c,tmp_d,c_val,d_val 
    378         SetDataFolder root: 
    379         Return(0) 
    380 End 
    381  
    382 //background only 
    383 // existence of data checked by dispatching routine 
    384 // data has already been copied to COR folder 
    385 xFunction V_CorrectMode_2() 
    386  
    387         //create the necessary wave references 
    388         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    389         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    390         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    391         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    392         WAVE bgd_data=$"root:Packages:NIST:BGD:linear_data" 
    393         WAVE bgd_reals=$"root:Packages:NIST:BGD:realsread" 
    394         WAVE bgd_ints=$"root:Packages:NIST:BGD:integersread" 
    395         WAVE/T bgd_text=$"root:Packages:NIST:BGD:textread" 
    396         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    397         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    398  
    399         // needed to propagate error 
    400         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    401         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    402         WAVE bgd_err =$"root:Packages:NIST:BGD:linear_data_error" 
    403         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    404          
    405         Variable sam_trans_err 
    406         sam_trans_err = sam_reals[41] 
    407  
    408          
    409         //get sam and bgd attenuation factors 
    410         String fileStr="" 
    411         Variable lambda,attenNo,sam_AttenFactor,bgd_attenFactor 
    412         Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
    413         Variable wcen=0.001 
    414         Variable sam_atten_err,bgd_atten_err 
    415         fileStr = sam_text[3] 
    416         lambda = sam_reals[26] 
    417         attenNo = sam_reals[3] 
    418         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    419         fileStr = bgd_text[3] 
    420         lambda = bgd_reals[26] 
    421         attenNo = bgd_reals[3] 
    422         bgd_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,bgd_atten_err) 
    423          
    424         //Print "atten = ",sam_attenFactor,bgd_attenFactor 
    425          
    426         //get relative monitor counts (should all be 10^8, since normalized in add step) 
    427         tmonsam = sam_reals[0]          //monitor count in SAM 
    428         csam = sam_reals[16]            //x center 
    429         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    430         tmonbgd = bgd_reals[0]          //monitor count in BGD 
    431         cbgd = bgd_reals[16] 
    432         rbgd = bgd_reals[17] 
    433  
    434         // set up beamcenter shift, relative to SAM 
    435         xshift = cbgd-csam 
    436         yshift = rbgd-rsam 
    437         if(abs(xshift) <= wcen) 
    438                 xshift = 0 
    439         Endif 
    440         if(abs(yshift) <= wcen) 
    441                 yshift = 0 
    442         Endif 
    443          
    444         NVAR pixelsX = root:myGlobals:gNPixelsX 
    445         NVAR pixelsY = root:myGlobals:gNPixelsY 
    446         //get shifted data arrays, relative to SAM 
    447         Make/D/O/N=(pixelsX,pixelsY) cor1,bgd_temp,noadd_bgd            //temp arrays 
    448         GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)              //bgd_temp is the BGD  
    449          
    450         //do the sam-bgd subtraction,  deposit result in cor1 
    451         fsam = 1 
    452         fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
    453          
    454         //print "fsam,fbgd = ",fsam,fbgd 
    455          
    456         cor1 = fsam*sam_data/sam_AttenFactor - fbgd*bgd_temp/bgd_AttenFactor 
    457         cor1 *= noadd_bgd               //zeros out regions where arrays do not overlap, one otherwise 
    458  
    459 // do the error propagation piecewise    
    460         Duplicate/O sam_err, tmp_a, tmp_b 
    461         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    462          
    463         tmp_b = (bgd_err/bgd_attenFactor)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2              //sig b ^2 
    464  
    465         cor_err = sqrt(tmp_a + tmp_b) 
    466  
    467  
    468         //we're done, get out w/no error 
    469         //set the COR_data to the result 
    470         cor_data = cor1 
    471         cor_data_display = cor1 
    472  
    473         //update COR header 
    474         cor_text[1] = date() + " " + time()             //date + time stamp 
    475  
    476         KillWaves/Z cor1,bgd_temp,noadd_bgd 
    477         Killwaves/Z tmp_a,tmp_b 
    478  
    479         SetDataFolder root: 
    480         Return(0) 
    481 End 
    482  
    483 // empty subtraction only 
    484 // data does exist, checked by dispatch routine 
    485 // 
    486 xFunction V_CorrectMode_3() 
    487         //create the necessary wave references 
    488         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    489         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    490         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    491         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    492         WAVE emp_data=$"root:Packages:NIST:EMP:linear_data" 
    493         WAVE emp_reals=$"root:Packages:NIST:EMP:realsread" 
    494         WAVE emp_ints=$"root:Packages:NIST:EMP:integersread" 
    495         WAVE/T emp_text=$"root:Packages:NIST:EMP:textread" 
    496         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    497         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    498          
    499         // needed to propagate error 
    500         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    501         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    502         WAVE emp_err =$"root:Packages:NIST:EMP:linear_data_error" 
    503         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    504          
    505         Variable sam_trans_err,emp_trans_err 
    506         sam_trans_err = sam_reals[41] 
    507         emp_trans_err = emp_reals[41]    
    508          
    509         //get sam and bgd attenuation factors 
    510         String fileStr="" 
    511         Variable lambda,attenNo,sam_AttenFactor,emp_attenFactor 
    512         Variable tmonsam,fsam,femp,xshift,yshift,rsam,csam,remp,cemp,tmonemp 
    513         Variable wcen=0.001,tsam,temp 
    514         Variable sam_atten_err,emp_atten_err 
    515         fileStr = sam_text[3] 
    516         lambda = sam_reals[26] 
    517         attenNo = sam_reals[3] 
    518         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    519         fileStr = emp_text[3] 
    520         lambda = emp_reals[26] 
    521         attenNo = emp_reals[3] 
    522         emp_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,emp_atten_err) 
    523          
    524         //get relative monitor counts (should all be 10^8, since normalized in add step) 
    525         tmonsam = sam_reals[0]          //monitor count in SAM 
    526         tsam = sam_reals[4]             //SAM transmission 
    527         csam = sam_reals[16]            //x center 
    528         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    529         tmonemp = emp_reals[0]          //monitor count in EMP 
    530         temp = emp_reals[4]                     //trans emp 
    531         cemp = emp_reals[16]            //beamcenter of EMP 
    532         remp = emp_reals[17] 
    533          
    534         if(temp==0) 
    535                 DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
    536                 temp=1 
    537         Endif 
    538          
    539         //Print "rbgd,cbgd = ",rbgd,cbgd 
    540         // set up beamcenter shift, relative to SAM 
    541         xshift = cemp-csam 
    542         yshift = remp-rsam 
    543         if(abs(xshift) <= wcen) 
    544                 xshift = 0 
    545         Endif 
    546         if(abs(yshift) <= wcen) 
    547                 yshift = 0 
    548         Endif 
    549          
    550         NVAR pixelsX = root:myGlobals:gNPixelsX 
    551         NVAR pixelsY = root:myGlobals:gNPixelsY 
    552         //get shifted data arrays, relative to SAM 
    553         Make/D/O/N=(pixelsX,pixelsY) cor1,emp_temp,noadd_emp            //temp arrays 
    554         GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)              //emp_temp is the EMP 
    555          
    556         //do the sam-bgd subtraction,  deposit result in cor1 
    557         fsam = 1 
    558         femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
    559          
    560         cor1 = fsam*sam_data/sam_AttenFactor - femp*(tsam/temp)*emp_temp/emp_AttenFactor 
    561         cor1 *= noadd_emp               //zeros out regions where arrays do not overlap, one otherwise 
    562  
    563 // do the error propagation piecewise    
    564         Duplicate/O sam_err, tmp_a, tmp_c ,c_val 
    565         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    566          
    567         tmp_c = (sam_trans_err*emp_data/(temp*emp_attenFactor))^2 + (emp_err*tsam/(temp*emp_attenFactor))^2 
    568         tmp_c += (tsam*emp_data*emp_trans_err/(temp*temp*emp_attenFactor))^2 + (tsam*emp_data*emp_atten_err/(temp*emp_attenFactor^2))^2//total of 6 terms 
    569  
    570         cor_err = sqrt(tmp_a + tmp_c) 
    571          
    572         //we're done, get out w/no error 
    573         //set the COR data to the result 
    574         cor_data = cor1 
    575         cor_data_display = cor1 
    576  
    577         //update COR header 
    578         cor_text[1] = date() + " " + time()             //date + time stamp 
    579  
    580         KillWaves/Z cor1,emp_temp,noadd_emp 
    581         Killwaves/Z tmp_a,tmp_c,c_val 
    582  
    583         SetDataFolder root: 
    584         Return(0) 
    585 End 
    586  
    587 // NO subtraction - simply rescales for attenuators 
    588 // SAM data does exist, checked by dispatch routine 
    589 // SAM data has already been copied to COR (both are the same at the start of the function) 
    590 // 
    591 // 
    592 // 
    593 Function V_CorrectMode_4() 
    594  
    595         //get SAM attenuation factor 
    596         Variable sam_AttenFactor,sam_atten_err,ii 
    597         String detStr 
    598          
     279         
     280        // these values apply to all of the detectors 
    599281        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
    600282        sam_atten_err = V_getAttenuator_trans_err("SAM") 
    601  
    602  
     283        bgd_AttenFactor = V_getAttenuator_transmission("BGD") 
     284        bgd_atten_err = V_getAttenuator_trans_err("BGD") 
     285        emp_AttenFactor = V_getAttenuator_transmission("EMP") 
     286        emp_atten_err = V_getAttenuator_trans_err("EMP") 
     287 
     288        //get relative monitor counts (should all be 10^8, since normalized in add step) 
     289        // get transmission and trans error for SAM, EMP 
     290        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     291 
     292        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     293        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     294        sam_trans_err = V_getSampleTransError("SAM") 
     295         
     296        tmonemp = V_getMonitorCount("EMP")              //monitor count in EMP 
     297        temp = V_getSampleTransmission("EMP")                   //trans emp 
     298        emp_trans_err = V_getSampleTransError("EMP") 
     299         
     300        tmonbgd = V_getMonitorCount("BGD")              //monitor count in BGD 
     301 
     302 
     303        // and now loop through all of the detectors 
    603304        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
    604305                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     
    607308                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
    608309                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     310                Wave bgd_data = V_getDetectorDataW("BGD",detStr) 
     311                Wave bgd_err = V_getDetectorDataErrW("BGD",detStr) 
     312                Wave emp_data = V_getDetectorDataW("EMP",detStr) 
     313                Wave emp_err = V_getDetectorDataErrW("EMP",detStr) 
     314                 
     315        // to check for beam center mismatch -- simply warn, but do no shift 
     316        // 
     317 
     318                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     319                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     320         
     321                cbgd = V_getDet_beam_center_x("BGD",detStr) 
     322                rbgd = V_getDet_beam_center_y("BGD",detStr) 
     323         
     324                cemp = V_getDet_beam_center_x("EMP",detStr)             //beamcenter of EMP 
     325                remp = V_getDet_beam_center_y("EMP",detStr) 
     326                 
     327                if(temp==0) 
     328                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     329                        temp=1 
     330                Endif 
     331         
     332         
     333                Duplicate/O cor_data cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
     334 
     335                // TODO -- document this, make a note, so everyone knows this is not done 
     336                // skip this part, but duplicate the results of no shift condition 
     337                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     338                 
     339                        //get the shifted data arrays, EMP and BGD, each relative to SAM 
     340         
     341                xshift = cbgd-csam 
     342                yshift = rbgd-rsam 
     343                if(abs(xshift) <= wcen) 
     344                        xshift = 0 
     345                Endif 
     346                if(abs(yshift) <= wcen) 
     347                        yshift = 0 
     348                Endif 
     349                // for the BGD file - alert if needed, generate dummy "pass-through" values 
     350                // 
     351                if(xshift != 0 || yshift != 0) 
     352                        DoAlert 0, "Beam center mismatch for BGD file. Data has NOT been corrected." 
     353                endif 
     354                bgd_temp = bgd_data             // no shift, no effect 
     355                noadd_bgd = 1 
     356                //GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)            //bgd_temp 
     357         
     358                xshift = cemp-csam 
     359                yshift = remp-rsam 
     360                if(abs(xshift) <= wcen) 
     361                        xshift = 0 
     362                Endif 
     363                if(abs(yshift) <= wcen) 
     364                        yshift = 0 
     365                Endif 
     366                // for the EMP file - alert if needed, generate dummy "pass-through" values 
     367                // 
     368                if(xshift != 0 || yshift != 0) 
     369                        DoAlert 0, "Beam center mismatch for EMP file. Data has NOT been corrected." 
     370                endif 
     371                emp_temp = emp_data // no shift, no effect 
     372                noadd_emp = 1 
     373                //GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)            //emp_temp 
     374 
     375 
     376                // ******* 
     377                //do the subtraction 
     378                fsam=1 
     379                femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
     380                fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
     381                cor1 = fsam*sam_data/sam_attenFactor - fbgd*bgd_temp/bgd_attenFactor 
     382                cor1 -= (tsam/temp)*(femp*emp_temp/emp_attenFactor - fbgd*bgd_temp/bgd_attenFactor) 
     383                cor1 *= noadd_bgd*noadd_emp             //zero out the array mismatch values 
     384         
     385                // do the error propagation piecewise    
     386                Duplicate/O sam_err, tmp_a, tmp_b, tmp_c, tmp_d,c_val,d_val 
     387                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     388                 
     389                tmp_b = (bgd_err/bgd_attenFactor)^2*(tsam/temp - 1)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2*(1-tsam/temp)^2            //sig b ^2 
     390         
     391                tmp_c = (sam_trans_err/temp)^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
     392                tmp_c += (tsam/temp^2)^2*emp_trans_err^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
     393                 
     394                tmp_d = (tsam/(temp*emp_attenFactor))^2*(emp_err)^2 + (tsam*emp_data/(temp*emp_attenFactor^2))^2*(emp_atten_err)^2 
     395         
     396                cor_err = sqrt(tmp_a + tmp_b + tmp_c + tmp_d) 
     397         
    609398        endfor 
    610399         
    611         cor_data = sam_data/sam_AttenFactor             //simply rescale the data 
    612  
    613 // do the error propagation piecewise 
    614         cor_err = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2            //sig a ^2 
    615         cor_err = sqrt(cor_err) 
    616  
     400        //we're done, get out w/no error 
     401 
     402         
     403        //TODO -- do I update COR header? 
     404//      cor_text[1] = date() + " " + time()             //date + time stamp 
     405 
     406        KillWaves/Z cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
     407        Killwaves/Z tmp_a,tmp_b,tmp_c,tmp_d,c_val,d_val 
     408         
     409        SetDataFolder root: 
     410        Return(0) 
     411End 
     412 
     413//background only 
     414// existence of data checked by dispatching routine 
     415// data has already been copied to COR folder 
     416Function V_CorrectMode_2() 
     417 
     418        //get SAM, BGD attenuation factor 
     419        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     420        Variable bgd_AttenFactor,bgd_atten_err 
     421        Variable ii 
     422        String detStr 
     423        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
     424        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     425         
     426        // these values apply to all of the detectors 
     427        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     428        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     429        bgd_AttenFactor = V_getAttenuator_transmission("BGD") 
     430        bgd_atten_err = V_getAttenuator_trans_err("BGD") 
     431 
     432        //get relative monitor counts (should all be 10^8, since normalized in add step) 
     433        // get transmission and trans error for SAM, EMP 
     434        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     435 
     436        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     437        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     438        sam_trans_err = V_getSampleTransError("SAM") 
     439         
     440        tmonbgd = V_getMonitorCount("BGD")              //monitor count in BGD 
     441 
     442 
     443        // and now loop through all of the detectors 
     444        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     445                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     446                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     447                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     448                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     449                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     450                Wave bgd_data = V_getDetectorDataW("BGD",detStr) 
     451                Wave bgd_err = V_getDetectorDataErrW("BGD",detStr) 
     452 
     453                 
     454        // to check for beam center mismatch -- simply warn, but do no shift 
     455        // 
     456 
     457                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     458                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     459         
     460                cbgd = V_getDet_beam_center_x("BGD",detStr) 
     461                rbgd = V_getDet_beam_center_y("BGD",detStr) 
     462 
     463         
     464                Duplicate/O cor_data cor1,bgd_temp,noadd_bgd 
     465 
     466                // TODO -- document this, make a note, so everyone knows this is not done 
     467                // skip this part, but duplicate the results of no shift condition 
     468                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     469                 
     470                        //get the shifted data array BGD, relative to SAM 
     471         
     472                xshift = cbgd-csam 
     473                yshift = rbgd-rsam 
     474                if(abs(xshift) <= wcen) 
     475                        xshift = 0 
     476                Endif 
     477                if(abs(yshift) <= wcen) 
     478                        yshift = 0 
     479                Endif 
     480                // for the BGD file - alert if needed, generate dummy "pass-through" values 
     481                // 
     482                if(xshift != 0 || yshift != 0) 
     483                        DoAlert 0, "Beam center mismatch for BGD file. Data has NOT been corrected." 
     484                endif 
     485                bgd_temp = bgd_data             // no shift, no effect 
     486                noadd_bgd = 1 
     487                //GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)            //bgd_temp 
     488         
     489 
     490                // **********    
     491                //do the sam-bgd subtraction,  deposit result in cor1 
     492                fsam = 1 
     493                fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
     494                 
     495                //print "fsam,fbgd = ",fsam,fbgd 
     496                 
     497                cor1 = fsam*sam_data/sam_AttenFactor - fbgd*bgd_temp/bgd_AttenFactor 
     498                cor1 *= noadd_bgd               //zeros out regions where arrays do not overlap, one otherwise 
     499         
     500        // do the error propagation piecewise    
     501                Duplicate/O sam_err, tmp_a, tmp_b 
     502                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     503                 
     504                tmp_b = (bgd_err/bgd_attenFactor)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2              //sig b ^2 
     505         
     506                cor_err = sqrt(tmp_a + tmp_b) 
     507 
     508        endfor 
     509 
     510        //we're done, get out w/no error 
     511 
     512        // TODO -- do I update COR header? 
     513        //cor_text[1] = date() + " " + time()           //date + time stamp 
     514 
     515        KillWaves/Z cor1,bgd_temp,noadd_bgd 
     516        Killwaves/Z tmp_a,tmp_b 
     517 
     518        SetDataFolder root: 
     519        Return(0) 
     520End 
     521 
     522// empty subtraction only 
     523// data does exist, checked by dispatch routine 
     524// 
     525Function V_CorrectMode_3() 
     526 
     527        //get SAM, EMP attenuation factor 
     528        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     529        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     530        Variable ii 
     531        String detStr 
     532        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
     533        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     534         
     535        // these values apply to all of the detectors 
     536        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     537        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     538        emp_AttenFactor = V_getAttenuator_transmission("EMP") 
     539        emp_atten_err = V_getAttenuator_trans_err("EMP") 
     540 
     541        //get relative monitor counts (should all be 10^8, since normalized in add step) 
     542        // get transmission and trans error for SAM, EMP 
     543        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     544 
     545        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     546        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     547        sam_trans_err = V_getSampleTransError("SAM") 
     548         
     549        tmonemp = V_getMonitorCount("EMP")              //monitor count in EMP 
     550        temp = V_getSampleTransmission("EMP")                   //trans emp 
     551        emp_trans_err = V_getSampleTransError("EMP") 
     552         
     553 
     554        // and now loop through all of the detectors 
     555        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     556                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     557                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     558                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     559                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     560                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     561                Wave emp_data = V_getDetectorDataW("EMP",detStr) 
     562                Wave emp_err = V_getDetectorDataErrW("EMP",detStr) 
     563                 
     564        // to check for beam center mismatch -- simply warn, but do no shift 
     565        // 
     566 
     567                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     568                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     569         
     570                cemp = V_getDet_beam_center_x("EMP",detStr)             //beamcenter of EMP 
     571                remp = V_getDet_beam_center_y("EMP",detStr) 
     572                 
     573                if(temp==0) 
     574                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     575                        temp=1 
     576                Endif 
     577         
     578         
     579                Duplicate/O cor_data cor1,emp_temp,noadd_emp 
     580 
     581                // TODO -- document this, make a note, so everyone knows this is not done 
     582                // skip this part, but duplicate the results of no shift condition 
     583                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     584                 
     585                        //get the shifted data array EMP, each relative to SAM 
     586         
     587                xshift = cemp-csam 
     588                yshift = remp-rsam 
     589                if(abs(xshift) <= wcen) 
     590                        xshift = 0 
     591                Endif 
     592                if(abs(yshift) <= wcen) 
     593                        yshift = 0 
     594                Endif 
     595                // for the EMP file - alert if needed, generate dummy "pass-through" values 
     596                // 
     597                if(xshift != 0 || yshift != 0) 
     598                        DoAlert 0, "Beam center mismatch for EMP file. Data has NOT been corrected." 
     599                endif 
     600                emp_temp = emp_data // no shift, no effect 
     601                noadd_emp = 1 
     602                //GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)            //emp_temp 
     603 
     604                // ********** 
     605         
     606                //do the sam-bgd subtraction,  deposit result in cor1 
     607                fsam = 1 
     608                femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
     609                 
     610                cor1 = fsam*sam_data/sam_AttenFactor - femp*(tsam/temp)*emp_temp/emp_AttenFactor 
     611                cor1 *= noadd_emp               //zeros out regions where arrays do not overlap, one otherwise 
     612         
     613        // do the error propagation piecewise    
     614                Duplicate/O sam_err, tmp_a, tmp_c ,c_val 
     615                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     616                 
     617                tmp_c = (sam_trans_err*emp_data/(temp*emp_attenFactor))^2 + (emp_err*tsam/(temp*emp_attenFactor))^2 
     618                tmp_c += (tsam*emp_data*emp_trans_err/(temp*temp*emp_attenFactor))^2 + (tsam*emp_data*emp_atten_err/(temp*emp_attenFactor^2))^2//total of 6 terms 
     619         
     620                cor_err = sqrt(tmp_a + tmp_c) 
     621         
     622        endfor 
     623         
     624        //we're done, get out w/no error 
     625 
     626        // TODO -- do I update COR header? 
     627        //cor_text[1] = date() + " " + time()           //date + time stamp 
     628 
     629        KillWaves/Z cor1,emp_temp,noadd_emp 
     630        Killwaves/Z tmp_a,tmp_c,c_val 
     631 
     632        SetDataFolder root: 
     633        Return(0) 
     634End 
     635 
     636// NO subtraction - simply rescales for attenuators 
     637// SAM data does exist, checked by dispatch routine 
     638// SAM data has already been copied to COR (both are the same at the start of the function) 
     639// 
     640//  TODO -- do I need to rescale to sam_trans here ?? 
     641// 
     642// 
     643Function V_CorrectMode_4() 
     644 
     645        //get SAM attenuation factor 
     646        Variable sam_AttenFactor,sam_atten_err,ii 
     647        String detStr 
     648         
     649        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     650        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     651 
     652 
     653        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     654                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     655                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     656                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     657                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     658                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     659         
     660         
     661                cor_data = sam_data/sam_AttenFactor             //simply rescale the data 
     662         
     663        // do the error propagation piecewise 
     664                cor_err = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2            //sig a ^2 
     665                cor_err = sqrt(cor_err) 
     666 
     667        endfor 
     668         
    617669        //TODO -- do I want to update COR header? 
    618670//      cor_text[1] = date() + " " + time()             //date + time stamp 
     
    622674End 
    623675 
    624 xFunction V_CorrectMode_11() 
    625         //create the necessary wave references 
    626         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    627         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    628         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    629         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    630         WAVE bgd_data=$"root:Packages:NIST:BGD:linear_data" 
    631         WAVE bgd_reals=$"root:Packages:NIST:BGD:realsread" 
    632         WAVE bgd_ints=$"root:Packages:NIST:BGD:integersread" 
    633         WAVE/T bgd_text=$"root:Packages:NIST:BGD:textread" 
    634         WAVE emp_data=$"root:Packages:NIST:EMP:linear_data" 
    635         WAVE emp_reals=$"root:Packages:NIST:EMP:realsread" 
    636         WAVE emp_ints=$"root:Packages:NIST:EMP:integersread" 
    637         WAVE/T emp_text=$"root:Packages:NIST:EMP:textread" 
    638         WAVE drk_data=$"root:Packages:NIST:DRK:linear_data" 
    639         WAVE drk_reals=$"root:Packages:NIST:DRK:realsread" 
    640         WAVE drk_ints=$"root:Packages:NIST:DRK:integersread" 
    641         WAVE/T drk_text=$"root:Packages:NIST:DRK:textread" 
    642         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    643         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    644  
    645         // needed to propagate error 
    646         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    647         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    648         WAVE bgd_err =$"root:Packages:NIST:BGD:linear_data_error" 
    649         WAVE emp_err =$"root:Packages:NIST:EMP:linear_data_error" 
    650         WAVE drk_err =$"root:Packages:NIST:DRK:linear_data_error" 
    651         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    652          
    653         Variable sam_trans_err,emp_trans_err 
    654         sam_trans_err = sam_reals[41] 
    655         emp_trans_err = emp_reals[41] 
    656          
    657         //get sam and bgd attenuation factors 
    658         String fileStr="" 
    659         Variable lambda,attenNo,sam_AttenFactor,bgd_attenFactor,emp_AttenFactor 
     676 
     677 
     678Function V_CorrectMode_11() 
     679 
     680        //get SAM, BGD, EMP attenuation factor 
     681        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     682        Variable bgd_AttenFactor,bgd_atten_err 
     683        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     684        Variable ii 
     685        String detStr 
    660686        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
    661         Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp,time_sam,time_drk,savmon_sam 
    662         Variable sam_atten_err,bgd_atten_err,emp_atten_err 
    663         fileStr = sam_text[3] 
    664         lambda = sam_reals[26] 
    665         attenNo = sam_reals[3] 
    666         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    667         fileStr = bgd_text[3] 
    668         lambda = bgd_reals[26] 
    669         attenNo = bgd_reals[3] 
    670         bgd_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,bgd_atten_err) 
    671         fileStr = emp_text[3] 
    672         lambda = emp_reals[26] 
    673         attenNo = emp_reals[3] 
    674         emp_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,emp_atten_err) 
    675          
     687        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     688        Variable savmon_sam,time_sam,time_drk 
     689         
     690        // these values apply to all of the detectors 
     691        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     692        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     693        bgd_AttenFactor = V_getAttenuator_transmission("BGD") 
     694        bgd_atten_err = V_getAttenuator_trans_err("BGD") 
     695        emp_AttenFactor = V_getAttenuator_transmission("EMP") 
     696        emp_atten_err = V_getAttenuator_trans_err("EMP") 
     697 
    676698        //get relative monitor counts (should all be 10^8, since normalized in add step) 
    677         tmonsam = sam_reals[0]          //monitor count in SAM 
    678         tsam = sam_reals[4]             //SAM transmission 
    679         csam = sam_reals[16]            //x center 
    680         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    681         tmonbgd = bgd_reals[0]          //monitor count in BGD 
    682         cbgd = bgd_reals[16] 
    683         rbgd = bgd_reals[17] 
    684         tmonemp = emp_reals[0]          //monitor count in EMP 
    685         temp = emp_reals[4]                     //trans emp 
    686         cemp = emp_reals[16]            //beamcenter of EMP 
    687         remp = emp_reals[17] 
    688         savmon_sam=sam_reals[1]         //true monitor count in SAM 
    689         time_sam = sam_ints[2]          //count time SAM 
    690         time_drk = drk_ints[2]          //drk count time 
    691          
    692         NVAR pixelsX = root:myGlobals:gNPixelsX 
    693         NVAR pixelsY = root:myGlobals:gNPixelsY 
    694         //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
    695         Make/D/O/N=(pixelsX,pixelsY) drk_temp, drk_tmp_err 
    696         drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
    697         drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK 
    698          
    699         if(temp==0) 
    700                 DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
    701                 temp=1 
    702         Endif 
    703          
    704         //get the shifted data arrays, EMP and BGD, each relative to SAM 
    705         Make/D/O/N=(pixelsX,pixelsY) cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
    706         xshift = cbgd-csam 
    707         yshift = rbgd-rsam 
    708         if(abs(xshift) <= wcen) 
    709                 xshift = 0 
    710         Endif 
    711         if(abs(yshift) <= wcen) 
    712                 yshift = 0 
    713         Endif 
    714         GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)              //bgd_temp 
    715          
    716         xshift = cemp-csam 
    717         yshift = remp-rsam 
    718         if(abs(xshift) <= wcen) 
    719                 xshift = 0 
    720         Endif 
    721         if(abs(yshift) <= wcen) 
    722                 yshift = 0 
    723         Endif 
    724         GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)              //emp_temp 
    725         //always ignore the DRK center shift 
    726          
    727         //do the subtraction 
    728         fsam=1 
    729         femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
    730         fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
    731         cor1 = fsam*sam_data/sam_attenFactor 
    732         cor1 -= (tsam/temp)*(femp*emp_temp/emp_attenFactor - fbgd*bgd_temp/bgd_attenFactor) 
    733         cor1 -= (fbgd*bgd_temp/bgd_attenFactor - drk_temp) 
    734         cor1 -= drk_temp/sam_attenFactor 
    735         cor1 *= noadd_bgd*noadd_emp             //zero out the array mismatch values 
    736          
    737 // do the error propagation piecewise    
    738         Duplicate/O sam_err, tmp_a, tmp_b, tmp_c, tmp_d,c_val,d_val 
    739         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    740          
    741         tmp_b = (bgd_err/bgd_attenFactor)^2*(tsam/temp - 1)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2*(1-tsam/temp)^2            //sig b ^2 
    742  
    743         tmp_c = (sam_trans_err/temp)^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
    744         tmp_c += (tsam/temp^2)^2*emp_trans_err^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
    745          
    746         tmp_d = (tsam/(temp*emp_attenFactor))^2*(emp_err)^2 + (tsam*emp_data/(temp*emp_attenFactor^2))^2*(emp_atten_err)^2 
    747  
    748         cor_err = sqrt(tmp_a + tmp_b + tmp_c + tmp_d + drk_tmp_err^2) 
     699        // get transmission and trans error for SAM, EMP 
     700        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     701         
     702        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     703        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     704        sam_trans_err = V_getSampleTransError("SAM") 
     705         
     706        tmonemp = V_getMonitorCount("EMP")              //monitor count in EMP 
     707        temp = V_getSampleTransmission("EMP")                   //trans emp 
     708        emp_trans_err = V_getSampleTransError("EMP") 
     709         
     710        tmonbgd = V_getMonitorCount("BGD")              //monitor count in BGD 
     711 
     712        // for proper scaling, get the time and actual monitor counts 
     713        // TODO -- make sure that these calls are reading the proper values 
     714        savmon_sam = V_getBeamMonNormSaved_count("SAM")         //true monitor count in SAM 
     715        time_sam = V_getCount_time("SAM")               //count time SAM 
     716        time_drk = V_getCount_time("DRK")       //drk count time 
     717 
     718 
     719        // and now loop through all of the detectors 
     720        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     721                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     722                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     723                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     724                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     725                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     726                Wave bgd_data = V_getDetectorDataW("BGD",detStr) 
     727                Wave bgd_err = V_getDetectorDataErrW("BGD",detStr) 
     728                Wave emp_data = V_getDetectorDataW("EMP",detStr) 
     729                Wave emp_err = V_getDetectorDataErrW("EMP",detStr) 
     730                Wave drk_data = V_getDetectorDataW("DRK",detStr) 
     731                Wave drk_err = V_getDetectorDataErrW("DRK",detStr) 
     732                 
     733        // to check for beam center mismatch -- simply warn, but do no shift 
     734        // 
     735 
     736                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     737                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     738         
     739                cbgd = V_getDet_beam_center_x("BGD",detStr) 
     740                rbgd = V_getDet_beam_center_y("BGD",detStr) 
     741         
     742                cemp = V_getDet_beam_center_x("EMP",detStr)             //beamcenter of EMP 
     743                remp = V_getDet_beam_center_y("EMP",detStr) 
     744                 
     745                //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
     746                Duplicate/O drk_data drk_temp, drk_tmp_err 
     747                drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
     748                drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK   
     749                 
     750                if(temp==0) 
     751                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     752                        temp=1 
     753                Endif 
     754         
     755                Duplicate/O cor_data cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp 
     756 
     757                // TODO -- document this, make a note, so everyone knows this is not done 
     758                // skip this part, but duplicate the results of no shift condition 
     759                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     760                 
     761                        //get the shifted data arrays, EMP and BGD, each relative to SAM 
     762         
     763                xshift = cbgd-csam 
     764                yshift = rbgd-rsam 
     765                if(abs(xshift) <= wcen) 
     766                        xshift = 0 
     767                Endif 
     768                if(abs(yshift) <= wcen) 
     769                        yshift = 0 
     770                Endif 
     771                // for the BGD file - alert if needed, generate dummy "pass-through" values 
     772                // 
     773                if(xshift != 0 || yshift != 0) 
     774                        DoAlert 0, "Beam center mismatch for BGD file. Data has NOT been corrected." 
     775                endif 
     776                bgd_temp = bgd_data             // no shift, no effect 
     777                noadd_bgd = 1 
     778                //GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)            //bgd_temp 
     779         
     780                xshift = cemp-csam 
     781                yshift = remp-rsam 
     782                if(abs(xshift) <= wcen) 
     783                        xshift = 0 
     784                Endif 
     785                if(abs(yshift) <= wcen) 
     786                        yshift = 0 
     787                Endif 
     788                // for the EMP file - alert if needed, generate dummy "pass-through" values 
     789                // 
     790                if(xshift != 0 || yshift != 0) 
     791                        DoAlert 0, "Beam center mismatch for EMP file. Data has NOT been corrected." 
     792                endif 
     793                emp_temp = emp_data // no shift, no effect 
     794                noadd_emp = 1 
     795                //GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)            //emp_temp 
     796 
     797 
     798                //always ignore the DRK center shift 
     799         
     800                // ************  
     801                //do the subtraction 
     802                fsam=1 
     803                femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
     804                fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
     805                cor1 = fsam*sam_data/sam_attenFactor 
     806                cor1 -= (tsam/temp)*(femp*emp_temp/emp_attenFactor - fbgd*bgd_temp/bgd_attenFactor) 
     807                cor1 -= (fbgd*bgd_temp/bgd_attenFactor - drk_temp) 
     808                cor1 -= drk_temp/sam_attenFactor 
     809                cor1 *= noadd_bgd*noadd_emp             //zero out the array mismatch values 
     810                 
     811        // do the error propagation piecewise    
     812                Duplicate/O sam_err, tmp_a, tmp_b, tmp_c, tmp_d,c_val,d_val 
     813                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     814                 
     815                tmp_b = (bgd_err/bgd_attenFactor)^2*(tsam/temp - 1)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2*(1-tsam/temp)^2            //sig b ^2 
     816         
     817                tmp_c = (sam_trans_err/temp)^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
     818                tmp_c += (tsam/temp^2)^2*emp_trans_err^2*(emp_data/emp_attenFactor-bgd_data/bgd_attenFactor)^2 
     819                 
     820                tmp_d = (tsam/(temp*emp_attenFactor))^2*(emp_err)^2 + (tsam*emp_data/(temp*emp_attenFactor^2))^2*(emp_atten_err)^2 
     821         
     822                cor_err = sqrt(tmp_a + tmp_b + tmp_c + tmp_d + drk_tmp_err^2) 
     823         
     824        endfor 
    749825         
    750826        //we're done, get out w/no error 
    751         //set the COR data to the result 
    752         cor_data = cor1 
    753         cor_data_display = cor1 
    754  
    755         //update COR header 
    756         cor_text[1] = date() + " " + time()             //date + time stamp 
     827 
     828 
     829        //TODO -- do I update COR header? 
     830//      cor_text[1] = date() + " " + time()             //date + time stamp 
    757831 
    758832        KillWaves/Z cor1,bgd_temp,noadd_bgd,emp_temp,noadd_emp,drk_temp 
     
    765839//bgd and drk subtraction 
    766840// 
    767 xFunction V_CorrectMode_12() 
    768         //create the necessary wave references 
    769         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    770         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    771         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    772         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    773         WAVE bgd_data=$"root:Packages:NIST:BGD:linear_data" 
    774         WAVE bgd_reals=$"root:Packages:NIST:BGD:realsread" 
    775         WAVE bgd_ints=$"root:Packages:NIST:BGD:integersread" 
    776         WAVE/T bgd_text=$"root:Packages:NIST:BGD:textread" 
    777         WAVE drk_data=$"root:Packages:NIST:DRK:linear_data" 
    778         WAVE drk_reals=$"root:Packages:NIST:DRK:realsread" 
    779         WAVE drk_ints=$"root:Packages:NIST:DRK:integersread" 
    780         WAVE/T drk_text=$"root:Packages:NIST:DRK:textread" 
    781         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    782         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    783  
    784         // needed to propagate error 
    785         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    786         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    787         WAVE bgd_err =$"root:Packages:NIST:BGD:linear_data_error" 
    788         WAVE drk_err =$"root:Packages:NIST:DRK:linear_data_error" 
    789         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    790          
    791         Variable sam_trans_err 
    792         sam_trans_err = sam_reals[41] 
    793          
    794          
    795         //get sam and bgd attenuation factors 
    796         String fileStr="" 
    797         Variable lambda,attenNo,sam_AttenFactor,bgd_attenFactor 
     841Function V_CorrectMode_12() 
     842 
     843        //get SAM, BGD, EMP attenuation factor 
     844        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     845        Variable bgd_AttenFactor,bgd_atten_err 
     846        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     847        Variable ii 
     848        String detStr 
    798849        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
    799         Variable wcen=0.001,time_drk,time_sam,savmon_sam,tsam 
    800         Variable sam_atten_err,bgd_atten_err 
    801         fileStr = sam_text[3] 
    802         lambda = sam_reals[26] 
    803         attenNo = sam_reals[3] 
    804         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    805         fileStr = bgd_text[3] 
    806         lambda = bgd_reals[26] 
    807         attenNo = bgd_reals[3] 
    808         bgd_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,bgd_atten_err) 
    809          
     850        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     851        Variable savmon_sam,time_sam,time_drk 
     852         
     853        // these values apply to all of the detectors 
     854        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     855        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     856        bgd_AttenFactor = V_getAttenuator_transmission("BGD") 
     857        bgd_atten_err = V_getAttenuator_trans_err("BGD") 
     858 
    810859        //get relative monitor counts (should all be 10^8, since normalized in add step) 
    811         tmonsam = sam_reals[0]          //monitor count in SAM 
    812         tsam = sam_reals[4]             //SAM transmission 
    813         csam = sam_reals[16]            //x center 
    814         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    815         tmonbgd = bgd_reals[0]          //monitor count in BGD 
    816         cbgd = bgd_reals[16] 
    817         rbgd = bgd_reals[17] 
    818         savmon_sam=sam_reals[1]         //true monitor count in SAM 
    819         time_sam = sam_ints[2]          //count time SAM 
    820         time_drk = drk_ints[2]          //drk count time 
    821          
    822         NVAR pixelsX = root:myGlobals:gNPixelsX 
    823         NVAR pixelsY = root:myGlobals:gNPixelsY 
    824         //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
    825         Make/D/O/N=(pixelsX,pixelsY) drk_temp,drk_tmp_err 
    826         drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
    827         drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK 
    828  
    829         // set up beamcenter shift, relative to SAM 
    830         xshift = cbgd-csam 
    831         yshift = rbgd-rsam 
    832         if(abs(xshift) <= wcen) 
    833                 xshift = 0 
    834         Endif 
    835         if(abs(yshift) <= wcen) 
    836                 yshift = 0 
    837         Endif 
    838         //get shifted data arrays, relative to SAM 
    839         Make/D/O/N=(pixelsX,pixelsY) cor1,bgd_temp,noadd_bgd            //temp arrays 
    840         GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)              //bgd_temp is the BGD  
    841         //always ignore the DRK center shift 
    842          
    843         //do the sam-bgd subtraction,  deposit result in cor1 
    844         fsam = 1 
    845         fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
    846          
    847         cor1 = fsam*sam_data/sam_AttenFactor + fbgd*tsam*bgd_temp/bgd_AttenFactor 
    848         cor1 += -1*(fbgd*bgd_temp/bgd_attenFactor - drk_temp) - drk_temp/sam_attenFactor 
    849         cor1 *= noadd_bgd               //zeros out regions where arrays do not overlap, one otherwise 
    850  
    851 // do the error propagation piecewise    
    852         Duplicate/O sam_err, tmp_a, tmp_b 
    853         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    854          
    855         tmp_b = (bgd_err/bgd_attenFactor)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2              //sig b ^2 
    856  
    857         cor_err = sqrt(tmp_a + tmp_b + drk_tmp_err^2) 
    858  
     860        // get transmission and trans error for SAM, EMP 
     861        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     862         
     863        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     864        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     865        sam_trans_err = V_getSampleTransError("SAM") 
     866         
     867        tmonbgd = V_getMonitorCount("BGD")              //monitor count in BGD 
     868 
     869        // for proper scaling, get the time and actual monitor counts 
     870        // TODO -- make sure that these calls are reading the proper values 
     871        savmon_sam = V_getBeamMonNormSaved_count("SAM")         //true monitor count in SAM 
     872        time_sam = V_getCount_time("SAM")               //count time SAM 
     873        time_drk = V_getCount_time("DRK")       //drk count time 
     874 
     875 
     876        // and now loop through all of the detectors 
     877        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     878                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     879                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     880                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     881                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     882                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     883                Wave bgd_data = V_getDetectorDataW("BGD",detStr) 
     884                Wave bgd_err = V_getDetectorDataErrW("BGD",detStr) 
     885                Wave drk_data = V_getDetectorDataW("DRK",detStr) 
     886                Wave drk_err = V_getDetectorDataErrW("DRK",detStr) 
     887                 
     888        // to check for beam center mismatch -- simply warn, but do no shift 
     889        // 
     890 
     891                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     892                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     893         
     894                cbgd = V_getDet_beam_center_x("BGD",detStr) 
     895                rbgd = V_getDet_beam_center_y("BGD",detStr) 
     896                 
     897                //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
     898                Duplicate/O drk_data drk_temp, drk_tmp_err 
     899                drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
     900                drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK   
     901                 
     902                if(temp==0) 
     903                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     904                        temp=1 
     905                Endif 
     906         
     907                Duplicate/O cor_data cor1,bgd_temp,noadd_bgd 
     908 
     909                // TODO -- document this, make a note, so everyone knows this is not done 
     910                // skip this part, but duplicate the results of no shift condition 
     911                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     912                 
     913                        //get the shifted data arrays, EMP and BGD, each relative to SAM 
     914         
     915                xshift = cbgd-csam 
     916                yshift = rbgd-rsam 
     917                if(abs(xshift) <= wcen) 
     918                        xshift = 0 
     919                Endif 
     920                if(abs(yshift) <= wcen) 
     921                        yshift = 0 
     922                Endif 
     923                // for the BGD file - alert if needed, generate dummy "pass-through" values 
     924                // 
     925                if(xshift != 0 || yshift != 0) 
     926                        DoAlert 0, "Beam center mismatch for BGD file. Data has NOT been corrected." 
     927                endif 
     928                bgd_temp = bgd_data             // no shift, no effect 
     929                noadd_bgd = 1 
     930                //GetShiftedArray(bgd_data,bgd_temp,noadd_bgd,xshift,yshift)            //bgd_temp 
     931 
     932 
     933                //always ignore the DRK center shift 
     934 
     935                // ************  
     936                //do the sam-bgd subtraction,  deposit result in cor1 
     937                fsam = 1 
     938                fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
     939                 
     940                cor1 = fsam*sam_data/sam_AttenFactor + fbgd*tsam*bgd_temp/bgd_AttenFactor 
     941                cor1 += -1*(fbgd*bgd_temp/bgd_attenFactor - drk_temp) - drk_temp/sam_attenFactor 
     942                cor1 *= noadd_bgd               //zeros out regions where arrays do not overlap, one otherwise 
     943         
     944        // do the error propagation piecewise    
     945                Duplicate/O sam_err, tmp_a, tmp_b 
     946                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     947                 
     948                tmp_b = (bgd_err/bgd_attenFactor)^2 + (bgd_atten_err*bgd_data/bgd_attenFactor^2)^2              //sig b ^2 
     949         
     950                cor_err = sqrt(tmp_a + tmp_b + drk_tmp_err^2) 
     951 
     952        endfor 
     953         
    859954        //we're done, get out w/no error 
    860         //set the COR_data to the result 
    861         cor_data = cor1 
    862         cor_data_display = cor1 
    863  
    864         //update COR header 
    865         cor_text[1] = date() + " " + time()             //date + time stamp 
     955         
     956        // TODO -- do I update COR header? 
     957//      cor_text[1] = date() + " " + time()             //date + time stamp 
    866958 
    867959        KillWaves/Z cor1,bgd_temp,noadd_bgd,drk_temp 
     
    876968//scale DRK by monitor count scaling factor and the ratio of couting times 
    877969//to place the DRK file on equal footing 
    878 xFunction V_CorrectMode_13() 
    879         //create the necessary wave references 
    880         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    881         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    882         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    883         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    884         WAVE emp_data=$"root:Packages:NIST:EMP:linear_data" 
    885         WAVE emp_reals=$"root:Packages:NIST:EMP:realsread" 
    886         WAVE emp_ints=$"root:Packages:NIST:EMP:integersread" 
    887         WAVE/T emp_text=$"root:Packages:NIST:EMP:textread" 
    888         WAVE drk_data=$"root:Packages:NIST:DRK:linear_data" 
    889         WAVE drk_reals=$"root:Packages:NIST:DRK:realsread" 
    890         WAVE drk_ints=$"root:Packages:NIST:DRK:integersread" 
    891         WAVE/T drk_text=$"root:Packages:NIST:DRK:textread" 
    892         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    893         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    894  
    895         // needed to propagate error 
    896         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    897         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    898         WAVE emp_err =$"root:Packages:NIST:EMP:linear_data_error" 
    899         WAVE drk_err =$"root:Packages:NIST:DRK:linear_data_error" 
    900         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    901          
    902         Variable sam_trans_err,emp_trans_err 
    903         sam_trans_err = sam_reals[41] 
    904         emp_trans_err = emp_reals[41] 
    905          
    906         //get sam and bgd attenuation factors (DRK irrelevant) 
    907         String fileStr="" 
    908         Variable lambda,attenNo,sam_AttenFactor,emp_attenFactor 
    909         Variable tmonsam,fsam,femp,xshift,yshift,rsam,csam,remp,cemp,tmonemp 
    910         Variable wcen=0.001,tsam,temp,savmon_sam,time_sam,time_drk 
    911         Variable sam_atten_err,emp_atten_err 
    912         fileStr = sam_text[3] 
    913         lambda = sam_reals[26] 
    914         attenNo = sam_reals[3] 
    915         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    916         fileStr = emp_text[3] 
    917         lambda = emp_reals[26] 
    918         attenNo = emp_reals[3] 
    919         emp_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,emp_atten_err) 
    920          
     970Function V_CorrectMode_13() 
     971 
     972        //get SAM, EMP attenuation factor 
     973        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     974        Variable bgd_AttenFactor,bgd_atten_err 
     975        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     976        Variable ii 
     977        String detStr 
     978        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
     979        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     980        Variable savmon_sam,time_sam,time_drk 
     981         
     982        // these values apply to all of the detectors 
     983        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     984        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     985        emp_AttenFactor = V_getAttenuator_transmission("EMP") 
     986        emp_atten_err = V_getAttenuator_trans_err("EMP") 
     987 
    921988        //get relative monitor counts (should all be 10^8, since normalized in add step) 
    922         tmonsam = sam_reals[0]          //monitor count in SAM 
    923         tsam = sam_reals[4]             //SAM transmission 
    924         csam = sam_reals[16]            //x center 
    925         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    926         tmonemp = emp_reals[0]          //monitor count in EMP 
    927         temp = emp_reals[4]                     //trans emp 
    928         cemp = emp_reals[16]            //beamcenter of EMP 
    929         remp = emp_reals[17] 
    930         savmon_sam=sam_reals[1]         //true monitor count in SAM 
    931         time_sam = sam_ints[2]          //count time SAM 
    932         time_drk = drk_ints[2]          //drk count time 
    933          
    934         NVAR pixelsX = root:myGlobals:gNPixelsX 
    935         NVAR pixelsY = root:myGlobals:gNPixelsY 
    936         //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
    937         Make/D/O/N=(pixelsX,pixelsY) drk_temp,drk_tmp_err 
    938         drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
    939         drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK 
    940  
    941          
    942         if(temp==0) 
    943                 DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
    944                 temp=1 
    945         Endif 
    946          
    947         //Print "rbgd,cbgd = ",rbgd,cbgd 
    948         // set up beamcenter shift, relative to SAM 
    949         xshift = cemp-csam 
    950         yshift = remp-rsam 
    951         if(abs(xshift) <= wcen) 
    952                 xshift = 0 
    953         Endif 
    954         if(abs(yshift) <= wcen) 
    955                 yshift = 0 
    956         Endif 
    957         //get shifted data arrays, relative to SAM 
    958         Make/D/O/N=(pixelsX,pixelsY) cor1,emp_temp,noadd_emp            //temp arrays 
    959         GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)              //emp_temp is the EMP 
    960         //always ignore beamcenter shift for DRK 
    961          
    962         //do the sam-bgd subtraction,  deposit result in cor1 
    963         fsam = 1 
    964         femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
    965          
    966         cor1 = fsam*sam_data/sam_AttenFactor - femp*(tsam/temp)*emp_temp/emp_AttenFactor 
    967         cor1 += drk_temp - drk_temp/sam_attenFactor 
    968         cor1 *= noadd_emp               //zeros out regions where arrays do not overlap, one otherwise 
    969  
    970 // do the error propagation piecewise    
    971         Duplicate/O sam_err, tmp_a, tmp_c, c_val 
    972         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    973          
    974         tmp_c = (sam_trans_err*emp_data/(temp*emp_attenFactor))^2 + (emp_err*tsam/(temp*emp_attenFactor))^2 
    975         tmp_c += (tsam*emp_data*emp_trans_err/(temp*temp*emp_attenFactor))^2 + (tsam*emp_data*emp_atten_err/(temp*emp_attenFactor^2))^2//total of 6 terms 
    976          
    977         cor_err = sqrt(tmp_a + tmp_c + drk_tmp_err^2) 
    978          
     989        // get transmission and trans error for SAM, EMP 
     990        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     991         
     992        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     993        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     994        sam_trans_err = V_getSampleTransError("SAM") 
     995         
     996        tmonemp = V_getMonitorCount("EMP")              //monitor count in EMP 
     997        temp = V_getSampleTransmission("EMP")                   //trans emp 
     998        emp_trans_err = V_getSampleTransError("EMP") 
     999 
     1000        // for proper scaling, get the time and actual monitor counts 
     1001        // TODO -- make sure that these calls are reading the proper values 
     1002        savmon_sam = V_getBeamMonNormSaved_count("SAM")         //true monitor count in SAM 
     1003        time_sam = V_getCount_time("SAM")               //count time SAM 
     1004        time_drk = V_getCount_time("DRK")       //drk count time 
     1005 
     1006 
     1007        // and now loop through all of the detectors 
     1008        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     1009                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     1010                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     1011                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     1012                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     1013                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     1014                Wave emp_data = V_getDetectorDataW("EMP",detStr) 
     1015                Wave emp_err = V_getDetectorDataErrW("EMP",detStr) 
     1016                Wave drk_data = V_getDetectorDataW("DRK",detStr) 
     1017                Wave drk_err = V_getDetectorDataErrW("DRK",detStr) 
     1018                 
     1019        // to check for beam center mismatch -- simply warn, but do no shift 
     1020        // 
     1021 
     1022                csam = V_getDet_beam_center_x("SAM",detStr)             //x center 
     1023                rsam = V_getDet_beam_center_y("SAM",detStr)             //beam (x,y) define center of corrected field 
     1024         
     1025                cemp = V_getDet_beam_center_x("EMP",detStr)             //beamcenter of EMP 
     1026                remp = V_getDet_beam_center_y("EMP",detStr) 
     1027                 
     1028                //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
     1029                Duplicate/O drk_data drk_temp, drk_tmp_err 
     1030                drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
     1031                drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK   
     1032                 
     1033                if(temp==0) 
     1034                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     1035                        temp=1 
     1036                Endif 
     1037         
     1038                Duplicate/O cor_data cor1,emp_temp,noadd_emp 
     1039 
     1040                // TODO -- document this, make a note, so everyone knows this is not done 
     1041                // skip this part, but duplicate the results of no shift condition 
     1042                //  where bgd_temp = input data, and noadd_bgd = 1 (so no data is zeroed out) 
     1043                 
     1044                        //get the shifted data arrays, EMP , each relative to SAM 
     1045         
     1046                xshift = cemp-csam 
     1047                yshift = remp-rsam 
     1048                if(abs(xshift) <= wcen) 
     1049                        xshift = 0 
     1050                Endif 
     1051                if(abs(yshift) <= wcen) 
     1052                        yshift = 0 
     1053                Endif 
     1054                // for the EMP file - alert if needed, generate dummy "pass-through" values 
     1055                // 
     1056                if(xshift != 0 || yshift != 0) 
     1057                        DoAlert 0, "Beam center mismatch for EMP file. Data has NOT been corrected." 
     1058                endif 
     1059                emp_temp = emp_data // no shift, no effect 
     1060                noadd_emp = 1 
     1061                //GetShiftedArray(emp_data,emp_temp,noadd_emp,xshift,yshift)            //emp_temp 
     1062 
     1063 
     1064                //always ignore the DRK center shift 
     1065                 
     1066                // ***************       
     1067                //do the sam-bgd subtraction,  deposit result in cor1 
     1068                fsam = 1 
     1069                femp = tmonsam/tmonemp          //this should be ==1 since normalized files 
     1070                 
     1071                cor1 = fsam*sam_data/sam_AttenFactor - femp*(tsam/temp)*emp_temp/emp_AttenFactor 
     1072                cor1 += drk_temp - drk_temp/sam_attenFactor 
     1073                cor1 *= noadd_emp               //zeros out regions where arrays do not overlap, one otherwise 
     1074         
     1075        // do the error propagation piecewise    
     1076                Duplicate/O sam_err, tmp_a, tmp_c, c_val 
     1077                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     1078                 
     1079                tmp_c = (sam_trans_err*emp_data/(temp*emp_attenFactor))^2 + (emp_err*tsam/(temp*emp_attenFactor))^2 
     1080                tmp_c += (tsam*emp_data*emp_trans_err/(temp*temp*emp_attenFactor))^2 + (tsam*emp_data*emp_atten_err/(temp*emp_attenFactor^2))^2//total of 6 terms 
     1081                 
     1082                cor_err = sqrt(tmp_a + tmp_c + drk_tmp_err^2) 
     1083 
     1084        endfor 
     1085                 
    9791086        //we're done, get out w/no error 
    980         //set the COR data to the result 
    981         cor_data = cor1 
    982         cor_data_display = cor1 
    983  
    984         //update COR header 
    985         cor_text[1] = date() + " " + time()             //date + time stamp 
     1087 
     1088        // TODO -- do I update COR header? 
     1089//      cor_text[1] = date() + " " + time()             //date + time stamp 
    9861090 
    9871091        KillWaves/Z cor1,emp_temp,noadd_emp,drk_temp 
     
    9941098// ONLY drk subtraction 
    9951099// 
    996 xFunction V_CorrectMode_14() 
    997         //create the necessary wave references 
    998         WAVE sam_data=$"root:Packages:NIST:SAM:linear_data" 
    999         WAVE sam_reals=$"root:Packages:NIST:SAM:realsread" 
    1000         WAVE sam_ints=$"root:Packages:NIST:SAM:integersread" 
    1001         WAVE/T sam_text=$"root:Packages:NIST:SAM:textread" 
    1002         WAVE drk_data=$"root:Packages:NIST:DRK:linear_data" 
    1003         WAVE drk_reals=$"root:Packages:NIST:DRK:realsread" 
    1004         WAVE drk_ints=$"root:Packages:NIST:DRK:integersread" 
    1005         WAVE/T drk_text=$"root:Packages:NIST:DRK:textread" 
    1006         WAVE cor_data=$"root:Packages:NIST:COR:linear_data" 
    1007         WAVE/T cor_text=$"root:Packages:NIST:COR:textread" 
    1008  
    1009         // needed to propagate error 
    1010         WAVE cor_data_display=$"root:Packages:NIST:COR:data"            //just for the final copy 
    1011         WAVE sam_err =$"root:Packages:NIST:SAM:linear_data_error" 
    1012         WAVE drk_err =$"root:Packages:NIST:DRK:linear_data_error" 
    1013         WAVE cor_err =$"root:Packages:NIST:COR:linear_data_error" 
    1014          
    1015         Variable sam_trans_err 
    1016         sam_trans_err = sam_reals[41] 
    1017          
    1018          
    1019         //get sam and bgd attenuation factors 
    1020         String fileStr="" 
    1021         Variable lambda,attenNo,sam_AttenFactor,bgd_attenFactor 
     1100Function V_CorrectMode_14() 
     1101 
     1102        //get SAM, EMP attenuation factor 
     1103        Variable sam_AttenFactor,sam_atten_err,sam_trans_err 
     1104        Variable bgd_AttenFactor,bgd_atten_err 
     1105        Variable emp_AttenFactor,emp_atten_err,emp_trans_err 
     1106        Variable ii 
     1107        String detStr 
    10221108        Variable tmonsam,fsam,fbgd,xshift,yshift,rsam,csam,rbgd,cbgd,tmonbgd 
    1023         Variable wcen=0.001,time_drk,time_sam,savmon_sam,tsam 
    1024         Variable sam_atten_err 
    1025         fileStr = sam_text[3] 
    1026         lambda = sam_reals[26] 
    1027         attenNo = sam_reals[3] 
    1028         sam_AttenFactor = AttenuationFactor(fileStr,lambda,AttenNo,sam_atten_err) 
    1029          
     1109        Variable wcen=0.001,tsam,temp,remp,cemp,tmonemp,femp 
     1110        Variable savmon_sam,time_sam,time_drk 
     1111         
     1112        // these values apply to all of the detectors 
     1113        sam_AttenFactor = V_getAttenuator_transmission("SAM") 
     1114        sam_atten_err = V_getAttenuator_trans_err("SAM") 
     1115 
     1116 
    10301117        //get relative monitor counts (should all be 10^8, since normalized in add step) 
    1031         tmonsam = sam_reals[0]          //monitor count in SAM 
    1032         tsam = sam_reals[4]             //SAM transmission 
    1033         csam = sam_reals[16]            //x center 
    1034         rsam = sam_reals[17]            //beam (x,y) define center of corrected field 
    1035  
    1036         savmon_sam=sam_reals[1]         //true monitor count in SAM 
    1037         time_sam = sam_ints[2]          //count time SAM 
    1038         time_drk = drk_ints[2]          //drk count time 
    1039          
    1040         NVAR pixelsX = root:myGlobals:gNPixelsX 
    1041         NVAR pixelsY = root:myGlobals:gNPixelsY 
    1042         //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
    1043         Make/D/O/N=(pixelsX,pixelsY) drk_temp,drk_tmp_err 
    1044         drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
    1045         drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK 
    1046  
    1047         Make/D/O/N=(pixelsX,pixelsY) cor1       //temp arrays 
    1048         //always ignore the DRK center shift 
    1049          
    1050         //do the subtraction,  deposit result in cor1 
    1051         fsam = 1 
    1052         fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
    1053          
    1054         //correct sam for attenuators, and do the same to drk, since it was scaled to sam count time 
    1055         cor1 = fsam*sam_data/sam_AttenFactor  - drk_temp/sam_attenFactor 
    1056  
    1057 // do the error propagation piecewise    
    1058         Duplicate/O sam_err, tmp_a 
    1059         tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
    1060  
    1061         cor_err = sqrt(tmp_a + drk_tmp_err^2) 
    1062          
     1118        // get transmission and trans error for SAM, EMP 
     1119        // TODO -- verify that the  call to V_getMonitorCount() is really rescaled to 10^8, and saved is the "true" count 
     1120         
     1121        tmonsam = V_getMonitorCount("SAM")              //monitor count in SAM 
     1122        tsam = V_getSampleTransmission("SAM")           //SAM transmission 
     1123        sam_trans_err = V_getSampleTransError("SAM") 
     1124         
     1125 
     1126        // for proper scaling, get the time and actual monitor counts 
     1127        // TODO -- make sure that these calls are reading the proper values 
     1128        savmon_sam = V_getBeamMonNormSaved_count("SAM")         //true monitor count in SAM 
     1129        time_sam = V_getCount_time("SAM")               //count time SAM 
     1130        time_drk = V_getCount_time("DRK")       //drk count time 
     1131 
     1132 
     1133        // and now loop through all of the detectors 
     1134        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1) 
     1135                detStr = StringFromList(ii, ksDetectorListAll, ";") 
     1136                Wave cor_data = V_getDetectorDataW("COR",detStr) 
     1137                Wave cor_err = V_getDetectorDataErrW("COR",detStr) 
     1138                Wave sam_data = V_getDetectorDataW("SAM",detStr) 
     1139                Wave sam_err = V_getDetectorDataErrW("SAM",detStr) 
     1140                Wave drk_data = V_getDetectorDataW("DRK",detStr) 
     1141                Wave drk_err = V_getDetectorDataErrW("DRK",detStr) 
     1142                 
     1143                 
     1144                //rescale drk to sam cnt time and then multiply by the same monitor scaling as SAM 
     1145                Duplicate/O drk_data drk_temp, drk_tmp_err 
     1146                drk_temp = drk_data*(time_sam/time_drk)*(tmonsam/savmon_sam) 
     1147                drk_tmp_err *= drk_err*(time_sam/time_drk)*(tmonsam/savmon_sam)                 //temporarily rescale the error of DRK   
     1148                 
     1149                if(temp==0) 
     1150                        DoAlert 0,"Empty Cell transmission was zero. It has been reset to one for the subtraction" 
     1151                        temp=1 
     1152                Endif 
     1153         
     1154                Duplicate/O cor_data cor1 
     1155 
     1156                //always ignore the DRK center shift 
     1157 
     1158                // ************ 
     1159                //do the subtraction,  deposit result in cor1 
     1160                fsam = 1 
     1161                fbgd = tmonsam/tmonbgd  //this should be ==1 since normalized files 
     1162                 
     1163                //correct sam for attenuators, and do the same to drk, since it was scaled to sam count time 
     1164                cor1 = fsam*sam_data/sam_AttenFactor  - drk_temp/sam_attenFactor 
     1165         
     1166        // do the error propagation piecewise    
     1167                Duplicate/O sam_err, tmp_a 
     1168                tmp_a = (sam_err/sam_attenFactor)^2 + (sam_atten_err*sam_data/sam_attenFactor^2)^2              //sig a ^2 
     1169         
     1170                cor_err = sqrt(tmp_a + drk_tmp_err^2) 
     1171 
     1172        endfor 
     1173                 
    10631174        //we're done, get out w/no error 
    1064         //set the COR_data to the result 
    1065         cor_data = cor1 
    1066         cor_data_display = cor1 
    1067  
    1068         //update COR header 
    1069         cor_text[1] = date() + " " + time()             //date + time stamp 
     1175         
     1176        //TODO -- do I update COR header? 
     1177//      cor_text[1] = date() + " " + time()             //date + time stamp 
    10701178 
    10711179        KillWaves/Z cor1,bgd_temp,noadd_bgd,drk_temp 
Note: See TracChangeset for help on using the changeset viewer.