CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch2/src/HLTrigger/HLTanalyzers/src/HLTJets.cc

Go to the documentation of this file.
00001 #include <iostream>
00002 #include <sstream>
00003 #include <istream>
00004 #include <fstream>
00005 #include <iomanip>
00006 #include <string>
00007 #include <cmath>
00008 #include <functional>
00009 #include <stdlib.h>
00010 #include <string.h>
00011 #include <set>
00012 
00013 #include "HLTrigger/HLTanalyzers/interface/HLTJets.h"
00014 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00015 #include "DataFormats/TrackReco/interface/Track.h"
00016 
00017 HLTJets::HLTJets() {
00018     evtCounter=0;
00019     
00020     //set parameter defaults 
00021     _Monte=false;
00022     _Debug=false;
00023     _CalJetMin=0.;
00024     _GenJetMin=0.;
00025 }
00026 
00027 /*  Setup the analysis to put the branch-variables into the tree. */
00028 void HLTJets::setup(const edm::ParameterSet& pSet, TTree* HltTree) {
00029     
00030     edm::ParameterSet myJetParams = pSet.getParameter<edm::ParameterSet>("RunParameters") ;
00031     std::vector<std::string> parameterNames = myJetParams.getParameterNames() ;
00032     
00033     for ( std::vector<std::string>::iterator iParam = parameterNames.begin();
00034          iParam != parameterNames.end(); iParam++ ){
00035         if  ( (*iParam) == "Monte" ) _Monte =  myJetParams.getParameter<bool>( *iParam );
00036         else if ( (*iParam) == "Debug" ) _Debug =  myJetParams.getParameter<bool>( *iParam );
00037         else if ( (*iParam) == "CalJetMin" ) _CalJetMin =  myJetParams.getParameter<double>( *iParam );
00038         else if ( (*iParam) == "GenJetMin" ) _GenJetMin =  myJetParams.getParameter<double>( *iParam );
00039     }
00040 
00041     jetID = new reco::helper::JetIDHelper(pSet.getParameter<edm::ParameterSet>("JetIDParams"));
00042 
00043     const int kMaxRecoPFJet = 10000;
00044     jpfrecopt=new float[kMaxRecoPFJet];
00045     jpfrecoe=new float[kMaxRecoPFJet];
00046     jpfrecophi=new float[kMaxRecoPFJet];
00047     jpfrecoeta=new float[kMaxRecoPFJet];
00048     jpfreconeutralHadronFraction=new float[kMaxRecoPFJet];
00049     jpfreconeutralEMFraction=new float[kMaxRecoPFJet];
00050     jpfrecochargedHadronFraction=new float[kMaxRecoPFJet];
00051     jpfrecochargedEMFraction=new float[kMaxRecoPFJet];
00052     jpfreconeutralMultiplicity=new int[kMaxRecoPFJet];
00053     jpfrecochargedMultiplicity=new int[kMaxRecoPFJet];
00054     
00055     const int kMaxJetCal = 10000;
00056     jhcalpt = new float[kMaxJetCal];
00057     jhcalphi = new float[kMaxJetCal];
00058     jhcaleta = new float[kMaxJetCal];
00059     jhcale = new float[kMaxJetCal];
00060     jhcalemf = new float[kMaxJetCal]; 
00061     jhcaln90 = new float[kMaxJetCal]; 
00062     jhcaln90hits = new float[kMaxJetCal]; 
00063     
00064     jhcorcalpt = new float[kMaxJetCal]; 
00065     jhcorcalphi = new float[kMaxJetCal]; 
00066     jhcorcaleta = new float[kMaxJetCal]; 
00067     jhcorcale = new float[kMaxJetCal]; 
00068     jhcorcalemf = new float[kMaxJetCal]; 
00069     jhcorcaln90 = new float[kMaxJetCal]; 
00070     jhcorcaln90hits = new float[kMaxJetCal]; 
00071 
00072     jhcorL1L2L3calpt = new float[kMaxJetCal]; 
00073     jhcorL1L2L3calphi = new float[kMaxJetCal]; 
00074     jhcorL1L2L3caleta = new float[kMaxJetCal]; 
00075     jhcorL1L2L3cale = new float[kMaxJetCal]; 
00076     jhcorL1L2L3calemf = new float[kMaxJetCal]; 
00077     jhcorL1L2L3caln90 = new float[kMaxJetCal]; 
00078     jhcorL1L2L3caln90hits = new float[kMaxJetCal]; 
00079    
00080     jrcalpt = new float[kMaxJetCal];
00081     jrcalphi = new float[kMaxJetCal];
00082     jrcaleta = new float[kMaxJetCal];
00083     jrcale = new float[kMaxJetCal];
00084     jrcalemf = new float[kMaxJetCal]; 
00085     jrcaln90 = new float[kMaxJetCal]; 
00086     jrcaln90hits = new float[kMaxJetCal]; 
00087 
00088     jrcorcalpt = new float[kMaxJetCal];
00089     jrcorcalphi = new float[kMaxJetCal];
00090     jrcorcaleta = new float[kMaxJetCal];
00091     jrcorcale = new float[kMaxJetCal];
00092     jrcorcalemf = new float[kMaxJetCal]; 
00093     jrcorcaln90 = new float[kMaxJetCal]; 
00094     jrcorcaln90hits = new float[kMaxJetCal]; 
00095 
00096     const int kMaxJetgen = 10000;
00097     jgenpt = new float[kMaxJetgen];
00098     jgenphi = new float[kMaxJetgen];
00099     jgeneta = new float[kMaxJetgen];
00100     jgene = new float[kMaxJetgen];
00101     const int kMaxTower = 10000;
00102     towet = new float[kMaxTower];
00103     toweta = new float[kMaxTower];
00104     towphi = new float[kMaxTower];
00105     towen = new float[kMaxTower];
00106     towem = new float[kMaxTower];
00107     towhd = new float[kMaxTower];
00108     towoe = new float[kMaxTower];
00109     towR45upper = new int[kMaxTower];
00110     towR45lower = new int[kMaxTower];
00111     towR45none = new int[kMaxTower];
00112     const int kMaxTau = 500;
00113     l2tauPt    = new float[kMaxTau];
00114     l2tauEta   = new float[kMaxTau];
00115     l2tauPhi   = new float[kMaxTau];
00116     l2tauemiso = new float[kMaxTau];
00117     l25tauPt = new float[kMaxTau];
00118     l3tautckiso = new int[kMaxTau];
00119     tauEta = new float[kMaxTau];
00120     tauPt = new float[kMaxTau];
00121     tauPhi = new float[kMaxTau];
00122     
00123     const int kMaxPFTau = 500;
00124     ohpfTauEta         =  new float[kMaxPFTau];
00125     ohpfTauProngs      =  new int[kMaxPFTau];
00126     ohpfTauPhi         =  new float[kMaxPFTau];
00127     ohpfTauPt          =  new float[kMaxPFTau];
00128     ohpfTauJetPt       =  new float[kMaxPFTau];
00129     ohpfTauLeadTrackPt =  new float[kMaxPFTau];
00130     ohpfTauLeadTrackVtxZ = new float[kMaxPFTau];
00131     ohpfTauLeadPionPt  =  new float[kMaxPFTau];
00132     ohpfTauTrkIso      =  new float[kMaxPFTau];
00133     ohpfTauGammaIso    =  new float[kMaxPFTau];
00134 
00135     ohpfTauTightConeProngs      =  new int[kMaxPFTau];
00136     ohpfTauTightConeEta         =  new float[kMaxPFTau];
00137     ohpfTauTightConePhi         =  new float[kMaxPFTau];
00138     ohpfTauTightConePt          =  new float[kMaxPFTau];
00139     ohpfTauTightConeJetPt       =  new float[kMaxPFTau];
00140     ohpfTauTightConeLeadTrackPt =  new float[kMaxPFTau];
00141     ohpfTauTightConeLeadPionPt  =  new float[kMaxPFTau];
00142     ohpfTauTightConeTrkIso      =  new float[kMaxPFTau];
00143     ohpfTauTightConeGammaIso    =  new float[kMaxPFTau];
00144     
00145     recopfTauEta         =  new float[kMaxPFTau];
00146     recopfTauPhi         =  new float[kMaxPFTau];
00147     recopfTauPt          =  new float[kMaxPFTau];
00148     recopfTauJetPt       =  new float[kMaxPFTau];
00149     recopfTauLeadTrackPt =  new float[kMaxPFTau];
00150     recopfTauLeadPionPt  =  new float[kMaxPFTau];
00151     recopfTauTrkIso      =  new int[kMaxPFTau];
00152     recopfTauGammaIso    =  new int[kMaxPFTau];
00153     recopfTauDiscrByTancOnePercent     =  new float[kMaxPFTau];
00154     recopfTauDiscrByTancHalfPercent    =  new float[kMaxPFTau];
00155     recopfTauDiscrByTancQuarterPercent =  new float[kMaxPFTau]; 
00156     recopfTauDiscrByTancTenthPercent   =  new float[kMaxPFTau];
00157     recopfTauDiscrByIso        =  new float[kMaxPFTau]; 
00158     recopfTauDiscrAgainstMuon  =  new float[kMaxPFTau];
00159     recopfTauDiscrAgainstElec  =  new float[kMaxPFTau];
00160     
00161     pfHT    = -100.;
00162     pfMHT   = -100.;    
00163     const int kMaxPFJet = 500;
00164     pfJetEta         = new float[kMaxPFJet];
00165     pfJetPhi         = new float[kMaxPFJet];
00166     pfJetPt         = new float[kMaxPFJet];
00167     pfJetE         = new float[kMaxPFJet];
00168     pfJetneutralHadronEnergyFraction = new float[kMaxPFJet];
00169     pfJetchargedHadronFraction = new float[kMaxPFJet];
00170     pfJetneutralMultiplicity = new float[kMaxPFJet];
00171     pfJetchargedMultiplicity = new float[kMaxPFJet];
00172     pfJetneutralEMFraction = new float[kMaxPFJet];
00173     pfJetchargedEMFraction = new float[kMaxPFJet];
00174 
00175 
00176     const int kMaxTauIso = 5000;
00177 
00178     // for offlineHPStau 
00179     signalTrToPFTauMatch = new int[kMaxTauIso];// index of reconstructed tau in tau collection
00180     recoPFTauSignalTrDz = new float[kMaxTauIso];
00181     recoPFTauSignalTrPt = new float[kMaxTauIso];
00182 
00183     isoTrToPFTauMatch = new int[kMaxTauIso]; // index of reconstructed tau in tau collection
00184     recoPFTauIsoTrDz = new float[kMaxTauIso];
00185     recoPFTauIsoTrPt = new float[kMaxTauIso];
00186 
00187     // HLT pf taus
00188     hltpftauSignalTrToPFTauMatch = new int[kMaxTauIso]; // index of HLTPF tau in tau collection
00189     HLTPFTauSignalTrDz = new float[kMaxTauIso];
00190     HLTPFTauSignalTrPt = new float[kMaxTauIso];
00191 
00192     hltpftauIsoTrToPFTauMatch = new int[kMaxTauIso]; // index of HLTPF tau in tau collection
00193     HLTPFTauIsoTrDz = new float[kMaxTauIso];
00194     HLTPFTauIsoTrPt = new float[kMaxTauIso];
00195 
00196     // offline pftau isolation and signal cands
00197     HltTree->Branch("NoRecoPFTausSignal",&noRecoPFTausSignal,"NoRecoPFTausSignal/I");
00198     HltTree->Branch("signalTrToPFTauMatch", signalTrToPFTauMatch,"signalTrToPFTauMatch[NoRecoPFTausSignal]/I");
00199     HltTree->Branch("recoPFTauSignalTrDz", recoPFTauSignalTrDz,"recoPFTauSignalTrDz[NoRecoPFTausSignal]/F");
00200     HltTree->Branch("recoPFTauSignalTrPt", recoPFTauSignalTrPt,"recoPFTauSignalTrPt[NoRecoPFTausSignal]/F");
00201 
00202     HltTree->Branch("NoRecoPFTausIso",&noRecoPFTausIso,"NoRecoPFTausIso/I");
00203     HltTree->Branch("isoTrToPFTauMatch", isoTrToPFTauMatch,"isoTrToPFTauMatch[NoRecoPFTausIso]/I");
00204     HltTree->Branch("recoPFTauIsoTrDz", recoPFTauIsoTrDz,"recoPFTauIsoTrDz[NoRecoPFTausIso]/F");
00205     HltTree->Branch("recoPFTauIsoTrPt", recoPFTauIsoTrPt,"recoPFTauIsoTrPt[NoRecoPFTausIso]/F");
00206 
00207     // HLT pftau isolation and signal cands
00208     HltTree->Branch("NoHLTPFTausSignal",&noHLTPFTausSignal,"NoHLTPFTausSignal/I");
00209     HltTree->Branch("hltpftauSignalTrToPFTauMatch",
00210                     hltpftauSignalTrToPFTauMatch,"hltpftauSignalTrToPFTauMatch[NoHLTPFTausSignal]/I");
00211     HltTree->Branch("HLTPFTauSignalTrDz", HLTPFTauSignalTrDz,"HLTPFTauSignalTrDz[NoHLTPFTausSignal]/F");
00212     HltTree->Branch("HLTPFTauSignalTrPt", HLTPFTauSignalTrPt,"HLTPFTauSignalTrPt[NoHLTPFTausSignal]/F");
00213 
00214     HltTree->Branch("NoHLTPFTausIso",&noHLTPFTausIso,"NoHLTPFTausIso/I");
00215     HltTree->Branch("hltpftauIsoTrToPFTauMatch",
00216                     hltpftauIsoTrToPFTauMatch,"hltpftauIsoTrToPFTauMatch[NoHLTPFTausIso]/I");
00217     HltTree->Branch("HLTPFTauIsoTrDz", HLTPFTauIsoTrDz,"HLTPFTauIsoTrDz[NoHLTPFTausIso]/F");
00218     HltTree->Branch("HLTPFTauIsoTrPt", HLTPFTauIsoTrPt,"HLTPFTauIsoTrPt[NoHLTPFTausIso]/F");
00219     
00220     // Jet- MEt-specific branches of the tree 
00221 
00222     HltTree->Branch("NrecoJetGen",&njetgen,"NrecoJetGen/I");
00223     HltTree->Branch("NrecoTowCal",&ntowcal,"NrecoTowCal/I");
00224 
00225     //ccla RECO JETs
00226     HltTree->Branch("NrecoJetCal",&nrjetcal,"NrecoJetCal/I");
00227     HltTree->Branch("recoJetCalPt",jrcalpt,"recoJetCalPt[NrecoJetCal]/F");
00228     HltTree->Branch("recoJetCalPhi",jrcalphi,"recoJetCalPhi[NrecoJetCal]/F");
00229     HltTree->Branch("recoJetCalEta",jrcaleta,"recoJetCalEta[NrecoJetCal]/F");
00230     HltTree->Branch("recoJetCalE",jrcale,"recoJetCalE[NrecoJetCal]/F");
00231     HltTree->Branch("recoJetCalEMF",jrcalemf,"recoJetCalEMF[NrecoJetCal]/F");
00232     HltTree->Branch("recoJetCalN90",jrcaln90,"recoJetCalN90[NrecoJetCal]/F");
00233     HltTree->Branch("recoJetCalN90hits",jrcaln90hits,"recoJetCalN90hits[NrecoJetCal]/F");
00234 
00235     HltTree->Branch("NrecoJetCorCal",&nrcorjetcal,"NrecoJetCorCal/I"); 
00236     HltTree->Branch("recoJetCorCalPt",jrcorcalpt,"recoJetCorCalPt[NrecoJetCorCal]/F"); 
00237     HltTree->Branch("recoJetCorCalPhi",jrcorcalphi,"recoJetCorCalPhi[NrecoJetCorCal]/F"); 
00238     HltTree->Branch("recoJetCorCalEta",jrcorcaleta,"recoJetCorCalEta[NrecoJetCorCal]/F"); 
00239     HltTree->Branch("recoJetCorCalE",jrcorcale,"recoJetCorCalE[NrecoJetCorCal]/F"); 
00240     HltTree->Branch("recoJetCorCalEMF",jrcorcalemf,"recoJetCorCalEMF[NrecoJetCorCal]/F");
00241     HltTree->Branch("recoJetCorCalN90",jrcorcaln90,"recoJetCorCalN90[NrecoJetCorCal]/F");
00242     HltTree->Branch("recoJetCorCalN90hits",jrcorcaln90hits,"recoJetCorCalN90hits[NrecoJetCorCal]/F");
00243  
00244     //ccla HLTJETS
00245     HltTree->Branch("NohJetCal",&nhjetcal,"NohJetCal/I");
00246     HltTree->Branch("ohJetCalPt",jhcalpt,"ohJetCalPt[NohJetCal]/F");
00247     HltTree->Branch("ohJetCalPhi",jhcalphi,"ohJetCalPhi[NohJetCal]/F");
00248     HltTree->Branch("ohJetCalEta",jhcaleta,"ohJetCalEta[NohJetCal]/F");
00249     HltTree->Branch("ohJetCalE",jhcale,"ohJetCalE[NohJetCal]/F");
00250     HltTree->Branch("ohJetCalEMF",jhcalemf,"ohJetCalEMF[NohJetCal]/F");
00251     HltTree->Branch("ohJetCalN90",jhcaln90,"ohJetCalN90[NohJetCal]/F");
00252     HltTree->Branch("ohJetCalN90hits",jhcaln90hits,"ohJetCalN90hits[NohJetCal]/F");
00253 
00254     HltTree->Branch("NohJetCorCal",&nhcorjetcal,"NohJetCorCal/I");
00255     HltTree->Branch("ohJetCorCalPt",jhcorcalpt,"ohJetCorCalPt[NohJetCorCal]/F");
00256     HltTree->Branch("ohJetCorCalPhi",jhcorcalphi,"ohJetCorCalPhi[NohJetCorCal]/F");
00257     HltTree->Branch("ohJetCorCalEta",jhcorcaleta,"ohJetCorCalEta[NohJetCorCal]/F");
00258     HltTree->Branch("ohJetCorCalE",jhcorcale,"ohJetCorCalE[NohJetCorCal]/F");
00259     HltTree->Branch("ohJetCorCalEMF",jhcorcalemf,"ohJetCorCalEMF[NohJetCorCal]/F");
00260     HltTree->Branch("ohJetCorCalN90",jhcorcaln90,"ohJetCorCalN90[NohJetCorCal]/F");
00261     HltTree->Branch("ohJetCorCalN90hits",jhcorcaln90hits,"ohJetCorCalN90hits[NohJetCorCal]/F");
00262 
00263     HltTree->Branch("NohJetCorL1L2L3Cal",&nhcorL1L2L3jetcal,"NohJetCorL1L2L3Cal/I");
00264     HltTree->Branch("ohJetCorL1L2L3CalPt",jhcorL1L2L3calpt,"ohJetCorL1L2L3CalPt[NohJetCorL1L2L3Cal]/F");
00265     HltTree->Branch("ohJetCorL1L2L3CalPhi",jhcorL1L2L3calphi,"ohJetCorL1L2L3CalPhi[NohJetCorL1L2L3Cal]/F");
00266     HltTree->Branch("ohJetCorL1L2L3CalEta",jhcorL1L2L3caleta,"ohJetCorL1L2L3CalEta[NohJetCorL1L2L3Cal]/F");
00267     HltTree->Branch("ohJetCorL1L2L3CalE",jhcorL1L2L3cale,"ohJetCorL1L2L3CalE[NohJetCorL1L2L3Cal]/F");
00268     HltTree->Branch("ohJetCorL1L2L3CalEMF",jhcorL1L2L3calemf,"ohJetCorL1L2L3CalEMF[NohJetCorL1L2L3Cal]/F");
00269     HltTree->Branch("ohJetCorL1L2L3CalN90",jhcorL1L2L3caln90,"ohJetCorL1L2L3CalN90[NohJetCorL1L2L3Cal]/F");
00270     HltTree->Branch("ohJetCorL1L2L3CalN90hits",jhcorL1L2L3caln90hits,"ohJetCorL1L2L3CalN90hits[NohJetCorL1L2L3Cal]/F");
00271     HltTree->Branch("rho",&jrho,"rho/D");
00272    
00273     //ccla GenJets
00274     HltTree->Branch("recoJetGenPt",jgenpt,"recoJetGenPt[NrecoJetGen]/F");
00275     HltTree->Branch("recoJetGenPhi",jgenphi,"recoJetGenPhi[NrecoJetGen]/F");
00276     HltTree->Branch("recoJetGenEta",jgeneta,"recoJetGenEta[NrecoJetGen]/F");
00277     HltTree->Branch("recoJetGenE",jgene,"recoJetGenE[NrecoJetGen]/F");
00278 
00279     HltTree->Branch("recoTowEt",towet,"recoTowEt[NrecoTowCal]/F");
00280     HltTree->Branch("recoTowEta",toweta,"recoTowEta[NrecoTowCal]/F");
00281     HltTree->Branch("recoTowPhi",towphi,"recoTowPhi[NrecoTowCal]/F");
00282     HltTree->Branch("recoTowE",towen,"recoTowE[NrecoTowCal]/F");
00283     HltTree->Branch("recoTowEm",towem,"recoTowEm[NrecoTowCal]/F");
00284     HltTree->Branch("recoTowHad",towhd,"recoTowHad[NrecoTowCal]/F");
00285     HltTree->Branch("recoTowOE",towoe,"recoTowOE[NrecoTowCal]/F");
00286     HltTree->Branch("recoTowHCalNoiseR45Upper",towR45upper,"recoTowHCalNoiseR45Upper[NrecoTowCal]/I");
00287     HltTree->Branch("recoTowHCalNoiseR45Lower",towR45lower,"recoTowHCalNoiseR45Lower[NrecoTowCal]/I");
00288     HltTree->Branch("recoTowHCalNoiseR45None",towR45none,"recoTowHCalNoiseR45None[NrecoTowCal]/I");
00289 
00290     HltTree->Branch("recoMetCal",&mcalmet,"recoMetCal/F");
00291     HltTree->Branch("recoMetCalPhi",&mcalphi,"recoMetCalPhi/F");
00292     HltTree->Branch("recoMetCalSum",&mcalsum,"recoMetCalSum/F");
00293     HltTree->Branch("recoMetGen",&mgenmet,"recoMetGen/F");
00294     HltTree->Branch("recoMetGenPhi",&mgenphi,"recoMetGenPhi/F");
00295     HltTree->Branch("recoMetGenSum",&mgensum,"recoMetGenSum/F");
00296     HltTree->Branch("recoHTCal",&htcalet,"recoHTCal/F");
00297     HltTree->Branch("recoHTCalPhi",&htcalphi,"recoHTCalPhi/F");
00298     HltTree->Branch("recoHTCalSum",&htcalsum,"recoHTCalSum/F");
00299     HltTree->Branch("recoMetPF", &pfmet, "recoMetPF/F");
00300     HltTree->Branch("recoMetPFSum", &pfsumet, "recoMetPFSum/F");
00301     HltTree->Branch("recoMetPFPhi", &pfmetphi, "recoMetPFPhi/F");
00302 
00303     //for(int ieta=0;ieta<NETA;ieta++){std::cout << " ieta " << ieta << " eta min " << CaloTowerEtaBoundries[ieta] <<std::endl;}
00304     
00305     
00306     // Taus
00307     nohl2tau = 0;
00308     HltTree->Branch("NohTauL2",&nohl2tau,"NohTauL2/I");
00309     HltTree->Branch("ohTauL2Pt",l2tauPt,"ohTauL2Pt[NohTauL2]/F");
00310     HltTree->Branch("ohTauL2Eta",l2tauEta,"ohTauL2Eta[NohTauL2]/F");
00311     HltTree->Branch("ohTauL2Phi",l2tauPhi,"ohTauL2Phi[NohTauL2]/F");
00312 
00313     nohtau = 0;
00314     HltTree->Branch("NohTau",&nohtau,"NohTau/I");
00315     HltTree->Branch("ohTauEta",tauEta,"ohTauEta[NohTau]/F");
00316     HltTree->Branch("ohTauPhi",tauPhi,"ohTauPhi[NohTau]/F");
00317     HltTree->Branch("ohTauPt",tauPt,"ohTauPt[NohTau]/F");
00318     HltTree->Branch("ohTauEiso",l2tauemiso,"ohTauEiso[NohTau]/F");
00319     HltTree->Branch("ohTauL25Tpt",l25tauPt,"ohTauL25Tpt[NohTau]/F");
00320     HltTree->Branch("ohTauL3Tiso",l3tautckiso,"ohTauL3Tiso[NohTau]/I");
00321 
00322     //ohpfTaus
00323     nohPFTau = 0;
00324     HltTree->Branch("NohpfTau",&nohPFTau,"NohpfTau/I");
00325     HltTree->Branch("ohpfTauPt",ohpfTauPt,"ohpfTauPt[NohpfTau]/F");
00326     HltTree->Branch("ohpfTauProngs",ohpfTauProngs,"ohpfTauProngs[NohpfTau]/I");
00327     HltTree->Branch("ohpfTauEta",ohpfTauEta,"ohpfTauEta[NohpfTau]/F");
00328     HltTree->Branch("ohpfTauPhi",ohpfTauPhi,"ohpfTauPhi[NohpfTau]/F");
00329     HltTree->Branch("ohpfTauLeadTrackPt",ohpfTauLeadTrackPt,"ohpfTauLeadTrackPt[NohpfTau]/F");
00330     HltTree->Branch("ohpfTauLeadTrackVtxZ",ohpfTauLeadTrackVtxZ,"ohpfTauLeadTrackVtxZ[NohpfTau]/F");
00331     HltTree->Branch("ohpfTauLeadPionPt",ohpfTauLeadPionPt,"ohpfTauLeadPionPt[NohpfTau]/F");
00332     HltTree->Branch("ohpfTauTrkIso",ohpfTauTrkIso,"ohpfTauTrkIso[NohpfTau]/F");
00333     HltTree->Branch("ohpfTauGammaIso",ohpfTauGammaIso,"ohpfTauGammaIso[NohpfTau]/F");
00334     HltTree->Branch("ohpfTauJetPt",ohpfTauJetPt,"ohpfTauJetPt[NohpfTau]/F");    
00335 
00336     //ohpfTaus tight cone
00337     nohPFTauTightCone = 0;
00338     HltTree->Branch("NohpfTauTightCone",&nohPFTauTightCone,"NohpfTauTightCone/I");
00339     HltTree->Branch("ohpfTauTightConePt",ohpfTauTightConePt,"ohpfTauTightConePt[NohpfTauTightCone]/F");
00340     HltTree->Branch("ohpfTauTightConeProngs",ohpfTauTightConeProngs,"ohpfTauProngs[NohpfTauTightCone]/I");
00341     HltTree->Branch("ohpfTauTightConeEta",ohpfTauTightConeEta,"ohpfTauEta[NohpfTauTightCone]/F");
00342     HltTree->Branch("ohpfTauTightConePhi",ohpfTauTightConePhi,"ohpfTauPhi[NohpfTauTightCone]/F");
00343     HltTree->Branch("ohpfTauTightConeLeadTrackPt",ohpfTauTightConeLeadTrackPt,"ohpfTauTightConeLeadTrackPt[NohpfTauTightCone]/F");
00344     HltTree->Branch("ohpfTauTightConeLeadPionPt",ohpfTauTightConeLeadPionPt,"ohpfTauTightConeLeadPionPt[NohpfTauTightCone]/F");
00345     HltTree->Branch("ohpfTauTightConeTrkIso",ohpfTauTightConeTrkIso,"ohpfTauTightConeTrkIso[NohpfTauTightCone]/F");
00346     HltTree->Branch("ohpfTauTightConeGammaIso",ohpfTauTightConeGammaIso,"ohpfTauTightConeGammaIso[NohpfTauTightCone]/F");
00347     HltTree->Branch("ohpfTauTightConeJetPt",ohpfTauTightConeJetPt,"ohpfTauTightConeJetPt[NohpfTauTightCone]/F");
00348    
00349    //Reco PFTaus
00350     nRecoPFTau = 0;
00351     HltTree->Branch("NRecoPFTau",&nRecoPFTau,"NRecoPFTau/I");
00352     HltTree->Branch("recopfTauPt",recopfTauPt,"recopfTauPt[NRecoPFTau]/F");
00353     HltTree->Branch("recopfTauEta",recopfTauEta,"recopfTauEta[NRecoPFTau]/F");
00354     HltTree->Branch("recopfTauPhi",recopfTauPhi,"recopfTauPhi[NRecoPFTau]/F");
00355     HltTree->Branch("recopfTauLeadTrackPt",recopfTauLeadTrackPt,"recopfTauLeadTrackPt[NRecoPFTau]/F");
00356     HltTree->Branch("recopfTauLeadPionPt",recopfTauLeadPionPt,"recopfTauLeadPionPt[NRecoPFTau]/F");
00357     HltTree->Branch("recopfTauTrkIso",recopfTauTrkIso,"recopfTauTrkIso[NRecoPFTau]/I");
00358     HltTree->Branch("recopfTauGammaIso",recopfTauGammaIso,"recopfTauGammaIso[NRecoPFTau]/I");
00359     HltTree->Branch("recopfTauJetPt",recopfTauJetPt,"recopfTauJetPt[NRecoPFTau]/F");   
00360     HltTree->Branch("recopfTauDiscrByTancOnePercent",recopfTauDiscrByTancOnePercent,"recopfTauDiscrByTancOnePercent[NRecoPFTau]/F");   
00361     HltTree->Branch("recopfTauDiscrByTancHalfPercent",recopfTauDiscrByTancHalfPercent,"recopfTauDiscrByTancHalfPercent[NRecoPFTau]/F");   
00362     HltTree->Branch("recopfTauDiscrByTancQuarterPercent",recopfTauDiscrByTancQuarterPercent,"recopfTauDiscrByTancQuarterPercent[NRecoPFTau]/F");   
00363     HltTree->Branch("recopfTauDiscrByTancTenthPercent",recopfTauDiscrByTancTenthPercent,"recopfTauDiscrByTancTenthPercent[NRecoPFTau]/F");       
00364     HltTree->Branch("recopfTauDiscrByIso",recopfTauDiscrByIso,"recopfTauDiscrByIso[NRecoPFTau]/F");   
00365     HltTree->Branch("recopfTauDiscrAgainstMuon",recopfTauDiscrAgainstMuon,"recopfTauDiscrAgainstMuon[NRecoPFTau]/F");
00366     HltTree->Branch("recopfTauDiscrAgainstElec",recopfTauDiscrAgainstElec,"recopfTauDiscrAgainstElec[NRecoPFTau]/F");
00367   
00368     //PFJets
00369     nohPFJet = 0;
00370     HltTree->Branch("pfHT",&pfHT,"pfHT/F");
00371     HltTree->Branch("pfMHT",&pfMHT,"pfMHT/F");
00372     HltTree->Branch("NohPFJet",&nohPFJet,"NohPFJet/I");
00373     HltTree->Branch("pfJetPt",pfJetPt,"pfJetPt[NohPFJet]/F");
00374     HltTree->Branch("pfJetE",pfJetE,"pfJetE[NohPFJet]/F");
00375     HltTree->Branch("pfJetEta",pfJetEta,"pfJetEta[NohPFJet]/F");
00376     HltTree->Branch("pfJetPhi",pfJetPhi,"pfJetPhi[NohPFJet]/F");
00377     HltTree->Branch("pfJetneutralHadronEnergyFraction",pfJetneutralHadronEnergyFraction,"pfJetneutralHadronEnergyFraction[NohPFJet]/F");
00378     HltTree->Branch("pfJetchargedHadronFraction",pfJetchargedHadronFraction,"pfJetchargedHadronFraction[NohPFJet]/F");
00379     HltTree->Branch("pfJetneutralMultiplicity",pfJetneutralMultiplicity,"pfJetneutralMultiplicity[NohPFJet]/F");
00380     HltTree->Branch("pfJetchargedMultiplicity",pfJetchargedMultiplicity,"pfJetchargedMultiplicity[NohPFJet]/F");
00381     HltTree->Branch("pfJetneutralEMFraction",pfJetneutralEMFraction,"pfJetneutralEMFraction[NohPFJet]/F");
00382     HltTree->Branch("pfJetchargedEMFraction",pfJetchargedEMFraction,"pfJetchargedEMFraction[NohPFJet]/F");
00383 
00384     //RECO PFJets
00385     HltTree->Branch("nrpj",&nrpj,"nrpj/I");
00386     HltTree->Branch("recopfJetpt",                    jpfrecopt,                     "recopfJetpt[nrpj]/F");
00387     HltTree->Branch("recopfJete",                     jpfrecoe,                      "recopfJete[nrpj]/F");
00388     HltTree->Branch("recopfJetphi",                   jpfrecophi,                    "recopfJetphi[nrpj]/F");
00389     HltTree->Branch("recopfJeteta",                   jpfrecoeta,                    "recopfJeteta[nrpj]/F");
00390     HltTree->Branch("recopfJetneutralHadronFraction", jpfreconeutralHadronFraction,  "recopfJetneutralHadronFraction[nrpj]/F");
00391     HltTree->Branch("recopfJetneutralEMFraction",     jpfreconeutralEMFraction,      "recopfJetneutralEMFraction[nrpj]/F");
00392     HltTree->Branch("recopfJetchargedHadronFraction", jpfrecochargedHadronFraction,  "recopfJetchargedHadronFraction[nrpj]/F");
00393     HltTree->Branch("recopfJetchargedEMFraction",     jpfrecochargedEMFraction,      "recopfJetchargedEMFraction[nrpj]/F");
00394     HltTree->Branch("recopfJetneutralMultiplicity",   jpfreconeutralMultiplicity,    "recopfJetneutralMultiplicity[nrpj]/I");
00395     HltTree->Branch("recopfJetchargedMultiplicity",   jpfrecochargedMultiplicity,    "recopfJetchargedMultiplicity[nrpj]/I"); 
00396     
00397 }
00398 
00399 /* **Analyze the event** */
00400 void HLTJets::analyze(edm::Event const& iEvent,
00401                       const edm::Handle<reco::CaloJetCollection>      & ohcalojets,
00402                       const edm::Handle<reco::CaloJetCollection>      & ohcalocorjets,
00403                       const edm::Handle<reco::CaloJetCollection>      & ohcalocorL1L2L3jets,
00404                       const edm::Handle< double >                     & rho,
00405                       const edm::Handle<reco::CaloJetCollection>      & rcalojets,
00406                       const edm::Handle<reco::CaloJetCollection>      & rcalocorjets,
00407                       const edm::Handle<reco::GenJetCollection>       & genjets,
00408                       const edm::Handle<reco::CaloMETCollection>      & recmets,
00409                       const edm::Handle<reco::GenMETCollection>       & genmets,
00410                       const edm::Handle<reco::METCollection>          & ht,
00411                       const edm::Handle<reco::CaloJetCollection>      & l2taujets,
00412                       const edm::Handle<reco::HLTTauCollection>       & taujets,
00413                       const edm::Handle<reco::PFTauCollection>        & pfTaus,
00414                       const edm::Handle<reco::PFTauCollection>        & pfTausTightCone,
00415                       const edm::Handle<reco::PFJetCollection>        & pfJets,
00416                       const edm::Handle<reco::PFTauCollection>        & recoPfTaus,  
00417                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrByTanCOnePercent,
00418                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrByTanCHalfPercent,
00419                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrByTanCQuarterPercent,
00420                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrByTanCTenthPercent,                      
00421                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrByIsolation,
00422                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrAgainstElec,
00423                       const edm::Handle<reco::PFTauDiscriminator>             & theRecoPFTauDiscrAgainstMuon,
00424                       const edm::Handle<reco::PFJetCollection>        & recoPFJets,
00425                       const edm::Handle<CaloTowerCollection>          & caloTowers,
00426                       const edm::Handle<CaloTowerCollection>          & caloTowersCleanerUpperR45,
00427                       const edm::Handle<CaloTowerCollection>          & caloTowersCleanerLowerR45,
00428                       const edm::Handle<CaloTowerCollection>          & caloTowersCleanerNoR45,
00429                       const edm::Handle<reco::PFMETCollection>        & pfmets, 
00430                       double thresholdForSavingTowers, 
00431                       double                minPtCH,
00432                       double               minPtGamma,
00433                       TTree * HltTree) {
00434     
00435     if (_Debug) std::cout << " Beginning HLTJets " << std::endl;
00436     
00437     //initialize branch variables
00438     nhjetcal=0; nhcorjetcal=0;nhcorL1L2L3jetcal=0; njetgen=0;ntowcal=0;
00439     jrho = 0;
00440     mcalmet=0.; mcalphi=0.;
00441     mgenmet=0.; mgenphi=0.;
00442     htcalet=0.,htcalphi=0.,htcalsum=0.;
00443 
00444     noRecoPFTausSignal = 0; noRecoPFTausIso =0;
00445     noHLTPFTausSignal = 0; noHLTPFTausIso = 0;
00446 
00447 
00448 
00449     if (rcalojets.isValid()) {
00450       reco::CaloJetCollection mycalojets;
00451       mycalojets=*rcalojets;
00452       std::sort(mycalojets.begin(),mycalojets.end(),PtGreater());
00453       typedef reco::CaloJetCollection::const_iterator cjiter;
00454       int jrcal=0;
00455       for ( cjiter i=mycalojets.begin(); i!=mycalojets.end(); i++) {
00456     
00457         if (i->pt()>_CalJetMin && i->energy()>0.){
00458           jrcalpt[jrcal] = i->pt();
00459           jrcalphi[jrcal] = i->phi();
00460           jrcaleta[jrcal] = i->eta();
00461           jrcale[jrcal] = i->energy();
00462           jrcalemf[jrcal] = i->emEnergyFraction();
00463           jrcaln90[jrcal] = i->n90();
00464           jetID->calculate( iEvent, *i );
00465           jrcaln90hits[jrcal] = jetID->n90Hits();
00466           jrcal++;
00467         }
00468       }
00469       nrjetcal = jrcal;
00470     }
00471     else {nrjetcal = 0;}
00472     
00473     if (rcalocorjets.isValid()) {
00474       reco::CaloJetCollection mycalojets;
00475       mycalojets=*rcalocorjets;
00476       std::sort(mycalojets.begin(),mycalojets.end(),PtGreater());
00477       typedef reco::CaloJetCollection::const_iterator cjiter;
00478       int jrcal=0;
00479       for ( cjiter i=mycalojets.begin(); i!=mycalojets.end(); i++) {
00480     
00481         if (i->pt()>_CalJetMin && i->energy()>0.){
00482           jrcorcalpt[jrcal] = i->pt();
00483           jrcorcalphi[jrcal] = i->phi();
00484           jrcorcaleta[jrcal] = i->eta();
00485           jrcorcale[jrcal] = i->energy();
00486           jrcorcalemf[jrcal] = i->emEnergyFraction();
00487           jrcorcaln90[jrcal] = i->n90();
00488           jetID->calculate( iEvent, *i );
00489           jrcorcaln90hits[jrcal] = jetID->n90Hits();
00490           jrcal++;
00491         }
00492       }
00493       nrcorjetcal = jrcal;
00494     }
00495     else {nrcorjetcal = 0;}
00496     
00497     if (ohcalojets.isValid()) {
00498         reco::CaloJetCollection mycalojets;
00499         mycalojets=*ohcalojets;
00500         std::sort(mycalojets.begin(),mycalojets.end(),PtGreater());
00501         typedef reco::CaloJetCollection::const_iterator cjiter;
00502         int jhcal=0;
00503         for ( cjiter i=mycalojets.begin(); i!=mycalojets.end(); i++) {
00504             
00505             if (i->pt()>_CalJetMin && i->energy()>0.){
00506                 jhcalpt[jhcal] = i->pt();
00507                 jhcalphi[jhcal] = i->phi();
00508                 jhcaleta[jhcal] = i->eta();
00509                 jhcale[jhcal] = i->energy();
00510                 jhcalemf[jhcal] = i->emEnergyFraction();
00511                 jhcaln90[jhcal] = i->n90();
00512                 jetID->calculate( iEvent, *i );
00513                 jhcaln90hits[jhcal] = jetID->n90Hits();
00514                 jhcal++;
00515             }
00516             
00517         }
00518         nhjetcal = jhcal;
00519     }
00520     else {nhjetcal = 0;}
00521     
00522     if (ohcalocorjets.isValid()) {
00523         reco::CaloJetCollection mycalocorjets;
00524         mycalocorjets=*ohcalocorjets;
00525         std::sort(mycalocorjets.begin(),mycalocorjets.end(),PtGreater());
00526         typedef reco::CaloJetCollection::const_iterator ccorjiter;
00527         int jhcorcal=0;
00528         for ( ccorjiter i=mycalocorjets.begin(); i!=mycalocorjets.end(); i++) {
00529             
00530             if (i->pt()>_CalJetMin && i->energy()>0.){
00531                 jhcorcalpt[jhcorcal] = i->pt();
00532                 jhcorcalphi[jhcorcal] = i->phi();
00533                 jhcorcaleta[jhcorcal] = i->eta();
00534                 jhcorcale[jhcorcal] = i->energy();
00535                 jhcorcalemf[jhcorcal] = i->emEnergyFraction();
00536                 jhcorcaln90[jhcorcal] = i->n90();
00537                 jetID->calculate( iEvent, *i );
00538                 jhcorcaln90hits[jhcorcal] = jetID->n90Hits();
00539                 jhcorcal++;
00540             }
00541             
00542         }
00543         nhcorjetcal = jhcorcal;
00544     }
00545     else {nhcorjetcal = 0;}
00546     
00547  if (ohcalocorL1L2L3jets.isValid()) {
00548         reco::CaloJetCollection mycalocorL1L2L3jets;
00549         mycalocorL1L2L3jets=*ohcalocorL1L2L3jets;
00550         std::sort(mycalocorL1L2L3jets.begin(),mycalocorL1L2L3jets.end(),PtGreater());
00551         typedef reco::CaloJetCollection::const_iterator ccorL1L2L3jiter;
00552         int jhcorL1L2L3cal=0;
00553         for ( ccorL1L2L3jiter i=mycalocorL1L2L3jets.begin(); i!=mycalocorL1L2L3jets.end(); i++) {
00554             
00555             if (i->pt()>_CalJetMin && i->energy()>0.){
00556                 jhcorL1L2L3calpt[jhcorL1L2L3cal] = i->pt();
00557                 jhcorL1L2L3calphi[jhcorL1L2L3cal] = i->phi();
00558                 jhcorL1L2L3caleta[jhcorL1L2L3cal] = i->eta();
00559                 jhcorL1L2L3cale[jhcorL1L2L3cal] = i->energy();
00560                 jhcorL1L2L3calemf[jhcorL1L2L3cal] = i->emEnergyFraction();
00561                 jhcorL1L2L3caln90[jhcorL1L2L3cal] = i->n90();
00562                 jetID->calculate( iEvent, *i );
00563                 jhcorL1L2L3caln90hits[jhcorL1L2L3cal] = jetID->n90Hits();
00564                 jhcorL1L2L3cal++;
00565             }
00566             
00567         }
00568         nhcorL1L2L3jetcal = jhcorL1L2L3cal;
00569    }
00570  else {nhcorL1L2L3jetcal = 0;}
00571 
00572  if (rho.isValid()){
00573    jrho = *rho;
00574  }
00575  else {
00576 
00577    if (_Debug) std::cout << "rho not found" << std::endl;
00578  }
00579 
00580  std::set<unsigned int> towersUpper;
00581  std::set<unsigned int> towersLower;
00582  std::set<unsigned int> towersNone;
00583 
00584  bool towersUpperValid=false;
00585  bool towersLowerValid=false;
00586  bool towersNoneValid=false;
00587  if( caloTowersCleanerUpperR45.isValid() ){
00588    towersUpperValid = true;
00589    for( CaloTowerCollection::const_iterator tow = caloTowersCleanerUpperR45->begin(); tow != caloTowersCleanerUpperR45->end(); tow++){
00590      towersUpper.insert(tow->id().denseIndex());
00591    }
00592  }
00593  if( caloTowersCleanerLowerR45.isValid() ){
00594    towersLowerValid = true;
00595    for( CaloTowerCollection::const_iterator tow = caloTowersCleanerLowerR45->begin(); tow != caloTowersCleanerLowerR45->end(); tow++){
00596      towersLower.insert(tow->id().denseIndex());
00597    }
00598  }
00599  if( caloTowersCleanerNoR45.isValid() ){
00600    towersNoneValid = true;
00601    for( CaloTowerCollection::const_iterator tow = caloTowersCleanerNoR45->begin(); tow != caloTowersCleanerNoR45->end(); tow++){
00602      towersNone.insert(tow->id().denseIndex());
00603    }
00604  }
00605  if (caloTowers.isValid()) {
00606    //    ntowcal = caloTowers->size();
00607    int jtow = 0;
00608    for ( CaloTowerCollection::const_iterator tower=caloTowers->begin(); tower!=caloTowers->end(); tower++) {
00609      if(tower->energy() > thresholdForSavingTowers)
00610        {
00611          towet[jtow] = tower->et();
00612          toweta[jtow] = tower->eta();
00613          towphi[jtow] = tower->phi();
00614          towen[jtow] = tower->energy();
00615          towem[jtow] = tower->emEnergy();
00616          towhd[jtow] = tower->hadEnergy();
00617          towoe[jtow] = tower->outerEnergy();
00618          // noise filters: true = no noise, false = noise
00619          if(towersUpperValid) {if(towersUpper.find(tower->id().denseIndex()) == towersUpper.end()) towR45upper[jtow]=true; else towR45upper[jtow]=false;}
00620          if(towersLowerValid) {if(towersLower.find(tower->id().denseIndex()) == towersLower.end()) towR45lower[jtow]=true; else towR45lower[jtow]=false;}
00621          if(towersNoneValid) {if(towersNone.find(tower->id().denseIndex()) == towersNone.end()) towR45none[jtow]=true; else towR45none[jtow]=false;}
00622          jtow++;
00623        }
00624    }
00625    ntowcal = jtow;
00626  }
00627  else {ntowcal = 0;}
00628     
00629     if (recmets.isValid()) {
00630         typedef reco::CaloMETCollection::const_iterator cmiter;
00631         for ( cmiter i=recmets->begin(); i!=recmets->end(); i++) {
00632             mcalmet = i->pt();
00633             mcalphi = i->phi();
00634             mcalsum = i->sumEt();
00635         }
00636     }
00637 
00638     if (pfmets.isValid()) {
00639       typedef reco::PFMETCollection::const_iterator pfmetiter;
00640       for( pfmetiter i=pfmets->begin(); i!=pfmets->end(); i++) {
00641         pfmet = i->pt();
00642         pfsumet = i->sumEt();
00643         pfmetphi = i->phi();
00644       }
00645     }
00646     
00647     if (ht.isValid()) {
00648         typedef reco::METCollection::const_iterator iter;
00649         for ( iter i=ht->begin(); i!=ht->end(); i++) {
00650             htcalet = i->pt();
00651             htcalphi = i->phi();
00652             htcalsum = i->sumEt();
00653         }
00654     }
00655     
00656     if (_Monte){
00657         
00658         if (genjets.isValid()) {
00659             reco::GenJetCollection mygenjets;
00660             mygenjets=*genjets;
00661             std::sort(mygenjets.begin(),mygenjets.end(),PtGreater());
00662             typedef reco::GenJetCollection::const_iterator gjiter;
00663             int jgen=0;
00664             for ( gjiter i=mygenjets.begin(); i!=mygenjets.end(); i++) {
00665                 
00666                 if (i->pt()>_GenJetMin){
00667                     jgenpt[jgen] = i->pt();
00668                     jgenphi[jgen] = i->phi();
00669                     jgeneta[jgen] = i->eta();
00670                     jgene[jgen] = i->energy();
00671                     jgen++;
00672                 }
00673                 
00674             }
00675             njetgen = jgen;
00676         }
00677         else {njetgen = 0;}
00678         
00679         if (genmets.isValid()) {
00680             typedef reco::GenMETCollection::const_iterator gmiter;
00681             for ( gmiter i=genmets->begin(); i!=genmets->end(); i++) {
00682                 mgenmet = i->pt();
00683                 mgenphi = i->phi();
00684                 mgensum = i->sumEt();
00685             }
00686         }
00687         
00688     }
00689     
00690     
00692     if (l2taujets.isValid()) {
00693         nohl2tau = l2taujets->size();
00694         reco::CaloJetCollection l2taus = *l2taujets;
00695         std::sort(l2taus.begin(),l2taus.end(),GetPFPtGreater());
00696         int itau=0;
00697         for(reco::CaloJetCollection::const_iterator i = l2taus.begin(); 
00698                                                    i!= l2taus.end(); ++i){
00699             l2tauPt[itau]  = i->pt();
00700             l2tauEta[itau] = i->eta();
00701             l2tauPhi[itau] = i->phi();
00702             itau++;
00703         }
00704     }else{
00705       nohl2tau = 0;
00706     }
00707     if (taujets.isValid()) {      
00708         nohtau = taujets->size();
00709         reco::HLTTauCollection mytaujets;
00710         mytaujets=*taujets;
00711         std::sort(mytaujets.begin(),mytaujets.end(),GetPtGreater());
00712         typedef reco::HLTTauCollection::const_iterator tauit;
00713         int itau=0;
00714         for(tauit i=mytaujets.begin(); i!=mytaujets.end(); i++){
00715             //Ask for Eta,Phi and Et of the tau:
00716             tauEta[itau] = i->getEta();
00717             tauPhi[itau] = i->getPhi();
00718             tauPt[itau] = i->getPt();
00719             //Ask for L2 EMIsolation cut: Nominal cut : < 5
00720             l2tauemiso[itau] = i->getEMIsolationValue();
00721             //Get L25 LeadTrackPt : Nominal cut : > 20 GeV
00722             l25tauPt[itau] = i->getL25LeadTrackPtValue();
00723             //Get TrackIsolation response (returns 0 = failed or 1= passed)
00724             l3tautckiso[itau] = i->getL3TrackIsolationResponse();
00725             //MET : > 65
00726             itau++;
00727         }      
00728     }
00729     else {nohtau = 0;}
00730 
00731     
00733     if(pfTaus.isValid()) {
00734         //float minTrkPt = minPtCH;
00735         //float minGammaPt = minPtGamma;
00736         nohPFTau  = pfTaus->size();
00737         reco::PFTauCollection taus = *pfTaus;
00738         std::sort(taus.begin(),taus.end(),GetPFPtGreater());
00739         typedef reco::PFTauCollection::const_iterator pftauit;
00740         int ipftau=0;
00741         for(pftauit i=taus.begin(); i!=taus.end(); i++){
00742             //Ask for Eta,Phi and Et of the tau:
00743             ohpfTauProngs[ipftau] = i->signalPFChargedHadrCands().size();
00744             ohpfTauEta[ipftau] = i->eta();
00745             ohpfTauPhi[ipftau] = i->phi();
00746             ohpfTauPt[ipftau] = i->pt();
00747             ohpfTauJetPt[ipftau] = i->pfTauTagInfoRef()->pfjetRef()->pt();
00748 
00749             
00750   /*
00751             if( (i->leadPFCand()).isNonnull())
00752                 pfTauLeadPionPt[ipftau] = i->leadPFCand()->pt();            
00753 */
00754             if( (i->leadPFNeutralCand()).isNonnull())
00755                 ohpfTauLeadPionPt[ipftau] = i->leadPFNeutralCand()->pt();        
00756             else 
00757                 ohpfTauLeadPionPt[ipftau] = -999.0;
00758 
00759             if((i->leadPFChargedHadrCand()).isNonnull()){
00760                 ohpfTauLeadTrackPt[ipftau] = i->leadPFChargedHadrCand()->pt();
00761                 ohpfTauLeadTrackVtxZ[ipftau] = i->leadPFChargedHadrCand()->vertex().z(); 
00762             }else{
00763                 ohpfTauLeadTrackPt[ipftau] = -999.0;
00764                 ohpfTauLeadTrackVtxZ[ipftau] = -999.0;  
00765             }
00766 
00767             float maxPtTrkIso = 0;
00768             for (unsigned int iTrk = 0; iTrk < i->isolationPFChargedHadrCands().size(); iTrk++)
00769             {
00770                 if(i->isolationPFChargedHadrCands()[iTrk]->pt() > maxPtTrkIso) maxPtTrkIso = i->isolationPFChargedHadrCands()[iTrk]->pt();
00771 
00772                 if (i->isolationPFChargedHadrCands()[iTrk]->trackRef().isNonnull()){
00773                   hltpftauIsoTrToPFTauMatch[noHLTPFTausIso]=ipftau;
00774                   HLTPFTauIsoTrDz[noHLTPFTausIso]=i->isolationPFChargedHadrCands()[iTrk]->trackRef()->dz(); // dz wrt (0,0,0), to compare offline with HLT 
00775                   HLTPFTauIsoTrPt[noHLTPFTausIso]=i->isolationPFChargedHadrCands()[iTrk]->pt();
00776                   /*
00777                   std::cout << "Adding isocand for hltpftau " << ipftau
00778                       << " pt " << HLTPFTauIsoTrPt[noHLTPFTausIso]
00779                       << " dz " << HLTPFTauIsoTrDz[noHLTPFTausIso]
00780                       << std::endl; // */
00781                   ++noHLTPFTausIso;
00782                 }
00783 
00784             }
00785                 
00786             ohpfTauTrkIso[ipftau] = maxPtTrkIso;
00787             float maxPtGammaIso = 0;
00788             for (unsigned int iGamma = 0; iGamma < i->isolationPFGammaCands().size(); iGamma++)
00789             {
00790                 if(i->isolationPFGammaCands()[iGamma]->pt() > maxPtGammaIso) maxPtGammaIso = i->isolationPFGammaCands()[iGamma]->pt();
00791             }                        
00792 
00793 
00794 
00795             for (unsigned int iTrk = 0; iTrk < i->signalPFChargedHadrCands().size(); iTrk++)
00796             {
00797               if (i->signalPFChargedHadrCands ()[iTrk]->trackRef().isNonnull()){
00798                 hltpftauSignalTrToPFTauMatch[noHLTPFTausSignal]=ipftau;
00799                 HLTPFTauSignalTrDz[noHLTPFTausSignal]=i->signalPFChargedHadrCands()[iTrk]->trackRef()->dz(); // dz wrt (0,0,0), to compare offline with HLT
00800                 HLTPFTauSignalTrPt[noHLTPFTausSignal]=i->signalPFChargedHadrCands()[iTrk]->pt();
00801                 /*
00802                   std::cout << "Adding sigcand for hltpftau " << ipftau
00803                       << " pt " << HLTPFTauSignalTrPt[noHLTPFTausSignal]
00804                       << " dz " << HLTPFTauSignalTrDz[noHLTPFTausSignal]
00805                       << std::endl; // */
00806                 ++noHLTPFTausSignal;
00807               }
00808             }
00809 
00810 
00811 
00812 
00813 
00814             ohpfTauGammaIso[ipftau] = maxPtGammaIso;
00815             ipftau++;
00816         } 
00817       
00818     }
00819 
00820     if(pfTausTightCone.isValid()) {
00821         //float minTrkPt = minPtCH;
00822         //float minGammaPt = minPtGamma;
00823         nohPFTauTightCone = pfTaus->size();
00824         reco::PFTauCollection taus = *pfTausTightCone;
00825         std::sort(taus.begin(),taus.end(),GetPFPtGreater());
00826         typedef reco::PFTauCollection::const_iterator pftauit;
00827         int ipftau=0;
00828         for(pftauit i=taus.begin(); i!=taus.end(); i++){
00829             //Ask for Eta,Phi and Et of the tau:
00830             ohpfTauTightConeProngs[ipftau] = i->signalPFChargedHadrCands().size();
00831             ohpfTauTightConeEta[ipftau] = i->eta();
00832             ohpfTauTightConePhi[ipftau] = i->phi();
00833             ohpfTauTightConePt[ipftau] = i->pt();
00834             ohpfTauTightConeJetPt[ipftau] = i->pfTauTagInfoRef()->pfjetRef()->pt();
00835 
00836     
00837             if( (i->leadPFNeutralCand()).isNonnull())
00838                 ohpfTauTightConeLeadPionPt[ipftau] = i->leadPFNeutralCand()->pt();
00839             else 
00840               ohpfTauTightConeLeadPionPt[ipftau] = -999.0;
00841 
00842 
00843             if((i->leadPFChargedHadrCand()).isNonnull())
00844                 ohpfTauTightConeLeadTrackPt[ipftau] = i->leadPFChargedHadrCand()->pt();
00845             else
00846               ohpfTauTightConeLeadTrackPt[ipftau] = -999.0;
00847 
00848             float maxPtTrkIso = 0;
00849             for (unsigned int iTrk = 0; iTrk < i->isolationPFChargedHadrCands().size(); iTrk++)
00850             {
00851                 if(i->isolationPFChargedHadrCands()[iTrk]->pt() > maxPtTrkIso) maxPtTrkIso = i->isolationPFChargedHadrCands()[iTrk]->pt();
00852             }
00853 
00854             ohpfTauTightConeTrkIso[ipftau] = maxPtTrkIso;
00855             float maxPtGammaIso = 0;
00856             for (unsigned int iGamma = 0; iGamma < i->isolationPFGammaCands().size(); iGamma++)
00857             {
00858                 if(i->isolationPFGammaCands()[iGamma]->pt() > maxPtGammaIso) maxPtGammaIso = i->isolationPFGammaCands()[iGamma]->pt();
00859             }
00860             ohpfTauTightConeGammaIso[ipftau] = maxPtGammaIso;
00861             ipftau++;
00862         }
00863 
00864     }
00865     
00867       
00868     if(recoPfTaus.isValid()) {
00869         float minTrkPt = minPtCH;
00870         float minGammaPt = minPtGamma;
00871         nRecoPFTau  = recoPfTaus->size();
00872         reco::PFTauCollection taus = *recoPfTaus;
00873 
00874         // disable sorting for proper access to discriminators
00875         //std::sort(taus.begin(),taus.end(),GetPFPtGreater());
00876         typedef reco::PFTauCollection::const_iterator pftauit;
00877         int ipftau=0;
00878         
00879         for(pftauit i=taus.begin(); i!=taus.end(); i++){
00880             //Ask for Eta,Phi and Et of the tau:
00881             recopfTauEta[ipftau] = i->eta();
00882             recopfTauPhi[ipftau] = i->phi();
00883             recopfTauPt[ipftau]  = i->pt();
00884 
00885             if( (i->leadPFNeutralCand()).isNonnull())
00886                 recopfTauLeadPionPt[ipftau] = i->leadPFNeutralCand()->pt();
00887             else
00888               recopfTauLeadPionPt[ipftau] = -999.0;  
00889 
00890 
00891             if((i->leadPFChargedHadrCand()).isNonnull())
00892                 recopfTauLeadTrackPt[ipftau] = i->leadPFChargedHadrCand()->pt();
00893             else
00894               recopfTauLeadTrackPt[ipftau]  = -999.0;
00895 
00896             int myTrks=0;
00897             for (unsigned int iTrk = 0; iTrk < i->isolationPFChargedHadrCands().size(); iTrk++)
00898             {
00899                 if(i->isolationPFChargedHadrCands()[iTrk]->pt() > minTrkPt) myTrks++;
00900                 if (i->isolationPFChargedHadrCands()[iTrk]->trackRef().isNonnull()){
00901                   isoTrToPFTauMatch[noRecoPFTausIso]=ipftau;
00902                   recoPFTauIsoTrDz[noRecoPFTausIso]=i->isolationPFChargedHadrCands()[iTrk]->trackRef()->dz(); // dz wrt (0,0,0), to compare offline with HLT
00903                   recoPFTauIsoTrPt[noRecoPFTausIso]=i->isolationPFChargedHadrCands()[iTrk]->pt();
00904                   /*
00905                   std::cout << "Adding isocand for tau " << ipftau
00906                             << " pt " << recoPFTauIsoTrPt[noRecoPFTausIso]
00907                             << " dz " << recoPFTauIsoTrDz[noRecoPFTausIso]
00908                             << std::endl;// */
00909                   ++noRecoPFTausIso;
00910                 }
00911 
00912             }
00913                
00914             recopfTauTrkIso[ipftau] = myTrks;
00915             int myGammas=0;
00916             for (unsigned int iGamma = 0; iGamma < i->isolationPFGammaCands().size(); iGamma++)
00917             {
00918                 if(i->isolationPFGammaCands()[iGamma]->pt() > minGammaPt) myGammas++;
00919             }                        
00920             recopfTauGammaIso[ipftau] = myGammas;
00921             
00922 
00923             for (unsigned int iTrk = 0; iTrk < i->signalPFChargedHadrCands().size(); iTrk++)
00924             {
00925                 if (i->signalPFChargedHadrCands ()[iTrk]->trackRef().isNonnull()){
00926                   signalTrToPFTauMatch[noRecoPFTausSignal]=ipftau;
00927                   recoPFTauSignalTrDz[noRecoPFTausSignal]=i->signalPFChargedHadrCands()[iTrk]->trackRef()->dz(); // dz wrt (0,0,0), to compare offline with HLT
00928                   recoPFTauSignalTrPt[noRecoPFTausSignal]=i->signalPFChargedHadrCands()[iTrk]->pt();
00929                   /*
00930                   std::cout << "Adding sigcand for tau " << ipftau
00931                             << " pt " << recoPFTauSignalTrPt[noRecoPFTausSignal]
00932                             << " dz " << recoPFTauSignalTrDz[noRecoPFTausSignal]
00933                             << std::endl;// */
00934                   ++noRecoPFTausSignal;
00935                 }
00936             }
00937 
00938             const reco::PFTauRef thisTauRef(recoPfTaus,ipftau);
00939             
00940             if(theRecoPFTauDiscrByTanCOnePercent.isValid()){
00941             recopfTauDiscrByTancOnePercent[ipftau] = (*theRecoPFTauDiscrByTanCOnePercent)[thisTauRef];}
00942             if(theRecoPFTauDiscrByIsolation.isValid()){ 
00943             recopfTauDiscrByIso[ipftau] = (*theRecoPFTauDiscrByIsolation)[thisTauRef];} 
00944             if(theRecoPFTauDiscrAgainstMuon.isValid()){
00945             recopfTauDiscrAgainstMuon[ipftau] = (*theRecoPFTauDiscrAgainstMuon)[thisTauRef];}
00946             if(theRecoPFTauDiscrAgainstElec.isValid()){
00947             recopfTauDiscrAgainstElec[ipftau] = (*theRecoPFTauDiscrAgainstElec)[thisTauRef];}
00948             if(theRecoPFTauDiscrByTanCHalfPercent.isValid()){
00949             recopfTauDiscrByTancHalfPercent[ipftau] = (*theRecoPFTauDiscrByTanCHalfPercent)[thisTauRef];}
00950             if(theRecoPFTauDiscrByTanCQuarterPercent.isValid()){
00951             recopfTauDiscrByTancQuarterPercent[ipftau] = (*theRecoPFTauDiscrByTanCQuarterPercent)[thisTauRef];}
00952             if(theRecoPFTauDiscrByTanCTenthPercent.isValid()){
00953             recopfTauDiscrByTancTenthPercent[ipftau] = (*theRecoPFTauDiscrByTanCTenthPercent)[thisTauRef];}
00954 
00955             ipftau++;
00956         }        
00957     }
00958    
00960     if(pfJets.isValid()) {
00961         nohPFJet  = pfJets->size();
00962         reco::PFJetCollection Jets = *pfJets;
00963         std::sort(Jets.begin(),Jets.end(),GetPFPtGreater());
00964         typedef reco::PFJetCollection::const_iterator pfJetit;
00965         int ipfJet=0;
00966         float pfMHTx = 0.;
00967         float pfMHTy = 0.;
00968         pfHT         = 0.;
00969 
00970         for(pfJetit i=Jets.begin(); i!=Jets.end(); i++){
00971             //Ask for Eta,Phi and Et of the Jet:
00972             pfJetEta[ipfJet] = i->eta();
00973             pfJetPhi[ipfJet] = i->phi();
00974             pfJetPt[ipfJet] = i->pt();           
00975             pfJetE[ipfJet] = i->energy();
00976             pfJetneutralHadronEnergyFraction[ipfJet]=i->neutralHadronEnergyFraction();
00977             pfJetchargedHadronFraction[ipfJet] = i->chargedHadronEnergyFraction ();
00978             pfJetneutralMultiplicity[ipfJet] = i->neutralMultiplicity ();
00979             pfJetchargedMultiplicity[ipfJet] = i->chargedMultiplicity ();
00980             pfJetneutralEMFraction[ipfJet] = i->neutralEmEnergyFraction ();
00981             pfJetchargedEMFraction[ipfJet] = i->chargedEmEnergyFraction ();
00982             //std::cout << "jet pT = " << i->pt() << " ; neutralHadronEnergyFraction = " << i->neutralHadronEnergyFraction() << std::endl;
00983 
00984             if (i->pt() > 40. && abs(i->eta())<3.0)
00985               pfHT  += i -> pt();
00986             if (i->pt() > 30.){
00987               pfMHTx = pfMHTx + i->px();
00988               pfMHTy = pfMHTy + i->py();
00989             }
00990             ipfJet++;   
00991         } 
00992         pfMHT = sqrt(pfMHTx*pfMHTx + pfMHTy*pfMHTy);
00993         
00994     }
00996     nrpj = 0;
00997     if(recoPFJets.isValid()){
00998             nrpj = recoPFJets->size();
00999             reco::PFJetCollection Jets = *recoPFJets;
01000             std::sort(Jets.begin(),Jets.end(),GetPFPtGreater());
01001             typedef reco::PFJetCollection::const_iterator pfJetit;
01002             int ipfJet=0;
01003             for(pfJetit i=Jets.begin(); i!=Jets.end(); i++){
01004                     //Ask for Eta,Phi and Et of the Jet:
01005                     jpfrecoeta[ipfJet] = i->eta();
01006                     jpfrecophi[ipfJet] = i->phi();
01007                     jpfrecopt[ipfJet] = i->pt();           
01008                     jpfrecoe[ipfJet] = i->energy();           
01009                     jpfreconeutralHadronFraction[ipfJet] = i->neutralHadronEnergyFraction ();
01010                     jpfrecochargedHadronFraction[ipfJet] = i->chargedHadronEnergyFraction ();
01011                     jpfreconeutralMultiplicity[ipfJet] = i->neutralMultiplicity ();
01012                     jpfrecochargedMultiplicity[ipfJet] = i->chargedMultiplicity ();
01013                     jpfreconeutralEMFraction[ipfJet] = i->neutralEmEnergyFraction ();
01014                     jpfrecochargedEMFraction[ipfJet] = i->chargedEmEnergyFraction ();
01015 
01016                     ipfJet++;   
01017             } 
01018 
01019     }
01020     
01021 }