CMS 3D CMS Logo

TruthTauDecayModeProducer.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    TruthTauDecayModeProducer
00004 // Class:      TruthTauDecayModeProducer
00005 // 
00015 //
00016 // Original Author:  Evan K. Friis, UC Davis (friis@physics.ucdavis.edu)
00017 //         Created:  Thu Sep 1 06:19:05 PST 2008
00018 // $Id: TruthTauDecayModeProducer.cc,v 1.2.2.1 2008/11/30 01:02:44 friis Exp $
00019 //
00020 //
00021 
00022 
00023 // system include files
00024 #include <memory>
00025 
00026 // user include files
00027 #include "FWCore/Framework/interface/Frameworkfwd.h"
00028 #include "FWCore/Framework/interface/EDProducer.h"
00029 
00030 #include "FWCore/Framework/interface/Event.h"
00031 #include "FWCore/Framework/interface/MakerMacros.h"
00032 #include "DataFormats/Candidate/interface/Candidate.h"
00033 #include "DataFormats/TauReco/interface/PFTauDecayMode.h"
00034 #include "RecoTauTag/TauTagTools/interface/GeneratorTau.h"
00035 #include "DataFormats/Candidate/interface/Candidate.h"
00036 #include "DataFormats/JetReco/interface/GenJet.h"
00037 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00038 #include "PhysicsTools/CandUtils/interface/AddFourMomenta.h"
00039 
00040 #include <vector>
00041 
00042 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00043 
00044 class TruthTauDecayModeProducer : public edm::EDProducer {
00045    public:
00046       struct tauObjectsHolder {
00047          std::vector<const reco::Candidate*> chargedObjects;
00048          std::vector<const reco::Candidate*> neutralObjects;
00049       };
00050       explicit TruthTauDecayModeProducer(const edm::ParameterSet&);
00051       ~TruthTauDecayModeProducer();
00052 
00053    private:
00054       virtual void beginJob(const edm::EventSetup&) ;
00055       virtual void produce(edm::Event&, const edm::EventSetup&);
00056       virtual void endJob() ;
00057 
00058       //for signal, the module takes input from a PdgIdAndStatusCandViewSelector 
00059       //for background, the module takes input from a collection of GenJets
00060       bool              iAmSignal_; 
00061       edm::InputTag     inputTag_;
00062       double            leadTrackPtCut_;
00063       double            leadTrackEtaCut_;
00064       double            totalPtCut_;
00065       double            totalEtaCut_;
00066       AddFourMomenta    addP4;
00067 };
00068 
00069 TruthTauDecayModeProducer::TruthTauDecayModeProducer(const edm::ParameterSet& iConfig)
00070 {
00071    edm::LogInfo("TruthTauDecayModeProducer") << "Initializing ctor of TruthTauDecayModeProducer";
00072    iAmSignal_           = iConfig.getParameter<bool>("iAmSignal");
00073    inputTag_            = iConfig.getParameter<edm::InputTag>("inputTag");
00074    leadTrackPtCut_      = iConfig.getParameter<double>("leadTrackPtCut");
00075    leadTrackEtaCut_     = iConfig.getParameter<double>("leadTrackEtaCut");
00076    totalPtCut_          = iConfig.getParameter<double>("totalPtCut");
00077    totalEtaCut_         = iConfig.getParameter<double>("totalEtaCut");
00078    //register your products
00079    edm::LogInfo("TruthTauDecayModeProducer") << "Registering products";
00080    produces<std::vector<reco::PFTauDecayMode> >();
00081    edm::LogInfo("TruthTauDecayModeProducer") << "TruthTauDecayModeProducer initialized";
00082 }
00083 
00084 
00085 TruthTauDecayModeProducer::~TruthTauDecayModeProducer()
00086 {
00087 }
00088 
00089 void
00090 TruthTauDecayModeProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
00091 {
00092    using namespace edm;
00093    using namespace std;
00094    using namespace reco;
00095 
00096    vector<tauObjectsHolder> tausToAdd_;
00097 
00098    /* **********************************************
00099     * **********  True Tau Case          ***********
00100     * ********************************************** */
00101    if (iAmSignal_)
00102    {
00103       Handle<RefToBaseVector<reco::Candidate> > decayedMCTaus;
00104       iEvent.getByLabel(inputTag_, decayedMCTaus);
00105       for(edm::RefToBaseVector<reco::Candidate>::const_iterator iterGen = decayedMCTaus->begin();
00106             iterGen != decayedMCTaus->end();
00107             ++iterGen)
00108       {
00109          //copy into custom format (this casting is bullshit)
00110          GeneratorTau tempTau = static_cast<const GenParticle&>(*(*iterGen));
00111          //LogInfo("MCTauCandidateProducer") << "Generator tau produced, initializing.. ";
00112          tempTau.init();
00113          //LogInfo("MCTauCandidateProducer") << "GenTau initialization done";
00114          if (tempTau.isFinalStateTau())
00115          {
00116             // Build a Tau Candidate from the information contained in the parsing class
00117             tauObjectsHolder tempTauHolder;
00118             tempTauHolder.chargedObjects = tempTau.getGenChargedPions();
00119             tempTauHolder.neutralObjects = tempTau.getGenNeutralPions();
00120             tausToAdd_.push_back(tempTauHolder);
00121          }
00122       }
00123    } else
00124    {
00125    /* **********************************************
00126     * **********  QCD Case               ***********
00127     * ********************************************** */
00128       Handle<GenJetCollection> genJets;
00129       iEvent.getByLabel(inputTag_, genJets);
00130       for(GenJetCollection::const_iterator aGenJet = genJets->begin(); aGenJet != genJets->end(); ++aGenJet)
00131       {
00132          // get all constituents
00133          vector<const GenParticle*> theJetConstituents = aGenJet->getGenConstituents();
00134 
00135          tauObjectsHolder tempTauHolder;
00136          // filter the constituents
00137          for( vector<const GenParticle*>::const_iterator aCandidate = theJetConstituents.begin();
00138                aCandidate != theJetConstituents.end();
00139                ++aCandidate)
00140          {
00141             int pdgId = abs((*aCandidate)->pdgId());
00142             const Candidate* theCandidate = static_cast<const Candidate*>(*aCandidate);
00143             //filter nus
00144             if (pdgId == 16 || pdgId == 12 || pdgId == 14)
00145             {
00146                //do nothing
00147             } else 
00148             {
00149                // call everything charged a pion
00150                // call everything neutral a neutral pion
00151                if (theCandidate->charge() != 0)
00152                   tempTauHolder.chargedObjects.push_back(theCandidate);
00153                else
00154                   tempTauHolder.neutralObjects.push_back(theCandidate);
00155             }
00156          }
00157          tausToAdd_.push_back(tempTauHolder);
00158       }
00159    }
00160 
00161    //output collection
00162    std::auto_ptr<vector<PFTauDecayMode> > pOut( new vector<PFTauDecayMode> );
00163    for(vector<tauObjectsHolder>::const_iterator iTempTau  = tausToAdd_.begin();
00164                                                 iTempTau != tausToAdd_.end();
00165                                               ++iTempTau)
00166    {
00167       double leadTrackPt  = 0.;
00168       double leadTrackEta = 0.;
00169       VertexCompositeCandidate chargedObjectsToAdd;
00170       const vector<const Candidate*>* chargedObjects = &(iTempTau->chargedObjects);
00171       for(vector<const Candidate*>::const_iterator iCharged  = chargedObjects->begin();
00172                                                    iCharged != chargedObjects->end();
00173                                                  ++iCharged)
00174       {
00175          chargedObjectsToAdd.addDaughter(**iCharged);
00176          double trackPt = (*iCharged)->pt();
00177          if (trackPt > leadTrackPt)
00178          {
00179             leadTrackPt  = trackPt;
00180             leadTrackEta = (*iCharged)->eta();
00181          }
00182       }
00183       //update the composite four vector
00184       addP4.set(chargedObjectsToAdd);
00185 
00186       CompositeCandidate neutralPionsToAdd;
00187       const vector<const Candidate*>* neutralObjects = &(iTempTau->neutralObjects);
00188       for(vector<const Candidate*>::const_iterator iNeutral  = neutralObjects->begin();
00189                                                    iNeutral != neutralObjects->end();
00190                                                  ++iNeutral)
00191       {
00192          neutralPionsToAdd.addDaughter(**iNeutral);
00193       }
00194       addP4.set(neutralPionsToAdd);
00195 
00196       Particle::LorentzVector myFourVector = chargedObjectsToAdd.p4();
00197       myFourVector += neutralPionsToAdd.p4();
00198 
00199       if(leadTrackPt > leadTrackPtCut_ && abs(leadTrackEta) < leadTrackEtaCut_ && myFourVector.pt() > totalPtCut_ && abs(myFourVector.eta()) < totalEtaCut_)
00200       {
00201          //TODO: add vertex fitting
00202          CompositeCandidate theOutliers;
00203          PFTauRef             noPFTau;     //empty REF to PFTau
00204          PFTauDecayMode decayModeToAdd(chargedObjectsToAdd, neutralPionsToAdd, theOutliers);
00205          decayModeToAdd.setPFTauRef(noPFTau);
00206          pOut->push_back(decayModeToAdd);
00207       }
00208    }
00209    iEvent.put(pOut);
00210 }
00211 
00212 void 
00213 TruthTauDecayModeProducer::beginJob(const edm::EventSetup&)
00214 {
00215 }
00216 
00217 // ------------ method called once each job just after ending the event loop  ------------
00218 void 
00219 TruthTauDecayModeProducer::endJob() {
00220 }
00221 
00222 //define this as a plug-in
00223 DEFINE_FWK_MODULE(TruthTauDecayModeProducer);

Generated on Tue Jun 9 17:45:09 2009 for CMSSW by  doxygen 1.5.4