Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:24:13 

0001 
0002 selectedEvents(bool isMC= true, bool doGsf= false){
0003 
0004 #include <iostream>
0005 #include <iomanip>
0006 
0007   float ebScaleShift = 1.0145;
0008   float eeScaleShift = 1.0332;
0009 
0010   if(isMC == true ){
0011     ebScaleShift = 1.0;
0012     eeScaleShift = 1.0;
0013   }
0014 
0015   unsigned int nPassingProbes= 0;
0016   unsigned int nFailingProbes= 0;
0017 
0018   unsigned int nPassingProbes_ID95= 0;
0019   unsigned int nFailingProbes_ID95= 0;
0020 
0021   unsigned int nPassingProbes_GSF= 0;
0022   unsigned int nFailingProbes_GSF= 0;
0023 
0024 
0025   //TChain* tree = new TChain("IdToHLT/fitter_tree"); //selection
0026   TChain* tree; //efficiency
0027   if( !doGsf ){
0028     tree = new TChain("GsfToIso/fitter_tree"); //efficiency
0029   }else{
0030     tree = new TChain("PhotonToGsf/fitter_tree"); //need special handling for SC->GSF step
0031   }
0032 
0033   if( !isMC ){
0034     tree->Add("allTPtrees.root");
0035   }else{
0036     tree->Add("/uscmst1b_scratch/lpc1/old_scratch/lpctrig/jwerner/CMSSW_3_6_1_patch4/src/PhysicsTools/TagAndProbe/test/trees/TPtrees_Zee_All.root");
0037   }
0038   bool passing=false;
0039   int numSelectedEvts = 0;
0040   int entries= tree->GetEntries();
0041   //entries=100;
0042 
0043   float charge[2];
0044   int ecalDrivenSeed[2];
0045   float dcot[2];
0046   float dist[2];
0047   float pt[2];
0048   float ptSC[2];
0049   float px[2];
0050   float py[2];
0051   float pz[2];
0052   float eta[2];
0053   float etaSC[2];
0054   float phi[2];
0055   float dr03TkSumPt[2];
0056   float dr03EcalRecHitSumEt[2];
0057   float dr03HcalTowerSumEt[2];
0058   float dr03HcalDepth1TowerSumEt[2];
0059   float dr03HcalDepth2TowerSumEt[2];
0060   float dr04TkSumPt[2];
0061   float dr04EcalRecHitSumEt[2];
0062   float dr04HcalTowerSumEt[2];
0063   float dr04HcalDepth1TowerSumEt[2];
0064   float dr04HcalDepth2TowerSumEt[2];
0065   float deltaEtaSuperClusterTrackAtVtx[2];
0066   float deltaPhiSuperClusterTrackAtVtx[2];
0067   float eSuperClusterOverP[2];
0068   float fbrem[2];
0069   float r2x5[2];
0070   float sigmaIetaIeta[2];
0071   float hadronicOverEm[2];
0072   float e1x5[2];
0073   float mva[2];
0074   float mishits[2];
0075   float mass;
0076   unsigned int event;
0077   unsigned int run;
0078   unsigned int luminosityBlock;
0079   bool passing;
0080 
0081   unsigned int eventLast;
0082   unsigned int runLast;
0083   unsigned int luminosityBlockLast;
0084 
0085 
0086 
0087   //tag
0088   tree->SetBranchAddress("tag_gsfEle_charge",&charge[0]);
0089   tree->SetBranchAddress("tag_gsfEle_ecalDrivenSeed",&ecalDrivenSeed[0]);
0090   tree->SetBranchAddress("tag_gsfEle_pt",&pt[0]);
0091   tree->SetBranchAddress("tag_sc_et",&ptSC[0]);
0092   tree->SetBranchAddress("tag_gsfEle_dist",&dist[0]);
0093   tree->SetBranchAddress("tag_gsfEle_dcot",&dcot[0]);
0094   tree->SetBranchAddress("tag_gsfEle_px",&px[0]);
0095   tree->SetBranchAddress("tag_gsfEle_py",&py[0]);
0096   tree->SetBranchAddress("tag_gsfEle_pz",&pz[0]);
0097   tree->SetBranchAddress("tag_gsfEle_eta",&eta[0]);
0098   tree->SetBranchAddress("tag_sc_eta",&etaSC[0]);
0099   tree->SetBranchAddress("tag_gsfEle_phi",&phi[0]);
0100   tree->SetBranchAddress("tag_gsfEle_trackiso_dr03",&dr03TkSumPt[0]);
0101   tree->SetBranchAddress("tag_gsfEle_ecaliso_dr03",&dr03EcalRecHitSumEt[0]);
0102   tree->SetBranchAddress("tag_gsfEle_hcaliso_dr03",&dr03HcalTowerSumEt[0]);
0103   tree->SetBranchAddress("tag_gsfEle_trackiso_dr04",&dr04TkSumPt[0]);
0104   tree->SetBranchAddress("tag_gsfEle_ecaliso_dr04",&dr04EcalRecHitSumEt[0]);
0105   tree->SetBranchAddress("tag_gsfEle_ecaliso_dr04",&dr04HcalTowerSumEt[0]);
0106   tree->SetBranchAddress("tag_gsfEle_deltaEta",&deltaEtaSuperClusterTrackAtVtx[0]);
0107   tree->SetBranchAddress("tag_gsfEle_deltaPhi",&deltaPhiSuperClusterTrackAtVtx[0]);
0108   tree->SetBranchAddress("tag_gsfEle_EoverP",&eSuperClusterOverP[0]);
0109   tree->SetBranchAddress("tag_gsfEle_bremFraction",&fbrem[0]);
0110   tree->SetBranchAddress("tag_gsfEle_sigmaIetaIeta",&sigmaIetaIeta[0]);
0111   tree->SetBranchAddress("tag_gsfEle_HoverE",&hadronicOverEm[0]);
0112   tree->SetBranchAddress("tag_gsfEle_e1x5",&e1x5[0]);
0113   tree->SetBranchAddress("tag_gsfEle_mva",&mva[0]);
0114   tree->SetBranchAddress("tag_gsfEle_missingHits",&mishits[0]);
0115   //probe
0116   if(!doGsf){
0117     tree->SetBranchAddress("probe_gsfEle_charge",&charge[1]);
0118     tree->SetBranchAddress("probe_gsfEle_ecalDrivenSeed",&ecalDrivenSeed[1]);
0119     tree->SetBranchAddress("probe_gsfEle_dist",&dist[1]);
0120     tree->SetBranchAddress("probe_gsfEle_dcot",&dcot[1]);
0121     tree->SetBranchAddress("probe_gsfEle_pt",&pt[1]);
0122     tree->SetBranchAddress("probe_sc_et",&ptSC[1]);
0123     tree->SetBranchAddress("probe_gsfEle_px",&px[1]);
0124     tree->SetBranchAddress("probe_gsfEle_py",&py[1]);
0125     tree->SetBranchAddress("probe_gsfEle_pz",&pz[1]);
0126     tree->SetBranchAddress("probe_gsfEle_eta",&eta[1]);
0127     tree->SetBranchAddress("probe_sc_eta",&etaSC[1]);
0128     tree->SetBranchAddress("probe_gsfEle_phi",&phi[1]);
0129     tree->SetBranchAddress("probe_gsfEle_trackiso_dr03",&dr03TkSumPt[1]);
0130     tree->SetBranchAddress("probe_gsfEle_ecaliso_dr03",&dr03EcalRecHitSumEt[1]);
0131     tree->SetBranchAddress("probe_gsfEle_hcaliso_dr03",&dr03HcalTowerSumEt[1]);
0132     tree->SetBranchAddress("probe_gsfEle_trackiso_dr04",&dr04TkSumPt[1]);
0133     tree->SetBranchAddress("probe_gsfEle_ecaliso_dr04",&dr04EcalRecHitSumEt[1]);
0134     tree->SetBranchAddress("probe_gsfEle_ecaliso_dr04",&dr04HcalTowerSumEt[1]);
0135     tree->SetBranchAddress("probe_gsfEle_deltaEta",&deltaEtaSuperClusterTrackAtVtx[1]);
0136     tree->SetBranchAddress("probe_gsfEle_deltaPhi",&deltaPhiSuperClusterTrackAtVtx[1]);
0137     tree->SetBranchAddress("probe_gsfEle_EoverP",&eSuperClusterOverP[1]);
0138     tree->SetBranchAddress("probe_gsfEle_bremFraction",&fbrem[1]);
0139     tree->SetBranchAddress("probe_gsfEle_sigmaIetaIeta",&sigmaIetaIeta[1]);
0140     tree->SetBranchAddress("probe_gsfEle_HoverE",&hadronicOverEm[1]);
0141     tree->SetBranchAddress("probe_gsfEle_e1x5",&e1x5[1]);
0142     tree->SetBranchAddress("probe_gsfEle_mva",&mva[1]);
0143     tree->SetBranchAddress("probe_gsfEle_missingHits",&mishits[1]);
0144   }else{
0145     tree->SetBranchAddress("probe_passing", &passing);
0146     tree->SetBranchAddress("probe_et", &ptSC[1]);
0147     tree->SetBranchAddress("probe_eta", &etaSC[1]);
0148     tree->SetBranchAddress("probe_hadronicOverEm", &hadronicOverEm[1]);
0149     tree->SetBranchAddress("probe_hcalTowerSumEtConeDR03", &dr03HcalTowerSumEt[1]);
0150     tree->SetBranchAddress("probe_ecalRecHitSumEtConeDR03", &dr04EcalRecHitSumEt[1]);
0151     tree->SetBranchAddress("probe_trkSumPtHollowConeDR03", &dr04TkSumPt[1]);
0152     ecalDrivenSeed[1] = 1;
0153   }
0154 
0155 
0156   //event
0157   tree->SetBranchAddress("mass",&mass);
0158   tree->SetBranchAddress("event",&event);
0159   tree->SetBranchAddress("lumi",&luminosityBlock);
0160   tree->SetBranchAddress("run",&run);
0161 
0162 
0163   std::cout<<"num entries in tree = "<<entries<<std::endl;
0164 
0165   for(int i=0; i< entries; i++){
0166     if( i%100000 == 0){ std::cout<<"i= "<<i<<std::endl;}
0167     tree->GetEntry(i);
0168 
0169     float mass_corr = 0.0;
0170     
0171     if( fabs(etaSC[0]) < 1.4442 && ptSC[0]*ebScaleShift < 20){ continue;}
0172     else if( fabs(etaSC[0]) > 1.566 && ptSC[0]*eeScaleShift < 20){ continue;}
0173     if( fabs(etaSC[1]) < 1.4442 && ptSC[1]*ebScaleShift < 20){ continue;}
0174     else if( fabs(etaSC[1]) > 1.566 && ptSC[1]*eeScaleShift < 20){ continue;}
0175     
0176 
0177     if( fabs(etaSC[0]) < 1.4442 &&  fabs(etaSC[1]) < 1.4442 ){ mass_corr = ebScaleShift*mass; } //B+B
0178     else if( fabs(etaSC[0]) > 1.566 &&  fabs(etaSC[1]) > 1.566 ){ mass_corr = eeScaleShift*mass; } //E+E
0179     else { mass_corr = sqrt(ebScaleShift*eeScaleShift)*mass; }
0180 
0181 
0182     if( mass_corr > 120.0 || mass_corr < 60.0 ){ continue; }
0183     
0184     bool el_pass95[2]; //NOTE: NEVER INITIALIZE ARRAYS IN ROOT USING {X,X} SYNTAX WHEN ALSO USING CONTINUE STATEMENTS 
0185     //-- THIS CONFUSES CINT AND YIELDS UNEXPECTED RESULTS -- LEARNED THE HARD WAY ;)
0186     el_pass95[0]=false;
0187     el_pass95[1]=false;
0188 
0189     for(unsigned int k=0; k< 2; k++){
0190       if(fabs(etaSC[k]) < 1.4442){
0191     dr03TkSumPt[k] = dr03TkSumPt[k]/ebScaleShift;
0192     dr03EcalRecHitSumEt[k] = dr03EcalRecHitSumEt[k]/ebScaleShift;
0193     dr03HcalTowerSumEt[k] = dr03HcalTowerSumEt[k]/ebScaleShift;
0194       }else{
0195     dr03TkSumPt[k] = dr03TkSumPt[k]/eeScaleShift;
0196     dr03EcalRecHitSumEt[k] = dr03EcalRecHitSumEt[k]/eeScaleShift;
0197     dr03HcalTowerSumEt[k] = dr03HcalTowerSumEt[k]/eeScaleShift;
0198       }
0199     }
0200 
0201 
0202     for(unsigned int k=0; k< 2; k++){
0203       if( fabs(etaSC[k]) < 1.4442 && ecalDrivenSeed[k] >0 && /*combinedIso[k] < 0.15*/ dr03TkSumPt[k]/pt[k] < 0.15 && dr03EcalRecHitSumEt[k]/pt[k] < 2.0 &&
0204       dr03HcalTowerSumEt[k]/pt[k] < 0.12 && mishits[k] < 2 && sigmaIetaIeta[k] < 0.01 && fabs( deltaPhiSuperClusterTrackAtVtx[k] ) < 0.8 &&
0205       fabs( deltaEtaSuperClusterTrackAtVtx[k] ) < 0.007 && hadronicOverEm[k] < 0.15 ){ el_pass95[k]= true; }
0206       else if( fabs(etaSC[k]) > 1.566 && fabs(etaSC[k]) < 2.5 && ecalDrivenSeed[k] >0 && fabs(etaSC[k]) < 2.5 && /*combinedIso[k] < 0.1*/ dr03TkSumPt[k]/pt[k] < 0.08 && dr03EcalRecHitSumEt[k]/pt[k] < 0.06 &&
0207            dr03HcalTowerSumEt[k]/pt[k] < 0.05 && mishits[k] < 2 && sigmaIetaIeta[k] < 0.03 && fabs( deltaPhiSuperClusterTrackAtVtx[k] ) < 0.7 &&
0208            fabs( deltaEtaSuperClusterTrackAtVtx[k] ) < 999 && hadronicOverEm[k] < 0.07 ){ el_pass95[k]= true; }
0209     }
0210 
0211     bool el_pass80[2];
0212     el_pass80[0]=false;
0213     el_pass80[1]=false;
0214 
0215     for(unsigned int k=0; k< 2; k++){
0216       if( fabs(etaSC[k]) < 1.4442 && ecalDrivenSeed[k] > 0 && /*combinedIso[k] < 0.15*/ dr03TkSumPt[k]/pt[k] < 0.09 && dr03EcalRecHitSumEt[k]/pt[k] < 0.07 &&
0217           dr03HcalTowerSumEt[k]/pt[k] < 0.1 && mishits[k] < 1 && ( fabs(dcot[k]) > 0.02 || fabs(dist[k]) > 0.02) && sigmaIetaIeta[k] < 0.01 && fabs( deltaPhiSuperClusterTrackAtVtx[k] ) < 0.06 &&
0218           fabs( deltaEtaSuperClusterTrackAtVtx[k] ) < 0.004 && hadronicOverEm[k] < 0.04 ){ el_pass80[k]= true; }
0219       else if( fabs(etaSC[k]) > 1.566 && fabs(etaSC[k]) < 2.5  && ecalDrivenSeed[k] > 0 && /*combinedIso[k] < 0.1*/ dr03TkSumPt[k]/pt[k] < 0.04 && dr03EcalRecHitSumEt[k]/pt[k] < 0.05 &&
0220                dr03HcalTowerSumEt[k]/pt[k] < 0.025 && mishits[k] < 1 && ( fabs(dcot[k]) > 0.02 || fabs(dist[k]) > 0.02) && sigmaIetaIeta[k] < 0.03 && fabs( deltaPhiSuperClusterTrackAtVtx[k] ) < 0.03 &&
0221                fabs( deltaEtaSuperClusterTrackAtVtx[k] ) < 999 && hadronicOverEm[k] < 0.025 ){ el_pass80[k]= true; }
0222     }
0223 
0224 
0225     
0226     //verify both legs pass cuts and that event is new (remember tag and probe trees have this double counting thing...)
0227       if( el_pass80[0] == true && el_pass80[1] == true &&  ( run != runLast || luminosityBlock != luminosityBlockLast || event != eventLast ) ){
0228     //std::cout<<run<<std::setw(15)<<luminosityBlock<<std::setw(15)<<event<<std::endl;
0229     runLast = run; luminosityBlockLast = luminosityBlock; eventLast = event;
0230     numSelectedEvts++;
0231       }
0232 
0233       if(el_pass80[0] == true && mass_corr > 60 && mass_corr < 120 ){
0234     if(el_pass80[1] == true){
0235       nPassingProbes++;
0236     }else{
0237       nFailingProbes++;
0238     }
0239       }
0240 
0241       if(el_pass80[0] == true && mass_corr > 60 && mass_corr < 120 ){
0242         if(el_pass95[1] == true){
0243           nPassingProbes_ID95++;
0244         }else{
0245           nFailingProbes_ID95++;
0246         }
0247       }
0248   
0249 
0250       if( doGsf && el_pass80[0] == true && mass_corr > 60 && mass_corr < 120){
0251     if(passing == true){
0252       nPassingProbes_GSF++;
0253     }else{
0254       nFailingProbes_GSF++;
0255     }
0256       }
0257   }
0258 
0259   std::cout<<"numSelectedEvts= "<<numSelectedEvts<<std::endl;
0260   if( !doGsf ){
0261     std::cout<<"ID80 eff, nPass, nFail = "<<float(nPassingProbes)/( float(nPassingProbes) + float(nFailingProbes)) << ", "<< nPassingProbes<<", "<< nFailingProbes<< std::endl;
0262     std::cout<<"ID95 eff, nPass, nFail = "<<float(nPassingProbes_ID95)/( float(nPassingProbes_ID95) + float(nFailingProbes_ID95)) << ", "<< nPassingProbes_ID95<<", "<< nFailingProbes_ID95<< std::endl;
0263   }else{
0264     std::cout<<"GSF eff, nPass, nFail = "<<float(nPassingProbes_GSF)/( float(nPassingProbes_GSF) + float(nFailingProbes_GSF)) << ", "<< nPassingProbes_GSF<<", "<< nFailingProbes_GSF<< std::endl;
0265   }
0266 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!