ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitProd/TreeFiller/src/FillerDecayParts.cc
Revision: 1.18
Committed: Fri Sep 25 08:42:50 2009 UTC (15 years, 7 months ago) by loizides
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_012f, Mit_012e, Mit_012d, Mit_012c, Mit_012b, Mit_012a, Mit_012, Mit_011a
Changes since 1.17: +2 -2 lines
Log Message:
Extended interface of BookDataBlock to contain event setup.

File Contents

# User Rev Content
1 loizides 1.18 // $Id: FillerDecayParts.cc,v 1.17 2009/06/15 15:00:25 loizides Exp $
2 paus 1.1
3 loizides 1.17 #include "MitProd/TreeFiller/interface/FillerDecayParts.h"
4 paus 1.1 #include "DataFormats/TrackReco/interface/Track.h"
5     #include "DataFormats/TrackReco/interface/TrackFwd.h"
6 bendavid 1.10 #include "DataFormats/TrackReco/interface/HitPattern.h"
7     #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
8 loizides 1.17 #include "DataFormats/RecoCandidate/interface/TrackAssociation.h"
9 paus 1.1 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
10     #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
11 loizides 1.17 #include "MitAna/DataTree/interface/StableDataCol.h"
12     #include "MitAna/DataTree/interface/DecayDataCol.h"
13     #include "MitAna/DataTree/interface/DecayParticleCol.h"
14     #include "MitAna/DataTree/interface/Names.h"
15 loizides 1.8 #include "MitEdm/DataFormats/interface/Collections.h"
16 loizides 1.17 #include "MitProd/ObjectService/interface/ObjectService.h"
17 paus 1.1
18     using namespace std;
19     using namespace edm;
20     using namespace mithep;
21    
22     //--------------------------------------------------------------------------------------------------
23 loizides 1.4 FillerDecayParts::FillerDecayParts(const ParameterSet &cfg, const char *name, bool active) :
24 paus 1.1 BaseFiller(cfg,name,active),
25 loizides 1.4 edmName_(Conf().getUntrackedParameter<string>("edmName","")),
26     mitName_(Conf().getUntrackedParameter<string>("mitName","")),
27 loizides 1.16 stableDataName_(mitName_ + "_StableDatas"),
28     decayDataName_(mitName_ + "_DecayDatas"),
29 bendavid 1.9 vertexMapName_(Conf().getUntrackedParameter<string>("vertexMap","")),
30 loizides 1.8 basePartMapNames_(Conf().exists("basePartMaps") ?
31     Conf().getUntrackedParameter<vector<string> >("basePartMaps") :
32     vector<string>()),
33 bendavid 1.9 vertexMap_(0),
34     decays_(new mithep::DecayParticleArr(250)),
35     stableData_(new mithep::StableDataArr(500)),
36     decayData_(new mithep::DecayDataArr(500))
37 paus 1.1 {
38 loizides 1.4 // Constructor.
39 paus 1.1 }
40    
41     //--------------------------------------------------------------------------------------------------
42     FillerDecayParts::~FillerDecayParts()
43     {
44 loizides 1.4 // Destructor.
45    
46 paus 1.1 delete decays_;
47 bendavid 1.9 delete stableData_;
48     delete decayData_;
49 paus 1.1 }
50    
51     //--------------------------------------------------------------------------------------------------
52 loizides 1.18 void FillerDecayParts::BookDataBlock(TreeWriter &tws, const edm::EventSetup &es)
53 paus 1.1 {
54 loizides 1.12 // Add our branches to tree and get our map.
55 bendavid 1.9
56 loizides 1.16 tws.AddBranch(mitName_, &decays_);
57     OS()->add<mithep::DecayParticleArr>(decays_,mitName_);
58     tws.AddBranch(stableDataName_, &stableData_);
59     OS()->add<mithep::StableDataArr>(stableData_,stableDataName_);
60     tws.AddBranch(decayDataName_, &decayData_);
61     OS()->add<mithep::DecayDataArr>(decayData_,decayDataName_);
62    
63     if (!vertexMapName_.empty()) {
64     vertexMap_ = OS()->get<VertexMap>(vertexMapName_);
65     if (vertexMap_)
66     AddBranchDep(mitName_,vertexMap_->GetBrName());
67     }
68 loizides 1.4
69 bendavid 1.6 for (std::vector<std::string>::const_iterator bmapName = basePartMapNames_.begin();
70     bmapName!=basePartMapNames_.end(); ++bmapName) {
71 loizides 1.16 if (!bmapName->empty()) {
72     const BasePartMap *map = OS()->get<BasePartMap>(*bmapName);
73     if (map) {
74     basePartMaps_.push_back(map);
75     AddBranchDep(stableDataName_,map->GetBrName());
76     AddBranchDep(decayDataName_,map->GetBrName());
77     }
78     }
79 bendavid 1.6 }
80 paus 1.1 }
81    
82     //--------------------------------------------------------------------------------------------------
83     void FillerDecayParts::FillDataBlock(const edm::Event &evt,
84     const edm::EventSetup &setup)
85     {
86 loizides 1.16 // Fill the EDM DecayPart collection into the MIT DecayParticle collection.
87 loizides 1.12
88 bendavid 1.15 decays_->Delete();
89     stableData_->Delete();
90     decayData_->Delete();
91 loizides 1.4
92 bendavid 1.2 Handle<mitedm::DecayPartCol> hParts;
93 paus 1.1 GetProduct(edmName_, hParts, evt);
94 bendavid 1.2 const mitedm::DecayPartCol *iParts = hParts.product();
95 paus 1.1
96 loizides 1.16 // loop through all DecayParts and fill the information
97 bendavid 1.2 for (UInt_t i=0; i<iParts->size(); ++i) {
98 loizides 1.16 const mitedm::DecayPart &p = iParts->at(i); // for convenience
99 bendavid 1.2 mithep::DecayParticle *d = decays_->Allocate();
100 bendavid 1.9 new (d) mithep::DecayParticle(p.pid());
101 bendavid 1.2
102 bendavid 1.9 d->SetMassError(p.fittedMassError());
103     d->SetLxy(p.lxyToPv());
104     d->SetLxyError(p.lxyToPvError());
105     d->SetDxy(p.dxyToPv());
106     d->SetDxyError(p.dxyToPvError());
107     d->SetLz(p.lzToPv());
108     d->SetLzError(p.lzToPvError());
109 bendavid 1.2 d->SetMom(p.fourMomentum());
110 mrudolph 1.5 d->SetChi2(p.chi2());
111     d->SetNdof(p.ndof());
112 bendavid 1.10
113     if (p.primaryVertex().isNonnull()) {
114 bendavid 1.9 d->SetPriVertex(vertexMap_->GetMit(p.primaryVertex()));
115 bendavid 1.10 }
116    
117 bendavid 1.6 if (basePartMaps_.size()) {
118 loizides 1.16 // loop through and add stable daughters
119 bendavid 1.9 for (Int_t j=0; j<p.nStableChild(); ++j) {
120     const mitedm::StableData &stable = p.getStableData(j);
121     mitedm::BasePartPtr thePtr = stable.originalPtr();
122     mithep::Particle *daughter = getMitParticle(thePtr);
123    
124     mithep::StableData *outStable = stableData_->Allocate();
125 loizides 1.16 new (outStable) mithep::StableData(daughter,
126     stable.p3().x(),stable.p3().y(),stable.p3().z());
127 bendavid 1.6
128 loizides 1.16 // fill bad layer mask information using corrected HitPattern
129 bendavid 1.11 const ChargedParticle *cDaughter = dynamic_cast<const ChargedParticle*>(daughter);
130 bendavid 1.14 if (stable.HitsFilled() && cDaughter) {
131 bendavid 1.11 const reco::HitPattern &hitPattern = stable.Hits();
132 bendavid 1.13 BitMask48 crossedLayers;
133 bendavid 1.11 UInt_t numCrossed=0;
134 loizides 1.16 // search for all good crossed layers (with or without hit)
135 bendavid 1.11 for (Int_t hi=0; hi<hitPattern.numberOfHits(); hi++) {
136     uint32_t hit = hitPattern.getHitPattern(hi);
137 loizides 1.16 if (hitPattern.getHitType(hit)<=1) {
138 bendavid 1.11 if (hitPattern.trackerHitFilter(hit)) {
139     numCrossed++;
140     crossedLayers.SetBit(hitReader_.Layer(hit));
141     }
142 loizides 1.16 }
143 bendavid 1.11 }
144 loizides 1.16
145 bendavid 1.13 BitMask48 badLayers = crossedLayers ^ cDaughter->Trk()->Hits();
146 bendavid 1.11 outStable->SetBadLayers(badLayers);
147    
148 loizides 1.16 if(verbose_>1) {
149 loizides 1.12 if (outStable->NWrongHits()) {
150     printf("numCrossed: %i\n",numCrossed);
151     printf("Hit Pattern Size: %i\n",hitPattern.numberOfHits());
152     printf("# of hits: %i\n",cDaughter->Trk()->NHits());
153     printf("# of crossed layers: %i\n",crossedLayers.NBitsSet());
154     printf("# of wrong/missing layers: %i\n",badLayers.NBitsSet());
155     printf("# of missed hits: %i\n",outStable->NMissedHits());
156     printf("# of wrong hits: %i\n",outStable->NWrongHits());
157     printf("Hits : ");
158     for (Int_t biti=63; biti >= 0; --biti) {
159     printf("%i",cDaughter->Trk()->Hits().TestBit(biti));
160     }
161     printf("\n");
162     printf("Crossed Layers : ");
163     for (Int_t biti=63; biti >= 0; --biti) {
164     printf("%i",crossedLayers.TestBit(biti));
165     }
166     printf("\n");
167     printf("WrongOrMissingHits: ");
168     for (Int_t biti=63; biti>=0; --biti) {
169     printf("%i",badLayers.TestBit(biti));
170     }
171     printf("\n");
172     }
173     }
174 bendavid 1.10 }
175 loizides 1.12
176 bendavid 1.9 d->AddDaughterData(outStable);
177     }
178 loizides 1.12
179 loizides 1.16 // loop through and add decay daughters
180 bendavid 1.9 for (Int_t j=0; j<p.nDecayChild(); ++j) {
181     const mitedm::DecayData &decay = p.getDecayData(j);
182     mitedm::BasePartPtr thePtr = decay.originalPtr();
183     mithep::Particle *daughter = getMitParticle(thePtr);
184    
185     mithep::DecayData *outDecay = decayData_->Allocate();
186     new (outDecay) mithep::DecayData(daughter,decay.p4());
187     outDecay->SetMassError(decay.massErr());
188     outDecay->SetLxy(decay.lxy());
189     outDecay->SetLxyError(decay.lxyErr());
190     outDecay->SetLz(decay.lz());
191     outDecay->SetLzError(decay.lzErr());
192     outDecay->SetDxy(decay.dxy());
193     outDecay->SetDxyError(decay.dxyErr());
194    
195     d->AddDaughterData(outDecay);
196 bendavid 1.2 }
197     }
198 paus 1.1 }
199     decays_->Trim();
200 bendavid 1.9 stableData_->Trim();
201     decayData_->Trim();
202     }
203    
204     //--------------------------------------------------------------------------------------------------
205     mithep::Particle *FillerDecayParts::getMitParticle(mitedm::BasePartPtr ptr) const
206     {
207 loizides 1.12 // Return our particle referenced by the edm pointer.
208    
209 bendavid 1.9 mithep::Particle *mitPart = 0;
210     for (std::vector<const mithep::BasePartMap*>::const_iterator bmap = basePartMaps_.begin();
211     bmap!=basePartMaps_.end(); ++bmap) {
212     const mithep::BasePartMap *theMap = *bmap;
213 loizides 1.16 if (theMap->HasMit(ptr)) {
214 bendavid 1.9 mitPart = theMap->GetMit(ptr);
215 loizides 1.16 return mitPart;
216     }
217 bendavid 1.9 }
218    
219     if (!mitPart)
220     throw edm::Exception(edm::errors::Configuration, "FillerDecayParts::FillDataBlock()\n")
221     << "Error! MITHEP Object "
222     << "not found in AssociationMaps (" << typeid(*this).name() << ")." << std::endl;
223    
224     return mitPart;
225 paus 1.1 }