CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
SVTagInfoValidationAnalyzer Class Reference
Inheritance diagram for SVTagInfoValidationAnalyzer:
edm::EDAnalyzer

Public Member Functions

 SVTagInfoValidationAnalyzer (const edm::ParameterSet &)
 
- Public Member Functions inherited from edm::EDAnalyzer
 EDAnalyzer ()
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 

Private Member Functions

virtual void analyze (const edm::Event &, const edm::EventSetup &)
 
void bookRecoToSim (std::string const &)
 
void bookSimToReco (std::string const &)
 
virtual void endJob ()
 
void fillRecoToSim (std::string const &, reco::Vertex const &, TrackingVertexRef const &)
 
void fillSimToReco (std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &)
 

Private Attributes

VertexClassifierByProxy
< reco::SecondaryVertexTagInfoCollection
classifier_
 
edm::Service< TFileServicefs_
 
Int_t n_event
 
Int_t numberVertexClassifier_
 
Int_t rs_total_nall
 
Int_t rs_total_nbsv
 
Int_t rs_total_nbv
 
Int_t rs_total_ncv
 
Int_t rs_total_nlv
 
Int_t rs_total_nsv
 
Int_t sr_total_nall
 
Int_t sr_total_nbsv
 
Int_t sr_total_nbv
 
Int_t sr_total_ncv
 
Int_t sr_total_nlv
 
Int_t sr_total_nsv
 
edm::InputTag svTagInfoProducer_
 
std::map< std::string, TH1D * > TH1Index_
 
Int_t total_nfake
 
Int_t total_nmiss
 
edm::InputTag trackingTruth_
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
typedef WorkerT< EDAnalyzerWorkerType
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 
- Protected Member Functions inherited from edm::EDAnalyzer
CurrentProcessingContext const * currentContext () const
 

Detailed Description

Definition at line 37 of file SVTagInfoValidationAnalyzer.cc.

Constructor & Destructor Documentation

SVTagInfoValidationAnalyzer::SVTagInfoValidationAnalyzer ( const edm::ParameterSet config)
explicit

Definition at line 91 of file SVTagInfoValidationAnalyzer.cc.

References bookRecoToSim(), bookSimToReco(), fs_, edm::ParameterSet::getUntrackedParameter(), i, TFileDirectory::make(), n_event, VertexCategories::Names, numberVertexClassifier_, rs_total_nall, rs_total_nbsv, rs_total_nbv, rs_total_ncv, rs_total_nlv, rs_total_nsv, sr_total_nall, sr_total_nbsv, sr_total_nbv, sr_total_ncv, sr_total_nlv, sr_total_nsv, svTagInfoProducer_, TH1Index_, total_nfake, total_nmiss, trackingTruth_, and VertexCategories::Unknown.

91  : classifier_(config)
92 {
93  //Initialize counters
94  n_event = 0;
95  rs_total_nall = 0;
96  rs_total_nsv = 0;
97  rs_total_nbv = 0;
98  rs_total_nbsv = 0;
99  rs_total_ncv = 0;
100  rs_total_nlv = 0;
101  total_nfake = 0;
102 
103  sr_total_nall = 0;
104  sr_total_nsv = 0;
105  sr_total_nbv = 0;
106  sr_total_nbsv = 0;
107  sr_total_ncv = 0;
108  sr_total_nlv = 0;
109  total_nmiss = 0;
110 
111  // Get the track collection
112  svTagInfoProducer_ = config.getUntrackedParameter<edm::InputTag> ( "svTagInfoProducer" );
113 
114  // Name of the traking pariticle collection
115  trackingTruth_ = config.getUntrackedParameter<edm::InputTag> ( "trackingTruth" );
116 
117  // Number of track categories
119 
120  // Define histogram for counting categories
121  TH1Index_["VertexClassifier"] = fs_->make<TH1D>(
122  "VertexClassifier",
123  "Frequency for the different track categories",
125  -0.5,
126  numberVertexClassifier_ - 0.5
127  );
128 
129  //--- RecoToSim
130  TH1Index_["rs_All_MatchQuality"]= fs_->make<TH1D>( "rs_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01 );
131  TH1Index_["rs_All_FlightDistance2d"]= fs_->make<TH1D>( "rs_All_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
132  TH1Index_["rs_SecondaryVertex_FlightDistance2d"]= fs_->make<TH1D>( "rs_SecondaryVertex_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
133  TH1Index_["rs_BSV_FlightDistance2d"]= fs_->make<TH1D>( "rs_BSV_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
134  TH1Index_["rs_BWeakDecay_FlightDistance2d"]= fs_->make<TH1D>( "rs_BWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
135  TH1Index_["rs_CWeakDecay_FlightDistance2d"]= fs_->make<TH1D>( "rs_CWeakDecay_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
136  TH1Index_["rs_Light_FlightDistance2d"]= fs_->make<TH1D>( "rs_Light_FlightDistance2d", "Transverse flight distance [cm]", 100, 0, 5 );
137 
138  TH1Index_["rs_All_nRecVtx"]= fs_->make<TH1D>( "rs_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
139  TH1Index_["rs_SecondaryVertex_nRecVtx"]= fs_->make<TH1D>( "rs_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
140  TH1Index_["rs_BSV_nRecVtx"]= fs_->make<TH1D>( "rs_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
141  TH1Index_["rs_BWeakDecay_nRecVtx"]= fs_->make<TH1D>( "rs_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
142  TH1Index_["rs_CWeakDecay_nRecVtx"]= fs_->make<TH1D>( "rs_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
143  TH1Index_["rs_Light_nRecVtx"]= fs_->make<TH1D>( "rs_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
144 
145 
146  //--- SimToReco
147  TH1Index_["sr_All_MatchQuality"]= fs_->make<TH1D>( "sr_All_MatchQuality", "Quality of Match", 51, -0.01, 1.01);
148  TH1Index_["sr_All_nRecVtx"]= fs_->make<TH1D>( "sr_All_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
149  TH1Index_["sr_SecondaryVertex_nRecVtx"]= fs_->make<TH1D>( "sr_SecondaryVertex_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
150  TH1Index_["sr_BSV_nRecVtx"]= fs_->make<TH1D>( "sr_BSV_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
151  TH1Index_["sr_BWeakDecay_nRecVtx"]= fs_->make<TH1D>( "sr_BWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
152  TH1Index_["sr_CWeakDecay_nRecVtx"]= fs_->make<TH1D>( "sr_CWeakDecay_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
153  TH1Index_["sr_Light_nRecVtx"]= fs_->make<TH1D>( "sr_Light_nRecVtx", "Number of Vertices per event", 11, -0.5, 10.5 );
154 
155 
156  // Set the proper categories names
157  for (Int_t i = 0; i < numberVertexClassifier_; ++i)
158  TH1Index_["VertexClassifier"]->GetXaxis()->SetBinLabel(i+1, VertexCategories::Names[i]);
159 
160  // book histograms
161  bookRecoToSim("rs_All");
162  bookRecoToSim("rs_SecondaryVertex");
163  bookRecoToSim("rs_BSV");
164  bookRecoToSim("rs_BWeakDecay");
165  bookRecoToSim("rs_CWeakDecay");
166  bookRecoToSim("rs_Light");
167 
168  bookSimToReco("sr_All");
169  bookSimToReco("sr_SecondaryVertex");
170  bookSimToReco("sr_BSV");
171  bookSimToReco("sr_BWeakDecay");
172  bookSimToReco("sr_CWeakDecay");
173  bookSimToReco("sr_Light");
174 }
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
edm::Service< TFileService > fs_
void bookSimToReco(std::string const &)
std::map< std::string, TH1D * > TH1Index_
void bookRecoToSim(std::string const &)
VertexClassifierByProxy< reco::SecondaryVertexTagInfoCollection > classifier_
T * make() const
make new ROOT object
static const char * Names[]
Name of the different categories.

Member Function Documentation

void SVTagInfoValidationAnalyzer::analyze ( const edm::Event event,
const edm::EventSetup setup 
)
privatevirtual

Implements edm::EDAnalyzer.

Definition at line 177 of file SVTagInfoValidationAnalyzer.cc.

References VertexCategories::BWeakDecay, classifier_, gather_cfg::cout, VertexCategories::CWeakDecay, VertexClassifierByProxy< Collection >::evaluate(), VertexCategories::Fake, fillRecoToSim(), fillSimToReco(), VertexClassifier::history(), i, getHLTprescales::index, VertexCategories::is(), n_event, VertexClassifierByProxy< Collection >::newEvent(), numberVertexClassifier_, VertexHistory::quality(), VertexCategories::Reconstructed, VertexHistory::recoVertex(), rs_total_nall, rs_total_nbsv, rs_total_nbv, rs_total_ncv, rs_total_nlv, rs_total_nsv, VertexCategories::SecondaryVertex, HistoryBase::simVertex(), sr_total_nall, sr_total_nbsv, sr_total_nbv, sr_total_ncv, sr_total_nlv, sr_total_nsv, svTagInfoProducer_, TH1Index_, total_nfake, total_nmiss, and trackingTruth_.

178 {
179  ++n_event;
180 
181  std::cout << "*** Analyzing " << event.id() << " n_event = " << n_event << std::endl << std::endl;
182 
183  // Set the classifier for a new event
184  classifier_.newEvent(event, setup);
185 
186 
187  // Vertex collection
189  event.getByLabel(svTagInfoProducer_, svTagInfoCollection);
190 
191  // Get a constant reference to the track history associated to the classifier
192  VertexHistory const & tracer = classifier_.history();
193 
194  cout << "* Event " << n_event << " ; svTagInfoCollection->size() = " << svTagInfoCollection->size() << endl;
195 
196  int rs_nall = 0;
197  int rs_nsv = 0;
198  int rs_nbv = 0;
199  int rs_nbsv = 0;
200  int rs_ncv = 0;
201  int rs_nlv = 0;
202  int nfake = 0;
203 
204  int sr_nall = 0;
205  int sr_nsv = 0;
206  int sr_nbv = 0;
207  int sr_nbsv = 0;
208  int sr_ncv = 0;
209  int sr_nlv = 0;
210  int nmiss = 0;
211 
212  // Loop over the svTagInfo collection.
213  for (std::size_t index = 0; index < svTagInfoCollection->size(); ++index){
214 
215  reco::SecondaryVertexTagInfoRef svTagInfo(svTagInfoCollection, index);
216 
217  // Loop over the vertexes in svTagInfo
218  for ( std::size_t vindex = 0; vindex < svTagInfo->nVertices(); ++vindex ){
219 
220 
221  // Classify the vertices
222  classifier_.evaluate(svTagInfo, vindex);
223 
224  //quality of the match
225  double rs_quality = tracer.quality();
226 
227 
228  // Fill the histogram with the categories
229  for (Int_t i = 0; i != numberVertexClassifier_; ++i) {
230 
232 
233  TH1Index_["VertexClassifier"]->Fill(i);
234  }
235  }
237 
238  cout << "R2S: MatchQuality = " << rs_quality << endl;
239 
240  TH1Index_["rs_All_MatchQuality"]->Fill( rs_quality );
241  fillRecoToSim("rs_All", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
242  TH1Index_["rs_All_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
243  rs_nall++;
244 
246 
247  //cout << " R2S VertexCategories::SecondaryVertex" << endl;
248  fillRecoToSim("rs_SecondaryVertex", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
249  TH1Index_["rs_SecondaryVertex_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
250  rs_nsv++;
251  }
252 
254 
255  //cout << " R2S VertexCategories::BWeakDecay" << endl;
256  fillRecoToSim("rs_BWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
257  TH1Index_["rs_BWeakDecay_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
258  rs_nbv++;
259 
261  //cout << " R2S VertexCategories::BWeakDecay SecondaryVertex" << endl;
262  fillRecoToSim("rs_BSV", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
263  TH1Index_["rs_BSV_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
264  rs_nbsv++;
265 
266  }
267  }//BWeakDecay
268 
270 
271  //cout << " R2S VertexCategories::CWeakDecay" << endl;
272  fillRecoToSim("rs_CWeakDecay", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
273  TH1Index_["rs_CWeakDecay_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
274  rs_ncv++;
275 
276  }
277  else {
278  //cout << " R2S Light (rest of categories)" << endl;
279  fillRecoToSim("rs_Light", svTagInfo->secondaryVertex(vindex), tracer.simVertex());
280  TH1Index_["rs_Light_FlightDistance2d"]->Fill( svTagInfo->flightDistance( vindex, true ).value() );
281  rs_nlv++;
282  }
283  }//end if classifier
284 
285  else {
286  cout << " VertexCategories::Fake!!" << endl;
287  nfake++;
288  }
289 
290  }//end loop over vertices in svTagInfo
291 
292  }//loop over svTagInfo
293 
294  TH1Index_["rs_All_nRecVtx"]->Fill( rs_nall );
295  TH1Index_["rs_SecondaryVertex_nRecVtx"]->Fill( rs_nsv );
296  TH1Index_["rs_BWeakDecay_nRecVtx"]->Fill( rs_nbv );
297  TH1Index_["rs_BSV_nRecVtx"]->Fill( rs_nbsv );
298  TH1Index_["rs_CWeakDecay_nRecVtx"]->Fill( rs_ncv );
299  TH1Index_["rs_Light_nRecVtx"]->Fill( rs_nlv );
300  cout << endl;
301 
302  //----------------------------------------------------------------
303  // SIM TO RECO!
304 
305  // Vertex collection
307  event.getByLabel(trackingTruth_, TVCollection);
308 
309  // Loop over the TV collection.
310  for (std::size_t index = 0; index < TVCollection->size(); ++index){
311 
312  TrackingVertexRef trackingVertex(TVCollection, index);
313 
314  classifier_.evaluate(trackingVertex);
315 
316  double sr_quality = tracer.quality();
317 
319 
320  cout << "S2R: MatchQuality = " << sr_quality << endl;
321 
322  //cout << " S2R VertexCategories::Reconstructed" << endl;
323  TH1Index_["sr_All_MatchQuality"]->Fill( sr_quality );
324  fillSimToReco("sr_All", tracer.recoVertex(), trackingVertex);
325  sr_nall++;
326 
328 
329  //cout << " S2R VertexCategories::Reconstructed::SecondaryVertex" << endl;
330  fillSimToReco("sr_SecondaryVertex", tracer.recoVertex(), trackingVertex);
331  sr_nsv++;
332  }
333 
335 
336  //cout << " S2R VertexCategories::Reconstructed::BWeakDecay" << endl;
337  fillSimToReco("sr_BWeakDecay", tracer.recoVertex(), trackingVertex);
338  sr_nbv++;
339 
341 
342  //cout << " S2R VertexCategories::Reconstructed::BWeakDecay SecondaryVertex" << endl;
343  fillSimToReco("sr_BSV", tracer.recoVertex(), trackingVertex);
344  sr_nbsv++;
345  }
346 
347  }//BWeakDecay
348 
350 
351  //cout << " S2R VertexCategories::CWeakDecay" << endl;
352  fillSimToReco("sr_CWeakDecay", tracer.recoVertex(), trackingVertex);
353  sr_ncv++;
354  }
355 
356  else {
357 
358  //cout << " S2R Light (rest of categories)" << endl;
359  fillSimToReco("sr_Light", tracer.recoVertex(), trackingVertex);
360  sr_nlv++;
361  }
362 
363  }//Reconstructed
364  else {
365  //cout << "##### Not reconstructed!" << endl;
366  nmiss++;
367  }
368 
369  }//TVCollection.size()
370 
371  TH1Index_["sr_All_nRecVtx"]->Fill( sr_nall );
372  TH1Index_["sr_SecondaryVertex_nRecVtx"]->Fill( sr_nsv );
373  TH1Index_["sr_BWeakDecay_nRecVtx"]->Fill( sr_nbv );
374  TH1Index_["sr_BSV_nRecVtx"]->Fill( sr_nbsv );
375  TH1Index_["sr_CWeakDecay_nRecVtx"]->Fill( sr_ncv );
376  TH1Index_["sr_Light_nRecVtx"]->Fill( rs_nlv );
377 
378  rs_total_nall += rs_nall;
379  rs_total_nsv += rs_nsv;
380  rs_total_nbv += rs_nbv;
381  rs_total_nbsv += rs_nbsv;
382  rs_total_ncv += rs_ncv ;
383  rs_total_nlv += rs_nlv;
384  total_nfake += nfake;
385 
386  sr_total_nall += sr_nall;
387  sr_total_nsv += sr_nsv;
388  sr_total_nbv += sr_nbv;
389  sr_total_nbsv += sr_nbsv;
390  sr_total_ncv += sr_ncv;
391  sr_total_nlv += sr_nlv;
392  total_nmiss += nmiss;
393 
394 
395 }
This class traces the simulated and generated history of a given track.
Definition: VertexHistory.h:17
int i
Definition: DBlmapReader.cc:9
VertexHistory const & history() const
Returns a reference to the vertex history used in the classification.
const reco::VertexBaseRef & recoVertex() const
Return a reference to the reconstructed track.
Definition: VertexHistory.h:62
void fillRecoToSim(std::string const &, reco::Vertex const &, TrackingVertexRef const &)
Category
Categories available to vertexes.
std::map< std::string, TH1D * > TH1Index_
bool is(Category category) const
Returns track flag for a given category.
const TrackingVertexRef & simVertex() const
Return the initial tracking vertex from the history.
Definition: HistoryBase.h:83
VertexClassifierByProxy< Collection > const & evaluate(TrackingVertexRef const &vertex)
Classify the TrackingVertex in categories.
void fillSimToReco(std::string const &, reco::VertexBaseRef const &, TrackingVertexRef const &)
VertexClassifierByProxy< reco::SecondaryVertexTagInfoCollection > classifier_
virtual void newEvent(edm::Event const &event, edm::EventSetup const &config)
Pre-process event information (for accessing reconstraction information).
tuple cout
Definition: gather_cfg.py:121
double quality() const
Return the quality of the match.
Definition: VertexHistory.h:68
void SVTagInfoValidationAnalyzer::bookRecoToSim ( std::string const &  prefix)
private

Definition at line 398 of file SVTagInfoValidationAnalyzer.cc.

References fs_, TFileDirectory::make(), mergeVDriftHistosByStation::name, and TH1Index_.

Referenced by SVTagInfoValidationAnalyzer().

398  {
399  // Book pull histograms
400 
401  std::string name = prefix + "_Pullx";
402  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
403  name = prefix + "_Pully";
404  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
405  name = prefix + "_Pullz";
406  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
407 
408  name = prefix + "_Resx";
409  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
410  name = prefix + "_Resy";
411  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
412  name = prefix + "_Resz";
413  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
414 
415  name = prefix + "_Chi2Norm";
416  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, 0, 10.);
417  name = prefix + "_Chi2Prob";
418  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, 0., 1.);
419 
420  name = prefix + "_nRecTrks";
421  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 501, -0.5, 500.5);
422 
423  name = prefix + "_AverageTrackWeight";
424  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.1, 1.1);
425 
426  name = prefix + "_Mass";
427  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 65, 0., 6.5);
428 
429  name = prefix + "_RecPt";
430  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
431 
432  name = prefix + "_RecEta";
433  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
434 
435  name = prefix + "_RecCharge";
436  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 21, -0.5, 20.5);
437 
438  name = prefix + "_RecTrackPt";
439  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
440 
441  name = prefix + "_RecTrackEta";
442  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
443 
444  name = prefix + "_nSimTrks";
445  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 501, -0.5, 500.5);
446 
447  name = prefix + "_SimPt";
448  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
449 
450  name = prefix + "_SimEta";
451  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
452 
453  name = prefix + "_SimCharge";
454  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 21, -0.5, 20.5);
455 
456  name = prefix + "_SimTrackPt";
457  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 500, 0., 500.);
458 
459  name = prefix + "_SimTrackEta";
460  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
461 
462 
463 }
edm::Service< TFileService > fs_
std::map< std::string, TH1D * > TH1Index_
T * make() const
make new ROOT object
void SVTagInfoValidationAnalyzer::bookSimToReco ( std::string const &  prefix)
private

Definition at line 466 of file SVTagInfoValidationAnalyzer.cc.

References fs_, TFileDirectory::make(), mergeVDriftHistosByStation::name, and TH1Index_.

Referenced by SVTagInfoValidationAnalyzer().

466  {
467  // Book pull histograms
468 
469  std::string name = prefix + "_Pullx";
470  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
471  name = prefix + "_Pully";
472  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
473  name = prefix + "_Pullz";
474  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -10., 10.);
475 
476  name = prefix + "_Resx";
477  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
478  name = prefix + "_Resy";
479  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
480  name = prefix + "_Resz";
481  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.05, 0.05);
482 
483  name = prefix + "_Chi2Norm";
484  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, 0, 10.);
485  name = prefix + "_Chi2Prob";
486  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, 0., 1.);
487 
488  name = prefix + "_nRecTrks";
489  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 501, -0.5, 500.5);
490 
491  name = prefix + "_AverageTrackWeight";
492  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 100, -0.1, 1.1);
493 
494  name = prefix + "_Mass";
495  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 65, 0., 6.5);
496 
497  name = prefix + "_RecPt";
498  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
499 
500  name = prefix + "_RecEta";
501  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
502 
503  name = prefix + "_RecCharge";
504  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 21, -0.5, 20.5);
505 
506  name = prefix + "_RecTrackPt";
507  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
508 
509  name = prefix + "_RecTrackEta";
510  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
511 
512  name = prefix + "_nSimTrks";
513  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 501, -0.5, 500.5);
514 
515  name = prefix + "_SimPt";
516  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 2000, 0., 1000.);
517 
518  name = prefix + "_SimEta";
519  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
520 
521  name = prefix + "_SimCharge";
522  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 21, -0.5, 20.5);
523 
524  name = prefix + "_SimTrackPt";
525  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 500, 0., 500.);
526 
527  name = prefix + "_SimTrackEta";
528  TH1Index_[name] = fs_->make<TH1D>(name.c_str(), name.c_str(), 200, -3., 3.);
529 
530 
531 }
edm::Service< TFileService > fs_
std::map< std::string, TH1D * > TH1Index_
T * make() const
make new ROOT object
void SVTagInfoValidationAnalyzer::endJob ( void  )
privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 717 of file SVTagInfoValidationAnalyzer.cc.

References gather_cfg::cout, n_event, rs_total_nall, rs_total_nbsv, rs_total_nbv, rs_total_ncv, rs_total_nlv, rs_total_nsv, sr_total_nall, sr_total_nbsv, sr_total_nbv, sr_total_ncv, sr_total_nlv, sr_total_nsv, and total_nfake.

717  {
718 
719  std::cout << std::endl;
720  std::cout << " ====== Total Number of analyzed events: " << n_event << " ====== " << std::endl;
721  std::cout << " ====== Total Number of R2S All: " << rs_total_nall << " ====== " << std::endl;
722  std::cout << " ====== Total Number of R2S SecondaryVertex: " << rs_total_nsv << " ====== " << std::endl;
723  std::cout << " ====== Total Number of R2S BWeakDecay: " << rs_total_nbv << " ====== " << std::endl;
724  std::cout << " ====== Total Number of R2S BWeakDecay::SecondaryVertex: " << rs_total_nbsv << " ====== " << std::endl;
725  std::cout << " ====== Total Number of R2S CWeakDecay: " << rs_total_ncv << " ====== " << std::endl;
726  std::cout << " ====== Total Number of R2S Light: " << rs_total_nlv << " ====== " << std::endl;
727  std::cout << std::endl;
728  std::cout << " ====== Total Number of S2R All: " << sr_total_nall << " ====== " << std::endl;
729  std::cout << " ====== Total Number of S2R SecondaryVertex: " << sr_total_nsv << " ====== " << std::endl;
730  std::cout << " ====== Total Number of S2R BWeakDecay: " << sr_total_nbv << " ====== " << std::endl;
731  std::cout << " ====== Total Number of S2R BWeakDecay::SecondaryVertex: " << sr_total_nbsv << " ====== " << std::endl;
732  std::cout << " ====== Total Number of S2R CWeakDecay: " << sr_total_ncv << " ====== " << std::endl;
733  std::cout << " ====== Total Number of S2R Light: " << sr_total_nlv << " ====== " << std::endl;
734  std::cout << std::endl;
735  std::cout << " ====== Total Number of Fake Vertices: " << total_nfake << " ====== " << std::endl;
736 
737 }
tuple cout
Definition: gather_cfg.py:121
void SVTagInfoValidationAnalyzer::fillRecoToSim ( std::string const &  prefix,
reco::Vertex const &  vertex,
TrackingVertexRef const &  simVertex 
)
private

Definition at line 534 of file SVTagInfoValidationAnalyzer.cc.

References DeDxDiscriminatorTools::charge(), reco::Vertex::chi2(), ChiSquaredProbability(), eta(), scaleCards::mass, reco::Vertex::ndof(), reco::Vertex::normalizedChi2(), AlCaHLTBitMon_ParallelJobs::p, mathSSE::sqrt(), TH1Index_, reco::Vertex::tracks_begin(), reco::Vertex::tracks_end(), reco::Vertex::tracksSize(), reco::Vertex::trackWeight(), CommonMethods::weight(), reco::Vertex::x(), reco::Vertex::xError(), reco::Vertex::y(), reco::Vertex::yError(), reco::Vertex::z(), and reco::Vertex::zError().

Referenced by analyze().

535 {
536 
537  double pullx = (vertex.x() - simVertex->position().x())/vertex.xError();
538  double pully = (vertex.y() - simVertex->position().y())/vertex.yError();
539  double pullz = (vertex.z() - simVertex->position().z())/vertex.zError();
540 
541  double resx = vertex.x() - simVertex->position().x();
542  double resy = vertex.y() - simVertex->position().y();
543  double resz = vertex.z() - simVertex->position().z();
544 
545  double chi2norm = vertex.normalizedChi2();
546  double chi2prob = ChiSquaredProbability( vertex.chi2(), vertex.ndof() );
547 
548  double sum_weight = 0.;
549  double weight = 0.;
550  double tracksize = vertex.tracksSize();
551  math::XYZVector momentum;
553  int charge = 0;
554  double thePiMass = 0.13957;
555  for ( reco::Vertex::trackRef_iterator recDaughter = vertex.tracks_begin() ; recDaughter != vertex.tracks_end(); ++recDaughter) {
556 
557  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double> > vec;
558 
559  vec.SetPx( (**recDaughter).px() );
560  vec.SetPy( (**recDaughter).py() );
561  vec.SetPz( (**recDaughter).pz() );
562  vec.SetM(thePiMass);
563 
564  sum += vec;
565 
566  weight = vertex.trackWeight(*recDaughter);
567  sum_weight += weight;
568 
570  p = (**recDaughter).momentum();
571  momentum += p;
572 
573  charge += (*recDaughter)->charge();
574 
575  TH1Index_[prefix + "_RecTrackPt"]->Fill( (*recDaughter)->pt() );
576  TH1Index_[prefix + "_RecTrackEta"]->Fill( (*recDaughter)->eta() );
577  }//end loop to recDaughters
578  //cout << " average sum of weights = " << sum_weight/tracksize << endl;
579 
580  double mass = sum.M();
581  double pt = sqrt( momentum.Perp2() );
582  double eta = momentum.Eta();
583 
584  math::XYZVector simmomentum;
585  int simcharge = 0;
586  for ( TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin(); simDaughter != simVertex->daughterTracks_end(); ++simDaughter ){
587 
589  p = (**simDaughter).momentum();
590  simmomentum += p;
591 
592  simcharge += (*simDaughter)->charge();
593 
594  TH1Index_[prefix + "_SimTrackPt"]->Fill( (*simDaughter)->pt() );
595  TH1Index_[prefix + "_SimTrackEta"]->Fill( (*simDaughter)->eta() );
596  }
597 
598  double simpt = sqrt( simmomentum.Perp2() );
599  double simeta = simmomentum.Eta();
600 
601  //cout << "[fillRecoToSim] vertex.tracksSize() = " << vertex.tracksSize() << " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() << endl;
602 
603  TH1Index_[prefix + "_nRecTrks"]->Fill( vertex.tracksSize() );
604  TH1Index_[prefix + "_nSimTrks"]->Fill( simVertex->nDaughterTracks() );
605  TH1Index_[prefix + "_Pullx"]->Fill(pullx);
606  TH1Index_[prefix + "_Pully"]->Fill(pully);
607  TH1Index_[prefix + "_Pullz"]->Fill(pullz);
608  TH1Index_[prefix + "_Resx"]->Fill(resx);
609  TH1Index_[prefix + "_Resy"]->Fill(resy);
610  TH1Index_[prefix + "_Resz"]->Fill(resz);
611  TH1Index_[prefix + "_AverageTrackWeight"]->Fill( sum_weight/tracksize );
612  TH1Index_[prefix + "_Chi2Norm"]->Fill(chi2norm);
613  TH1Index_[prefix + "_Chi2Prob"]->Fill(chi2prob);
614  TH1Index_[prefix + "_RecPt"]->Fill(pt);
615  TH1Index_[prefix + "_RecEta"]->Fill(eta);
616  TH1Index_[prefix + "_RecCharge"]->Fill(charge);
617  TH1Index_[prefix + "_Mass"]->Fill(mass);
618  TH1Index_[prefix + "_SimPt"]->Fill(simpt);
619  TH1Index_[prefix + "_SimEta"]->Fill(simeta);
620  TH1Index_[prefix + "_SimCharge"]->Fill(simcharge);
621 }
T eta() const
double charge(const std::vector< uint8_t > &Ampls)
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
T sqrt(T t)
Definition: SSEVec.h:46
std::map< std::string, TH1D * > TH1Index_
float ChiSquaredProbability(double chiSquared, double nrDOF)
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
tuple mass
Definition: scaleCards.py:27
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector&lt;TrackRef&gt;
Definition: Vertex.h:38
void SVTagInfoValidationAnalyzer::fillSimToReco ( std::string const &  prefix,
reco::VertexBaseRef const &  vertex,
TrackingVertexRef const &  simVertex 
)
private

Definition at line 624 of file SVTagInfoValidationAnalyzer.cc.

References DeDxDiscriminatorTools::charge(), reco::Vertex::chi2(), ChiSquaredProbability(), eta(), scaleCards::mass, reco::Vertex::ndof(), reco::Vertex::normalizedChi2(), AlCaHLTBitMon_ParallelJobs::p, mathSSE::sqrt(), TH1Index_, reco::Vertex::tracks_begin(), reco::Vertex::tracks_end(), reco::Vertex::tracksSize(), reco::Vertex::trackWeight(), CommonMethods::weight(), reco::Vertex::x(), reco::Vertex::xError(), reco::Vertex::y(), reco::Vertex::yError(), reco::Vertex::z(), and reco::Vertex::zError().

Referenced by analyze().

625 {
626 
627  double pullx = (vertex->x() - simVertex->position().x())/vertex->xError();
628  double pully = (vertex->y() - simVertex->position().y())/vertex->yError();
629  double pullz = (vertex->z() - simVertex->position().z())/vertex->zError();
630 
631  double resx = vertex->x() - simVertex->position().x();
632  double resy = vertex->y() - simVertex->position().y();
633  double resz = vertex->z() - simVertex->position().z();
634 
635  double chi2norm = vertex->normalizedChi2();
636  double chi2prob = ChiSquaredProbability(vertex->chi2(), vertex->ndof());
637 
638  double sum_weight = 0.;
639  double weight = 0.;
640  double tracksize = vertex->tracksSize();
641  math::XYZVector momentum;
643  int charge = 0;
644  double thePiMass = 0.13957;
645  for ( reco::Vertex::trackRef_iterator recDaughter = vertex->tracks_begin() ; recDaughter != vertex->tracks_end(); ++recDaughter ) {
646 
647  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double> > vec;
648 
649  vec.SetPx( (**recDaughter).px() );
650  vec.SetPy( (**recDaughter).py() );
651  vec.SetPz( (**recDaughter).pz() );
652  vec.SetM(thePiMass);
653 
654  sum += vec;
655 
656  weight = vertex->trackWeight(*recDaughter);
657  sum_weight += weight;
658 
660  p = (**recDaughter).momentum();
661  momentum += p;
662 
663  charge += (*recDaughter)->charge();
664 
665  TH1Index_[prefix + "_RecTrackPt"]->Fill( (*recDaughter)->pt() );
666  TH1Index_[prefix + "_RecTrackEta"]->Fill( (*recDaughter)->eta() );
667 
668  }
669  //cout << " average sum of weights = " << sum_weight/tracksize << endl;
670 
671  double mass = sum.M();
672  double pt = sqrt( momentum.Perp2() );
673  double eta = momentum.Eta();
674 
675  math::XYZVector simmomentum;
676  int simcharge = 0;
677  for ( TrackingVertex::tp_iterator simDaughter = simVertex->daughterTracks_begin(); simDaughter != simVertex->daughterTracks_end(); ++simDaughter ){
678 
680  p = (**simDaughter).momentum();
681  simmomentum += p;
682 
683  simcharge += (*simDaughter)->charge();
684 
685  TH1Index_[prefix + "_SimTrackPt"]->Fill( (*simDaughter)->pt() );
686  TH1Index_[prefix + "_SimTrackEta"]->Fill( (*simDaughter)->eta() );
687  }
688 
689  double simpt = sqrt( simmomentum.Perp2() );
690  double simeta = simmomentum.Eta();
691 
692  //cout << "[fillSimToReco] vertex->tracksSize() = " << vertex->tracksSize() << " ; simVertex->nDaughterTracks() = " << simVertex->nDaughterTracks() << endl;
693 
694  TH1Index_[prefix + "_nRecTrks"]->Fill( vertex->tracksSize() );
695  TH1Index_[prefix + "_nSimTrks"]->Fill( simVertex->nDaughterTracks() );
696  TH1Index_[prefix + "_Pullx"]->Fill(pullx);
697  TH1Index_[prefix + "_Pully"]->Fill(pully);
698  TH1Index_[prefix + "_Pullz"]->Fill(pullz);
699  TH1Index_[prefix + "_Resx"]->Fill(resx);
700  TH1Index_[prefix + "_Resy"]->Fill(resy);
701  TH1Index_[prefix + "_Resz"]->Fill(resz);
702  TH1Index_[prefix + "_AverageTrackWeight"]->Fill( sum_weight/tracksize );
703  TH1Index_[prefix + "_Chi2Norm"]->Fill(chi2norm);
704  TH1Index_[prefix + "_Chi2Prob"]->Fill(chi2prob);
705  TH1Index_[prefix + "_RecPt"]->Fill(pt);
706  TH1Index_[prefix + "_RecEta"]->Fill(eta);
707  TH1Index_[prefix + "_RecCharge"]->Fill(charge);
708  TH1Index_[prefix + "_Mass"]->Fill(mass);
709  TH1Index_[prefix + "_SimPt"]->Fill(simpt);
710  TH1Index_[prefix + "_SimEta"]->Fill(simeta);
711  TH1Index_[prefix + "_SimCharge"]->Fill(simcharge);
712 
713 }
T eta() const
double charge(const std::vector< uint8_t > &Ampls)
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
T sqrt(T t)
Definition: SSEVec.h:46
std::map< std::string, TH1D * > TH1Index_
float ChiSquaredProbability(double chiSquared, double nrDOF)
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
tuple mass
Definition: scaleCards.py:27
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector&lt;TrackRef&gt;
Definition: Vertex.h:38

Member Data Documentation

VertexClassifierByProxy<reco::SecondaryVertexTagInfoCollection> SVTagInfoValidationAnalyzer::classifier_
private

Definition at line 50 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze().

edm::Service<TFileService> SVTagInfoValidationAnalyzer::fs_
private
Int_t SVTagInfoValidationAnalyzer::n_event
private

Definition at line 59 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::numberVertexClassifier_
private

Definition at line 52 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_nall
private

Definition at line 60 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_nbsv
private

Definition at line 63 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_nbv
private

Definition at line 62 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_ncv
private

Definition at line 64 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_nlv
private

Definition at line 65 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::rs_total_nsv
private

Definition at line 61 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_nall
private

Definition at line 68 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_nbsv
private

Definition at line 71 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_nbv
private

Definition at line 70 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_ncv
private

Definition at line 72 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_nlv
private

Definition at line 73 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::sr_total_nsv
private

Definition at line 69 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

edm::InputTag SVTagInfoValidationAnalyzer::svTagInfoProducer_
private

Definition at line 54 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), and SVTagInfoValidationAnalyzer().

std::map<std::string, TH1D *> SVTagInfoValidationAnalyzer::TH1Index_
private
Int_t SVTagInfoValidationAnalyzer::total_nfake
private

Definition at line 66 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), endJob(), and SVTagInfoValidationAnalyzer().

Int_t SVTagInfoValidationAnalyzer::total_nmiss
private

Definition at line 74 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), and SVTagInfoValidationAnalyzer().

edm::InputTag SVTagInfoValidationAnalyzer::trackingTruth_
private

Definition at line 53 of file SVTagInfoValidationAnalyzer.cc.

Referenced by analyze(), and SVTagInfoValidationAnalyzer().