ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/CmsHi/HiGenTools/plugins/HiEventTuner.cc
Revision: 1.1
Committed: Fri Apr 2 21:12:50 2010 UTC (15 years, 1 month ago) by bazterra
Content type: text/plain
Branch: MAIN
Log Message:
Seeding the HiEventTuner.

File Contents

# User Rev Content
1 bazterra 1.1
2     #include <map>
3     #include <memory>
4     #include <string>
5     #include <utility>
6     #include <vector>
7    
8     #include "DataFormats/Common/interface/Handle.h"
9    
10     #include "FWCore/Framework/interface/EDProducer.h"
11     #include "FWCore/Framework/interface/Frameworkfwd.h"
12     #include "FWCore/Framework/interface/MakerMacros.h"
13     #include "FWCore/ParameterSet/interface/ParameterSet.h"
14    
15     #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
16     #include "SimDataFormats/Track/interface/SimTrack.h"
17     #include "SimDataFormats/Track/interface/SimTrackContainer.h"
18     #include "SimDataFormats/Vertex/interface/SimVertex.h"
19     #include "SimDataFormats/Vertex/interface/SimVertexContainer.h"
20    
21     //! Producer to reduce the number of charge in simulation
22     class HiEventTuner : public edm::EDProducer
23     {
24     public:
25    
26     //! Contructor
27     explicit HiEventTuner(const edm::ParameterSet&);
28    
29     //! Destructor
30     ~HiEventTuner(){};
31    
32     private:
33    
34     static std::string MessageCategory_;
35    
36     virtual void produce(edm::Event&, const edm::EventSetup&);
37    
38     bool useMultipleHepMCLabels_;
39     std::vector<std::string> hepMCLabels_;
40     std::vector<edm::Handle<edm::HepMCProduct> > hepMCProducts_;
41    
42     std::string simHitLabel_;
43    
44     edm::Handle<edm::SimTrackContainer> simTracks_;
45     edm::Handle<edm::SimVertexContainer> simVertexes_;
46    
47     std::auto_ptr<edm::SimTrackContainer> tunedSimTracks_;
48     std::auto_ptr<edm::SimVertexContainer> tunedSimVertexes_;
49    
50     typedef std::multimap<std::size_t, std::size_t> IndexToIndexes;
51    
52     IndexToIndexes trackIndexToDecayTrackIndexes_;
53    
54     void mapTrackIdToDecayTrackIndexes();
55    
56     void createTunedCollections();
57    
58     void addSimulatedHistory(SimTrack const &, int parentSimTrackIndex = -1);
59    
60     bool selectSimulatedTrack(SimTrack const &);
61    
62     HepMC::GenParticle * getGeneratorParticle(SimTrack const &);
63     };
64    
65    
66     std::string HiEventTuner::MessageCategory_ = "HiEventTuner";
67    
68    
69     HiEventTuner::HiEventTuner(const edm::ParameterSet& config)
70     {
71     // Initialize global parameters
72     useMultipleHepMCLabels_ = config.getParameter<bool>("useMultipleHepMCLabels");
73     hepMCLabels_ = config.getParameter<std::vector<std::string> >("hepMCLabels");
74     simHitLabel_ = config.getParameter<std::string>("simHitLabel");
75    
76     // Declaration of the products
77     produces<edm::SimTrackContainer>();
78     produces<edm::SimVertexContainer>();
79     }
80    
81    
82     void HiEventTuner::produce(edm::Event& event, const edm::EventSetup& setup)
83     {
84     // Clean the list of hepmc products
85     hepMCProducts_.clear();
86    
87     // Collect all the HepMCProducts
88     edm::Handle<edm::HepMCProduct> hepMCHandle;
89    
90     for (std::vector<std::string>::const_iterator source = hepMCLabels_.begin(); source != hepMCLabels_.end(); ++source)
91     {
92     if ( event.getByLabel(*source, hepMCHandle) )
93     {
94     hepMCProducts_.push_back(hepMCHandle);
95     edm::LogInfo (MessageCategory_) << "Using HepMC source " << *source;
96     if (!useMultipleHepMCLabels_) break;
97     }
98     }
99    
100     // Warning and checking preconditions
101     if ( hepMCProducts_.empty() )
102     {
103     edm::LogWarning (MessageCategory_) << "No HepMC source found";
104     return;
105     }
106     else if ( hepMCProducts_.size() > 1 || useMultipleHepMCLabels_ )
107     {
108     edm::LogInfo (MessageCategory_) << "You are using more than one HepMC source.";
109     edm::LogInfo (MessageCategory_) << "If the labels are not in the same order as the events in the crossing frame (i.e. signal, pileup(s) ) ";
110     edm::LogInfo (MessageCategory_) << "or there are fewer labels than events in the crossing frame";
111     edm::LogInfo (MessageCategory_) << MessageCategory_ << " may try to access data in the wrong HepMCProduct and crash.";
112     }
113    
114     // Collect all the SimTracks
115     event.getByLabel("mix", simHitLabel_, simTracks_);
116    
117     // Collect all the simvertex from the crossing frame
118     event.getByLabel("mix", simHitLabel_, simVertexes_);
119    
120     // Create collections of things we will put in event
121     tunedSimTracks_ = std::auto_ptr<edm::SimTrackContainer>( new edm::SimTrackContainer );
122     tunedSimVertexes_ = std::auto_ptr<edm::SimVertexContainer>( new edm::SimVertexContainer );
123    
124     // Create a map between trackId and vertex (decay) index
125     mapTrackIdToDecayTrackIndexes();
126    
127     // Create the actual tuned collection
128     createTunedCollections();
129    
130     // Put TrackingParticles and TrackingVertices in event
131     event.put(tunedSimTracks_);
132     event.put(tunedSimVertexes_);
133     }
134    
135    
136     void HiEventTuner::mapTrackIdToDecayTrackIndexes()
137     {
138     std::size_t index = 0;
139    
140     // Clear the association map
141     trackIndexToDecayTrackIndexes_.clear();
142    
143     // Loop over the track collection
144     for (edm::SimTrackContainer::const_iterator simTrack = simTracks_->begin(); simTrack != simTracks_->end(); ++simTrack, ++index)
145     {
146     // Check for a source vertex
147     if (simTrack->noVertex()) continue;
148    
149     // Get the source vertex
150     SimVertex const & sourceVertex = simVertexes_->at( simTrack->vertIndex() );
151    
152     // Check if the source vertex has a parent track
153     if (sourceVertex.noParent()) continue;
154    
155     // Get the source vertex
156     SimTrack const & parentTrack = simTracks_->at( sourceVertex.parentIndex() );
157    
158     // Making the association between parentSimTrack and decay index
159     trackIndexToDecayTrackIndexes_.insert( std::make_pair(parentTrack.trackId(), index) );
160     }
161     }
162    
163    
164     void HiEventTuner::createTunedCollections()
165     {
166     // Loop over the tracks looking for the seeding tracks
167     for (edm::SimTrackContainer::const_iterator simTrack = simTracks_->begin(); simTrack != simTracks_->end(); ++simTrack)
168     if ( selectSimulatedTrack(*simTrack) ) addSimulatedHistory(*simTrack);
169     }
170    
171    
172     void HiEventTuner::addSimulatedHistory(SimTrack const & simTrack, int parentSimTrackIndex)
173     {
174     // Add the simTrack to the tuned collection
175     tunedSimTracks_->push_back(simTrack);
176    
177     // Reference to the tunedSimTrack
178     SimTrack & tunedSimTrack = tunedSimTracks_->back();
179    
180     // Add the origin vertex of a selected track
181     if ( !simTrack.noVertex() )
182     {
183     // Add the parent simVertex to the tunedSimVertexes collection
184     tunedSimVertexes_->push_back(
185     SimVertex(simVertexes_->at(simTrack.vertIndex()), parentSimTrackIndex)
186     );
187    
188     // Update the tunedSimTrack vertex index
189     tunedSimTrack.setVertexIndex((int)tunedSimVertexes_->size());
190     }
191    
192     // Define trackId
193     std::size_t simTrackIndex = simTrack.trackId();
194    
195     // Loop over the decay tracks
196     for (
197     IndexToIndexes::const_iterator decaySimTrackIndex = trackIndexToDecayTrackIndexes_.lower_bound(simTrackIndex);
198     decaySimTrackIndex != trackIndexToDecayTrackIndexes_.upper_bound(simTrackIndex);
199     ++decaySimTrackIndex
200     )
201     addSimulatedHistory(simTracks_->at(decaySimTrackIndex->second), (int)tunedSimTracks_->size());
202     }
203    
204    
205     bool HiEventTuner::selectSimulatedTrack(SimTrack const & simTrack)
206     {
207     // This is the basic rule, select SimTracks associate with a GenParticle
208     if (getGeneratorParticle(simTrack)) return true;
209     return false;
210     }
211    
212    
213     HepMC::GenParticle * HiEventTuner::getGeneratorParticle(SimTrack const & simTrack)
214     {
215     // HepMC holder
216     edm::Handle<edm::HepMCProduct> hepmc;
217    
218     // Index to the generator particle
219     int genParticleIndex = simTrack.genpartIndex();
220    
221     // Check in case there is a GenParticle associated to SimTrack
222     if (genParticleIndex >= 0)
223     {
224     // Get the correct HepMC product (for signal and pileup if any)
225     hepmc = (useMultipleHepMCLabels_) ? hepMCProducts_.at(simTrack.eventId().rawId()) : hepmc = hepMCProducts_.at(0);
226    
227     // Return the pointer to the GenParticle if there is any
228     return hepmc->GetEvent()->barcode_to_particle(genParticleIndex);
229     }
230    
231     return 0;
232     }
233    
234    
235     //define this as a plug-in
236     DEFINE_FWK_MODULE(HiEventTuner);