Morgan/src/MuonAnalyzer.cc

#include "../interface/MuonAnalyzer.h"

using namespace std;
using namespace reco;
using namespace edm;

MuonAnalyzer::MuonAnalyzer(const edm::ParameterSet& producersNames):LeptonAnalyzer(producersNames),useMC_(false),doPrimaryVertex_(false),doBeamSpot_(false)
{
        //dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
        muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
}

MuonAnalyzer::MuonAnalyzer(const edm::ParameterSet& producersNames, const edm::ParameterSet& myConfig, int verbosity):LeptonAnalyzer(producersNames, verbosity)
{
        //dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
        muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
        useMC_ = myConfig.getUntrackedParameter<bool>("doMuonMC");
        doPrimaryVertex_ = myConfig.getUntrackedParameter<bool>("doPrimaryVertex");
        doBeamSpot_ = myConfig.getUntrackedParameter<bool>("doBeamSpot");
}

MuonAnalyzer::~MuonAnalyzer()
{
}

void MuonAnalyzer::Process(const edm::Event& iEvent, TRootBeamSpot* rootBeamSpot, TClonesArray* rootMuons)
{

        Float_t sintheta = 0.;
        unsigned int nMuons=0;

        edm::Handle < std::vector <reco::Muon> > recoMuons;
        if( dataType_=="RECO" || dataType_=="AOD" )
        {
                iEvent.getByLabel(muonProducer_, recoMuons);
                nMuons = recoMuons->size();
        }

        edm::Handle < std::vector <pat::Muon> > patMuons;
        if( dataType_=="PAT" || dataType_=="PATAOD" )
        {
                iEvent.getByLabel(muonProducer_, patMuons);
                nMuons = patMuons->size();
        }

        if(verbosity_>1) std::cout << "   Number of muons = " << nMuons << "   Label: " << muonProducer_.label() << "   Instance: " << muonProducer_.instance() << std::endl;

        for (unsigned int j=0; j<nMuons; j++)
        {
                const reco::Muon* muon = 0;     
                if( dataType_=="RECO" || dataType_=="AOD" ) muon =  &((*recoMuons)[j]);
                if( dataType_=="PAT" || dataType_=="PATAOD" ) muon = (const reco::Muon*) ( & ((*patMuons)[j]) );

                TRootMuon localMuon(
                        muon->px()
                        ,muon->py()
                        ,muon->pz()
                        ,muon->energy()
                        ,muon->vx()
                        ,muon->vy()
                        ,muon->vz()
                        ,muon->pdgId()
                        ,muon->charge()
                ); 

                // Variables from reco::Muon
                sintheta = sin( localMuon.Theta() );

                localMuon.setCaloEnergy(
                        muon->calEnergy().em * sintheta
                        ,muon->calEnergy().emS9 * sintheta
                        ,muon->calEnergy().had * sintheta
                        ,muon->calEnergy().hadS9 * sintheta
                        ,muon->calEnergy().ho * sintheta
                        ,muon->calEnergy().hoS9 * sintheta
                        ,muon->caloCompatibility()
                );

                localMuon.setIsoR03(
                        muon->isolationR03().emEt
                        ,muon->isolationR03().hadEt
                        ,muon->isolationR03().hoEt
                        ,muon->isolationR03().sumPt
                        ,muon->isolationR03().nTracks
                        ,muon->isolationR03().nJets
                );

                localMuon.setIsoR05(
                        muon->isolationR05().emEt
                        ,muon->isolationR05().hadEt
                        ,muon->isolationR05().hoEt
                        ,muon->isolationR05().sumPt
                        ,muon->isolationR05().nTracks
                        ,muon->isolationR05().nJets
                );

                localMuon.setValidity(
                        muon->isEnergyValid()
                        ,muon->isMatchesValid()
                        ,muon->isIsolationValid()
                );

                localMuon.setDirection( muon->time().direction() );
                localMuon.setAlgo( muon->type() );
                localMuon.setID(
                        int(muon->isGood(reco::Muon::TrackerMuonArbitrated))
                        ,int(muon->isGood(reco::Muon::AllArbitrated))
                        ,int(muon->isGood(reco::Muon::GlobalMuonPromptTight))
                        ,int(muon->isGood(reco::Muon::TMLastStationLoose))
                        ,int(muon->isGood(reco::Muon::TMLastStationTight))
                        ,int(muon->isGood(reco::Muon::TM2DCompatibilityLoose))
                        ,int(muon->isGood(reco::Muon::TM2DCompatibilityTight))
                );

                // Variables from reco::Track
                reco::TrackRef track = muon->innerTrack();
                if ( track.isNonnull() )
                {
                        const reco::HitPattern& hit = track->hitPattern();
                        localMuon.setNumberOfValidPixelHits(hit.numberOfValidPixelHits());
                        localMuon.setNumberOfValidTrackerHits(hit.numberOfValidTrackerHits());
                        localMuon.setPixelLayersWithMeasurement(hit.pixelLayersWithMeasurement());
                        localMuon.setStripLayersWithMeasurement(hit.stripLayersWithMeasurement());

                        if (doBeamSpot_)
                        {
                                if ( rootBeamSpot->sigmaZ() > 0. ) // Check beam spot was correctly initialised
                                {
                                        const reco::TrackBase::Point point( rootBeamSpot->x(), rootBeamSpot->y(), rootBeamSpot->z() );
                                        localMuon.setD0( -1.*(track->dxy(point)) );
                                        localMuon.setDsz( track->dsz(point) );
                                }
                                else
                                {
                                        localMuon.setD0( -1.*(track->dxy()) );
                                        localMuon.setDsz( track->dsz() );
                                }
                        }
                        else
                        {
                                localMuon.setD0( -1.*(track->dxy()) );
                                localMuon.setDsz( track->dsz() );
                        }
                        
                        localMuon.setD0Error(track->d0Error());
                        localMuon.setDszError(track->dszError());
        
                        localMuon.setNormalizedChi2(track->normalizedChi2());
                        localMuon.setPtError(track->ptError());
                        localMuon.setEtaError(track->etaError());
                        localMuon.setPhiError(track->phiError());
                        
                        if(doPrimaryVertex_)
                        {
                                // FIXME - Should be Tracker track... what if STA muon ?
                                const reco::TransientTrack transtrack = trackBuilder_->build( track ) ;
                                float sig3d = IP3DSignificance(transtrack);
                                localMuon.setIP3DSignificance(sig3d);
                        }
                }


                if( dataType_=="RECO" || dataType_=="AOD" )
                {
                        // Some specific methods requiring  RECO / AOD format
                        // Do association to genParticle ?
                        // Add InnerTrack, OuterTrack, GlobalTrack infos ?
                }

                if( dataType_=="PATAOD" || dataType_=="PAT" )
                {
                        // Some specific methods to pat::Muon
                        const pat::Muon *patMuon = dynamic_cast<const pat::Muon*>(&*muon);
                        // Use existing reference to genParticle [ pat::PATObject::genParticleRef() ] ?
                        // Alternative methode for isolation (isoDeposit) ?
                        //
                        // leptonID apparently not initialised in PAT...
                        // cout << "Valeur pourrie du leptonID=" << patMuon->leptonID() << endl;

                        if(useMC_)
                        {
                                // MC truth associator index
                                if ((patMuon->genParticleRef()).isNonnull()) {
                                        localMuon.setGenParticleIndex((patMuon->genParticleRef()).index());
                                } else {
                                        localMuon.setGenParticleIndex(-1);
                                }
                        }
                }

                new( (*rootMuons)[j] ) TRootMuon(localMuon);
                if(verbosity_>2) cout << "   ["<< setw(3) << j << "] " << localMuon << endl;
        }
}
Revision:	1.11
Committed:	Tue Apr 21 10:42:23 2009 UTC (16 years ago) by lethuill
Content type:	text/plain
Branch:	MAIN
CVS Tags:	pat_2_2_7_01
Changes since 1.10:	+28 -6 lines
Log Message:	Use beam spot in impact parameter calculation
#	User	Rev	Content
1	lethuill	1.2	#include "../interface/MuonAnalyzer.h"
2	lethuill	1.1
3			using namespace std;
4			using namespace reco;
5			using namespace edm;
6
7	lethuill	1.11	MuonAnalyzer::MuonAnalyzer(const edm::ParameterSet& producersNames):LeptonAnalyzer(producersNames),useMC_(false),doPrimaryVertex_(false),doBeamSpot_(false)
8	lethuill	1.1	{
9	lethuill	1.10	//dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
10	lethuill	1.1	muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
11			}
12
13	lethuill	1.10	MuonAnalyzer::MuonAnalyzer(const edm::ParameterSet& producersNames, const edm::ParameterSet& myConfig, int verbosity):LeptonAnalyzer(producersNames, verbosity)
14	lethuill	1.1	{
15	lethuill	1.10	//dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
16	lethuill	1.1	muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
17	lethuill	1.4	useMC_ = myConfig.getUntrackedParameter<bool>("doMuonMC");
18	lethuill	1.10	doPrimaryVertex_ = myConfig.getUntrackedParameter<bool>("doPrimaryVertex");
19	lethuill	1.11	doBeamSpot_ = myConfig.getUntrackedParameter<bool>("doBeamSpot");
20	lethuill	1.1	}
21
22			MuonAnalyzer::~MuonAnalyzer()
23			{
24			}
25
26	lethuill	1.11	void MuonAnalyzer::Process(const edm::Event& iEvent, TRootBeamSpot* rootBeamSpot, TClonesArray* rootMuons)
27	lethuill	1.1	{
28
29			Float_t sintheta = 0.;
30	lethuill	1.2	unsigned int nMuons=0;
31
32			edm::Handle < std::vector <reco::Muon> > recoMuons;
33			if( dataType_=="RECO" \|\| dataType_=="AOD" )
34			{
35			iEvent.getByLabel(muonProducer_, recoMuons);
36			nMuons = recoMuons->size();
37			}
38
39			edm::Handle < std::vector <pat::Muon> > patMuons;
40			if( dataType_=="PAT" \|\| dataType_=="PATAOD" )
41			{
42			iEvent.getByLabel(muonProducer_, patMuons);
43			nMuons = patMuons->size();
44			}
45
46			if(verbosity_>1) std::cout << " Number of muons = " << nMuons << " Label: " << muonProducer_.label() << " Instance: " << muonProducer_.instance() << std::endl;
47
48			for (unsigned int j=0; j<nMuons; j++)
49	lethuill	1.1	{
50	lethuill	1.3	const reco::Muon* muon = 0;
51			if( dataType_=="RECO" \|\| dataType_=="AOD" ) muon = &((*recoMuons)[j]);
52			if( dataType_=="PAT" \|\| dataType_=="PATAOD" ) muon = (const reco::Muon) ( & ((patMuons)[j]) );
53	lethuill	1.2
54	lethuill	1.1	TRootMuon localMuon(
55	lethuill	1.3	muon->px()
56			,muon->py()
57			,muon->pz()
58			,muon->energy()
59			,muon->vx()
60			,muon->vy()
61			,muon->vz()
62	lethuill	1.6	,muon->pdgId()
63	lethuill	1.3	,muon->charge()
64	lethuill	1.1	);
65
66	lethuill	1.4	// Variables from reco::Muon
67	lethuill	1.1	sintheta = sin( localMuon.Theta() );
68	lethuill	1.2
69	lethuill	1.3	localMuon.setCaloEnergy(
70			muon->calEnergy().em * sintheta
71			,muon->calEnergy().emS9 * sintheta
72			,muon->calEnergy().had * sintheta
73			,muon->calEnergy().hadS9 * sintheta
74			,muon->calEnergy().ho * sintheta
75			,muon->calEnergy().hoS9 * sintheta
76			,muon->caloCompatibility()
77			);
78
79			localMuon.setIsoR03(
80			muon->isolationR03().emEt
81			,muon->isolationR03().hadEt
82			,muon->isolationR03().hoEt
83			,muon->isolationR03().sumPt
84			,muon->isolationR03().nTracks
85			,muon->isolationR03().nJets
86			);
87
88			localMuon.setIsoR05(
89			muon->isolationR05().emEt
90			,muon->isolationR05().hadEt
91			,muon->isolationR05().hoEt
92			,muon->isolationR05().sumPt
93			,muon->isolationR05().nTracks
94			,muon->isolationR05().nJets
95			);
96
97			localMuon.setValidity(
98			muon->isEnergyValid()
99			,muon->isMatchesValid()
100			,muon->isIsolationValid()
101			);
102
103			localMuon.setDirection( muon->time().direction() );
104			localMuon.setAlgo( muon->type() );
105			localMuon.setID(
106			int(muon->isGood(reco::Muon::TrackerMuonArbitrated))
107			,int(muon->isGood(reco::Muon::AllArbitrated))
108			,int(muon->isGood(reco::Muon::GlobalMuonPromptTight))
109			,int(muon->isGood(reco::Muon::TMLastStationLoose))
110			,int(muon->isGood(reco::Muon::TMLastStationTight))
111			,int(muon->isGood(reco::Muon::TM2DCompatibilityLoose))
112			,int(muon->isGood(reco::Muon::TM2DCompatibilityTight))
113			);
114
115	lethuill	1.7	// Variables from reco::Track
116			reco::TrackRef track = muon->innerTrack();
117			if ( track.isNonnull() )
118			{
119			const reco::HitPattern& hit = track->hitPattern();
120	lethuill	1.11	localMuon.setNumberOfValidPixelHits(hit.numberOfValidPixelHits());
121			localMuon.setNumberOfValidTrackerHits(hit.numberOfValidTrackerHits());
122	lethuill	1.8	localMuon.setPixelLayersWithMeasurement(hit.pixelLayersWithMeasurement());
123			localMuon.setStripLayersWithMeasurement(hit.stripLayersWithMeasurement());
124	lethuill	1.11
125			if (doBeamSpot_)
126			{
127			if ( rootBeamSpot->sigmaZ() > 0. ) // Check beam spot was correctly initialised
128			{
129			const reco::TrackBase::Point point( rootBeamSpot->x(), rootBeamSpot->y(), rootBeamSpot->z() );
130			localMuon.setD0( -1.*(track->dxy(point)) );
131			localMuon.setDsz( track->dsz(point) );
132			}
133			else
134			{
135			localMuon.setD0( -1.*(track->dxy()) );
136			localMuon.setDsz( track->dsz() );
137			}
138			}
139			else
140			{
141			localMuon.setD0( -1.*(track->dxy()) );
142			localMuon.setDsz( track->dsz() );
143			}
144
145	lethuill	1.7	localMuon.setD0Error(track->d0Error());
146			localMuon.setDszError(track->dszError());
147	lethuill	1.11
148	lethuill	1.7	localMuon.setNormalizedChi2(track->normalizedChi2());
149			localMuon.setPtError(track->ptError());
150			localMuon.setEtaError(track->etaError());
151			localMuon.setPhiError(track->phiError());
152	lethuill	1.10
153			if(doPrimaryVertex_)
154			{
155	lethuill	1.11	// FIXME - Should be Tracker track... what if STA muon ?
156	lethuill	1.10	const reco::TransientTrack transtrack = trackBuilder_->build( track ) ;
157			float sig3d = IP3DSignificance(transtrack);
158			localMuon.setIP3DSignificance(sig3d);
159			}
160	lethuill	1.7	}
161
162
163	lethuill	1.3	if( dataType_=="RECO" \|\| dataType_=="AOD" )
164	lethuill	1.2	{
165	lethuill	1.3	// Some specific methods requiring RECO / AOD format
166			// Do association to genParticle ?
167			// Add InnerTrack, OuterTrack, GlobalTrack infos ?
168	lethuill	1.2	}
169
170	lethuill	1.3	if( dataType_=="PATAOD" \|\| dataType_=="PAT" )
171	lethuill	1.2	{
172			// Some specific methods to pat::Muon
173			const pat::Muon patMuon = dynamic_cast<const pat::Muon>(&*muon);
174	lethuill	1.3	// Use existing reference to genParticle [ pat::PATObject::genParticleRef() ] ?
175			// Alternative methode for isolation (isoDeposit) ?
176			//
177			// leptonID apparently not initialised in PAT...
178			// cout << "Valeur pourrie du leptonID=" << patMuon->leptonID() << endl;
179	lethuill	1.4
180			if(useMC_)
181			{
182			// MC truth associator index
183			if ((patMuon->genParticleRef()).isNonnull()) {
184	lethuill	1.5	localMuon.setGenParticleIndex((patMuon->genParticleRef()).index());
185	lethuill	1.4	} else {
186	lethuill	1.5	localMuon.setGenParticleIndex(-1);
187	lethuill	1.4	}
188			}
189	lethuill	1.2	}
190	lethuill	1.6
191	lethuill	1.1	new( (*rootMuons)[j] ) TRootMuon(localMuon);
192			if(verbosity_>2) cout << " ["<< setw(3) << j << "] " << localMuon << endl;
193			}
194			}