#include "AliEMCALRawUtils.h"
#include "AliEMCALDigit.h"
#include "AliEMCALClusterizerv1.h"
+#include "AliEMCALClusterizerv2.h"
#include "AliEMCALClusterizerNxN.h"
#include "AliEMCALRecPoint.h"
#include "AliEMCALPID.h"
-#include "AliEMCALTrigger.h"
+#include "AliEMCALRecoUtils.h"
#include "AliRawReader.h"
#include "AliCDBEntry.h"
#include "AliCDBManager.h"
: fGeom(0),fCalibData(0),fPedestalData(0),fTriggerData(0x0), fMatches(0x0)
{
// ctor
-
+
+ // AliDebug(2, "Mark.");
+
fgRawUtils = new AliEMCALRawUtils;
//To make sure we match with the geometry in a simulation file,
if(!fCalibData)
{
AliCDBEntry *entry = (AliCDBEntry*)
- AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
+ AliCDBManager::Instance()->Get("EMCAL/Calib/Data");
if (entry) fCalibData = (AliEMCALCalibData*) entry->GetObject();
}
if(!fPedestalData)
{
AliCDBEntry *entry = (AliCDBEntry*)
- AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals");
+ AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals");
if (entry) fPedestalData = (AliCaloCalibPedestal*) entry->GetObject();
}
//Init temporary list of digits
fgDigitsArr = new TClonesArray("AliEMCALDigit",1000);
fgClustersArr = new TObjArray(1000);
- fgTriggerDigits = new TClonesArray("AliEMCALTriggerRawDigit",1000);
+ const int kNTRU = fGeom->GetNTotalTRU();
+ fgTriggerDigits = new TClonesArray("AliEMCALTriggerRawDigit", kNTRU * 96);
+
//Track matching
fMatches = new TList();
fMatches->SetOwner(kTRUE);
AliEMCALReconstructor::~AliEMCALReconstructor()
{
// dtor
-
- if(fGeom) delete fGeom;
+
+ //AliDebug(2, "Mark.");
+
+ ////RS if(fGeom) delete fGeom;
//No need to delete, recovered from OCDB
//if(fCalibData) delete fCalibData;
//if(fPedestalData) delete fPedestalData;
- if(fgDigitsArr){
- fgDigitsArr->Clear("C");
- delete fgDigitsArr;
- }
+ if(fgDigitsArr) fgDigitsArr->Clear("C");
+ delete fgDigitsArr;
+ fgDigitsArr = 0;
- if(fgClustersArr){
- fgClustersArr->Clear();
- delete fgClustersArr;
- }
+ if(fgClustersArr) fgClustersArr->Clear();
+ delete fgClustersArr;
+ fgClustersArr = 0;
- if(fgTriggerDigits){
- fgTriggerDigits->Clear();
- delete fgTriggerDigits;
- }
+ if(fgTriggerDigits) fgTriggerDigits->Clear();
+ delete fgTriggerDigits;
+ fgTriggerDigits = 0;
+
+ delete fgRawUtils;
+ fgRawUtils = 0;
+ delete fgClusterizer;
+ fgClusterizer = 0;
- if(fgRawUtils) delete fgRawUtils;
- if(fgClusterizer) delete fgClusterizer;
- if(fgTriggerProcessor) delete fgTriggerProcessor;
+ delete fgTriggerProcessor;
+ fgTriggerProcessor = 0;
if(fMatches) { fMatches->Delete(); delete fMatches; fMatches = 0;}
//printf("ReCreate clusterizer? Clusterizer set <%d>, Clusterizer in use <%s>\n",
// clusterizerType, fgClusterizer->Version());
- if (clusterizerType == AliEMCALRecParam::kClusterizerv1 && !strcmp(fgClusterizer->Version(),"clu-v1")) return;
+ if (clusterizerType == AliEMCALRecParam::kClusterizerv1 && !strcmp(fgClusterizer->Version(),"clu-v1")) return;
else if(clusterizerType == AliEMCALRecParam::kClusterizerNxN && !strcmp(fgClusterizer->Version(),"clu-NxN")) return;
+ else if(clusterizerType == AliEMCALRecParam::kClusterizerv2 && !strcmp(fgClusterizer->Version(),"clu-v2")) return;
+
//Need to create new clusterizer, the one set previously is not the correct one
delete fgClusterizer;
}
else return;
}
- if (clusterizerType == AliEMCALRecParam::kClusterizerv1)
+ if (clusterizerType == AliEMCALRecParam::kClusterizerv1)
{
- fgClusterizer = new AliEMCALClusterizerv1(fGeom, fCalibData,fPedestalData);
+ fgClusterizer = new AliEMCALClusterizerv1 (fGeom, fCalibData,fPedestalData);
}
- else
+ else if (clusterizerType == AliEMCALRecParam::kClusterizerNxN)
{
fgClusterizer = new AliEMCALClusterizerNxN(fGeom, fCalibData,fPedestalData);
}
+ else if (clusterizerType == AliEMCALRecParam::kClusterizerv2)
+ {
+ fgClusterizer = new AliEMCALClusterizerv2 (fGeom, fCalibData,fPedestalData);
+ }
+ else
+ {
+ AliFatal(Form("Unknown clusterizer %d ", clusterizerType));
+ }
}
//____________________________________________________________________________
}//not a LED event
clustersTree->Fill();
+
+ // Deleting the recpoints at the end of the reconstruction call
+ fgClusterizer->DeleteRecPoints();
}
//____________________________________________________________________________
if(fgDigitsArr) fgDigitsArr->Clear("C");
- TClonesArray *digitsTrg = new TClonesArray("AliEMCALTriggerRawDigit", 32 * 96);
+ const int kNTRU = fGeom->GetNTotalTRU();
+ TClonesArray *digitsTrg = new TClonesArray("AliEMCALTriggerRawDigit", kNTRU * 96);
Int_t bufsize = 32000;
digitsTree->Branch("EMCAL", &fgDigitsArr, bufsize);
fgRawUtils->SetNoiseThreshold(GetRecParam()->GetNoiseThreshold());
fgRawUtils->SetNPedSamples(GetRecParam()->GetNPedSamples());
fgRawUtils->SetRemoveBadChannels(GetRecParam()->GetRemoveBadChannels());
- fgRawUtils->SetFittingAlgorithm(GetRecParam()->GetFittingAlgorithm());
+ if (!fgRawUtils->GetFittingAlgorithm()) fgRawUtils->SetFittingAlgorithm(GetRecParam()->GetFittingAlgorithm());
fgRawUtils->SetFALTROUsage(GetRecParam()->UseFALTRO());
-
+ // fgRawUtils->SetFALTROUsage(0);
+
//fgRawUtils->SetTimeMin(GetRecParam()->GetTimeMin());
//fgRawUtils->SetTimeMax(GetRecParam()->GetTimeMax());
if (esdV0)
{
- for (Int_t i = 0; i < 32; i++)
- {
- v0M[0] += (Int_t)esdV0->GetAdcV0C(i);
- v0M[1] += (Int_t)esdV0->GetAdcV0A(i);
- }
+ v0M[0] = esdV0->GetTriggerChargeC();
+ v0M[1] = esdV0->GetTriggerChargeA();
}
else
{
- AliWarning("Cannot retrieve V0 ESD! Run w/ null V0 charges");
+ AliWarning("No V0 ESD! Run trigger processor w/ null V0 charges");
}
- if (fgTriggerDigits) fgTriggerDigits->Clear();
+ if (fgTriggerDigits && fgTriggerDigits->GetEntriesFast()) fgTriggerDigits->Delete();
TBranch *branchtrg = digitsTree->GetBranch("EMTRG");
for (Int_t i = 0; i < fgTriggerDigits->GetEntriesFast(); i++)
{
AliEMCALTriggerRawDigit* rdig = (AliEMCALTriggerRawDigit*)fgTriggerDigits->At(i);
-
+ if (AliDebugLevel() > 999) rdig->Print("");
+
Int_t px, py;
if (fGeom->GetPositionInEMCALFromAbsFastORIndex(rdig->GetId(), px, py))
{
rdig->GetMaximum(a, t);
rdig->GetL0Times(times);
-
+
trgESD->Add(px, py, a, t, times, rdig->GetNL0Times(), rdig->GetL1TimeSum(), rdig->GetTriggerBits());
}
}
- trgESD->SetL1Threshold(0, fTriggerData->GetL1GammaThreshold());
-
- trgESD->SetL1Threshold(1, fTriggerData->GetL1JetThreshold() );
+ for (int i = 0; i < 2; i++) {
+ trgESD->SetL1Threshold(2 * i , fTriggerData->GetL1JetThreshold( i));
+ trgESD->SetL1Threshold(2 * i + 1, fTriggerData->GetL1GammaThreshold(i));
+ }
Int_t v0[2];
fTriggerData->GetL1V0(v0);
if (!saveOnce && fTriggerData->GetL1DataDecoded())
{
- int type[8] = {0};
+ int type[15] = {0};
fTriggerData->GetL1TriggerType(type);
esd->SetCaloTriggerType(type);
AliESDCaloCells &emcCells = *(esd->GetEMCALCells());
emcCells.CreateContainer(nDigits);
emcCells.SetType(AliVCaloCells::kEMCALCell);
+
Float_t energy = 0;
- for (Int_t idig = 0 ; idig < nDigits ; idig++) {
+ Float_t time = 0;
+ for (Int_t idig = 0 ; idig < nDigits ; idig++)
+ {
const AliEMCALDigit * dig = (const AliEMCALDigit*)fgDigitsArr->At(idig);
- if(dig->GetAmplitude() > 0 ){
- energy = fgClusterizer->Calibrate(dig->GetAmplitude(),dig->GetTime(),dig->GetId()); //TimeR or Time?
- if(energy > 0){ //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them
- emcCells.SetCell(idignew,dig->GetId(),energy, dig->GetTime());
+ time = dig->GetTime(); // Time already calibrated in clusterizer
+ energy = dig->GetAmplitude(); // energy calibrated in clusterizer
+
+ if(energy > 0 )
+ {
+ fgClusterizer->Calibrate(energy,time,dig->GetId()); //Digits already calibrated in clusterizers
+
+ if(energy > 0) //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them
+ {
+ // Only for MC
+ // Get the label of the primary particle that generated the cell
+ // Assign the particle that deposited more energy
+ Int_t nprimaries = dig->GetNprimary() ;
+ Int_t digLabel =-1 ;
+ Float_t edep =-1.;
+ if ( nprimaries > 0 )
+ {
+ Int_t jndex ;
+ for ( jndex = 0 ; jndex < nprimaries ; jndex++ ) { // all primaries in digit
+
+ if(edep < dig->GetDEParent(jndex+1))
+ {
+ digLabel = dig->GetIparent (jndex+1);
+ edep = dig->GetDEParent(jndex+1);
+ }
+
+ } // all primaries in digit
+ } // select primary label
+
+ Bool_t highGain = kFALSE;
+ if( dig->GetType() == AliEMCALDigit::kHG ) highGain = kTRUE;
+
+ emcCells.SetCell(idignew,dig->GetId(),energy, time,digLabel,0.,highGain);
idignew++;
}
}
}
+
emcCells.SetNumberOfCells(idignew);
emcCells.Sort();
AliESDtrack * track = esd->GetTrack(itrack) ; // retrieve track
if(track->GetEMCALcluster()==iClust)
{
- Double_t dEta=-999, dPhi=-999;
+ Float_t dEta=-999, dPhi=-999;
Bool_t isMatch = CalculateResidual(track, ec, dEta, dPhi);
if(!isMatch)
{
+ // AliDebug(10, "Not good");
continue;
- cout<<"Not good"<<endl;
}
AliEMCALMatch *match = new AliEMCALMatch();
match->SetIndexT(itrack);
const Int_t bufsize = 255;
char path[bufsize] ;
TGeoHMatrix * m = 0x0;
+ Int_t tmpType = -1;
+ Int_t SMOrder = 0;
+ TString SMName;
for(Int_t sm = 0; sm < fGeom->GetNumberOfSuperModules(); sm++){
- snprintf(path,bufsize,"/ALIC_1/XEN1_1/SMOD_%d",sm+1) ; //In Geometry modules numbered 1,2,.,5
- if(sm >= 10) snprintf(path,bufsize,"/ALIC_1/XEN1_1/SM10_%d",sm-10+1) ;
-
+ if(fGeom->GetSMType(sm) == AliEMCALGeometry::kEMCAL_Standard ) SMName = "SMOD";
+ else if(fGeom->GetSMType(sm) == AliEMCALGeometry::kEMCAL_Half ) SMName = "SM10";
+ else if(fGeom->GetSMType(sm) == AliEMCALGeometry::kEMCAL_3rd ) SMName = "SM3rd";
+ else if( fGeom->GetSMType(sm) == AliEMCALGeometry::kDCAL_Standard ) SMName = "DCSM";
+ else if( fGeom->GetSMType(sm) == AliEMCALGeometry::kDCAL_Ext ) SMName = "DCEXT";
+ else AliError("Unkown SM Type!!");
+
+ if(fGeom->GetSMType(sm) == tmpType) {
+ SMOrder++;
+ } else {
+ tmpType = fGeom->GetSMType(sm);
+ SMOrder = 1;
+ }
+ snprintf(path,bufsize,"/ALIC_1/XEN1_1/%s_%d", SMName.Data(), SMOrder) ;
+
if (gGeoManager->CheckPath(path)){
gGeoManager->cd(path);
m = gGeoManager->GetCurrentMatrix() ;
}
//==================================================================================
-Bool_t AliEMCALReconstructor::CalculateResidual(AliESDtrack *track, AliESDCaloCluster *cluster, Double_t &dEta, Double_t &dPhi)const
+Bool_t AliEMCALReconstructor::CalculateResidual(AliESDtrack *track, AliESDCaloCluster *cluster, Float_t &dEta, Float_t &dPhi)const
{
//
// calculate the residual between track and cluster
// If the esdFriend is available, use the TPCOuter point as the starting point of extrapolation
// Otherwise use the TPCInner point
- AliExternalTrackParam *trkParam;
+
+ dEta = -999, dPhi = -999;
+ Bool_t ITSTrackSA = 0;
+
+ AliExternalTrackParam *trkParam = 0;
+
const AliESDfriendTrack* friendTrack = track->GetFriendTrack();
if(friendTrack && friendTrack->GetTPCOut())
trkParam = const_cast<AliExternalTrackParam*>(friendTrack->GetTPCOut());
- else
+ else if(track->GetInnerParam())
trkParam = const_cast<AliExternalTrackParam*>(track->GetInnerParam());
+ else{
+ trkParam = new AliExternalTrackParam(*track); //If there is ITSSa track
+ ITSTrackSA = 1;
+ }
if(!trkParam) return kFALSE;
+
+ AliExternalTrackParam trkParamTmp (*trkParam);
+ if(!AliEMCALRecoUtils::ExtrapolateTrackToCluster(&trkParamTmp, cluster, track->GetMass(kTRUE), GetRecParam()->GetExtrapolateStep(), dEta, dPhi)){
+ if(ITSTrackSA) delete trkParam;
+ return kFALSE;
+ }
- //Perform extrapolation
- Double_t trkPos[3];
- Float_t clsPos[3];
-
- AliExternalTrackParam *trkParamTmp = new AliExternalTrackParam(*trkParam);
- cluster->GetPosition(clsPos);
- TVector3 vec(clsPos[0],clsPos[1],clsPos[2]);
- Double_t alpha = ((int)(vec.Phi()*TMath::RadToDeg()/20)+0.5)*20*TMath::DegToRad();
- //Rotate to proper local coordinate
- vec.RotateZ(-alpha);
- trkParamTmp->Rotate(alpha);
- //extrapolation is done here
- if(!AliTrackerBase::PropagateTrackToBxByBz(trkParamTmp, vec.X(), track->GetMass(), GetRecParam()->GetExtrapolateStep(), kFALSE))
- return kFALSE;
-
- //Calculate the residuals
- trkParamTmp->GetXYZ(trkPos);
- delete trkParamTmp;
-
- TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
- TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
-
- Double_t clsPhi = clsPosVec.Phi();
- if(clsPhi<0) clsPhi+=2*TMath::Pi();
- Double_t trkPhi = trkPosVec.Phi();
- if(trkPhi<0) trkPhi+=2*TMath::Pi();
-
- dPhi = clsPhi-trkPhi;
- dEta = clsPosVec.Eta()-trkPosVec.Eta();
-
+ if(ITSTrackSA) delete trkParam;
return kTRUE;
}
{
//copy ctor
}
+//_____________________________________________________________________
+AliEMCALReconstructor::AliEMCALMatch& AliEMCALReconstructor::AliEMCALMatch::AliEMCALMatch::operator = (const AliEMCALMatch &source)
+{ // assignment operator; use copy ctor
+ if (&source == this) return *this;
+ new (this) AliEMCALMatch(source);
+ return *this;
+}
//
//==================================================================================
//