1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 1 October 2000. Yuri Kharlov:
20 PPSD upper layer is considered if number of layers>1
21 18 October 2000. Yuri Kharlov:
22 AliEMCALClusterizerv1()
23 CPV clusterizing parameters added
25 After first PPSD digit remove EMC digits only once
27 //*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (SUBATECH & Kurchatov Institute)
28 // August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
29 // of new IO (à la PHOS)
30 //////////////////////////////////////////////////////////////////////////////
31 // Clusterization class. Performs clusterization (collects neighbouring active cells) and
32 // unfolds the clusters having several local maxima.
33 // Results are stored in TreeR#, branches EMCALTowerRP (EMC recPoints),
34 // EMCALPreShoRP (CPV RecPoints) and AliEMCALClusterizer (Clusterizer with all
35 // parameters including input digits branch title, thresholds etc.)
36 // This TTask is normally called from Reconstructioner, but can as well be used in
39 // root [0] AliEMCALClusterizerv1 * cl = new AliEMCALClusterizerv1("galice.root")
40 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
41 // //reads gAlice from header file "..."
42 // root [1] cl->ExecuteTask()
43 // //finds RecPoints in all events stored in galice.root
44 // root [2] cl->SetDigitsBranch("digits2")
45 // //sets another title for Digitis (input) branch
46 // root [3] cl->SetRecPointsBranch("recp2")
47 // //sets another title four output branches
48 // root [4] cl->SetTowerLocalMaxCut(0.03)
49 // //set clusterization parameters
50 // root [5] cl->ExecuteTask("deb all time")
51 // //once more finds RecPoints options are
52 // // deb - print number of found rec points
53 // // deb all - print number of found RecPoints and some their characteristics
54 // // time - print benchmarking results
56 // --- ROOT system ---
65 #include "TBenchmark.h"
67 // --- Standard library ---
70 // --- AliRoot header files ---
72 #include "AliEMCALClusterizerv1.h"
73 #include "AliEMCALDigit.h"
74 #include "AliEMCALDigitizer.h"
75 #include "AliEMCALTowerRecPoint.h"
77 #include "AliEMCALGetter.h"
78 #include "AliEMCALGeometry.h"
81 ClassImp(AliEMCALClusterizerv1)
83 //____________________________________________________________________________
84 AliEMCALClusterizerv1::AliEMCALClusterizerv1() : AliEMCALClusterizer()
86 // default ctor (to be used mainly by Streamer)
89 fDefaultInit = kTRUE ;
92 //____________________________________________________________________________
93 AliEMCALClusterizerv1::AliEMCALClusterizerv1(const char* headerFile, const char* name, const Bool_t toSplit)
94 :AliEMCALClusterizer(headerFile, name, toSplit)
96 // ctor with the indication of the file where header Tree and digits Tree are stored
100 fDefaultInit = kFALSE ;
104 //____________________________________________________________________________
105 AliEMCALClusterizerv1::~AliEMCALClusterizerv1()
112 //____________________________________________________________________________
113 const TString AliEMCALClusterizerv1::BranchName() const
115 TString branchName(GetName() ) ;
116 branchName.Remove(branchName.Index(Version())-1) ;
120 //____________________________________________________________________________
121 Float_t AliEMCALClusterizerv1::Calibrate(Int_t amp, Bool_t inpresho) const
123 //To be replased later by the method, reading individual parameters from the database
124 if ( inpresho ) // calibrate as pre shower
125 return -fADCpedestalPreSho + amp * fADCchannelPreSho ;
126 else //calibrate as tower
127 return -fADCpedestalTower + amp * fADCchannelTower ;
130 //____________________________________________________________________________
131 void AliEMCALClusterizerv1::Exec(Option_t * option)
135 if( strcmp(GetName(), "")== 0 )
138 if(strstr(option,"tim"))
139 gBenchmark->Start("EMCALClusterizer");
141 if(strstr(option,"print"))
144 AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ;
145 if(gime->BranchExists("RecPoints"))
147 Int_t nevents = gime->MaxEvent() ;
150 for(ievent = 0; ievent < nevents; ievent++){
152 gime->Event(ievent,"D") ;
155 GetCalibrationParameters() ;
157 fNumberOfTowerClusters = fNumberOfPreShoClusters = 0 ;
164 WriteRecPoints(ievent) ;
166 if(strstr(option,"deb"))
167 PrintRecPoints(option) ;
169 //increment the total number of digits per run
170 fRecPointsInRun += gime->TowerRecPoints()->GetEntriesFast() ;
171 fRecPointsInRun += gime->PreShowerRecPoints()->GetEntriesFast() ;
174 if(strstr(option,"tim")){
175 gBenchmark->Stop("EMCALClusterizer");
176 Info("Exec", "took %f seconds for Clusterizing %f seconds per event",
177 gBenchmark->GetCpuTime("EMCALClusterizer"), gBenchmark->GetCpuTime("EMCALClusterizer")/nevents ) ;
182 //____________________________________________________________________________
183 Bool_t AliEMCALClusterizerv1::FindFit(AliEMCALTowerRecPoint * emcRP, AliEMCALDigit ** maxAt, Float_t * maxAtEnergy,
184 Int_t nPar, Float_t * fitparameters) const
186 // Calls TMinuit to fit the energy distribution of a cluster with several maxima
187 // The initial values for fitting procedure are set equal to the positions of local maxima.
188 // Cluster will be fitted as a superposition of nPar/3 electromagnetic showers
190 AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ;
191 TClonesArray * digits = gime->Digits() ;
194 gMinuit->mncler(); // Reset Minuit's list of paramters
195 gMinuit->SetPrintLevel(-1) ; // No Printout
196 gMinuit->SetFCN(AliEMCALClusterizerv1::UnfoldingChiSquare) ;
197 // To set the address of the minimization function
198 TList * toMinuit = new TList();
199 toMinuit->AddAt(emcRP,0) ;
200 toMinuit->AddAt(digits,1) ;
202 gMinuit->SetObjectFit(toMinuit) ; // To tranfer pointer to UnfoldingChiSquare
204 // filling initial values for fit parameters
205 AliEMCALDigit * digit ;
209 Int_t nDigits = (Int_t) nPar / 3 ;
213 AliEMCALGeometry * geom = gime->EMCALGeometry() ;
215 for(iDigit = 0; iDigit < nDigits; iDigit++){
216 digit = maxAt[iDigit];
221 geom->AbsToRelNumbering(digit->GetId(), relid) ;
222 geom->PosInAlice(relid, x, z) ;
224 Float_t energy = maxAtEnergy[iDigit] ;
226 gMinuit->mnparm(index, "x", x, 0.1, 0, 0, ierflg) ;
229 Error("FindFit", "EMCAL Unfolding Unable to set initial value for fit procedure : x = %f", x ) ;
232 gMinuit->mnparm(index, "z", z, 0.1, 0, 0, ierflg) ;
235 Error("FindFit", "EMCAL Unfolding Unable to set initial value for fit procedure : z = %f", z) ;
238 gMinuit->mnparm(index, "Energy", energy , 0.05*energy, 0., 4.*energy, ierflg) ;
241 Error("FindFit", "EMCAL Unfolding Unable to set initial value for fit procedure : energy = %f", energy) ;
246 Double_t p0 = 0.1 ; // "Tolerance" Evaluation stops when EDM = 0.0001*p0 ; The number of function call slightly
251 gMinuit->mnexcm("SET STR", &p2, 0, ierflg) ; // force TMinuit to reduce function calls
252 gMinuit->mnexcm("SET GRA", &p1, 1, ierflg) ; // force TMinuit to use my gradient
253 gMinuit->SetMaxIterations(5);
254 gMinuit->mnexcm("SET NOW", &p2 , 0, ierflg) ; // No Warnings
256 gMinuit->mnexcm("MIGRAD", &p0, 0, ierflg) ; // minimize
258 if(ierflg == 4){ // Minimum not found
259 Error("FindFit", "EMCAL Unfolding Fit not converged, cluster abandoned " ) ;
262 for(index = 0; index < nPar; index++){
265 gMinuit->GetParameter(index, val, err) ; // Returns value and error of parameter index
266 fitparameters[index] = val ;
274 //____________________________________________________________________________
275 void AliEMCALClusterizerv1::GetCalibrationParameters()
277 AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ;
278 const AliEMCALDigitizer * dig = gime->Digitizer(BranchName()) ;
280 fADCchannelTower = dig->GetTowerchannel() ;
281 fADCpedestalTower = dig->GetTowerpedestal();
283 fADCchannelPreSho = dig->GetPreShochannel() ;
284 fADCpedestalPreSho = dig->GetPreShopedestal() ;
287 //____________________________________________________________________________
288 void AliEMCALClusterizerv1::Init()
290 // Make all memory allocations which can not be done in default constructor.
291 // Attach the Clusterizer task to the list of EMCAL tasks
293 if ( strcmp(GetTitle(), "") == 0 )
294 SetTitle("galice.root") ;
296 TString branchname = GetName() ;
297 branchname.Remove(branchname.Index(Version())-1) ;
299 AliEMCALGetter * gime = AliEMCALGetter::GetInstance(GetTitle(), branchname.Data(), fToSplit ) ;
301 Error("Init", "Could not obtain the Getter object !" ) ;
307 // construct the name of the file as /path/EMCAL.SDigits.root
308 //First - extract full path if necessary
309 TString fileName(GetTitle()) ;
310 Ssiz_t islash = fileName.Last('/') ;
311 if(islash<fileName.Length())
312 fileName.Remove(islash+1,fileName.Length()) ;
315 // Next - append the file name
316 fileName+="EMCAL.RecData." ;
317 if((strcmp(branchname.Data(),"Default")!=0)&&(strcmp(branchname.Data(),"")!=0)){
318 fileName+=branchname ;
322 // Finally - check if the file already opened or open the file
323 fSplitFile = static_cast<TFile*>(gROOT->GetFile(fileName.Data()));
325 fSplitFile = TFile::Open(fileName.Data(),"update") ;
328 const AliEMCALGeometry * geom = gime->EMCALGeometry() ;
329 fNTowers = geom->GetNZ() * geom->GetNPhi() ;
331 gMinuit = new TMinuit(100) ;
333 gime->PostClusterizer(this) ;
334 gime->PostRecPoints(branchname ) ;
338 //____________________________________________________________________________
339 void AliEMCALClusterizerv1::InitParameters()
341 fNumberOfPreShoClusters = fNumberOfTowerClusters = 0 ;
342 fPreShoClusteringThreshold = 0.0001;
343 fTowerClusteringThreshold = 0.2;
344 fTowerLocMaxCut = 0.03 ;
345 fPreShoLocMaxCut = 0.03 ;
354 TString clusterizerName( GetName()) ;
355 if (clusterizerName.IsNull() )
356 clusterizerName = "Default" ;
357 clusterizerName.Append(":") ;
358 clusterizerName.Append(Version()) ;
359 SetName(clusterizerName) ;
360 fRecPointsInRun = 0 ;
364 //____________________________________________________________________________
365 Int_t AliEMCALClusterizerv1::AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2)const
367 // Gives the neighbourness of two digits = 0 are not neighbour but continue searching
369 // = 2 are not neighbour but do not continue searching
370 // neighbours are defined as digits having at least a common vertex
371 // The order of d1 and d2 is important: first (d1) should be a digit already in a cluster
372 // which is compared to a digit (d2) not yet in a cluster
374 AliEMCALGeometry * geom = AliEMCALGetter::GetInstance()->EMCALGeometry() ;
379 geom->AbsToRelNumbering(d1->GetId(), relid1) ;
382 geom->AbsToRelNumbering(d2->GetId(), relid2) ;
384 if ( (relid1[0] == relid2[0]) && (relid1[1]==relid2[1]) ) { // inside the same EMCAL Arm
385 Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
386 Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
388 if (( coldiff <= 1 ) && ( rowdiff <= 1 )){
389 if((relid1[1] != 0) || (TMath::Abs(d1->GetTime() - d2->GetTime() ) < fTimeGate))
393 if((relid2[2] > relid1[2]) && (relid2[3] > relid1[3]+1))
394 rv = 2; // Difference in row numbers is too large to look further
400 if( (relid1[0] < relid2[0]) || (relid1[1] != relid2[1]) )
406 //____________________________________________________________________________
407 Bool_t AliEMCALClusterizerv1::IsInTower(AliEMCALDigit * digit) const
409 // Tells if (true) or not (false) the digit is in a EMCAL-Tower
412 if (!digit->IsInPreShower())
417 //____________________________________________________________________________
418 Bool_t AliEMCALClusterizerv1::IsInPreShower(AliEMCALDigit * digit) const
420 // Tells if (true) or not (false) the digit is in a EMCAL-PreShower
423 if (digit->IsInPreShower())
428 //____________________________________________________________________________
429 void AliEMCALClusterizerv1::WriteRecPoints(Int_t event)
432 // Creates new branches with given title
433 // fills and writes into TreeR.
435 AliEMCALGetter *gime = AliEMCALGetter::GetInstance() ;
436 TObjArray * towerRecPoints = gime->TowerRecPoints() ;
437 TObjArray * preshoRecPoints = gime->PreShowerRecPoints() ;
438 TClonesArray * digits = gime->Digits() ;
445 TString name("TreeR") ;
447 treeR = dynamic_cast<TTree*>(fSplitFile->Get(name));
450 treeR = gAlice->TreeR();
454 gAlice->MakeTree("R", fSplitFile);
455 treeR = gAlice->TreeR() ;
459 //Evaluate position, dispersion and other RecPoint properties...
460 for(index = 0; index < towerRecPoints->GetEntries(); index++)
461 (dynamic_cast<AliEMCALTowerRecPoint *>(towerRecPoints->At(index)))->EvalAll(fW0,digits) ;
463 towerRecPoints->Sort() ;
465 for(index = 0; index < towerRecPoints->GetEntries(); index++)
466 (dynamic_cast<AliEMCALTowerRecPoint *>(towerRecPoints->At(index)))->SetIndexInList(index) ;
468 towerRecPoints->Expand(towerRecPoints->GetEntriesFast()) ;
469 //Now the same for pre shower
470 for(index = 0; index < preshoRecPoints->GetEntries(); index++)
471 (dynamic_cast<AliEMCALRecPoint *>(preshoRecPoints->At(index)))->EvalAll(fW0CPV,digits) ;
472 preshoRecPoints->Sort() ;
474 for(index = 0; index < preshoRecPoints->GetEntries(); index++)
475 (dynamic_cast<AliEMCALRecPoint *>(preshoRecPoints->At(index)))->SetIndexInList(index) ;
477 preshoRecPoints->Expand(preshoRecPoints->GetEntriesFast()) ;
479 Int_t bufferSize = 32000 ;
480 Int_t splitlevel = 0 ;
483 TBranch * towerBranch = treeR->Branch("EMCALTowerRP","TObjArray",&towerRecPoints,bufferSize,splitlevel);
484 towerBranch->SetTitle(BranchName());
486 //Now Pre Shower branch
487 TBranch * preshoBranch = treeR->Branch("EMCALPreShoRP","TObjArray",&preshoRecPoints,bufferSize,splitlevel);
488 preshoBranch->SetTitle(BranchName());
490 //And Finally clusterizer branch
491 AliEMCALClusterizerv1 * cl = (AliEMCALClusterizerv1*)gime->Clusterizer(BranchName()) ;
492 TBranch * clusterizerBranch = treeR->Branch("AliEMCALClusterizer","AliEMCALClusterizerv1",
493 &cl,bufferSize,splitlevel);
494 clusterizerBranch->SetTitle(BranchName());
496 towerBranch ->Fill() ;
497 preshoBranch ->Fill() ;
498 clusterizerBranch->Fill() ;
500 treeR->AutoSave() ; //Write(0,kOverwrite) ;
501 if(gAlice->TreeR()!=treeR)
505 //____________________________________________________________________________
506 void AliEMCALClusterizerv1::MakeClusters()
508 // Steering method to construct the clusters stored in a list of Reconstructed Points
509 // A cluster is defined as a list of neighbour digits
511 AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ;
513 TObjArray * towerRecPoints = gime->TowerRecPoints(BranchName()) ;
514 TObjArray * preshoRecPoints = gime->PreShowerRecPoints(BranchName()) ;
515 towerRecPoints->Delete() ;
516 preshoRecPoints->Delete() ;
518 TClonesArray * digits = gime->Digits() ;
520 Fatal("MakeClusters -> Digits with name %s not found", GetName() ) ;
522 TClonesArray * digitsC = dynamic_cast<TClonesArray*>(digits->Clone()) ;
525 // Clusterization starts
527 TIter nextdigit(digitsC) ;
528 AliEMCALDigit * digit ;
529 Bool_t notremoved = kTRUE ;
531 while ( (digit = dynamic_cast<AliEMCALDigit *>(nextdigit())) ) { // scan over the list of digitsC
532 AliEMCALRecPoint * clu = 0 ;
534 TArrayI clusterdigitslist(1500) ;
536 if (( IsInTower (digit) && Calibrate(digit->GetAmp(),digit->IsInPreShower()) > fTowerClusteringThreshold ) ||
537 ( IsInPreShower (digit) && Calibrate(digit->GetAmp(),digit->IsInPreShower()) > fPreShoClusteringThreshold ) ) {
539 Int_t iDigitInCluster = 0 ;
541 if ( IsInTower(digit) ) {
542 // start a new Tower RecPoint
543 if(fNumberOfTowerClusters >= towerRecPoints->GetSize())
544 towerRecPoints->Expand(2*fNumberOfTowerClusters+1) ;
546 towerRecPoints->AddAt(new AliEMCALTowerRecPoint(""), fNumberOfTowerClusters) ;
547 clu = dynamic_cast<AliEMCALTowerRecPoint *>(towerRecPoints->At(fNumberOfTowerClusters)) ;
548 fNumberOfTowerClusters++ ;
549 clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->IsInPreShower())) ;
550 clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;
552 digitsC->Remove(digit) ;
556 // start a new Pre Shower cluster
557 if(fNumberOfPreShoClusters >= preshoRecPoints->GetSize())
558 preshoRecPoints->Expand(2*fNumberOfPreShoClusters+1);
560 preshoRecPoints->AddAt(new AliEMCALTowerRecPoint(""), fNumberOfPreShoClusters) ;
562 clu = dynamic_cast<AliEMCALTowerRecPoint *>(preshoRecPoints->At(fNumberOfPreShoClusters)) ;
563 fNumberOfPreShoClusters++ ;
564 clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->IsInPreShower() ) );
565 clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;
567 digitsC->Remove(digit) ;
570 // Here we remove remaining Tower digits, which cannot make a cluster
573 while( ( digit = dynamic_cast<AliEMCALDigit *>(nextdigit()) ) ) {
574 if( IsInTower(digit) )
575 digitsC->Remove(digit) ;
579 notremoved = kFALSE ;
586 AliEMCALDigit * digitN ;
588 while (index < iDigitInCluster){ // scan over digits already in cluster
589 digit = (AliEMCALDigit*)digits->At(clusterdigitslist[index]) ;
591 while ( (digitN = (AliEMCALDigit *)nextdigit()) ) { // scan over the reduced list of digits
592 Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
594 case 0 : // not a neighbour
596 case 1 : // are neighbours
597 clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), digitN->IsInPreShower() ) ) ;
598 clusterdigitslist[iDigitInCluster] = digitN->GetIndexInList() ;
600 digitsC->Remove(digitN) ;
602 case 2 : // too far from each other
610 } // loop over cluster
616 //____________________________________________________________________________
617 void AliEMCALClusterizerv1::MakeUnfolding()
619 Fatal("AliEMCALClusterizerv1::MakeUnfolding", "--> Unfolding not implemented") ;
623 //____________________________________________________________________________
624 Double_t AliEMCALClusterizerv1::ShowerShape(Double_t r)
626 // Shape of the shower (see EMCAL TDR)
627 // If you change this function, change also the gradient evaluation in ChiSquare()
629 Double_t r4 = r*r*r*r ;
630 Double_t r295 = TMath::Power(r, 2.95) ;
631 Double_t shape = TMath::Exp( -r4 * (1. / (2.32 + 0.26 * r4) + 0.0316 / (1 + 0.0652 * r295) ) ) ;
635 //____________________________________________________________________________
636 void AliEMCALClusterizerv1::UnfoldCluster(AliEMCALTowerRecPoint * iniTower,
638 AliEMCALDigit ** maxAt,
639 Float_t * maxAtEnergy)
641 // Performs the unfolding of a cluster with nMax overlapping showers
643 Fatal("UnfoldCluster", "--> Unfolding not implemented") ;
647 //_____________________________________________________________________________
648 void AliEMCALClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Double_t & fret, Double_t * x, Int_t iflag)
650 // Calculates the Chi square for the cluster unfolding minimization
651 // Number of parameters, Gradient, Chi squared, parameters, what to do
653 ::Fatal("UnfoldingChiSquare","Unfolding not implemented") ;
655 //____________________________________________________________________________
656 void AliEMCALClusterizerv1::Print(Option_t * option)const
658 // Print clusterizer parameters
660 TString message("\n") ;
662 if( strcmp(GetName(), "") !=0 ){
666 TString taskName(GetName()) ;
667 taskName.ReplaceAll(Version(), "") ;
669 message += "--------------- " ;
670 message += taskName.Data() ;
672 message += GetTitle() ;
673 message += "-----------\n" ;
674 message += "Clusterizing digits from the file: " ;
675 message += taskName.Data() ;
676 message += "\n Branch: " ;
677 message += GetName() ;
678 message += "\n EMC Clustering threshold = " ;
679 message += fTowerClusteringThreshold ;
680 message += "\n EMC Local Maximum cut = " ;
681 message += fTowerLocMaxCut ;
682 message += "\n EMC Logarothmic weight = " ;
684 message += "\n CPV Clustering threshold = " ;
685 message += fPreShoClusteringThreshold ;
686 message += "\n CPV Local Maximum cut = " ;
687 message += fPreShoLocMaxCut ;
688 message += "\n CPV Logarothmic weight = " ;
691 message +="\nUnfolding on\n" ;
693 message += "\nUnfolding off\n";
695 message += "------------------------------------------------------------------" ;
698 message += "AliEMCALClusterizerv1 not initialized " ;
700 Info("Print", message.Data() ) ;
703 //____________________________________________________________________________
704 void AliEMCALClusterizerv1::PrintRecPoints(Option_t * option)
706 // Prints list of RecPoints produced at the current pass of AliEMCALClusterizer
708 TObjArray * towerRecPoints = AliEMCALGetter::GetInstance()->TowerRecPoints() ;
709 TObjArray * preshoRecPoints = AliEMCALGetter::GetInstance()->PreShowerRecPoints() ;
711 TString message("\n") ;
713 message += "event " ;
714 message += gAlice->GetEvNumber() ;
715 message += "\n Found " ;
716 message += towerRecPoints->GetEntriesFast() ;
717 message += " TOWER Rec Points and " ;
718 message += preshoRecPoints->GetEntriesFast() ;
719 message += " PRE SHOWER RecPoints\n" ;
721 fRecPointsInRun += towerRecPoints->GetEntriesFast() ;
722 fRecPointsInRun += preshoRecPoints->GetEntriesFast() ;
724 char * tempo = new char[8192];
726 if(strstr(option,"all")) {
728 message += "Tower clusters\n" ;
729 message += "Index Ene(MeV) Multi Module phi r theta Lambda 1 Lambda 2 # of prim Primaries list\n" ;
732 for (index = 0 ; index < towerRecPoints->GetEntries() ; index++) {
733 AliEMCALTowerRecPoint * rp = dynamic_cast<AliEMCALTowerRecPoint * >(towerRecPoints->At(index)) ;
735 rp->GetGlobalPosition(globalpos);
737 rp->GetElipsAxis(lambda);
740 primaries = rp->GetPrimaries(nprimaries);
741 sprintf(tempo, "\n%6d %8.2f %3d %2d %4.1f %4.1f %4.1f %4f %4f %2d : ",
742 rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(), rp->GetEMCALArm(),
743 globalpos.X(), globalpos.Y(), globalpos.Z(), lambda[0], lambda[1], nprimaries) ;
746 for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
747 sprintf(tempo, "%d ", primaries[iprimary] ) ;
752 //Now plot Pre shower recPoints
754 message += "\n-----------------------------------------------------------------------\n" ;
756 message += "PreShower clusters\n" ;
757 message += "Index Ene(MeV) Multi Module phi r theta Lambda 1 Lambda 2 # of prim Primaries list\n" ;
759 for (index = 0 ; index < preshoRecPoints->GetEntries() ; index++) {
760 AliEMCALTowerRecPoint * rp = dynamic_cast<AliEMCALTowerRecPoint *>(preshoRecPoints->At(index)) ;
762 rp->GetGlobalPosition(globalpos);
764 rp->GetElipsAxis(lambda);
767 primaries = rp->GetPrimaries(nprimaries);
768 sprintf(tempo, "\n%6d %8.2f %3d %2d %4.1f %4.1f %4.1f %4f %4f %2d : ",
769 rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(), rp->GetEMCALArm(),
770 globalpos.X(), globalpos.Y(), globalpos.Z(), lambda[0], lambda[1], nprimaries) ;
773 for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
774 sprintf(tempo, "%d ", primaries[iprimary] ) ;
779 message += "\n-----------------------------------------------------------------------" ;
782 Info("PrintRecPoints", message.Data() ) ;