All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups
ZFinder.cc
Go to the documentation of this file.
1 #include "EVENT/LCIO.h"
2 #include "EVENT/LCRunHeader.h"
3 #include "EVENT/LCCollection.h"
4 #include "EVENT/LCParameters.h"
5 #include "EVENT/ReconstructedParticle.h"
6 #include "IMPL/ReconstructedParticleImpl.h"
7 #include "ZFinder.h"
8 //fg: for now MarlinReco should not depend on ROOT - use CLHEP instead
9 //#include "TLorentzVector.h"
10 //#include "TVector3.h"
11 #include "CLHEP/Vector/LorentzVector.h"
12 #include "CLHEP/Vector/ThreeVector.h"
13 typedef CLHEP::HepLorentzVector LorentzVector ;
14 typedef CLHEP::Hep3Vector Vector3D ;
15 
16 
17 // Marlin stuff
18 #include <marlin/Global.h>
19 // the event display
20 
21 #include <cstdlib>
22 
23 using namespace lcio;
24 
25 ZFinder aZFinder;
26 
27 ZFinder::ZFinder() : marlin::Processor("ZFinder") {
28 
29  registerProcessorParameter( "Printing" ,
30  "Print certain messages" ,
31  _printing,
32  (int)1 ) ;
33 
34 
35  std::string zdecay = "ee";
36  registerProcessorParameter( "ZDecay" ,
37  "Z decay mode" ,
38  _zdecay,
39  zdecay) ;
40 
41 
42  std::string inputParticleCollectionName = "PandoraPFOs";
43  registerInputCollection( LCIO::RECONSTRUCTEDPARTICLE,
44  "InputParticleCollectionName" ,
45  "Input Particle Collection Name " ,
46  _inputParticleCollectionName,
47  inputParticleCollectionName);
48 
49 
50  std::string outputParticleCollectionName = "eeX";
51  registerOutputCollection( LCIO::RECONSTRUCTEDPARTICLE,
52  "OutputParticleCollectionName" ,
53  "Output Particle Collection Name " ,
54  _outputParticleCollectionName,
55  outputParticleCollectionName);
56 
57 
58  registerProcessorParameter( "AddPhotons" ,
59  "Include photons with reconstructed Z" ,
60  _addPhotons,
61  (int)1 ) ;
62 
63  registerProcessorParameter( "CanUseClusterEnergyForElectrons" ,
64  "Allow possibilit of cluster not track used for electrons",
65  _canUseClusterEnergyForElectrons,
66  (int)1 );
67 
68 
69  registerProcessorParameter( "FermionMomentumCut" ,
70  "Momentum Cut for fermion from Z decay" ,
71  _momentumCut,
72  (float)10.) ;
73 
74  registerProcessorParameter( "CosTrackGammaCut" ,
75  "Minimum cosine of track-photon angle" ,
76  _cosTrackGammaCut,
77  (float)0.999);
78 
79  registerProcessorParameter( "MaxDeltaMz" ,
80  "Maximum difference between candidate and Z mass" ,
81  _dmzcut,
82  (float)50.);
83 
84  registerProcessorParameter( "MuonEcalEnergyCut" ,
85  "Cut on muon ECAL energy" ,
86  _muonEcalEnergyCut,
87  (float)2.5);
88 
89  registerProcessorParameter( "MuonHcalEnergyCut" ,
90  "Cut on Muon HCAL energy" ,
91  _muonHcalEnergyCut,
92  (float)10.);
93 
94 
95  registerProcessorParameter( "MuonHcalEnergyCut1" ,
96  "Cut on Muon HCAL energy (alt)" ,
97  _muonHcalEnergyCut1,
98  (float)5.);
99 
100  registerProcessorParameter( "ElectronEcalEnergyCut" ,
101  "Cut on electron ECAL energy" ,
102  _electronEcalEnergyCut,
103  (float)10.);
104 
105  registerProcessorParameter( "ElectronHcalEnergyCut" ,
106  "Cut on Electron HCAL energy" ,
107  _electronHcalEnergyCut,
108  (float)10.);
109 
110 
111  registerProcessorParameter( "ElectronEoverPCutLow" ,
112  "Cut on Electron E/p cut" ,
113  _electronEoPCutLow,
114  (float)0.5);
115 
116  registerProcessorParameter( "ElectronEoverPCutHigh" ,
117  "Cut on Electron E/p cut" ,
118  _electronEoPCutHigh,
119  (float)1.3);
120 
121 
122  return;
123 
124 }
125 
126 //===================================================================================
127 
128 void ZFinder::init() {
129  return;
130 }
131 
132 //===================================================================================
133 
134 void ZFinder::processRunHeader( LCRunHeader* /*run*/) {
135  return;
136 }
137 
138 //===================================================================================
139 
140 void ZFinder::processEvent( LCEvent * evt ) {
141 
142  // Make a new vector of particles
143  LCCollectionVec * recparcol = new LCCollectionVec(LCIO::RECONSTRUCTEDPARTICLE);
144 
145  // Access PFO collection
146  bool found = this->FindPFOs(evt);
147  if(found){
148  if(_printing>1)std::cout << "Analysis : " << _zdecay << std::endl;
149  if(_zdecay.find("ee")!=std::string::npos)this->FindZee(recparcol);
150  if(_zdecay.find("mm")!=std::string::npos)this->FindZmumu(recparcol);
151  if(_zdecay.find("mumu")!=std::string::npos)this->FindZmumu(recparcol);
152  if(_zdecay.find("tautau")!=std::string::npos)streamlog_out(ERROR) << "ZFinder tautau decay mode not implemented" << std::endl;
153  if(_zdecay.find("qq")!=std::string::npos)streamlog_out(ERROR) << "ZFinder qq decay mode not implemented" << std::endl;
154  }
155 
156  // Add new collection to event
157  evt->addCollection( recparcol , _outputParticleCollectionName.c_str() );
158 
159  return;
160 
161 }
162 
163 //===================================================================================
164 
165 void ZFinder::end(){
166  return;
167 }
168 
169 //===================================================================================
170 
171 bool ZFinder::FindPFOs( LCEvent* evt ) {
172 
173  bool tf = false;
174 
175  // clear old vector
176  _pfovec.clear();
177  typedef const std::vector<std::string> StringVec ;
178  StringVec* strVec = evt->getCollectionNames() ;
179  for(StringVec::const_iterator name=strVec->begin(); name!=strVec->end(); name++){
180  if(*name==_inputParticleCollectionName){
181  LCCollection* col = evt->getCollection(*name);
182  unsigned int nelem = col->getNumberOfElements();
183  tf = true;
184  for(unsigned int i=0;i<nelem;i++){
185  ReconstructedParticle* recoPart = dynamic_cast<ReconstructedParticle*>(col->getElementAt(i));
186  _pfovec.push_back(recoPart);
187  }
188  }
189  }
190 
191  if(_printing>1)std::cout << "Find PFOs : " << tf << std::endl;
192 
193  return tf;
194 
195 }
196 
197 //===================================================================================
198 
199 void ZFinder::FindZmumu(LCCollectionVec * recparcol) {
200 
201  if(_printing>1)std::cout << "FindZmumu : " << _pfovec.size() << std::endl;
202 
203  // Look for muon candidates (loose cuts)
204  std::vector<LorentzVector>pmuplus;
205  std::vector<LorentzVector>pmuminus;
206  std::vector<ReconstructedParticle*>pMplusPfoVec;
207  std::vector<ReconstructedParticle*>pMminusPfoVec;
208  for(unsigned int i=0;i<_pfovec.size();i++){
209  if(_pfovec[i]->getCharge()!=0){
210  ClusterVec c1 = _pfovec[i]->getClusters();
211  float ecal1 = 0;
212  float hcal1 = 0;
213  float muon1 = 0;
214  if(c1.size()==1){
215  ecal1 = c1[0]->getSubdetectorEnergies()[0] + c1[0]->getSubdetectorEnergies()[3];
216  hcal1 = c1[0]->getSubdetectorEnergies()[1] + c1[0]->getSubdetectorEnergies()[4];
217  muon1 = c1[0]->getSubdetectorEnergies()[2];
218  }
219  bool muon = true;
220  if(_pfovec[i]->getEnergy()<_momentumCut)muon=false;
221  if(ecal1>_muonEcalEnergyCut)muon=false;
222  if(hcal1>_muonHcalEnergyCut)muon=false;
223  if(hcal1>_muonHcalEnergyCut1&& muon1<0.01)muon=false;
224  if(muon){
225  LorentzVector pmu(_pfovec[i]->getMomentum()[0] ,_pfovec[i]->getMomentum()[1] ,_pfovec[i]->getMomentum()[2], _pfovec[i]->getEnergy() );
226  if(_pfovec[i]->getCharge()>0.5){
227  pmuplus.push_back(pmu);
228  pMplusPfoVec.push_back(_pfovec[i]);
229  }
230  if(_pfovec[i]->getCharge()<-0.5){
231  pmuminus.push_back(pmu);
232  pMminusPfoVec.push_back(_pfovec[i]);
233  }
234  }
235  }
236  }
237 
238  // loop over pairs of candidate muons and look for best Z candidates
239 
240  ReconstructedParticle* pMplus = NULL;
241  ReconstructedParticle* pMminus = NULL;
242  float bestdmz = 999.;
243  unsigned int besti=0;
244  unsigned int bestj=0;
245  for(unsigned int i=0;i<pmuplus.size();i++){
246  for(unsigned int j=0;j<pmuminus.size();j++){
247  LorentzVector pmumu = pmuplus[i]+pmuminus[j];
248  float massmumu = pmumu.m();
249  if( fabs(massmumu-91.2) < bestdmz){
250  bestdmz = fabs(massmumu-91.2);
251  besti=i; bestj=j;
252  pMplus = pMplusPfoVec[i];
253  pMminus = pMminusPfoVec[j];
254  }
255  }
256  }
257 
258  // Make the reconstructed particle and add any photons
259  float Energy = 0;
260  float Mom[3] = {0.,0.,0.};
261  if(bestdmz<_dmzcut && pMplus!=NULL && pMminus!=NULL){
262  ReconstructedParticleImpl * recoPart = new ReconstructedParticleImpl();
263  recoPart->addParticle(pMplus);
264  recoPart->addParticle(pMminus);
265  float pxem = pMminus->getMomentum()[0];
266  float pyem = pMminus->getMomentum()[1];
267  float pzem = pMminus->getMomentum()[2];
268  float eem = pMminus->getEnergy();
269  float pxep = pMplus->getMomentum()[0];
270  float pyep = pMplus->getMomentum()[1];
271  float pzep = pMplus->getMomentum()[2];
272  float eep = pMplus->getEnergy();
273  Energy = eep+eem;
274  Mom[0] = Mom[0] + pxem+pxep;
275  Mom[1] = Mom[1] + pyem+pyep;
276  Mom[2] = Mom[2] + pzem+pzep;
277 
278  // look for FSR photon candidates
279  if(_addPhotons>0){
280  for(unsigned int i=0;i<_pfovec.size();i++){
281  if(_pfovec[i]->getType()==22){
282  float pxg = _pfovec[i]->getMomentum()[0];
283  float pyg = _pfovec[i]->getMomentum()[1];
284  float pzg = _pfovec[i]->getMomentum()[2];
285  float eg = _pfovec[i]->getEnergy();
286  float cosm = 0.;
287  if(eg>0&&eem>0)cosm=(pxg*pxem+pyg*pyem+pzg*pzem)/eg/eem;
288  float cosp = 0.;
289  if(eg>0&&eep>0)cosp=(pxg*pxep+pyg*pyep+pzg*pzep)/eg/eep;
290  if(cosp>_cosTrackGammaCut||cosm>_cosTrackGammaCut){
291  recoPart->addParticle(_pfovec[i]);
292  Energy+=eg;
293  Mom[0]+=pxg;
294  Mom[1]+=pyg;
295  Mom[2]+=pzg;
296  }
297  }
298  }
299  }
300 
301  // set reconstructed particle parameters
302  float Mass = (Energy*Energy-Mom[0]*Mom[0]-Mom[1]*Mom[1]-Mom[2]*Mom[2]);
303  if(Mass>0)Mass=sqrt(Mass);
304  recoPart->setMomentum( Mom );
305  recoPart->setEnergy( Energy );
306  recoPart->setMass( Mass );
307  recoPart->setCharge( 0. );
308  recoPart->setType( 94 );
309  if(_printing>1)std::cout << "Found mmX : " << Mass << std::endl;
310 
311  // add it to the collection
312  recparcol->addElement( recoPart );
313  }
314 
315  return;
316 }
317 
318 //===================================================================================
319 
320 void ZFinder::FindZee(LCCollectionVec * recparcol) {
321 
322  if(_printing>1)std::cout << "FindZee : " << _pfovec.size() << std::endl;
323 
324  // Look for electron candidates (loose cuts)
325  std::vector<LorentzVector>peplus;
326  std::vector<LorentzVector>peminus;
327  std::vector<ReconstructedParticle*>pEplusPfoVec;
328  std::vector<ReconstructedParticle*>pEminusPfoVec;
329  for(unsigned int i=0;i<_pfovec.size();i++){
330  if(_pfovec[i]->getCharge()!=0 &&_pfovec[i]->getEnergy()>0){
331  ClusterVec c1 = _pfovec[i]->getClusters();
332  float ecal1 = 0;
333  float hcal1 = 0;
334  float muon1 = 0;
335  if(c1.size()==1){
336  ecal1 = c1[0]->getSubdetectorEnergies()[0] + c1[0]->getSubdetectorEnergies()[3];
337  hcal1 = c1[0]->getSubdetectorEnergies()[1] + c1[0]->getSubdetectorEnergies()[4];
338  muon1 = c1[0]->getSubdetectorEnergies()[2];
339  }
340  bool electron = true;
341  if(_pfovec[i]->getEnergy()<_momentumCut)electron=false;
342  if(ecal1<_electronEcalEnergyCut)electron=false;
343  if(hcal1>_electronHcalEnergyCut)electron=false;
344  float eop = ecal1/_pfovec[i]->getEnergy();
345  if(eop < _electronEoPCutLow)electron=false;
346  if(eop > _electronEoPCutHigh)electron=false;
347  if(muon1>0.01)electron=false;
348  if(abs(_pfovec[i]->getType())==11 && _pfovec[i]->getEnergy()>_momentumCut/2.0)electron=true;
349  if(electron){
350  LorentzVector pe(_pfovec[i]->getMomentum()[0] ,_pfovec[i]->getMomentum()[1] ,_pfovec[i]->getMomentum()[2], _pfovec[i]->getEnergy() );
351  if(_pfovec[i]->getCharge()>0.5){
352  peplus.push_back(pe);
353  pEplusPfoVec.push_back(_pfovec[i]);
354  }
355  if(_pfovec[i]->getCharge()<-0.5){
356  peminus.push_back(pe);
357  pEminusPfoVec.push_back(_pfovec[i]);
358  }
359  }
360  }
361  }
362 
363  // loop over pairs of candidate electrons and look for best Z candidates
364 
365  ReconstructedParticle* pEplus = NULL;
366  ReconstructedParticle* pEminus = NULL;
367  float bestdmz = 999.;
368  unsigned int besti = 999;
369  unsigned int bestj = 999;
370  for(unsigned int i=0;i<peplus.size();i++){
371  for(unsigned int j=0;j<peminus.size();j++){
372  LorentzVector pee = peplus[i]+peminus[j];
373  float massee = pee.m();
374  if( fabs(massee-91.2) < bestdmz){
375  bestdmz = fabs(massee-91.2);
376  besti=i; bestj=j;
377  pEplus = pEplusPfoVec[i];
378  pEminus = pEminusPfoVec[j];
379  }
380  }
381  }
382 
383  // Make the reconstructed particle and add any photons
384  float Energy = 0;
385  float Mom[3] = {0.,0.,0.};
386  if(bestdmz<_dmzcut && pEplus!=NULL && pEminus!=NULL){
387  ReconstructedParticleImpl * recoPart = new ReconstructedParticleImpl();
388  recoPart->addParticle(pEplus);
389  recoPart->addParticle(pEminus);
390  float pxem = pEminus->getMomentum()[0];
391  float pyem = pEminus->getMomentum()[1];
392  float pzem = pEminus->getMomentum()[2];
393  float eem = pEminus->getEnergy();
394  float pxep = pEplus->getMomentum()[0];
395  float pyep = pEplus->getMomentum()[1];
396  float pzep = pEplus->getMomentum()[2];
397  float eep = pEplus->getEnergy();
398 
399  // look for FSR and Bremstahlung photon candidates
400  std::vector<ReconstructedParticle*>photons;
401  float trackEPlus = 0;
402  float trackEMinus = 0;
403  float clusterEPlus = 0;
404  float clusterEMinus = 0;
405  float sigmaEPlus = 999.;
406  float sigmaEMinus = 999.;
407  float chiEPlus = 999.;
408  float chiEMinus = 999.;
409  float sigpPlus = 0.;
410  float sigpMinus = 0.;
411 
412 
413  Vector3D vecEPlus;
414  Vector3D vecEMinus;
415  const EVENT::ClusterVec ci = pEplus->getClusters();
416  const EVENT::TrackVec ti = pEplus->getTracks();
417  trackEPlus = pEplus->getEnergy();
418  if(ti.size()==1)sigpPlus = trackEPlus*sqrt(ti[0]->getCovMatrix()[5])/fabs(ti[0]->getOmega());
419  if(ci.size()>0){
420  clusterEPlus = ci[0]->getEnergy();
421  sigmaEPlus = 0.18*sqrt(clusterEPlus);
422  chiEPlus = (trackEPlus-clusterEPlus)/sqrt(sigmaEPlus*sigmaEPlus+sigpPlus*sigpPlus);
423  vecEPlus.set(ci[0]->getPosition()[0],ci[0]->getPosition()[1],ci[0]->getPosition()[2]);
424  }
425  const EVENT::ClusterVec cj = pEminus->getClusters();
426  const EVENT::TrackVec tj = pEminus->getTracks();
427  trackEMinus = pEminusPfoVec[bestj]->getEnergy();
428  if(tj.size()==1)sigpMinus = trackEMinus*sqrt(tj[0]->getCovMatrix()[5])/fabs(tj[0]->getOmega());
429  if(cj.size()>0){
430  clusterEMinus = cj[0]->getEnergy();
431  sigmaEMinus = 0.18*sqrt(clusterEMinus);
432  chiEMinus = (trackEMinus-clusterEMinus)/sqrt(sigmaEMinus*sigmaEMinus+sigpMinus*sigpMinus);
433  vecEMinus.set(cj[0]->getPosition()[0],cj[0]->getPosition()[1],cj[0]->getPosition()[2]);
434  }
435 
436  if(_addPhotons>0){
437  // loop over particles
438  for(unsigned int i=0;i<_pfovec.size();i++){
439  if(_pfovec[i]->getType()==22){
440  float pxg = _pfovec[i]->getMomentum()[0];
441  float pyg = _pfovec[i]->getMomentum()[1];
442  float pzg = _pfovec[i]->getMomentum()[2];
443  float eg = _pfovec[i]->getEnergy();
444  float cosm = 0.;
445  if(eg>0&&eem>0)cosm=(pxg*pxem+pyg*pyem+pzg*pzem)/eg/eem;
446  float cosp = 0.;
447  if(eg>0&&eep>0)cosp=(pxg*pxep+pyg*pyep+pzg*pzep)/eg/eep;
448  // found a photon candidate near to electron/positron
449  if(cosp>_cosTrackGammaCut||cosm>_cosTrackGammaCut){
450 
451  // check that this isn't a split-off cluster
452  float drp = 999.;
453  float drm = 999.;
454  const EVENT::ClusterVec c = _pfovec[i]->getClusters();
455  if(c.size()==1){
456  // distance between cluster and the electron/positron cluster
457  Vector3D vecg(c[0]->getPosition()[0],c[0]->getPosition()[1],c[0]->getPosition()[2]);
458  Vector3D v = vecg.cross(vecEPlus);
459  float magg = vecg.mag();
460  if(magg>0)drp = v.mag()/magg;
461  v = vecg.cross(vecEMinus);
462  if(magg>0)drm = v.mag()/vecg.mag();
463  }
464  if(cosp>cosm){
465  bool merge = false;
466  float chiNew = (trackEPlus-clusterEPlus-eg)/sigmaEPlus;
467  if(drp<20.0)merge = true;
468  // if fairly close merge if chi2 for matching improves greatly
469  if(drp<30.0 && chiEPlus>4.0 && fabs(chiNew)<chiEPlus)merge = true;
470  if(drp<40.0 && chiEPlus>5.0 && fabs(chiNew)<chiEPlus)merge = true;
471  if(drp<50.0 && chiEPlus>7.0 && fabs(chiNew)<chiEPlus)merge = true;
472  // sanity check
473  if( fabs(chiEPlus)<2.0 && chiNew*chiNew>chiEPlus*chiEPlus+5.0)merge = false;
474  // always merge if very close - can't expect reconstruction to work
475  if(drp<10.0)merge = true;
476  // if we merge cluster do nothing except add photon to Z
477  if(merge){
478  recoPart->addParticle(_pfovec[i]);
479  clusterEPlus += eg;
480  sigmaEPlus = 0.18*sqrt(clusterEPlus);
481  chiEPlus = (trackEPlus-clusterEPlus)/sqrt(sigmaEPlus*sigmaEPlus+sigpPlus*sigpPlus);
482  }
483  // if not merged add photon and energy to Z
484  if(!merge){
485  recoPart->addParticle(_pfovec[i]);
486  photons.push_back(_pfovec[i]);
487  }
488  }else{
489  float chiNew = (trackEMinus-clusterEMinus-eg)/sigmaEMinus;
490  bool merge = false;
491  // always merge if very close
492  if(drm<20.0)merge = true;
493  // if fairly close merge if chi2 for matching improves greatly
494  if(drm<30.0 && chiEMinus>4.0 && fabs(chiNew)<chiEMinus)merge = true;
495  if(drm<40.0 && chiEMinus>5.0 && fabs(chiNew)<chiEMinus)merge = true;
496  if(drm<50.0 && chiEMinus>7.0 && fabs(chiNew)<chiEMinus)merge = true;
497  // sanity check
498  if( fabs(chiEMinus)<2.0 && chiNew*chiNew>chiEMinus*chiEMinus+5.0)merge = false;
499  // always merge if very close - can't expect reconstruction to work
500  if(drm<10.0)merge = true;
501  // if we merge cluster do nothing
502  if(merge){
503  recoPart->addParticle(_pfovec[i]);
504  clusterEMinus += eg;
505  sigmaEMinus = 0.18*sqrt(clusterEMinus);
506  chiEMinus = (trackEMinus-clusterEMinus)/sqrt(sigmaEMinus*sigmaEMinus+sigpMinus*sigpMinus);
507  }
508  // if not merged then add photon to Z
509  if(!merge){
510  recoPart->addParticle(_pfovec[i]);
511  photons.push_back(_pfovec[i]);
512  }
513  }
514  }
515  }
516  } // end of loop over particles
517  } // end of adding photons
518 
519 
520  // first the electrons and positrons
521  Energy = eep+eem;
522  Mom[0] = Mom[0] + pxem+pxep;
523  Mom[1] = Mom[1] + pyem+pyep;
524  Mom[2] = Mom[2] + pzem+pzep;
525  // any photons
526  for(unsigned int i=0;i<photons.size();i++){
527  Energy += photons[i]->getEnergy();
528  Mom[0] += photons[i]->getMomentum()[0];
529  Mom[1] += photons[i]->getMomentum()[1];
530  Mom[2] += photons[i]->getMomentum()[2];
531  }
532  float Mass = sqrt(Energy*Energy-Mom[0]*Mom[0]-Mom[1]*Mom[1]-Mom[2]*Mom[2]);
533 
534 
535  // Use cluster for electron/positron energy ?
536  if(_canUseClusterEnergyForElectrons && (chiEPlus<-3.5 || chiEMinus < -3.5)){
537  // look for evidence of merged electrons/photons
538  float MomX[3]={0.,0.,0.};
539  float EnergyX=0;
540  float scaleP = 1.0;
541  if(chiEPlus<-3.5)scaleP = clusterEPlus/trackEPlus;
542  float scaleM = 1.0;
543  if(chiEMinus<-3.5)scaleM = clusterEMinus/trackEMinus;
544  EnergyX = Energy + eep*(scaleP-1.) +eem*(scaleM-1.);
545  MomX[0] = Mom[0] + pxem*(scaleM-1.)+pxep*(scaleP-1.);
546  MomX[1] = Mom[1] + pyem*(scaleM-1.)+pyep*(scaleP-1.);
547  MomX[2] = Mom[2] + pzem*(scaleM-1.)+pzep*(scaleP-1.);
548  float MassX = sqrt(EnergyX*EnergyX-MomX[0]*MomX[0]-MomX[1]*MomX[1]-MomX[2]*MomX[2]);
549  if(fabs(MassX-91.2)<fabs(Mass-91.2)){
550  Energy = EnergyX;
551  Mom[0] = MomX[0];
552  Mom[1] = MomX[1];
553  Mom[2] = MomX[2];
554  Mass = sqrt(Energy*Energy-Mom[0]*Mom[0]-Mom[1]*Mom[1]-Mom[2]*Mom[2]);
555  }
556  }
557 
558  // Set parameters of reconstructed Z
559  if(Mass>0)Mass=sqrt(Mass);
560  recoPart->setMomentum( Mom );
561  recoPart->setEnergy( Energy );
562  recoPart->setMass( Mass );
563  recoPart->setCharge( 0. );
564  recoPart->setType( 94 );
565  if(_printing>1)std::cout << "Found eeX : " << Mass << std::endl;
566 
567  // add it to the collection
568  recparcol->addElement( recoPart );
569  }
570 
571  return;
572 }
573 
ZFinder::_cosTrackGammaCut
float _cosTrackGammaCut
Definition: ZFinder.h:67
ZFinder::_electronEoPCutLow
float _electronEoPCutLow
Definition: ZFinder.h:62
ZFinder::_pfovec
std::vector< ReconstructedParticle * > _pfovec
Definition: ZFinder.h:51
ZFinder::FindZee
void FindZee(LCCollectionVec *)
Definition: ZFinder.cc:320
ZFinder::_zdecay
std::string _zdecay
Definition: ZFinder.h:53
Vector3D
CLHEP::Hep3Vector Vector3D
Definition: DistilledPFOCreator.cc:11
ZFinder::end
virtual void end()
Called after data processing for clean up.
Definition: ZFinder.cc:165
ZFinder::_electronHcalEnergyCut
float _electronHcalEnergyCut
Definition: ZFinder.h:61
ZFinder::_inputParticleCollectionName
std::string _inputParticleCollectionName
Definition: ZFinder.h:54
ZFinder::init
virtual void init()
Called at the begin of the job before anything is read.
Definition: ZFinder.cc:128
ZFinder::FindZmumu
void FindZmumu(LCCollectionVec *)
Definition: ZFinder.cc:199
ZFinder::processRunHeader
virtual void processRunHeader(LCRunHeader *run)
Called for every run.
Definition: ZFinder.cc:134
ZFinder::_momentumCut
float _momentumCut
Definition: ZFinder.h:56
ZFinder::_printing
int _printing
Definition: ZFinder.h:52
ZFinder::_electronEcalEnergyCut
float _electronEcalEnergyCut
Definition: ZFinder.h:60
ZFinder::_electronEoPCutHigh
float _electronEoPCutHigh
Definition: ZFinder.h:63
ZFinder::_canUseClusterEnergyForElectrons
int _canUseClusterEnergyForElectrons
Definition: ZFinder.h:66
LorentzVector
CLHEP::HepLorentzVector LorentzVector
Definition: DistilledPFOCreator.cc:10
ZFinder::processEvent
virtual void processEvent(LCEvent *evt)
Called for every event - the working horse.
Definition: ZFinder.cc:140
sistrip::LCCollectionVec
std::vector< LCCollection * > LCCollectionVec
Definition: SiStripClus.h:55
aZFinder
ZFinder aZFinder
Definition: ZFinder.cc:25
ZFinder::_muonHcalEnergyCut
float _muonHcalEnergyCut
Definition: ZFinder.h:58
ZFinder::_dmzcut
float _dmzcut
Definition: ZFinder.h:64
ZFinder::ZFinder
ZFinder()
Definition: ZFinder.cc:27
ZFinder::_outputParticleCollectionName
std::string _outputParticleCollectionName
Definition: ZFinder.h:55
ZFinder::_muonEcalEnergyCut
float _muonEcalEnergyCut
Definition: ZFinder.h:57
ZFinder::_addPhotons
int _addPhotons
Definition: ZFinder.h:65
ZFinder::_muonHcalEnergyCut1
float _muonHcalEnergyCut1
Definition: ZFinder.h:59
ZFinder::FindPFOs
bool FindPFOs(LCEvent *evt)
Definition: ZFinder.cc:171
ZFinder
ZFinder: Returns the best Z-&gt;ee/Z-&gt;mm candidate in the event.
Definition: ZFinder.h:19
sistrip::StringVec
std::vector< std::string > StringVec
Definition: SiStripClus.h:56