TreeFiller/src/FillerTracks.cc

// $Id: FillerTracks.cc,v 1.23 2008/11/03 11:22:35 bendavid Exp $

#include "MitProd/TreeFiller/interface/FillerTracks.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "DataFormats/Common/interface/Handle.h"
#include "DataFormats/TrackReco/interface/HitPattern.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
#include "DataFormats/RecoCandidate/interface/TrackAssociation.h"
#include "TrackingTools/TrackAssociator/interface/TrackDetectorAssociator.h"

using namespace std;
using namespace edm;
using namespace mithep;

//--------------------------------------------------------------------------------------------------
FillerTracks::FillerTracks(const ParameterSet &cfg, const char *name, bool active, bool ecalActive,
                           const char *edmName, const char *mitName) :
  BaseFiller(cfg,name,active),
  ecalAssocActive_(Conf().getUntrackedParameter<bool>("ecalAssocActive",ecalActive)),
  edmName_(Conf().getUntrackedParameter<string>("edmName",edmName)),
  mitName_(Conf().getUntrackedParameter<string>("mitName",mitName)),
  edmSimAssocName_(Conf().getUntrackedParameter<string>("edmSimAssociationName","")),
  trackingMapName_(Conf().getUntrackedParameter<string>("trackingMapName","TrackingMap")),
  barrelSuperClusterIdMapName_(Conf().getUntrackedParameter<string>("superClusterIdMapName","barrelSuperClusterIdMap")),
  endcapSuperClusterIdMapName_(Conf().getUntrackedParameter<string>("endcapClusterIdMapName","endcapSuperClusterIdMap")),
  trackMapName_(Conf().getUntrackedParameter<string>("trackMapName",
                                                     Form("%sMapName",mitName_.c_str()))),
  trackingMap_(0),
  tracks_(new mithep::TrackArr(250)), 
  trackMap_(new mithep::TrackMap)
{
  // Constructor.
  
  //initialize track associator configuration if needed
  if (ecalAssocActive_)
    assocParams_.loadParameters(cfg.getUntrackedParameter<ParameterSet>("TrackAssociatorParameters"));

}

//--------------------------------------------------------------------------------------------------
FillerTracks::~FillerTracks()
{
  // Destructor.

  delete tracks_;
  delete trackMap_;
}

//--------------------------------------------------------------------------------------------------
void FillerTracks::BookDataBlock(TreeWriter &tws)
{
  // Add tracks branch to tree, publish and get our objects.

  tws.AddBranch(mitName_.c_str(),&tracks_);

  trackingMap_ = OS()->get<TrackingParticleMap>(trackingMapName_.c_str());
  barrelSuperClusterIdMap_ = OS()->get<SuperClusterIdMap>(barrelSuperClusterIdMapName_.c_str());
  endcapSuperClusterIdMap_ = OS()->get<SuperClusterIdMap>(endcapSuperClusterIdMapName_.c_str());
  OS()->add<TrackMap>(trackMap_,trackMapName_.c_str());
  OS()->add<TrackArr>(tracks_,mitName_.c_str());
}

//--------------------------------------------------------------------------------------------------
void FillerTracks::FillDataBlock(const edm::Event      &event, 
                                 const edm::EventSetup &setup)
{
  // Fill tracks from edm collection into our collection.

  tracks_  ->Reset();
  trackMap_->Reset();

  Handle<View<reco::Track> > hTrackProduct;
  GetProduct(edmName_, hTrackProduct, event);  
        
  trackMap_->SetEdmProductId(hTrackProduct.id().id());
  const View<reco::Track> inTracks = *(hTrackProduct.product());  
  
  // for MC SimParticle association (reco->sim mappings)
  reco::RecoToSimCollection simAssociation;
  if (trackingMap_ && !edmSimAssocName_.empty()) {
    Handle<reco::RecoToSimCollection> simAssociationProduct;
    GetProduct(edmSimAssocName_, simAssociationProduct, event);  
    simAssociation = *(simAssociationProduct.product());
    //printf("SimAssociation Map Size = %i\n",simAssociation.size());
  }
  
  //set up associator for Track-Ecal associations
  TrackDetectorAssociator trackAssociator;
  trackAssociator.useDefaultPropagator();

  // loop through all tracks and fill the information
  for (View<reco::Track>::const_iterator it = inTracks.begin();
         it != inTracks.end(); ++it) {
         
    mithep::Track *outTrack = tracks_->Allocate();
    // create track and set the core parameters
    new (outTrack) mithep::Track(it->qoverp(),it->lambda(),
                                 it->phi(),it->dxy(),it->dsz());
    outTrack->SetErrors(it->qoverpError(),it->lambdaError(),
                        it->phiError(),it->dxyError(),it->dszError());
        
    // fill track quality information
    outTrack->SetChi2(it->chi2());
    outTrack->SetNdof(static_cast<Int_t>(it->ndof()));
    
    //fill hit layer map
    const reco::HitPattern &hits = it->hitPattern();
    for (Int_t i=0; i<hits.numberOfHits(); i++) {
      uint32_t hit = hits.getHitPattern(i);
      if (hits.validHitFilter(hit))
        if (hits.trackerHitFilter(hit))
          outTrack->SetHit(hitReader_.Layer(hit));
                
      if (0) {
        if (hits.muonDTHitFilter(hit))
          printf("Muon DT Layer = %i\n", hits.getLayer(hit));
        if (hits.muonCSCHitFilter(hit))
          printf("Muon CSC Layer = %i\n", hits.getLayer(hit));        
        if (hits.muonRPCHitFilter(hit))
          printf("Muon RPC Layer = %i\n", hits.getLayer(hit));
      }
    }

    //make ecal associations
    if (ecalAssocActive_) {
      TrackDetMatchInfo matchInfo = trackAssociator.associate(event,setup,
                                                             *it,
                                                              assocParams_);
      outTrack->SetEtaEcal(matchInfo.trkGlobPosAtEcal.eta());
      outTrack->SetPhiEcal(matchInfo.trkGlobPosAtEcal.phi());

      //fill supercluster link
      if (barrelSuperClusterIdMap_ || endcapSuperClusterIdMap_) {
        mithep::SuperCluster *cluster = 0;
        for (std::vector<const ::CaloTower*>::const_iterator iTower = matchInfo.crossedTowers.begin();
            iTower<matchInfo.crossedTowers.end() && !cluster; iTower++) {

          for (uint ihit=0; ihit<(*iTower)->constituentsSize() && !cluster; ihit++) {
            DetId hit = (*iTower)->constituent(ihit);
            if (barrelSuperClusterIdMap_ && barrelSuperClusterIdMap_->HasMit(hit))
              cluster = barrelSuperClusterIdMap_->GetMit(hit);
            else if (endcapSuperClusterIdMap_ && endcapSuperClusterIdMap_->HasMit(hit))
              cluster = endcapSuperClusterIdMap_->GetMit(hit);
          }
        }
        if (cluster)
          outTrack->SetSCluster(cluster);

      }
    }
    
    // add reference between mithep and edm object
    edm::Ptr<reco::Track> thePtr = inTracks.ptrAt(it - inTracks.begin());
    trackMap_->Add(thePtr, outTrack);
        
    //do dim associations
    if (trackingMap_ && !edmSimAssocName_.empty()) {
      //printf("Trying Track-Sim association\n");
      reco::TrackBaseRef theBaseRef = inTracks.refAt(it - inTracks.begin());
      vector<pair<TrackingParticleRef, double> > simRefs;
      Bool_t noSimParticle = 0;
      try {
        simRefs = simAssociation[theBaseRef]; //try to get the sim references if existing
      }
      catch (edm::Exception &ex) {
        noSimParticle = 1;
      }

      if (!noSimParticle) { //loop through sim match candidates
        //printf("Applying track-sim association\n");
        for (vector<pair<TrackingParticleRef, double> >::const_iterator simRefPair=simRefs.begin(); 
             simRefPair != simRefs.end(); ++simRefPair) 

          if (simRefPair->second > 0.5) // require more than 50% shared hits between reco and sim
            outTrack->SetMCPart(trackingMap_->GetMit(simRefPair->first)); //add reco->sim reference
      }
    }
  }
  tracks_->Trim();
}

Revision:	1.24
Committed:	Mon Nov 3 18:11:10 2008 UTC (16 years, 6 months ago) by bendavid
Content type:	text/plain
Branch:	MAIN
Changes since 1.23:	+9 -8 lines
Log Message:	Changed to base class usage of reco::Track throughout and removed FillerGsfTracks
#	User	Rev	Content
1	bendavid	1.24	// $Id: FillerTracks.cc,v 1.23 2008/11/03 11:22:35 bendavid Exp $
2	loizides	1.1
3	bendavid	1.11	#include "MitProd/TreeFiller/interface/FillerTracks.h"
4	loizides	1.1	#include "FWCore/MessageLogger/interface/MessageLogger.h"
5			#include "DataFormats/Common/interface/Handle.h"
6	bendavid	1.14	#include "DataFormats/TrackReco/interface/HitPattern.h"
7	loizides	1.1	#include "DataFormats/TrackReco/interface/Track.h"
8			#include "DataFormats/TrackReco/interface/TrackFwd.h"
9			#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
10			#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
11			#include "DataFormats/RecoCandidate/interface/TrackAssociation.h"
12	bendavid	1.23	#include "TrackingTools/TrackAssociator/interface/TrackDetectorAssociator.h"
13	loizides	1.1
14			using namespace std;
15			using namespace edm;
16			using namespace mithep;
17
18	loizides	1.4	//--------------------------------------------------------------------------------------------------
19	bendavid	1.23	FillerTracks::FillerTracks(const ParameterSet &cfg, const char *name, bool active, bool ecalActive,
20	loizides	1.16	const char edmName, const char mitName) :
21	loizides	1.7	BaseFiller(cfg,name,active),
22	bendavid	1.23	ecalAssocActive_(Conf().getUntrackedParameter<bool>("ecalAssocActive",ecalActive)),
23	loizides	1.16	edmName_(Conf().getUntrackedParameter<string>("edmName",edmName)),
24			mitName_(Conf().getUntrackedParameter<string>("mitName",mitName)),
25	loizides	1.19	edmSimAssocName_(Conf().getUntrackedParameter<string>("edmSimAssociationName","")),
26	bendavid	1.20	trackingMapName_(Conf().getUntrackedParameter<string>("trackingMapName","TrackingMap")),
27	bendavid	1.23	barrelSuperClusterIdMapName_(Conf().getUntrackedParameter<string>("superClusterIdMapName","barrelSuperClusterIdMap")),
28			endcapSuperClusterIdMapName_(Conf().getUntrackedParameter<string>("endcapClusterIdMapName","endcapSuperClusterIdMap")),
29	loizides	1.13	trackMapName_(Conf().getUntrackedParameter<string>("trackMapName",
30			Form("%sMapName",mitName_.c_str()))),
31	bendavid	1.17	trackingMap_(0),
32	loizides	1.7	tracks_(new mithep::TrackArr(250)),
33	loizides	1.1	trackMap_(new mithep::TrackMap)
34			{
35	loizides	1.13	// Constructor.
36	bendavid	1.23
37			//initialize track associator configuration if needed
38			if (ecalAssocActive_)
39			assocParams_.loadParameters(cfg.getUntrackedParameter<ParameterSet>("TrackAssociatorParameters"));
40
41	loizides	1.1	}
42
43	loizides	1.4	//--------------------------------------------------------------------------------------------------
44	loizides	1.1	FillerTracks::~FillerTracks()
45			{
46	loizides	1.13	// Destructor.
47
48	loizides	1.5	delete tracks_;
49	loizides	1.1	delete trackMap_;
50			}
51
52	loizides	1.4	//--------------------------------------------------------------------------------------------------
53	loizides	1.1	void FillerTracks::BookDataBlock(TreeWriter &tws)
54			{
55	loizides	1.13	// Add tracks branch to tree, publish and get our objects.
56
57	loizides	1.1	tws.AddBranch(mitName_.c_str(),&tracks_);
58	loizides	1.12
59	bendavid	1.17	trackingMap_ = OS()->get<TrackingParticleMap>(trackingMapName_.c_str());
60	bendavid	1.23	barrelSuperClusterIdMap_ = OS()->get<SuperClusterIdMap>(barrelSuperClusterIdMapName_.c_str());
61			endcapSuperClusterIdMap_ = OS()->get<SuperClusterIdMap>(endcapSuperClusterIdMapName_.c_str());
62	loizides	1.13	OS()->add<TrackMap>(trackMap_,trackMapName_.c_str());
63			OS()->add<TrackArr>(tracks_,mitName_.c_str());
64	loizides	1.1	}
65
66	loizides	1.4	//--------------------------------------------------------------------------------------------------
67	loizides	1.1	void FillerTracks::FillDataBlock(const edm::Event &event,
68			const edm::EventSetup &setup)
69			{
70	loizides	1.13	// Fill tracks from edm collection into our collection.
71
72	paus	1.10	tracks_ ->Reset();
73	loizides	1.1	trackMap_->Reset();
74	loizides	1.6
75	bendavid	1.24	Handle<View<reco::Track> > hTrackProduct;
76	loizides	1.6	GetProduct(edmName_, hTrackProduct, event);
77	loizides	1.1
78	loizides	1.6	trackMap_->SetEdmProductId(hTrackProduct.id().id());
79	bendavid	1.24	const View<reco::Track> inTracks = *(hTrackProduct.product());
80	loizides	1.1
81	paus	1.10	// for MC SimParticle association (reco->sim mappings)
82	loizides	1.1	reco::RecoToSimCollection simAssociation;
83	loizides	1.19	if (trackingMap_ && !edmSimAssocName_.empty()) {
84	loizides	1.1	Handle<reco::RecoToSimCollection> simAssociationProduct;
85	loizides	1.19	GetProduct(edmSimAssocName_, simAssociationProduct, event);
86	loizides	1.1	simAssociation = *(simAssociationProduct.product());
87	bendavid	1.17	//printf("SimAssociation Map Size = %i\n",simAssociation.size());
88	loizides	1.1	}
89
90	bendavid	1.23	//set up associator for Track-Ecal associations
91			TrackDetectorAssociator trackAssociator;
92			trackAssociator.useDefaultPropagator();
93
94	paus	1.10	// loop through all tracks and fill the information
95	bendavid	1.24	for (View<reco::Track>::const_iterator it = inTracks.begin();
96			it != inTracks.end(); ++it) {
97
98	loizides	1.5	mithep::Track *outTrack = tracks_->Allocate();
99	paus	1.10	// create track and set the core parameters
100	loizides	1.19	new (outTrack) mithep::Track(it->qoverp(),it->lambda(),
101			it->phi(),it->dxy(),it->dsz());
102			outTrack->SetErrors(it->qoverpError(),it->lambdaError(),
103			it->phiError(),it->dxyError(),it->dszError());
104	loizides	1.1
105	loizides	1.19	// fill track quality information
106	bendavid	1.14	outTrack->SetChi2(it->chi2());
107	loizides	1.16	outTrack->SetNdof(static_cast<Int_t>(it->ndof()));
108	bendavid	1.14
109	loizides	1.19	//fill hit layer map
110	bendavid	1.14	const reco::HitPattern &hits = it->hitPattern();
111			for (Int_t i=0; i<hits.numberOfHits(); i++) {
112			uint32_t hit = hits.getHitPattern(i);
113	loizides	1.19	if (hits.validHitFilter(hit))
114			if (hits.trackerHitFilter(hit))
115	bendavid	1.22	outTrack->SetHit(hitReader_.Layer(hit));
116	bendavid	1.14
117	loizides	1.19	if (0) {
118			if (hits.muonDTHitFilter(hit))
119			printf("Muon DT Layer = %i\n", hits.getLayer(hit));
120			if (hits.muonCSCHitFilter(hit))
121			printf("Muon CSC Layer = %i\n", hits.getLayer(hit));
122			if (hits.muonRPCHitFilter(hit))
123			printf("Muon RPC Layer = %i\n", hits.getLayer(hit));
124			}
125	bendavid	1.14	}
126	bendavid	1.23
127			//make ecal associations
128			if (ecalAssocActive_) {
129			TrackDetMatchInfo matchInfo = trackAssociator.associate(event,setup,
130			*it,
131			assocParams_);
132			outTrack->SetEtaEcal(matchInfo.trkGlobPosAtEcal.eta());
133			outTrack->SetPhiEcal(matchInfo.trkGlobPosAtEcal.phi());
134
135			//fill supercluster link
136			if (barrelSuperClusterIdMap_ \|\| endcapSuperClusterIdMap_) {
137			mithep::SuperCluster *cluster = 0;
138			for (std::vector<const ::CaloTower*>::const_iterator iTower = matchInfo.crossedTowers.begin();
139			iTower<matchInfo.crossedTowers.end() && !cluster; iTower++) {
140
141			for (uint ihit=0; ihit<(*iTower)->constituentsSize() && !cluster; ihit++) {
142			DetId hit = (*iTower)->constituent(ihit);
143			if (barrelSuperClusterIdMap_ && barrelSuperClusterIdMap_->HasMit(hit))
144			cluster = barrelSuperClusterIdMap_->GetMit(hit);
145			else if (endcapSuperClusterIdMap_ && endcapSuperClusterIdMap_->HasMit(hit))
146			cluster = endcapSuperClusterIdMap_->GetMit(hit);
147			}
148			}
149			if (cluster)
150			outTrack->SetSCluster(cluster);
151
152			}
153			}
154	bendavid	1.14
155	loizides	1.1	// add reference between mithep and edm object
156	bendavid	1.24	edm::Ptr<reco::Track> thePtr = inTracks.ptrAt(it - inTracks.begin());
157			trackMap_->Add(thePtr, outTrack);
158	loizides	1.1
159	bendavid	1.23	//do dim associations
160	loizides	1.19	if (trackingMap_ && !edmSimAssocName_.empty()) {
161	bendavid	1.17	//printf("Trying Track-Sim association\n");
162	bendavid	1.24	reco::TrackBaseRef theBaseRef = inTracks.refAt(it - inTracks.begin());
163	loizides	1.1	vector<pair<TrackingParticleRef, double> > simRefs;
164	paus	1.10	Bool_t noSimParticle = 0;
165	loizides	1.1	try {
166			simRefs = simAssociation[theBaseRef]; //try to get the sim references if existing
167			}
168	loizides	1.5	catch (edm::Exception &ex) {
169	paus	1.10	noSimParticle = 1;
170	loizides	1.1	}
171	paus	1.10
172	loizides	1.1	if (!noSimParticle) { //loop through sim match candidates
173	bendavid	1.17	//printf("Applying track-sim association\n");
174	loizides	1.1	for (vector<pair<TrackingParticleRef, double> >::const_iterator simRefPair=simRefs.begin();
175			simRefPair != simRefs.end(); ++simRefPair)
176
177	paus	1.10	if (simRefPair->second > 0.5) // require more than 50% shared hits between reco and sim
178	bendavid	1.17	outTrack->SetMCPart(trackingMap_->GetMit(simRefPair->first)); //add reco->sim reference
179	loizides	1.1	}
180			}
181			}
182			tracks_->Trim();
183			}
184	bendavid	1.14