Morgan/src/VertexAnalyzer.cc

#include "../interface/VertexAnalyzer.h"

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

VertexAnalyzer::VertexAnalyzer(const edm::ParameterSet& producersNames):verbosity_(0)
{
        primaryVertexProducer_ = producersNames.getParameter<edm::InputTag>("primaryVertexProducer");
        beamSpotProducer_ = producersNames.getParameter<edm::InputTag>("beamSpotProducer");
}


VertexAnalyzer::VertexAnalyzer(const edm::ParameterSet& producersNames, int verbosity):verbosity_(verbosity)
{
        primaryVertexProducer_ = producersNames.getParameter<edm::InputTag>("primaryVertexProducer");
        beamSpotProducer_ = producersNames.getParameter<edm::InputTag>("beamSpotProducer");
}


VertexAnalyzer::~VertexAnalyzer()
{
}


void VertexAnalyzer::getBeamSpot(const edm::Event& iEvent, TRootBeamSpot* rootBeamSpot)
{

                edm::Handle<reco::BeamSpot> recoBeamSpot;
                iEvent.getByLabel(beamSpotProducer_, recoBeamSpot);

                if ( recoBeamSpot.isValid() )
                {
                        if(verbosity_>1) std::cout  << "   Beam Spot - Label: " << beamSpotProducer_.label() << "   Instance: " << beamSpotProducer_.instance() << std::endl;
                        rootBeamSpot->fill(
                                recoBeamSpot->x0()
                                ,recoBeamSpot->y0()
                                ,recoBeamSpot->z0()
                                ,recoBeamSpot->x0Error()
                                ,recoBeamSpot->y0Error()
                                ,recoBeamSpot->z0Error()
                                ,recoBeamSpot->sigmaZ()
                                ,recoBeamSpot->sigmaZ0Error()
                                ,recoBeamSpot->BeamWidth()
                                ,recoBeamSpot->BeamWidthError()
                        );                      
                        if(verbosity_>2) cout << "   "<< *rootBeamSpot << endl << endl;
                }               
                else
                {
                        if(verbosity_>1) std::cout  << "   ***** ERROR in VertexAnalyzer::getBeamSpot - No beam spot available from EventSetup *****" << std::endl;
                }

}


void VertexAnalyzer::getVertices(const edm::Event& iEvent, TClonesArray* rootVertices)
{

        edm::Handle< reco::VertexCollection > recoVertices;
        iEvent.getByLabel(primaryVertexProducer_, recoVertices);
        if(verbosity_>1) std::cout  << "   Number of primary vertices = " << recoVertices->size() << "   Label: " << primaryVertexProducer_.label() << "   Instance: " << primaryVertexProducer_.instance() << std::endl;


        int iRootVertex = 0;
        for (unsigned int j=0; j<recoVertices->size(); j++)
        {
                const reco::Vertex* vertex = & ((*recoVertices)[j]);

                // Put your vertex selection here....
                if (! vertex->isValid() ) continue;
                if ( vertex->isFake() ) continue;
                        
                TRootVertex localVertex(
                        vertex->x()
                        ,vertex->y()
                        ,vertex->z()
                        ,vertex->xError()
                        ,vertex->yError()
                        ,vertex->zError()
                );

                localVertex.setChi2( vertex->chi2() );
                localVertex.setNdof( vertex->ndof() );

                Int_t ntracks = 0;
                Float_t higherPt = 0.;
                Float_t scalarSumPt = 0.;

                for( std::vector< reco::TrackBaseRef >::const_iterator it = vertex->tracks_begin(); it != vertex->tracks_end(); it++)
                {
                        scalarSumPt += (**it).pt();
                        if( (**it).pt()>higherPt ) higherPt=(**it).pt();
                        ntracks++;
                }

                // No refitted tracks embeded in reco::Vertex....
                //cout  << "vertex->refittedTracks().size()=" << vertex->refittedTracks().size() << endl;
                
                localVertex.setNtracks( ntracks );
                localVertex.setHigherTrackPt( higherPt );
                localVertex.setScalarSumPt( scalarSumPt );

                new( (*rootVertices)[iRootVertex] ) TRootVertex(localVertex);
                if(verbosity_>2) cout << "   ["<< setw(3) << iRootVertex << "] " << localVertex << endl;
                iRootVertex++;
        }

}


void VertexAnalyzer::selectPrimary(TRootEvent* rootEvent, TClonesArray* rootVertices)
{

        int ibestVertex = -1;
        float highestScalarSumPt = -1.;
        TRootVertex* localVertex = 0;
        
        for (int ivtx=0; ivtx<rootVertices->GetEntriesFast(); ivtx++)
        {
                localVertex = (TRootVertex*) rootVertices->At(ivtx);
                // Put you primary vertex selection here...
                if ( localVertex->scalarSumPt()>highestScalarSumPt )
                {
                        highestScalarSumPt = localVertex->scalarSumPt();
                        ibestVertex = ivtx;
                }
        }

        if ( ibestVertex >= 0)
        {
                localVertex = (TRootVertex*) rootVertices->At(ibestVertex);
                rootEvent->setPrimaryVertexIndex( ibestVertex );
                rootEvent->setPrimaryVertex( localVertex );
        }
        
        if(verbosity_>1) std::cout  << "   Selected primary vertex: TRootVertex[" << ibestVertex << "]" << std::endl;
}
Revision:	1.6
Committed:	Tue Apr 21 10:42:23 2009 UTC (16 years ago) by lethuill
Content type:	text/plain
Branch:	MAIN
CVS Tags:	JeSuisBeaucoupPlusGrosQunReco_2_2_7_01, RecoPhoton_2_2_7_02, pat_2_2_7_01
Changes since 1.5:	+35 -2 lines
Log Message:	Use beam spot in impact parameter calculation
#	User	Rev	Content
1	lethuill	1.2	#include "../interface/VertexAnalyzer.h"
2	lethuill	1.1
3			using namespace std;
4			using namespace reco;
5			using namespace edm;
6
7			VertexAnalyzer::VertexAnalyzer(const edm::ParameterSet& producersNames):verbosity_(0)
8			{
9			primaryVertexProducer_ = producersNames.getParameter<edm::InputTag>("primaryVertexProducer");
10	lethuill	1.6	beamSpotProducer_ = producersNames.getParameter<edm::InputTag>("beamSpotProducer");
11	lethuill	1.1	}
12
13	lethuill	1.4
14	lethuill	1.1	VertexAnalyzer::VertexAnalyzer(const edm::ParameterSet& producersNames, int verbosity):verbosity_(verbosity)
15			{
16			primaryVertexProducer_ = producersNames.getParameter<edm::InputTag>("primaryVertexProducer");
17	lethuill	1.6	beamSpotProducer_ = producersNames.getParameter<edm::InputTag>("beamSpotProducer");
18	lethuill	1.1	}
19
20	lethuill	1.4
21	lethuill	1.1	VertexAnalyzer::~VertexAnalyzer()
22			{
23			}
24
25	lethuill	1.4
26	lethuill	1.6	void VertexAnalyzer::getBeamSpot(const edm::Event& iEvent, TRootBeamSpot* rootBeamSpot)
27			{
28
29			edm::Handle<reco::BeamSpot> recoBeamSpot;
30			iEvent.getByLabel(beamSpotProducer_, recoBeamSpot);
31
32			if ( recoBeamSpot.isValid() )
33			{
34			if(verbosity_>1) std::cout << " Beam Spot - Label: " << beamSpotProducer_.label() << " Instance: " << beamSpotProducer_.instance() << std::endl;
35			rootBeamSpot->fill(
36			recoBeamSpot->x0()
37			,recoBeamSpot->y0()
38			,recoBeamSpot->z0()
39			,recoBeamSpot->x0Error()
40			,recoBeamSpot->y0Error()
41			,recoBeamSpot->z0Error()
42			,recoBeamSpot->sigmaZ()
43			,recoBeamSpot->sigmaZ0Error()
44			,recoBeamSpot->BeamWidth()
45			,recoBeamSpot->BeamWidthError()
46			);
47			if(verbosity_>2) cout << " "<< *rootBeamSpot << endl << endl;
48			}
49			else
50			{
51			if(verbosity_>1) std::cout << " *** ERROR in VertexAnalyzer::getBeamSpot - No beam spot available from EventSetup ***" << std::endl;
52			}
53
54			}
55
56
57			void VertexAnalyzer::getVertices(const edm::Event& iEvent, TClonesArray* rootVertices)
58	lethuill	1.1	{
59
60	lethuill	1.3	edm::Handle< reco::VertexCollection > recoVertices;
61			iEvent.getByLabel(primaryVertexProducer_, recoVertices);
62			if(verbosity_>1) std::cout << " Number of primary vertices = " << recoVertices->size() << " Label: " << primaryVertexProducer_.label() << " Instance: " << primaryVertexProducer_.instance() << std::endl;
63	lethuill	1.4
64
65			int iRootVertex = 0;
66	lethuill	1.3	for (unsigned int j=0; j<recoVertices->size(); j++)
67	lethuill	1.1	{
68	lethuill	1.3	const reco::Vertex* vertex = & ((*recoVertices)[j]);
69
70			// Put your vertex selection here....
71			if (! vertex->isValid() ) continue;
72			if ( vertex->isFake() ) continue;
73
74			TRootVertex localVertex(
75			vertex->x()
76			,vertex->y()
77			,vertex->z()
78			,vertex->xError()
79			,vertex->yError()
80			,vertex->zError()
81			);
82
83			localVertex.setChi2( vertex->chi2() );
84			localVertex.setNdof( vertex->ndof() );
85
86			Int_t ntracks = 0;
87			Float_t higherPt = 0.;
88	lethuill	1.5	Float_t scalarSumPt = 0.;
89	lethuill	1.3
90			for( std::vector< reco::TrackBaseRef >::const_iterator it = vertex->tracks_begin(); it != vertex->tracks_end(); it++)
91			{
92	lethuill	1.5	scalarSumPt += (**it).pt();
93			if( (it).pt()>higherPt ) higherPt=(it).pt();
94	lethuill	1.3	ntracks++;
95			}
96	lethuill	1.5
97			// No refitted tracks embeded in reco::Vertex....
98			//cout << "vertex->refittedTracks().size()=" << vertex->refittedTracks().size() << endl;
99	lethuill	1.3
100			localVertex.setNtracks( ntracks );
101			localVertex.setHigherTrackPt( higherPt );
102	lethuill	1.5	localVertex.setScalarSumPt( scalarSumPt );
103	lethuill	1.3
104	lethuill	1.4	new( (*rootVertices)[iRootVertex] ) TRootVertex(localVertex);
105			if(verbosity_>2) cout << " ["<< setw(3) << iRootVertex << "] " << localVertex << endl;
106			iRootVertex++;
107	lethuill	1.1	}
108	lethuill	1.4
109			}
110
111
112	lethuill	1.6	void VertexAnalyzer::selectPrimary(TRootEvent* rootEvent, TClonesArray* rootVertices)
113	lethuill	1.4	{
114
115			int ibestVertex = -1;
116	lethuill	1.5	float highestScalarSumPt = -1.;
117	lethuill	1.4	TRootVertex* localVertex = 0;
118	lethuill	1.3
119	lethuill	1.4	for (int ivtx=0; ivtx<rootVertices->GetEntriesFast(); ivtx++)
120			{
121			localVertex = (TRootVertex*) rootVertices->At(ivtx);
122			// Put you primary vertex selection here...
123	lethuill	1.5	if ( localVertex->scalarSumPt()>highestScalarSumPt )
124	lethuill	1.4	{
125	lethuill	1.5	highestScalarSumPt = localVertex->scalarSumPt();
126	lethuill	1.4	ibestVertex = ivtx;
127			}
128			}
129	lethuill	1.1
130	lethuill	1.4	if ( ibestVertex >= 0)
131			{
132			localVertex = (TRootVertex*) rootVertices->At(ibestVertex);
133			rootEvent->setPrimaryVertexIndex( ibestVertex );
134			rootEvent->setPrimaryVertex( localVertex );
135			}
136
137			if(verbosity_>1) std::cout << " Selected primary vertex: TRootVertex[" << ibestVertex << "]" << std::endl;
138	lethuill	1.1	}