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)
{
        muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
}

MuonAnalyzer::MuonAnalyzer(const edm::ParameterSet& producersNames, const edm::ParameterSet& myConfig, int verbosity):LeptonAnalyzer(producersNames, myConfig, verbosity)
{
        useMC_ = myConfig.getUntrackedParameter<bool>("doMuonMC");
        muonProducer_ = producersNames.getParameter<edm::InputTag>("muonProducer");
}

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