Advertisement
Guest User

C++ Source HEP Analysis

a guest
Aug 29th, 2014
215
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 32.72 KB | None | 0 0
  1. #include "NtupleTools3.h"
  2. #include <fstream>
  3. #include <iostream>
  4. #include <vector>
  5. #include <cmath>
  6. #include <TLorentzVector.h>
  7. #include <TCanvas.h>
  8. #include <TH1D.h>
  9. #include <TFile.h>
  10. #include <TLatex.h>
  11. #include <TMath.h>
  12.  
  13. using namespace std;
  14.  
  15. void readerSummerStd(TString list, TString outname,bool useW=true){
  16.  
  17.   TObjArray* arr = list.Tokenize(" ");
  18.   int size=arr->GetEntries();
  19.   if(size%2!=0) {
  20.     cout<<"unbalance file/weight list: "<<list<<endl;
  21.     exit(0);
  22.   }
  23.   vector<TString> files;
  24.   vector<double> weights;
  25.   for(int i=0;i<size;i+=2){
  26.     files.push_back( arr->At(i)->GetName() );
  27.     weights.push_back( atof( arr->At(i+1)->GetName() ) );
  28.   }
  29.  
  30.  
  31.   // here is the place for histograms --- learn how to determine ranges >> Learned
  32.   TH1D * hNel        = new TH1D("hNel","Number of electrons",5,0.5,5.5);
  33.   TH1D * hNel2       = new TH1D("hNel2","Number of electrons",5,0.5,5.5);
  34.   TH1D * hNel3       = new TH1D("hNel3","Number of electrons",5,0.5,5.5);
  35.   TH1D * hNel4       = new TH1D("CuthNel","Number of electrons",5,0.5,5.5);
  36.   TH1D * hNmu        = new TH1D("hNmu","Number of muons",5,0.5,5.5);
  37.   TH1D * hNmu2       = new TH1D("hNmu2","Number of muons",5,0.5,5.5);
  38.   TH1D * hNmu3       = new TH1D("hNmu3","Number of muons",5,0.5,5.5);
  39.   TH1D * hNmu4       = new TH1D("CuthNmu","Number of muons",5,0.5,5.5);
  40.   TH1D * hNlep       = new TH1D("hNlep","Number of leptons",5,0.5,5.5);
  41.   TH1D * hNlep2      = new TH1D("hNlep2","Number of leptons",5,0.5,5.5);
  42.   TH1D * hNlep3      = new TH1D("hNlep3","Number of leptons",5,0.5,5.5);
  43.   TH1D * hNlep4      = new TH1D("CuthNlep","Number of leptons",5,0.5,5.5);
  44.  
  45.   TH1D * hNj         = new TH1D("hNj","Number of jets",20,0.5,20.5);
  46.   TH1D * hNj2        = new TH1D("hNj2","Number of jets",20,0.5,20.5);
  47.   TH1D * hNj3        = new TH1D("hNj3","Number of jets",20,0.5,20.5);
  48.   TH1D * hNj4        = new TH1D("CuthNj","Number of jets",20,0.5,20.5);
  49.   TH1D * hNbj        = new TH1D("hNbj", "Number of bJets", 10, 0.5, 10.5);
  50.   TH1D * hNbj2       = new TH1D("hNbj2", "Number of bJets", 10, 0.5, 10.5);
  51.   TH1D * hNbj3       = new TH1D("hNbj3", "Number of bJets", 10, 0.5, 10.5);
  52.   TH1D * hNbj4       = new TH1D("CuthNbj", "Number of bJets", 10, 0.5, 10.5);
  53.   TH1D * hj1Pt       = new TH1D("Jet1Pt", "Jet1 Pt", 50, 0, 1000);
  54.   TH1D * hj1Pt2      = new TH1D("Jet1Pt2", "Jet1 Pt", 50, 0, 1000);
  55.   TH1D * hj1Pt3      = new TH1D("Jet1Pt3", "Jet1 Pt", 50, 0, 1000);
  56.   TH1D * hj1Pt4      = new TH1D("CutJet1Pt", "Jet1 Pt", 50, 0, 1000);
  57.   TH1D * hj2Pt       = new TH1D("Jet2Pt", "Jet2 Pt", 50, 0, 1000);
  58.   TH1D * hj2Pt2      = new TH1D("Jet2Pt2", "Jet2 Pt", 50, 0, 1000);
  59.   TH1D * hj2Pt3      = new TH1D("Jet2Pt3", "Jet2 Pt", 50, 0, 1000);
  60.   TH1D * hj2Pt4      = new TH1D("CutJet2Pt", "Jet2 Pt", 50, 0, 1000);
  61.   TH1D * hj3Pt       = new TH1D("Jet3Pt", "Jet3 Pt", 50, 0, 1000);
  62.   TH1D * hj3Pt2      = new TH1D("Jet3Pt2", "Jet3 Pt", 50, 0, 1000);
  63.   TH1D * hj3Pt3      = new TH1D("Jet3Pt3", "Jet3 Pt", 50, 0, 1000);
  64.   TH1D * hj3Pt4      = new TH1D("CutJet3Pt", "Jet3 Pt", 50, 0, 1000);
  65.   TH1D * hj4Pt       = new TH1D("Jet4Pt", "Jet4 Pt", 50, 0, 1000);
  66.   TH1D * hj4Pt2      = new TH1D("Jet4Pt2", "Jet4 Pt", 50, 0, 1000);
  67.   TH1D * hj4Pt3      = new TH1D("Jet4Pt3", "Jet4 Pt", 50, 0, 1000);
  68.   TH1D * hj4Pt4      = new TH1D("CutJet4Pt", "Jet4 Pt", 50, 0, 1000);
  69.   TH1D * hj1Eta      = new TH1D("Jet1Eta", "Jet 1 Eta", 20, -5, 5);
  70.   TH1D * hj1Eta2     = new TH1D("Jet1Eta2", "Jet 1 Eta", 20, -5, 5);
  71.   TH1D * hj1Eta3     = new TH1D("Jet1Eta3", "Jet 1 Eta", 20, -5, 5);
  72.   TH1D * hj1Eta4     = new TH1D("CutJet1Eta", "Jet 1 Eta", 20, -5, 5);
  73.   TH1D * hj2Eta      = new TH1D("Jet2Eta", "Jet 2 Eta", 20, -5, 5);
  74.   TH1D * hj2Eta2     = new TH1D("Jet2Eta2", "Jet 2 Eta", 20, -5, 5);
  75.   TH1D * hj2Eta3     = new TH1D("Jet2Eta3", "Jet 2 Eta", 20, -5, 5);
  76.   TH1D * hj2Eta4     = new TH1D("CutJet2Eta", "Jet 2 Eta", 20, -5, 5);
  77.   TH1D * hj3Eta      = new TH1D("Jet3Eta", "Jet 3 Eta", 20, -5, 5);
  78.   TH1D * hj3Eta2     = new TH1D("Jet3Eta2", "Jet 3 Eta", 20, -5, 5);
  79.   TH1D * hj3Eta3     = new TH1D("Jet3Eta3", "Jet 3 Eta", 20, -5, 5);
  80.   TH1D * hj3Eta4     = new TH1D("CutJet3Eta", "Jet 3 Eta", 20, -5, 5);
  81.   TH1D * hj4Eta      = new TH1D("Jet4Eta", "Jet 4 Eta", 20, -5, 5);
  82.   TH1D * hj4Eta2     = new TH1D("Jet4Eta2", "Jet 4 Eta", 20, -5, 5);
  83.   TH1D * hj4Eta3     = new TH1D("Jet4Eta3", "Jet 4 Eta", 20, -5, 5);
  84.   TH1D * hj4Eta4     = new TH1D("CutJet4Eta", "Jet 4 Eta", 20, -5, 5);
  85.  
  86.   TH1D * hHT         = new TH1D("hHT","HT",100,0,5000);
  87.   TH1D * hHT2        = new TH1D("hHT2","HT",100,0,5000);
  88.   TH1D * hHT3        = new TH1D("hHT3","HT",100,0,5000);
  89.   TH1D * hHT4        = new TH1D("CuthHT", "HT", 100, 0, 5000);
  90.   TH1D * hMET        = new TH1D("hMET","MET",40,0,2000);
  91.   TH1D * hMET2       = new TH1D("hMET2","MET",40,0,2000);
  92.   TH1D * hMET3       = new TH1D("hMET3","MET",40,0,2000);
  93.   TH1D * hMET4       = new TH1D("CuthMET","MET",40,0,2000);
  94.   TH1D * hST         = new TH1D("hST", "ST", 20, 0, 4000);
  95.   TH1D * hST2        = new TH1D("hST2", "ST", 20, 0, 4000);
  96.   TH1D * hST3        = new TH1D("hST3", "ST", 20, 0, 4000);
  97.   TH1D * hST4        = new TH1D("CuthST", "ST", 20, 0, 4000);
  98.   TH1D * hdPhi1      = new TH1D("hdPhi1", "#Delta#phi", 40,0,4);
  99.   TH1D * hdPhi12     = new TH1D("hdPhi12", "#Delta#phi", 40,0,4);
  100.   TH1D * hdPhi13     = new TH1D("hdPhi13", "#Delta#phi", 40,0,4);
  101.   TH1D * hdPhi14     = new TH1D("CuthdPhi1", "#Delta#phi",40,0,4);
  102.   TH1D * hdPhi2      = new TH1D("hdPhi2", "#Delta#phi ",40,0,4);
  103.   TH1D * hdPhi22     = new TH1D("hdPhi22", "#Delta#phi ",40,0,4);
  104.   TH1D * hdPhi23     = new TH1D("hdPhi23", "#Delta#phi ",40,0,4);
  105.   TH1D * hdPhi24     = new TH1D("CuthdPhi2", "#Delta#phi ",40,0,4);
  106.   TH1D * hMeff       = new TH1D("hMeff","Meff",20,0,5000);
  107.   TH1D * hMeff2      = new TH1D("hMeff2","Meff",20,0,5000);
  108.   TH1D * hMeff3      = new TH1D("hMeff3","Meff",20,0,5000);
  109.   TH1D * hMeff4      = new TH1D("CuthMeff","Cut Meff",20,0,5000);
  110.  
  111.   TH1D * hDiLMassOS  = new TH1D("hDiLMassOS","OS DiLepton mass",30,0,450);
  112.   TH1D * hDiLMassOS2 = new TH1D("hDiLMassOS2","OS DiLepton mass",30,0,450);
  113.   TH1D * hDiLMassOS3 = new TH1D("hDiLMassOS3","OS DiLepton mass",30,0,450);
  114.   TH1D * hDiLMassOS4 = new TH1D("CuthDiLMassOS","OS DiLepton mass",30,0,450);
  115.   TH1D * hDiLMassSS  = new TH1D("hDiLMassSS","SS DiLepton mass",30,0,450);
  116.   TH1D * hDiLMassSS2 = new TH1D("hDiLMassSS2","SS DiLepton mass",30,0,450);
  117.   TH1D * hDiLMassSS3 = new TH1D("hDiLMassSS3","SS DiLepton mass",30,0,450);
  118.   TH1D * hDiLMassSS4 = new TH1D("CuthDiLMassSS","SS DiLepton mass after cuts",30,0,450);
  119.   TH1D * hMETMeff    = new TH1D("hMETMeff", "MET / Meff", 100, 0, 1);
  120.   TH1D * hMETMeff2   = new TH1D("hMETMeff2", "MET / Meff", 100, 0, 1);
  121.   TH1D * hMETMeff3   = new TH1D("hMETMeff3", "MET / Meff", 100, 0, 1);
  122.   TH1D * hMETMeff4   = new TH1D("CuthMETMeff", "MET / Meff", 100, 0, 1);
  123.   TH1D * hLeadLepPt  = new TH1D("hLeadLepPt","1 Leading Lepton Pt ",60,0,600);
  124.   TH1D * hLeadLepPt2 = new TH1D("hLeadLepPt2","1 Leading Lepton Pt ",60,0,600);
  125.   TH1D * hLeadLepPt3 = new TH1D("hLeadLepPt3","1 Leading Lepton Pt ",60,0,600);
  126.   TH1D * hLeadLepPt4 = new TH1D("CuthLeadLepPt","1 Leading Lepton Pt",60,0,600);
  127.   TH1D * hLeadLep2Pt = new TH1D("hLeadLep2Pt","2 Leading Lepton Pt",60,0,600);
  128.   TH1D * hLeadLep2Pt2= new TH1D("hLeadLep2Pt2","2 Leading Lepton Pt",60,0,600);
  129.   TH1D * hLeadLep2Pt3= new TH1D("hLeadLep2Pt3","2 Leading Lepton Pt",60,0,600);
  130.   TH1D * hLeadLep2Pt4= new TH1D("CuthLeadLep2Pt","2 Leading Lepton Pt",60,0,600);
  131.   TH1D * hLeadLep3Pt = new TH1D("hLeadLep3Pt","3 Leading Lepton Pt",60,0,600);
  132.   TH1D * hLeadLep3Pt2= new TH1D("hLeadLep3Pt2","3 Leading Lepton Pt",60,0,600);
  133.   TH1D * hLeadLep3Pt3= new TH1D("hLeadLep3Pt3","3 Leading Lepton Pt",60,0,600);
  134.   TH1D * hLeadLep3Pt4= new TH1D("CuthLeadLep3Pt","3 Leading Lepton Pt",60,0,600);
  135.   TH1D * hLeadLepEta  = new TH1D("hLeadLepEta","1 Leading Lepton Eta ",20, -5, 5);
  136.   TH1D * hLeadLepEta2 = new TH1D("hLeadLepEta2","1 Leading Lepton Eta ",20, -5, 5);
  137.   TH1D * hLeadLepEta3 = new TH1D("hLeadLepEta3","1 Leading Lepton Eta ",20, -5, 5);
  138.   TH1D * hLeadLepEta4 = new TH1D("CuthLeadLepEta","1 Leading Lepton Eta",20, -5, 5);
  139.   TH1D * hLeadLep2Eta = new TH1D("hLeadLep2Eta","2 Leading Lepton Eta",20, -5, 5);
  140.   TH1D * hLeadLep2Eta2= new TH1D("hLeadLep2Eta2","2 Leading Lepton Eta",20, -5, 5);
  141.   TH1D * hLeadLep2Eta3= new TH1D("hLeadLep2Eta3","2 Leading Lepton Eta",20, -5, 5);
  142.   TH1D * hLeadLep2Eta4= new TH1D("CuthLeadLep2Eta","2 Leading Lepton Eta",20, -5, 5);
  143.   TH1D * hLeadLep3Eta = new TH1D("hLeadLep3Eta","3 Leading Lepton Eta",20, -5, 5);
  144.   TH1D * hLeadLep3Eta2= new TH1D("hLeadLep3Eta2","3 Leading Lepton Eta",20, -5, 5);
  145.   TH1D * hLeadLep3Eta3= new TH1D("hLeadLep3Eta3","3 Leading Lepton Eta",20, -5, 5);
  146.   TH1D * hLeadLep3Eta4= new TH1D("CuthLeadLep3Eta","3 Leading Lepton Eta",20, -5, 5);
  147.   TH1D * hLeadLepIso  = new TH1D("hLeadLepIso","1 Leading Lepton Iso ",20, 0, 1);
  148.   TH1D * hLeadLepIso2 = new TH1D("hLeadLepIso2","1 Leading Lepton Iso ",20, 0, 1);
  149.   TH1D * hLeadLepIso3 = new TH1D("hLeadLepIso3","1 Leading Lepton Iso ",20, 0, 1);
  150.   TH1D * hLeadLepIso4 = new TH1D("CuthLeadLepIso","1 Leading Lepton Iso",20, 0, 1);
  151.   TH1D * hLeadLep2Iso = new TH1D("hLeadLep2Iso","2 Leading Lepton Iso",20, 0, 1);
  152.   TH1D * hLeadLep2Iso2= new TH1D("hLeadLep2Iso2","2 Leading Lepton Iso",20, 0, 1);
  153.   TH1D * hLeadLep2Iso3= new TH1D("hLeadLep2Iso3","2 Leading Lepton Iso",20, 0, 1);
  154.   TH1D * hLeadLep2Iso4= new TH1D("CuthLeadLep2Iso","2 Leading Lepton Iso",20, 0, 1);
  155.   TH1D * hLeadLep3Iso = new TH1D("hLeadLep3Iso","3 Leading Lepton Iso",20, 0, 1);
  156.   TH1D * hLeadLep3Iso2= new TH1D("hLeadLep3Iso2","3 Leading Lepton Iso",20, 0, 1);
  157.   TH1D * hLeadLep3Iso3= new TH1D("hLeadLep3Iso3","3 Leading Lepton Iso",20, 0, 1);
  158.   TH1D * hLeadLep3Iso4= new TH1D("CuthLeadLep3Iso","3 Leading Lepton Iso",20, 0, 1);
  159.  
  160.   TH1D * hHasZ      = new TH1D("hHasZ", "ZZ count", 5, -0.5, 4.5);
  161.   TH1D * hHasZ2     = new TH1D("hHasZ2", "ZZ count", 5, -0.5, 4.5);
  162.   TH1D * hHasZ3     = new TH1D("hHasZ3", "ZZ count", 5, -0.5, 4.5);
  163.   TH1D * hHasZ4     = new TH1D("CuthHasZ", "ZZ count", 5, -0.5, 4.5);
  164.   //TH1D * hMT      = new TH1D("hMT", "Transverse Mass", 50, 0, 500);
  165.   //TH1D * hMT2     = new TH1D("hMT2", "Transverse Mass", 50, 0, 500);
  166.   //TH1D * hMT3     = new TH1D("hMT3", "Transverse Mass", 50, 0, 500);
  167.   //TH1D * hMT4     = new TH1D("CuthMT", "Transverse Mass", 50, 0, 500);
  168.  
  169.   EasyChain* tree   = new EasyChain("delphTree");
  170.  
  171.   for(unsigned i=0;i<files.size();i++){
  172.     tree->AddSmartW(files[i],weights[i]);
  173.     cout<<"add: "<<files[i]<<" "<<weights[i]<<endl;
  174.   }
  175.  
  176.   int Nevents=tree->GetEntries();
  177.   cout<<">>>>>>>>>>>>>>>>>>>>>>> total number of events:\t" << Nevents <<endl;
  178.  
  179.   // CutFlow variables:
  180.   const int CutNumb = 13;
  181.   const char * CutList[CutNumb] = {"PreSel","Elec", "Jets$\\geq$4", "Bjets$\\geq$3", "MET$>$100", "MET$>$150", "MET$>$200", "MET$>$250", "MET$>$300", "MET$>$350", "MET$>$400", "MET$>$450", "MET$>$500"};
  182.   //const char * CutList[CutNumb] = {"Before","Electrons","SMZ", "Jets>4", "Bjets=0", "Iso<0.20", "Iso<0.15", "Iso<0.10", "dPhi", "MET>1000", "HT>1000\\_MET>250", "HT>1000\\_MET>500", "MET/Meff>0.2"};
  183.   const int CutFlow = 5;
  184.   const char * Cutter[CutFlow] = {"Bef. sel.", "3 Elec. req.", "Jets > 4", "Bjets > 3", "MET>250"};
  185.   //const char * Cutter[CutFlow] = {"Bef. sel.", "3 Elec. req.", "No Z", "4 Jets", "0 Bjets", "Iso > 0.15", "dPhi > 0.3", "MET>500", "HT>1000", "MET/Meff>0.2"};
  186.   // Counter start here:
  187.     double CFCounter[CutNumb];
  188.     int   iCFCounter[CutNumb];
  189.     int  iCFC[CutFlow];
  190.     double CFC[CutFlow];
  191.  
  192.    
  193.     for (int i=0;i < CutNumb ;i++){
  194.         CFCounter[i] = 0;
  195.         iCFCounter[i] = 0;
  196.     }
  197.    
  198.     for (int i=0;i < CutFlow ;i++){
  199.         CFC[i] = 0;
  200.         iCFC[i] = 0;
  201.     }
  202.     TH1D* CF  = new TH1D("CF","Cut Flow",CutFlow,0.5,CutFlow+0.5);
  203.     // Label bins in CF histogram
  204.     for(int cj = 0; cj < CutFlow; cj++)
  205.       CF->GetXaxis()->SetBinLabel(cj+1,Cutter[cj]);
  206.     ////////////////////////////////////////////
  207.     //LOOP STARTS HERE
  208.     ////////////////////////////////////////////
  209.     for(int entry=0; entry < Nevents; entry+=1){
  210.     ////////////////////////////////////////////
  211.     ////////////////////////////////////////////
  212.       progressT();
  213.       double fw = tree->GetEntryW(entry); // the applied with AddSmartW for the current file/dir
  214.      
  215.       double EvWeight = 1;
  216.       if(useW) EvWeight = tree->Get(EvWeight,"EventWeight");
  217.       EvWeight *= fw * 1000; // 1000 for pb->fb
  218.      
  219.      
  220.       // New Vector
  221.       vector<TLorentzVector> Test;// = new vector<TLorentzVector>;
  222.       vector<TLorentzVector> Lepton;// = new vector<TLorentzVector>;
  223.       vector<TLorentzVector> GoodJets;
  224.       vector<TLorentzVector> CleanJets;
  225.       vector<TLorentzVector> GoodElectrons;
  226.       vector<TLorentzVector> GoodMuons;
  227.       vector<TLorentzVector> GoodLepton;
  228.       // From ntupler
  229.       vector<TLorentzVector> &Electrons = tree->Get(&Electrons,"Electrons");                
  230.       vector<TLorentzVector> &Muons = tree->Get(&Muons,"Muons");
  231.       vector<TLorentzVector> &Jets = tree->Get(&Jets,"Jets");
  232.       TLorentzVector lv_sum;
  233.  
  234.       vector<int> &ElCh = tree->Get(&ElCh,"ElectronCh");
  235.       vector<int> &MuCh = tree->Get(&MuCh,"MuonCh");
  236.       vector<int> &JetB = tree->Get(&JetB,"JetB");
  237.       vector<int> LepFlavour;
  238.       vector<int> CleanJetsInd;
  239.       vector<int> GoodJetsInd;
  240.       vector<int> GoodElectronsInd;
  241.       vector<int> GoodMuonsInd;
  242.  
  243.       vector<double> &ElIso = tree->Get(&ElIso,"ElectronIso");
  244.       vector<double> &MuIso = tree->Get(&MuIso,"MuonIso");
  245.       vector<double> &JetMETdPhi = tree->Get(&JetMETdPhi,"JetMETdPhi");    
  246.       vector<double> LeptonIso;
  247.       vector<double> GoodLeptonIso;
  248.  
  249.       int Nel = tree->Get(Nel,   "Nel"); // Get number of electrons from root file
  250.       int Nmu = tree->Get(Nmu,   "Nmu"); // Get number of muons fromm root file
  251.       double MET  = tree->Get(MET,"DelphMET");
  252.       double HT40   = tree->Get(HT40,"HT40");
  253.       double HT = 0;
  254.       int hasSS = 0; // has di-leptons
  255.       int hasZ = 0; // has Z like pair
  256.       int hasOS = 0;
  257.      
  258.  
  259.       // The derived variables are defined below
  260.       int Nlept = Nel + Nmu;
  261.      
  262.       bool  vetoE[100];for(int ie=0;ie<Nel;ie++)vetoE[ie]=0;
  263.       bool  vetoM[100];for(int im=0;im<Nmu;im++)vetoM[im]=0;
  264.      
  265.     // EWKino leptons
  266.     for(int ie=0;ie<Nel;ie++) {
  267.         // Previous line vetoE[ie]|=(Electrons[ie].Pt()<20)||(fabs(Electrons[ie].Eta())>2.4);
  268.         vetoE[ie]|=(Electrons[ie].Pt()<10.)||(fabs(Electrons[ie].Eta())>2.4) || (ElIso[ie]  > 0.15); // Edited in accordance with p.3 on article
  269.         if(!vetoE[ie])
  270.         {
  271.             GoodElectrons.push_back(Electrons[ie]);
  272.             GoodElectronsInd.push_back(ie);
  273.             GoodLepton.push_back(Electrons[ie]);
  274.             GoodLeptonIso.push_back(ElIso[ie]);
  275.         }
  276.     }
  277.     for(int im=0;im<Nmu;im++) {
  278.         // Previous line vetoM[im]|=(Muons[im].Pt()<20)||(fabs(Muons[im].Eta())>2.4);
  279.         vetoM[im]|=(Muons[im].Pt()<10.)||(fabs(Muons[im].Eta())>2.4) || (MuIso[im] > 0.15); // Edited in accordance with p.3 on article iso
  280.         if(!vetoM[im])
  281.         {
  282.             GoodMuons.push_back(Muons[im]);
  283.             GoodMuonsInd.push_back(im);
  284.             GoodLepton.push_back(Muons[im]);
  285.             GoodLeptonIso.push_back(MuIso[im]);
  286.         }
  287.     }
  288.     int NGel, NGmu, NGlept;
  289.     NGel = GoodElectrons.size();
  290.     NGmu = GoodMuons.size();
  291.     NGlept = NGel + NGmu;
  292.     vector<bool> veto;
  293.     // Look at the particle handbook for numbers 11, 13
  294.     for(int elid = 0; elid < Nel; elid++)            
  295.     {
  296.         veto.push_back(vetoE[elid]);
  297.         Lepton.push_back(Electrons[elid]);
  298.         LepFlavour.push_back(11*ElCh[elid]);
  299.         LeptonIso.push_back(ElIso[elid]);
  300.     }
  301.     for(int muid = 0; muid < Nmu; muid++)
  302.     {
  303.         veto.push_back(vetoM[muid]);
  304.         Lepton.push_back(Muons[muid]);
  305.         LepFlavour.push_back(13*MuCh[muid]);
  306.         LeptonIso.push_back(MuIso[muid]);
  307.     }
  308.     ////////////////////////////////////
  309.     ////////////// VETO ////////////////
  310.     ////////////////////////////////////
  311.     bool abort = false;
  312.     for(int lid = 0; lid < Nlept; lid++)
  313.     {
  314.         if(veto[lid])
  315.             abort = true;
  316.     }
  317.     //cout << " \v NEW SETS" << endl;
  318.     //cout << "Lepton.size(): " << Lepton.size() << "\n" << "GoodLepton.size(): " << GoodLepton.size() << endl;
  319.     if(abort)
  320.     {
  321.         //cout << "ABORTED!!!" << endl;
  322.         continue;
  323.     }
  324.     // 0. CF presel
  325.     CFCounter[0]+= EvWeight;
  326.     iCFCounter[0]++;
  327.     CFC[0]+= EvWeight;
  328.     iCFC[0]++;
  329.     //cout << "Lepton.size(): " << Lepton.size() << "\n" << "GoodLepton.size(): " << GoodLepton.size() << endl;
  330.     ////
  331.     ////
  332.     for(int lid = 0; lid < Nlept; lid++)
  333.     {
  334.         if(veto[lid]) continue;
  335.         for(int sid = 0; (sid < Nlept) && (lid != sid); sid++)
  336.         {
  337.             if(veto[sid]) continue;
  338.             int sumflvr = LepFlavour[lid] + LepFlavour[sid]; // sum of lepton flavour, possible comb: SS: 22,24,26 OS: 0,2
  339.             if(abs(sumflvr) > 0)
  340.             {
  341.                 hasSS++;
  342.                 double m=lv_sum.M();
  343.                 if(m>0)hDiLMassSS->Fill(m,EvWeight);
  344.             }
  345.  
  346.             if(sumflvr == 0) // may be Z decay
  347.             {
  348.                    lv_sum = Lepton[lid] + Lepton[sid];//lv_lep1 + lv_lep2;
  349.                    hasOS++; // NEW ADDITION AS OF 24.8.2014
  350.                    if(abs(lv_sum.M()-90.) < 15.) hasZ++; // For book keeping purposes
  351.                    if(lv_sum.M() > 0) hDiLMassOS->Fill(lv_sum.M(),EvWeight);
  352.             }
  353.         }
  354.     }
  355.        
  356.     // These need to be defined after preselection in order to fill our histograms before the cuts are applied
  357.     /////////////////////////
  358.     // GOOD JETS SELECTION //
  359.     /////////////////////////
  360.     int Nbjet = 0;
  361.     bool checkCLEAN = true;
  362.     for(int ij = 0; ij < (int)Jets.size(); ij++)
  363.     {
  364.         if(Jets[ij].Pt() > 40. && fabs(Jets[ij].Eta()) < 2.4)
  365.         {
  366.             GoodJets.push_back(Jets[ij]);
  367.             GoodJetsInd.push_back(ij);
  368.         }
  369.     }
  370.     //////////////////////////
  371.     // CLEAN JETS SELECTION //
  372.     //////////////////////////
  373.     for(int gj = 0; gj < (int)GoodJets.size(); gj++)
  374.     {
  375.         checkCLEAN = true;
  376.         for(int l = 0; l < NGlept; l++)
  377.             if(GoodJets[gj].DeltaR(GoodLepton[l]) <= 0.4) checkCLEAN = false;
  378.         if(checkCLEAN)
  379.         {
  380.             CleanJets.push_back(GoodJets[gj]);
  381.             CleanJetsInd.push_back(GoodJetsInd[gj]);
  382.         }
  383.     }
  384.     ////////////////////////
  385.     //////// NEW HT ////////
  386.     ///////////////////////
  387.     for(int cj = 0; cj < (int)CleanJets.size(); cj++)
  388.         HT += CleanJets[cj].Pt();
  389.     //cout << "New HT:  " << HT << "\n  Old HT: " << HT40 << endl;
  390.     //////////////////
  391.     // BJET TAGGING //
  392.     //////////////////
  393.     for(int i=0; i < (int)CleanJets.size();i++) // Two lines of for loop
  394.         if((JetB[CleanJetsInd[i]] & 2) && (CleanJets[i].Pt() > 50.) && (fabs(CleanJets[i].Eta()) < 1.8)) Nbjet++; // Batool
  395.     //////////////////
  396.     double LepPtAll = 0;
  397.     for(int lid = 0; lid < NGlept; lid++)
  398.         LepPtAll += GoodLepton[lid].Pt();
  399.     double ST = MET + LepPtAll;
  400.     double MET_Phi = tree->Get(MET_Phi,"DelphMET_Phi");
  401.     double Meff = MET + HT;
  402.     //double MT = 0;
  403.     //if(Nlept > 2) MT = sqrt(2*MET*Lepton[2].Pt()*(1-cos(fabs(MET_Phi - Lepton[2].Phi()))));
  404.     //////////////////////////////////////////////////////
  405.     //HISTOGRAMS BEFORE THE CUTS
  406.     //////////////////////////////////////////////////////
  407.     //cout << (MET / Meff) << endl;
  408.     hNel->Fill(NGel,EvWeight);
  409.     hNmu->Fill(NGmu,EvWeight);
  410.     hNlep->Fill(NGlept, EvWeight);
  411.     hHT->Fill(HT, EvWeight);
  412.     hMET->Fill(MET, EvWeight);
  413.     hNj->Fill(CleanJets.size(), EvWeight);
  414.     hMETMeff->Fill((MET / Meff), EvWeight);
  415.     // The following checks are necessary to avoid segmentation faults
  416.     if(CleanJets.size() > 0) hdPhi1->Fill(JetMETdPhi[CleanJetsInd[0]], EvWeight);
  417.     if(CleanJets.size() > 1) hdPhi2->Fill(JetMETdPhi[CleanJetsInd[1]], EvWeight);
  418.     if(CleanJets.size() > 0) hj1Pt->Fill(CleanJets[0].Pt(), EvWeight);
  419.     if(CleanJets.size() > 1) hj2Pt->Fill(CleanJets[1].Pt(), EvWeight);
  420.     if(CleanJets.size() > 2) hj3Pt->Fill(CleanJets[2].Pt(), EvWeight);
  421.     if(CleanJets.size() > 3) hj4Pt->Fill(CleanJets[3].Pt(), EvWeight);
  422.     hNbj->Fill(Nbjet, EvWeight);
  423.     hHasZ->Fill(hasZ);
  424.     if(GoodLepton.size() > 0 && GoodLepton[0].Pt() > 0) hLeadLepPt->Fill(GoodLepton[0].Pt(), EvWeight);
  425.     if(GoodLepton.size() > 1 && GoodLepton[1].Pt() > 0) hLeadLep2Pt->Fill(GoodLepton[1].Pt(), EvWeight);
  426.     if(GoodLepton.size() > 2 && GoodLepton[2].Pt() > 0) hLeadLep3Pt->Fill(GoodLepton[2].Pt(), EvWeight);
  427.     hMeff->Fill(Meff, EvWeight);
  428.     if(CleanJets.size() > 0) hj1Eta->Fill(CleanJets[0].Eta(), EvWeight);
  429.     if(CleanJets.size() > 1) hj2Eta->Fill(CleanJets[1].Eta(), EvWeight);
  430.     if(CleanJets.size() > 2) hj3Eta->Fill(CleanJets[2].Eta(), EvWeight);
  431.     if(CleanJets.size() > 3) hj4Eta->Fill(CleanJets[3].Eta(), EvWeight);
  432.     if(GoodLepton.size() > 0) hLeadLepEta->Fill(GoodLepton[0].Eta(), EvWeight);
  433.     if(GoodLepton.size() > 1) hLeadLep2Eta->Fill(GoodLepton[1].Eta(), EvWeight);
  434.     if(GoodLepton.size() > 2) hLeadLep3Eta->Fill(GoodLepton[2].Eta(), EvWeight);
  435.     if(GoodLepton.size() > 0) hLeadLepIso->Fill(GoodLeptonIso[0], EvWeight);
  436.     if(GoodLepton.size() > 1) hLeadLep2Iso->Fill(GoodLeptonIso[1], EvWeight);
  437.     if(GoodLepton.size() > 2) hLeadLep3Iso->Fill(GoodLeptonIso[2], EvWeight);
  438.     hST->Fill(ST, EvWeight);
  439.     //if(MT > 0) hMT->Fill(MT, EvWeight);
  440.     //////////////////////////////////////////////////////
  441.     //////////////////////
  442.     // Start of CF //
  443.     //////////////////////--------->>>>>>> here is the place, which needs to be filled by the students, Jack or Ongun, or both ---
  444.     // 1. nElectrons
  445.     if((NGel != 3) || (NGmu != 0))continue; // Select three electrons only
  446.     if((GoodLepton[0].Pt() < 25.) || (GoodLepton[1].Pt() < 15.) || (GoodLepton[2].Pt() < 10.)) continue;
  447.     CFCounter[1]+= EvWeight;
  448.     iCFCounter[1]++;
  449.     CFC[1]+= EvWeight;
  450.     iCFC[1]++;
  451.  
  452.     // 2. Z-WINDOW and 15 GeV
  453.     //if(lv_sum.M() < 15) continue;
  454.     //if(hasZ > 0) continue; // For OS lepton pair
  455.     //CFCounter[2]+= EvWeight;
  456.     //iCFCounter[2]++; 
  457.     //CFC[2]+= EvWeight;
  458.     //iCFC[2]++;   
  459.    
  460.     // 3. nJets >= 4
  461.     //The construct is as follows for putting cuts on leading jets IMO(Ongun)
  462.     if(CleanJets.size() < 4) continue;
  463.     CFCounter[2]+= EvWeight;
  464.     iCFCounter[2]++;   
  465.     CFC[2]+= EvWeight;
  466.     iCFC[2]++; 
  467.  
  468.     // 4. NBjets => 3
  469.     if(Nbjet < 3) continue;
  470.     CFCounter[3]+= EvWeight;
  471.     iCFCounter[3]++;
  472.     CFC[3]+= EvWeight;
  473.     iCFC[3]++;
  474.     //////////////////////////////////////////////////////
  475.     //////////////////////////////////////////////////////
  476.     //////// HISTOGRAMS IN THE MIDDLE OF CUTS SET 2 //////
  477.     hNlep2->Fill(NGlept, EvWeight);
  478.     hNel2->Fill(NGel,EvWeight);
  479.     hNmu2->Fill(NGmu,EvWeight);
  480.     hHT2->Fill(HT,EvWeight);
  481.     hMET2->Fill(MET,EvWeight);
  482.     hNj2->Fill(CleanJets.size(), EvWeight);
  483.     hdPhi12->Fill(JetMETdPhi[CleanJetsInd[0]], EvWeight);
  484.     hdPhi22->Fill(JetMETdPhi[CleanJetsInd[1]], EvWeight);
  485.     hNbj2->Fill(Nbjet, EvWeight);
  486.     hMETMeff2->Fill((MET / Meff), EvWeight);
  487.     hj1Pt2->Fill(CleanJets[0].Pt(), EvWeight);
  488.     hj2Pt2->Fill(CleanJets[1].Pt(), EvWeight);
  489.     hj3Pt2->Fill(CleanJets[2].Pt(), EvWeight);
  490.     hj4Pt2->Fill(CleanJets[3].Pt(), EvWeight);
  491.     if((hasSS > 0) && (lv_sum.M()  > 0)) hDiLMassSS2->Fill(lv_sum.M(),EvWeight);
  492.     if(lv_sum.M() > 0) hDiLMassOS2->Fill(lv_sum.M(),EvWeight);    
  493.     hHasZ2->Fill(hasZ);
  494.     if(GoodLepton.size() > 0 && GoodLepton[0].Pt() > 0) hLeadLepPt2->Fill(GoodLepton[0].Pt(), EvWeight);
  495.     if(GoodLepton.size() > 1 && GoodLepton[1].Pt() > 0) hLeadLep2Pt2->Fill(GoodLepton[1].Pt(), EvWeight);
  496.     if(GoodLepton.size() > 2 && GoodLepton[2].Pt() > 0) hLeadLep3Pt2->Fill(GoodLepton[2].Pt(), EvWeight);
  497.     hMeff2->Fill(Meff, EvWeight);
  498.     if(CleanJets.size() > 0) hj1Eta2->Fill(CleanJets[0].Eta(), EvWeight);
  499.     if(CleanJets.size() > 1) hj2Eta2->Fill(CleanJets[1].Eta(), EvWeight);
  500.     if(CleanJets.size() > 2) hj3Eta2->Fill(CleanJets[2].Eta(), EvWeight);
  501.     if(CleanJets.size() > 3) hj4Eta2->Fill(CleanJets[3].Eta(), EvWeight);
  502.     if(GoodLepton.size() > 0) hLeadLepEta2->Fill(GoodLepton[0].Eta(), EvWeight);
  503.     if(GoodLepton.size() > 1) hLeadLep2Eta2->Fill(GoodLepton[1].Eta(), EvWeight);
  504.     if(GoodLepton.size() > 2) hLeadLep3Eta2->Fill(GoodLepton[2].Eta(), EvWeight);
  505.     if(GoodLepton.size() > 0) hLeadLepIso2->Fill(GoodLeptonIso[0], EvWeight);
  506.     if(GoodLepton.size() > 1) hLeadLep2Iso2->Fill(GoodLeptonIso[1], EvWeight);
  507.     if(GoodLepton.size() > 2) hLeadLep3Iso2->Fill(GoodLeptonIso[2], EvWeight);
  508.     hST2->Fill(ST, EvWeight);
  509.     //if(MT > 0) hMT2->Fill(MT, EvWeight);
  510.     //////////////////////////////////////////////////////
  511.     //////////////////////////////////////////////////////
  512.    
  513.     //5. dPhi > 0.3
  514.     //double dPhi1 = TMath::Min(JetMETdPhi[CleanJetsInd[0]],JetMETdPhi[CleanJetsInd[1]]);
  515.     //double dPhi2 = TMath::Min(JetMETdPhi[CleanJetsInd[2]],JetMETdPhi[CleanJetsInd[3]]);
  516.     //double dPhi = TMath::Min(dPhi1, dPhi2);
  517.     //if(dPhi < 0.3) continue;
  518.     //CFCounter[5]+= EvWeight;
  519.     //iCFCounter[5]++;
  520.     //CFC[5]+= EvWeight;
  521.     //iCFC[5]++;
  522.    
  523.     //////////////////////////////////////////////////////
  524.     //////// HISTOGRAMS IN THE MIDDLE OF CUTS SET 3 //////
  525.     //////////////////////////////////////////////////////
  526.     hNlep3->Fill(NGlept, EvWeight);
  527.     hNel3->Fill(NGel,EvWeight);
  528.     hNmu3->Fill(NGmu,EvWeight);
  529.     hHT3->Fill(HT,EvWeight);
  530.     hMET3->Fill(MET,EvWeight);
  531.     hNj3->Fill(CleanJets.size(), EvWeight);
  532.     hdPhi13->Fill(JetMETdPhi[CleanJetsInd[0]], EvWeight);
  533.     hdPhi23->Fill(JetMETdPhi[CleanJetsInd[1]], EvWeight);
  534.     hNbj3->Fill(Nbjet, EvWeight);
  535.     hMETMeff3->Fill((MET / Meff), EvWeight);
  536.     hj1Pt3->Fill(CleanJets[0].Pt(), EvWeight);
  537.     hj2Pt3->Fill(CleanJets[1].Pt(), EvWeight);
  538.     hj3Pt3->Fill(CleanJets[2].Pt(), EvWeight);
  539.     hj4Pt3->Fill(CleanJets[3].Pt(), EvWeight);
  540.     if((hasSS > 0) && (lv_sum.M()  > 0)) hDiLMassSS3->Fill(lv_sum.M(),EvWeight);
  541.     if(lv_sum.M() > 0) hDiLMassOS3->Fill(lv_sum.M(),EvWeight);    
  542.     hHasZ3->Fill(hasZ);
  543.     if(GoodLepton.size() > 0 && GoodLepton[0].Pt() > 0) hLeadLepPt3->Fill(GoodLepton[0].Pt(), EvWeight);
  544.     if(GoodLepton.size() > 1 && GoodLepton[1].Pt() > 0) hLeadLep2Pt3->Fill(GoodLepton[1].Pt(), EvWeight);
  545.     if(GoodLepton.size() > 2 && GoodLepton[2].Pt() > 0) hLeadLep3Pt3->Fill(GoodLepton[2].Pt(), EvWeight);
  546.     hMeff3->Fill(Meff, EvWeight);
  547.     if(CleanJets.size() > 0) hj1Eta3->Fill(CleanJets[0].Eta(), EvWeight);
  548.     if(CleanJets.size() > 1) hj2Eta3->Fill(CleanJets[1].Eta(), EvWeight);
  549.     if(CleanJets.size() > 2) hj3Eta3->Fill(CleanJets[2].Eta(), EvWeight);
  550.     if(CleanJets.size() > 3) hj4Eta3->Fill(CleanJets[3].Eta(), EvWeight);
  551.     if(GoodLepton.size() > 0) hLeadLepEta3->Fill(GoodLepton[0].Eta(), EvWeight);
  552.     if(GoodLepton.size() > 1) hLeadLep2Eta3->Fill(GoodLepton[1].Eta(), EvWeight);
  553.     if(GoodLepton.size() > 2) hLeadLep3Eta3->Fill(GoodLepton[2].Eta(), EvWeight);
  554.     if(GoodLepton.size() > 0) hLeadLepIso3->Fill(GoodLeptonIso[0], EvWeight);
  555.     if(GoodLepton.size() > 1) hLeadLep2Iso3->Fill(GoodLeptonIso[1], EvWeight);
  556.     if(GoodLepton.size() > 2) hLeadLep3Iso3->Fill(GoodLeptonIso[2], EvWeight);
  557.     hST3->Fill(ST, EvWeight);
  558.     //if(MT > 0) hMT3->Fill(MT, EvWeight);
  559.  
  560.     // 6. HT > 1000
  561.     //if(HT < 500) continue;
  562.     ////if( MET < 500) continue;
  563.     //CFCounter[6]+= EvWeight;
  564.     //iCFCounter[6]++;
  565.     //CFC[6]+= EvWeight;
  566.     //iCFC[6]++;
  567.    
  568.     if(MET > 100)
  569.     {
  570.         CFCounter[4] += EvWeight;
  571.         iCFCounter[4]++;
  572.     }
  573.     if(MET > 150)
  574.     {
  575.         CFCounter[5] += EvWeight;
  576.         iCFCounter[5]++;
  577.     }
  578.  
  579.     if(MET > 200)
  580.     {
  581.         CFCounter[6] += EvWeight;
  582.         iCFCounter[6]++;
  583.     }
  584.     // 7. MET > 500 HT
  585.     if( MET < 250) continue;
  586.     //if(HT < 1000) continue;
  587.     CFCounter[7]+= EvWeight;
  588.     iCFCounter[7]++;
  589.     CFC[4]+= EvWeight;
  590.     iCFC[4]++;
  591.  
  592.     if(MET > 300)  
  593.     {
  594.         CFCounter[8] += EvWeight;
  595.         iCFCounter[8]++;
  596.     }
  597.     if(MET > 350)
  598.     {
  599.         CFCounter[9] += EvWeight;
  600.         iCFCounter[9]++;
  601.     }
  602.     if(MET > 400)
  603.     {
  604.         CFCounter[10] += EvWeight;
  605.         iCFCounter[10]++;
  606.     }
  607.     if(MET > 450)
  608.     {
  609.         CFCounter[11] += EvWeight;
  610.         iCFCounter[11]++;
  611.     }
  612.     if(MET > 500)
  613.     {
  614.         CFCounter[12] += EvWeight;
  615.         iCFCounter[12]++;
  616.     }
  617.  
  618.     //////////////////////////////////////////////////////
  619.     //////////////////////////////////////////////////////
  620.     ///////////////////////
  621.     //  END of Selection //
  622.     ///////////////////////
  623.    
  624.     ///////////////////////////
  625.     // HISTOGRAMS AFTER CUTS //
  626.     ///////////////////////////
  627.     hNlep4->Fill(NGlept, EvWeight);
  628.     hNel4->Fill(NGel,EvWeight);
  629.     hNmu4->Fill(NGmu,EvWeight);
  630.     hHT4->Fill(HT,EvWeight);
  631.     hMET4->Fill(MET,EvWeight);
  632.     hNj4->Fill(CleanJets.size(), EvWeight);
  633.     hdPhi14->Fill(JetMETdPhi[CleanJetsInd[0]], EvWeight);
  634.     hdPhi24->Fill(JetMETdPhi[CleanJetsInd[1]], EvWeight);
  635.     hNbj4->Fill(Nbjet, EvWeight);
  636.     hMETMeff4->Fill((MET / Meff), EvWeight);
  637.     hj1Pt4->Fill(CleanJets[0].Pt(), EvWeight);
  638.     hj2Pt4->Fill(CleanJets[1].Pt(), EvWeight);
  639.     hj3Pt4->Fill(CleanJets[2].Pt(), EvWeight);
  640.     hj4Pt4->Fill(CleanJets[3].Pt(), EvWeight);
  641.     if((hasSS > 0) && (lv_sum.M()  > 0)) hDiLMassSS4->Fill(lv_sum.M(),EvWeight);
  642.     if(lv_sum.M() > 0) hDiLMassOS4->Fill(lv_sum.M(),EvWeight);    
  643.     hHasZ4->Fill(hasZ);
  644.     if(GoodLepton.size() > 0 && GoodLepton[0].Pt() > 0) hLeadLepPt4->Fill(GoodLepton[0].Pt(), EvWeight);
  645.     if(GoodLepton.size() > 1 && GoodLepton[1].Pt() > 0) hLeadLep2Pt4->Fill(GoodLepton[1].Pt(), EvWeight);
  646.     if(GoodLepton.size() > 2 && GoodLepton[2].Pt() > 0) hLeadLep3Pt4->Fill(GoodLepton[2].Pt(), EvWeight);
  647.     hMeff4->Fill(Meff, EvWeight);
  648.     if(CleanJets.size() > 0) hj1Eta4->Fill(CleanJets[0].Eta(), EvWeight);
  649.     if(CleanJets.size() > 1) hj2Eta4->Fill(CleanJets[1].Eta(), EvWeight);
  650.     if(CleanJets.size() > 2) hj3Eta4->Fill(CleanJets[2].Eta(), EvWeight);
  651.     if(CleanJets.size() > 3) hj4Eta4->Fill(CleanJets[3].Eta(), EvWeight);
  652.     if(GoodLepton.size() > 0) hLeadLepEta4->Fill(GoodLepton[0].Eta(), EvWeight);
  653.     if(GoodLepton.size() > 1) hLeadLep2Eta4->Fill(GoodLepton[1].Eta(), EvWeight);
  654.     if(GoodLepton.size() > 2) hLeadLep3Eta4->Fill(GoodLepton[2].Eta(), EvWeight);
  655.     if(GoodLepton.size() > 0) hLeadLepIso4->Fill(GoodLeptonIso[0], EvWeight);
  656.     if(GoodLepton.size() > 1) hLeadLep2Iso4->Fill(GoodLeptonIso[1], EvWeight);
  657.     if(GoodLepton.size() > 2) hLeadLep3Iso4->Fill(GoodLeptonIso[2], EvWeight);
  658.     hST4->Fill(ST, EvWeight);
  659.     //if(MT > 0) hMT4->Fill(MT, EvWeight);
  660.     }
  661.     /////////////////////////////
  662.     //LOOP ENDS HERE
  663.     /////////////////////////////
  664.  
  665.     ofstream tfile;
  666.     tfile.open("SummerStd_"+outname+".txt");
  667.     tfile << "########################################" << endl;
  668.     tfile << "Cut efficiency numbers:" << endl;
  669.     for(int ci = 0; ci < CutNumb; ci++)
  670.     {
  671.         tfile << "After cut " << CutList[ci] <<  "\t\t\t"
  672.           << CFCounter[ci]  << "\t events left\t"<< iCFCounter[ci] <<" cnt\t"<< endl;
  673.         //CF->SetBinContent(1+ci,CFCounter[ci]);
  674.     }
  675.     // Extreme care needed here
  676.     for(int cf = 0; cf < CutFlow; cf++)
  677.     {
  678.         CF->SetBinContent(1+cf, CFC[cf]);
  679.     }
  680.    
  681.     TFile * outf = new TFile("SummerStd_"+outname+"_his.root","RECREATE"); // The root file is generated
  682.    
  683.     hNmu->Write();
  684.     hNmu2->Write();
  685.     hNmu3->Write();
  686.     hNmu4->Write();
  687.     hNel->Write();
  688.     hNel2->Write();
  689.     hNel3->Write();
  690.     hNel4->Write();
  691.     hNlep->Write();
  692.     hNlep2->Write();
  693.     hNlep3->Write();
  694.     hNlep4->Write();
  695.     hHT->Write();
  696.     hHT2->Write();
  697.     hHT3->Write();
  698.     hHT4->Write();
  699.     hMET->Write();
  700.     hMET2->Write();
  701.     hMET3->Write();
  702.     hMET4->Write();
  703.     hNj->Write();
  704.     hNj2->Write();
  705.     hNj3->Write();
  706.     hNj4->Write();
  707.     hdPhi1->Write();
  708.     hdPhi12->Write();
  709.     hdPhi13->Write();
  710.     hdPhi14->Write();
  711.     hdPhi2->Write();
  712.     hdPhi22->Write();
  713.     hdPhi23->Write();
  714.     hdPhi24->Write();
  715.     hNbj->Write();
  716.     hNbj2->Write();
  717.     hNbj3->Write();
  718.     hNbj4->Write();
  719.     //hMETMeff->Write();
  720.     //hMETMeff2->Write();
  721.     //hMETMeff3->Write();
  722.     //hMETMeff4->Write();
  723.     hDiLMassOS->Write();
  724.     hDiLMassOS2->Write();
  725.     hDiLMassOS3->Write();
  726.     hDiLMassOS4->Write();
  727.     hDiLMassSS->Write();
  728.     hDiLMassSS2->Write();
  729.     hDiLMassSS3->Write();
  730.     hDiLMassSS4->Write();
  731.     hj1Pt->Write();
  732.     hj1Pt2->Write();
  733.     hj1Pt3->Write();
  734.     hj1Pt4->Write();
  735.     hj2Pt->Write();
  736.     hj2Pt2->Write();
  737.     hj2Pt3->Write();
  738.     hj2Pt4->Write();
  739.     hj3Pt->Write();
  740.     hj3Pt2->Write();
  741.     hj3Pt3->Write();
  742.     hj3Pt4->Write();
  743.     hj4Pt->Write();
  744.     hj4Pt2->Write();
  745.     hj4Pt3->Write();
  746.     hj4Pt4->Write();
  747.     hj1Eta->Write();
  748.     hj1Eta2->Write();
  749.     hj1Eta3->Write();
  750.     hj1Eta4->Write();
  751.     hj2Eta->Write();
  752.     hj2Eta2->Write();
  753.     hj2Eta3->Write();
  754.     hj2Eta4->Write();
  755.     hj3Eta->Write();
  756.     hj3Eta2->Write();
  757.     hj3Eta3->Write();
  758.     hj3Eta4->Write();
  759.     hj4Eta->Write();
  760.     hj4Eta2->Write();
  761.     hj4Eta3->Write();
  762.     hj4Eta4->Write();
  763.     hHasZ->Write();
  764.     hHasZ2->Write();
  765.     hHasZ3->Write();
  766.     hHasZ4->Write();
  767.     hMeff->Write();
  768.     hMeff2->Write();
  769.     hMeff3->Write();
  770.     hMeff4->Write();
  771.     hLeadLepPt->Write();
  772.     hLeadLepPt2->Write();
  773.     hLeadLepPt3->Write();
  774.     hLeadLepPt4->Write();
  775.     hLeadLep2Pt->Write();
  776.     hLeadLep2Pt2->Write();
  777.     hLeadLep2Pt3->Write();
  778.     hLeadLep2Pt4->Write();
  779.     hLeadLep3Pt->Write();
  780.     hLeadLep3Pt2->Write();
  781.     hLeadLep3Pt3->Write();
  782.     hLeadLep3Pt4->Write();
  783.     hLeadLepEta->Write();
  784.     hLeadLepEta2->Write();
  785.     hLeadLepEta3->Write();
  786.     hLeadLepEta4->Write();
  787.     hLeadLep2Eta->Write();
  788.     hLeadLep2Eta2->Write();
  789.     hLeadLep2Eta3->Write();
  790.     hLeadLep2Eta4->Write();
  791.     hLeadLep3Eta->Write();
  792.     hLeadLep3Eta2->Write();
  793.     hLeadLep3Eta3->Write();
  794.     hLeadLep3Eta4->Write();
  795.     hLeadLepIso->Write();
  796.     hLeadLepIso2->Write();
  797.     hLeadLepIso3->Write();
  798.     hLeadLepIso4->Write();
  799.     hLeadLep2Iso->Write();
  800.     hLeadLep2Iso2->Write();
  801.     hLeadLep2Iso3->Write();
  802.     hLeadLep2Iso4->Write();
  803.     hLeadLep3Iso->Write();
  804.     hLeadLep3Iso2->Write();
  805.     hLeadLep3Iso3->Write();
  806.     hLeadLep3Iso4->Write();
  807.     hST->Write();
  808.     hST2->Write();
  809.     hST3->Write();
  810.     hST4->Write();
  811.     //hMT->Write();
  812.     //hMT2->Write();
  813.     //hMT3->Write();
  814.     //hMT4->Write();
  815.     CF->Write(); // This allows us to visualize the cuts
  816.     outf->Close();
  817. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement