/* $Id$ */
+/* History of cvs commits:
+ *
+ * $Log$
+ * Revision 1.88 2006/01/11 08:54:52 hristov
+ * Additional protection in case no calibration entry was found
+ *
+ * Revision 1.87 2005/11/22 08:46:43 kharlov
+ * Updated to new CDB (Boris Polichtchouk)
+ *
+ * Revision 1.86 2005/11/14 21:52:43 hristov
+ * Coding conventions
+ *
+ * Revision 1.85 2005/09/27 16:08:08 hristov
+ * New version of CDB storage framework (A.Colla)
+ *
+ * Revision 1.84 2005/09/21 10:02:47 kharlov
+ * Reading calibration from CDB (Boris Polichtchouk)
+ *
+ * Revision 1.82 2005/09/02 15:43:13 kharlov
+ * Add comments in GetCalibrationParameters and Calibrate
+ *
+ * Revision 1.81 2005/09/02 14:32:07 kharlov
+ * Calibration of raw data
+ *
+ * Revision 1.80 2005/08/24 15:31:36 kharlov
+ * Setting raw digits flag
+ *
+ * Revision 1.79 2005/07/25 15:53:53 kharlov
+ * Read raw data
+ *
+ * Revision 1.78 2005/05/28 14:19:04 schutz
+ * Compilation warnings fixed by T.P.
+ *
+ */
+
//*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (SUBATECH & Kurchatov Institute)
//////////////////////////////////////////////////////////////////////////////
// Clusterization class. Performs clusterization (collects neighbouring active cells) and
// --- Standard library ---
// --- AliRoot header files ---
+#include "AliLog.h"
#include "AliPHOSGetter.h"
#include "AliPHOSGeometry.h"
#include "AliPHOSClusterizerv1.h"
#include "AliPHOSCpvRecPoint.h"
#include "AliPHOSDigit.h"
#include "AliPHOSDigitizer.h"
+#include "AliPHOSCalibrationDB.h"
+#include "AliCDBManager.h"
+#include "AliCDBStorage.h"
+#include "AliCDBEntry.h"
ClassImp(AliPHOSClusterizerv1)
{
// dtor
- // delete fPedestals ;
-// delete fGains ;
}
-
-
//____________________________________________________________________________
const TString AliPHOSClusterizerv1::BranchName() const
{
}
//____________________________________________________________________________
-Float_t AliPHOSClusterizerv1::Calibrate(Int_t amp, Int_t absId) const
+Float_t AliPHOSClusterizerv1::Calibrate(Int_t amp, Int_t absId)
{
- //To be replaced later by the method, reading individual parameters from the database
- // if(fCalibrationDB)
- // return fCalibrationDB->Calibrate(amp,absId) ;
- // else{ //simulation
- if(absId <= fEmcCrystals) //calibrate as EMC
- return fADCpedestalEmc + amp*fADCchanelEmc ;
- else //calibrate as CPV
- return fADCpedestalCpv+ amp*fADCchanelCpv ;
- // }
+ // Convert digitized amplitude into energy.
+ // Calibration parameters are taken from calibration data base for raw data,
+ // or from digitizer parameters for simulated data.
+
+ if(fCalibData){
+ Int_t relId[4];
+ AliPHOSGetter *gime = AliPHOSGetter::Instance();
+ gime->PHOSGeometry()->AbsToRelNumbering(absId,relId) ;
+ Int_t module = relId[0];
+ Int_t column = relId[3];
+ Int_t row = relId[2];
+ if(absId <= fEmcCrystals) { //calibrate as EMC
+ fADCchanelEmc = fCalibData->GetADCchannelEmc (module,column,row);
+ fADCpedestalEmc = fCalibData->GetADCpedestalEmc(module,column,row);
+ return fADCpedestalEmc + amp*fADCchanelEmc ;
+ }
+ else { //calibrate as CPV
+ fADCchanelCpv = fCalibData->GetADCchannelCpv (module,column,row);
+ fADCpedestalCpv = fCalibData->GetADCpedestalCpv(module,column,row);
+ return fADCpedestalCpv + amp*fADCchanelCpv ;
+ }
+ }
+ else{ //simulation
+ if(absId <= fEmcCrystals) //calibrate as EMC
+ return fADCpedestalEmc + amp*fADCchanelEmc ;
+ else //calibrate as CPV
+ return fADCpedestalCpv+ amp*fADCchanelCpv ;
+ }
}
//____________________________________________________________________________
return ;
}
+ GetCalibrationParameters() ;
+
AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
+ if (fRawReader == 0)
+ gime->SetRawDigits(kFALSE);
+ else
+ gime->SetRawDigits(kTRUE);
- if (fLastEvent == -1) fLastEvent = gime->MaxEvent() - 1 ;
- else fLastEvent = TMath::Min(fLastEvent,gime->MaxEvent());
+ if (fLastEvent == -1)
+ fLastEvent = gime->MaxEvent() - 1 ;
+ else
+ fLastEvent = TMath::Min(fFirstEvent, gime->MaxEvent()); // one event at the time
Int_t nEvents = fLastEvent - fFirstEvent + 1;
-
- for (Int_t ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
- gime->Event(ievent, "D");
-
- GetCalibrationParameters() ;
+ Int_t ievent ;
+ for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
+ if (fRawReader == 0)
+ gime->Event(ievent ,"D"); // Read digits from simulated data
+ else
+ gime->Event(fRawReader,"W"); // Read digits from raw data
+
fNumberOfEmcClusters = fNumberOfCpvClusters = 0 ;
MakeClusters() ;
fRecPointsInRun += gime->CpvRecPoints()->GetEntriesFast() ;
}
- Unload();
+ if(fWrite) //do not unload in "on flight" mode
+ Unload();
if(strstr(option,"tim")){
gBenchmark->Stop("PHOSClusterizer");
- Info("Exec", " took %f seconds for Clusterizing %f seconds per event \n",
+ AliInfo(Form("took %f seconds for Clusterizing %f seconds per event \n",
gBenchmark->GetCpuTime("PHOSClusterizer"),
- gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents ) ;
+ gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents )) ;
}
}
//____________________________________________________________________________
void AliPHOSClusterizerv1::GetCalibrationParameters()
{
+ // Set calibration parameters:
+ // if calibration database exists, they are read from database,
+ // otherwise, they are taken from digitizer.
+ //
+ // It is a user responsilibity to open CDB before reconstruction, for example:
+ // AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage("local://CalibDB");
AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- if ( !gime->Digitizer() )
- gime->LoadDigitizer();
- AliPHOSDigitizer * dig = gime->Digitizer();
-
- fADCchanelEmc = dig->GetEMCchannel() ;
- fADCpedestalEmc = dig->GetEMCpedestal();
+ fCalibData = new AliPHOSCalibData(gAlice->GetRunNumber());
- fADCchanelCpv = dig->GetCPVchannel() ;
- fADCpedestalCpv = dig->GetCPVpedestal() ;
-
- // fCalibrationDB = gime->CalibrationDB();
+ if(!fCalibData)
+ {
+ if ( !gime->Digitizer() )
+ gime->LoadDigitizer();
+ AliPHOSDigitizer * dig = gime->Digitizer();
+ fADCchanelEmc = dig->GetEMCchannel() ;
+ fADCpedestalEmc = dig->GetEMCpedestal();
+
+ fADCchanelCpv = dig->GetCPVchannel() ;
+ fADCpedestalCpv = dig->GetCPVpedestal() ;
+ }
}
//____________________________________________________________________________
{
// Make all memory allocations which can not be done in default constructor.
// Attach the Clusterizer task to the list of PHOS tasks
-
- AliPHOSGetter* gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
+
+ AliPHOSGetter* gime = AliPHOSGetter::Instance() ;
+ if(!gime)
+ gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
AliPHOSGeometry * geom = gime->PHOSGeometry();
-
+
fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
if(!gMinuit)
gMinuit = new TMinuit(100);
- if ( !gime->Clusterizer() )
+ if ( !gime->Clusterizer() ) {
gime->PostClusterizer(this);
+ }
}
//____________________________________________________________________________
fEmcLocMaxCut = 0.03 ;
fCpvLocMaxCut = 0.03 ;
-
+
+ fEmcMinE = 0.01 ;
+ fCpvMinE = 0.0 ;
+
fW0 = 4.5 ;
fW0CPV = 4.0 ;
fRecPointsInRun = 0 ;
+ fWrite = kTRUE ;
+
+ fCalibData = 0 ;
+
SetEventRange(0,-1) ;
}
return rv ;
}
+//____________________________________________________________________________
+void AliPHOSClusterizerv1::CleanDigits(TClonesArray * digits){
+ for(Int_t i=0; i<digits->GetEntriesFast(); i++){
+ AliPHOSDigit * digit = static_cast<AliPHOSDigit*>(digits->At(i)) ;
+ Float_t cut = IsInEmc(digit) ? fEmcMinE : fCpvMinE ;
+ if(Calibrate(digit->GetAmp(),digit->GetId()) < cut)
+ digits->RemoveAt(i) ;
+ }
+ digits->Compress() ;
+ for (Int_t i = 0 ; i < digits->GetEntriesFast() ; i++) {
+ AliPHOSDigit *digit = static_cast<AliPHOSDigit*>( digits->At(i) ) ;
+ digit->SetIndexInList(i) ;
+ }
-
+}
//____________________________________________________________________________
Bool_t AliPHOSClusterizerv1::IsInEmc(AliPHOSDigit * digit) const
{
TObjArray * emcRecPoints = gime->EmcRecPoints() ;
TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
TClonesArray * digits = gime->Digits() ;
-
- TTree * treeR = gime->TreeR();
-
+
+
Int_t index ;
- //Evaluate position, dispersion and other RecPoint properties...
- for(index = 0; index < emcRecPoints->GetEntries(); index++)
- dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->EvalAll(fW0,digits) ;
-
+ //Evaluate position, dispersion and other RecPoint properties..
+ Int_t nEmc = emcRecPoints->GetEntriesFast();
+ for(index = 0; index < nEmc; index++){
+ AliPHOSEmcRecPoint * rp = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) );
+ rp->Purify(fEmcMinE) ;
+ if(rp->GetMultiplicity()>0.) //If this RP is not empty
+ rp->EvalAll(fW0,digits) ;
+ else{
+ emcRecPoints->RemoveAt(index) ;
+ delete rp ;
+ }
+ }
+ emcRecPoints->Compress() ;
emcRecPoints->Sort() ;
- for(index = 0; index < emcRecPoints->GetEntries(); index++)
+ // emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ;
+ for(index = 0; index < emcRecPoints->GetEntries(); index++){
dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->SetIndexInList(index) ;
-
- emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ;
+ }
//Now the same for CPV
- for(index = 0; index < cpvRecPoints->GetEntries(); index++)
- dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) )->EvalAll(digits) ;
-
+ for(index = 0; index < cpvRecPoints->GetEntries(); index++){
+ AliPHOSCpvRecPoint * rp = dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) );
+ rp->EvalAll(fW0CPV,digits) ;
+ }
cpvRecPoints->Sort() ;
for(index = 0; index < cpvRecPoints->GetEntries(); index++)
- dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
+ dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ;
- Int_t bufferSize = 32000 ;
- Int_t splitlevel = 0 ;
-
- //First EMC
- TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
- emcBranch->SetTitle(BranchName());
+ if(fWrite){ //We write TreeR
+ TTree * treeR = gime->TreeR();
- //Now CPV branch
- TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
- cpvBranch->SetTitle(BranchName());
+ Int_t bufferSize = 32000 ;
+ Int_t splitlevel = 0 ;
- emcBranch ->Fill() ;
- cpvBranch ->Fill() ;
-
- gime->WriteRecPoints("OVERWRITE");
- gime->WriteClusterizer("OVERWRITE");
+ //First EMC
+ TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
+ emcBranch->SetTitle(BranchName());
+
+ //Now CPV branch
+ TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
+ cpvBranch->SetTitle(BranchName());
+
+ emcBranch ->Fill() ;
+ cpvBranch ->Fill() ;
+
+ gime->WriteRecPoints("OVERWRITE");
+ gime->WriteClusterizer("OVERWRITE");
+ }
}
//____________________________________________________________________________
TClonesArray * digits = gime->Digits() ;
+ //Remove digits below threshold
+ CleanDigits(digits) ;
+
+
TClonesArray * digitsC = static_cast<TClonesArray*>( digits->Clone() ) ;
fNumberOfEmcClusters -- ;
numberofNotUnfolded-- ;
}
+ else{
+ emcRecPoint->SetNExMax(1) ; //Only one local maximum
+ }
delete[] maxAt ;
delete[] maxAtEnergy ;
Bool_t rv = FindFit(iniEmc, maxAt, maxAtEnergy, nPar, fitparameters) ;
if( !rv ) {
// Fit failed, return and remove cluster
+ iniEmc->SetNExMax(-1) ;
delete[] fitparameters ;
return ;
}
(*emcRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ;
emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(fNumberOfEmcClusters) ) ;
fNumberOfEmcClusters++ ;
+ emcRP->SetNExMax((Int_t)nPar/3) ;
}
else{//create new entries in fCpvRecPoints
if(fNumberOfCpvClusters >= cpvRecPoints->GetSize())
}
//____________________________________________________________________________
-void AliPHOSClusterizerv1::Print()const
+void AliPHOSClusterizerv1::Print(const Option_t *)const
{
// Print clusterizer parameters
else
message = " AliPHOSClusterizerv1 not initialized " ;
- Info("Print", message.Data(),
+ AliInfo(Form("%s, %s %s %s %s %s %s %s %s %s %s", message.Data(),
taskName.Data(),
GetTitle(),
taskName.Data(),
fW0,
fCpvClusteringThreshold,
fCpvLocMaxCut,
- fW0CPV ) ;
+ fW0CPV )) ;
}
TObjArray * emcRecPoints = gime->EmcRecPoints() ;
TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
- printf("\nevent %d \n", gAlice->GetEvNumber()) ;
- printf(" Found %d EMC RecPoints and %d CPV RecPoints \n",
- emcRecPoints->GetEntriesFast(),cpvRecPoints->GetEntriesFast()) ;
+ AliInfo(Form("\nevent %d \n Found %d EMC RecPoints and %d CPV RecPoints",
+ gAlice->GetEvNumber(),
+ emcRecPoints->GetEntriesFast(),
+ cpvRecPoints->GetEntriesFast() )) ;
fRecPointsInRun += emcRecPoints->GetEntriesFast() ;
fRecPointsInRun += cpvRecPoints->GetEntriesFast() ;
printf("\n CPV clusters \n") ;
printf("Index Ene(MeV) Module X Y Z \n") ;
for (index = 0 ; index < cpvRecPoints->GetEntries() ; index++) {
- AliPHOSRecPoint * rp = (AliPHOSRecPoint * )cpvRecPoints->At(index) ;
-
- TVector3 locpos;
- rp->GetLocalPosition(locpos);
-
- printf("\n%6d %8.2f %2d %4.1f %4.1f %4.1f \n",
- rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(),
- locpos.X(), locpos.Y(), locpos.Z()) ;
+ AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )cpvRecPoints->At(index) ;
+
+ TVector3 locpos;
+ rp->GetLocalPosition(locpos);
+
+ printf("\n%6d %8.2f %2d %4.1f %4.1f %4.1f \n",
+ rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(),
+ locpos.X(), locpos.Y(), locpos.Z()) ;
}
}
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff