EDAnalyzer to exercise and demonstrate usage of Calibratable tree. More...

#include <CalibratableTest.h>

Inheritance diagram for CalibratableTest:

Public Member Functions
	CalibratableTest (const edm::ParameterSet &)

template<class T >
void	getCollection (edm::Handle< T > &c, const edm::InputTag &tag, const edm::Event &event) const

	~CalibratableTest ()

Public Member Functions inherited from edm::EDAnalyzer
	EDAnalyzer ()

std::string	workerType () const

virtual	~EDAnalyzer ()

Private Member Functions
virtual void	analyze (const edm::Event &, const edm::EventSetup &)

virtual void	beginJob ()

double	deltaR (const double &eta1, const double &eta2, const double &phi1, const double &phi2)

virtual void	endJob ()

void	extractCandidate (const reco::PFCandidate &cand)

virtual void	fillTreeAndReset ()

std::vector< unsigned >	findCandidatesInDeltaR (const reco::PFSimParticle &pft, const std::vector< reco::PFCandidate > &cands, const double &deltaR)

std::vector< unsigned >	findPrimarySimParticles (const std::vector< reco::PFSimParticle > &sims)

Private Attributes
pftools::Calibratable *	calib_

edm::Handle < reco::PFClusterCollection > *	clustersEcal_

edm::Handle < reco::PFClusterCollection > *	clustersHcal_

int	debug_

double	deltaRCandToSim_

edm::Service< TFileService >	fileservice_

edm::InputTag	inputTagCandidates_

edm::InputTag	inputTagClustersEcal_

edm::InputTag	inputTagClustersHcal_

edm::InputTag	inputTagSimParticles_

unsigned	nEventFails_

unsigned	nEventWrites_

unsigned	nParticleFails_

unsigned	nParticleWrites_

edm::Handle < reco::PFCandidateCollection > *	pfCandidates_

edm::Handle < reco::PFSimParticleCollection > *	simParticles_

bool	thisEventPasses_

bool	thisParticlePasses_

TTree *	tree_

Additional Inherited Members
Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer	ModuleType

typedef WorkerT< EDAnalyzer >	WorkerType

Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string &	baseType ()

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &)

Protected Member Functions inherited from edm::EDAnalyzer
CurrentProcessingContext const *	currentContext () const

Detailed Description

EDAnalyzer to exercise and demonstrate usage of Calibratable tree.

Author: Jamie Ballin, Imperial College London

Date: November 2008

OBJECTIVE: this analyzer will create a tree of PFClusterTools/Calibratable instances. Each entry of the tree represents particle flow information relating to one pion in a multi-pion event. The use of monte carlo is assumed as sim particle information is extracted too.

USAGE: This is an analyzer, not a producer: the tree it produces is stored in a seperate ROOT file using the CMSSW TFileService.

Consult test/CalibratableTest_FastSim.py for a sample CMSSW configuration. Run the example from that directory with, cmsRun CalibratableTest_FastSim.py This will produce a file with a tree inside. Browse the contents in bare root with a TBrowser.

NOTE: This analyzer does not exercise the complete functionality of the Calibratable class. For instance, tracks and rechits are not stored. Contact me (the author) for more details, and for an analyzer that we use to extract particle information to calibrate the particle flow algorithm. A MUCH more involved example from my private analysis is to be found in, UserCode/JamieBallin/interface/DipionDelegate.h

General details about the usage of Calibratable may be found at, https://twiki.cern.ch/twiki/bin/view/CMS/PFClusterToolsPackage

Definition at line 57 of file CalibratableTest.h.

Constructor & Destructor Documentation

CalibratableTest::CalibratableTest ( const edm::ParameterSet & parameters )

explicit

Definition at line 28 of file CalibratableTest.cc.

References calib_, gather_cfg::cout, debug_, deltaRCandToSim_, fileservice_, edm::ParameterSet::getParameter(), inputTagCandidates_, inputTagClustersEcal_, inputTagClustersHcal_, inputTagSimParticles_, TFileDirectory::make(), and tree_.

                                                                     :
     debug_(0), nParticleWrites_(0), nParticleFails_(0), nEventWrites_(0), nEventFails_(0),
              deltaRCandToSim_(0.4) {
 
     std::cout << __PRETTY_FUNCTION__ << std::endl;
     
     /* This procedure is GENERIC to storing any dictionary enable class in a ROOT tree. */
     tree_ = fileservice_->make<TTree>("CalibratableTest", "");
     calib_ = new Calibratable();
     tree_->Branch("Calibratable", "pftools::Calibratable", &calib_, 32000, 2);
     
     inputTagCandidates_= parameters.getParameter<InputTag>("PFCandidates");
     inputTagSimParticles_= parameters.getParameter<InputTag>("PFSimParticles");
     inputTagClustersEcal_= parameters.getParameter<InputTag>("PFClustersEcal");
     inputTagClustersHcal_= parameters.getParameter<InputTag>("PFClustersHcal");
     deltaRCandToSim_ = parameters.getParameter<double>("deltaRCandToSim");
     debug_= parameters.getParameter<int>("debug");
 }

CalibratableTest::~CalibratableTest ( )

Definition at line 47 of file CalibratableTest.cc.

                                     {
 
 }

Member Function Documentation

void CalibratableTest::analyze	(	const edm::Event &	event,
		const edm::EventSetup &	iSetup
	)

privatevirtual

Implements edm::EDAnalyzer.

Definition at line 57 of file CalibratableTest.cc.

                                      {
     if (debug_ > 1)
         std::cout << __PRETTY_FUNCTION__ << "\n";
 
     //Extract new collection references
     pfCandidates_ = new Handle<PFCandidateCollection>;
     simParticles_ = new Handle<PFSimParticleCollection>;
     clustersEcal_ = new Handle<PFClusterCollection>;
     clustersHcal_ = new Handle<PFClusterCollection>;
 
     getCollection(*pfCandidates_, inputTagCandidates_, event);
     getCollection(*simParticles_, inputTagSimParticles_, event);
     getCollection(*clustersEcal_, inputTagClustersEcal_, event);
     getCollection(*clustersHcal_, inputTagClustersHcal_, event);
 
     //Reset calibratable branch
     thisEventPasses_ = true;
     thisParticlePasses_ = true;
     calib_->reset();
 
     if (debug_ > 1)
         std::cout << "\tStarting event..."<< std::endl;
 
     PFSimParticleCollection sims = **simParticles_;
     PFCandidateCollection candidates = **pfCandidates_;
     PFClusterCollection clustersEcal = **clustersEcal_;
     PFClusterCollection clustersHcal = **clustersHcal_;
 
     if (sims.size() == 0) {
         std::cout << "\tNo sim particles found!" << std::endl;
         thisEventPasses_ = false;
     }
 
     //Find primary pions in the event
     std::vector<unsigned> primarySims = findPrimarySimParticles(sims);
     if (debug_) {
         std::cout << "\tFound "<< primarySims.size()
                 << " primary sim particles, "<< (**pfCandidates_).size() << " pfCandidates\n";
     }
     for (std::vector<unsigned>::const_iterator cit = primarySims.begin(); cit
             != primarySims.end(); ++cit) {
         //There will be one write to the tree for each pion found.
         if (debug_ > 1)
             std::cout << "\t**Starting particle...**\n";
         const PFSimParticle& sim = sims[*cit];
         //One sim per calib =>
         calib_->sim_numEvent_ = 1;
         calib_->sim_isMC_ = true;
         calib_->sim_energyEvent_ = sim.trajectoryPoint(PFTrajectoryPoint::ClosestApproach).momentum().E();
         calib_->sim_phi_ = sim.trajectoryPoint(PFTrajectoryPoint::ClosestApproach).momentum().Phi();
         calib_->sim_eta_ = sim.trajectoryPoint(PFTrajectoryPoint::ClosestApproach).momentum().Eta();
 
         if (sim.nTrajectoryPoints() > PFTrajectoryPoint::ECALEntrance) {
             calib_->sim_etaEcal_ = sim.trajectoryPoint(PFTrajectoryPoint::ECALEntrance).positionREP().Eta();
             calib_->sim_phiEcal_ = sim.trajectoryPoint(PFTrajectoryPoint::ECALEntrance).positionREP().Phi();
         }
         if (sim.nTrajectoryPoints() > PFTrajectoryPoint::HCALEntrance) {
             calib_->sim_etaHcal_ = sim.trajectoryPoint(PFTrajectoryPoint::HCALEntrance).positionREP().Eta();
             calib_->sim_phiHcal_ = sim.trajectoryPoint(PFTrajectoryPoint::HCALEntrance).positionREP().Phi();
         }
 
         // Find candidates near this sim particle
         std::vector<unsigned> matchingCands = findCandidatesInDeltaR(sim,
                 candidates, deltaRCandToSim_);
         if (debug_ > 3)
             std::cout << "\t\tFound candidates near sim, found "
                     << matchingCands.size()<< " of them.\n";
         if (matchingCands.size() == 0)
             thisParticlePasses_ = false;
         for (std::vector<unsigned>::const_iterator mcIt = matchingCands.begin(); mcIt
                 != matchingCands.end(); ++mcIt) {
             const PFCandidate& theCand = candidates[*mcIt];
             extractCandidate(theCand);
         }
         //Finally,
         fillTreeAndReset();
         
     }
 
     delete pfCandidates_;
     delete simParticles_;
     delete clustersEcal_;
     delete clustersHcal_;
 
     if (thisEventPasses_)
         ++nEventWrites_;
     else
         ++nEventFails_;
 
 }

void CalibratableTest::beginJob ( void )

privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 51 of file CalibratableTest.cc.

References gather_cfg::cout, and debug_.

                                 {
     if (debug_ > 1)
         std::cout << __PRETTY_FUNCTION__ << "\n";
 
 }

double CalibratableTest::deltaR	(	const double &	eta1,
		const double &	eta2,
		const double &	phi1,
		const double &	phi2
	)

private

Definition at line 305 of file CalibratableTest.cc.

References HLTFastRecoForTau_cff::deltaEta, SiPixelRawToDigiRegional_cfi::deltaPhi, M_PI, funct::pow(), and mathSSE::sqrt().

Referenced by findCandidatesInDeltaR().

                                                 {
     double deltaEta = fabs(eta1 - eta2);
     double deltaPhi = fabs(phi1 - phi2);
     if (deltaPhi > M_PI) {
         deltaPhi = 2 * M_PI- deltaPhi;
     }
     return sqrt(pow(deltaEta, 2) + pow(deltaPhi, 2));
 }

void CalibratableTest::endJob ( void )

privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 291 of file CalibratableTest.cc.

References gather_cfg::cout, debug_, nEventFails_, nEventWrites_, nParticleFails_, and nParticleWrites_.

                               {
 
     if (debug_> 0) {
         std::cout << __PRETTY_FUNCTION__ << std::endl;
 
         std::cout << "\tnParticleWrites: "<< nParticleWrites_
                 << ", nParticleFails: "<< nParticleFails_ << "\n";
         std::cout << "\tnEventWrites: "<< nEventWrites_ << ", nEventFails: "
                 << nEventFails_ << "\n";
         std::cout << "Leaving "<< __PRETTY_FUNCTION__ << "\n";
     }
 
 }

void CalibratableTest::extractCandidate ( const reco::PFCandidate & cand )

private

Definition at line 171 of file CalibratableTest.cc.

References createPayload::block, calib_, pftools::Calibratable::cands_, pftools::Calibratable::cluster_ecal_, pftools::Calibratable::cluster_hcal_, gather_cfg::cout, debug_, ECAL, reco::PFBlock::elements(), asciidump::elements, reco::PFCluster::energy(), pftools::CandidateWrapper::energy_, pftools::CandidateWrapper::energyEcal_, pftools::CandidateWrapper::energyHcal_, pftools::CandidateWrapper::eta_, HCAL, reco::PFCluster::layer(), pftools::CandidateWrapper::phi_, reco::PFCluster::positionREP(), and pftools::CandidateWrapper::type_.

Referenced by analyze().

                                                                {
     if (debug_ > 3)
         std::cout << "\tCandidate: "<< cand << "\n";
 
     //There may be several PFCandiates per sim particle. So we create a mini-class
     //to represent each one. CandidateWrapper is defined in Calibratable.
     //It's very easy to use, as we shall see...
     CandidateWrapper cw;
     cw.energy_ = cand.energy();
     cw.eta_ = cand.eta();
     cw.phi_ = cand.phi();
     cw.type_ = cand.particleId();
     cw.energyEcal_ = cand.ecalEnergy();
     cw.energyHcal_ = cand.hcalEnergy();
     if (debug_ > 4)
         std::cout << "\t\tECAL energy = " << cand.ecalEnergy()
                 << ", HCAL energy = " << cand.hcalEnergy() << "\n";
 
     //Now, extract block elements from the pfCandidate:
     PFCandidate::ElementsInBlocks eleInBlocks = cand.elementsInBlocks();
     if (debug_ > 3)
         std::cout << "\tLooping over elements in blocks, "
                 << eleInBlocks.size() << " of them."<< std::endl;
     for (PFCandidate::ElementsInBlocks::iterator bit = eleInBlocks.begin(); bit
             != eleInBlocks.end(); ++bit) {
 
         /* 
          * Find PFClusters associated with this candidate.
          */
         
         //Extract block reference
         PFBlockRef blockRef((*bit).first);
         //Extract index
         unsigned indexInBlock((*bit).second);
         //Dereference the block (what a palava!)
         const PFBlock& block = *blockRef;
         //And finally get a handle on the elements
         const edm::OwnVector<reco::PFBlockElement> & elements = block.elements();
         //get references to the candidate's track, ecal clusters and hcal clusters
         switch (elements[indexInBlock].type()) {
         case PFBlockElement::ECAL: {
             reco::PFClusterRef clusterRef = elements[indexInBlock].clusterRef();
             const PFCluster theRealCluster = *clusterRef;
             CalibratableElement d(theRealCluster.energy(),
                     theRealCluster.positionREP().eta(), theRealCluster.positionREP().phi(), theRealCluster.layer() );
             calib_->cluster_ecal_.push_back(d);
             if (debug_ > 4)
                 std::cout << "\t\tECAL cluster: "<< theRealCluster << "\n";
 
             break;
         }
 
         case PFBlockElement::HCAL: {
             reco::PFClusterRef clusterRef = elements[indexInBlock].clusterRef();
             const PFCluster theRealCluster = *clusterRef;
             CalibratableElement d(theRealCluster.energy(),
                     theRealCluster.positionREP().eta(), theRealCluster.positionREP().phi(), theRealCluster.layer() );
             calib_->cluster_hcal_.push_back(d);
             if (debug_ > 4)
                 std::cout << "\t\tHCAL cluster: "<< theRealCluster << "\n";
 
             break;
         }
 
         default:
             if (debug_ > 3)
                 std::cout << "\t\tOther block type: "
                         << elements[indexInBlock].type() << "\n";
             break;
         }
 
     }
     //For each candidate found,
     calib_->cands_.push_back(cw);
 }

void CalibratableTest::fillTreeAndReset ( )

privatevirtual

Definition at line 274 of file CalibratableTest.cc.

References calib_, gather_cfg::cout, debug_, nParticleFails_, nParticleWrites_, pftools::Calibratable::recompute(), pftools::Calibratable::reset(), thisEventPasses_, thisParticlePasses_, and tree_.

Referenced by analyze().

                                         {
     if (thisEventPasses_ && thisParticlePasses_) {
         ++nParticleWrites_;
         calib_->recompute();
         if (debug_> 4) {
             //print a summary
             std::cout << *calib_;
         }
         tree_->Fill();
     } else {
         ++nParticleFails_;
     }
     if (debug_ > 1)
         std::cout << "\t**Finished particle.**\n";
     calib_->reset();
 }

std::vector< unsigned > CalibratableTest::findCandidatesInDeltaR	(	const reco::PFSimParticle &	pft,
		const std::vector< reco::PFCandidate > &	cands,
		const double &	deltaR
	)

private

Definition at line 247 of file CalibratableTest.cc.

References deltaR(), reco::LeafCandidate::eta(), getHLTprescales::index, reco::LeafCandidate::phi(), reco::PFTrajectoryPoint::positionREP(), and reco::PFTrack::trajectoryPoint().

Referenced by analyze().

                                  {
 
     unsigned index(0);
     std::vector<unsigned> answers;
 
     double trEta = pft.trajectoryPoint(PFTrajectoryPoint::ECALEntrance).positionREP().Eta();
     double trPhi = pft.trajectoryPoint(PFTrajectoryPoint::ECALEntrance).positionREP().Phi();
 
     for (std::vector<PFCandidate>::const_iterator cit = cands.begin(); cit
             != cands.end(); ++cit) {
 
         PFCandidate cand = *cit;
         double cEta = cand.eta();
         double cPhi = cand.phi();
 
         if (deltaR(cEta, trEta, cPhi, trPhi) < deltaRCut) {
             //accept
             answers.push_back(index);
         }
 
         ++index;
     }
     return answers;
 }

std::vector< unsigned > CalibratableTest::findPrimarySimParticles ( const std::vector< reco::PFSimParticle > & sims )

private

Definition at line 149 of file CalibratableTest.cc.

References abs, reco::PFSimParticle::daughterIds(), getHLTprescales::index, reco::PFSimParticle::motherId(), and reco::PFSimParticle::pdgCode().

Referenced by analyze().

                                               {
     std::vector<unsigned> answers;
     unsigned index(0);
     for (std::vector<PFSimParticle>::const_iterator cit = sims.begin(); cit
             != sims.end(); ++cit) {
         PFSimParticle theSim = *cit;
         //TODO: what about rejected events?
         if (theSim.motherId() >= 0)
             continue;
         int particleId = abs(theSim.pdgCode());
         if (particleId != 211)
             continue;
         //TODO: ...particularly interacting pions?
         if (theSim.daughterIds().size() > 0)
             continue;
         answers.push_back(index);
         ++index;
     }
     return answers;
 }

template<class T >

void CalibratableTest::getCollection	(	edm::Handle< T > &	c,
		const edm::InputTag &	tag,
		const edm::Event &	event
	)		const

Definition at line 147 of file CalibratableTest.h.

References gather_cfg::cout, and edm::HandleBase::isValid().

Referenced by analyze().

                                                              {
 
     try {
         event.getByLabel(tag, c);
         if(!c.isValid()) {
             std::cout << "Warning! Collection for label " << tag << " is not valid!" << std::endl;
         }
     }
     catch (cms::Exception& err) {
         std::cout << "Couldn't get collection\n";
         std::ostringstream err;
         //LogError("Error getting collection!") << err;
     }
 }