LCFIVertex  0.7.2
VertexChargeProcessor.cc
1 #include "VertexChargeProcessor.h"
2 #include <iostream>
3 #include <sstream>
4 
5 #include <EVENT/LCParameters.h>
6 #include <EVENT/LCCollection.h>
7 #include <EVENT/ReconstructedParticle.h>
8 #include <EVENT/Track.h>
9 #include <EVENT/LCRelation.h>
10 #include <IMPL/LCCollectionVec.h>
11 #include <UTIL/LCRelationNavigator.h>
12 #include <IMPL/ParticleIDImpl.h>
13 #include <EVENT/LCFloatVec.h>
14 
15 #include <inc/event.h>
16 #include <util/inc/memorymanager.h>
17 #include <inc/jet.h>
18 #include <inc/track.h>
19 #include <inc/algo.h>
20 #include <inc/vertex.h>
21 #include <algo/inc/vertexcharge.h>
22 #include <algo/inc/trackattach.h>
23 #include <util/inc/helixrep.h>
24 #include <util/inc/projection.h>
25 #include <inc/track.h>
26 #include <inc/lciointerface.h>
27 #include <algo/inc/twotrackpid.h>
28 
29 #include <vector>
30 #include <string>
31 #include <map>
32 
33 using namespace marlin ;
34 using namespace lcio;
35 using namespace vertex_lcfi;
37 using std::vector;
38 using std::string;
39 
40 VertexChargeProcessor aVertexChargeProcessor ;
41 
42 VertexChargeProcessor::VertexChargeProcessor() : Processor("VertexChargeProcessor") {
43  // modify processor description
44  _description = "VertexChargeProcessor - takes a set of vertices as a decay chain with its associated jet and calculates the vertex charge." ;
45 
46  // register steering parameters: name, description, class-variable, default value
47  registerInputCollection( lcio::LCIO::RECONSTRUCTEDPARTICLE,
48  "JetRPCollection" ,
49  "Name of the ReconstructedParticle collection that represents jets",
50  _JetRPColName ,
51  std::string("Jets") ) ;
52  registerInputCollection( lcio::LCIO::VERTEX,
53  "IPVertexCollection" ,
54  "Name of the Vertex collection that contains the primary vertex (Optional)" ,
55  _IPVertexCollectionName ,
56  std::string("IPVertex") ) ;
57  registerInputCollection( lcio::LCIO::RECONSTRUCTEDPARTICLE,
58  "DecayChainRPCollection" ,
59  "Name of the ReconstructedParticle collection that represents decay chains" ,
60  _DecayChainRPColName ,
61  std::string("DecayChains") ) ;
62  registerOutputCollection( lcio::LCIO::LCFLOATVEC,
63  "VertexChargeCollection" ,
64  "Name of the LCFloatVec Collection that will be created to contain the flavour tag inputs" ,
65  _VertexChargeCollectionName,
66  "BCharge" ) ;
67  registerOptionalParameter( "ChargeAllSecondaryTracks",
68  "Parameter determining whether all tracks from secondary are included in the B-Charge" ,
69  _ChargeAddAllTracksFromSecondary,
70  bool(true));
71  registerOptionalParameter( "ChargeLoDCutmin",
72  "Cut determining the minimum L/D for the B-Charge" ,
73  _ChargeLoDCutmin,
74  double(0.18 )) ;
75  registerOptionalParameter( "ChargeLoDCutmax",
76  "Cut determining the maximum L/D for the B-Charge" ,
77  _ChargeLoDCutmax,
78  double(2.5)) ;
79  registerOptionalParameter( "ChargeCloseapproachCut",
80  "Upper cut on track distance of closest approach to the seed axis for the B-Charge " ,
81  _ChargeCloseapproachCut,
82  double(1.0)) ;
83 
84 }
85 
86 void VertexChargeProcessor::init()
87 {
88  // usually a good idea to
89  printParameters() ;
90 
91  _nRun = 0 ;
92  _nEvt = 0 ;
93 
94  _Attach = new TrackAttach();
96  _Attach->setDoubleParameter("AddAllTracksFromSecondary",_ChargeAddAllTracksFromSecondary );
97  _Attach->setDoubleParameter("LoDCutmin",_ChargeLoDCutmin );
98  _Attach->setDoubleParameter("LoDCutmax",_ChargeLoDCutmax );
99  _Attach->setDoubleParameter("CloseapproachCut",_ChargeCloseapproachCut );
100 
101  _VertexCharge = new VertexCharge();
103 
104 
105 }
106 
107 void VertexChargeProcessor::processRunHeader( LCRunHeader* run) {
108 
109  _JetVariableNames.clear();
110  _JetVariableNames.push_back("Charge");
111  run->parameters().setValues(_VertexChargeCollectionName, _JetVariableNames);
112 
113  _nRun++ ;
114 }
115 
116 void VertexChargeProcessor::processEvent( LCEvent * evt ) {
117 
118 
119  // this gets called for every event
120  // usually the working horse ...
121 if( isFirstEvent() )
122 {
123  std::cout << "Collections in Event" << std::endl;
124  const std::vector<std::string>* names = evt->getCollectionNames();
125  for (unsigned int i = 0 ; i<names->size() ; ++i) std::cout << (*names)[i] << std::endl;
126 }
127 
128  LCCollection* JetRPCol = evt->getCollection( _JetRPColName );
129  LCCollection* DecayChainRPCol = evt->getCollection( _DecayChainRPColName );
130 
131  //Find the primary vertex in the event
132  LCCollection* VertexCol;
133  VertexCol = evt->getCollection( _IPVertexCollectionName );
134 
135  //Search throught the vertices in this colection to find the primary
136  Vector3 IPPos;
137  Matrix3x3 IPErr;
138  int nVerts = VertexCol->getNumberOfElements() ;
139  bool done = 0;
140  for(int i=0; i< nVerts ; i++)
141  {
142  lcio::Vertex* iVertex = dynamic_cast<lcio::Vertex*>(VertexCol->getElementAt(i));
143  if (iVertex->isPrimary())
144  {
145  IPPos.x() = iVertex->getPosition()[0];
146  IPPos.y() = iVertex->getPosition()[1];
147  IPPos.z() = iVertex->getPosition()[2];
148  IPErr(0,0) = iVertex->getCovMatrix()[0];
149  IPErr(1,0) = iVertex->getCovMatrix()[1];
150  IPErr(1,1) = iVertex->getCovMatrix()[2];
151  IPErr(2,0) = iVertex->getCovMatrix()[3];
152  IPErr(2,1) = iVertex->getCovMatrix()[4];
153  IPErr(2,2) = iVertex->getCovMatrix()[5];
154  done = 1;
155  }
156  if (done) break;
157  }
158  if (!done) done = 1;//TODO Throw something
159 
160  Event* MyEvent = new Event(IPPos,IPErr);
161  MemoryManager<Event>::Event()->registerObject(MyEvent);
162 
163  std::map<Jet*,DecayChain*> DecayChainOf;
164  std::map<Jet*,ReconstructedParticle*> LCIORPOf;
165  //LCRelationNavigator RelNav(evt->getCollection(_RelationColName));
166  //std::cout << evt->getCollection(_RelationColName)->getNumberOfElements() << std::endl;
167  //Jet Loop, for each jet add it to the event and its corresponding decay chain to the map
168  int nRCP = JetRPCol->getNumberOfElements() ;
169 
170  if(nRCP ==0 ) std::cerr<<"Warning: VertexChargeProcessor.cc:336 : NO jets present "<<std::endl;
171 
172 
173  for(int i=0; i< nRCP ; i++)
174  {
175  ReconstructedParticle* JetRP = dynamic_cast<ReconstructedParticle*>(JetRPCol->getElementAt(i));
176  Jet* ThisJet = jetFromLCIORP(MyEvent,JetRP);
177  LCIORPOf[ThisJet] = JetRP;
178  //Assume Jets and DecayChains in same order in LCIO
179  DecayChainOf[ThisJet] = decayChainFromLCIORP(ThisJet,dynamic_cast<ReconstructedParticle*>(DecayChainRPCol->getElementAt(i)));
180  //Commented Out as we rely on the order of decay chains and jets being the same
181  /*//Find the Decay chain RP associated with this jet
182  std::cout << JetRPCol->getElementAt(i) <<std::endl;
183  LCObjectVec RelatedDecayChains = RelNav.getRelatedFromObjects(JetRPCol->getElementAt(i));
184 
185  LCRelation* rel = dynamic_cast<LCRelation*>( evt->getCollection(_RelationColName)->getElementAt(0) ) ;
186  std::cout << rel->getFrom() << " "<< rel->getTo() << " " << rel->getWeight() << std::endl;
187  LCRelation* rel2 = dynamic_cast<LCRelation*>( evt->getCollection(_RelationColName)->getElementAt(1) ) ;
188  std::cout << rel2->getFrom() << " "<< rel2->getTo() << " " << rel2->getWeight() << std::endl;
189  if (RelatedDecayChains.size() != 1)
190  {
191  std::cout << RelatedDecayChains.size() << " relations found!!" << std::endl;
192  //TODO Throw Something
193  }
194  */
195  }
196 
197  //Create the collection to store the result
198  LCCollectionVec* OutCollection = new LCCollectionVec("LCFloatVec");
199  evt->addCollection(OutCollection,_VertexChargeCollectionName);
200 
201 
202  //Loop over the jets
203  for (vector<Jet*>::const_iterator iJet=MyEvent->jets().begin();iJet != MyEvent->jets().end();++iJet)
204  {
205 
206  LCFloatVec VertexCharge;
207 
208 
209  VertexCharge.push_back(_VertexCharge->calculateFor(_Attach->calculateFor(DecayChainOf[*iJet])));
210  LCFloatVec* OutVec = new LCFloatVec(VertexCharge);
211  OutCollection->addElement(OutVec);
212 
213  }//End iJet Loop
214 
215  //std::cout << ",";std::cout.flush();
216  //Clear all objects created for this event
218  _nEvt ++ ;
219 }
220 
221 
222 
223 void VertexChargeProcessor::check( LCEvent * evt ) {
224  // nothing to check here - could be used to fill checkplots in reconstruction processor
225 }
226 
227 
228 void VertexChargeProcessor::end(){
229 
231  std::cout << "VertexChargeProcessor::end() " << name()
232  << " processed " << _nEvt << " events in " << _nRun << " runs "
233  << std::endl ;
234 
235 }
void registerObject(T *pointer)
Register an object for memory management.
static MetaMemoryManager * Event()
Returns the Event duration singleton instance of the controller.
const std::vector< Jet * > & jets() const
Get Jets.
Definition: event.cpp:35
void delAllObjects()
Delete all objects of all types held by this instance.
Calculates the Vertex Charge.
static MetaMemoryManager * Run()
Returns the Run duration singleton instance of the controller.
static MemoryManager< T > * Event()
Returns the Event duration singleton instance of the MemoryManager for type T.
Calculation of the charge of the vertex.
Track attachment algorithm. Calculates the seed vertex and the Tracks attached to it...