{
// Default ctor
fName="EMCAL";
- //fQATask = 0;
- fTreeQA = 0;
fGeom = 0 ;
}
AliEMCAL::AliEMCAL(const char* name, const char* title): AliDetector(name,title)
{
// ctor : title is used to identify the layout
-
- //fQATask = 0;
- fTreeQA = 0;
fGeom = 0;
}
AliMaterial(3, "Al$", 26.98, 13., 2.7, 8.9, 999., 0, 0) ;
// --- Absorption length is ignored ^
- // --- Copper ---
- AliMaterial(4, "Cu$", 63.546, 29, 8.96, 1.43, 14.8, 0, 0) ;
- // --- Absorption length is ignored ^
-
-
// DEFINITION OF THE TRACKING MEDIA
// for EMCAL: idtmed[1599->1698] equivalent to fIdtmed[0->100]
AliMedium(1, "Lead $", 1, 0,
isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
-
// The scintillator of the CPV made of Polystyrene scintillator -> idtmed[1601]
AliMedium(2, "CPV scint. $", 2, 1,
isxfld, sxmgmx, 10.0, 0.001, 0.1, 0.001, 0.001, 0, 0) ;
AliMedium(3, "Al parts $", 3, 0,
isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
- // Copper for HCal (post shower) -> idtmed[1603]
- AliMedium(4, "Copper $", 4, 0,
- isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
-
-
// --- Set decent energy thresholds for gamma and electron tracking
gMC->Gstpar(idtmed[1602], "DCUTE",0.00001) ;
gMC->Gstpar(idtmed[1602], "DCUTM",0.00001) ;
-// --- in copper parts ---
- gMC->Gstpar(idtmed[1603], "LOSS",3.) ;
- gMC->Gstpar(idtmed[1603], "DRAY",1.) ;
- gMC->Gstpar(idtmed[1603], "DCUTE",0.00001) ;
- gMC->Gstpar(idtmed[1603], "DCUTM",0.00001) ;
-
-
-
// --- and finally thresholds for photons and electrons in the scintillator ---
gMC->Gstpar(idtmed[1601],"CUTGAM",0.00008) ;
gMC->Gstpar(idtmed[1601],"CUTELE",0.001) ;
gMC->Gstpar(idtmed[1601],"BCUTE",0.0001) ;
-
}
//____________________________________________________________________________
{
if (fHits == 0x0)
fHits= new TClonesArray("AliEMCALHit",1000);
- //Info("SetTreeAddress","<%s> Setting Hits Address",GetName());
branch->SetAddress(&fHits);
}
else
}
}
}
-//____________________________________________________________________________
-void AliEMCAL::WriteQA()
-{
-
- // Make TreeQA in the output file.
-
- if(fTreeQA == 0)
- fTreeQA = new TTree("TreeQA", "QA Alarms") ;
- // Create Alarms branches
- Int_t bufferSize = 32000 ;
- Int_t splitlevel = 0 ;
-
- TFolder* topfold = GetLoader()->GetTopFolder(); //get top aliroot folder; skowron
- TString emcalqafn(AliConfig::Instance()->GetQAFolderName()+"/"); //get name of QAaut folder relative to top event; skowron
- emcalqafn+=GetName(); //hard wired string!!! add the detector name to the pathname; skowron
- TFolder * alarmsF = (TFolder*)topfold->FindObjectAny(emcalqafn); //get the folder
-
- if (alarmsF == 0x0)
- {
- Error("WriteQA","Can not find folder with qa alarms");
- return;
- }
- TString branchName(alarmsF->GetName());
- TBranch * alarmsBranch = fTreeQA->Branch(branchName,"TFolder", &alarmsF, bufferSize, splitlevel);
- TString branchTitle = branchName + " QA alarms" ;
- alarmsBranch->SetTitle(branchTitle);
- alarmsBranch->Fill() ;
-
- //fTreeQA
-}
//____________________________________________________________________________
AliLoader* AliEMCAL::MakeLoader(const char* topfoldername)
Fatal("cpy ctor", "not implemented") ;
}
virtual ~AliEMCAL() ;
- virtual void AddHit(Int_t, Int_t*, Float_t *) {
+ virtual void AddHit(Int_t, Int_t*, Float_t *) const{
Fatal("AddHit(Int_t, Int_t*, Float_t *", "not to be used: use AddHit( Int_t shunt, Int_t primary, Int_t track,Int_t id, Float_t *hits )") ;
}
virtual void CreateMaterials() ;
- virtual void FinishRun() {WriteQA();}
+ virtual void FinishRun() {}
virtual AliEMCALGeometry * GetGeometry() const ;
- virtual Int_t IsVersion(void) const = 0 ;
- //AliEMCALQAChecker * QAChecker() const {return fQATask;}
- virtual void SetTreeAddress() ;
- virtual TTree * TreeQA() const {return fTreeQA; }
- virtual const TString Version() const {return TString(" ") ; }
- virtual void WriteQA() ;
+ virtual Int_t IsVersion(void) const = 0 ;
+ virtual void SetTreeAddress() ;
+ virtual const TString Version() const {return TString(" ") ; }
AliEMCAL & operator = (const AliEMCAL & /*rvalue*/) {
Fatal("operator =", "not implemented") ; return *this ; }
virtual AliDigitizer* CreateDigitizer(AliRunDigitizer* manager) const;
protected:
-
- //AliEMCALQAChecker * fQATask ; //! PHOS checkers container
- TTree * fTreeQA ; // the QA tree that contains the alarms
AliEMCALGeometry * fGeom ; // the geometry object
- ClassDef(AliEMCAL,4) // Electromagnetic calorimeter (base class)
+ ClassDef(AliEMCAL,5) // Electromagnetic calorimeter (base class)
} ;
virtual Float_t GetTowerLocalMaxCut()const {Warning("GetTowerLocalMaxCut", "Not Defined") ; return 0. ; }
virtual Float_t GetTowerLogWeight()const {Warning("GetTowerLogWeight", "Not Defined") ; return 0. ; }
virtual Float_t GetTimeGate() const {Warning("GetTimeGate", "Not Defined") ; return 0. ; }
- virtual Float_t GetPreShoClusteringThreshold()const {Warning("GetPreShoClusteringThreshold", "Not Defined") ; return 0. ; }
- virtual Float_t GetPreShoLocalMaxCut()const {Warning("GetPreShoLocalMaxCut", "Not Defined") ; return 0. ; }
- virtual Float_t GetPreShoLogWeight()const {Warning("GetPreShoLogWeight", "Not Defined") ; return 0. ; }
virtual const char * GetRecPointsBranch() const {Warning("GetRecPointsBranch", "Not Defined") ; return 0 ; }
virtual const Int_t GetRecPointsInRun() const {Warning("GetRecPointsInRun", "Not Defined") ; return 0 ; }
virtual const char * GetDigitsBranch() const {Warning("GetDigitsBranch", "Not Defined") ; return 0 ; }
- virtual void MakeClusters() {Warning("MakeClusters", "Not Defined") ; }
+ virtual void MakeClusters() const {Warning("MakeClusters", "Not Defined") ; }
virtual void Print(Option_t * /*option*/)const {Warning("Print", "Not Defined") ; }
virtual void SetECAClusteringThreshold(Float_t) = 0;
virtual void SetECALocalMaxCut(Float_t) = 0;
virtual void SetECALogWeight(Float_t) = 0;
- virtual void SetHCAClusteringThreshold(Float_t) = 0;
- virtual void SetHCALocalMaxCut(Float_t) = 0;
- virtual void SetHCALogWeight(Float_t) = 0;
virtual void SetTimeGate(Float_t) = 0;
- virtual void SetPREClusteringThreshold(Float_t) = 0;
- virtual void SetPRELocalMaxCut(Float_t) = 0;
- virtual void SetPRELogWeight(Float_t) = 0;
virtual void SetUnfolding(Bool_t) = 0;
virtual const char * Version() const {Warning("Version", "Not Defined") ; return 0 ; }
protected:
TString fEventFolderName ; // event folder name
- ClassDef(AliEMCALClusterizer,3) // Clusterization algorithm class
+ ClassDef(AliEMCALClusterizer,4) // Clusterization algorithm class
} ;
}
//____________________________________________________________________________
-Float_t AliEMCALClusterizerv1::Calibrate(Int_t amp, Int_t where) const
+Float_t AliEMCALClusterizerv1::Calibrate(Int_t amp) const
{
- //To be replased later by the method, reading individual parameters from the database
- // where = 0 == PRE ; where = 1 == ECAL ; where = 2 == HCAL
- if ( where == 0 ) // calibrate as PRE section
- return -fADCpedestalPRE + amp * fADCchannelPRE ;
- else if (where == 1) //calibrate as ECA section
- return -fADCpedestalECA + amp * fADCchannelECA ;
- else if (where == 2) //calibrate as HCA section
- return -fADCpedestalHCA + amp * fADCchannelHCA ;
- else
- Fatal("Calibrate", "Something went wrong!") ;
- return -9999999. ;
+ //To be replased later by the method, reading individual parameters from the database
+ return -fADCpedestalECA + amp * fADCchannelECA ;
}
//____________________________________________________________________________
GetCalibrationParameters() ;
- fNumberOfPREClusters = fNumberOfECAClusters = fNumberOfHCAClusters = 0 ;
+ fNumberOfECAClusters = 0;
MakeClusters() ;
if(strstr(option,"deb"))
PrintRecPoints(option) ;
- //increment the total number of recpoints per run
- fRecPointsInRun += gime->PRERecPoints()->GetEntriesFast() ;
+ //increment the total number of recpoints per run
fRecPointsInRun += gime->ECARecPoints()->GetEntriesFast() ;
- fRecPointsInRun += gime->HCARecPoints()->GetEntriesFast() ;
}
Unload();
if(strstr(option,"tim")){
gBenchmark->Stop("EMCALClusterizer");
- Info("Exec", "took %f seconds for Clusterizing %f seconds per event",
+ printf("Exec took %f seconds for Clusterizing %f seconds per event",
gBenchmark->GetCpuTime("EMCALClusterizer"), gBenchmark->GetCpuTime("EMCALClusterizer")/nevents ) ;
}
}
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
TClonesArray * digits = gime->Digits() ;
-
gMinuit->mncler(); // Reset Minuit's list of paramters
gMinuit->SetPrintLevel(-1) ; // No Printout
gMinuit->SetFCN(AliEMCALClusterizerv1::UnfoldingChiSquare) ;
for(iDigit = 0; iDigit < nDigits; iDigit++){
digit = maxAt[iDigit];
- Int_t relid[4] ;
+ Int_t relid[3] ;
Float_t x = 0.;
Float_t z = 0.;
geom->AbsToRelNumbering(digit->GetId(), relid) ;
//____________________________________________________________________________
void AliEMCALClusterizerv1::GetCalibrationParameters()
{
+ // Gets the parameters for the calibration from the digitizer
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
if ( !gime->Digitizer() )
gime->LoadDigitizer();
- AliEMCALDigitizer * dig = gime->Digitizer();
-
- fADCchannelPRE = dig->GetPREchannel() ;
- fADCpedestalPRE = dig->GetPREpedestal() ;
+ AliEMCALDigitizer * dig = gime->Digitizer();
fADCchannelECA = dig->GetECAchannel() ;
fADCpedestalECA = dig->GetECApedestal();
-
- fADCchannelHCA = dig->GetHCAchannel() ;
- fADCpedestalHCA = dig->GetHCApedestal();
}
//____________________________________________________________________________
//____________________________________________________________________________
void AliEMCALClusterizerv1::InitParameters()
-{
- fNumberOfPREClusters = fNumberOfECAClusters = fNumberOfHCAClusters = 0 ;
- fPREClusteringThreshold = 0.0001; // must be adjusted according to the noise leve set by digitizer
+{
+ // Initializes the parameters for the Clusterizer
+ fNumberOfECAClusters = 0;
fECAClusteringThreshold = 0.0045; // must be adjusted according to the noise leve set by digitizer
- fHCAClusteringThreshold = 0.001; // must be adjusted according to the noise leve set by digitizer
- fPRELocMaxCut = 0.03 ;
fECALocMaxCut = 0.03 ;
- fHCALocMaxCut = 0.03 ;
-
- fPREW0 = 4.0 ;
fECAW0 = 4.5 ;
- fHCAW0 = 4.5 ;
-
fTimeGate = 1.e-8 ;
-
fToUnfold = kFALSE ;
-
- fRecPointsInRun = 0 ;
-
+ fRecPointsInRun = 0 ;
}
//____________________________________________________________________________
Int_t rv = 0 ;
- Int_t relid1[4] ;
+ Int_t relid1[3] ;
geom->AbsToRelNumbering(d1->GetId(), relid1) ;
- Int_t relid2[4] ;
+ Int_t relid2[3] ;
geom->AbsToRelNumbering(d2->GetId(), relid2) ;
- if ( (relid1[0] == relid2[0]) && // inside the same EMCAL Arm
- (relid1[1]==relid2[1]) ) { // and same tower section
- Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
- Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
+ if ( (relid1[0] == relid2[0])){ // inside the same EMCAL Arm
+ Int_t rowdiff = TMath::Abs( relid1[1] - relid2[1] ) ;
+ Int_t coldiff = TMath::Abs( relid1[2] - relid2[2] ) ;
if (( coldiff <= 1 ) && ( rowdiff <= 1 )){
- if((relid1[1] != 0) || (TMath::Abs(d1->GetTime() - d2->GetTime() ) < fTimeGate))
- rv = 1 ;
+ if(TMath::Abs(d1->GetTime() - d2->GetTime() ) < fTimeGate)
+ rv = 1 ;
}
else {
- if((relid2[2] > relid1[2]) && (relid2[3] > relid1[3]+1))
+ if((relid2[1] > relid1[1]) && (relid2[2] > relid1[2]+1))
rv = 2; // Difference in row numbers is too large to look further
}
}
else {
- if( (relid1[0] < relid2[0]) || (relid1[1] != relid2[1]) )
+ if(relid1[0] < relid2[0])
rv=0 ;
}
if (gDebug == 2 )
- Info("AreNeighbours", "neighbours=%d, id1=%d, relid1=%d,%d,%d,%d \n id2=%d, relid2=%d,%d,%d,%d",
- rv, d1->GetId(), relid1[0], relid1[1], relid1[2], relid1[3],
- d2->GetId(), relid2[0], relid2[1], relid2[2], relid2[3]) ;
+ printf("AreNeighbours: neighbours=%d, id1=%d, relid1=%d,%d,%d \n id2=%d, relid2=%d,%d,%d ",
+ rv, d1->GetId(), relid1[0], relid1[1], relid1[2],
+ d2->GetId(), relid2[0], relid2[1], relid2[2]) ;
return rv ;
}
//____________________________________________________________________________
void AliEMCALClusterizerv1::Unload()
{
+ // Unloads the Digits and RecPoints
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
gime->EmcalLoader()->UnloadDigits() ;
gime->EmcalLoader()->UnloadRecPoints() ;
AliEMCALGetter *gime = AliEMCALGetter::Instance() ;
- TObjArray * aPRERecPoints = gime->PRERecPoints() ;
TObjArray * aECARecPoints = gime->ECARecPoints() ;
- TObjArray * aHCARecPoints = gime->HCARecPoints() ;
TClonesArray * digits = gime->Digits() ;
TTree * treeR = gime->TreeR(); ;
Int_t index ;
- //Evaluate position, dispersion and other RecPoint properties for PRE section
- for(index = 0; index < aPRERecPoints->GetEntries(); index++)
- (dynamic_cast<AliEMCALRecPoint *>(aPRERecPoints->At(index)))->EvalAll(fPREW0,digits) ;
- aPRERecPoints->Sort() ;
-
- for(index = 0; index < aPRERecPoints->GetEntries(); index++)
- (dynamic_cast<AliEMCALRecPoint *>(aPRERecPoints->At(index)))->SetIndexInList(index) ;
-
- aPRERecPoints->Expand(aPRERecPoints->GetEntriesFast()) ;
-
//Evaluate position, dispersion and other RecPoint properties for EC section
for(index = 0; index < aECARecPoints->GetEntries(); index++)
(dynamic_cast<AliEMCALTowerRecPoint *>(aECARecPoints->At(index)))->EvalAll(fECAW0,digits) ;
aECARecPoints->Expand(aECARecPoints->GetEntriesFast()) ;
- //Evaluate position, dispersion and other RecPoint properties for HCA section
- for(index = 0; index < aHCARecPoints->GetEntries(); index++)
- (dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(index)))->EvalAll(fHCAW0,digits) ;
-
- aHCARecPoints->Sort() ;
-
- for(index = 0; index < aHCARecPoints->GetEntries(); index++)
- (dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(index)))->SetIndexInList(index) ;
-
- aHCARecPoints->Expand(aHCARecPoints->GetEntriesFast()) ;
-
Int_t bufferSize = 32000 ;
Int_t splitlevel = 0 ;
-
- //PRE section branch
- TBranch * branchPRE = treeR->Branch("EMCALPRERP","TObjArray",&aPRERecPoints,bufferSize,splitlevel);
- branchPRE->SetTitle(BranchName());
//EC section branch
TBranch * branchECA = treeR->Branch("EMCALECARP","TObjArray",&aECARecPoints,bufferSize,splitlevel);
branchECA->SetTitle(BranchName());
- //HCA section branch
- TBranch * branchHCA = treeR->Branch("EMCALHCARP","TObjArray",&aHCARecPoints,bufferSize,splitlevel);
- branchHCA->SetTitle(BranchName());
-
- branchPRE->Fill() ;
branchECA->Fill() ;
- branchHCA->Fill() ;
gime->WriteRecPoints("OVERWRITE");
gime->WriteClusterizer("OVERWRITE");
AliEMCALGeometry * geom = gime->EMCALGeometry() ;
-
- TObjArray * aPRERecPoints = gime->PRERecPoints() ;
TObjArray * aECARecPoints = gime->ECARecPoints() ;
- TObjArray * aHCARecPoints = gime->HCARecPoints() ;
- aPRERecPoints->Delete() ;
aECARecPoints->Delete() ;
- aHCARecPoints->Delete() ;
TClonesArray * digits = gime->Digits() ;
TIter next(digits) ;
AliEMCALDigit * digit ;
- Int_t ndigECA=0, ndigPRE=0, ndigHCA=0 ;
+ Int_t ndigECA=0 ;
// count the number of digits in ECA section
while ( (digit = dynamic_cast<AliEMCALDigit *>(next())) ) { // scan over the list of digits
if (geom->IsInECA(digit->GetId()))
ndigECA++ ;
- else if (geom->IsInPRE(digit->GetId()))
- ndigPRE++;
- else if (geom->IsInHCA(digit->GetId()))
- ndigHCA++;
else {
Error("MakeClusters", "id = %d is a wrong ID!", digit->GetId()) ;
abort() ;
}
}
-
- // add amplitude of PRE and ECA sections
- Int_t digECA ;
- for (digECA = 0 ; digECA < ndigECA ; digECA++) {
- digit = dynamic_cast<AliEMCALDigit *>(digits->At(digECA)) ;
- Int_t digPRE ;
- for (digPRE = ndigECA ; digPRE < ndigECA+ndigPRE ; digPRE++) {
- AliEMCALDigit * digitPRE = dynamic_cast<AliEMCALDigit *>(digits->At(digPRE)) ;
- if ( geom->AreInSameTower(digit->GetId(), digitPRE->GetId()) ){
- Float_t amp = static_cast<Float_t>(digit->GetAmp()) + geom->GetSummationFraction() * static_cast<Float_t>(digitPRE->GetAmp()) + 0.5 ;
- digit->SetAmp(static_cast<Int_t>(amp)) ;
- break ;
- }
- }
- if (gDebug)
- Info("MakeClusters","id = %d amp = %d", digit->GetId(), digit->GetAmp()) ;
- }
-
TClonesArray * digitsC = dynamic_cast<TClonesArray*>(digits->Clone()) ;
-
-
// Clusterization starts
TIter nextdigit(digitsC) ;
- Bool_t notremovedECA = kTRUE, notremovedPRE = kTRUE ;
while ( (digit = dynamic_cast<AliEMCALDigit *>(nextdigit())) ) { // scan over the list of digitsC
AliEMCALRecPoint * clu = 0 ;
- TArrayI clusterPREdigitslist(50), clusterECAdigitslist(50), clusterHCAdigitslist(50);
+ TArrayI clusterECAdigitslist(50);
- Bool_t inPRE = kFALSE, inECA = kFALSE, inHCA = kFALSE ;
- if( geom->IsInPRE(digit->GetId()) ) {
- inPRE = kTRUE ;
- }
- else if( geom->IsInECA(digit->GetId()) ) {
- inECA = kTRUE ;
- }
- else if( geom->IsInHCA(digit->GetId()) ) {
- inHCA = kTRUE ;
- }
-
+ Bool_t inECA = kFALSE;
+ if( geom->IsInECA(digit->GetId()) ) {
+ inECA = kTRUE ;
+ }
if (gDebug == 2) {
- if (inPRE)
- Info("MakeClusters","id = %d, ene = %f , thre = %f ",
- digit->GetId(),Calibrate(digit->GetAmp(), 0), fPREClusteringThreshold) ;
if (inECA)
- Info("MakeClusters","id = %d, ene = %f , thre = %f",
- digit->GetId(),Calibrate(digit->GetAmp(), 1), fECAClusteringThreshold) ;
- if (inHCA)
- Info("MakeClusters","id = %d, ene = %f , thre = %f",
- digit->GetId(),Calibrate(digit->GetAmp(), 2), fHCAClusteringThreshold ) ;
+ printf("MakeClusters: id = %d, ene = %f , thre = %f",
+ digit->GetId(),Calibrate(digit->GetAmp()), fECAClusteringThreshold) ;
}
-
- if ( (inPRE && (Calibrate(digit->GetAmp(), 0) > fPREClusteringThreshold )) ||
- (inECA && (Calibrate(digit->GetAmp(), 1) > fECAClusteringThreshold )) ||
- (inHCA && (Calibrate(digit->GetAmp(), 2) > fHCAClusteringThreshold )) ) {
-
- Int_t iDigitInPRECluster = 0, iDigitInECACluster = 0, iDigitInHCACluster = 0;
- Int_t where ; // PRE = 0, ECAl = 1, HCAL = 2
-
- // Find the seed in each of the section ECAL/PRE/HCAL
+ if (inECA && (Calibrate(digit->GetAmp()) > fECAClusteringThreshold ) ){
+ Int_t iDigitInECACluster = 0;
+ // Find the seed
if( geom->IsInECA(digit->GetId()) ) {
- where = 1 ; // to tell we are in ECAL
// start a new Tower RecPoint
if(fNumberOfECAClusters >= aECARecPoints->GetSize())
aECARecPoints->Expand(2*fNumberOfECAClusters+1) ;
aECARecPoints->AddAt(rp, fNumberOfECAClusters) ;
clu = dynamic_cast<AliEMCALTowerRecPoint *>(aECARecPoints->At(fNumberOfECAClusters)) ;
fNumberOfECAClusters++ ;
- clu->AddDigit(*digit, Calibrate(digit->GetAmp(), where)) ;
+ clu->AddDigit(*digit, Calibrate(digit->GetAmp())) ;
clusterECAdigitslist[iDigitInECACluster] = digit->GetIndexInList() ;
iDigitInECACluster++ ;
digitsC->Remove(digit) ;
if (gDebug == 2 )
- Info("MakeClusters","OK id = %d, ene = %f , thre = %f ", digit->GetId(),Calibrate(digit->GetAmp(), 1), fECAClusteringThreshold) ;
-
- }
- else if( geom->IsInPRE(digit->GetId()) ) {
- where = 0 ; // to tell we are in PRE
- // start a new Pre Shower cluster
- if(fNumberOfPREClusters >= aPRERecPoints->GetSize())
- aPRERecPoints->Expand(2*fNumberOfPREClusters+1);
- AliEMCALTowerRecPoint * rp = new AliEMCALTowerRecPoint("") ;
- rp->SetPRE() ;
- aPRERecPoints->AddAt(rp, fNumberOfPREClusters) ;
- clu = dynamic_cast<AliEMCALTowerRecPoint *>(aPRERecPoints->At(fNumberOfPREClusters)) ;
- fNumberOfPREClusters++ ;
- clu->AddDigit(*digit, Calibrate(digit->GetAmp(), where));
- clusterPREdigitslist[iDigitInPRECluster] = digit->GetIndexInList() ;
- iDigitInPRECluster++ ;
- digitsC->Remove(digit) ;
- if (gDebug == 2 )
- Info("MakeClusters","OK id = %d, ene = %f , thre = %f ", digit->GetId(),Calibrate(digit->GetAmp(), 0), fPREClusteringThreshold) ;
-
- nextdigit.Reset() ;
-
- // Here we remove remaining ECA digits, which cannot make a cluster
-
- if( notremovedECA ) {
- while( ( digit = dynamic_cast<AliEMCALDigit *>(nextdigit()) ) ) {
- if( geom->IsInECA(digit->GetId()) )
- digitsC->Remove(digit) ;
- else
- break ;
- }
- notremovedECA = kFALSE ;
- }
-
- }
- else if( geom->IsInHCA(digit->GetId()) ) {
- where = 2 ; // to tell we are in HCAL
- // start a new HCAL cluster
- if(fNumberOfHCAClusters >= aHCARecPoints->GetSize())
- aHCARecPoints->Expand(2*fNumberOfHCAClusters+1);
- AliEMCALTowerRecPoint * rp = new AliEMCALTowerRecPoint("") ;
- rp->SetHCA() ;
- aHCARecPoints->AddAt(rp, fNumberOfHCAClusters) ;
- clu = dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(fNumberOfHCAClusters)) ;
- fNumberOfHCAClusters++ ;
- clu->AddDigit(*digit, Calibrate(digit->GetAmp(), where));
- clusterHCAdigitslist[iDigitInHCACluster] = digit->GetIndexInList() ;
- iDigitInHCACluster++ ;
- digitsC->Remove(digit) ;
- if (gDebug == 2 )
- Info("MakeClusters","OK id = %d, ene = %f , thre = %f ", digit->GetId(),Calibrate(digit->GetAmp(), 2), fHCAClusteringThreshold) ;
-
- nextdigit.Reset() ;
-
- // Here we remove remaining PRE digits, which cannot make a cluster
+ printf("MakeClusters: OK id = %d, ene = %f , thre = %f ", digit->GetId(),Calibrate(digit->GetAmp()), fECAClusteringThreshold) ;
- if( notremovedPRE ) {
- while( ( digit = dynamic_cast<AliEMCALDigit *>(nextdigit()) ) ) {
- if( geom->IsInPRE(digit->GetId()) )
- digitsC->Remove(digit) ;
- else
- break ;
- }
- notremovedPRE = kFALSE ;
- }
- }
-
+ }
nextdigit.Reset() ;
AliEMCALDigit * digitN ;
Int_t index = 0 ;
- // Do the Clustering in each of the three section ECAL/PRE/HCAL
+ // Do the Clustering
while (index < iDigitInECACluster){ // scan over digits already in cluster
digit = (AliEMCALDigit*)digits->At(clusterECAdigitslist[index]) ;
index++ ;
while ( (digitN = (AliEMCALDigit *)nextdigit()) ) { // scan over the reduced list of digits
- Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
- // Info("MakeClusters","id1 = %d, id2 = %d , neighbours = %d", digit->GetId(), digitN->GetId(), ineb) ;
+ Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
switch (ineb ) {
case 0 : // not a neighbour
break ;
case 1 : // are neighbours
- clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), 1) ) ;
+ clu->AddDigit(*digitN, Calibrate( digitN->GetAmp()) ) ;
clusterECAdigitslist[iDigitInECACluster] = digitN->GetIndexInList() ;
iDigitInECACluster++ ;
digitsC->Remove(digitN) ;
endofloop1: ;
nextdigit.Reset() ;
} // loop over ECA cluster
-
- index = 0 ;
- while (index < iDigitInPRECluster){ // scan over digits already in cluster
- digit = (AliEMCALDigit*)digits->At(clusterPREdigitslist[index]) ;
- index++ ;
- while ( (digitN = (AliEMCALDigit *)nextdigit()) ) { // scan over the reduced list of digits
- Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
- // Info("MakeClusters","id1 = %d, id2 = %d , neighbours = %d", digit->GetId(), digitN->GetId(), ineb) ;
- switch (ineb ) {
- case 0 : // not a neighbour
- break ;
- case 1 : // are neighbours
- clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), 0) ) ;
- clusterPREdigitslist[iDigitInPRECluster] = digitN->GetIndexInList() ;
- iDigitInPRECluster++ ;
- digitsC->Remove(digitN) ;
- break ;
- case 2 : // too far from each other
- goto endofloop2;
- } // switch
-
- } // while digitN
-
- endofloop2: ;
- nextdigit.Reset() ;
- } // loop over PRE cluster
-
- index = 0 ;
- while (index < iDigitInHCACluster){ // scan over digits already in cluster
- digit = (AliEMCALDigit*)digits->At(clusterHCAdigitslist[index]) ;
- index++ ;
- while ( (digitN = (AliEMCALDigit *)nextdigit()) ) { // scan over the reduced list of digits
- Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
- //Info("MakeClusters","id1 = %d, id2 = %d , neighbours = %d", digit->GetId(), digitN->GetId(), ineb) ;
- switch (ineb ) {
- case 0 : // not a neighbour
- break ;
- case 1 : // are neighbours
- clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), 2) ) ;
- clusterHCAdigitslist[iDigitInHCACluster] = digitN->GetIndexInList() ;
- iDigitInHCACluster++ ;
- digitsC->Remove(digitN) ;
- break ;
- case 2 : // too far from each other
- goto endofloop3;
- } // switch
- } // while digitN
-
- endofloop3: ;
- nextdigit.Reset() ;
- } // loop over HCA cluster
-
} // energy theshold
} // while digit
delete digitsC ;
}
//____________________________________________________________________________
-void AliEMCALClusterizerv1::MakeUnfolding()
+void AliEMCALClusterizerv1::MakeUnfolding() const
{
Fatal("AliEMCALClusterizerv1::MakeUnfolding", "--> Unfolding not implemented") ;
void AliEMCALClusterizerv1::UnfoldCluster(AliEMCALTowerRecPoint * /*iniTower*/,
Int_t /*nMax*/,
AliEMCALDigit ** /*maxAt*/,
- Float_t * /*maxAtEnergy*/)
+ Float_t * /*maxAtEnergy*/) const
{
// Performs the unfolding of a cluster with nMax overlapping showers
// Print parameters
- TString taskName(GetName()) ;
+ TString taskName(GetName()) ;
taskName.ReplaceAll(Version(), "") ;
- message += "--------------- " ;
- message += taskName.Data() ;
- message += " " ;
- message += GetTitle() ;
- message += "-----------\n" ;
- message += "Clusterizing digits from the file: " ;
- message += taskName.Data() ;
- message += "\n Branch: " ;
- message += GetName() ;
- message += "\n Pre Shower Clustering threshold = " ;
- message += fPREClusteringThreshold ;
- message += "\n Pre Shower Local Maximum cut = " ;
- message += fPRELocMaxCut ;
- message += "\n Pre Shower Logarothmic weight = " ;
- message += fPREW0 ;
- message += "\n ECA Clustering threshold = " ;
- message += fECAClusteringThreshold ;
- message += "\n ECA Local Maximum cut = " ;
- message += fECALocMaxCut ;
- message += "\n ECA Logarothmic weight = " ;
- message += fECAW0 ;
- message += "\n Pre Shower Clustering threshold = " ;
- message += fHCAClusteringThreshold ;
- message += "\n HCA Local Maximum cut = " ;
- message += fHCALocMaxCut ;
- message += "\n HCA Logarothmic weight = " ;
- message += fHCAW0 ;
+ printf("--------------- ");
+ printf(taskName.Data()) ;
+ printf(" ");
+ printf(GetTitle()) ;
+ printf("-----------\n");
+ printf("Clusterizing digits from the file: ");
+ printf(taskName.Data());
+ printf("\n Branch: ");
+ printf(GetName());
+ printf("\n ECA Local Maximum cut = %f", fECALocMaxCut);
+ printf("\n ECA Logarothmic weight = %f", fECAW0);
if(fToUnfold)
- message +="\nUnfolding on\n" ;
+ printf("\nUnfolding on\n");
else
- message += "\nUnfolding off\n";
+ printf("\nUnfolding off\n");
- message += "------------------------------------------------------------------" ;
+ printf("------------------------------------------------------------------");
}
else
- message += "AliEMCALClusterizerv1 not initialized " ;
-
- Info("Print", message.Data() ) ;
+ printf("AliEMCALClusterizerv1 not initialized ") ;
}
//____________________________________________________________________________
void AliEMCALClusterizerv1::PrintRecPoints(Option_t * option)
{
- // Prints list of RecPoints produced at the current pass of AliEMCALClusterizer
-
- TObjArray * aPRERecPoints = AliEMCALGetter::Instance()->PRERecPoints() ;
+ // Prints list of RecPoints produced at the current pass of AliEMCALClusterizer
TObjArray * aECARecPoints = AliEMCALGetter::Instance()->ECARecPoints() ;
- TObjArray * aHCARecPoints = AliEMCALGetter::Instance()->HCARecPoints() ;
-
- Info("PrintRecPoints", "Clusterization result:") ;
+ printf("PrintRecPoints: Clusterization result:") ;
printf("event # %d\n", gAlice->GetEvNumber() ) ;
- printf(" Found %d PRE SHOWER RecPoints, %d ECA Rec Points and %d HCA Rec Points\n ",
- aPRERecPoints->GetEntriesFast(), aECARecPoints->GetEntriesFast(), aHCARecPoints->GetEntriesFast() ) ;
+ printf(" Found %d ECA Rec Points\n ",
+ aECARecPoints->GetEntriesFast()) ;
- fRecPointsInRun += aPRERecPoints->GetEntriesFast() ;
fRecPointsInRun += aECARecPoints->GetEntriesFast() ;
- fRecPointsInRun += aHCARecPoints->GetEntriesFast() ;
if(strstr(option,"all")) {
-
- //Pre shower recPoints
-
- printf("-----------------------------------------------------------------------\n") ;
- printf("Clusters in PRE section\n") ;
- printf("Index Ene(GeV) Multi Module phi r theta X Y Z Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
-
- Int_t index ;
-
- for (index = 0 ; index < aPRERecPoints->GetEntries() ; index++) {
- AliEMCALTowerRecPoint * rp = dynamic_cast<AliEMCALTowerRecPoint *>(aPRERecPoints->At(index)) ;
- TVector3 globalpos;
- rp->GetGlobalPosition(globalpos);
- TVector3 localpos;
- rp->GetLocalPosition(localpos);
- Float_t lambda[2];
- rp->GetElipsAxis(lambda);
- Int_t * primaries;
- Int_t nprimaries;
- primaries = rp->GetPrimaries(nprimaries);
- printf("\n%6d %8.4f %3d %2d %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4f %4f %2d : ",
- rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(), rp->GetEMCALArm(),
- globalpos.X(), globalpos.Y(), globalpos.Z(), localpos.X(), localpos.Y(), localpos.Z(),
- rp->GetDispersion(), lambda[0], lambda[1], nprimaries) ;
- for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
- printf("%d ", primaries[iprimary] ) ;
- }
- }
-
+ Int_t index =0;
printf("\n-----------------------------------------------------------------------\n") ;
printf("Clusters in ECAL section\n") ;
printf("Index Ene(GeV) Multi Module phi r theta X Y Z Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
printf("%d ", primaries[iprimary] ) ;
}
}
-
- printf("\n-----------------------------------------------------------------------\n") ;
- printf("Clusters in HCAL section\n") ;
- printf("Index Ene(GeV) Multi Module phi r theta X Y Z Dispersion Lambda 1 Lambda 2 # of prim Primaries list\n") ;
-
- for (index = 0 ; index < aHCARecPoints->GetEntries() ; index++) {
- AliEMCALTowerRecPoint * rp = dynamic_cast<AliEMCALTowerRecPoint * >(aHCARecPoints->At(index)) ;
- TVector3 globalpos;
- rp->GetGlobalPosition(globalpos);
- TVector3 localpos;
- rp->GetLocalPosition(localpos);
- Float_t lambda[2];
- rp->GetElipsAxis(lambda);
- Int_t * primaries;
- Int_t nprimaries;
- primaries = rp->GetPrimaries(nprimaries);
- printf("\n%6d %8.4f %3d %2d %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4.1f %4f %4f %2d : ",
- rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(), rp->GetEMCALArm(),
- globalpos.X(), globalpos.Y(), globalpos.Z(), localpos.X(), localpos.Y(), localpos.Z(),
- rp->GetDispersion(), lambda[0], lambda[1], nprimaries) ;
- for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
- printf("%d ", primaries[iprimary] ) ;
- }
- }
-
printf("\n-----------------------------------------------------------------------\n");
}
}
virtual Int_t AreNeighbours(AliEMCALDigit * d1, AliEMCALDigit * d2)const ;
// Checks if digits are in neighbour cells
- virtual Float_t Calibrate(Int_t amp, Int_t where)const ; // Tranforms Amp to energy
+ virtual Float_t Calibrate(Int_t amp)const ; // Tranforms Amp to energy
- virtual void GetNumberOfClustersFound(int * numb )const{ numb[0] = fNumberOfPREClusters ;
- numb[1] = fNumberOfECAClusters ;
- numb[2] = fNumberOfHCAClusters ; }
-
- virtual Float_t GetPREClusteringThreshold()const{ return fPREClusteringThreshold; }
- virtual Float_t GetECAClusteringThreshold()const{ return fECAClusteringThreshold;}
- virtual Float_t GetHCAClusteringThreshold()const{ return fHCAClusteringThreshold;}
-
- virtual Float_t GetPRELocalMaxCut()const { return fPRELocMaxCut;}
- virtual Float_t GetPREShoLogWeight()const { return fPREW0;}
+ virtual void GetNumberOfClustersFound(int numb )const{ numb = fNumberOfECAClusters ;}
+ virtual Float_t GetECAClusteringThreshold()const{ return fECAClusteringThreshold;}
virtual Float_t GetECALocalMaxCut()const { return fECALocMaxCut;}
- virtual Float_t GetECALogWeight()const { return fECAW0;}
- virtual Float_t GetHCALocalMaxCut()const { return fHCALocMaxCut;}
- virtual Float_t GetHCALogWeight()const { return fHCAW0;}
+ virtual Float_t GetECALogWeight()const { return fECAW0;}
virtual Float_t GetTimeGate() const { return fTimeGate ; }
virtual const char * GetRecPointsBranch() const{ return GetName() ;}
virtual void SetECAClusteringThreshold(Float_t cluth) { fECAClusteringThreshold = cluth ; }
virtual void SetECALocalMaxCut(Float_t cut) { fECALocMaxCut = cut ; }
virtual void SetECALogWeight(Float_t w) { fECAW0 = w ; }
- virtual void SetHCAClusteringThreshold(Float_t cluth) { fHCAClusteringThreshold = cluth ; }
- virtual void SetHCALocalMaxCut(Float_t cut) { fHCALocMaxCut = cut ; }
- virtual void SetHCALogWeight(Float_t w) { fHCAW0 = w ; }
virtual void SetTimeGate(Float_t gate) { fTimeGate = gate ;}
- virtual void SetPREClusteringThreshold(Float_t cluth) { fPREClusteringThreshold = cluth ; }
- virtual void SetPRELocalMaxCut(Float_t cut) { fPRELocMaxCut = cut ; }
- virtual void SetPRELogWeight(Float_t w) { fPREW0 = w ; }
virtual void SetUnfolding(Bool_t toUnfold = kTRUE ) {fToUnfold = toUnfold ;}
static Double_t ShowerShape(Double_t r) ; // Shape of EM shower used in unfolding;
//class member function (not object member function)
void Init() ;
void InitParameters() ;
- virtual void MakeUnfolding() ;
+ virtual void MakeUnfolding() const;
void UnfoldCluster(AliEMCALTowerRecPoint * /*iniEmc*/, Int_t /*Nmax*/,
AliEMCALDigit ** /*maxAt*/,
- Float_t * /*maxAtEnergy*/ ) ; //Unfolds cluster using TMinuit package
- void PrintRecPoints(Option_t * /*option*/) ;
+ Float_t * /*maxAtEnergy*/ ) const; //Unfolds cluster using TMinuit package
+ void PrintRecPoints(Option_t * option) ;
private:
Int_t fNTowers ; // number of Towers in EMCAL
Bool_t fToUnfold ; // To perform unfolding
-
- Int_t fNumberOfPREClusters ; // number of clusters found in PRE section
Int_t fNumberOfECAClusters ; // number of clusters found in EC section
- Int_t fNumberOfHCAClusters ; // number of clusters found in HC section
//Calibration parameters... to be replaced by database
- Float_t fADCchannelPRE ; // width of one ADC channel for PRE section (GeV)
- Float_t fADCpedestalPRE ; // pedestal of ADC for PRE section (GeV)
Float_t fADCchannelECA ; // width of one ADC channel for EC section (GeV)
Float_t fADCpedestalECA ; // pedestal of ADC for EC section (GeV)
- Float_t fADCchannelHCA ; // width of one ADC channel for HC section (GeV)
- Float_t fADCpedestalHCA ; // pedestal of ADC for HC section (GeV)
Float_t fECAClusteringThreshold ; // minimum energy to include a EC digit in a cluster
- Float_t fHCAClusteringThreshold ; // minimum energy to include a HC digit in a cluster
- Float_t fPREClusteringThreshold ; // minimum energy to include a PRE digit in a cluster
Float_t fECALocMaxCut ; // minimum energy difference to distinguish local maxima in a cluster
Float_t fECAW0 ; // logarithmic weight for the cluster center of gravity calculation
- Float_t fHCALocMaxCut ; // minimum energy difference to distinguish local maxima in a cluster
- Float_t fHCAW0 ; // logarithmic weight for the cluster center of gravity calculation
- Float_t fPRELocMaxCut ; // minimum energy difference to distinguish local maxima in a CPV cluster
- Float_t fPREW0 ; // logarithmic weight for the CPV cluster center of gravity calculation
Int_t fRecPointsInRun ; //! Total number of recpoints in one run
Float_t fTimeGate ; // Maximum time difference between the digits in ont EMC cluster
{
// retrieves the primary particle number given its index in the list
Int_t rv = -1 ;
- if ( index <= fNprimary && index > 0){
+ if ( (index <= fNprimary) && (index > 0)){
rv = fPrimary[index-1] ;
}
//____________________________________________________________________________
void AliEMCALDigit::ShiftPrimary(Int_t shift){
- //shifts primary nimber to BIG offset, to separate primary in different TreeK
+ //shifts primary number to BIG offset, to separate primary in different TreeK
Int_t index ;
for(index = 0; index <fNprimary; index++ ){
fPrimary[index] = fPrimary[index]+ shift * 10000000 ;}
// EMCAL digit:
//
// A Digit is the sum of energy in a Tower (Hit sum) and stores information, about primaries
-// and enterring particle contributing to a Digit
+// and entering particle contributing to a Digit
//
//*-- Author: Sahal Yacoob (LBL)
// based on : AliPHOSDigit
void SetAmp(Int_t amp) { fAmp= amp ; }
void SetId(Int_t id) {fId = id ;}
void SetTime(Float_t time) {fTime = time ;}
- void ShiftPrimary(Int_t shift); // shift to semarate different TreeK in merging
+ void ShiftPrimary(Int_t shift); // shift to separate different TreeK in merging
private:
fADCchannelEC = d.fADCchannelEC ;
fADCpedestalEC = d.fADCpedestalEC ;
fNADCEC = d.fNADCEC ;
- fADCchannelHC = d.fADCchannelHC ;
- fADCpedestalHC = d.fADCpedestalHC ;
- fNADCHC = d.fNADCHC ;
- fADCchannelPRE = d.fADCchannelPRE ;
- fADCpedestalPRE = d.fADCpedestalPRE ;
- fNADCPRE = d.fNADCPRE ;
fEventFolderName = d.fEventFolderName;
}
digits->Clear() ;
const AliEMCALGeometry *geom = gime->EMCALGeometry() ;
- //Making digits with noise, first EMC
- Int_t nEMC = 0 ;
- if (geom->GetNHCLayers() > 0 )
- nEMC = 3*geom->GetNPhi()*geom->GetNZ(); //max number of digits possible (Preshower, ECAL, HCAL)
- else
- nEMC = 2*geom->GetNPhi()*geom->GetNZ(); //max number of digits possible (Preshower, ECAL)
+
+ //Making digits from noise first
+ Int_t nEMC = 0 ;
+ nEMC = geom->GetNPhi()*geom->GetNZ(); //max number of digits possible
Int_t absID ;
Float_t b = TMath::Abs( a /fTimeSignalLength) ;
//Mark the beginning of the signal
new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime(),0, b);
- //Mark the end of the ignal
+ //Mark the end of the signal
new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime()+fTimeSignalLength, -a, -b);
// loop over input
else
curSDigit = 0 ;
//May be several digits will contribute from the same input
- while(curSDigit && curSDigit->GetId() == absID){
+ while(curSDigit && (curSDigit->GetId() == absID)){
//Shift primary to separate primaries belonging different inputs
Int_t primaryoffset ;
if(fManager)
}
}
// add the noise now
-
- if (geom->IsInECA(digit->GetId()))
- amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
- else if (geom->IsInPRE(digit->GetId()))
- amp += TMath::Abs(gRandom->Gaus(0., fPinNoise/100.)) ; // arbitrarely divide by 100, assuming that the gain of APD will be higher
- else if (geom->IsInHCA(digit->GetId()))
- amp += TMath::Abs(gRandom->Gaus(0., fPinNoise/10.)) ; // arbitrarely divide by 10, assuming that the gain of APD will be higher
+ amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
digit->SetAmp(sDigitizer->Digitize(amp)) ;
}
//remove digits below thresholds
for(i = 0 ; i < nEMC ; i++){
digit = dynamic_cast<AliEMCALDigit*>( digits->At(i) ) ;
- Float_t threshold = 0 ;
-
- if (geom->IsInECA(digit->GetId()))
- threshold = fDigitThreshold ;
- else if (geom->IsInPRE(digit->GetId()))
- threshold = fDigitThreshold / 100. ; // arbitrary see before when noise is added
- else if (geom->IsInHCA(digit->GetId()))
- threshold = fDigitThreshold / 10. ; // arbitrary see before when noise is added
-
+ Float_t threshold = fDigitThreshold ;
if(sDigitizer->Calibrate( digit->GetAmp() ) < threshold)
digits->RemoveAt(i) ;
else
digit = dynamic_cast<AliEMCALDigit *>( digits->At(i) ) ;
digit->SetIndexInList(i) ;
Float_t energy = sDigitizer->Calibrate(digit->GetAmp()) ;
- digit->SetAmp(DigitizeEnergy(energy,digit->GetId()) ) ;
+ digit->SetAmp(DigitizeEnergy(energy) ) ;
}
}
//____________________________________________________________________________
-Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t absId)
+Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy)
{
- Int_t channel = -999;
- AliEMCALGeometry * geom = AliEMCALGetter::Instance()->EMCALGeometry() ;
-
- if(geom->IsInPRE(absId)){ //digitize as PRE section
- channel = static_cast<Int_t>(TMath::Ceil( (energy + fADCpedestalPRE)/fADCchannelPRE )) ;
- if(channel > fNADCPRE )
- channel = fNADCPRE ;
- }
- else if(geom->IsInECA(absId)){ //digitize as ECAL section
- channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC )) ;
- if(channel > fNADCEC )
- channel = fNADCEC ;
- }
- else if(geom->IsInHCA(absId)){ //digitize as HCAL section
- channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalHC)/fADCchannelHC )) ;
- if(channel > fNADCHC )
- channel = fNADCHC ;
- }
-
+ Int_t channel = -999;
+ channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC )) ;
+ if(channel > fNADCEC )
+ channel = fNADCEC ;
return channel ;
}
if(strstr(option,"tim")){
gBenchmark->Stop("EMCALDigitizer");
- Info("Exec", "took %f seconds for Digitizing %f seconds per event",
+ printf("Exec: took %f seconds for Digitizing %f seconds per event",
gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nevents ) ;
}
}
fADCpedestalEC = 0.005 ; // GeV
fNADCEC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC
- fADCchannelHC = 0.000220; // width of one ADC channel in GeV
- fADCpedestalHC = 0.005 ; // GeV
- fNADCHC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC
-
- fADCchannelPRE = 0.0000300; // width of one ADC channel in Pre Shower
- fADCpedestalPRE = 0.005 ; // GeV
- fNADCPRE = (Int_t) TMath::Power(2,12); // number of channels in Pre ShowerADC
-
fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
}
void AliEMCALDigitizer::Print()const
{
// Print Digitizer's parameters
- Info("Print", "\n------------------- %s -------------", GetName() ) ;
+ printf("Print: \n------------------- %s -------------", GetName() ) ;
if( strcmp(fEventFolderName.Data(), "") != 0 ){
printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
printf("---------------------------------------------------\n") ;
}
else
- Info("Print", "AliEMCALDigitizer not initialized") ;
+ printf("Print: AliEMCALDigitizer not initialized") ;
}
//__________________________________________________________________
AliEMCALGetter * gime = AliEMCALGetter::Instance(GetTitle(), fEventFolderName) ;
TClonesArray * digits = gime->Digits() ;
- Info("PrintDigits", "%d", digits->GetEntriesFast()) ;
+ printf("PrintDigits: %d", digits->GetEntriesFast()) ;
printf("\nevent %d", gAlice->GetEvNumber()) ;
printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
//__________________________________________________________________
void AliEMCALDigitizer::Unload()
{
-
+ // Unloads the SDigits and Digits
Int_t i ;
for(i = 1 ; i < fInput ; i++){
TString tempo(fEventNames[i]) ;
const Float_t GetTimeResolution() const { return fTimeResolution ; }
const Float_t GetECAchannel() const { return fADCchannelEC ; }
const Float_t GetECApedestal() const { return fADCpedestalEC ; }
- const Float_t GetHCAchannel() const { return fADCchannelHC ; }
- const Float_t GetHCApedestal() const { return fADCpedestalHC ; }
- const Float_t GetPREchannel() const { return fADCchannelPRE ; }
- const Float_t GetPREpedestal() const { return fADCpedestalPRE ; }
void SetDigitThreshold(Float_t EMCThreshold) {fDigitThreshold = EMCThreshold;}
void SetPinNoise(Float_t PinNoise ) {fPinNoise = PinNoise;}
//Calculate the time of crossing of the threshold by front edge
Float_t FrontEdgeTime(TClonesArray * ticks) ;
- Int_t DigitizeEnergy(Float_t energy, Int_t absId) ;
+ Int_t DigitizeEnergy(Float_t energy) ;
private:
Float_t fADCchannelEC ; // width of one ADC channel in EC section (GeV)
Float_t fADCpedestalEC ; //
Int_t fNADCEC ; // number of channels in EC section ADC
- Float_t fADCchannelHC ; // width of one ADC channel in HC section (GeV)
- Float_t fADCpedestalHC ; //
- Int_t fNADCHC ; // number of channels in HC section ADC
- Float_t fADCchannelPRE ; // width of one ADC channel in PRE section (GeV)
- Float_t fADCpedestalPRE ; //
- Int_t fNADCPRE ; // number of channels in PRE section ADC
TString fEventFolderName; // skowron: name of EFN to read data from in stand alone mode
/* $Id$ */
-
+//____________________________________________________________________________
+//*--
//*-- Author: Andreas Morsch (CERN)
+//*--
+//*--
+////////////////////////////////////////////////////////////////////////////
#include "TMath.h"
#include <TRandom.h>
//_________________________________________________________________________
// A Particle modified by EMCAL response and produced by AliEMCALvFast
+//*--
// To become a general class of AliRoot ?
-//
+//*--
//*-- Author: Yves Schutz (SUBATECH)
+//*--
+/////////////////////////////////////////////////////////////////////////////
// --- ROOT system ---
{
// Print the type, energy and momentum of the reconstructed particle
- Info("Print", "Summary:") ;
+ printf("Print Summary:") ;
printf("AliEMCALFastRecParticle > type is %s\n", Name().Data()) ;
printf(" Energy = %f\n", fE) ;
printf(" Px = %f\n", fPx) ;
// and : Yves Schutz (SUBATECH)
// and : Jennifer Klay (LBL)
-// --- ROOT system ---
-
-// --- Standard library ---
-//#include <stdlib.h>
-
// --- AliRoot header files ---
-//#include <TError.h>
#include <TMath.h>
#include <TVector3.h>
//______________________________________________________________________
const Bool_t AliEMCALGeometry::AreInSameTower(Int_t id1, Int_t id2) const {
+ // Find out whether two hits are in the same tower
Int_t idmax = TMath::Max(id1, id2) ;
Int_t idmin = TMath::Min(id1, id2) ;
if ( ((idmax - GetNZ() * GetNPhi()) == idmin ) ||
//______________________________________________________________________
void AliEMCALGeometry::Init(void){
// Initializes the EMCAL parameters
- // naming convention : GUV_L_WX_N_YZ_M gives the composition of a tower
- // UV inform about the compsition of the pre-shower section:
- // thickness in mm of Pb radiator (U) and of scintillator (V), and number of scintillator layers (L)
+ // naming convention : GUV_WX_N_ gives the composition of a tower
// WX inform about the composition of the EM calorimeter section:
- // thickness in mm of Pb radiator (W) and of scintillator (X), and number of scintillator layers (N)
- // YZ inform about the composition of the hadron calorimeter section:
- // thickness in mm of Cu radiator (Y) and of scintillator (Z), and number of scintillator layers (M)
- // Valid geometries are G56_2_55_19_104_14
- // G56_2_55_19 or EMCAL_5655_21
- // G65_2_64_19 or EMCAL_6564_21
-
- fgInit = kFALSE; // Assume failer untill proven otherwise.
- TString name(GetName()) ;
-
- if ( name == "G56_2_55_19_104_14" ) {
- fPRPbRadThickness = 0.5; // cm, Thickness of the Pb radiators for the preshower section
- fPRScintThick = 0.6; // cm, Thickness of the sintilator for the preshower section of the tower
- fNPRLayers = 2; // number of scintillator layers in the preshower section
-
- fECPbRadThickness = 0.5; // cm, Thickness of the Pb radiators for the EM calorimeter section
- fECScintThick = 0.5; // cm, Thickness of the sintilator for the EM alorimeter section of the tower
- fNECLayers = 19; // number of scintillator layers in the EM calorimeter section
-
- fHCCuRadThickness = 1.0; // cm, Thickness of the Cu radiators.
- fHCScintThick = 0.4; // cm, Thickness of the sintilator for the hadronic alorimeter section of the tower
- fNHCLayers = 14; // number of scintillator layers in the hadronic calorimeter section
-
- fSampling = 11.3 ;
- fSummationFraction = 0.8 ;
+ // thickness in mm of Pb radiator (W) and of scintillator (X), and number of scintillator layers (N)
+ // New geometry: EMCAL_55_25
- fAlFrontThick = 3.0; // cm, Thickness of front Al layer
- fGap2Active = 1.0; // cm, Gap between Al and 1st Scintillator
- }
- else if ( name == "G56_2_55_19" || name == "EMCAL_5655_21" ) {
- fPRPbRadThickness = 0.5; // cm, Thickness of the Pb radiators for the preshower section
- fPRScintThick = 0.6; // cm, Thickness of the sintilator for the preshower section of the tower
- fNPRLayers = 2; // number of scintillator layers in the preshower section
-
- fECPbRadThickness = 0.5; // cm, Thickness of the Pb radiators for the EM calorimeter section
- fECScintThick = 0.5; // cm, Thickness of the sintilator for the EM alorimeter section of the tower
- fNECLayers = 19; // number of scintillator layers in the EM calorimeter section
-
- fHCCuRadThickness = 0.0; // cm, Thickness of the Cu radiators.
- fHCScintThick = 0.0; // cm, Thickness of the sintilator for the hadronic alorimeter section of the tower
- fNHCLayers = 0; // number of scintillator layers in the hadronic calorimeter section
+ fgInit = kFALSE; // Assume failed until proven otherwise.
+ TString name(GetName()) ;
+ if (name == "EMCAL_55_25") {
+ fECPbRadThickness = 0.5; // cm, Thickness of the Pb radiators
+ fECScintThick = 0.5; // cm, Thickness of the scintillator
+ fNECLayers = 25; // number of scintillator layers
- fSampling = 11.3 ;
- fSummationFraction = 0.8 ;
+ fSampling = 11.8;
- fAlFrontThick = 3.0; // cm, Thickness of front Al layer
+ fAlFrontThick = 3.5; // cm, Thickness of front Al layer
fGap2Active = 1.0; // cm, Gap between Al and 1st Scintillator
}
- else if ( name == "G65_2_64_19" || name == "EMCAL_6564_21" ) {
- fPRPbRadThickness = 0.6; // cm, Thickness of the Pb radiators for the preshower section
- fPRScintThick = 0.5; // cm, Thickness of the sintilator for the preshower section of the tower
- fNPRLayers = 2; // number of scintillator layers in the preshower section
-
- fECPbRadThickness = 0.6; // cm, Thickness of the Pb radiators for the EM calorimeter section
- fECScintThick = 0.4; // cm, Thickness of the sintilator for the EM alorimeter section of the tower
- fNECLayers = 19; // number of scintillator layers in the EM calorimeter section
-
- fHCCuRadThickness = 0.0; // cm, Thickness of the Cu radiators.
- fHCScintThick = 0.0; // cm, Thickness of the sintilator for the hadronic alorimeter section of the tower
- fNHCLayers = 0; // number of scintillator layers in the hadronic calorimeter section
-
- fSampling = 16. ;
- fSummationFraction = 0.8 ;
-
- fAlFrontThick = 3.0; // cm, Thickness of front Al layer
- fGap2Active = 1.0; // cm, Gap between Al and 1st Scintillator
+ else if( name == "G56_2_55_19" || name == "EMCAL_5655_21" || name == "G56_2_55_19_104_14"|| name == "G65_2_64_19" || name == "EMCAL_6564_21"){
+ Fatal("Init", "%s is an old geometry! Please update your Config file", name.Data()) ;
}
else
Fatal("Init", "%s is an undefined geometry!", name.Data()) ;
- // if( name != "EMCALArch1a" &&
-// name != "EMCALArch1b" &&
-// name != "EMCALArch2a" &&
-// name != "EMCALArch2b" &&
-// name != "EMCALArch1aN" ){
-// Fatal("Init", "%s is not a known geometry (choose among EMCALArch1a, EMCALArch1b, EMCALArch2a and EMCALArch2b, EMCALArch1aN)", name.Data()) ;
-// } // end if
-// //
-// if ( name == "EMCALArch1a" ||
-// name == "EMCALArch1b" ||
-// name == "EMCALArch1aN") {
-// fNZ = 96;
-// fNPhi = 144;
-// } // end if
-// if ( name == "EMCALArch2a" ||
-// name == "EMCALArch2b" ) {
-// fNZ = 112;
-// fNPhi = 168;
-// } // end if
-// if ( name == "EMCALArch1a" ||
-// name == "EMCALArch2a" ) {
-// fNPRLayers = 2;
-// fNECLayers = 19;
-// fNHCLayers = 0;
-// } // end if
-// if ( name == "EMCALArch1b" ||
-// name == "EMCALArch2b" ) {
-// fNPRLayers = 2;
-// fNECLayers = 23;
-// fNHCLayers = 0;
-// } // end if
-// if ( name == "EMCALArch1aN") {
-// fNPRLayers = 2;
-// fNECLayers = 19;
-// fNHCLayers = 14;
-// }
-
// geometry
- fNZ = 96; // granularity along Z (eta)
- fNPhi = 144; // granularity in phi (azimuth)
- fArm1PhiMin = 60.0; // degrees, Starting EMCAL Phi position
- fArm1PhiMax = 180.0; // degrees, Ending EMCAL Phi position
- fArm1EtaMin = -0.7; // pseudorapidity, Starting EMCAL Eta position
- fArm1EtaMax = +0.7; // pseudorapidity, Ending EMCAL Eta position
+ fNZ = 114; // granularity along Z (eta)
+ fNPhi = 168; // granularity in phi (azimuth)
+ fArm1PhiMin = 60.0; // degrees, Starting EMCAL Phi position
+ fArm1PhiMax = 180.0; // degrees, Ending EMCAL Phi position
+ fArm1EtaMin = -0.7; // pseudorapidity, Starting EMCAL Eta position
+ fArm1EtaMax = +0.7; // pseudorapidity, Ending EMCAL Eta position
fIPDistance = 454.0; // cm, Radial distance to inner surface of EMCAL
- fShellThickness = fAlFrontThick + fGap2Active + 2.*(GetPRScintThick() + GetPRPbRadThick()) + // pre shower
- (fNECLayers-1)*(GetECScintThick()+ GetECPbRadThick()) + // E cal -1 because the last element is a scintillator
- fNHCLayers*(GetHCScintThick()+ GetHCCuRadThick()) + // H cal
- GetHCScintThick() ; // last scintillator
+
+ //There is always one more scintillator than radiator layer because of the first block of aluminium
+ fShellThickness = fAlFrontThick + fGap2Active + fNECLayers*GetECScintThick()+(fNECLayers-1)*GetECPbRadThick();
+
fZLength = 2.*ZFromEtaR(fIPDistance+fShellThickness,fArm1EtaMax); // Z coverage
fEnvelop[0] = fIPDistance; // mother volume inner radius
fEnvelop[1] = fIPDistance + fShellThickness; // mother volume outer r.
fgInit = kTRUE;
if (gDebug) {
- Info("Init", "geometry of EMCAL named %s is as follows:", name.Data());
- printf( "Tower geometry pre-shower: %d x (%f mm Pb, %f mm Sc) \n", GetNPRLayers(), GetPRPbRadThick(), GetPRScintThick() ) ;
+ printf("Init: geometry of EMCAL named %s is as follows:", name.Data());
printf( " ECAL : %d x (%f mm Pb, %f mm Sc) \n", GetNECLayers(), GetECPbRadThick(), GetECScintThick() ) ;
- if ( GetNHCLayers() > 0 )
- printf( " HCAL : %d x (%f mm Pb, %f mm Sc) \n", GetNHCLayers(), GetHCCuRadThick(), GetHCScintThick() ) ;
printf("Granularity: %d in eta and %d in phi\n", GetNZ(), GetNPhi()) ;
printf("Layout: phi = (%f, %f), eta = (%f, %f), y = %f\n",
GetArm1PhiMin(), GetArm1PhiMax(),GetArm1EtaMin(), GetArm1EtaMax(), GetIPDistance() ) ;
} // end if strcmp(name,"")
}else{
if ( strcmp(fgGeom->GetName(), name) != 0 ) {
- TString message("\n") ;
- message += "current geometry is " ;
- message += fgGeom->GetName() ;
- message += "\n you cannot call " ;
- message += name ;
- ::Info("GetGeometry", message.Data() ) ;
+ printf("\ncurrent geometry is ") ;
+ printf(fgGeom->GetName());
+ printf("\n you cannot call ");
+ printf(name);
}else{
rv = (AliEMCALGeometry *) fgGeom;
} // end if
//______________________________________________________________________
Int_t AliEMCALGeometry::TowerIndex(Int_t ieta,Int_t iphi) const {
// Returns the tower index number from the based on the Z and Phi
- // index numbers. There are 2 times the number of towers to separate
- // out the full towers from the pre-showers.
+ // index numbers.
// Inputs:
- // Int_t ieta // index allong z axis [1-fNZ]
- // Int_t iphi // index allong phi axis [1-fNPhi]
- // Int_t where // 1 = PRE section, 0 = EC section, 2 = HC section
+ // Int_t ieta // index along z axis [1-fNZ]
+ // Int_t iphi // index along phi axis [1-fNPhi]
// Outputs:
// none.
// Returned
}
//______________________________________________________________________
-void AliEMCALGeometry::TowerIndexes(Int_t index,Int_t &ieta,Int_t &iphi,
- Int_t &ipre) const {
+void AliEMCALGeometry::TowerIndexes(Int_t index,Int_t &ieta,Int_t &iphi) const {
// Inputs:
- // Int_t index // Tower index number [1-i*fNZ*fNPhi] PRE(i=1)/ECAL(i=2)/HCAL(i=3)
+ // Int_t index // Tower index number [1-fNZ*fNPhi]
// Outputs:
// Int_t ieta // index allong z axis [1-fNZ]
// Int_t iphi // index allong phi axis [1-fNPhi]
- // Int_t ipre // 0 = ECAL section, 1 = Pre-shower section, 2 = HCAL section
// Returned
// none.
-
- Int_t nindex = 0, itowers = GetNEta() * GetNPhi();
+ Int_t nindex = 0;
- if ( IsInPRE(index) ) { // PRE index
- nindex = index - itowers;
- ipre = 1 ;
- }
- else if ( IsInECA(index) ) { // ECAL index
+ if ( IsInECA(index) ) { // ECAL index
nindex = index ;
- ipre = 0 ;
- }
- else if ( IsInHCA(index) ) { // HCAL index
- nindex = index - 2*itowers;
- ipre = 2 ;
}
else
Fatal("TowerIndexes", "Unexpected Id number!") ;
ieta = nindex - (iphi - 1) * GetNZ() ;
if (gDebug==2)
- Info("TowerIndexes", "index=%d,%d, ieta=%d, iphi = %d", index, nindex,ieta, iphi) ;
+ printf("TowerIndexes: index=%d,%d, ieta=%d, iphi = %d", index, nindex,ieta, iphi) ;
return;
}
// given the tower index number it returns the based on the eta and phi
// of the tower.
// Inputs:
- // Int_t index // Tower index number [1-i*fNZ*fNPhi] PRE(i=1)/ECAL(i=2)/HCAL(i=3)
+ // Int_t index // Tower index number [1-fNZ*fNPhi]
// Outputs:
// Float_t eta // eta of center of tower in pseudorapidity
// Float_t phi // phi of center of tower in degrees
// Returned
// none.
- Int_t ieta, iphi, ipre ;
+ Int_t ieta, iphi;
Float_t deta, dphi ;
- TowerIndexes(index,ieta,iphi,ipre);
+ TowerIndexes(index,ieta,iphi);
if (gDebug == 2)
- Info("EtaPhiFromIndex","index = %d, ieta = %d, iphi = %d", index, ieta, iphi) ;
+ printf("EtaPhiFromIndex: index = %d, ieta = %d, iphi = %d", index, ieta, iphi) ;
deta = (GetArm1EtaMax()-GetArm1EtaMin())/(static_cast<Float_t>(GetNEta()));
eta = GetArm1EtaMin() + ((static_cast<Float_t>(ieta) - 0.5 ))*deta;
Bool_t AliEMCALGeometry::AbsToRelNumbering(Int_t AbsId, Int_t *relid) const {
// Converts the absolute numbering into the following array/
// relid[0] = EMCAL Arm number 1:1
- // relid[1] = 0 ECAL section ; = 1 PRE section; = 2 HCA section
- // relid[2] = Row number inside EMCAL
- // relid[3] = Column number inside EMCAL
+ // relid[1] = Row number inside EMCAL
+ // relid[2] = Column number inside EMCAL
// Input:
// Int_t AbsId // Tower index number [1-2*fNZ*fNPhi]
// Outputs:
- // Int_t *relid // array of 5. Discribed above.
+ // Int_t *relid // array of 3. Discribed above.
Bool_t rv = kTRUE ;
- Int_t ieta=0,iphi=0,ipre=0,index=AbsId;
+ Int_t ieta=0,iphi=0,index=AbsId;
- TowerIndexes(index,ieta,iphi,ipre);
+ TowerIndexes(index,ieta,iphi);
relid[0] = 1;
- relid[1] = ipre;
- relid[2] = ieta;
- relid[3] = iphi;
+ relid[1] = ieta;
+ relid[2] = iphi;
return rv;
}
{
// Converts the relative numbering into the local EMCAL-module (x, z)
// coordinates
- Int_t sect = relid[1]; // PRE/ECAL/HCAL section 1/0/2
- Int_t ieta = relid[2]; // offset along x axis
- Int_t iphi = relid[3]; // offset along z axis
+ Int_t ieta = relid[1]; // offset along x axis
+ Int_t iphi = relid[2]; // offset along z axis
Int_t index;
Float_t eta;
index = TowerIndex(ieta,iphi);
EtaPhiFromIndex(index,eta,phi);
- theta = 180.*(2.0*TMath::ATan(TMath::Exp(-eta)))/TMath::Pi();
+ //theta = 180.*(2.0*TMath::ATan(TMath::Exp(-eta)))/TMath::Pi();
+ theta = 2.0*TMath::ATan(TMath::Exp(-eta));
- // correct for distance to IP different in PRE/ECAL/HCAL
- Float_t d = 0. ;
- if (sect == 1)
- d = GetIP2PRESection() - GetIPDistance() ;
- else if (sect == 0)
- d = GetIP2ECASection() - GetIPDistance() ;
- else if (sect == 2)
- d = GetIP2HCASection() - GetIPDistance() ;
- else
- Fatal("PosInAlice", "Unexpected tower section!") ;
+ // correct for distance to IP
+ Float_t d = GetIP2ECASection() - GetIPDistance() ;
Float_t correction = 1 + d/GetIPDistance() ;
Float_t tantheta = TMath::Tan(theta) * correction ;
// Converts the relative numbering into the local EMCAL-module (x, z)
// coordinates
- Int_t relid[4] ;
+ Int_t relid[3] ;
AbsToRelNumbering(absid, relid) ;
- Int_t ieta = relid[2]; // offset along x axis
- Int_t iphi = relid[3]; // offset along z axis
+ Int_t ieta = relid[1]; // offset along x axis
+ Int_t iphi = relid[2]; // offset along z axis
Int_t index;
Float_t eta;
EtaPhiFromIndex(index,eta,phi);
theta = 2.0*TMath::ATan(TMath::Exp(-eta)) ;
- // correct for distance to IP different in PRE/ECAL/HCAL
+ // correct for distance to IP
Float_t d = 0. ;
- if (IsInPRE(absid))
- d = GetIP2PRESection() - GetIPDistance() ;
- else if (IsInECA(absid))
+ if (IsInECA(absid))
d = GetIP2ECASection() - GetIPDistance() ;
- else if (IsInHCA(absid))
- d = GetIP2HCASection() - GetIPDistance() ;
else
Fatal("PosInAlice", "Unexpected id # %d!", absid) ;
// Returned
// none.
- Float_t eta,theta, phi,cyl_radius=0. ;
+ Float_t eta,theta, phi,cylradius=0. ;
- Int_t ieta = relid[2]; // offset along x axis
- Int_t iphi = relid[3]; // offset along z axis
- Int_t ipre = relid[1]; // indicates 0 ECAL section, 1 PRE section, 2 HCAL section.
+ Int_t ieta = relid[1]; // offset along x axis
+ Int_t iphi = relid[2]; // offset along z axis.
Int_t index;
index = TowerIndex(ieta,iphi);
EtaPhiFromIndex(index,eta,phi);
theta = 180.*(2.0*TMath::ATan(TMath::Exp(-eta)))/TMath::Pi();
- if ( ipre == 0 )
- cyl_radius = GetIP2ECASection() ;
- else if ( ipre == 1 )
- cyl_radius = GetIP2PRESection() ;
- else if ( ipre == 2 )
- cyl_radius = GetIP2HCASection() ;
- else
- Fatal("XYZFromIndex", "Unexpected Tower section # %d", ipre) ;
+ cylradius = GetIP2ECASection() ;
Double_t kDeg2Rad = TMath::DegToRad() ;
- x = cyl_radius * TMath::Cos(phi * kDeg2Rad ) ;
- y = cyl_radius * TMath::Sin(phi * kDeg2Rad ) ;
- z = cyl_radius / TMath::Tan(theta * kDeg2Rad ) ;
+ x = cylradius * TMath::Cos(phi * kDeg2Rad ) ;
+ y = cylradius * TMath::Sin(phi * kDeg2Rad ) ;
+ z = cylradius / TMath::Tan(theta * kDeg2Rad ) ;
return;
}
// Returned
// none.
- Float_t theta, phi,cyl_radius=0. ;
+ Float_t theta, phi,cylradius=0. ;
PosInAlice(absid, theta, phi) ;
if ( IsInECA(absid) )
- cyl_radius = GetIP2ECASection() ;
- else if ( IsInPRE(absid) )
- cyl_radius = GetIP2PRESection() ;
- else if ( IsInHCA(absid) )
- cyl_radius = GetIP2HCASection() ;
+ cylradius = GetIP2ECASection() ;
else
Fatal("XYZFromIndex", "Unexpected Tower section") ;
Double_t kDeg2Rad = TMath::DegToRad() ;
- v.SetX(cyl_radius * TMath::Cos(phi * kDeg2Rad ) );
- v.SetY(cyl_radius * TMath::Sin(phi * kDeg2Rad ) );
- v.SetZ(cyl_radius / TMath::Tan(theta * kDeg2Rad ) ) ;
+ v.SetX(cylradius * TMath::Cos(phi * kDeg2Rad ) );
+ v.SetY(cylradius * TMath::Sin(phi * kDeg2Rad ) );
+ v.SetZ(cylradius / TMath::Tan(theta * kDeg2Rad ) ) ;
return;
}
// Returned
// Boot_t kTRUE if the towers are neighbours otherwise false.
Boot_t anb = kFALSE;
- Int_t ieta1 = 0, ieta2 = 0, iphi1 = 0, iphi2 = 0, ipre1 = 0, ipre2 = 0;
+ Int_t ieta1 = 0, ieta2 = 0, iphi1 = 0, iphi2 = 0;
- TowerIndexes(index1,ieta1,iphi1,ipre1);
- TowerIndexes(index2,ieta2,iphi2,ipre2);
- if(ipre1!=ipre2) return anb;
+ TowerIndexes(index1,ieta1,iphi1);
+ TowerIndexes(index2,ieta2,iphi2);
if((ieta1>=ieta2-1 && ieta1<=ieta2+1) && (iphi1>=iphi2-1 &&iphi1<=iphi2+1))
anb = kTRUE;
return anb;
//*-- Author: Sahal Yacoob (LBL / UCT)
//*-- and : Yves Schutz (Subatech)
- //#include <assert.h>
+//#include <assert.h>
// --- ROOT system ---
class TString ;
static AliEMCALGeometry * GetInstance() ;
AliEMCALGeometry & operator = (const AliEMCALGeometry & /*rvalue*/) const {
// assignement operator requested by coding convention but not needed
- Fatal("operator =", "not implemented") ;
+ Fatal("operator =", "not implemented");
return *(GetInstance()) ;
};
virtual Bool_t Impact(const TParticle *) const {return kTRUE;}
// General
Bool_t IsInitialized(void) const { return fgInit ; }
- // Return EMCA geometrical parameters
+ // Return EMCA geometrical parameters
// geometry
const Float_t GetAlFrontThickness() const { return fAlFrontThick;}
const Float_t GetArm1PhiMin() const { return fArm1PhiMin ; }
const Float_t GetArm1PhiMax() const { return fArm1PhiMax ; }
const Float_t GetArm1EtaMin() const { return fArm1EtaMin;}
const Float_t GetArm1EtaMax() const { return fArm1EtaMax;}
- const Float_t GetIPDistance() const { return fIPDistance ; }
- const Float_t GetIP2PRESection() const { return (GetIPDistance() + GetAlFrontThickness() + GetGap2Active() ) ;}
- const Float_t GetIP2ECASection() const { return ( GetIP2PRESection() + GetNPRLayers() * ( GetPRScintThick() + GetPRPbRadThick() ) ) ; }
- const Float_t GetIP2HCASection() const { return ( GetIP2ECASection() + GetNECLayers() * ( GetECScintThick() + GetECPbRadThick() ) ) ; }
+ const Float_t GetIPDistance() const { return fIPDistance;}
+ const Float_t GetIP2ECASection() const { return ( GetIPDistance() + GetAlFrontThickness() + GetGap2Active() ) ; }
const Float_t GetEnvelop(Int_t index) const { return fEnvelop[index] ; }
const Float_t GetShellThickness() const { return fShellThickness ; }
const Float_t GetZLength() const { return fZLength ; }
const Float_t GetDeltaPhi() const {return (fArm1PhiMax-fArm1PhiMin)/
((Float_t)fNPhi);}
const Int_t GetNECLayers() const {return fNECLayers ;}
- const Int_t GetNHCLayers() const {return fNHCLayers ;}
- const Int_t GetNPRLayers() const {return fNPRLayers;}
const Int_t GetNZ() const {return fNZ ;}
const Int_t GetNEta() const {return fNZ ;}
const Int_t GetNPhi() const {return fNPhi ;}
const Int_t GetNTowers() const {return fNPhi * fNZ ;}
- const Float_t GetPRPbRadThick()const {return fPRPbRadThickness;}
const Float_t GetECPbRadThick()const {return fECPbRadThickness;}
- const Float_t GetHCCuRadThick()const {return fHCCuRadThickness;}
- const Float_t GetPRScintThick() const {return fPRScintThick;}
const Float_t GetECScintThick() const {return fECScintThick;}
- const Float_t GetHCScintThick() const {return fECScintThick;}
const Float_t GetSampling() const {return fSampling ; }
- const Float_t GetSummationFraction() const {return fSummationFraction ; }
-
- const Bool_t IsInPRE(Int_t index) const { if ( (index > (GetNZ() * GetNPhi()) && (index <= 2 * (GetNZ() * GetNPhi())))) return kTRUE; else return kFALSE ;}
const Bool_t IsInECA(Int_t index) const { if ( (index > 0 && (index <= GetNZ() * GetNPhi()))) return kTRUE; else return kFALSE ;}
- const Bool_t IsInHCA(Int_t index) const { if ( (index > 2*(GetNZ() * GetNPhi()) && (index <= 3 * (GetNZ() * GetNPhi())))) return kTRUE; else return kFALSE ;} ;
-
- Float_t AngleFromEta(Float_t eta){ // returns angle in radians for a given
- // pseudorapidity.
+
+ Float_t AngleFromEta(Float_t eta){ // returns theta in radians for a given pseudorapidity
return 2.0*TMath::ATan(TMath::Exp(-eta));
}
Float_t ZFromEtaR(Float_t r,Float_t eta){ // returns z in for a given
return r/TMath::Tan(AngleFromEta(eta));
}
Int_t TowerIndex(Int_t iz,Int_t iphi) const; // returns tower index
- // returns tower indexs iz, iphi.
- void TowerIndexes(Int_t index,Int_t &iz,Int_t &iphi,Int_t &ipre) const;
- // for a given tower index it returns eta and phi of center of that tower.
+ // returns tower indexs iz, iphi.
+ void TowerIndexes(Int_t index,Int_t &iz,Int_t &iphi) const;
+ // for a given tower index it returns eta and phi of center of that tower.
void EtaPhiFromIndex(Int_t index,Float_t &eta,Float_t &phi) const;
- // returns x, y, and z (cm) on the inner surface of a given EMCAL Cell specified by relid.
+ // returns x, y, and z (cm) on the inner surface of a given EMCAL Cell specified by relid.
void XYZFromIndex(const Int_t *relid,Float_t &x,Float_t &y, Float_t &z) const;
void XYZFromIndex(const Int_t absid, TVector3 &v) const;
- // for a given eta and phi in the EMCAL it returns the tower index.
+ // for a given eta and phi in the EMCAL it returns the tower index.
Int_t TowerIndexFromEtaPhi(Float_t eta,Float_t phi) const;
- // for a given eta and phi in the EMCAL it returns the pretower index.
+ // for a given eta and phi in the EMCAL it returns the pretower index.
Int_t PreTowerIndexFromEtaPhi(Float_t eta,Float_t phi) const;
- // Returns theta and phi (degree) for a given EMCAL cell indicated by relid or absid
+ // Returns theta and phi (degree) for a given EMCAL cell indicated by relid or absid
void PosInAlice(const Int_t *relid, Float_t &theta, Float_t &phi) const ;
void PosInAlice(const Int_t absid, Float_t &theta, Float_t &phi) const ;
Bool_t AbsToRelNumbering(Int_t AbsId, Int_t *relid) const;
- /*
- // Returns kTRUE if the two indexs are neighboring towers or preshowers.
- Boot_t AliEMCALGeometry::AreNeighbours(Int_t index1,Int_t index2) const;
- */
-
- void SetNZ(Int_t nz) { fNZ= nz ; Info("SetNZ", "Number of modules in Z set to %d", fNZ) ; }
- void SetNPhi(Int_t nphi) { fNPhi= nphi ; Info("SetNPhi", "Number of modules in Phi set to %d", fNPhi) ; }
- void SetSampling(Float_t samp) { fSampling = samp; Info("SetSampling", "Sampling factor set to %f", fSampling) ; }
+ void SetNZ(Int_t nz) { fNZ= nz ; printf("SetNZ: Number of modules in Z set to %d", fNZ) ; }
+ void SetNPhi(Int_t nphi) { fNPhi= nphi ; printf("SetNPhi: Number of modules in Phi set to %d", fNPhi) ; }
+ void SetSampling(Float_t samp) { fSampling = samp; printf("SetSampling: Sampling factor set to %f", fSampling) ; }
protected:
AliEMCALGeometry(const Text_t* name, const Text_t* title="") :
AliGeometry(name, title) {// ctor only for internal usage (singleton)
Init();
};
- void Init(void) ; // initializes the parameters of EMCAL
+ void Init(void); // initializes the parameters of EMCAL
private:
- static AliEMCALGeometry * fgGeom ; // pointer to the unique instance
- // of the singleton
- static Bool_t fgInit;// Tells if geometry has been succesfully set up.
- Float_t fAlFrontThick; // Thickness of the front Al face of the support box
-
- Float_t fPRPbRadThickness ; // cm, Thickness of the Pb radiators for the preshower section
- Float_t fPRScintThick ; // cm, Thickness of the sintilator for the preshower section of the tower
- Int_t fNPRLayers ; // number of scintillator layers in the preshower section
-
- Float_t fECPbRadThickness ; // cm, Thickness of the Pb radiators for the EM calorimeter section
- Float_t fECScintThick ; // cm, Thickness of the sintilator for the EM alorimeter section of the tower
- Int_t fNECLayers ; // number of scintillator layers in the EM calorimeter section
+ static AliEMCALGeometry * fgGeom; // pointer to the unique instance of the singleton
+ static Bool_t fgInit; // Tells if geometry has been succesfully set up.
+ Float_t fAlFrontThick; // Thickness of the front Al face of the support box
- Float_t fHCCuRadThickness ; // cm, Thickness of the Cu radiators.
- Float_t fHCScintThick ; // cm, Thickness of the sintilator for the hadronic alorimeter section of the tower
- Int_t fNHCLayers ; // number of scintillator layers in the hadronic calorimeter section
+ Float_t fECPbRadThickness; // cm, Thickness of the Pb radiators
+ Float_t fECScintThick; // cm, Thickness of the scintillators
+ Int_t fNECLayers; // number of scintillator layers
- Float_t fArm1PhiMin; // Minimum angular position of EMCAL in Phi (degrees)
- Float_t fArm1PhiMax; // Maximum angular position of EMCAL in Phi (degrees)
- Float_t fArm1EtaMin; // Minimum pseudorapidity position of EMCAL in Eta
- Float_t fArm1EtaMax; // Maximum pseudorapidity position of EMCAL in Eta
+ Float_t fArm1PhiMin; // Minimum angular position of EMCAL in Phi (degrees)
+ Float_t fArm1PhiMax; // Maximum angular position of EMCAL in Phi (degrees)
+ Float_t fArm1EtaMin; // Minimum pseudorapidity position of EMCAL in Eta
+ Float_t fArm1EtaMax; // Maximum pseudorapidity position of EMCAL in Eta
- // It is assumed that Arm1 and Arm2 have the same following parameters
- Float_t fEnvelop[3]; // the GEANT TUB for the detector
- Float_t fIPDistance; // Radial Distance of the inner surface of the EMCAL
- Float_t fShellThickness; // Total thickness in (x,y) direction
- Float_t fZLength; // Total length in z direction
- Float_t fGap2Active; // Gap between the envelop and the active material
- Int_t fNZ; // Number of Towers in the Z direction
- Int_t fNPhi; // Number of Towers in the Phi Direction
- Float_t fSampling; // Sampling factor
- Float_t fSummationFraction; // Fraction of the energy collected in the PRE section to be added to the EC section
+ // Geometry Parameters
+ Float_t fEnvelop[3]; // the GEANT TUB for the detector
+ Float_t fIPDistance; // Radial Distance of the inner surface of the EMCAL
+ Float_t fShellThickness; // Total thickness in (x,y) direction
+ Float_t fZLength; // Total length in z direction
+ Float_t fGap2Active; // Gap between the envelop and the active material
+ Int_t fNZ; // Number of Towers in the Z direction
+ Int_t fNPhi; // Number of Towers in the Phi Direction
+ Float_t fSampling; // Sampling factor
- ClassDef(AliEMCALGeometry,6) // EMCAL geometry class
+ ClassDef(AliEMCALGeometry,7) // EMCAL geometry class
};
return rv ;
}
-
-//____________________________________________________________________________
-TObjArray * AliEMCALGetter::PRERecPoints() const
-{
- // asks the Loader to return the EMC RecPoints container
-
- TObjArray * rv = 0 ;
-
- rv = EmcalLoader()->PRERecPoints() ;
- if (!rv) {
- EmcalLoader()->MakeRecPointsArray() ;
- rv = EmcalLoader()->PRERecPoints() ;
- }
- return rv ;
-}
-
//____________________________________________________________________________
TObjArray * AliEMCALGetter::ECARecPoints() const
{
return rv ;
}
-//____________________________________________________________________________
-TObjArray * AliEMCALGetter::HCARecPoints() const
-{
- // asks the Loader to return the EMC RecPoints container
-
- TObjArray * rv = 0 ;
-
- rv = EmcalLoader()->HCARecPoints() ;
- if (!rv) {
- EmcalLoader()->MakeRecPointsArray() ;
- rv = EmcalLoader()->HCARecPoints() ;
- }
- return rv ;
-}
-
//____________________________________________________________________________
TClonesArray * AliEMCALGetter::TrackSegments() const
{
if( strstr(opt,"P") )
ReadTreeP() ;
-
-// if( strstr(opt,"Q") )
-// ReadTreeQA() ;
}
// Creates and returns the pointer of the unique instance
// Must be called only when the environment has changed
- //::Info("Instance","alirunFileName=%s version=%s openingOption=%s",alirunFileName,version,openingOption);
-
if(!fgObjGetter){ // first time the getter is called
fgObjGetter = new AliEMCALGetter(alirunFileName, version, openingOption) ;
}
// creates the Primaries container if needed
if ( !fPrimaries ) {
if (fgDebug)
- Info("Primaries", "Creating a new TClonesArray for primaries") ;
+ printf("Primaries: Creating a new TClonesArray for primaries") ;
fPrimaries = new TClonesArray("TParticle", 1000) ;
}
return fPrimaries ;
// Print usefull information about the getter
AliRunLoader * rl = AliRunLoader::GetRunLoader(fgEmcalLoader->GetTitle());
- ::Info( "Print", "gAlice file is %s -- version name is %s", (rl->GetFileName()).Data(), rl->GetEventFolder()->GetName() ) ;
+ printf("Print: gAlice file is %s -- version name is %s", (rl->GetFileName()).Data(), rl->GetEventFolder()->GetName() ) ;
}
//____________________________________________________________________________
fNPrimaries = rl->Stack()->GetNtrack() ;
if (fgDebug)
- Info( "ReadTreeK", "Found %d particles in event # %d", fNPrimaries, EventNumber() ) ;
+ printf("ReadTreeK: Found %d particles in event # %d", fNPrimaries, EventNumber() ) ;
// first time creates the container
return EmcalLoader()->WriteDigitizer(opt) ; }
//========== RecPoints =============
- TObjArray * PRERecPoints() const ;
- AliEMCALRecPoint * PRERecPoint(const Int_t index) const { return static_cast<AliEMCALRecPoint *>(PRERecPoints()->At(index)) ;}
- TObjArray * ECARecPoints() const ;
- AliEMCALTowerRecPoint * ECARecPoint(const Int_t index) const { return static_cast<AliEMCALTowerRecPoint *>(ECARecPoints()->At(index)) ;}
- TObjArray * HCARecPoints() const ;
- AliEMCALTowerRecPoint * HCARecPoint(const Int_t index) const { return static_cast<AliEMCALTowerRecPoint *>(HCARecPoints()->At(index)) ;}
+ TObjArray * ECARecPoints() const;
+ AliEMCALTowerRecPoint * ECARecPoint(const Int_t index) const{ return static_cast<AliEMCALTowerRecPoint *>(ECARecPoints()->At(index)) ;}
TTree * TreeR() const ;
AliEMCALClusterizer * Clusterizer() ;
TString GetRecPointsFileName() const { return EmcalLoader()->GetRecPointsFileName() ; }
-/**************************************************************************
+ /**************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
* *
* Author: The ALICE Off-line Project. *
// from the same enterring Particle
Bool_t rv = kFALSE;
- if ( (fId == rValue.GetId()) && ( fIparent == rValue.GetIparent()) )
+ if ( (fId == rValue.GetId()) && ( fIparent == rValue.GetIparent()))
rv = kTRUE;
return rv;
friend ostream& operator << (ostream&,AliEMCALHit&);
public:
-
AliEMCALHit(); // default ctor
AliEMCALHit(const AliEMCALHit & hit);
- AliEMCALHit(Int_t shunt, Int_t primary, Int_t tracknumber, Int_t iparent, Float_t ienergy, Int_t id,
- Float_t *hits,Float_t *p);
+ AliEMCALHit(Int_t shunt, Int_t primary, Int_t tracknumber, Int_t iparent, Float_t ienergy, Int_t id, Float_t *hits,Float_t *p);
virtual ~AliEMCALHit(void) {}// dtor
//returns the energy loss for this hit
Float_t GetEnergy(void) const{return fELOS;}
Float_t GetTime(void) const {
// returns the time of the first energy deposition
return fTime ;}
-
+
Float_t GetPx(void) const{return fPx;}
Float_t GetPy(void) const{return fPy;}
Float_t GetPz(void) const{return fPz;}
Float_t GetPe(void) const{return fPe;}
+
Bool_t operator == (AliEMCALHit const &rValue) const;
AliEMCALHit operator + (const AliEMCALHit& rValue);
Int_t fId; // Absolute Id number EMCAL segment
Float_t fELOS; // Energy deposited
Int_t fPrimary; // Primary particles at the origin of the hit
- Float_t fPx; // Primary partical enetrence momentum/energy
- Float_t fPy; // Primary partical enetrence momentum/energy
- Float_t fPz; // Primary partical enetrence momentum/energy
- Float_t fPe; // Primary partical enetrence momentum/energy
+ Float_t fPx; // Primary particle entrance momentum/energy
+ Float_t fPy; // Primary particle entrance momentum/energy
+ Float_t fPz; // Primary particle entrance momentum/energy
+ Float_t fPe; // Primary particle entrance momentum/energy
Int_t fIparent; // Parent particle that entered emcal
Float_t fIenergy; // Initial energy of parent particle that enterred the emcal
Float_t fTime ; // Time of the energy deposition
nbytes += branchDr->GetEntry(0);
//
// Get digitizer parameters
- Float_t preADCped = digr->GetPREpedestal();
- Float_t preADCcha = digr->GetPREchannel();
Float_t ecADCped = digr->GetECApedestal();
Float_t ecADCcha = digr->GetECAchannel();
- Float_t hcADCped = digr->GetHCApedestal();
- Float_t hcADCcha = digr->GetHCAchannel();
AliEMCAL* pEMCAL = (AliEMCAL*) gAlice->GetModule("EMCAL");
AliEMCALGeometry* geom =
if (fDebug) {
Int_t ndig = digs->GetEntries();
- Info("FillFromDigits","Number of Digits: %d %d\n Parameters: PRE : %f %f EC: %f %f HC: %f %f\n Geometry: %d %d",
- ndig, nent, preADCped, preADCcha, ecADCped, ecADCcha, hcADCped, hcADCcha, geom->GetNEta(), geom->GetNPhi());
+ Info("FillFromDigits","Number of Digits: %d %d\n Parameters: EC: %f %f\n Geometry: %d %d",
+ ndig, nent, ecADCped, ecADCcha, geom->GetNEta(), geom->GetNPhi());
}
//
{
Double_t pedestal = 0.;
Double_t channel = 0.;
- if (geom->IsInPRE(sdg->GetId())) {
- pedestal = preADCped;
- channel = preADCcha;
- }
- else if (geom->IsInECA(sdg->GetId())) {
+ if (geom->IsInECA(sdg->GetId())) {
pedestal = ecADCped;
channel = ecADCcha;
- }
- else if (geom->IsInHCA(sdg->GetId())) {
- pedestal = hcADCped;
- channel = hcADCcha;
- }
+ }
else
Fatal("FillFromDigits", "unexpected digit is number!") ;
ClassImp(AliEMCALLink)
//____________________________________________________________________________
- AliEMCALLink::AliEMCALLink(Float_t prod, Int_t eca, Int_t rp, unsigned int what)
+ AliEMCALLink::AliEMCALLink(Float_t prod, Int_t eca, Int_t rp)
{
// ctor
if (gDebug == 2 )
- Info("ctor", "prod = %f, ec=%d , rp=%d, what=%d", prod, eca, rp, what) ;
+ printf("ctor: prod = %f, ec=%d , rp=%d", prod, eca, rp) ;
fProd = prod ;
fECAN = eca ;
fOtherN = rp ;
- fWhat = what ;
}
//____________________________________________________________________________
public:
- AliEMCALLink( Float_t prod, Int_t ec, Int_t rp, unsigned int what) ; // ctor
+ AliEMCALLink( Float_t prod, Int_t ec, Int_t rp) ; // ctor
virtual ~AliEMCALLink(){
// dtor
}
const Int_t GetECA(void) const { return fECAN ; }
const Int_t GetOther(void) const { return fOtherN ; }
const Float_t GetProd(void) const { return fProd ; }
- const Bool_t IsLinkToPRE(void) const {if (fWhat) return kFALSE; else return kTRUE;}
- const Bool_t IsLinkToHCA(void) const {if (fWhat) return kTRUE; else return kFALSE;}
Bool_t IsSortable() const{ return kTRUE; }
private:
Int_t fECAN ; // ECAL index
Int_t fOtherN ; // index of the linked recpoint
Float_t fProd ; // Scalar produc of the direction of the 2 recpoints
- unsigned int fWhat ; // PRE (=0) or HCAL (=1)
ClassDef(AliEMCALLink,1) // Auxilliary algorithm class used by AliEMCALTrackSegmentMaker
const TString AliEMCALLoader::fgkHitsName("HITS");//Name for TClonesArray with hits from one event
const TString AliEMCALLoader::fgkSDigitsName("SDIGITS");//Name for TClonesArray
const TString AliEMCALLoader::fgkDigitsName("DIGITS");//Name for TClonesArray
-const TString AliEMCALLoader::fgkPRERecPointsName("PRERECPOINTS");//Name for TClonesArray
const TString AliEMCALLoader::fgkECARecPointsName("ECARECPOINTS");//Name for TClonesArray
-const TString AliEMCALLoader::fgkHCARecPointsName("HCARECPOINTS");//Name for TClonesArray
const TString AliEMCALLoader::fgkTracksName("TRACKS");//Name for TClonesArray
const TString AliEMCALLoader::fgkRecParticlesName("RECPARTICLES");//Name for TClonesArray
-const TString AliEMCALLoader::fgkPRERecPointsBranchName("EMCALPRERP");//Name for branch with PreShower Reconstructed Points
const TString AliEMCALLoader::fgkECARecPointsBranchName("EMCALECARP");//Name for branch with ECA Reconstructed Points
-const TString AliEMCALLoader::fgkHCARecPointsBranchName("EMCALHCARP");//Name for branch with HCA Reconstructed Points
const TString AliEMCALLoader::fgkTrackSegmentsBranchName("EMCALTS");//Name for branch with TrackSegments
const TString AliEMCALLoader::fgkRecParticlesBranchName("EMCALRP");//Name for branch with Reconstructed Particles
Clean(fgkHitsName);
Clean(fgkSDigitsName);
Clean(fgkDigitsName);
- Clean(fgkPRERecPointsName);
Clean(fgkECARecPointsName);
- Clean(fgkHCARecPointsName);
Clean(fgkTracksName);
Clean(fgkRecParticlesName);
}
if (Hits()) Hits()->Clear();
if (SDigits()) SDigits()->Clear();
if (Digits()) Digits()->Clear();
- if (PRERecPoints()) PRERecPoints()->Clear();
if (ECARecPoints()) ECARecPoints()->Clear();
- if (HCARecPoints()) HCARecPoints()->Clear();
if (TrackSegments()) TrackSegments()->Clear();
if (RecParticles()) RecParticles()->Clear();
{
//Loads Tracks: Open File, Reads Tree and posts, Read Data and Posts
if (GetDebug())
- Info("LoadTracks","opt = %s",opt);
+ printf("LoadTracks: opt = %s",opt);
if (fTracksLoaded)
{
Warning("LoadTracks","Tracks are already loaded");
}
if (GetDebug())
- Info("ReadHits","Reading Hits");
+ printf("ReadHits: Reading Hits");
if (hitsbranch->GetEntries() > 1)
{
}
-//____________________________________________________________________________
-void AliEMCALLoader::ReadTreeQA()
-{
- // Read the digit tree gAlice->TreeQA()
- // so far only EMCAL knows about this Tree
-
- if(EMCAL()->TreeQA()== 0){
- cerr << "ERROR: AliEMCALLoader::ReadTreeQA: can not read TreeQA " << endl ;
- return ;
- }
-
- TBranch * qabranch = EMCAL()->TreeQA()->GetBranch("EMCAL");
- if (!qabranch) {
- if (fDebug)
- cout << "WARNING: AliEMCALLoader::ReadTreeQA -> Cannot find QA Alarms for EMCAL" << endl ;
- return ;
- }
-
- // if(!Alarms()) PostQA();
-
- qabranch->SetAddress(AlarmsRef()) ;
-
- qabranch->GetEntry(0) ;
-
-}
-
//____________________________________________________________________________
Int_t AliEMCALLoader::ReadRecPoints()
{
//connects branch in tree (if exists), and reads data to array
MakeRecPointsArray();
-
- TObjArray * pre = 0x0 ;
- TObjArray * eca = 0x0 ;
- TObjArray * hca = 0x0 ;
+
+ TObjArray * eca = 0x0 ;
TTree * treeR = TreeR();
}
Int_t retval = 0;
- TBranch * prebranch = treeR->GetBranch(fgkPRERecPointsBranchName);
-
- if (prebranch == 0x0)
- {
- Error("ReadRecPoints","Can not get branch with PRE Rec. Points named %s",fgkPRERecPointsBranchName.Data());
- retval = 1;
- }
- else
- {
- prebranch->SetAddress(&pre) ;
- prebranch->GetEntry(0) ;
- }
+
TBranch * ecabranch = treeR->GetBranch(fgkECARecPointsBranchName);
if (ecabranch == 0x0)
{
retval = 2;
}
else
- {
- ecabranch->SetAddress(&eca);
- ecabranch->GetEntry(0) ;
- }
- TBranch * hcabranch = treeR->GetBranch(fgkHCARecPointsBranchName);
- if (hcabranch == 0x0)
- {
- Error("ReadRecPoints","Can not get branch with HCA Rec. Points named %s",fgkHCARecPointsBranchName.Data());
- retval = 2;
- }
- else
- {
- hcabranch->SetAddress(&hca);
- hcabranch->GetEntry(0) ;
- }
-
+ {
+ ecabranch->SetAddress(&eca);
+ ecabranch->GetEntry(0) ;
+ }
+
+
Int_t ii ;
- Int_t maxpre = pre->GetEntries() ;
- for ( ii= 0 ; ii < maxpre ; ii++ )
- PRERecPoints()->Add(pre->At(ii)) ;
-
+
Int_t maxeca = eca->GetEntries() ;
for ( ii= 0 ; ii < maxeca ; ii++ )
ECARecPoints()->Add(eca->At(ii)) ;
-
- Int_t maxhca = hca->GetEntries() ;
- for ( ii= 0 ; ii < maxhca ; ii++ )
- HCARecPoints()->Add(hca->At(ii)) ;
-
+
return retval;
}
dataname = GetDetectorName();
zername = "AliEMCALDigitizer" ;
}
- else if(recName == "PRERecPoints"){
- tree = TreeR();
- dataname = fgkPRERecPointsBranchName;
- zername = "AliEMCALClusterizer" ;
- }
else if(recName == "ECARecPoints"){
tree = TreeR();
dataname = fgkECARecPointsBranchName;
zername = "AliEMCALClusterizer" ;
}
- else if(recName == "HCARecPoints"){
- tree = TreeR();
- dataname = fgkHCARecPointsBranchName;
- zername = "AliEMCALClusterizer" ;
- }
else if(recName == "TrackSegments"){
tree = TreeT();
dataname = fgkTrackSegmentsBranchName;
{
// Cleans the RecPoints array in the folder structure
AliLoader::CleanRecPoints();
- TObjArray* recpoints = PRERecPoints();
- if (recpoints)
- recpoints->Clear();
- recpoints = ECARecPoints();
- if (recpoints)
- recpoints->Clear();
- recpoints = HCARecPoints();
- if (recpoints)
- recpoints->Clear();
+ TObjArray* recpoints = ECARecPoints();
+ if (recpoints) recpoints->Clear();
+
}
//____________________________________________________________________________
//____________________________________________________________________________
void AliEMCALLoader::MakeRecPointsArray()
-{
- // Create the array for RecPoints
- if ( PRERecPoints() == 0x0) {
- if (GetDebug()>9)
- Info("MakeRecPointsArray","Making array for PRE");
- TObjArray* pre = new TObjArray(100) ;
- pre->SetName(fgkPRERecPointsName) ;
- GetDetectorDataFolder()->Add(pre);
- }
+{
if ( ECARecPoints() == 0x0) {
if (GetDebug()>9)
- Info("MakeRecPointsArray","Making array for ECA");
+ printf("MakeRecPointsArray: Making array for ECA");
TObjArray* eca = new TObjArray(100) ;
eca->SetName(fgkECARecPointsName) ;
GetDetectorDataFolder()->Add(eca);
- }
- if ( HCARecPoints() == 0x0) {
- if (GetDebug()>9)
- Info("MakeRecPointsArray","Making array for HCA");
- TObjArray* hca = new TObjArray(100) ;
- hca->SetName(fgkHCARecPointsName) ;
- GetDetectorDataFolder()->Add(hca);
- }
+ }
}
//____________________________________________________________________________
// Int_t WriteRecParticles(Option_t* opt="");//writes the reconstructed particles
// Int_t WritePID(Option_t* opt="");//writes the task for PID to file
// Bool_t PostPID (AliEMCALPID * pid) const {return kTRUE;}
-// Bool_t PostQA (void) const ; //it was empty anyway
-
+
/*******************************************************************/
/*******************************************************************/
/*******************************************************************/
TObject** HitsRef(){return GetDetectorDataRef(Hits());}
TObject** SDigitsRef(){return GetDetectorDataRef(SDigits());}
TObject** DigitsRef(){return GetDetectorDataRef(Digits());}
- TObject** PRERecPointsRef(){return GetDetectorDataRef(PRERecPoints());}
TObject** ECARecPointsRef(){return GetDetectorDataRef(ECARecPoints());}
- TObject** HCARecPointsRef(){return GetDetectorDataRef(HCARecPoints());}
TObject** TracksRef(){return GetDetectorDataRef(TrackSegments());}
TObject** RecParticlesRef(){return GetDetectorDataRef(RecParticles());}
- TObject** AlarmsRef(){return GetDetectorDataRef(Alarms());}
+
void Track(Int_t itrack) ;
static AliEMCALGeometry* GetEMCALGeometry();
const AliEMCAL * EMCAL();
const AliEMCALGeometry *EMCALGeometry() ;
- // Alarms
- // TFolder * Alarms() const { return (TFolder*)(ReturnO("Alarms", 0)); }
- TObjArray * Alarms();
/*********************************************/
/************ TClonesArrays ***********/
const AliEMCALDigit * Digit(Int_t index);
void MakeDigitsArray();
/**** R e c P o i n t s ****/
- TObjArray * PRERecPoints();
TObjArray * ECARecPoints();
- TObjArray * HCARecPoints();
- const AliEMCALRecPoint * PRERecPoint(Int_t index) ;
const AliEMCALTowerRecPoint * ECARecPoint(Int_t index) ;
- const AliEMCALTowerRecPoint * HCARecPoint(Int_t index) ;
void MakeRecPointsArray();
/**** T r a c k S e g m e n t s ****/
TClonesArray * TrackSegments();
static TString HitsName() { return fgkHitsName ; } //Name for TClonesArray with hits from one event
static TString SDigitsName() { return fgkSDigitsName ;} //Name for TClonesArray
static TString DigitsName() { return fgkDigitsName ;} //Name for TClonesArray
- static TString PreRecPointsName() { return fgkPRERecPointsName ;} //Name for TClonesArray
- static TString ECARecPointsName() { return fgkECARecPointsName ;} //Name for TClonesArray
- static TString HCARecPointsName() { return fgkHCARecPointsName ;} //Name for TClonesArray
+ static TString ECARecPointsName() { return fgkECARecPointsName ;} //Name for TClonesArray y
static TString TracksName() { return fgkTracksName ;} //Name for TClonesArray
static TString RecParticlesName() { return fgkRecParticlesName ;} //Name for TClonesArray
- static TString PRERecPointsBranchName() { return fgkPRERecPointsBranchName ;} //Name for branch
static TString ECARecPointsBranchName() { return fgkECARecPointsBranchName ;} //Name for branch
- static TString HCARecPointsBranchName() { return fgkHCARecPointsBranchName ;} //Name for branch
static TString TrackSegmentsBranchName() { return fgkTrackSegmentsBranchName ;} //Name for branch
static TString RecParticlesBranchName() { return fgkRecParticlesBranchName ;} //Name for branch
Int_t ReadRecPoints();
Int_t ReadTracks();
Int_t ReadRecParticles();
-
- void ReadTreeQA() ;
static const TString fgkHitsName;//Name for TClonesArray with hits from one event
static const TString fgkSDigitsName;//Name for TClonesArray
static const TString fgkDigitsName;//Name for TClonesArray
- static const TString fgkPRERecPointsName;//Name for TClonesArray
static const TString fgkECARecPointsName;//Name for TClonesArray
- static const TString fgkHCARecPointsName;//Name for TClonesArray
static const TString fgkTracksName;//Name for TClonesArray
static const TString fgkRecParticlesName;//Name for TClonesArray
- static const TString fgkPRERecPointsBranchName;//Name for branch
static const TString fgkECARecPointsBranchName;//Name for branch
- static const TString fgkHCARecPointsBranchName;//Name for branch
static const TString fgkTrackSegmentsBranchName;//Name for branch
static const TString fgkRecParticlesBranchName;//Name for branch
Int_t fDebug ; // Debug level
/******************************************************************************/
-inline TObjArray * AliEMCALLoader::PRERecPoints()
-{
- return dynamic_cast<TObjArray*>(GetDetectorData(fgkPRERecPointsName));
-}
-/******************************************************************************/
-
-inline const AliEMCALRecPoint * AliEMCALLoader::PRERecPoint(Int_t index)
-{
- TObjArray* tcarr = PRERecPoints();
- if (tcarr)
- return dynamic_cast<const AliEMCALRecPoint*>(tcarr->At(index));
- return 0x0;
-}
-
-/******************************************************************************/
-
inline TObjArray * AliEMCALLoader::ECARecPoints()
{
return dynamic_cast<TObjArray*>(GetDetectorData(fgkECARecPointsName));
/******************************************************************************/
-inline TObjArray * AliEMCALLoader::HCARecPoints()
-{
- return dynamic_cast<TObjArray*>(GetDetectorData(fgkHCARecPointsName));
-}
-
-/******************************************************************************/
-
-inline const AliEMCALTowerRecPoint * AliEMCALLoader::HCARecPoint(Int_t index)
-{
- TObjArray* tcarr = HCARecPoints();
- if (tcarr)
- return dynamic_cast<const AliEMCALTowerRecPoint*>(tcarr->At(index));
- return 0x0;
-}
-
-/******************************************************************************/
-
inline TClonesArray * AliEMCALLoader::TrackSegments()
{
return dynamic_cast<TClonesArray*>(GetDetectorData(fgkTracksName));
return (const AliEMCALRecParticle*) tcarr->At(index);
return 0x0;
}
-/******************************************************************************/
-inline TObjArray * AliEMCALLoader::Alarms()
-{ return (TObjArray*)(GetQAFolder()->FindObject(fDetectorName));}
#endif // AliEMCALLOADER_H
}
if(strstr(option,"tim")){
gBenchmark->Stop("EMCALPID");
- Info("Exec", "took %f seconds for PID %f seconds per event",
+ printf("Exec: took %f seconds for PID %f seconds per event",
gBenchmark->GetCpuTime("EMCALPID"),
gBenchmark->GetCpuTime("EMCALPID")/nevents) ;
}
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
TObjArray * aECARecPoints = gime->ECARecPoints() ;
- TObjArray * aPRERecPoints = gime->PRERecPoints() ;
- TObjArray * aHCARecPoints = gime->HCARecPoints() ;
TClonesArray * trackSegments = gime->TrackSegments() ;
- if ( !aECARecPoints || !aPRERecPoints || !aHCARecPoints || !trackSegments ) {
+ if ( !aECARecPoints || !trackSegments ) {
Fatal("MakeRecParticles", "RecPoints or TrackSegments not found !") ;
}
TClonesArray * recParticles = gime->RecParticles() ;
AliEMCALTowerRecPoint * eca = 0 ;
if(ts->GetECAIndex()>=0)
eca = dynamic_cast<AliEMCALTowerRecPoint *>(aECARecPoints->At(ts->GetECAIndex())) ;
-
- AliEMCALTowerRecPoint * pre = 0 ;
- if(ts->GetPREIndex()>=0)
- pre = dynamic_cast<AliEMCALTowerRecPoint *>(aPRERecPoints->At(ts->GetPREIndex())) ;
-
- AliEMCALTowerRecPoint * hca = 0 ;
- if(ts->GetHCAIndex()>=0)
- hca = dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(ts->GetHCAIndex())) ;
// Now set type (reconstructed) of the particle
{
// Print the parameters used for the particle type identification
- Info("Print", "=============== AliEMCALPID1 ================") ;
+ printf("Print: =============== AliEMCALPID1 ================") ;
printf("Making PID\n");
printf(" Pricipal analysis file from 0.5 to 100 %s\n", fFileName.Data() ) ;
printf(" Name of parameters file %s\n", fFileNamePar.Data() ) ;
TClonesArray * recParticles = gime->RecParticles() ;
- TString message ;
- message = "\nevent " ;
- message += gAlice->GetEvNumber() ;
- message += " found " ;
- message += recParticles->GetEntriesFast();
- message += " RecParticles\n" ;
+ printf("\nevent %i", gAlice->GetEvNumber());
+ printf(" found %i", recParticles->GetEntriesFast());
+ printf(" RecParticles\n");
if(strstr(option,"all")) { // printing found TS
- message += "\n PARTICLE Index \n" ;
+ printf("\n PARTICLE Index \n");
Int_t index ;
for (index = 0 ; index < recParticles->GetEntries() ; index++) {
AliEMCALRecParticle * rp = (AliEMCALRecParticle * ) recParticles->At(index) ;
- message += "\n" ;
- message += rp->Name().Data() ;
- message += " " ;
- message += rp->GetIndexInList() ;
- message += " " ;
- message += rp->GetType() ;
+ printf("\n");
+ printf(rp->Name().Data());
+ printf(" %i", rp->GetIndexInList());
+ printf(" %i", rp->GetType());
}
}
- Info("Print", message.Data() ) ;
}
//____________________________________________________________________________
fMaxTrack = 0 ;
fTheta = fPhi = 0. ;
fEMCALArm = 0;
- fPRESection = fECASection = fHCASection = kFALSE ;
+ fECASection = kFALSE ;
}
//____________________________________________________________________________
AliEMCALGeometry * emcalgeom = (AliEMCALGetter::Instance())->EMCALGeometry() ;
Int_t iDigit;
- Int_t relid[4] ;
+ Int_t relid[3] ;
const Int_t kMulDigit=AliEMCALRecPoint::GetDigitsMultiplicity() ;
Float_t * xi = new Float_t [kMulDigit] ;
if( fEMCALArm == 0){
- Int_t relid[4] ;
+ Int_t relid[3] ;
AliEMCALGeometry * emcalgeom = (AliEMCALGetter::Instance())->EMCALGeometry();
gpos.SetX(fPhi) ;
if ( IsInECA() )
gpos.SetY(emcalgeom->GetIP2ECASection()) ;
- else if ( IsInPRE() )
- gpos.SetY(emcalgeom->GetIP2PRESection()) ;
- else if ( IsInHCA() )
- gpos.SetY(emcalgeom->GetIP2HCASection()) ;
else
Fatal("GetGlobalPosition", "Unexpected tower section") ;
gpos.SetZ(fTheta) ;
return fTracksList ; }
virtual Bool_t IsEmc(void)const { return kTRUE ; }
const Bool_t IsInECA(void) const { return fECASection ; }
- const Bool_t IsInHCA(void) const { return fHCASection ; }
- const Bool_t IsInPRE(void) const { return fPRESection ; }
virtual Bool_t IsSortable() const {
// tells that this is a sortable object
return kTRUE ;
}
void SetECA() { fECASection = kTRUE ; }
- void SetHCA() { fHCASection = kTRUE ; }
- void SetPRE() { fPRESection = kTRUE ; }
AliEMCALRecPoint & operator = (const AliEMCALRecPoint & ) {
Fatal("operator =", "not implemented") ;
return *this ;
Float_t fTheta ; // theta angle in Alice
Float_t fPhi ; // phi angle in Alice
Bool_t fECASection ; // tells if the recpoint is in ECAL section
- Bool_t fHCASection ; // tells if the recpoint is in HCAL section
- Bool_t fPRESection ; // tells if the recpoint is in PRE section
ClassDef(AliEMCALRecPoint,3) // RecPoint for EMCAL (Base Class)
void AliEMCALReconstructioner::Print(Option_t * option)const {
// Print reconstructioner data
- TString message("\n") ;
+ printf("\n") ;
- message += " Reconstruction of the header file " ;
- message += fHeaderFileName.Data() ;
- message += "\n with the following modules:\n" ;
+ printf(" Reconstruction of the header file ");
+ printf(fHeaderFileName.Data());
+ printf("\n with the following modules:\n");
if(fSDigitizer->IsActive()){
- message += " (+) " ;
- message += fSDigitizer->GetName() ;
- message +=" to branch : " ;
- message += fSDigitsBranch.Data() ;
+ printf(" (+) ");
+ printf(fSDigitizer->GetName());
+ printf(" to branch : ");
+ printf(fSDigitsBranch.Data());
}
if(fDigitizer->IsActive()){
- message += "\n (+) " ;
- message += fDigitizer->GetName() ;
- message += " to branch : " ;
- message += fDigitsBranch.Data() ;
+ printf("\n (+) ");
+ printf(fDigitizer->GetName());
+ printf(" to branch : ");
+ printf(fDigitsBranch.Data());
}
if(fClusterizer->IsActive()){
- message += "\n (+) " ;
- message += fClusterizer->GetName() ;
- message += " to branch : " ;
- message += fRecPointBranch.Data() ;
- }
- Info("Print", message.Data() ) ;
+ printf("\n (+) ");
+ printf(fClusterizer->GetName());
+ printf(" to branch : ");
+ printf(fRecPointBranch.Data());
+ }
}
#include "TTask.h"
class AliEMCALDigitizer ;
-class AliEMCALClusterizer ;
+class AliEMCALClusterizerv1 ;
//class AliEMCALTrackSegmentMaker ;
class AliEMCALPID ;
class AliEMCALSDigitizer ;
virtual void Exec(Option_t * option) ;
AliEMCALDigitizer * GetDigitizer() const { return fDigitizer ; }
- AliEMCALClusterizer * GetClusterizer()const { return fClusterizer ; }
+ AliEMCALClusterizerv1 * GetClusterizer()const { return fClusterizer ; }
//AliEMCALPID * GetPID() const { return fPID; }
//AliEMCALTrackSegmentMaker * GetTSMaker() const { return fTSMaker ; }
AliEMCALSDigitizer * GetSDigitizer() const { return fSDigitizer ; }
AliEMCALDigitizer * fDigitizer ; //! Pointer to AliEMCALDigitizer
- AliEMCALClusterizer * fClusterizer ; //! Pointer to AliEMCALClusterizer
+ AliEMCALClusterizerv1 * fClusterizer ; //! Pointer to AliEMCALClusterizer
//AliEMCALPID * fPID ; //! Pointer to AliEMCALPID
//AliEMCALTrackSegmentMaker * fTSMaker ; //! Pointer to AliEMCALTrackSegmentMaker
AliEMCALSDigitizer * fSDigitizer ; //! Pointer to AliEMCALSDigitizer
Bool_t fIsInitialized ; // kTRUE if reconstructioner is initialized
- ClassDef(AliEMCALReconstructioner,2) // Reconstruction algorithm class (Base Class)
+ ClassDef(AliEMCALReconstructioner,3) // Reconstruction algorithm class (Base Class)
};
fA = sd.fA ;
fB = sd.fB ;
- fPREPrimThreshold = sd.fPREPrimThreshold ;
fECPrimThreshold = sd.fECPrimThreshold ;
- fHCPrimThreshold = sd.fHCPrimThreshold ;
fSDigitsInRun = sd.fSDigitsInRun ;
SetName(sd.GetName()) ;
SetTitle(sd.GetTitle()) ;
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
const AliEMCALGeometry * geom = gime->EMCALGeometry() ;
if (geom->GetSampling() == 0.) {
- Error("InitParameters", "Sampling factor not set !") ;
+ printf("InitParameters: Sampling factor not set !") ;
abort() ;
}
else
- Info("InitParameters", "Sampling factor set to %f\n", geom->GetSampling()) ;
+ printf("InitParameters: Sampling factor set to %f\n", geom->GetSampling()) ;
// this threshold corresponds approximately to 100 MeV
- fECPrimThreshold = 100E-3 / ( geom->GetSampling() * ( geom->GetNPRLayers() + geom->GetNECLayers()) ) * geom->GetNECLayers() ;
- fPREPrimThreshold = 100E-3 / ( geom->GetSampling() * ( geom->GetNPRLayers() + geom->GetNECLayers()) ) * geom->GetNPRLayers() ;
- fHCPrimThreshold = fECPrimThreshold/5. ; // 5 is totally arbitrary
-
+ fECPrimThreshold = 100E-3 / ( geom->GetSampling() * geom->GetNECLayers()) * geom->GetNECLayers() ;
}
//____________________________________________________________________________
void AliEMCALSDigitizer::Exec(Option_t *option)
{
- // Collects all hits in the section (PRE/ECAL/HCAL) of the same tower into digit
-
+ // Collects all hit of the same tower into digits
if (strstr(option, "print") ) {
Print() ;
return ;
if (!fInit) { // to prevent overwrite existing file
Error( "Exec", "Give a version name different from %s", fEventFolderName.Data() ) ;
return ;
- }
-
+ }
+
Int_t nevents = gime->MaxEvent() ;
Int_t ievent ;
for(ievent = 0; ievent < nevents; ievent++){
-
gime->Event(ievent,"H") ;
TTree * treeS = gime->TreeS();
// Attention nPrim is the number of primaries tracked by Geant
// and this number could be different to the number of Primaries in TreeK;
Int_t iprim ;
-
for ( iprim = 0 ; iprim < nPrim ; iprim++ ) {
//=========== Get the EMCAL branch from Hits Tree for the Primary iprim
gime->Track(iprim) ;
Bool_t newsdigit = kTRUE;
// Assign primary number only if deposited energy is significant
-
AliEMCALGeometry * geom = gime->EMCALGeometry() ;
-
- if( geom->IsInPRE(hit->GetId()) )
- if( hit->GetEnergy() > fPREPrimThreshold )
- curSDigit = new AliEMCALDigit( hit->GetPrimary(),
- hit->GetIparent(), hit->GetId(),
- Digitize(hit->GetEnergy()), hit->GetTime() ) ;
- else
- curSDigit = new AliEMCALDigit( -1 ,
- -1 ,
- hit->GetId(),
- Digitize(hit->GetEnergy()), hit->GetTime() ) ;
- else if( geom->IsInECA(hit->GetId()) )
+ if( geom->IsInECA(hit->GetId()) ){
if( hit->GetEnergy() > fECPrimThreshold )
curSDigit = new AliEMCALDigit( hit->GetPrimary(),
hit->GetIparent(), hit->GetId(),
-1 ,
hit->GetId(),
Digitize(hit->GetEnergy()), hit->GetTime() ) ;
- else if( geom->IsInHCA(hit->GetId()) )
- if( hit->GetEnergy() > fHCPrimThreshold )
-
- curSDigit = new AliEMCALDigit( hit->GetPrimary(),
- hit->GetIparent(), hit->GetId(),
- Digitize(hit->GetEnergy()), hit->GetTime() ) ;
- else
- curSDigit = new AliEMCALDigit( -1 ,
- -1 ,
- hit->GetId(),
- Digitize(hit->GetEnergy()), hit->GetTime() ) ;
-
+ }
+ else{
+ newsdigit = kFALSE;
+ }
+
Int_t check = 0 ;
+
for(check= 0; check < nSdigits ; check++) {
sdigit = dynamic_cast<AliEMCALDigit *>(sdigits->At(check)) ;
- if( sdigit->GetId() == curSDigit->GetId()) { // Are we in the same ECAL/HCAL/preshower tower ?
+ if(curSDigit != 0){
+ if( sdigit->GetId() == curSDigit->GetId()) { // Are we in the same ECAL tower ?
*sdigit = *sdigit + *curSDigit;
newsdigit = kFALSE;
+ }
}
}
- if (newsdigit) {
+ if (newsdigit) {
new((*sdigits)[nSdigits]) AliEMCALDigit(*curSDigit);
nSdigits++ ;
}
delete curSDigit ;
- } // loop over all hits (hit = deposited energy/layer/entering particle)
+ } // loop over all hits (hit = deposited energy/entering particle)
} // loop over iprim
-
sdigits->Sort() ;
nSdigits = sdigits->GetEntriesFast() ;
}
// Now write SDigits
-
//First list of sdigits
Int_t bufferSize = 32000 ;
if(strstr(option,"tim")){
gBenchmark->Stop("EMCALSDigitizer");
- Info("Exec", "took %f seconds for SDigitizing %f seconds per event",
+ printf("Exec: took %f seconds for SDigitizing %f seconds per event",
gBenchmark->GetCpuTime("EMCALSDigitizer"), gBenchmark->GetCpuTime("EMCALSDigitizer") ) ;
}
}
void AliEMCALSDigitizer::Print()const
{
// Prints parameters of SDigitizer
- Info("Print", "\n------------------- %s -------------", GetName() ) ;
+ printf("Print: \n------------------- %s -------------", GetName() ) ;
printf(" Writing SDigits to branch with title %s\n", fEventFolderName.Data()) ;
printf(" with digitization parameters A = %f\n", fA) ;
printf(" B = %f\n", fB) ;
- printf(" Threshold for PRE Primary assignment= %f\n", fPREPrimThreshold) ;
printf(" Threshold for EC Primary assignment= %f\n", fECPrimThreshold) ;
- printf(" Threshold for HC Primary assignment= %f\n", fHCPrimThreshold) ;
printf("---------------------------------------------------\n") ;
}
{
// Equal operator.
// SDititizers are equal if their pedestal, slope and threshold are equal
-
if( (fA==sd.fA)&&(fB==sd.fB)&&
- (fECPrimThreshold==sd.fECPrimThreshold) &&
- (fHCPrimThreshold==sd.fHCPrimThreshold) &&
- (fPREPrimThreshold==sd.fPREPrimThreshold))
+ (fECPrimThreshold==sd.fECPrimThreshold))
return kTRUE ;
else
return kFALSE ;
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
const TClonesArray * sdigits = gime->SDigits() ;
- TString message("\n") ;
- message += "event " ;
- message += gAlice->GetEvNumber() ;
- message += "\n Number of entries in SDigits list " ;
- message += sdigits->GetEntriesFast() ;
+ printf("\n") ;
+ printf("event %i", gAlice->GetEvNumber()) ;
+ printf("\n Number of entries in SDigits list %i", sdigits->GetEntriesFast());
if(strstr(option,"all")||strstr(option,"EMC")){
//loop over digits
AliEMCALDigit * digit;
- message += "\n Id Amplitude Time Index Nprim: Primaries list \n" ;
+ printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
Int_t index ;
char * tempo = new char[8192];
for (index = 0 ; index < sdigits->GetEntries() ; index++) {
digit = dynamic_cast<AliEMCALDigit *>( sdigits->At(index) ) ;
sprintf(tempo, "\n%6d %8d %6.5e %4d %2d :",
digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
- message += tempo ;
+ printf(tempo);
Int_t iprimary;
for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
sprintf(tempo, "%d ",digit->GetPrimary(iprimary+1) ) ;
- message += tempo ;
+ printf(tempo);
}
}
delete tempo ;
}
- Info("PrintSDigits", message.Data() ) ;
}
//____________________________________________________________________________
private:
Float_t fA ; // Pedestal parameter
Float_t fB ; // Slope Digitizition parameters
- Float_t fPREPrimThreshold ; // To store primary if Pre Shower Elos > threshold
Float_t fECPrimThreshold ; // To store primary if EC Shower Elos > threshold
- Float_t fHCPrimThreshold ; // To store primary if HC Shower Elos > threshold
Bool_t fDefaultInit; //! Says if the task was created by defaut ctor (only parameters are initialized)
TString fEventFolderName; // event folder name
Bool_t fInit ; //! tells if initialisation wennt OK, will revent exec if not
AliEMCALGeometry * phosgeom = (AliEMCALGetter::Instance())->EMCALGeometry();
- Int_t relid1[4] ;
+ Int_t relid1[3] ;
phosgeom->AbsToRelNumbering(digit1->GetId(), relid1) ;
- Int_t relid2[4] ;
+ Int_t relid2[3] ;
phosgeom->AbsToRelNumbering(digit2->GetId(), relid2) ;
- Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
- Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
+ Int_t rowdiff = TMath::Abs( relid1[1] - relid2[1] ) ;
+ Int_t coldiff = TMath::Abs( relid1[2] - relid2[2] ) ;
if (( coldiff <= 1 ) && ( rowdiff <= 1 ) && (coldiff + rowdiff > 0))
aren = kTRUE ;
// case kButton1Down: {
// AliEMCALDigit * digit ;
// Int_t iDigit;
-// Int_t relid[4] ;
+// Int_t relid[3] ;
// const Int_t kMulDigit = AliEMCALTowerRecPoint::GetDigitsMultiplicity() ;
// Float_t * xi = new Float_t[kMulDigit] ;
const Float_t kDeg2Rad = TMath::DegToRad() ;
Float_t cyl_radius = 0 ;
-
- if (IsInPRE())
- cyl_radius = emcalgeom->GetIP2PRESection() ;
- else if (IsInECA())
+
+ if (IsInECA())
cyl_radius = emcalgeom->GetIP2ECASection() ;
- else if (IsInHCA())
- cyl_radius = emcalgeom->GetIP2HCASection() ;
else
Fatal("EvalDispersion", "Unexpected tower section!") ;
Float_t z = fLocPos.Z() ;
if (gDebug == 2)
- Info("EvalDispersion", "x,y,z = %f,%f,%f", x, y, z) ;
+ printf("EvalDispersion: x,y,z = %f,%f,%f", x, y, z) ;
// Calculates the dispersion in coordinates
wtot = 0.;
Float_t zi = cyl_radius / TMath::Tan(thetai * kDeg2Rad ) ;
if (gDebug == 2)
- Info("EvalDispersion", "id = %d, xi,yi,zi = %f,%f,%f", digit->GetId(), xi, yi, zi) ;
+ printf("EvalDispersion: id = %d, xi,yi,zi = %f,%f,%f", digit->GetId(), xi, yi, zi) ;
Float_t w = TMath::Max(0.,logWeight+TMath::Log(fEnergyList[iDigit]/fAmp ) ) ;
d += w * ( (xi-x)*(xi-x) + (zi-z)*(zi-z) ) ;
Float_t cyl_radius = 0 ;
- if (IsInPRE())
- cyl_radius = emcalgeom->GetIP2PRESection() ;
- else if (IsInECA())
+ if (IsInECA())
cyl_radius = emcalgeom->GetIP2ECASection() ;
- else if (IsInHCA())
- cyl_radius = emcalgeom->GetIP2HCASection() ;
else
Fatal("EvalDispersion", "Unexpected tower section!") ;
// Calculates the center of gravity in the local EMCAL-module coordinates
Float_t wtot = 0. ;
- // Int_t relid[4] ;
+ // Int_t relid[3] ;
AliEMCALDigit * digit ;
AliEMCALGeometry * emcalgeom = (AliEMCALGetter::Instance())->EMCALGeometry();
Float_t cyl_radius = 0 ;
- if (IsInPRE())
- cyl_radius = emcalgeom->GetIP2PRESection() ;
- else if (IsInECA())
+ if (IsInECA())
cyl_radius = emcalgeom->GetIP2ECASection() ;
- else if (IsInHCA())
- cyl_radius = emcalgeom->GetIP2HCASection() ;
else
Fatal("EvalGlobalPosition", "Unexpected tower section!") ;
fLocPos.SetZ(z) ;
if (gDebug==2)
- Info("EvalGlobalPosition", "x,y,z = %f,%f,%f", fLocPos.X(), fLocPos.Y(), fLocPos.Z()) ;
-
-
+ printf("EvalGlobalPosition: x,y,z = %f,%f,%f", fLocPos.X(), fLocPos.Y(), fLocPos.Z()) ;
fLocPosM = 0 ;
}
AliEMCALDigit * digit ;
AliEMCALDigit * digitN ;
-
Int_t iDigitN ;
Int_t iDigit ;
for(iDigit = 0; iDigit < fMulDigit; iDigit++)
maxAt[iDigit] = (AliEMCALDigit*) digits->At(fDigitsList[iDigit]) ;
-
for(iDigit = 0 ; iDigit < fMulDigit; iDigit++) {
if(maxAt[iDigit]) {
{
// Print the list of digits belonging to the cluster
- TString message("\n") ;
+ printf("\n") ;
Int_t iDigit;
- message += "digits # = " ;
+ printf("digits # = ");
for(iDigit=0; iDigit<fMulDigit; iDigit++) {
- message += fDigitsList[iDigit] ;
- message += " " ;
+ printf("%i ", fDigitsList[iDigit]);
}
- message += "\nEnergies = " ;
+ printf("\nEnergies = ");
for(iDigit=0; iDigit<fMulDigit; iDigit++) {
- message += fEnergyList[iDigit] ;
- message += " " ;
+ printf("%f ", fEnergyList[iDigit]);
}
- message += "\nPrimaries " ;
+ printf("\nPrimaries ");
for(iDigit = 0;iDigit < fMulTrack; iDigit++) {
- message += fTracksList[iDigit] ;
- message += " " ;
+ printf("%i ", fTracksList[iDigit]);
}
- message += "\n Multiplicity = " ;
- message += fMulDigit ;
- message += "\n Cluster Energy = " ;
- message += fAmp ;
- message += "\n Number of primaries " ;
- message += fMulTrack ;
- message += "\n Stored at position " ;
- message += GetIndexInList() ;
-
- Info("Print", message.Data() ) ;
+ printf("\n Multiplicity = %i", fMulDigit);
+ printf("\n Cluster Energy = %f", fAmp);
+ printf("\n Number of primaries %i", fMulTrack);
+ printf("\n Stored at position: %i", GetIndexInList());
}
//____________________________________________________________________________
// spherical coordinates of recpoint in Alice reference frame
if (gDebug == 2)
- Info("XYZInAlice", "this= %d , r = %f, theta = %f, phi = %f", this, r, theta, phi) ;
+ printf("XYZInAlice: r = %f, theta = %f, phi = %f", r, theta, phi) ;
if (theta == 9999. || phi == 9999. || r == 9999.) {
TVector3 globalpos;
// --- Standard library ---
// --- AliRoot header files ---
-
-#include "AliEMCALTrackSegment.h"
+#include "AliEMCALTrackSegment.h"
ClassImp(AliEMCALTrackSegment)
//____________________________________________________________________________
-AliEMCALTrackSegment::AliEMCALTrackSegment( AliEMCALTowerRecPoint * eca, AliEMCALTowerRecPoint * pre, AliEMCALTowerRecPoint * hca)
+AliEMCALTrackSegment::AliEMCALTrackSegment( AliEMCALTowerRecPoint * eca)
{
// ctor
-
- if( pre )
- fPRERecPoint = pre->GetIndexInList() ;
- else
- fPRERecPoint = -1 ;
-
if( eca )
fECARecPoint = eca->GetIndexInList() ;
else
fECARecPoint = -1 ;
-
- if( hca )
- fHCARecPoint = hca->GetIndexInList() ;
- else
- fHCARecPoint = -1 ;
-
fIndexInList = -1 ;
}
// Copy of a track segment into another track segment
TObject::Copy(obj) ;
- ( (AliEMCALTrackSegment &)obj ).fPRERecPoint = fPRERecPoint ;
( (AliEMCALTrackSegment &)obj ).fECARecPoint = fECARecPoint ;
- ( (AliEMCALTrackSegment &)obj ).fHCARecPoint = fHCARecPoint ;
( (AliEMCALTrackSegment &)obj ).fIndexInList = fIndexInList ;
}
void AliEMCALTrackSegment::Print(Option_t *) const
{
// Print all information on this track Segment
-
-
- Info("Print", "TrackSegment information:") ;
+ printf("Print: TrackSegment information:") ;
printf("--------AliEMCALTrackSegment-------- \n");
printf("Stored at position %d\n", fIndexInList) ;
- if (fPRERecPoint)
- printf("PRE RecPoint # %d\n", fPRERecPoint) ;
if (fECARecPoint)
printf("EC RecPoint # %d\n", fECARecPoint) ;
- if (fHCARecPoint)
- printf("HC RecPoint # %d\n", fHCARecPoint) ;
-
- printf("------------------------------------ \n") ;
-
-}
-
-//____________________________________________________________________________
-void AliEMCALTrackSegment::SetPRERecPoint(AliEMCALRecPoint * pre)
-{
- // gives an id from its position in the list
- if( pre )
- fPRERecPoint = pre->GetIndexInList() ;
- else
- fPRERecPoint = -1 ;
-}
-
-//____________________________________________________________________________
-void AliEMCALTrackSegment::SetHCARecPoint(AliEMCALRecPoint * hca)
-{
- // gives an id from its position in the list
- if( hca )
- fHCARecPoint = hca->GetIndexInList() ;
- else
- fHCARecPoint = -1 ;
+ printf("------------------------------------ \n") ;
}
/* $Id$ */
//_________________________________________________________________________
-// Track segment in EMCAL
-// Can be any combination of : 1 PRERecPoint, ECRecPoint and HCRecPoint
+// Track segment in EMCAL
//
//*-- Author: Dmitri Peressounko (RRC KI & SUBATECH)
// Adapted from PHOS by Y. Schutz (SUBATECH)
public:
AliEMCALTrackSegment() {}
- AliEMCALTrackSegment(AliEMCALTowerRecPoint * ec, AliEMCALTowerRecPoint * pre, AliEMCALTowerRecPoint * hc) ;
+ AliEMCALTrackSegment(AliEMCALTowerRecPoint * ec) ;
AliEMCALTrackSegment(const AliEMCALTrackSegment & ts) ; // ctor
virtual ~AliEMCALTrackSegment() { }
void Copy(TObject & obj) ;
Int_t GetIndexInList() const { return fIndexInList ; }
- Int_t GetPREIndex() const { return fPRERecPoint ; }
Int_t GetECAIndex() const { return fECARecPoint; }
- Int_t GetHCAIndex() const { return fHCARecPoint; }
virtual void Print(Option_t * option) const;
void SetIndexInList(Int_t val){ fIndexInList = val ; }
- void SetPRERecPoint(AliEMCALRecPoint * pre ) ;
- void SetHCARecPoint(AliEMCALRecPoint * hc ) ;
typedef TClonesArray TrackSegmentsList ;
private:
-
- Int_t fPRERecPoint ; // The PRE reconstructed point index in array stored in TreeR/EMCALPRERP
Int_t fECARecPoint ; // The EC reconstructed point index in array stored in TreeR/EMCALECRP
- Int_t fHCARecPoint ; // The HC reconstructed point index in array stored in TreeR/EMCALHCRP
Int_t fIndexInList ; // The index of this TrackSegment in the list stored in TreeR (to be set by analysis)
ClassDef(AliEMCALTrackSegment,1) // Track segment in EMCAL
{
// dtor
// fDefaultInit = kTRUE if TrackSegmentMaker created by default ctor (to get just the parameters)
-
- if (!fDefaultInit) {
- delete fPRELinkArray ;
- delete fHCALinkArray ;
- }
+
}
//____________________________________________________________________________
toofar = kTRUE ;
if (gDebug == 2 )
- Info("HowClose", "ec = %d, rp = %d pro = %f, toofar=%d", ec->GetIndexInList(), rp->GetIndexInList(), pro, toofar ) ;
+ printf("HowClose: ec = %d, rp = %d pro = %f, toofar=%d", ec->GetIndexInList(), rp->GetIndexInList(), pro, toofar ) ;
return r ;
}
// Make all memory allocations that are not possible in default constructor
AliEMCALGetter* gime = AliEMCALGetter::Instance(GetTitle(), fEventFolderName.Data());
-
- fPRELinkArray = new TClonesArray("AliEMCALLink", 1000);
- fHCALinkArray = new TClonesArray("AliEMCALLink", 1000);
+
if ( !gime->TrackSegmentMaker() ) {
gime->PostTrackSegmentMaker(this);
}
void AliEMCALTrackSegmentMakerv1::InitParameters()
{
fClose = 10e-3 ;
- fPRELinkArray = 0 ;
- fHCALinkArray = 0 ;
fTrackSegmentsInRun = 0 ;
}
// which are not further apart from each other than fDangle
// and sort them in accordance with this distance
- AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
+ /* AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
TObjArray * aECARecPoints = gime->ECARecPoints() ;
- TObjArray * aPRERecPoints = gime->PRERecPoints() ;
- TObjArray * aHCARecPoints = gime->HCARecPoints() ;
+ // TObjArray * aPRERecPoints = gime->PRERecPoints() ;
+ //TObjArray * aHCARecPoints = gime->HCARecPoints() ;
fPRELinkArray->Clear() ;
fHCALinkArray->Clear() ;
fPRELinkArray->Sort() ; //first links with largest scalar product
fHCALinkArray->Sort() ; //first links with largest scalar product
+ */
}
//____________________________________________________________________________
// unassigned RecParticles. We assign these RecPoints to TrackSegment and
// remove them from the list of "unassigned".
- AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
- TObjArray * aECARecPoints = gime->ECARecPoints() ;
+ /*AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
+ TObjArray * aECARecPoints = gime->ECARecPoints() ;
TObjArray * aPRERecPoints = gime->PRERecPoints() ;
TObjArray * aHCARecPoints = gime->HCARecPoints() ;
TClonesArray * trackSegments = gime->TrackSegments() ;
(dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
fNTrackSegments++ ;
if (gDebug == 2 )
- Info("MakePairs", "ECAL section with PRE section") ;
+ printf("MakePairs: ECAL section with PRE section") ;
ecaExist[linkPRE->GetECA()] = -1 ; //Mark ecal that pre was found
//mark PRE recpoint as already used
preExist[linkPRE->GetOther()] = kFALSE ;
if (found){
ts->SetHCARecPoint( dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(linkHCA->GetOther())) ) ;
if (gDebug == 2 )
- Info("MakePairs", "ECAL section with PRE and HCAL sections") ;
+ printf("MakePairs: ECAL section with PRE and HCAL sections") ;
}
if (!found) {
new ((* trackSegments)[fNTrackSegments])
(dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
fNTrackSegments++ ;
if (gDebug == 2 )
- Info("MakePairs", "ECAL section with HCAL section") ;
+ printf("MakePairs: ECAL section with HCAL section") ;
}
ecaExist[linkHCA->GetECA()] = -2 ; //Mark ecal that hcal was found
//mark HCAL recpoint as already used
(dynamic_cast<AliEMCALTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
fNTrackSegments++;
if( gDebug == 2 )
- Info("MakePairs", "ECAL section alone") ;
+ printf("MakePairs: ECAL section alone") ;
}
}
}
delete [] ecaExist ;
delete [] preExist ;
delete [] hcaExist ;
+ */
}
//____________________________________________________________________________
if(strstr(option,"tim")){
gBenchmark->Stop("EMCALTSMaker");
- Info("Exec", "took %f seconds for making TS %f seconds per event",
+ printf("Exec: took %f seconds for making TS %f seconds per event",
gBenchmark->GetCpuTime("EMCALTSMaker"), gBenchmark->GetCpuTime("EMCALTSMaker")/nevents) ;
}
Unload();
//____________________________________________________________________________
void AliEMCALTrackSegmentMakerv1::Unload()
{
+ // Unloads the RecPoints and Tracks
AliEMCALGetter * gime = AliEMCALGetter::Instance() ;
gime->EmcalLoader()->UnloadRecPoints() ;
gime->EmcalLoader()->UnloadTracks() ;
{
// Print TrackSegmentMaker parameters
- Info("Print", "TrackSegmentMakerv1 parameters:") ;
+ printf("Print: TrackSegmentMakerv1 parameters:") ;
if( strcmp(GetName(), "") != 0 ) {
printf("Making Track segments with parameters:\n") ;
printf(" Allowed spred on the scalar product of two recpoints with same direction: %f\n", fClose) ;
TClonesArray * trackSegments = AliEMCALGetter::Instance()->TrackSegments() ;
- Info("PrintTrackSegments", "Results from TrackSegmentMaker:") ;
+ printf("PrintTrackSegments: Results from TrackSegmentMaker:") ;
printf("nevent: %d\n", gAlice->GetEvNumber()) ;
printf(" Found %d TrackSegments\n", trackSegments->GetEntriesFast() );
if(strstr(option,"all")) { // printing found TS
- printf("TrackSegment# ECAL RP# PRE RP# HCAL RP# \n") ;
+ printf("TrackSegment# ECAL RP# \n") ;
Int_t index;
for (index = 0 ; index < fNTrackSegments ; index++) {
AliEMCALTrackSegment * ts = (AliEMCALTrackSegment * )trackSegments->At(index) ;
- printf(" %d %d %d %d \n",
- ts->GetIndexInList(), ts->GetECAIndex(), ts->GetPREIndex(), ts->GetHCAIndex() );
+ printf(" %d %d \n",
+ ts->GetIndexInList(), ts->GetECAIndex());
}
}
}
virtual void Exec(Option_t * option) ;
Float_t HowClose(AliEMCALTowerRecPoint * ec, AliEMCALTowerRecPoint * rp, Bool_t &toofar) const ;
- void MakeLinks() const; //Evaluates distances(links) between PRE/EC/HC recpoints
+ void MakeLinks() const; //Evaluates distances(links) between recpoints
void MakePairs() ; //Finds pairs(triplets) with smallest link
virtual void Print(Option_t * option) const ;
virtual const char * Version() const { return "tsm-v1" ; }
Float_t fClose ; // Spread within which 2 recpoints are declared to have the same direction
Bool_t fDefaultInit ; //! Says if the task was created by defaut ctor (only parameters are initialized)
Int_t fNTrackSegments ; // number of track segments found
- TClonesArray * fPRELinkArray ;//! Contains the links ECAL-PRE
- TClonesArray * fHCALinkArray ;//! Contains the links ECAL-HCAL
Int_t fTrackSegmentsInRun ; //! Total number of track segments in one run
ClassDef( AliEMCALTrackSegmentMakerv1,3) // Implementation version 1 of algorithm class to make EMCAL track segments
// This Version of AliEMCALv0 reduces the number of volumes placed in XEN1 (the envelope) to less than five hundred
// The Envelope is Placed in Alice, And the Aluminium layer. Mini envelopes (XU) are then placed in XEN1.
-// Each mini envelope contains 2 scintillator, and 2 lead layers, except the last one which contains just one scintillator layer.
+// Each mini envelope contains 1 scintillator, and 1 lead layer, except the last one which contains just one scintillator layer.
// At the moment I cannot place the 36 and above layers in the mini envelopes so all layers are still placed in XEN1
TNode * top = gAlice->GetGeometry()->GetNode("alice") ;
new TTUBS("Envelop1", "Tubs that contains arm 1", "void",
geom->GetEnvelop(0), // rmin
- geom->GetEnvelop(1) +30 , // rmax
+ geom->GetEnvelop(1) +30 ,// rmax
geom->GetEnvelop(2)/2.0, // half length in Z
- geom->GetArm1PhiMin(), // minimun phi angle
- geom->GetArm1PhiMax() // maximun phi angle
+ geom->GetArm1PhiMin(), // minimum phi angle
+ geom->GetArm1PhiMax() // maximum phi angle
);
// Place the Node
//| | Air Gap = GetGap2Active() | |
//| | | |
//| ------------------------------------------------- |
- //| | XU0 : XPST (PreShower e = GetPRScintThick() )| |
- //| ------------------------------------------------- |
- //| | XU0 : XPBX (PreShower e = GetPRPbRadThick() )| |
- //| ------------------------------------------------- |
- //| | XU0 : XPST (PreShower e = GetPEScintThick() )| |
- //| ------------------------------------------------- |
- //| | XU0 : XPBX (PreShower e = GetPRPbRadThick() )| |
- //| ------------------------------------------------- |
- //| etc ..... GetNPRLayers() times |
- //| ------------------------------------------------- |
//| | XU1 : XPST (ECAL e = GetECScintThick() ) | |
//| ------------------------------------------------- |
//| | XU1 : XPBX (ECAL e = GetECPbRadThick() ) | |
//| ------------------------------------------------- |
//| | XU1 : XPBX (ECAL e = GetECPbRadThick() ) | |
//| ------------------------------------------------- |
- //| etc ..... GetNECLayers() times |
- //| ------------------------------------------------- |
- //| | XU1 : XPST (ECAL e = GetHCScintThick() ) | |
+ //| etc ..... GetNECLayers() - 1 times |
//| ------------------------------------------------- |
- //| | XU1 : XPBX (ECAL e = GetHCPbRadThick() ) | |
+ //| | XUNLayer : XPST (ECAL e = GetECScintThick() | |
//| ------------------------------------------------- |
- //| | XU1 : XPST (ECAL e = GetHCScintThick() | |
- //| ------------------------------------------------- |
- //| | XU1 : XPBX (ECAL e = GetHCPbRadThick() ) | |
- //| ------------------------------------------------- |
- //| etc ..... GetNHCLayers() times |
- //| ------------------------------------------------- |
- //| | XU10 : XPST (HCAL e = GetHCScintThick() ) | |
- //|-----------------------------------------------------|
-
+
Float_t etamin,etamax;
Float_t *dum=0;
Int_t idrotm = 1;
AliMatrix(idrotm, 90.0, 0., 90.0, 90.0, 0.0, 0.0) ;
-
-
// Create the EMCAL Mother Volume (a polygone) within which to place the Detector and named XEN1
Float_t envelopA[10];
envelopA[3] = 2; // 2 z coordinates
envelopA[4] = geom->ZFromEtaR(geom->GetEnvelop(1),
geom->GetArm1EtaMin()); // z coordinate 1
+ //add some padding for mother volume
envelopA[5] = geom->GetEnvelop(0) ; // rmin at z1
envelopA[6] = geom->GetEnvelop(1) ; // rmax at z1
envelopA[7] = geom->ZFromEtaR(geom->GetEnvelop(1),
gMC->Gspos("XEN1", 1, "ALIC", 0.0, 0.0, 0.0, idrotm, "ONLY") ;
- if (gDebug==2)
- Info("CreateGeometry","rXEN1 = %f, %f\n", envelopA[5], envelopA[6]);
-
- // Create mini-envelopes which will contain the PreShower scintillator-Lead-scintillator-lead
-
- if (gDebug==2)
- Info("CreateGeometry","XU0 = %f, %f\n", envelopA[5], envelopA[6]);
-
- // Create mini-envelopes which will contain the Tower scintillator-radiator-scintillator-radiator
+ if (gDebug==2) {
+ printf("CreateGeometry: XEN1 = %f, %f\n", envelopA[5], envelopA[6]);
+ printf("CreateGeometry: XU0 = %f, %f\n", envelopA[5], envelopA[6]);
+ }
+ // Create mini-envelopes which will contain the Tower scintillator-radiator
TString label ;
-
+
envelopA[5] = envelopA[5] + geom->GetGap2Active() // we are at the first scintllator
+ geom->GetAlFrontThickness(); // rmin at z1
envelopA[6] = envelopA[5] ;
Int_t i ;
- Int_t nLayers = geom->GetNPRLayers() + geom->GetNECLayers() + geom->GetNHCLayers() ;
+ Int_t nLayers = geom->GetNECLayers();
- for (i = 0; i < nLayers/2 ; i++ ){
+ for (i = 0; i < (nLayers-1); i++ ){
label = "XU" ;
label += i ;
Float_t tseg ;
- if (i == 0 )
- tseg = 2 *(geom->GetPRScintThick()+geom->GetPRPbRadThick()); // thickness of 2 * scintillator+Pb in pre shower
- else if ( i <= geom->GetNECLayers()/2)
- tseg = 2* (geom->GetECScintThick()+geom->GetECPbRadThick()); // thickness of 2 * scintillator+Pb in E Cal
- else
- tseg = 2* (geom->GetHCScintThick()+geom->GetHCCuRadThick()); // thickness of 2 * scintillator+Cu in H Cal
+ tseg = geom->GetECScintThick()+geom->GetECPbRadThick(); // thickness of scintillator+Pb in E Cal
envelopA[5] = envelopA[6] ; // rmin at z1
envelopA[4] = geom->ZFromEtaR(envelopA[5] + tseg,
geom->GetArm1EtaMin()); // z coordinate 1
gMC->Gspos(label.Data(), 1, "XEN1", 0.0, 0.0, 0.0, idrotm, "ONLY") ;
if (gDebug == 2)
- Info("CreateGeometry","XU%d = %f, %f\n", i, envelopA[5], envelopA[6]);
+ printf("CreateGeometry: XU%d = %f, %f\n", i, envelopA[5], envelopA[6]);
} // end i
- // Create one mini-envelope which will contain the last Tower scintillator (XU(nlayers-1)/2)
+ // Create one mini-envelope which will contain the last scintillator XU(nlayers-1) because there is one more scintillator than Pb layer XU(nlayers-1)
label = "XU" ;
label += i ;
envelopA[5] = envelopA[6] ; // rmin at z1
- envelopA[4] = geom->ZFromEtaR(envelopA[5],
+ envelopA[4] = geom->ZFromEtaR(envelopA[5] + geom->GetECScintThick(),
geom->GetArm1EtaMin()); // z coordinate 1
- envelopA[7] = geom->ZFromEtaR(envelopA[5],
+ envelopA[7] = geom->ZFromEtaR(envelopA[5] + geom->GetECScintThick(),
geom->GetArm1EtaMax()); // z coordinate 2
- if (geom->GetNHCLayers() == 0) // last scintillator is in ECAL
- envelopA[6] = envelopA[5] + geom->GetECScintThick() ; // rmax at z1
- else // last scintillator is in HCAL
- envelopA[6] = envelopA[5] + geom->GetECScintThick() ; // rmax at z1
+ envelopA[6] = envelopA[5] + geom->GetECScintThick() ; // rmax at z1
envelopA[8] = envelopA[5] ; // radii are the same.
envelopA[9] = envelopA[6] ; // radii are the same.
// Position the last minienvelope in XEN1
gMC->Gspos(label.Data(), 1, "XEN1", 0.0, 0.0, 0.0, idrotm, "ONLY") ;
-
+
if(gDebug == 2)
- Info("CreateGeometry","XEN%d = %f, %f\n", i, envelopA[5], envelopA[6]);
+ printf("CreateGeometry: XEN%d = %f, %f\n", i, envelopA[5], envelopA[6]);
// Create the shapes of active material (LEAD/Aluminium/Scintillator)
// to be placed
envelopC[0] = envelopD[0] = envelopB[0] = envelopA[0] ; // starting position in Phi
envelopC[1] = envelopD[1] = envelopB[1] = envelopA[1] ; // angular range in phi
envelopC[2] = envelopD[2] = envelopB[2] = envelopA[2] ; // number of sections in Phi
- envelopD[3] = envelopC[3] = envelopB[3] = envelopA[3] ; // 2 z coordinates
+ envelopC[3] = envelopD[3] = envelopB[3] = envelopA[3] ; // 2 z coordinates
Float_t dist = geom->GetEnvelop(0) + geom->GetAlFrontThickness() + geom->GetGap2Active() ;
envelopB[4] = geom->ZFromEtaR(dist,
gMC->Gsvolu("XPBX", "PGON", idtmed[1600], dum, 0); // PGON filled with Lead (shape to be defined by GSPOSP)
- gMC->Gsvolu("XCUX", "PGON", idtmed[1603], dum, 0); // PGON filled with Copper (shape to be defined by GSPOSP)
+ //gMC->Gsvolu("XCUX", "PGON", idtmed[1603], dum, 0); // PGON filled with Copper (shape to be defined by GSPOSP)
gMC->Gsdvn("XPHI", "XPST", geom->GetNPhi(), 2); // Divide eta section of scintillators into phi segments.
for (int i = 0; i < nLayers; i++ ){
label = "XU" ;
- label += static_cast<Int_t> (i/2) ; // we will place two layers (i = one layer) in each mini envelope)
+ label += i ; // we will place one layer in each mini envelope)
Float_t scthick ; // scintillator thickness
- if ( i < geom->GetNPRLayers() ) // its a preshower
- scthick = geom->GetPRScintThick() ;
- else if( i < geom->GetNPRLayers() + geom->GetNECLayers() ) // its an EMCAL section
- scthick = geom->GetECScintThick() ;
- else // its an HCAL section
- scthick = geom->GetHCScintThick() ;
+ scthick = geom->GetECScintThick() ;
envelopC[5] = envelopD[6] ; //rmin
envelopC[6] = envelopC[5] + scthick ; //rmax
envelopC[8] = envelopC[5] ; //rmin
envelopC[9] = envelopC[6] ; //rmax
-
- // envelopC[6] = envelopD[6] + ((i > 1) ? geom->GetFullSintThick() : geom->GetPreSintThick());//rmax larger for first two layers (preshower)
- // envelopC[9] = envelopD[6] + ((i > 1 ) ? geom->GetFullSintThick() :geom->GetPreSintThick());//rmax larger for first two layers (preshower)
-
if(gDebug == 2 )
- Info("CreateGeometry", "volume = %s, name = XPST thickness = %f deb = %f/%f fin = %f/%f", label.Data(), scthick, envelopC[5], envelopC[8], envelopC[6], envelopC[9]) ;
+ printf("CreateGeometry: volume = %s, name = XPST thickness = %f deb = %f/%f fin = %f/%f", label.Data(), scthick, envelopC[5], envelopC[8], envelopC[6], envelopC[9]) ;
for (int j =0; j < (geom->GetNEta()) ; j++){
etamin = geom->GetArm1EtaMin()+
0.0, 0.0, 0.0 , idrotm, "ONLY", envelopC, 10); // Position and define layer
} // end for j
- if (i < nLayers){
Float_t radthick ; // radiator thickness
TString radname ; // radiator name
- if ( i <= 1 ) { // its a preshower
- radthick = geom->GetPRPbRadThick();
- radname = "XPBX" ;
- }
- else if( i <= geom->GetNECLayers()) {// its an EMCAL section
- radthick = geom->GetECPbRadThick();
- radname = "XPBX" ;
- }
- else { // its an HCAL section
- radthick = geom->GetHCCuRadThick();
- radname = "XCUX" ;
- }
+ radthick = geom->GetECPbRadThick();
+ radname = "XPBX" ;
if ( i < nLayers -1 ) { // except for the last XU which contains only one scintillator layer
envelopD[5] = envelopC[6] ; //rmin
envelopD[8] = envelopD[5] ; //rmin
envelopD[6] = envelopD[5] + radthick ; // rmax
- // envelopD[6] = envelopC[6] + geom->GetPbRadThick(); //rmax
envelopD[9] = envelopD[6] ; //rmax
- // envelopD[9] = envelopC[6] + geom->GetPbRadThick(); //rmax
if(gDebug == 2 )
- Info("CreateGeometry", "volume = %s, name = %s thickness = %f deb = %f/%f fin = %f/%f", label.Data(), radname.Data(), radthick, envelopD[5], envelopD[8], envelopD[6], envelopD[9]) ;
+ printf("CreateGeometry: volume = %s, name = %s thickness = %f deb = %f/%f fin = %f/%f", label.Data(), radname.Data(), radthick, envelopD[5], envelopD[8], envelopD[6], envelopD[9]) ;
for (int j =0; j < (geom->GetNEta()) ; j++){
etamin = geom->GetArm1EtaMin()+
gMC->Gsposp(radname.Data(),1+j+i*(geom->GetNEta()), label.Data(),
0.0, 0.0, 0.0 , idrotm, "ONLY", envelopD, 10);
- } // end for j
- } // if not last layer
- } // end if i
+ } // end for j
+ } // if not last layer
} // for i
}
message += "EMCAL geometry initialization failed !" ;
}
message += "\n*****************************************" ;
- Info("Init", message.Data() ) ;
+ printf(message.Data() ) ;
}
}
/* $Id$ */
//_________________________________________________________________________
-// Implementation version v1 of EMCAL Manager class
-// An object of this class does not produce digits
-// It is the one to use if you do want to produce outputs in TREEH
-//
+//*-- Implementation version v1 of EMCAL Manager class
+//*-- An object of this class does not produce digits
+//*-- It is the one to use if you do want to produce outputs in TREEH
+//*--
//*-- Author: Sahal Yacoob (LBL /UCT)
//*-- : Jennifer Klay (LBL)
-
// This Class not stores information on all particles prior to EMCAL entry - in order to facilitate analysis.
// This is done by setting fIShunt =2, and flagging all parents of particles entering the EMCAL.
void AliEMCALv1::AddHit(Int_t shunt, Int_t primary, Int_t tracknumber, Int_t iparent, Float_t ienergy,
Int_t id, Float_t * hits,Float_t * p){
// Add a hit to the hit list.
- // An EMCAL hit is the sum of all hits in a tower section (PRE, ECAL, HCAL)
+ // An EMCAL hit is the sum of all hits in a tower section
// originating from the same entering particle
Int_t hitCounter;
Bool_t deja = kFALSE;
newHit = new AliEMCALHit(shunt, primary, tracknumber, iparent, ienergy, id, hits, p);
- for ( hitCounter = fNhits-1; hitCounter >= 0 && !deja; hitCounter-- ) {
+ for ( hitCounter = fNhits-1; hitCounter >= 0 && !deja; hitCounter-- ) {
curHit = (AliEMCALHit*) (*fHits)[hitCounter];
// We add hits with the same tracknumber, while GEANT treats
// primaries succesively
if( *curHit == *newHit ) {
*curHit = *curHit + *newHit;
deja = kTRUE;
- } // end if
+ } // end if
} // end for hitCounter
-
+
if ( !deja ) {
new((*fHits)[fNhits]) AliEMCALHit(*newHit);
fNhits++;
void AliEMCALv1::StepManager(void){
// Accumulates hits as long as the track stays in a tower
- Int_t id[2]; // (layer, phi, Eta) indices
+ Int_t id[2]; // (phi, Eta) indices
// position wrt MRS and energy deposited
Float_t xyzte[5]={0.,0.,0.,0.,0.};// position wrt MRS, time and energy deposited
Float_t pmom[4]={0.,0.,0.,0.};
Int_t tower = (id[0]-1) % geom->GetNZ() + 1 + (id[1] - 1) * geom->GetNZ() ;
Int_t layer = static_cast<Int_t>((id[0]-1)/(geom->GetNZ())) + 1 ;
- Int_t absid = tower ;
- if (layer <= geom->GetNPRLayers() )
- absid += geom->GetNZ() * geom->GetNPhi() ;
- else if (layer > geom->GetNECLayers() + geom->GetNPRLayers() )
- absid += 2 * geom->GetNZ() * geom->GetNPhi() ;
- else {
- Int_t nlayers = geom->GetNPRLayers()+ geom->GetNECLayers()+ geom->GetNHCLayers() ;
- if (layer > nlayers)
- Fatal("StepManager", "Wrong calculation of layer number: layer = %d > %d\n", layer, nlayers) ;
- }
+ Int_t absid = tower;
+ Int_t nlayers = geom->GetNECLayers();
+ if ((layer > nlayers)||(layer<1))
+ Fatal("StepManager", "Wrong calculation of layer number: layer = %d > %d\n", layer, nlayers) ;
+
Float_t lightYield = depositedEnergy ;
- ;
xyzte[4] = lightYield ;
primary = gAlice->GetMCApp()->GetPrimary(tracknumber);
if (gDebug == 2)
- Info("StepManager", "id0 = %d, id1 = %d, absid = %d tower = %d layer = %d energy = %f\n", id[0], id[1], absid, tower, layer, xyzte[4]) ;
+ printf("StepManager: id0 = %d, id1 = %d, absid = %d tower = %d layer = %d energy = %f\n", id[0], id[1], absid, tower, layer, xyzte[4]) ;
- AddHit(fIshunt, primary,tracknumber, iparent, ienergy, absid, xyzte, pmom);
+ AddHit(fIshunt, primary,tracknumber, iparent, ienergy, absid, xyzte, pmom);
} // there is deposited energy
}
}
// gener->Init();
AliGenBox *gener = new AliGenBox(1);
- gener->SetMomentumRange(10,11.);
+ gener->SetMomentumRange(50.,100.);
gener->SetPhiRange(60.0,180.0);
gener->SetThetaRange(EtaToTheta(-0.7), EtaToTheta(0.7));
gener->SetOrigin(0,0,0); //vertex position
gener->SetSigma(0,0,0); //Sigma in (X,Y,Z) (cm) on IP position
- gener->SetPart(22);
+ gener->SetPart(kGamma);
gener->Init();
//