All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups
EMShowerFinder.cc
Go to the documentation of this file.
1 #include "EMShowerFinder.h"
2 
3 using namespace lcio ;
4 using namespace marlin ;
5 
6 EMShowerFinder aEMShowerFinder ;
7 
8 
9 
10 EMShowerFinder::EMShowerFinder() : Processor("EMShowerFinder") {
11 
12  // modify processor description
13  _description = "a photon finder processor based on the KIT and KITutil 'classes and fuctions'" ;
14 
15  registerProcessorParameter( "colNameECAL" ,
16  "ECAL Collection Name" ,
17  _colNameECAL,
18  std::string("ECAL") );
19 
20  registerProcessorParameter( "collectionNameOfEMShowerCandidates",
21  "Name of the collection of EM shower candidates",
22  _collectionNameOfEMShowerCandidates,
23  std::string("EMShowerCandidates"));
24 
25  registerProcessorParameter( "Cleaning",
26  "To do the cleaning on hits or not ",
27  _ToClean,
28  std::string("YES"));
29  registerProcessorParameter( "TopologicalCut",
30  "At which number of neighbors to put the threshold, condition is < so you need to put N+1 ",
31  _CleanCut,
32  (int)5);
33 
34  registerProcessorParameter( "NumberOfLevels",
35  "Number of levels for central loop ",
36  _N,
37  (int)10);
38  vector<float> miipstep;
39  miipstep.push_back(0.1);
40  miipstep.push_back(1.5);
41  miipstep.push_back(2.5);
42  miipstep.push_back(4.0);
43  miipstep.push_back(6.0);
44  miipstep.push_back(9.0);
45  miipstep.push_back(16.0);
46  miipstep.push_back(26.0);
47  miipstep.push_back(41.0);
48  miipstep.push_back(65.0);
49 
50  registerProcessorParameter( "Levels",
51  "Levels for central loop in MIP ",
52  _miipstep,
53  miipstep);
54 
55  registerProcessorParameter( "MinHit0",
56  "Minimal Number of hits for ground level cluster ",
57  _MinHit0,
58  (int)4);
59 
60  registerProcessorParameter( "MinHitSplit",
61  "Minimal Number of hits for i-th level cluster ",
62  _MinHitSplit,
63  (int)2);
64 
65  registerProcessorParameter( "Rcut",
66  "Fluctuation suprresion cut",
67  _Rcut,
68  (double)0.4);
69  registerProcessorParameter( "Distcut",
70  "Square of distance cut for merging ",
71  _Distcut,
72  (double)35.0);
73  registerProcessorParameter( "Coscut",
74  "Cosine of the angle for merging ",
75  _Coscut,
76  (double)0.95);
77 
78  registerProcessorParameter( "energyDeviationCut",
79  "cut on energy deviation of em shower candidates from estimated energy ( abs( (Ecluster - Eestimated)/Eestimated ) < cut )",
80  _energyDeviationCut,
81  (double)0.25);
82 
83  registerProcessorParameter( "probabilityDensityCut",
84  "cut on the probability density to assign hits to shower cores",
85  _probabilityDensityCut,
86  (double)1E-3);
87 
88  registerProcessorParameter( "DebugLevel",
89  "limits the amount of information written to std out (0 - none, 9 - maximal information)",
90  _debugLevel,
91  int(0) );
92 
93  registerProcessorParameter( "DrawOnCED",
94  "draw objects on CED",
95  _drawOnCED,
96  int(0) );
97 
98 
99 }
100 
101 
102 void EMShowerFinder::init() {
103 
104  // usually a good idea to
105  printParameters();
106 
107  // FIXME: hard coded cell id's for old Mokka (e.g. Mokka v5.4) versions)
108  CellIDDecoder<CalorimeterHit>::setDefaultEncoding("M:3,S-1:3,I:9,J:9,K-1:6");
109 
110 
111  // debug
112  if (_drawOnCED) MarlinCED::init(this);
113 
114  _nRun = 0 ;
115  _nEvt = 0 ;
116 
117 
118 }
119 
120 
121 void EMShowerFinder::processRunHeader( LCRunHeader* /*run*/) {
122 
123  ++_nRun;
124 
125 }
126 
127 
128 void EMShowerFinder::processEvent( LCEvent * evt ) {
129 
130 
131  // debug
132  if (_drawOnCED) {
133 
134  MarlinCED::newEvent(this,0);
135  MarlinCED::drawMCParticleTree(evt,"MCParticle",0.05,4.0,15.5,50.0,1626.0,2500.0);
136 
137  }
138 
139 
140  try {
141 
142  LCCollection* colt = evt->getCollection(_colNameECAL.c_str()) ;
143  CellIDDecoder<CalorimeterHit> CDECAL(colt);
144 
145  LCCollectionVec* clscol = new LCCollectionVec(LCIO::CLUSTER);
146 
147  if( colt!=0 ) {
148 
149  unsigned int nelem=colt->getNumberOfElements();
150  // MAIN CONTAINER OF SHITS
151  vector<Superhit2*> calo[10];
152 
153  // creating all superhits
154  CreateAllShits2(colt,CDECAL,calo);
155 
156  // precalculation
157  TotalPrecalc2(calo,CDECAL,nelem,_CleanCut);
158 
159  // setting the parameters of the alghorithm
160  vector <PROTSEED2> prs2;
161  CoreCut2 Ccut;
162  Ccut.Rcut=_Rcut;
163  Ccut.Distcut=_Distcut;
164  Ccut.Coscut=_Coscut;
165  Ccut.MinHit0=(unsigned int) _MinHit0;
166  Ccut.MinHitSplit=(unsigned int) _MinHitSplit;
167  const unsigned int N=_N;
168  vector< vector<Tmpcl2*> > bbb(N);
169 
170  // finding cores.
171  if( _ToClean=="YES" || _ToClean=="yes")
172  {
173  // cout << " da koristim cat " << endl;
174  FindCores2(&(calo[4]), bbb , &prs2,_N,_miipstep,Ccut);
175  }else{
176  FindCores2(&(calo[0]), bbb , &prs2,_N,_miipstep,Ccut);
177  }
178 
179 
180  // container to store information to assign hits to photon candidates
181  std::vector<ECALHitWithAttributes> ECALHitsWithAttributes;
182 
183  std::vector<CoreCalib2> coreCalibrationLDC00;
184  CreateCalibrationLDC00(&coreCalibrationLDC00);
185 
186 
187  for(unsigned int i=0;i<prs2.size();i++)
188  {
189  if(prs2[i].active==true)
190  {
191 
192  double centerPosition[3];
193  centerPosition[0]=(double)prs2[i].cl->getCenter()[0];
194  centerPosition[1]=(double)prs2[i].cl->getCenter()[1];
195  centerPosition[2]=(double)prs2[i].cl->getCenter()[2];
196 
197  double directionOfCluster[3];
198  prs2[i].cl->findInertia();
199  directionOfCluster[0] = prs2[i].cl->direction[0]; // already normalised
200  directionOfCluster[1] = prs2[i].cl->direction[1]; // already normalised
201  directionOfCluster[2] = prs2[i].cl->direction[2]; // already normalised
202 
203 
204  double coreEnergy = 0.0;
205  double startPosition[3];
206  double startPositionMinProjection[3];
207  double startPositionMaxProjection[3];
208  double minProjectionAlongMainPrincipalAxis = DBL_MAX;
209  double maxProjectionAlongMainPrincipalAxis = 0.0;
210 
211 
212  for(unsigned int j=0;j<prs2[i].cl->hits.size();j++) {
213 
214  // sum up core energy
215  coreEnergy += prs2[i].cl->hits[j]->chit->getEnergy();
216 
217  // debug
218  if (_drawOnCED) ced_hit( prs2[i].cl->hits[j]->chit->getPosition()[0], prs2[i].cl->hits[j]->chit->getPosition()[1],
219  prs2[i].cl->hits[j]->chit->getPosition()[2], 0 | 4 << CED_LAYER_SHIFT, 2, 0xff47e6);
220 
221 
222 
223  double vCenterToHit[3];
224  for(unsigned int k=0;k<3;k++) vCenterToHit[k] = centerPosition[k] - prs2[i].cl->hits[j]->chit->getPosition()[k];
225 
226  double projectionAlongMainPrincipalAxis = 0.0;
227  for(unsigned int k=0;k<3;k++) projectionAlongMainPrincipalAxis += directionOfCluster[k]*vCenterToHit[k];
228 
229  if ( projectionAlongMainPrincipalAxis < minProjectionAlongMainPrincipalAxis ) {
230 
231  minProjectionAlongMainPrincipalAxis = projectionAlongMainPrincipalAxis;
232  for(unsigned int k=0;k<3;k++) startPositionMinProjection[k] = centerPosition[k] - 2.0*minProjectionAlongMainPrincipalAxis*directionOfCluster[k];
233  //for(unsigned int k=0;k<3;k++) startPositionMinProjection[k] = prs2[i].cl->hits[j]->chit->getPosition()[k];
234 
235  }
236 
237  if ( projectionAlongMainPrincipalAxis > maxProjectionAlongMainPrincipalAxis ) {
238 
239  maxProjectionAlongMainPrincipalAxis = projectionAlongMainPrincipalAxis;
240  for(unsigned int k=0;k<3;k++) startPositionMaxProjection[k] = centerPosition[k] - 2.0*maxProjectionAlongMainPrincipalAxis*directionOfCluster[k];
241  //for(unsigned int k=0;k<3;k++) startPositionMaxProjection[k] = prs2[i].cl->hits[j]->chit->getPosition()[k];
242 
243  }
244 
245 
246  double rStartPositionMinProjection = sqrt( pow(startPositionMinProjection[0],2) + pow(startPositionMinProjection[1],2) + pow(startPositionMinProjection[2],2));
247  double rStartPositionMaxProjection = sqrt( pow(startPositionMaxProjection[0],2) + pow(startPositionMaxProjection[1],2) + pow(startPositionMaxProjection[2],2));
248 
249  if ( rStartPositionMinProjection < rStartPositionMaxProjection ) {
250 
251  startPosition[0] = startPositionMinProjection[0];
252  startPosition[1] = startPositionMinProjection[1];
253  startPosition[2] = startPositionMinProjection[2];
254 
255  }
256  else {
257 
258  startPosition[0] = startPositionMaxProjection[0];
259  startPosition[1] = startPositionMaxProjection[1];
260  startPosition[2] = startPositionMaxProjection[2];
261 
262  }
263 
264 
265 
266  }
267 
268 
269  double photonE = giveMeEEstimate2(prs2[i].level,prs2[i].cl->getEnergy(),coreCalibrationLDC00);
270 
271 
272 
273 
274 
275  // debug
276  if ( _debugLevel > 5 ) {
277  std::cout << "EM Shower Core found by the KIT package: " << std::endl
278  << "CoreEnergy: " << coreEnergy << " " << " estimated 'em' particle energy: " << photonE << std::endl
279  << "centerPos: " << "(" << centerPosition[0] << ", " << centerPosition[1] << ", " << centerPosition[2] << ")"
280  << " " << "startPos: " << "(" << startPosition[0] << ", " << startPosition[1] << ", " << startPosition[2] << ")"
281  << std::endl << std::endl;
282  }
283 
284  // debug
285  if (_drawOnCED) {
286  ced_hit(centerPosition[0],centerPosition[1],centerPosition[2], 2 | 4 << CED_LAYER_SHIFT, 6, 0xffffff);
287  ced_hit(startPosition[0],startPosition[1],startPosition[2], 2 | 4 << CED_LAYER_SHIFT, 6, 0xff0004);
288  }
289 
290 
291 
292 
293 
294  if ( photonE > 0.0 ) {
295 
296 
297  Photon2* photonFinder = new Photon2(photonE,centerPosition,startPosition);
298 
299  // loop over all ECAL hits
300  for(unsigned int j = 0; j < nelem; ++j) {
301 
302  double probabilityAndDistance[2];
303  CalorimeterHit* ECALHit = dynamic_cast<CalorimeterHit*>(colt->getElementAt(j));
304 
305  photonFinder->Prob(ECALHit,_probabilityDensityCut,probabilityAndDistance);
306 
307  bool isECALHitAlreadyAssigned = false;
308  int indexOfECALHitAlreadyAssigned = 0;
309  for(unsigned int k = 0; k < ECALHitsWithAttributes.size(); ++k) {
310 
311  if ( ECALHitsWithAttributes.at(k).ECALHit == ECALHit) {
312 
313  isECALHitAlreadyAssigned = true;
314  indexOfECALHitAlreadyAssigned = k;
315  break;
316 
317  }
318 
319  }
320 
321 
322 
323  if ( !isECALHitAlreadyAssigned ) {
324 
325  if ( probabilityAndDistance[0] > 0.0 ) {
326 
327  ECALHitWithAttributes hitWithAttributes;
328 
329  hitWithAttributes.ECALHit = ECALHit;
330  hitWithAttributes.relatedCores.push_back(&(prs2[i]));
331  hitWithAttributes.probabilitiesForThisECALHit.push_back(probabilityAndDistance[0]);
332  hitWithAttributes.distancesToCoresForThisECALHit.push_back(probabilityAndDistance[1]);
333  hitWithAttributes.estimatedEnergyPerCore.push_back(photonE);
334 
335  ECALHitsWithAttributes.push_back(hitWithAttributes);
336 
337  }
338 
339  }
340  else {
341 
342  if ( probabilityAndDistance[0] > 0.0 ) {
343 
344  ECALHitsWithAttributes.at(indexOfECALHitAlreadyAssigned).relatedCores.push_back(&(prs2[i]));
345  ECALHitsWithAttributes.at(indexOfECALHitAlreadyAssigned).probabilitiesForThisECALHit.push_back(probabilityAndDistance[0]);
346  ECALHitsWithAttributes.at(indexOfECALHitAlreadyAssigned).distancesToCoresForThisECALHit.push_back(probabilityAndDistance[1]);
347  ECALHitsWithAttributes.at(indexOfECALHitAlreadyAssigned).estimatedEnergyPerCore.push_back(photonE);
348 
349  }
350 
351  }
352 
353  }
354 
355  delete photonFinder;
356  photonFinder = 0;
357 
358  }
359  }
360  }
361 
362 
363 
364 
365  // debug
366 
367  if ( _debugLevel > 5 ) {
368  std::cout << "CalorimeterHits in ECAL and their probability to belong to EM showers: " << std::endl;
369  for(unsigned int i = 0; i < ECALHitsWithAttributes.size(); ++i) {
370  std::cout << "ECAL Hit: " << ECALHitsWithAttributes.at(i).ECALHit << " (" << ECALHitsWithAttributes.at(i).ECALHit->getPosition()[0] << ", "
371  << ECALHitsWithAttributes.at(i).ECALHit->getPosition()[1] << ", " << ECALHitsWithAttributes.at(i).ECALHit->getPosition()[2] << ")" << " "
372  <<"E = " << ECALHitsWithAttributes.at(i).ECALHit->getEnergy() << " " << std::endl;
373  std::cout << "Cores: " << "(";
374  for(unsigned int j = 0; j < ECALHitsWithAttributes.at(i).relatedCores.size(); ++j) {
375  std::cout << ECALHitsWithAttributes.at(i).relatedCores.at(j);
376  if (j < ECALHitsWithAttributes.at(i).relatedCores.size() - 1) std::cout << ",";
377  }
378  std::cout << ")" << " " << std::endl;
379  std::cout << "Probs: " << "(";
380  for(unsigned int j = 0; j < ECALHitsWithAttributes.at(i).probabilitiesForThisECALHit.size(); ++j) {
381  std::cout << ECALHitsWithAttributes.at(i).probabilitiesForThisECALHit.at(j);
382  if (j < ECALHitsWithAttributes.at(i).probabilitiesForThisECALHit.size() - 1) std::cout << ",";
383  }
384  std::cout << ")" << " " << std::endl;
385  std::cout << "Dists: " << "(";
386  for(unsigned int j = 0; j < ECALHitsWithAttributes.at(i).distancesToCoresForThisECALHit.size(); ++j) {
387  std::cout << ECALHitsWithAttributes.at(i).distancesToCoresForThisECALHit.at(j);
388  if (j < ECALHitsWithAttributes.at(i).distancesToCoresForThisECALHit.size() - 1) std::cout << ",";
389  }
390  std::cout << ")" << " " << std::endl;
391  std::cout << "Estimated EM Energy: " << "(";
392  for(unsigned int j = 0; j < ECALHitsWithAttributes.at(i).estimatedEnergyPerCore.size(); ++j) {
393  std::cout << ECALHitsWithAttributes.at(i).estimatedEnergyPerCore.at(j);
394  if (j < ECALHitsWithAttributes.at(i).estimatedEnergyPerCore.size() - 1) std::cout << ",";
395  }
396  std::cout << ")" << std::endl << std::endl;
397  }
398  }
399 
400 
401 
402 
403 
404  std::map<PROTSEED2*,ClusterImpl*> mapCoreToCluster;
405  std::map<PROTSEED2*,double> mapCoreToEstimatedEnergy;
406 
407  for(unsigned int i = 0; i < ECALHitsWithAttributes.size(); ++i) {
408 
409  PROTSEED2* coreToAddHitTo;
410  double estimatedEnergyOfCoreToAddHitTo = 0.0;
411 
412  if (ECALHitsWithAttributes.at(i).relatedCores.size() == 1) {
413 
414  coreToAddHitTo = ECALHitsWithAttributes.at(i).relatedCores.at(0);
415  estimatedEnergyOfCoreToAddHitTo = ECALHitsWithAttributes.at(i).estimatedEnergyPerCore.at(0);
416 
417  }
418  else {
419 
420  unsigned int indexOfMaxProbability = 0;
421  double maxProbability = 0.0;
422 
423  for(unsigned int j = 0; j < ECALHitsWithAttributes.at(i).relatedCores.size(); ++j){
424 
425  if ( ECALHitsWithAttributes.at(i).probabilitiesForThisECALHit.at(j) > maxProbability ) {
426 
427  indexOfMaxProbability = j;
428  maxProbability = ECALHitsWithAttributes.at(i).probabilitiesForThisECALHit.at(j);
429 
430  }
431 
432  }
433 
434  coreToAddHitTo = ECALHitsWithAttributes.at(i).relatedCores.at(indexOfMaxProbability);
435  estimatedEnergyOfCoreToAddHitTo = ECALHitsWithAttributes.at(i).estimatedEnergyPerCore.at(indexOfMaxProbability);
436 
437  }
438 
439 
440  std::map<PROTSEED2*,ClusterImpl*>::iterator position = mapCoreToCluster.find(coreToAddHitTo);
441 
442  if ( position == mapCoreToCluster.end() ) {
443 
444  mapCoreToCluster[coreToAddHitTo] = new ClusterImpl();
445  mapCoreToEstimatedEnergy[coreToAddHitTo] = estimatedEnergyOfCoreToAddHitTo;
446  position = mapCoreToCluster.find(coreToAddHitTo);
447 
448  }
449 
450  (*position).second->addHit(ECALHitsWithAttributes.at(i).ECALHit,(float)1.0);
451 
452  }
453 
454 
455  int iCluster = 0;
456 
457  for(std::map<PROTSEED2*,ClusterImpl*>::iterator i = mapCoreToCluster.begin(); i != mapCoreToCluster.end(); ++i) {
458 
459  ClusterImpl* cluster = (*i).second;
460 
461  double estimatedEnergyOfCluster = mapCoreToEstimatedEnergy[(*i).first];
462 
463  // set energy and position of cluster
464  unsigned int n = cluster->getCalorimeterHits().size();
465 
466  float* x = new float[n];
467  float* y = new float[n];
468  float* z = new float[n];
469  float* a = new float[n];
470  float energy = 0.0;
471 
472  for (unsigned int j = 0; j < n; ++j) {
473 
474  x[j] = cluster->getCalorimeterHits().at(j)->getPosition()[0];
475  y[j] = cluster->getCalorimeterHits().at(j)->getPosition()[1];
476  z[j] = cluster->getCalorimeterHits().at(j)->getPosition()[2];
477  a[j] = cluster->getCalorimeterHits().at(j)->getEnergy();
478 
479  energy += cluster->getCalorimeterHits().at(j)->getEnergy();
480 
481  }
482 
483 
484  // cut of EM candidtate clusters which differ more than _energyDeviationCut from enery estimate
485  double energyDeviation = fabs( (energy-estimatedEnergyOfCluster)/estimatedEnergyOfCluster );
486 
487 
488  // simple cut on the starting layer of the shower
489  // FIXME: this is not working properly for the edges of the ECAL and HCAL => need pseudo layer !!!
490  const unsigned int cutOnLayer = 10; // FIXME: remove hard-coded cut
491  const unsigned int contOnNumberOfHitsInLayerCut = 4; // FIXME: remove hard-coded cut
492 
493  unsigned int nHitsInLayerCut = 0;
494  bool hitWithinCutRangeFound = false;
495 
496  for (unsigned int j = 0; j < cluster->getCalorimeterHits().size(); ++j) {
497 
498  unsigned int layer = CDECAL(cluster->getCalorimeterHits().at(j))["K-1"];
499 
500  if ( layer <= cutOnLayer ) ++nHitsInLayerCut;
501 
502  if ( nHitsInLayerCut >= contOnNumberOfHitsInLayerCut ) {
503  hitWithinCutRangeFound = true;
504  break;
505 
506  }
507 
508  }
509 
510 
511 
512 
513  // debug
514  if ( _debugLevel > 5 ) {
515  std::cout << std::endl << "Cluster Size:" << cluster->getCalorimeterHits().size() << std::endl
516  << "energy of cluster: " << energy << " " << "estimated energy: " << estimatedEnergyOfCluster << " " << "deviation: " << energyDeviation << std::endl
517  << "Calorimeter Hist in this Cluster:" << std::endl;
518  for (unsigned int j = 0; j < cluster->getCalorimeterHits().size(); ++j) {
519 
520  std::cout << "Position: " << "(" << cluster->getCalorimeterHits().at(j)->getPosition()[0] << "," << cluster->getCalorimeterHits().at(j)->getPosition()[1] << ","
521  << cluster->getCalorimeterHits().at(j)->getPosition()[2] << ")" << " " << "E = " << cluster->getCalorimeterHits().at(j)->getEnergy() << " "
522  << "layer: " << CDECAL(cluster->getCalorimeterHits().at(j))["K-1"] << std::endl;
523 
524 
525  }
526  }
527 
528 
529  if ( (energyDeviation < _energyDeviationCut) && hitWithinCutRangeFound ) {
530 
531 
532 
533  // flag hits
534  for (unsigned int j = 0; j < cluster->getCalorimeterHits().size(); ++j) cluster->getCalorimeterHits().at(j)->ext<isPartOfEMShowerCandidate>() = 1;
535 
536  ClusterShapes* clusterShape = new ClusterShapes(n,a,x,y,z);
537 
538  float position[3];
539  position[0]=clusterShape->getCentreOfGravity()[0];
540  position[1]=clusterShape->getCentreOfGravity()[1];
541  position[2]=clusterShape->getCentreOfGravity()[2];
542 
543 
544  cluster->setPosition(position);
545  cluster->setEnergy(energy);
546  clscol->addElement(cluster);
547 
548 
549  // debug
550  if ( _debugLevel > 5 ) {
551  std::cout << "EM Shower Candidate FOUND : " << "cluster " << cluster << " with energy E = " << cluster->getEnergy() << " " << "n of hits = "
552  << cluster->getCalorimeterHits().size() << std::endl;
553  }
554 
555 
556  if (_drawOnCED) {
557 
558  int color = 0x8f57ff * ( iCluster + 32);
559  MarlinCED::drawClusterImpl(cluster,2,2,color,5);
560  ced_send_event();
561  ++iCluster;
562  getchar();
563 
564  }
565 
566 
567 
568  delete clusterShape;
569  clusterShape = 0;
570 
571  }
572 
573 
574 
575  delete[] x;
576  x = 0;
577  delete[] y;
578  y = 0;
579  delete[] z;
580  z = 0;
581  delete[] a;
582  a = 0;
583 
584  }
585 
586  mapCoreToCluster.clear();
587  ECALHitsWithAttributes.clear();
588 
589 
590 
591 
592  // for strong memory and nice dreams ..
593  for(unsigned int i=0;i<N;i++)
594  {
595  if( bbb[i].size()!=0)
596  for(unsigned int im=0;im<bbb[i].size();im++)
597  {
598  delete bbb[i][im];
599  }
600  }
601 
602  for(unsigned int im=0;im<2;im++)
603  {
604  if(calo[im].size()!=0)
605  for( unsigned int iij=0;iij<calo[im].size();iij++)
606  delete (calo[im])[iij];
607  }
608 
609  }
610 
611  evt->addCollection(clscol,_collectionNameOfEMShowerCandidates.c_str());
612 
613  }catch(DataNotAvailableException &e) {std::cout << "no valid ECAL collection in event " << _nEvt << std::endl; }
614 
615 
616  // debug
617  if (_drawOnCED) MarlinCED::draw(this,true);
618 
619 
620  ++_nEvt;
621 
622 }
623 
624 
625 void EMShowerFinder::check( LCEvent * /*evt*/ ) {
626 
627 }
628 
629 
630 void EMShowerFinder::end() {
631 
632 }
ECALHitWithAttributes::distancesToCoresForThisECALHit
std::vector< double > distancesToCoresForThisECALHit
Definition: EMShowerFinder.h:42
Photon2::Prob
void Prob(CalorimeterHit *ch, double cut, double *out)
Definition: KITutil.cc:1361
EMShowerFinder::processEvent
virtual void processEvent(LCEvent *evt)
Definition: EMShowerFinder.cc:128
ECALHitWithAttributes::relatedCores
std::vector< PROTSEED2 * > relatedCores
Definition: EMShowerFinder.h:40
EMShowerFinder::_Distcut
double _Distcut
Definition: EMShowerFinder.h:95
CoreCut2
container for keeping the parameters of the core fineder together
Definition: KITutil.h:289
CoreCut2::MinHitSplit
unsigned int MinHitSplit
minimal number of hits in i-th level cluster to be accepted for splitting of the core ...
Definition: KITutil.h:309
EMShowerFinder::_debugLevel
int _debugLevel
Definition: EMShowerFinder.h:100
EMShowerFinder::_CleanCut
int _CleanCut
Definition: EMShowerFinder.h:89
EMShowerFinder::_collectionNameOfEMShowerCandidates
std::string _collectionNameOfEMShowerCandidates
Definition: EMShowerFinder.h:87
sistrip::k
static const float k
Definition: PhysicalConstants.h:73
CoreCut2::MinHit0
unsigned int MinHit0
minimal number of hits needed for 0-th level cluster to be accepted as a core candidate ...
Definition: KITutil.h:305
EMShowerFinder::_ToClean
std::string _ToClean
Definition: EMShowerFinder.h:88
EMShowerFinder::_Rcut
double _Rcut
Definition: EMShowerFinder.h:94
EMShowerFinder::processRunHeader
virtual void processRunHeader(LCRunHeader *run)
Definition: EMShowerFinder.cc:121
EMShowerFinder::_drawOnCED
int _drawOnCED
Definition: EMShowerFinder.h:101
EMShowerFinder::_miipstep
vector< float > _miipstep
Definition: EMShowerFinder.h:91
FindCores2
void FindCores2(Shitvec2 *secal1, vector< Tmpclvec2 > &bbb, vector< PROTSEED2 > *prs2, unsigned int N, vector< float > miipstep, CoreCut2 Ccut)
Global EM core finding function , iput is vector of superhits, array of Tmpcl vectors bbb for interna...
Definition: KITutil.cc:286
EMShowerFinder::_Coscut
double _Coscut
Definition: EMShowerFinder.h:96
CoreCut2::Distcut
double Distcut
distance cut for core merging
Definition: KITutil.h:297
EMShowerFinder::_colNameECAL
std::string _colNameECAL
Definition: EMShowerFinder.h:86
EMShowerFinder::_MinHit0
int _MinHit0
Definition: EMShowerFinder.h:92
ECALHitWithAttributes::probabilitiesForThisECALHit
std::vector< double > probabilitiesForThisECALHit
Definition: EMShowerFinder.h:41
EMShowerFinder::EMShowerFinder
EMShowerFinder()
Definition: EMShowerFinder.cc:10
EMShowerFinder::check
virtual void check(LCEvent *evt)
Definition: EMShowerFinder.cc:625
TotalPrecalc2
void TotalPrecalc2(vector< Superhit2 * > *calo, CellIDDecoder< CalorimeterHit > &id, unsigned int nelem, int Ncut)
Global precalculation function , iput is vector of superhits, ECAL decoder, number of hits...
Definition: KITutil.cc:29
aEMShowerFinder
EMShowerFinder aEMShowerFinder
Definition: EMShowerFinder.cc:6
CoreCut2::Coscut
double Coscut
angular cut for core merging (value of the cosine is stored not the angle!)
Definition: KITutil.h:301
PROTSEED2
container for storing the EM shower core candidates
Definition: KITutil.h:267
sistrip::e
static const float e
Definition: PhysicalConstants.h:26
ECALHitWithAttributes::estimatedEnergyPerCore
std::vector< double > estimatedEnergyPerCore
Definition: EMShowerFinder.h:43
EMShowerFinder::end
virtual void end()
Definition: EMShowerFinder.cc:630
EMShowerFinder::_energyDeviationCut
double _energyDeviationCut
Definition: EMShowerFinder.h:97
sistrip::LCCollectionVec
std::vector< LCCollection * > LCCollectionVec
Definition: SiStripClus.h:55
EMShowerFinder::_nEvt
int _nEvt
Definition: EMShowerFinder.h:104
EMShowerFinder::_probabilityDensityCut
double _probabilityDensityCut
Definition: EMShowerFinder.h:98
EMShowerFinder::_nRun
int _nRun
Definition: EMShowerFinder.h:103
CoreCut2::Rcut
double Rcut
fluctuation suppresion cut
Definition: KITutil.h:293
giveMeEEstimate2
double giveMeEEstimate2(int nivo, double Ecore, vector< CoreCalib2 > cc)
returns energy estimate for a given core , input int level, double core energy in GeV...
Definition: KITutil.cc:1450
EMShowerFinder::init
virtual void init()
Definition: EMShowerFinder.cc:102
EMShowerFinder::_MinHitSplit
int _MinHitSplit
Definition: EMShowerFinder.h:93
EMShowerFinder::_N
int _N
Definition: EMShowerFinder.h:90
CreateCalibrationLDC00
void CreateCalibrationLDC00(vector< CoreCalib2 > *cc)
Example function for creation of the energy estimaiton table.
Definition: KITutil.cc:1463
ECALHitWithAttributes::ECALHit
CalorimeterHit * ECALHit
Definition: EMShowerFinder.h:39
CreateAllShits2
void CreateAllShits2(LCCollection *colt, CellIDDecoder< CalorimeterHit > &id, vector< Superhit2 * > *calo)
Creation of superhits, input is ECAL collection ,it&#39;s decoded and pointer to resulting container for ...
Definition: KITutil.cc:7