LCFIVertex  0.7.2
DSTCollectionProcessor.cc
1 #include "DSTCollectionProcessor.h"
2 #include <iostream>
3 #include <fstream>
4 #include <sstream>
5 #include <cmath>
6 
7 #include <EVENT/LCParameters.h>
8 #include <EVENT/LCCollection.h>
9 #include <EVENT/ReconstructedParticle.h>
10 #include <EVENT/Track.h>
11 #include <EVENT/LCRelation.h>
12 #include <EVENT/MCParticle.h>
13 #include <EVENT/LCIntVec.h>
14 #include "EVENT/LCFloatVec.h"
15 
16 
17 #include <UTIL/LCRelationNavigator.h>
18 #include <UTIL/PIDHandler.h>
19 #include <IMPL/ParticleIDImpl.h>
20 #include <IMPL/LCCollectionVec.h>
21 
22 #include <vector>
23 #include <string>
24 #include <map>
25 #include <set>
26 
27 using namespace marlin ;
28 using namespace lcio;
29 using namespace std;
30 
31 using std::vector;
32 using std::string;
33 using std::map;
34 using EVENT::Track;
35 
36 
37 DSTCollectionProcessor aDSTCollectionProcessor ;
38 
39 DSTCollectionProcessor::DSTCollectionProcessor() : Processor("DSTCollectionProcessor") {
40 
41  // modify processor description
42  _description = "DSTCollectionProcessor - Takes the flavour tag info, mc truth info, and some NN input info, and adds it as Particle ID info to the jets, for the DSTs" ;
43 
44 
45  // register steering parameters: name, description, class-variable, default value
46 
47 
48 
49  registerInputCollection( LCIO::RECONSTRUCTEDPARTICLE,
50  "JetCollectionName" ,
51  "Name of the collection of ReconstructedParticles that is the jet" ,
53  std::string("Jets") );
54 
55 
56  registerInputCollection( LCIO::LCFLOATVEC,
57  "FlavourTagCollection" ,
58  "Name of the LCFloatVec Collection containing the flavour tags (in same order as jet collection)" ,
59  _FlavourTagCollectionName,
60  "FlavourTag" ) ;
61 
62  registerInputCollection( LCIO::LCFLOATVEC,
63  "FlavourTagInputsCollection" ,
64  "Name of the LCFloatVec Collection that contains the flavour tag inputs (in same order as jet collection)" ,
65  _FlavourTagInputsCollectionName,
66  "FlavourTagInputs" ) ;
67 
68 
69 
70  registerInputCollection( LCIO::LCFLOATVEC,
71  "TrueJetFlavourCollection" ,
72  "Name of the LCIntVec collection containing the true flavour of the jets (same order as jets)" ,
73  _TrueJetFlavourColName ,
74  std::string("TrueJetFlavour") ) ;
75 
76  registerProcessorParameter( "Debug",
77  "Debugging option",
78  _debug,
79  int(0));
80 
81 
82  }
83 
84 void DSTCollectionProcessor::init()
85 {
86  // usually a good idea to
87  printParameters() ;
88 
89  _nRun = 0 ;
90  _nEvt = 0 ;
91 }
92 
93 void DSTCollectionProcessor::processRunHeader( LCRunHeader* pRun) {
94 
95 
96  // Marlin doesn't necessarily process the run header, e.g. if you use the "SkipNEvents"
97  // parameter in the steering file. The flavour tag variable/tag value names are held in
98  // the run header though, so this processor has to have access to it. Set this variable
99  // so that "processEvent" can tell if "processRunHeader" has been called for the run
100  // it's in.
101  _lastRunHeaderProcessed=pRun->getRunNumber();
102 
103 
104  //map the parameter names to the indexes, for the flavout tags, NN input vales, and mc truth info
105 
106  std::vector<std::string> VarNames;
107  (pRun->parameters()).getStringVals(_FlavourTagCollectionName,VarNames);
108  for (size_t i = 0;i < VarNames.size();++i)
109  _IndexOfForEachTag[VarNames[i]] = i;
110 
111 
112  //do the same for the true jet values
113  std::vector<std::string> trueJetFlavourVarNames;
114  (pRun->parameters()).getStringVals(_TrueJetFlavourColName,trueJetFlavourVarNames);
115  for (size_t i = 0;i < trueJetFlavourVarNames.size();++i)
116  _FlavourIndex[trueJetFlavourVarNames[i]] = i;
117 
118 
119  //do the same for the FT inputs
120  std::vector<std::string> InputVarNames;
121  (pRun->parameters()).getStringVals(_FlavourTagInputsCollectionName,InputVarNames);
122  for (size_t i = 0;i < InputVarNames.size();++i)
123  _InputsIndex[InputVarNames[i]] = i;
124 
125 
126 
127  _nRun++ ;
128 
129 }
130 
131 
132 void DSTCollectionProcessor::processEvent( LCEvent* pEvent) {
133 
134  try{
135  LCCollection* JetCollection=pEvent->getCollection( _JetCollectionName );
136 
137  try{
138  LCCollection* trueJetFlavourCollection = pEvent->getCollection(_TrueJetFlavourColName);
139 
140  try{
141  LCCollection* flavourTagInputsCollection = pEvent->getCollection( _FlavourTagInputsCollectionName);
142 
143  try{
144  LCCollection* TagCollection=pEvent->getCollection( _FlavourTagCollectionName );
145 
146 
147  //make a PIDHandler to add the info with
148  PIDHandler jetPID ( JetCollection ) ;
149 
150  //the string vec containing the parameter names for FT info
151  StringVec pNames ;
152  pNames.push_back( "BTag" );
153  pNames.push_back( "CTag" );
154  pNames.push_back( "BCTag" ) ;
155 
156  pNames.push_back("NumVertices");
157  pNames.push_back("JointProbRPhi");
158  pNames.push_back("JointProbZ");
159  pNames.push_back("NumTracksInVertices");
160  pNames.push_back("DecayLength");
161  pNames.push_back("DecayLengthSignificance");
162  pNames.push_back("RawMomentum");
163  pNames.push_back("PTCorrectedMass");
164  pNames.push_back("SecondaryVertexProbability");
165 
166 
167  //the string vec containing the parameter names for MC info
168  StringVec mcNames ;
169 
170  mcNames.push_back("TruePDGCode");
171  mcNames.push_back("TruePartonCharge");
172  mcNames.push_back("TrueHadronCharge");
173  mcNames.push_back("TrueJetFlavour");
174 
175  //make two new PID algorithms, one for the FT info, one for the MC info
176  int ftagID = jetPID.addAlgorithm( "LCFIFlavourTag",pNames ) ;
177  int mcID = jetPID.addAlgorithm( "MCTruth",mcNames ) ;
178 
179  //fill the float vec that will contain the paramters with dummy values = 9999
180  FloatVec fv(pNames.size());
181  //fill with dummy values
182  for(int i = 0; i < (int)pNames.size(); i++)
183  fv[i] = 9999;
184 
185  FloatVec mcv(mcNames.size());
186  //fill with dummy values
187  for(int i = 0; i < (int)mcNames.size(); i++)
188  mcv[i] = 9999;
189 
190  if (_debug>0)
191  cout<<"Number of jets in the event "<<JetCollection->getNumberOfElements()<<endl;
192 
193  //loop over all jets
194  for( int jet=0; jet<JetCollection->getNumberOfElements(); jet++ )
195  {
196 
197 
198  ReconstructedParticle* Jet = dynamic_cast<ReconstructedParticle*>( JetCollection->getElementAt(jet) );
199 
200  //get the FT info, add it to the parameter float vec
201  LCFloatVec* flavourTags = dynamic_cast<LCFloatVec*>(TagCollection->getElementAt( jet ));
202  fv[0]= (*flavourTags)[_IndexOfForEachTag["BTag"]];
203  fv[1]= (*flavourTags)[_IndexOfForEachTag["CTag"]];
204  fv[2]=(*flavourTags)[_IndexOfForEachTag["BCTag"]];
205 
206  //get the NN input info, add it to the parameter float vec
207  LCFloatVec* tagInputs = dynamic_cast<LCFloatVec*>(flavourTagInputsCollection->getElementAt( jet ));
208  int NumVertices = int((*tagInputs)[_InputsIndex["NumVertices"]]);
209  fv[3]= NumVertices;
210  //only add this info if the no of vertices >=2
211  if(NumVertices > 1)
212  {
213  fv[4]=(*tagInputs)[_InputsIndex["JointProbRPhi"]];
214  fv[5]=(*tagInputs)[_InputsIndex["JointProbZ"]];
215  fv[6]=(*tagInputs)[_InputsIndex["NumTracksInVertices"]];
216  fv[7]=(*tagInputs)[_InputsIndex["DecayLength"]];
217  fv[8]=(*tagInputs)[_InputsIndex["DecayLengthSignificance"]];
218  fv[9]=(*tagInputs)[_InputsIndex["RawMomentum"]];
219  fv[10]=(*tagInputs)[_InputsIndex["PTCorrectedMass"]];
220  fv[11]=(*tagInputs)[_InputsIndex["SecondaryVertexProbability"]];
221  }
222 
223 
224  //add the FT particle id to the jet collection
225  jetPID.setParticleID( Jet , 42 , // user type
226  99999 , // PDG
227  0.0, // likelihood
228  ftagID ,
229  fv ) ;
230 
231  //get the mc truth info, add it to the mc parameter float vec
232  LCFloatVec* pJetFlavour = dynamic_cast<LCFloatVec*>(trueJetFlavourCollection->getElementAt( jet ));
233  if(pJetFlavour==0) std::cerr << "The wrong type of true jet flavour collection was found, dynamic cast failed" << std::endl;
234  mcv[0] = (*pJetFlavour)[_FlavourIndex["TruePDGCode"]];
235  mcv[1] =(*pJetFlavour)[_FlavourIndex["TruePartonCharge"]];
236  mcv[2]= (*pJetFlavour)[_FlavourIndex["TrueHadronCharge"]];
237  mcv[3]= (*pJetFlavour)[_FlavourIndex["TrueJetFlavour"]];
238 
239  //add the MC particle id to the jet collection
240  jetPID.setParticleID( Jet , 43 , // user type
241  99998 , // PDG
242  0.0, // likelihood
243  mcID ,
244  mcv ) ;
245  }
246 
247 
248  }
249  catch(DataNotAvailableException &e){
250  if (_debug == 1)
251  std::cout << "Collection " <<_FlavourTagCollectionName << " is unavailable in event " << _nEvt << std::endl;
252  }
253  }
254  catch(DataNotAvailableException &e){
255  if (_debug == 1)
256  std::cout << "Collection " <<_FlavourTagInputsCollectionName << " is unavailable in event " << _nEvt << std::endl;
257  }
258  }
259  catch(DataNotAvailableException &e){
260  if (_debug == 1)
261  std::cout << "Collection " << _TrueJetFlavourColName << " is unavailable in event " << _nEvt << std::endl;
262  }
263  }
264  catch(DataNotAvailableException &e){
265  if (_debug == 1)
266  std::cout << "Collection " <<_JetCollectionName << " is unavailable in event " << _nEvt << std::endl;
267  }
268 
269 
270 
271 
272  _nEvt ++ ;
273 
274 }
275 
276 
277 
278 
279 void DSTCollectionProcessor::check( LCEvent * evt ) {
280  // nothing to check here - could be used to fill checkplots in reconstruction processor
281 }
282 
283 
284 void DSTCollectionProcessor::end(){
285 
286  std::cout << "DSTCollectionProcessor::end() " << name()
287  << " processed " << _nEvt << " events in " << _nRun << " runs "
288  << std::endl ;
289 
290 }
291