//*-- Author: Heather Gray (LBL) merged AliEMCALRecPoint and AliEMCALTowerRecPoint 02/04
// --- ROOT system ---
-#include <Riostream.h>
+class Riostream;
#include <TPad.h>
-#include <TGraph.h>
-#include <TPaveText.h>
+class TGraph;
+class TPaveText;
#include <TClonesArray.h>
#include <TMath.h>
// --- Standard library ---
// --- AliRoot header files ---
-#include "AliGenerator.h"
+//#include "AliGenerator.h"
+class AliGenerator;
#include "AliRunLoader.h"
#include "AliRun.h"
-#include "AliEMCAL.h"
+class AliEMCAL;
#include "AliEMCALLoader.h"
#include "AliEMCALGeometry.h"
#include "AliEMCALHit.h"
//____________________________________________________________________________
AliEMCALRecPoint::AliEMCALRecPoint()
- : AliRecPoint()
+ : AliRecPoint(),
+ fGeomPtr(0),
+ fClusterType(-1),
+ fCoreEnergy(0),
+ fDispersion(0),
+ fEnergyList(0),
+ fTimeList(0),
+ fAbsIdList(0),
+ fTime(0.),
+ fCoreRadius(10), //HG check this
+ fDETracksList(0),
+ fMulParent(0),
+ fMaxParent(0),
+ fParentsList(0),
+ fDEParentsList(0),
+ fSuperModuleNumber(0)
{
// ctor
- fClusterType = -1;
fMaxTrack = 0 ;
- fMulDigit = 0 ;
- fMaxParent = 0;
- fMulParent = 0;
+ fMulDigit = 0 ;
fAmp = 0. ;
- fCoreEnergy = 0 ;
- fEnergyList = 0 ;
- fTimeList = 0 ;
- fAbsIdList = 0;
- fParentsList = 0;
- fTime = 0. ;
- // fLocPos.SetX(1.e+6) ; //Local position should be evaluated
- fCoreRadius = 10; //HG Check this
AliRunLoader *rl = AliRunLoader::GetRunLoader();
- fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
- //fGeom = AliEMCALGeometry::GetInstance();
- fGeom->GetTransformationForSM(); // Global <-> Local
+ if (rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL"))
+ fGeomPtr = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
+ else
+ fGeomPtr = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaulGeometryName());
+ fGeomPtr->GetTransformationForSM(); // Global <-> Local
}
//____________________________________________________________________________
-AliEMCALRecPoint::AliEMCALRecPoint(const char * opt) : AliRecPoint(opt)
+AliEMCALRecPoint::AliEMCALRecPoint(const char * opt)
+ : AliRecPoint(opt),
+ fGeomPtr(0),
+ fClusterType(-1),
+ fCoreEnergy(0),
+ fDispersion(0),
+ fEnergyList(0),
+ fTimeList(0),
+ fAbsIdList(0),
+ fTime(-1.),
+ fCoreRadius(10), //HG check this
+ fDETracksList(0),
+ fMulParent(0),
+ fMaxParent(1000),
+ fParentsList(0),
+ fDEParentsList(0),
+ fSuperModuleNumber(0)
{
// ctor
- fClusterType = -1;
fMaxTrack = 1000 ;
- fMaxParent = 1000;
fMulDigit = 0 ;
- fMulParent = 0;
fAmp = 0. ;
- fCoreEnergy = 0 ;
- fEnergyList = 0 ;
- fTimeList = 0 ;
- fAbsIdList = 0;
+ fDETracksList = new Float_t[fMaxTrack];
fParentsList = new Int_t[fMaxParent];
- fTime = -1. ;
- //fLocPos.SetX(1.e+6) ; //Local position should be evaluated
- fCoreRadius = 10; //HG Check this
- //fGeom = AliEMCALGeometry::GetInstance();
+ fDEParentsList = new Float_t[fMaxParent];
+ for (Int_t i = 0; i < fMaxTrack; i++)
+ fDETracksList[i] = 0;
+ for (Int_t i = 0; i < fMaxParent; i++)
+ fDEParentsList[i] = 0;
+
AliRunLoader *rl = AliRunLoader::GetRunLoader();
- fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
- fGeom->GetTransformationForSM(); // Global <-> Local
+ if (rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL"))
+ fGeomPtr = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
+ else
+ fGeomPtr = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaulGeometryName());
+ fGeomPtr->GetTransformationForSM(); // Global <-> Local
+}
+
+//____________________________________________________________________________
+AliEMCALRecPoint::AliEMCALRecPoint(const AliEMCALRecPoint & rp)
+ : AliRecPoint(rp),
+ fGeomPtr(rp.fGeomPtr),
+ fClusterType(rp.fClusterType),
+ fCoreEnergy(rp.fCoreEnergy),
+ fDispersion(rp.fDispersion),
+ fEnergyList(0),
+ fTimeList(0),
+ fAbsIdList(0),
+ fTime(rp.fTime),
+ fCoreRadius(rp.fCoreRadius),
+ fDETracksList(0),
+ fMulParent(rp.fMulParent),
+ fMaxParent(rp.fMaxParent),
+ fParentsList(0),
+ fDEParentsList(0),
+ fSuperModuleNumber(rp.fSuperModuleNumber)
+{
+ //copy ctor
+ fLambda[0] = rp.fLambda[0];
+ fLambda[1] = rp.fLambda[1];
+
+ fEnergyList = new Float_t[rp.fMaxDigit];
+ fTimeList = new Float_t[rp.fMaxDigit];
+ fAbsIdList = new Int_t[rp.fMaxDigit];
+ for(Int_t i = 0; i < rp.fMulDigit; i++) {
+ fEnergyList[i] = rp.fEnergyList[i];
+ fTimeList[i] = rp.fTimeList[i];
+ fAbsIdList[i] = rp.fAbsIdList[i];
+ }
+ fDETracksList = new Float_t[rp.fMaxTrack];
+ for(Int_t i = 0; i < rp.fMulTrack; i++) fDETracksList[i] = rp.fDETracksList[i];
+ fParentsList = new Int_t[rp.fMaxParent];
+ for(Int_t i = 0; i < rp.fMulParent; i++) fParentsList[i] = rp.fParentsList[i];
+ fDEParentsList = new Float_t[rp.fMaxParent];
+ for(Int_t i = 0; i < rp.fMulParent; i++) fDEParentsList[i] = rp.fDEParentsList[i];
+
}
//____________________________________________________________________________
AliEMCALRecPoint::~AliEMCALRecPoint()
delete[] fTimeList ;
if ( fAbsIdList )
delete[] fAbsIdList ;
+ if ( fDETracksList)
+ delete[] fDETracksList;
if ( fParentsList)
delete[] fParentsList;
+ if ( fDEParentsList)
+ delete[] fDEParentsList;
}
//____________________________________________________________________________
fTimeList = new Float_t[fMaxDigit];
if(fAbsIdList == 0) {
fAbsIdList = new Int_t[fMaxDigit];
- fSuperModuleNumber = fGeom->GetSuperModuleNumber(digit.GetId());
+ fSuperModuleNumber = fGeomPtr->GetSuperModuleNumber(digit.GetId());
}
if ( fMulDigit >= fMaxDigit ) { // increase the size of the lists
// A neighbour is defined as being two digits which share a corner
static Bool_t areNeighbours = kFALSE ;
- static Int_t nSupMod=0, nTower=0, nIphi=0, nIeta=0;
- static int nSupMod1=0, nTower1=0, nIphi1=0, nIeta1=0;
+ static Int_t nSupMod=0, nModule=0, nIphi=0, nIeta=0;
+ static int nSupMod1=0, nModule1=0, nIphi1=0, nIeta1=0;
static Int_t relid1[2] , relid2[2] ; // ieta, iphi
static Int_t rowdiff=0, coldiff=0;
areNeighbours = kFALSE ;
- fGeom->GetCellIndex(digit1->GetId(), nSupMod,nTower,nIphi,nIeta);
- fGeom->GetCellPhiEtaIndexInSModule(nSupMod,nTower,nIphi,nIeta, relid1[0],relid1[1]);
+ fGeomPtr->GetCellIndex(digit1->GetId(), nSupMod,nModule,nIphi,nIeta);
+ fGeomPtr->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, relid1[0],relid1[1]);
- fGeom->GetCellIndex(digit2->GetId(), nSupMod1,nTower1,nIphi1,nIeta1);
- fGeom->GetCellPhiEtaIndexInSModule(nSupMod1,nTower1,nIphi1,nIeta1, relid2[0],relid2[1]);
+ fGeomPtr->GetCellIndex(digit2->GetId(), nSupMod1,nModule1,nIphi1,nIeta1);
+ fGeomPtr->GetCellPhiEtaIndexInSModule(nSupMod1,nModule1,nIphi1,nIeta1, relid2[0],relid2[1]);
rowdiff = TMath::Abs( relid1[0] - relid2[0] ) ;
coldiff = TMath::Abs( relid1[1] - relid2[1] ) ;
void AliEMCALRecPoint::EvalDispersion(Float_t logWeight, TClonesArray * digits)
{
// Calculates the dispersion of the shower at the origin of the RecPoint
+ // in cell units - Nov 16,2006
- Double_t d = 0., wtot = 0., w = 0., xyzi[3], diff=0.;
- Int_t iDigit=0, nstat=0, i=0;
+ Double_t d = 0., wtot = 0., w = 0.;
+ Int_t iDigit=0, nstat=0;
AliEMCALDigit * digit ;
-
- // Calculates the centre of gravity in the local EMCAL-module coordinates
- if (!fLocPos.Mag())
- EvalLocalPosition(logWeight, digits) ;
- // Calculates the dispersion in coordinates
+ // Calculates the dispersion in cell units
+ Double_t etai, phii, etaMean=0.0, phiMean=0.0;
+ int nSupMod=0, nModule=0, nIphi=0, nIeta=0;
+ int iphi=0, ieta=0;
+ // Calculate mean values
for(iDigit=0; iDigit < fMulDigit; iDigit++) {
digit = (AliEMCALDigit *) digits->At(fDigitsList[iDigit]) ;
- fGeom->RelPosCellInSModule(digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
- w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
+ if (fAmp>0 && fEnergyList[iDigit]>0) {
+ fGeomPtr->GetCellIndex(digit->GetId(), nSupMod,nModule,nIphi,nIeta);
+ fGeomPtr->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
+ etai=(Double_t)ieta;
+ phii=(Double_t)iphi;
+ w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
+
+ if(w>0.0) {
+ phiMean += phii*w;
+ etaMean += etai*w;
+ wtot += w;
+ }
+ }
+ }
+ if (wtot>0) {
+ phiMean /= wtot ;
+ etaMean /= wtot ;
+ } else AliError(Form("Wrong weight %f\n", wtot));
- if(w>0.0) {
- wtot += w;
- nstat++;
- for(i=0; i<3; i++ ) {
- diff = xyzi[i] - double(fLocPos[i]);
- d += w * diff*diff;
+ // Calculate dispersion
+ for(iDigit=0; iDigit < fMulDigit; iDigit++) {
+ digit = (AliEMCALDigit *) digits->At(fDigitsList[iDigit]) ;
+
+ if (fAmp>0 && fEnergyList[iDigit]>0) {
+ fGeomPtr->GetCellIndex(digit->GetId(), nSupMod,nModule,nIphi,nIeta);
+ fGeomPtr->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
+ etai=(Double_t)ieta;
+ phii=(Double_t)iphi;
+ w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
+
+ if(w>0.0) {
+ nstat++;
+ d += w*((etai-etaMean)*(etai-etaMean)+(phii-phiMean)*(phii-phiMean));
}
}
}
for(Int_t iDigit=0; iDigit<fMulDigit; iDigit++) {
digit = dynamic_cast<AliEMCALDigit *>(digits->At(fDigitsList[iDigit])) ;
- fGeom->RelPosCellInSModule(digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
+ fGeomPtr->RelPosCellInSModule(digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
// printf(" Id %i : Local x,y,z %f %f %f \n", digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
- if(logWeight > 0.0) w = TMath::Max( 0., logWeight + TMath::Log( fEnergyList[iDigit] / fAmp ));
- else w = fEnergyList[iDigit]; // just energy
+ if(logWeight > 0.0) w = TMath::Max( 0., logWeight + TMath::Log( fEnergyList[iDigit] / fAmp ));
+ else w = fEnergyList[iDigit]; // just energy
if(w>0.0) {
wtot += w ;
// i.e. within a radius rad = fCoreEnergy around the center. Beyond this radius
// in accordance with shower profile the energy deposition
// should be less than 2%
+ // Unfinished - Nov 15,2006
+ // Distance is calculate in (phi,eta) units
AliEMCALDigit * digit ;
- const Float_t kDeg2Rad = TMath::DegToRad() ;
Int_t iDigit;
EvalLocalPosition(logWeight, digits);
}
+ Double_t phiPoint = fLocPos.Phi(), etaPoint = fLocPos.Eta();
+ Double_t eta, phi, distance;
for(iDigit=0; iDigit < fMulDigit; iDigit++) {
digit = (AliEMCALDigit *) ( digits->At(fDigitsList[iDigit]) ) ;
- Float_t ieta = 0.0;
- Float_t iphi = 0.0;
- //fGeom->PosInAlice(digit->GetId(), ieta, iphi);
- fGeom->EtaPhiFromIndex(digit->GetId(),ieta, iphi) ;
- iphi = iphi * kDeg2Rad;
+ eta = phi = 0.0;
+ fGeomPtr->EtaPhiFromIndex(digit->GetId(),eta, phi) ;
+ phi = phi * TMath::DegToRad();
- Float_t distance = TMath::Sqrt((ieta-fLocPos.X())*(ieta-fLocPos.X())+(iphi-fLocPos.Y())*(iphi-fLocPos.Y())) ;
+ distance = TMath::Sqrt((eta-etaPoint)*(eta-etaPoint)+(phi-phiPoint)*(phi-phiPoint));
if(distance < fCoreRadius)
fCoreEnergy += fEnergyList[iDigit] ;
}
void AliEMCALRecPoint::EvalElipsAxis(Float_t logWeight,TClonesArray * digits)
{
// Calculates the axis of the shower ellipsoid in eta and phi
+ // in cell units
+
+ static TString gn(fGeomPtr->GetName());
- Double_t wtot = 0. ;
+ Double_t wtot = 0.;
Double_t x = 0.;
Double_t z = 0.;
Double_t dxx = 0.;
Double_t dzz = 0.;
Double_t dxz = 0.;
- const Float_t kDeg2Rad = TMath::DegToRad();
- AliEMCALDigit * digit ;
+ AliEMCALDigit * digit = 0;
- TString gn(fGeom->GetName());
-
- Int_t iDigit;
-
- for(iDigit=0; iDigit<fMulDigit; iDigit++) {
+ Double_t etai , phii, w;
+ int nSupMod=0, nModule=0, nIphi=0, nIeta=0;
+ int iphi=0, ieta=0;
+ for(Int_t iDigit=0; iDigit<fMulDigit; iDigit++) {
digit = (AliEMCALDigit *) digits->At(fDigitsList[iDigit]) ;
- Float_t etai = 0. ;
- Float_t phii = 0. ;
- if(gn.Contains("SHISH")) { // have to be change - Feb 28, 2006
- //copied for shish-kebab geometry, ieta,iphi is cast as float as eta,phi conversion
- // for this geometry does not exist
- int nSupMod=0, nTower=0, nIphi=0, nIeta=0;
- int iphi=0, ieta=0;
- fGeom->GetCellIndex(digit->GetId(), nSupMod,nTower,nIphi,nIeta);
- fGeom->GetCellPhiEtaIndexInSModule(nSupMod,nTower,nIphi,nIeta, iphi,ieta);
- etai=(Float_t)ieta;
- phii=(Float_t)iphi;
- } else {
- fGeom->EtaPhiFromIndex(digit->GetId(), etai, phii);
- phii = phii * kDeg2Rad;
+ etai = phii = 0.;
+ if(gn.Contains("SHISH")) {
+ // Nov 15,2006 - use cell numbers as coordinates
+ // Copied for shish-kebab geometry, ieta,iphi is cast as double as eta,phi
+ // We can use the eta,phi(or coordinates) of cell
+ nSupMod = nModule = nIphi = nIeta = iphi = ieta = 0;
+
+ fGeomPtr->GetCellIndex(digit->GetId(), nSupMod,nModule,nIphi,nIeta);
+ fGeomPtr->GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
+ etai=(Double_t)ieta;
+ phii=(Double_t)iphi;
+ } else { //
+ fGeomPtr->EtaPhiFromIndex(digit->GetId(), etai, phii);
+ phii = phii * TMath::DegToRad();
}
- Double_t w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
+ w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
+ // fAmp summed amplitude of digits, i.e. energy of recpoint
+ // Gives smaller value of lambda than log weight
+ // w = fEnergyList[iDigit] / fAmp; // Nov 16, 2006 - try just energy
dxx += w * etai * etai ;
x += w * etai ;
//______________________________________________________________________________
void AliEMCALRecPoint::EvalPrimaries(TClonesArray * digits)
{
- // Constructs the list of primary particles (tracks) which have contributed to this RecPoint
+ // Constructs the list of primary particles (tracks) which
+ // have contributed to this RecPoint and calculate deposited energy
+ // for each track
AliEMCALDigit * digit ;
- Int_t * tempo = new Int_t[fMaxTrack] ;
+ Int_t * primArray = new Int_t[fMaxTrack] ;
+ Float_t * dEPrimArray = new Float_t[fMaxTrack] ;
Int_t index ;
for ( index = 0 ; index < GetDigitsMultiplicity() ; index++ ) { // all digits
digit = dynamic_cast<AliEMCALDigit *>(digits->At( fDigitsList[index] )) ;
Int_t nprimaries = digit->GetNprimary() ;
if ( nprimaries == 0 ) continue ;
- Int_t * newprimaryarray = new Int_t[nprimaries] ;
- Int_t ii ;
- for ( ii = 0 ; ii < nprimaries ; ii++)
- newprimaryarray[ii] = digit->GetPrimary(ii+1) ;
-
Int_t jndex ;
for ( jndex = 0 ; jndex < nprimaries ; jndex++ ) { // all primaries in digit
if ( fMulTrack > fMaxTrack ) {
Error("GetNprimaries", "increase fMaxTrack ") ;
break ;
}
- Int_t newprimary = newprimaryarray[jndex] ;
+ Int_t newPrimary = digit->GetPrimary(jndex+1);
+ Float_t dEPrimary = digit->GetDEPrimary(jndex+1);
Int_t kndex ;
Bool_t already = kFALSE ;
for ( kndex = 0 ; kndex < fMulTrack ; kndex++ ) { //check if not already stored
- if ( newprimary == tempo[kndex] ){
+ if ( newPrimary == primArray[kndex] ){
already = kTRUE ;
+ dEPrimArray[kndex] += dEPrimary;
break ;
}
} // end of check
if ( !already && (fMulTrack < fMaxTrack)) { // store it
- tempo[fMulTrack] = newprimary ;
+ primArray[fMulTrack] = newPrimary ;
+ dEPrimArray[fMulTrack] = dEPrimary ;
fMulTrack++ ;
} // store it
} // all primaries in digit
- delete [] newprimaryarray ;
} // all digits
-
- fTracksList = new Int_t[fMulTrack] ;
- for(index = 0; index < fMulTrack; index++)
- fTracksList[index] = tempo[index] ;
-
- delete [] tempo ;
+ Int_t *sortIdx = new Int_t[fMulTrack];
+ TMath::Sort(fMulTrack,dEPrimArray,sortIdx);
+ for(index = 0; index < fMulTrack; index++) {
+ fTracksList[index] = primArray[sortIdx[index]] ;
+ fDETracksList[index] = dEPrimArray[sortIdx[index]] ;
+ }
+ delete [] sortIdx;
+ delete [] primArray ;
+ delete [] dEPrimArray ;
}
// Constructs the list of parent particles (tracks) which have contributed to this RecPoint
AliEMCALDigit * digit ;
- Int_t * tempo = new Int_t[fMaxParent] ;
+ Int_t * parentArray = new Int_t[fMaxTrack] ;
+ Float_t * dEParentArray = new Float_t[fMaxTrack] ;
Int_t index ;
for ( index = 0 ; index < GetDigitsMultiplicity() ; index++ ) { // all digits
digit = dynamic_cast<AliEMCALDigit *>(digits->At( fDigitsList[index] )) ;
Int_t nparents = digit->GetNiparent() ;
if ( nparents == 0 ) continue ;
- Int_t * newparentarray = new Int_t[nparents] ;
- Int_t ii ;
- for ( ii = 0 ; ii < nparents ; ii++)
- newparentarray[ii] = digit->GetIparent(ii+1) ;
Int_t jndex ;
for ( jndex = 0 ; jndex < nparents ; jndex++ ) { // all primaries in digit
Error("GetNiparent", "increase fMaxParent") ;
break ;
}
- Int_t newparent = newparentarray[jndex] ;
+ Int_t newParent = digit->GetIparent(jndex+1) ;
+ Float_t newdEParent = digit->GetDEParent(jndex+1) ;
Int_t kndex ;
Bool_t already = kFALSE ;
for ( kndex = 0 ; kndex < fMulParent ; kndex++ ) { //check if not already stored
- if ( newparent == tempo[kndex] ){
+ if ( newParent == parentArray[kndex] ){
+ dEParentArray[kndex] += newdEParent;
already = kTRUE ;
break ;
}
} // end of check
if ( !already && (fMulTrack < fMaxTrack)) { // store it
- tempo[fMulParent] = newparent ;
+ parentArray[fMulParent] = newParent ;
+ dEParentArray[fMulParent] = newdEParent ;
fMulParent++ ;
} // store it
} // all parents in digit
- delete [] newparentarray ;
} // all digits
if (fMulParent>0) {
- fParentsList = new Int_t[fMulParent] ;
- for(index = 0; index < fMulParent; index++)
- fParentsList[index] = tempo[index] ;
+ Int_t *sortIdx = new Int_t[fMulParent];
+ TMath::Sort(fMulParent,dEParentArray,sortIdx);
+ for(index = 0; index < fMulParent; index++) {
+ fParentsList[index] = parentArray[sortIdx[index]] ;
+ fDEParentsList[index] = dEParentArray[sortIdx[index]] ;
+ }
+ delete [] sortIdx;
}
- delete [] tempo ;
-
+ delete [] parentArray;
+ delete [] dEParentArray;
}
//____________________________________________________________________________
// returns the position of the cluster in the global reference system of ALICE
// These are now the Cartesian X, Y and Z
// cout<<" geom "<<geom<<endl;
- fGeom->GetGlobal(fLocPos, gpos, fSuperModuleNumber);
+ fGeomPtr->GetGlobal(fLocPos, gpos, fSuperModuleNumber);
}
//____________________________________________________________________________
Int_t AliEMCALRecPoint::GetPrimaryIndex() const
{
// Get the primary track index in TreeK which deposits the most energy
- // in Digits which forms RecPoint. Kinematics, Hits and Digits must be
- // loaded before the call of the method.
-
- AliRunLoader *rl = AliRunLoader::GetRunLoader();
- if (!rl)
- AliError(Form(" No Runloader ")) ;
-
- AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>
- (rl->GetDetectorLoader("EMCAL"));
-
- // Get the list of digits forming this RecPoint
- Int_t nDigits = fMulDigit ;
- Int_t *digitList = fDigitsList ;
-
- // Find the digit with maximum amplitude
- AliEMCALDigit *digit = 0;
- TClonesArray *digits = emcalLoader->Digits();
- Int_t maxAmp = 0;
- Int_t bestDigitIndex = -1;
- for (Int_t iDigit=0; iDigit<nDigits; iDigit++) {
- digit = static_cast<AliEMCALDigit *>(digits->At(digitList[iDigit]));
- if (digit->GetAmp() > maxAmp) {
- maxAmp = digit->GetAmp();
- bestDigitIndex = iDigit;
- }
- }
+ // in Digits which forms RecPoint.
- digit = static_cast<AliEMCALDigit *>(digits->At(digitList[bestDigitIndex]));
-
- // Get the list of hits producing this digit,
- // find which hit has deposited more energy
- // and find the primary track.
-
- AliEMCALHit *hit = 0;
- TClonesArray *hits = emcalLoader->Hits();
-
- Double_t maxedep = 0;
- Int_t maxtrack = -1;
- Int_t nHits = hits ->GetEntries();
- Int_t id = digit->GetId();
- for (Int_t iHit=0; iHit<nHits; iHit++) {
- hit = static_cast<AliEMCALHit*> (hits->At(iHit)) ;
- if(hit->GetId() == id){
- Double_t edep = hit->GetEnergy();
- Int_t track = hit->GetIparent();//Primary();
- if(edep > maxedep){
- maxedep = edep;
- maxtrack = track;
- }
- }
- }
- if (maxtrack != -1) return maxtrack;
- return -12345; // no track found :(
+ if (fMulTrack)
+ return fTracksList[0];
+ return -12345;
}
//____________________________________________________________________________
message += " Stored at position %d" ;
Info("Print", message.Data(), fClusterType, fMulDigit, fAmp, fCoreEnergy, fCoreRadius, fMulTrack, GetIndexInList() ) ;
}
+
+Double_t AliEMCALRecPoint::GetPointEnergy() const
+{
+ static double e;
+ e=0.0;
+ for(int ic=0; ic<GetMultiplicity(); ic++) e += double(fEnergyList[ic]);
+ return e;
+}