LCFIVertex  0.7.2
DSTPlotProcessor.cc
1 
2 #include "DSTPlotProcessor.h"
3 #include <iostream>
4 #include <fstream>
5 #include <string>
6 #include <sstream>
7 #include <vector>
8 #include <cmath>
9 #include <set>
10 #ifdef USEROOT
11 #include "TFile.h"
13 #include "TString.h"
14 #include "TGraph.h"
15 #include "TH1.h"
16 #include "TCanvas.h"
17 #include "TMultiGraph.h"
19 #endif // of #ifdef USEROOT
20 
21 #include "EVENT/LCCollection.h"
22 #include "UTIL/PIDHandler.h"
23 #include "IMPL/ReconstructedParticleImpl.h"
24 #include "EVENT/LCParameters.h"
25 #include "EVENT/LCFloatVec.h"
26 #include "EVENT/LCIntVec.h"
27 
28 #include "util/inc/vector3.h"
29 #include "util/inc/util.h"
30 
31 using std::includes;
32 using std::map;
33 using std::string;
34 using std::endl;
35 using std::cout;
36 using std::vector;
37 using std::stringstream;
38 using std::abs;
39 using std::pair;
40 using std::ofstream;
41 using std::cerr;
42 using std::set;
43 using namespace lcio;
44 
45 //Needs to be instantiated for Marlin to know about it (I think)
46 DSTPlotProcessor aDSTPlotProcessor;
47 
48 DSTPlotProcessor::DSTPlotProcessor() : marlin::Processor("DSTPlotProcessor")
49 {
50  _description = "Plots various outputs from the flavour tag" ;
51 
52  // register steering parameters: name, description, class-variable, default value
53  //The name of the collection of ReconstructedParticles that is the jet
54  registerInputCollection( LCIO::RECONSTRUCTEDPARTICLE,
55  "JetCollectionName" ,
56  "Name of the collection of ReconstructedParticles that is the jet" ,
57  _JetCollectionName ,
58  string("FTSelectedJets") ) ;
59 
60  //The output filename
61  registerProcessorParameter( "OutputFilename" ,
62  "Filename for the output" ,
63  _OutputFilename ,
64  string("DSTPlotProcessorOutput") ) ;
65 
66  registerProcessorParameter( "CheckDSTParameters" ,
67  "Checks the DST parameters against the full parameters - gives a lot of print out" ,
68  _checkDST ,
69  int(0) ) ;
70 
71 
72 
73  //Add these collections id you want to check the DST collections against the full collections
74  registerOptionalParameter("FlavourTagCollection" ,
75  "Name of the LCFloatVec Collection containing the flavour tags, only need this if you want to check DST parameters" ,
76  _FlavourTagCollectionName,
77  std::string("FlavourTag") ) ;
78 
79  registerOptionalParameter("FlavourTagInputsCollection" ,
80  "Name of the LCFloatVec Collection that contains the flavour tag inputs, only need this if you want to check DST parameters" ,
81  _FlavourTagInputsCollectionName,
82  std::string("FlavourTagInputs") ) ;
83 
84 
85 
86  registerOptionalParameter("TrueJetFlavourCollection" ,
87  "Name of the LCIntVec collection containing the true flavour of the jets, only need this if you want to check DST parameters" ,
88  _TrueJetFlavourColName ,
89  std::string("TrueJetFlavour") ) ;
90 
91 
92 }
93 
94 DSTPlotProcessor::~DSTPlotProcessor()
95 {
96 }
97 
98 void DSTPlotProcessor::init()
99 {
100  printParameters();
101  cout << _description << endl
102  << "-------------------------------------------------" << endl
103  << endl;
104 
105  _nRun=0;
106  _jetEMax=0;
107 
108 }
109 
110 void DSTPlotProcessor::processRunHeader( LCRunHeader* pRun )
111 {
112  //
113 
114  //map the parameter names to the indexes, for the flavout tags, NN input vales, and mc truth info, if the DST parameters are to be checked
115  if(_checkDST==1)
116  {
117  std::vector<std::string> VarNames;
118  (pRun->parameters()).getStringVals(_FlavourTagCollectionName,VarNames);
119  for (size_t i = 0;i < VarNames.size();++i)
120  _IndexOfForEachTag[VarNames[i]] = i;
121 
122 
123  //do the same for the true jet values
124  std::vector<std::string> trueJetFlavourVarNames;
125  (pRun->parameters()).getStringVals(_TrueJetFlavourColName,trueJetFlavourVarNames);
126  for (size_t i = 0;i < trueJetFlavourVarNames.size();++i)
127  _FlavourIndex[trueJetFlavourVarNames[i]] = i;
128 
129 
130  //do the same for the FT inputs
131  std::vector<std::string> InputVarNames;
132  (pRun->parameters()).getStringVals(_FlavourTagInputsCollectionName,InputVarNames);
133  for (size_t i = 0;i < InputVarNames.size();++i)
134  _InputsIndex[InputVarNames[i]] = i;
135  }
136 
137 
138 }
139 
140 void DSTPlotProcessor::processEvent( LCEvent* pEvent )
141 {
142 
143  if(_checkDST==1)
144  _checkDSTParameters( pEvent);
145 
146  //Get the collection of jets. Can't do anything if the collection isn't there
147  //so don't bother catching the exception and terminate.
148  LCCollection* pJetCollection=pEvent->getCollection( _JetCollectionName );
149 
150  //make sure the collection is of the right type
151  if( pJetCollection->getTypeName()!=LCIO::RECONSTRUCTEDPARTICLE )
152  {
153  stringstream message;
154  message << endl
155  << "########################################################################################\n"
156  << "# DSTPlotProcessor - #\n"
157  << "# The jet collection requested (\"" << _JetCollectionName << "\") is not of the type \"" << LCIO::RECONSTRUCTEDPARTICLE << "\" #\n"
158  << "########################################################################################" << endl;
159  throw EventException( message.str() );
160  }
161 
162  //apply any cuts on the event here
163  if( _passesEventCuts(pEvent) )
164  {
165  ReconstructedParticle* pJet;
166  //loop over the jets
167  for( int a=0; a<pJetCollection->getNumberOfElements(); ++a )
168  {
169  //Dynamic casts are not the best programming practice in the world, but I can't see another way of doing this
170  //in the LCIO framework. This cast should be safe though because we've already tested the type.
171  pJet=dynamic_cast<ReconstructedParticle*>( pJetCollection->getElementAt(a) );
172 
173  if( _passesJetCuts(pJet) )
174  {
175  _fillPlots(pEvent, a );
176  _jetEnergy.add(pJet->getEnergy());
177  if (pJet->getEnergy() > _jetEMax) _jetEMax = pJet->getEnergy();
178  }
179  }
180  }
181 }
182 
183 
184 void DSTPlotProcessor::end()
185 {
187 
188  cout << "The largest jet energy was " << _jetEMax << endl;
189 }
190 
191 // IMPORTANT - If you change the cuts make sure you change the line below to show the changes in the docs
193 bool DSTPlotProcessor::_passesEventCuts( LCEvent* pEvent )
194 {
195  //
196  // No event cuts at present
197  //
198 
199  return true;
200 }
201 
202 // IMPORTANT - If you change the cuts make sure you change the line below to show the changes in the docs
205 bool DSTPlotProcessor::_passesJetCuts( ReconstructedParticle* pJet )
206 {
207  //
208  // This cut added on the suggestion of Sonja Hillert 12/Jan/07.
209  //
210  // Selects jets for which the cosine of the jet polar
211  // angle theta for all jets is not too large.
212  //
213  // Make something that's easy to search for to track down erroneous cuts:
214  // GREPTAG_CUT : Jet cut on abs(cos(theta)) of jet axis
215  //
216 
217  double CThJ_lower=0; //lower cut on cos(theta) of the jet axis
218  double CThJ_upper=0.95; //upper cut (Sho Ryu Ken!)
219 
220  vertex_lcfi::util::Vector3 zAxis( 0, 0, 1 ); //work out theta from dot product with jet axis
221 
222  const double* mom=pJet->getMomentum();
223  vertex_lcfi::util::Vector3 jetMomentum( mom[0], mom[1], mom[2] );
224 
225  jetMomentum.makeUnit();
226  double cosTheta=jetMomentum.dot( zAxis );
227  if( fabs(cosTheta)<=CThJ_lower || fabs(cosTheta)>=CThJ_upper ) return false;
228 
229 
230  // If control gets to this point then the jet has passed
231  return true;
232 }
233 
234 void DSTPlotProcessor::_fillPlots( LCEvent* pEvent, unsigned int jet)
235 {
236 
237  //jet the jet collection
238  LCCollection* pJetCollection=pEvent->getCollection( _JetCollectionName );
239 
240  //get the jet
241  ReconstructedParticle* pJet=dynamic_cast<ReconstructedParticle*>( pJetCollection->getElementAt(jet) );
242 
243  //get the PIDHandler fot this jet collection
244  PIDHandler pidh( pJetCollection ) ;
245 
246  //get the algorithm id for the FT info
247  int alid = pidh.getAlgorithmID("LCFIFlavourTag");
248 
249  //get the Particle id object containing the FT info
250  const ParticleID& pid = pidh.getParticleID(pJet,alid);
251 
252  //get the paramters for the FT info
253  FloatVec params = pid.getParameters();
254 
255  //get the actual flavour tags
256  double bTag= params[pidh.getParameterIndex(alid,"BTag")];
257  double cTag= params[pidh.getParameterIndex(alid,"CTag")];
258  double cTagBBack= params[pidh.getParameterIndex(alid,"BCTag")];
259 
260  //get the algorithm id for the MC info
261  int mcid = pidh.getAlgorithmID("MCTruth");
262 
263  //get the Particle id object containing the MC info
264  const ParticleID& mcpid = pidh.getParticleID(pJet,mcid);
265 
266  //get the paramters for the MC info
267  FloatVec mcparams = mcpid.getParameters();
268 
269  //get the true jet flavour
270  float jetType = mcparams[pidh.getParameterIndex(mcid,"TrueJetFlavour")];
271 
272 
273  //cout<<"JetType = "<<jetType<<endl;
274 
275  if( jetType==B_JET )
276  {
277  if( bTag > -10 ) _BTagEfficiencyPurity.add_signal( bTag );
278  if( cTag > -10 ) _CTagEfficiencyPurity.add_background( cTag );
279  if( cTagBBack > -10 ) _BCTagEfficiencyPurity.add_background( cTagBBack );
280  }
281  else if( jetType==C_JET )
282  {
283  if( bTag > -10 ) _BTagEfficiencyPurity.add_background( bTag );
284  if( cTag > -10 ) _CTagEfficiencyPurity.add_signal( cTag );
285  if( cTagBBack > -10 ) _BCTagEfficiencyPurity.add_signal( cTagBBack );
286  }
287  else
288  {
289  if( bTag > -10 ) _BTagEfficiencyPurity.add_background( bTag );
290  if( cTag > -10 ) _CTagEfficiencyPurity.add_background( cTag );
291  //don't add background for the bcnet result because we're only considering the b background
292  }
293 
294 }
295 
296 void DSTPlotProcessor::_outputDataToFile( string filename )
297 {
298 #ifdef USEROOT
299  stringstream filenameStream;
300  filenameStream << filename << ".root";
301 
302  TFile rootFile( filenameStream.str().c_str(), "RECREATE", "Various plots from the plot processor" );
303 
304  //have to put all the results into normal c arrays for root to deal with
305  //there may be an easier way of doing this conversion but I don't know what it is
306 
307 
308  int numberOfPoints=100;
309 
310  //get the results in STL form
311  vector< pair<double,double> > BTagefficiencyPurityPairs=_BTagEfficiencyPurity.eff_pur( numberOfPoints );
312  vector< pair<double,double> > CTagefficiencyPurityPairs=_CTagEfficiencyPurity.eff_pur( numberOfPoints );
313  vector< pair<double,double> > BCTagefficiencyPurityPairs=_BCTagEfficiencyPurity.eff_pur( numberOfPoints );
314 
315  //have to put all the results into normal c arrays for root to deal with
316  //there may be an easier way of doing this conversion but I don't know what it is
317  double BNetEfficiency[100];
318  double BNetPurity[100];
319  double CNetEfficiency[100];
320  double CNetPurity[100];
321  double BCNetEfficiency[100];
322  double BCNetPurity[100];
323 
324  for( int a=0; a<numberOfPoints; ++a )
325  {
326  BNetEfficiency[a]=BTagefficiencyPurityPairs[a].first;
327  BNetPurity[a]=BTagefficiencyPurityPairs[a].second;
328  CNetEfficiency[a]=CTagefficiencyPurityPairs[a].first;
329  CNetPurity[a]=CTagefficiencyPurityPairs[a].second;
330  BCNetEfficiency[a]=BCTagefficiencyPurityPairs[a].first;
331  BCNetPurity[a]=BCTagefficiencyPurityPairs[a].second;
332  }
333 
334  TGraph BNetGraph( 99, BNetEfficiency, BNetPurity );
335  BNetGraph.SetName( (" B Tag") );
336  TGraph CNetGraph( 99, CNetEfficiency, CNetPurity );
337  CNetGraph.SetName( ( " C Tag") );
338  TGraph BCNetGraph( 99, BCNetEfficiency, BCNetPurity );
339  BCNetGraph.SetName( ( " BC Tag") );
340 
341  TMultiGraph TagGraphs;
342  TagGraphs.SetName( "TagGraphs" );
343  TagGraphs.Add( &BNetGraph );
344  TagGraphs.Add( &CNetGraph );
345  TagGraphs.Add( &BCNetGraph );
346 
347  TAxis* pAxis;
348  TagGraphs.SetTitle( (string("Efficiency-purity for ") ).c_str() );
349 
350  //why the hell doesn't this work?
351  pAxis=TagGraphs.GetXaxis();
352  if(pAxis) pAxis->SetTitle( "Efficiency" );
353 
354  pAxis=TagGraphs.GetYaxis();
355  if(pAxis) pAxis->SetTitle( "Purity" );
356 
357  rootFile.Add( &TagGraphs );
358 
359  rootFile.Write();
360 
361 
362  // now add in the jet energies
363  TH1F jetEnergyHistogram( "jetEnergyHistogram", "Jet energies", 200, 0, 120 );
364  const vector<double> jetEnergyData=_jetEnergy.sorted_data();
365  for( vector<double>::const_iterator i=jetEnergyData.begin(); i<jetEnergyData.end(); ++i ) jetEnergyHistogram.Fill( (*i) );
366 
367  rootFile.Write();
368  rootFile.Close();
369 #endif //of
370 
371 }
372 
374 {
375  const vector<string>* pCollectionNames=pEvent->getCollectionNames();
376 
377  cout << "The available collections are: (name - type)" << endl;
378  for( vector<string>::const_iterator i=pCollectionNames->begin(); i<pCollectionNames->end(); ++i )
379  {
380  LCCollection* pCollection=pEvent->getCollection( (*i) );
381  const string typeName=pCollection->getTypeName();
382  cout << " " << (*i) << " - " << typeName << endl;
383  }
384  cout << endl;
385 }
386 
387 
388 void DSTPlotProcessor::_checkDSTParameters( LCEvent* pEvent)
389 {
390  try{
391  LCCollection* JetCollection=pEvent->getCollection( _JetCollectionName );
392 
393  try{
394  LCCollection* trueJetFlavourCollection = pEvent->getCollection(_TrueJetFlavourColName);
395 
396  try{
397  LCCollection* flavourTagInputsCollection = pEvent->getCollection( _FlavourTagInputsCollectionName);
398 
399  try{
400  LCCollection* TagCollection=pEvent->getCollection( _FlavourTagCollectionName );
401 
402 
403 
404  for( int jet=0; jet<JetCollection->getNumberOfElements(); jet++ )
405  {
406 
407 
408  //get all the info from the full collections
409  ReconstructedParticle* pJet = dynamic_cast<ReconstructedParticle*>( JetCollection->getElementAt(jet) );
410 
411  //get the flavour tags
412  LCFloatVec* flavourTags = dynamic_cast<LCFloatVec*>(TagCollection->getElementAt( jet ));
413  float fullBTag = (*flavourTags)[_IndexOfForEachTag["BTag"]];
414  float fullCTag = (*flavourTags)[_IndexOfForEachTag["CTag"]];
415  float fullBCTag =(*flavourTags)[_IndexOfForEachTag["BCTag"]];
416 
417  LCFloatVec* pJetFlavour = dynamic_cast<LCFloatVec*>(trueJetFlavourCollection->getElementAt( jet ));
418  float fullTruePDGCode = (*pJetFlavour)[_FlavourIndex["TruePDGCode"]];
419  float fullTruePartonCharge =(*pJetFlavour)[_FlavourIndex["TruePartonCharge"]];
420  float fullTrueHadronCharge= (*pJetFlavour)[_FlavourIndex["TrueHadronCharge"]];
421  float fullTrueJetFlavour = (*pJetFlavour)[_FlavourIndex["TrueJetFlavour"]];
422 
423  LCFloatVec* tagInputs = dynamic_cast<LCFloatVec*>(flavourTagInputsCollection->getElementAt( jet ));
424  int fullNumVertices = int((*tagInputs)[_InputsIndex["NumVertices"]]);
425 
426  float fullJointProbRPhi=999;
427  float fullJointProbZ=999;
428  float fullNumTracksInVertices=999;
429  float fullDecayLength=999;
430  float fullDecayLengthSignificance=999;
431  float fullRawMomentum=999;
432  float fullPTCorrectedMass=999;
433  float fullSecondaryVertexProbability=999;
434 
435  if(fullNumVertices > 1)
436  {
437  fullJointProbRPhi=(*tagInputs)[_InputsIndex["JointProbRPhi"]];
438  fullJointProbZ=(*tagInputs)[_InputsIndex["JointProbZ"]];
439  fullNumTracksInVertices=(*tagInputs)[_InputsIndex["NumTracksInVertices"]];
440  fullDecayLength=(*tagInputs)[_InputsIndex["DecayLength"]];
441  fullDecayLengthSignificance=(*tagInputs)[_InputsIndex["DecayLengthSignificance"]];
442  fullRawMomentum=(*tagInputs)[_InputsIndex["RawMomentum"]];
443  fullPTCorrectedMass=(*tagInputs)[_InputsIndex["PTCorrectedMass"]];
444  fullSecondaryVertexProbability=(*tagInputs)[_InputsIndex["SecondaryVertexProbability"]];
445  }
446 
447 
448  //get all the info from the DST collections
449  PIDHandler pidh( JetCollection ) ;
450  int alid = pidh.getAlgorithmID("LCFIFlavourTag");
451  const ParticleID& pid = pidh.getParticleID(pJet,alid);
452  FloatVec params = pid.getParameters();
453 
454  float dstBTag = params[pidh.getParameterIndex(alid,"BTag")];
455  float dstCTag = params[pidh.getParameterIndex(alid,"CTag")];
456  float dstBCTag =params[pidh.getParameterIndex(alid,"BCTag")];
457 
458  int mcid = pidh.getAlgorithmID("MCTruth");
459  const ParticleID& mcpid = pidh.getParticleID(pJet,mcid);
460  FloatVec mcparams = mcpid.getParameters();
461 
462  float dstTruePDGCode = mcparams[pidh.getParameterIndex(mcid,"TruePDGCode")];
463  float dstTruePartonCharge =mcparams[pidh.getParameterIndex(mcid,"TruePartonCharge")];
464  float dstTrueHadronChange= mcparams[pidh.getParameterIndex(mcid,"TrueHadronCharge")];
465  float dstTrueJetFlavour = mcparams[pidh.getParameterIndex(mcid,"TrueJetFlavour")];
466 
467  int dstNumVertices = int(params[pidh.getParameterIndex(alid,"NumVertices")]);
468 
469  float dstJointProbRPhi=999;
470  float dstJointProbZ=999;
471  float dstNumTracksInVertices=999;
472  float dstDecayLength=999;
473  float dstDecayLengthSignificance=999;
474  float dstRawMomentum=999;
475  float dstPTCorrectedMass=999;
476  float dstSecondaryVertexProbability=999;
477 
478  if(dstNumVertices > 1)
479  {
480  dstJointProbRPhi=params[pidh.getParameterIndex(alid,"JointProbRPhi")];
481  dstJointProbZ=params[pidh.getParameterIndex(alid,"JointProbZ")];
482  dstNumTracksInVertices=params[pidh.getParameterIndex(alid,"NumTracksInVertices")];
483  dstDecayLength=params[pidh.getParameterIndex(alid,"DecayLength")];
484  dstDecayLengthSignificance=params[pidh.getParameterIndex(alid,"DecayLengthSignificance")];
485  dstRawMomentum=params[pidh.getParameterIndex(alid,"RawMomentum")];
486  dstPTCorrectedMass=params[pidh.getParameterIndex(alid,"PTCorrectedMass")];
487  dstSecondaryVertexProbability=params[pidh.getParameterIndex(alid,"SecondaryVertexProbability")];
488  }
489 
490 
491  //compare the values
492  cout<< dstBTag <<" "<<fullBTag<<endl;;
493  cout<< dstCTag <<" "<<fullCTag<<endl;;
494  cout<< dstBCTag <<" "<<fullBCTag<<endl;;
495 
496 
497  cout<< dstTruePDGCode <<" "<<fullTruePDGCode<<endl;;
498  cout<< dstTruePartonCharge <<" "<<fullTruePartonCharge<<endl;;
499  cout<< dstTrueHadronChange<<" "<<fullTrueHadronCharge<<endl;;
500  cout<< dstTrueJetFlavour <<" "<<fullTrueJetFlavour<<endl;;
501  cout<< dstNumVertices <<" "<<fullNumVertices<<endl;;
502 
503  if(dstNumVertices > 1)
504  {
505  cout<< dstJointProbRPhi <<" "<<fullJointProbRPhi<<endl;;
506  cout<< dstJointProbZ <<" "<<fullJointProbZ<<endl;;
507  cout<< dstNumTracksInVertices <<" "<<fullNumTracksInVertices<<endl;;
508  cout<< dstDecayLength <<" "<<fullDecayLength<<endl;;
509  cout<< dstDecayLengthSignificance <<" "<<fullDecayLengthSignificance<<endl;;
510  cout<< dstRawMomentum <<" "<<fullRawMomentum<<endl;;
511  cout<< dstPTCorrectedMass <<" "<<fullPTCorrectedMass<<endl;;
512  cout<< dstSecondaryVertexProbability <<" "<<fullSecondaryVertexProbability<<endl;;
513  }
514 
515 
516  float ep = 0.00001;
517  if(fabs( dstBTag -fullBTag) > ep)
518  cout<<"B Parameters not equal"<<endl;
519  if(fabs( dstCTag -fullCTag) > ep)
520  cout<<"C Parameters not equal"<<endl;
521  if(fabs( dstBCTag -fullBCTag) > ep)
522  cout<<"BC Parameters not equal"<<endl;
523 
524 
525  if(fabs( dstTruePDGCode -fullTruePDGCode) > ep)
526  cout<<"TPDG Parameters not equal"<<endl;
527  if(fabs( dstTruePartonCharge -fullTruePartonCharge) > ep)
528  cout<<"TP Parameters not equal"<<endl;
529  if(fabs( dstTrueHadronChange-fullTrueHadronCharge) > ep)
530  cout<<"TH Parameters not equal"<<endl;
531  if(fabs( dstTrueJetFlavour -fullTrueJetFlavour) > ep)
532  cout<<"TJ Parameters not equal"<<endl;
533  if(fabs( dstNumVertices -fullNumVertices) > ep)
534  cout<<"NV Parameters not equal"<<endl;
535 
536  if(dstNumVertices > 1)
537  {
538  if(fabs( dstJointProbRPhi -fullJointProbRPhi) > ep)
539  cout<<"JPRP Parameters not equal"<<endl;
540  if(fabs( dstJointProbZ -fullJointProbZ) > ep)
541  cout<<"JPZ Parameters not equal"<<endl;
542  if(fabs( dstNumTracksInVertices -fullNumTracksInVertices) > ep)
543  cout<<"NT Parameters not equal"<<endl;
544  if(fabs( dstDecayLength -fullDecayLength) > ep)
545  cout<<"DL Parameters not equal"<<endl;
546  if(fabs( dstDecayLengthSignificance -fullDecayLengthSignificance) > ep)
547  cout<<"DLS Parameters not equal"<< dstDecayLengthSignificance -fullDecayLengthSignificance<<endl;
548  if(fabs( dstRawMomentum -fullRawMomentum) > ep)
549  cout<<"raw p Parameters not equal"<<endl;
550  if(fabs( dstPTCorrectedMass -fullPTCorrectedMass) > ep)
551  cout<<"pt cor Parameters not equal"<<endl;
552  if(fabs( dstSecondaryVertexProbability -fullSecondaryVertexProbability) > ep)
553  cout<<"sv Parameters not equal "<<endl;
554  }
555 
556 
557 
558  }
559 
560  }
561 
562  catch(DataNotAvailableException &e){
563  std::cout << "Collection " <<_JetCollectionName << " is unavailable " << std::endl;
564  }
565  }
566  catch(DataNotAvailableException &e){
567  std::cout << "Collection " <<_TrueJetFlavourColName << " is unavailable " << std::endl;
568  }
569  }
570  catch(DataNotAvailableException &e){
571  std::cout << "Collection " << _FlavourTagInputsCollectionName << " is unavailable " << std::endl;
572  }
573  }
574  catch(DataNotAvailableException &e){
575  std::cout << "Collection " <<_FlavourTagCollectionName << " is unavailable " << std::endl;
576  }
577 }
578 
579 
580 
void _fillPlots(LCEvent *pEvent, unsigned int jet)
bool _passesJetCuts(lcio::ReconstructedParticle *pJet)
static const int C_JET
std::string _JetCollectionName
std::string _TrueJetFlavourColName
std::string _OutputFilename
bool _passesEventCuts(lcio::LCEvent *pEvent)
histogram_data< double > _jetEnergy
void _outputDataToFile(std::string filename)
static const int B_JET
void _displayCollectionNames(lcio::LCEvent *pEvent)
Creates some sample plots from the data calculated by the LCFI vertex package.