/* History of cvs commits:
*
- * $Log$
+ * $Log: AliPHOSClusterizerv1.cxx,v $
+ * Revision 1.118 2007/12/11 21:23:26 kharlov
+ * Added possibility to swith off unfolding
+ *
+ * Revision 1.117 2007/10/18 08:42:05 kharlov
+ * Bad channels cleaned before clusterization
+ *
+ * Revision 1.116 2007/10/01 20:24:08 kharlov
+ * Memory leaks fixed
+ *
+ * Revision 1.115 2007/09/26 14:22:17 cvetan
+ * Important changes to the reconstructor classes. Complete elimination of the run-loaders, which are now steered only from AliReconstruction. Removal of the corresponding Reconstruct() and FillESD() methods.
+ *
+ * Revision 1.114 2007/09/06 16:06:44 kharlov
+ * Absence of sorting results in loose of all unfolded clusters
+ *
+ * Revision 1.113 2007/08/28 12:55:07 policheh
+ * Loaders removed from the reconstruction code (C.Cheshkov)
+ *
+ * Revision 1.112 2007/08/22 09:20:50 hristov
+ * Updated QA classes (Yves)
+ *
+ * Revision 1.111 2007/08/08 12:11:28 kharlov
+ * Protection against uninitialized fQADM
+ *
+ * Revision 1.110 2007/08/07 14:16:00 kharlov
+ * Quality assurance added (Yves Schutz)
+ *
+ * Revision 1.109 2007/07/24 17:20:35 policheh
+ * Usage of RecoParam objects instead of hardcoded parameters in reconstruction.
+ * (See $ALICE_ROOT/PHOS/macros/BeamTest2006/RawReconstruction.C).
+ *
+ * Revision 1.108 2007/06/18 07:00:51 kharlov
+ * Bug fix for attempt to use AliPHOSEmcRecPoint after its deletion
+ *
+ * Revision 1.107 2007/05/25 14:12:26 policheh
+ * Local to tracking CS transformation added for CPV rec. points
+ *
+ * Revision 1.106 2007/05/24 13:01:22 policheh
+ * Local to tracking CS transformation invoked for each EMC rec.point
+ *
+ * Revision 1.105 2007/05/02 13:41:22 kharlov
+ * Mode protection against absence of calib.data from AliPHOSCalibData to AliPHOSClusterizerv1::GetCalibrationParameters()
+ *
+ * Revision 1.104 2007/04/27 16:55:53 kharlov
+ * Calibration stops if PHOS CDB objects do not exist
+ *
+ * Revision 1.103 2007/04/11 11:55:45 policheh
+ * SetDistancesToBadChannels() added.
+ *
+ * Revision 1.102 2007/03/28 19:18:15 kharlov
+ * RecPoints recalculation in TSM removed
+ *
+ * Revision 1.101 2007/03/06 06:51:27 kharlov
+ * Calculation of cluster properties dep. on vertex posponed to TrackSegmentMaker
+ *
+ * Revision 1.100 2007/01/10 11:07:26 kharlov
+ * Raw digits writing to file (B.Polichtchouk)
+ *
+ * Revision 1.99 2006/11/07 16:49:51 kharlov
+ * Corrections for next event switch in case of raw data (B.Polichtchouk)
+ *
+ * Revision 1.98 2006/10/27 17:14:27 kharlov
+ * Introduce AliDebug and AliLog (B.Polichtchouk)
+ *
+ * Revision 1.97 2006/08/29 11:41:19 kharlov
+ * Missing implementation of ctors and = operator are added
+ *
+ * Revision 1.96 2006/08/25 16:56:30 kharlov
+ * Compliance with Effective C++
+ *
* Revision 1.95 2006/08/11 12:36:26 cvetan
* Update of the PHOS code needed in order to read and reconstruct the beam test raw data (i.e. without an existing galice.root)
*
// This TTask is normally called from Reconstructioner, but can as well be used in
// standalone mode.
// Use Case:
-// root [0] AliPHOSClusterizerv1 * cl = new AliPHOSClusterizerv1("galice.root", "recpointsname", "digitsname")
-// Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
-// // reads gAlice from header file "galice.root", uses digits stored in the branch names "digitsname" (default = "Default")
-// // and saves recpoints in branch named "recpointsname" (default = "digitsname")
-// root [1] cl->ExecuteTask()
-// //finds RecPoints in all events stored in galice.root
+// root [0] AliPHOSClusterizerv1 * cl = new AliPHOSClusterizerv1(<pointer_to_phos_geometry_onject>)
+// root [1] cl->Digits2Clusters(digitsTree,clusterTree)
+// //finds RecPoints in the current event
// root [2] cl->SetDigitsBranch("digits2")
// //sets another title for Digitis (input) branch
// root [3] cl->SetRecPointsBranch("recp2")
// //sets another title four output branches
// root [4] cl->SetEmcLocalMaxCut(0.03)
// //set clusterization parameters
-// root [5] cl->ExecuteTask("deb all time")
-// //once more finds RecPoints options are
-// // deb - print number of found rec points
-// // deb all - print number of found RecPoints and some their characteristics
-// // time - print benchmarking results
// --- ROOT system ---
#include "TMinuit.h"
#include "TTree.h"
#include "TBenchmark.h"
+#include "TClonesArray.h"
// --- Standard library ---
// --- AliRoot header files ---
#include "AliLog.h"
-#include "AliPHOSGetter.h"
+#include "AliConfig.h"
#include "AliPHOSGeometry.h"
#include "AliPHOSClusterizerv1.h"
#include "AliPHOSEmcRecPoint.h"
#include "AliCDBManager.h"
#include "AliCDBStorage.h"
#include "AliCDBEntry.h"
+#include "AliPHOSRecoParam.h"
+#include "AliPHOSReconstructor.h"
+#include "AliPHOSCalibData.h"
ClassImp(AliPHOSClusterizerv1)
AliPHOSClusterizer(),
fDefaultInit(0), fEmcCrystals(0), fToUnfold(0),
fWrite(0), fNumberOfEmcClusters(0), fNumberOfCpvClusters(0),
- fCalibData(0), fADCchanelEmc(0), fADCpedestalEmc(0),
- fADCchanelCpv(0), fADCpedestalCpv(0), fEmcClusteringThreshold(0),
- fCpvClusteringThreshold(0), fEmcMinE(0), fCpvMinE(0),
+ fEmcClusteringThreshold(0), fCpvClusteringThreshold(0),
fEmcLocMaxCut(0), fW0(0), fCpvLocMaxCut(0),
- fW0CPV(0), fRecPointsInRun(0), fEmcTimeGate(0),
- fIsOldRCUFormat(0)
+ fW0CPV(0), fEmcTimeGate(0)
{
// default ctor (to be used mainly by Streamer)
}
//____________________________________________________________________________
-AliPHOSClusterizerv1::AliPHOSClusterizerv1(const TString alirunFileName, const TString eventFolderName) :
- AliPHOSClusterizer(alirunFileName, eventFolderName),
+AliPHOSClusterizerv1::AliPHOSClusterizerv1(AliPHOSGeometry *geom) :
+ AliPHOSClusterizer(geom),
fDefaultInit(0), fEmcCrystals(0), fToUnfold(0),
fWrite(0), fNumberOfEmcClusters(0), fNumberOfCpvClusters(0),
- fCalibData(0), fADCchanelEmc(0), fADCpedestalEmc(0),
- fADCchanelCpv(0), fADCpedestalCpv(0), fEmcClusteringThreshold(0),
- fCpvClusteringThreshold(0), fEmcMinE(0), fCpvMinE(0),
+ fEmcClusteringThreshold(0), fCpvClusteringThreshold(0),
fEmcLocMaxCut(0), fW0(0), fCpvLocMaxCut(0),
- fW0CPV(0), fRecPointsInRun(0), fEmcTimeGate(0),
- fIsOldRCUFormat(0)
+ fW0CPV(0), fEmcTimeGate(0)
{
// ctor with the indication of the file where header Tree and digits Tree are stored
}
//____________________________________________________________________________
-const TString AliPHOSClusterizerv1::BranchName() const
-{
- return GetName();
-}
-
-//____________________________________________________________________________
-Float_t AliPHOSClusterizerv1::CalibrateEMC(Float_t amp, Int_t absId)
-{
- // Convert EMC measured amplitude into real 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) { // this is EMC
- fADCchanelEmc = fCalibData->GetADCchannelEmc (module,column,row);
- return amp*fADCchanelEmc ;
- }
- }
- else{ //simulation
- if(absId <= fEmcCrystals) // this is EMC
- return fADCpedestalEmc + amp*fADCchanelEmc ;
- }
- return 0;
-}
-
-//____________________________________________________________________________
-Float_t AliPHOSClusterizerv1::CalibrateCPV(Int_t amp, Int_t absId)
-{
- // Convert digitized CPV amplitude into charge.
- // 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) { // this is CPV
- fADCchanelCpv = fCalibData->GetADCchannelCpv (module,column,row);
- fADCpedestalCpv = fCalibData->GetADCpedestalCpv(module,column,row);
- return fADCpedestalCpv + amp*fADCchanelCpv ;
- }
- }
- else{ //simulation
- if(absId > fEmcCrystals) // this is CPV
- return fADCpedestalCpv+ amp*fADCchanelCpv ;
- }
- return 0;
-}
-
-//____________________________________________________________________________
-void AliPHOSClusterizerv1::Exec(Option_t *option)
+void AliPHOSClusterizerv1::Digits2Clusters(Option_t *option)
{
- // Steering method to perform clusterization for events
- // in the range from fFirstEvent to fLastEvent.
- // This range is optionally set by SetEventRange().
- // if fLastEvent=-1 (by default), then process events until the end.
+ // Steering method to perform clusterization for one event
+ // The input is the tree with digits.
+ // The output is the tree with clusters.
if(strstr(option,"tim"))
gBenchmark->Start("PHOSClusterizer");
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(fFirstEvent, gime->MaxEvent()); // one event at the time
- Int_t nEvents = fLastEvent - fFirstEvent + 1;
-
- 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",fIsOldRCUFormat); // Read digits from raw data
- }
- fNumberOfEmcClusters = fNumberOfCpvClusters = 0 ;
+ MakeClusters() ;
- MakeClusters() ;
+ AliDebug(2,Form(" ---- Printing clusters (%d)\n",
+ fEMCRecPoints->GetEntries()));
+ if(AliLog::GetGlobalDebugLevel()>1)
+ fEMCRecPoints->Print();
- if(fToUnfold)
- MakeUnfolding() ;
+ if(fToUnfold)
+ MakeUnfolding();
- WriteRecPoints();
+ WriteRecPoints();
- if(strstr(option,"deb"))
- PrintRecPoints(option) ;
+ if(strstr(option,"deb"))
+ PrintRecPoints(option) ;
- //increment the total number of recpoints per run
- fRecPointsInRun += gime->EmcRecPoints()->GetEntriesFast() ;
- fRecPointsInRun += gime->CpvRecPoints()->GetEntriesFast() ;
- }
-
- if(fWrite) //do not unload in "on flight" mode
- Unload();
-
if(strstr(option,"tim")){
gBenchmark->Stop("PHOSClusterizer");
- AliInfo(Form("took %f seconds for Clusterizing %f seconds per event \n",
- gBenchmark->GetCpuTime("PHOSClusterizer"),
- gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents )) ;
+ AliInfo(Form("took %f seconds for Clusterizing\n",
+ gBenchmark->GetCpuTime("PHOSClusterizer")));
}
}
// Cluster will be fitted as a superposition of nPar/3 electromagnetic showers
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
- TClonesArray * digits = gime->Digits();
-
-
gMinuit->mncler(); // Reset Minuit's list of paramters
gMinuit->SetPrintLevel(-1) ; // No Printout
gMinuit->SetFCN(AliPHOSClusterizerv1::UnfoldingChiSquare) ;
TList * toMinuit = new TList();
toMinuit->AddAt(emcRP,0) ;
- toMinuit->AddAt(digits,1) ;
+ toMinuit->AddAt(fDigitsArr,1) ;
+ toMinuit->AddAt(fGeom,2) ;
gMinuit->SetObjectFit(toMinuit) ; // To tranfer pointer to UnfoldingChiSquare
Int_t iDigit ;
- const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
-
for(iDigit = 0; iDigit < nDigits; iDigit++){
digit = maxAt[iDigit];
Int_t relid[4] ;
Float_t x = 0.;
Float_t z = 0.;
- geom->AbsToRelNumbering(digit->GetId(), relid) ;
- geom->RelPosInModule(relid, x, z) ;
+ fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
+ fGeom->RelPosInModule(relid, x, z) ;
Float_t energy = maxAtEnergy[iDigit] ;
}
-//____________________________________________________________________________
-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();
- // fCalibData = new AliPHOSCalibData(gAlice->GetRunNumber()); //original
-
- fCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number
-
- if(!fCalibData)
- {
- if ( !gime->Digitizer() )
- gime->LoadDigitizer();
- AliPHOSDigitizer * dig = gime->Digitizer();
- fADCchanelEmc = dig->GetEMCchannel() ;
- fADCpedestalEmc = dig->GetEMCpedestal();
-
- fADCchanelCpv = dig->GetCPVchannel() ;
- fADCpedestalCpv = dig->GetCPVpedestal() ;
- }
-}
//____________________________________________________________________________
void AliPHOSClusterizerv1::Init()
// 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() ;
- if(!gime)
- gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
-
- AliPHOSGeometry * geom = gime->PHOSGeometry();
-
- fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
+ fEmcCrystals = fGeom->GetNModules() * fGeom->GetNCristalsInModule() ;
if(!gMinuit)
gMinuit = new TMinuit(100);
- if ( !gime->Clusterizer() ) {
- gime->PostClusterizer(this);
- }
+ if (!fgCalibData)
+ fgCalibData = new AliPHOSCalibData(-1); //use AliCDBManager's run number
+ if (fgCalibData->GetCalibDataEmc() == 0)
+ AliFatal("Calibration parameters for PHOS EMC not found. Stop reconstruction.\n");
+
}
//____________________________________________________________________________
fNumberOfCpvClusters = 0 ;
fNumberOfEmcClusters = 0 ;
-
- fCpvClusteringThreshold = 0.0;
- fEmcClusteringThreshold = 0.2;
-
- fEmcLocMaxCut = 0.03 ;
- fCpvLocMaxCut = 0.03 ;
- fEmcMinE = 0.01 ;
- fCpvMinE = 0.0 ;
+ const AliPHOSRecoParam* parEmc = AliPHOSReconstructor::GetRecoParamEmc();
+ if(!parEmc) AliFatal("Reconstruction parameters for EMC not set!");
- fW0 = 4.5 ;
- fW0CPV = 4.0 ;
+ const AliPHOSRecoParam* parCpv = AliPHOSReconstructor::GetRecoParamCpv();
+ if(!parCpv) AliFatal("Reconstruction parameters for CPV not set!");
- fEmcTimeGate = 1.e-8 ;
+ fCpvClusteringThreshold = parCpv->GetClusteringThreshold();
+ fEmcClusteringThreshold = parEmc->GetClusteringThreshold();
- fToUnfold = kTRUE ;
-
- fRecPointsInRun = 0 ;
-
- fWrite = kTRUE ;
-
- fCalibData = 0 ;
+ fEmcLocMaxCut = parEmc->GetLocalMaxCut();
+ fCpvLocMaxCut = parCpv->GetLocalMaxCut();
- SetEventRange(0,-1) ;
+ fW0 = parEmc->GetLogWeight();
+ fW0CPV = parCpv->GetLogWeight();
- fIsOldRCUFormat = kFALSE;
+ fEmcTimeGate = 1.e-6 ;
+
+ fToUnfold = parEmc->ToUnfold() ;
+
+ fWrite = kTRUE ;
}
//____________________________________________________________________________
// The order of d1 and d2 is important: first (d1) should be a digit already in a cluster
// which is compared to a digit (d2) not yet in a cluster
- AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
-
Int_t rv = 0 ;
Int_t relid1[4] ;
- geom->AbsToRelNumbering(d1->GetId(), relid1) ;
+ fGeom->AbsToRelNumbering(d1->GetId(), relid1) ;
Int_t relid2[4] ;
- geom->AbsToRelNumbering(d2->GetId(), relid2) ;
+ fGeom->AbsToRelNumbering(d2->GetId(), relid2) ;
if ( (relid1[0] == relid2[0]) && (relid1[1]==relid2[1]) ) { // inside the same PHOS module
Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
- if (( coldiff <= 1 ) && ( rowdiff <= 1 )){
+ if (( coldiff <= 1 ) && ( rowdiff <= 1 )){ //At least common vertex
+ // if (( relid1[2]==relid2[2] && coldiff <= 1 ) || ( relid1[3]==relid2[3] && rowdiff <= 1 )){ //common side
if((relid1[1] != 0) || (TMath::Abs(d1->GetTime() - d2->GetTime() ) < fEmcTimeGate))
rv = 1 ;
}
return rv ;
}
//____________________________________________________________________________
-void AliPHOSClusterizerv1::CleanDigits(TClonesArray * digits)
-{
- // Remove digits with amplitudes below threshold
-
- for(Int_t i=0; i<digits->GetEntriesFast(); i++){
- AliPHOSDigit * digit = static_cast<AliPHOSDigit*>(digits->At(i)) ;
- if ( (IsInEmc(digit) && CalibrateEMC(digit->GetEnergy(),digit->GetId()) < fEmcMinE) ||
- (IsInCpv(digit) && CalibrateCPV(digit->GetAmp() ,digit->GetId()) < fCpvMinE) )
- 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
{
// Tells if (true) or not (false) the digit is in a PHOS-EMC module
Bool_t rv = kFALSE ;
- AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
- Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
+ Int_t nEMC = fGeom->GetNModules()*fGeom->GetNPhi()*fGeom->GetNZ();
if(digit->GetId() <= nEMC ) rv = kTRUE;
Bool_t rv = kFALSE ;
- AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
-
- Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
+ Int_t nEMC = fGeom->GetNModules()*fGeom->GetNPhi()*fGeom->GetNZ();
if(digit->GetId() > nEMC ) rv = kTRUE;
return rv ;
}
-//____________________________________________________________________________
-void AliPHOSClusterizerv1::Unload()
-{
- AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
- gime->PhosLoader()->UnloadDigits() ;
- gime->PhosLoader()->UnloadRecPoints() ;
-}
-
//____________________________________________________________________________
void AliPHOSClusterizerv1::WriteRecPoints()
{
// Creates new branches with given title
// fills and writes into TreeR.
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- TObjArray * emcRecPoints = gime->EmcRecPoints() ;
- TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
- TClonesArray * digits = gime->Digits() ;
-
-
Int_t index ;
//Evaluate position, dispersion and other RecPoint properties..
- Int_t nEmc = emcRecPoints->GetEntriesFast();
+ Int_t nEmc = fEMCRecPoints->GetEntriesFast();
+ Float_t emcMinE= AliPHOSReconstructor::GetRecoParamEmc()->GetMinE(); //Minimal digit energy
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) ;
+ AliPHOSEmcRecPoint * rp =
+ dynamic_cast<AliPHOSEmcRecPoint *>( fEMCRecPoints->At(index) );
+ rp->Purify(emcMinE) ;
+ if(rp->GetMultiplicity()==0){
+ fEMCRecPoints->RemoveAt(index) ;
delete rp ;
+ continue;
}
+
+ // No vertex is available now, calculate corrections in PID
+ rp->EvalAll(fW0,fDigitsArr) ;
+ TVector3 fakeVtx(0.,0.,0.) ;
+ rp->EvalAll(fW0,fakeVtx,fDigitsArr) ;
+ rp->EvalLocal2TrackingCSTransform();
}
- emcRecPoints->Compress() ;
- emcRecPoints->Sort() ;
- // emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ;
- for(index = 0; index < emcRecPoints->GetEntries(); index++){
- dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->SetIndexInList(index) ;
+ fEMCRecPoints->Compress() ;
+ fEMCRecPoints->Sort() ;
+ // fEMCRecPoints->Expand(fEMCRecPoints->GetEntriesFast()) ;
+ for(index = 0; index < fEMCRecPoints->GetEntries(); index++){
+ dynamic_cast<AliPHOSEmcRecPoint *>( fEMCRecPoints->At(index) )->SetIndexInList(index) ;
}
+ //For each rec.point set the distance to the nearest bad crystal (BVP)
+ SetDistancesToBadChannels();
+
//Now the same for CPV
- for(index = 0; index < cpvRecPoints->GetEntries(); index++){
- AliPHOSCpvRecPoint * rp = dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) );
- rp->EvalAll(fW0CPV,digits) ;
+ for(index = 0; index < fCPVRecPoints->GetEntries(); index++){
+ AliPHOSCpvRecPoint * rp = dynamic_cast<AliPHOSCpvRecPoint *>( fCPVRecPoints->At(index) );
+ rp->EvalAll(fW0CPV,fDigitsArr) ;
+ rp->EvalLocal2TrackingCSTransform();
}
- cpvRecPoints->Sort() ;
+ fCPVRecPoints->Sort() ;
- for(index = 0; index < cpvRecPoints->GetEntries(); index++)
- dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
+ for(index = 0; index < fCPVRecPoints->GetEntries(); index++)
+ dynamic_cast<AliPHOSCpvRecPoint *>( fCPVRecPoints->At(index) )->SetIndexInList(index) ;
- cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ;
+ fCPVRecPoints->Expand(fCPVRecPoints->GetEntriesFast()) ;
if(fWrite){ //We write TreeR
- TTree * treeR = gime->TreeR();
-
- Int_t bufferSize = 32000 ;
- Int_t splitlevel = 0 ;
-
- //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");
+ fTreeR->Fill();
}
}
// Steering method to construct the clusters stored in a list of Reconstructed Points
// A cluster is defined as a list of neighbour digits
-
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- TObjArray * emcRecPoints = gime->EmcRecPoints() ;
- TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
-
- emcRecPoints->Delete() ;
- cpvRecPoints->Delete() ;
-
- TClonesArray * digits = gime->Digits() ;
-
- //Remove digits below threshold
- CleanDigits(digits) ;
-
- TClonesArray * digitsC = static_cast<TClonesArray*>( digits->Clone() ) ;
+ TClonesArray * digitsC = static_cast<TClonesArray*>( fDigitsArr->Clone() ) ;
// Clusterization starts
TArrayI clusterdigitslist(1500) ;
Int_t index ;
- if (( IsInEmc (digit) &&
- CalibrateEMC(digit->GetEnergy(),digit->GetId()) > fEmcClusteringThreshold ) ||
- ( IsInCpv (digit) &&
- CalibrateCPV(digit->GetAmp() ,digit->GetId()) > fCpvClusteringThreshold ) ) {
+ if (( IsInEmc(digit) && digit->GetEnergy() > fEmcClusteringThreshold ) ||
+ ( IsInCpv(digit) && digit->GetEnergy() > fCpvClusteringThreshold ) ) {
Int_t iDigitInCluster = 0 ;
if ( IsInEmc(digit) ) {
// start a new EMC RecPoint
- if(fNumberOfEmcClusters >= emcRecPoints->GetSize())
- emcRecPoints->Expand(2*fNumberOfEmcClusters+1) ;
+ if(fNumberOfEmcClusters >= fEMCRecPoints->GetSize())
+ fEMCRecPoints->Expand(2*fNumberOfEmcClusters+1) ;
- emcRecPoints->AddAt(new AliPHOSEmcRecPoint(""), fNumberOfEmcClusters) ;
- clu = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(fNumberOfEmcClusters) ) ;
+ fEMCRecPoints->AddAt(new AliPHOSEmcRecPoint(""), fNumberOfEmcClusters) ;
+ clu = dynamic_cast<AliPHOSEmcRecPoint *>( fEMCRecPoints->At(fNumberOfEmcClusters) ) ;
fNumberOfEmcClusters++ ;
- clu->AddDigit(*digit, CalibrateEMC(digit->GetEnergy(),digit->GetId())) ;
+ clu->AddDigit(*digit, digit->GetEnergy()) ;
clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;
iDigitInCluster++ ;
digitsC->Remove(digit) ;
} else {
// start a new CPV cluster
- if(fNumberOfCpvClusters >= cpvRecPoints->GetSize())
- cpvRecPoints->Expand(2*fNumberOfCpvClusters+1);
+ if(fNumberOfCpvClusters >= fCPVRecPoints->GetSize())
+ fCPVRecPoints->Expand(2*fNumberOfCpvClusters+1);
- cpvRecPoints->AddAt(new AliPHOSCpvRecPoint(""), fNumberOfCpvClusters) ;
+ fCPVRecPoints->AddAt(new AliPHOSCpvRecPoint(""), fNumberOfCpvClusters) ;
- clu = dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(fNumberOfCpvClusters) ) ;
+ clu = dynamic_cast<AliPHOSCpvRecPoint *>( fCPVRecPoints->At(fNumberOfCpvClusters) ) ;
fNumberOfCpvClusters++ ;
- clu->AddDigit(*digit, CalibrateCPV(digit->GetAmp(),digit->GetId()) ) ;
+ clu->AddDigit(*digit, digit->GetEnergy()) ;
clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;
iDigitInCluster++ ;
digitsC->Remove(digit) ;
AliPHOSDigit * digitN ;
index = 0 ;
while (index < iDigitInCluster){ // scan over digits already in cluster
- digit = dynamic_cast<AliPHOSDigit*>( digits->At(clusterdigitslist[index]) ) ;
+ digit = dynamic_cast<AliPHOSDigit*>( fDigitsArr->At(clusterdigitslist[index]) ) ;
index++ ;
while ( (digitN = dynamic_cast<AliPHOSDigit *>( nextdigit() ) ) ) { // scan over the reduced list of digits
Int_t ineb = AreNeighbours(digit, digitN); // call (digit,digitN) in THAT oder !!!!!
case 0 : // not a neighbour
break ;
case 1 : // are neighbours
- if (IsInEmc (digitN))
- clu->AddDigit(*digitN, CalibrateEMC( digitN->GetEnergy(), digitN->GetId() ) );
- else
- clu->AddDigit(*digitN, CalibrateCPV( digitN->GetAmp() , digitN->GetId() ) );
-
+ clu->AddDigit(*digitN, digitN->GetEnergy());
clusterdigitslist[iDigitInCluster] = digitN->GetIndexInList() ;
iDigitInCluster++ ;
digitsC->Remove(digitN) ;
// Unfolds clusters using the shape of an ElectroMagnetic shower
// Performs unfolding of all EMC/CPV clusters
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
-
- TObjArray * emcRecPoints = gime->EmcRecPoints() ;
- TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
- TClonesArray * digits = gime->Digits() ;
-
// Unfold first EMC clusters
if(fNumberOfEmcClusters > 0){
- Int_t nModulesToUnfold = geom->GetNModules() ;
+ Int_t nModulesToUnfold = fGeom->GetNModules() ;
Int_t numberofNotUnfolded = fNumberOfEmcClusters ;
Int_t index ;
for(index = 0 ; index < numberofNotUnfolded ; index++){
- AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) ) ;
+ AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast<AliPHOSEmcRecPoint *>( fEMCRecPoints->At(index) ) ;
if(emcRecPoint->GetPHOSMod()> nModulesToUnfold)
break ;
Int_t nMultipl = emcRecPoint->GetMultiplicity() ;
AliPHOSDigit ** maxAt = new AliPHOSDigit*[nMultipl] ;
Float_t * maxAtEnergy = new Float_t[nMultipl] ;
- Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fEmcLocMaxCut,digits) ;
+ Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fEmcLocMaxCut,fDigitsArr) ;
if( nMax > 1 ) { // if cluster is very flat (no pronounced maximum) then nMax = 0
UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
- emcRecPoints->Remove(emcRecPoint);
- emcRecPoints->Compress() ;
+
+ fEMCRecPoints->Remove(emcRecPoint);
+ fEMCRecPoints->Compress() ;
index-- ;
fNumberOfEmcClusters -- ;
numberofNotUnfolded-- ;
// Unfold now CPV clusters
if(fNumberOfCpvClusters > 0){
- Int_t nModulesToUnfold = geom->GetNModules() ;
+ Int_t nModulesToUnfold = fGeom->GetNModules() ;
Int_t numberofCpvNotUnfolded = fNumberOfCpvClusters ;
Int_t index ;
for(index = 0 ; index < numberofCpvNotUnfolded ; index++){
- AliPHOSRecPoint * recPoint = dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) ) ;
+ AliPHOSRecPoint * recPoint = dynamic_cast<AliPHOSRecPoint *>( fCPVRecPoints->At(index) ) ;
if(recPoint->GetPHOSMod()> nModulesToUnfold)
break ;
Int_t nMultipl = emcRecPoint->GetMultiplicity() ;
AliPHOSDigit ** maxAt = new AliPHOSDigit*[nMultipl] ;
Float_t * maxAtEnergy = new Float_t[nMultipl] ;
- Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fCpvLocMaxCut,digits) ;
+ Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fCpvLocMaxCut,fDigitsArr) ;
if( nMax > 1 ) { // if cluster is very flat (no pronounced maximum) then nMax = 0
UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
- cpvRecPoints->Remove(emcRecPoint);
- cpvRecPoints->Compress() ;
+ fCPVRecPoints->Remove(emcRecPoint);
+ fCPVRecPoints->Compress() ;
index-- ;
numberofCpvNotUnfolded-- ;
fNumberOfCpvClusters-- ;
}
//____________________________________________________________________________
-Double_t AliPHOSClusterizerv1::ShowerShape(Double_t r)
+Double_t AliPHOSClusterizerv1::ShowerShape(Double_t x, Double_t z)
{
// Shape of the shower (see PHOS TDR)
// If you change this function, change also the gradient evaluation in ChiSquare()
- Double_t r4 = r*r*r*r ;
- Double_t r295 = TMath::Power(r, 2.95) ;
+ //for the moment we neglect dependence on the incident angle.
+
+ Double_t r2 = x*x + z*z ;
+ Double_t r4 = r2*r2 ;
+ Double_t r295 = TMath::Power(r2, 2.95/2.) ;
Double_t shape = TMath::Exp( -r4 * (1. / (2.32 + 0.26 * r4) + 0.0316 / (1 + 0.0652 * r295) ) ) ;
return shape ;
}
{
// Performs the unfolding of a cluster with nMax overlapping showers
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
-
- const TClonesArray * digits = gime->Digits() ;
- TObjArray * emcRecPoints = gime->EmcRecPoints() ;
- TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
-
Int_t nPar = 3 * nMax ;
Float_t * fitparameters = new Float_t[nPar] ;
Bool_t rv = FindFit(iniEmc, maxAt, maxAtEnergy, nPar, fitparameters) ;
+
if( !rv ) {
// Fit failed, return and remove cluster
iniEmc->SetNExMax(-1) ;
Int_t nDigits = iniEmc->GetMultiplicity() ;
Float_t * efit = new Float_t[nDigits] ;
- Float_t xDigit=0.,zDigit=0.,distance=0. ;
+ Float_t xDigit=0.,zDigit=0. ;
Float_t xpar=0.,zpar=0.,epar=0. ;
Int_t relid[4] ;
AliPHOSDigit * digit = 0 ;
Int_t * emcDigits = iniEmc->GetDigitsList() ;
+ TVector3 vIncid ;
+
Int_t iparam ;
Int_t iDigit ;
for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
- digit = dynamic_cast<AliPHOSDigit*>( digits->At(emcDigits[iDigit] ) ) ;
- geom->AbsToRelNumbering(digit->GetId(), relid) ;
- geom->RelPosInModule(relid, xDigit, zDigit) ;
+ digit = dynamic_cast<AliPHOSDigit*>( fDigitsArr->At(emcDigits[iDigit] ) ) ;
+ fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
+ fGeom->RelPosInModule(relid, xDigit, zDigit) ;
efit[iDigit] = 0;
iparam = 0 ;
zpar = fitparameters[iparam+1] ;
epar = fitparameters[iparam+2] ;
iparam += 3 ;
- distance = (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) ;
- distance = TMath::Sqrt(distance) ;
- efit[iDigit] += epar * ShowerShape(distance) ;
+// fGeom->GetIncidentVector(fVtx,relid[0],xpar,zpar,vIncid) ;
+// efit[iDigit] += epar * ShowerShape(xDigit - xpar,zDigit - zpar,vIncid) ;
+ efit[iDigit] += epar * ShowerShape(xDigit - xpar,zDigit - zpar) ;
}
}
-
// Now create new RecPoints and fill energy lists with efit corrected to fluctuations
// so that energy deposited in each cell is distributed betwin new clusters proportionally
// to its contribution to efit
zpar = fitparameters[iparam+1] ;
epar = fitparameters[iparam+2] ;
iparam += 3 ;
-
+// fGeom->GetIncidentVector(fVtx,iniEmc->GetPHOSMod(),xpar,zpar,vIncid) ;
+
AliPHOSEmcRecPoint * emcRP = 0 ;
- if(iniEmc->IsEmc()){ //create new entries in fEmcRecPoints...
+ if(iniEmc->IsEmc()){ //create new entries in fEMCRecPoints...
- if(fNumberOfEmcClusters >= emcRecPoints->GetSize())
- emcRecPoints->Expand(2*fNumberOfEmcClusters) ;
+ if(fNumberOfEmcClusters >= fEMCRecPoints->GetSize())
+ fEMCRecPoints->Expand(2*fNumberOfEmcClusters) ;
- (*emcRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ;
- emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(fNumberOfEmcClusters) ) ;
+ (*fEMCRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ;
+ emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( fEMCRecPoints->At(fNumberOfEmcClusters) ) ;
fNumberOfEmcClusters++ ;
emcRP->SetNExMax((Int_t)nPar/3) ;
}
- else{//create new entries in fCpvRecPoints
- if(fNumberOfCpvClusters >= cpvRecPoints->GetSize())
- cpvRecPoints->Expand(2*fNumberOfCpvClusters) ;
+ else{//create new entries in fCPVRecPoints
+ if(fNumberOfCpvClusters >= fCPVRecPoints->GetSize())
+ fCPVRecPoints->Expand(2*fNumberOfCpvClusters) ;
- (*cpvRecPoints)[fNumberOfCpvClusters] = new AliPHOSCpvRecPoint("") ;
- emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( cpvRecPoints->At(fNumberOfCpvClusters) ) ;
+ (*fCPVRecPoints)[fNumberOfCpvClusters] = new AliPHOSCpvRecPoint("") ;
+ emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( fCPVRecPoints->At(fNumberOfCpvClusters) ) ;
fNumberOfCpvClusters++ ;
}
Float_t eDigit ;
for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
- digit = dynamic_cast<AliPHOSDigit*>( digits->At( emcDigits[iDigit] ) ) ;
- geom->AbsToRelNumbering(digit->GetId(), relid) ;
- geom->RelPosInModule(relid, xDigit, zDigit) ;
- distance = (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) ;
- distance = TMath::Sqrt(distance) ;
- ratio = epar * ShowerShape(distance) / efit[iDigit] ;
+ digit = dynamic_cast<AliPHOSDigit*>( fDigitsArr->At( emcDigits[iDigit] ) ) ;
+ fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
+ fGeom->RelPosInModule(relid, xDigit, zDigit) ;
+// ratio = epar * ShowerShape(xDigit - xpar,zDigit - zpar,vIncid) / efit[iDigit] ;
+ ratio = epar * ShowerShape(xDigit - xpar,zDigit - zpar) / efit[iDigit] ;
eDigit = emcEnergies[iDigit] * ratio ;
emcRP->AddDigit( *digit, eDigit ) ;
}
AliPHOSEmcRecPoint * emcRP = dynamic_cast<AliPHOSEmcRecPoint*>( toMinuit->At(0) ) ;
TClonesArray * digits = dynamic_cast<TClonesArray*>( toMinuit->At(1) ) ;
+ // A bit buggy way to get an access to the geometry
+ // To be revised!
+ AliPHOSGeometry *geom = dynamic_cast<AliPHOSGeometry *>(toMinuit->At(2));
-
+// TVector3 * vtx = dynamic_cast<TVector3*>(toMinuit->At(3)) ; //Vertex position
// AliPHOSEmcRecPoint * emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( gMinuit->GetObjectFit() ) ; // EmcRecPoint to fit
Float_t * emcEnergies = emcRP->GetEnergiesList() ;
- const AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
+// TVector3 vInc ;
fret = 0. ;
Int_t iparam ;
Int_t iParam = 0 ;
efit = 0 ;
while(iParam < nPar ){
- Double_t distance = (xDigit - x[iParam]) * (xDigit - x[iParam]) ;
+ Double_t dx = (xDigit - x[iParam]) ;
iParam++ ;
- distance += (zDigit - x[iParam]) * (zDigit - x[iParam]) ;
- distance = TMath::Sqrt( distance ) ;
+ Double_t dz = (zDigit - x[iParam]) ;
iParam++ ;
- efit += x[iParam] * ShowerShape(distance) ;
+// fGeom->GetIncidentVector(*vtx,emcRP->GetPHOSMod(),x[iParam-2],x[iParam-1],vInc) ;
+// efit += x[iParam] * ShowerShape(dx,dz,vInc) ;
+ efit += x[iParam] * ShowerShape(dx,dz) ;
iParam++ ;
}
Double_t sum = 2. * (efit - emcEnergies[iDigit]) / emcEnergies[iDigit] ; // Here we assume, that sigma = sqrt(E)
Double_t xpar = x[iParam] ;
Double_t zpar = x[iParam+1] ;
Double_t epar = x[iParam+2] ;
+// fGeom->GetIncidentVector(*vtx,emcRP->GetPHOSMod(),xpar,zpar,vInc) ;
Double_t dr = TMath::Sqrt( (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) );
- Double_t shape = sum * ShowerShape(dr) ;
+// Double_t shape = sum * ShowerShape(xDigit - xpar,zDigit - zpar,vInc) ;
+ Double_t shape = sum * ShowerShape(xDigit - xpar,zDigit - zpar) ;
+//DP: No incident angle dependence in gradient yet!!!!!!
Double_t r4 = dr*dr*dr*dr ;
Double_t r295 = TMath::Power(dr,2.95) ;
Double_t deriv =-4. * dr*dr * ( 2.32 / ( (2.32 + 0.26 * r4) * (2.32 + 0.26 * r4) ) +
Double_t zpar = x[iparam+1] ;
Double_t epar = x[iparam+2] ;
iparam += 3 ;
- Double_t distance = (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) ;
- distance = TMath::Sqrt(distance) ;
- efit += epar * ShowerShape(distance) ;
+// fGeom->GetIncidentVector(*vtx,emcRP->GetPHOSMod(),xpar,zpar,vInc) ;
+// efit += epar * ShowerShape(xDigit - xpar,zDigit - zpar,vInc) ;
+ efit += epar * ShowerShape(xDigit - xpar,zDigit - zpar) ;
}
fret += (efit-emcEnergies[iDigit])*(efit-emcEnergies[iDigit])/emcEnergies[iDigit] ;
fCpvLocMaxCut,
fW0CPV )) ;
}
-
-
//____________________________________________________________________________
void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option)
{
// Prints list of RecPoints produced at the current pass of AliPHOSClusterizer
- AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
- TObjArray * emcRecPoints = gime->EmcRecPoints() ;
- TObjArray * cpvRecPoints = gime->CpvRecPoints() ;
-
- AliInfo(Form("\nevent %d \n Found %d EMC RecPoints and %d CPV RecPoints",
- gAlice->GetEvNumber(),
- emcRecPoints->GetEntriesFast(),
- cpvRecPoints->GetEntriesFast() )) ;
+ AliInfo(Form("\nFound %d EMC RecPoints and %d CPV RecPoints",
+ fEMCRecPoints->GetEntriesFast(),
+ fCPVRecPoints->GetEntriesFast() )) ;
- fRecPointsInRun += emcRecPoints->GetEntriesFast() ;
- fRecPointsInRun += cpvRecPoints->GetEntriesFast() ;
-
-
if(strstr(option,"all")) {
printf("\n EMC clusters \n") ;
printf("Index Ene(MeV) Multi Module X Y Z Lambdas_1 Lambda_2 # of prim Primaries list\n") ;
Int_t index ;
- for (index = 0 ; index < emcRecPoints->GetEntries() ; index++) {
- AliPHOSEmcRecPoint * rp = (AliPHOSEmcRecPoint * )emcRecPoints->At(index) ;
+ for (index = 0 ; index < fEMCRecPoints->GetEntries() ; index++) {
+ AliPHOSEmcRecPoint * rp = (AliPHOSEmcRecPoint * )fEMCRecPoints->At(index) ;
TVector3 locpos;
rp->GetLocalPosition(locpos);
Float_t lambda[2];
//Now plot CPV recPoints
printf("\n CPV clusters \n") ;
printf("Index Ene(MeV) Module X Y Z \n") ;
- for (index = 0 ; index < cpvRecPoints->GetEntries() ; index++) {
- AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )cpvRecPoints->At(index) ;
+ for (index = 0 ; index < fCPVRecPoints->GetEntries() ; index++) {
+ AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )fCPVRecPoints->At(index) ;
TVector3 locpos;
rp->GetLocalPosition(locpos);
}
}
+
+//____________________________________________________________________________
+void AliPHOSClusterizerv1::SetDistancesToBadChannels()
+{
+ //For each EMC rec. point set the distance to the nearest bad crystal.
+ //Author: Boris Polichtchouk
+
+ if(!fgCalibData->GetNumOfEmcBadChannels()) return;
+ AliInfo(Form("%d bad channel(s) found.\n",fgCalibData->GetNumOfEmcBadChannels()));
+
+ Int_t badIds[8000];
+ fgCalibData->EmcBadChannelIds(badIds);
+
+ AliPHOSEmcRecPoint* rp;
+
+ TMatrixF gmat;
+ TVector3 gposRecPoint; // global (in ALICE frame) position of rec. point
+ TVector3 gposBadChannel; // global position of bad crystal
+ TVector3 dR;
+
+ Float_t dist,minDist;
+
+ for(Int_t iRP=0; iRP<fEMCRecPoints->GetEntries(); iRP++){
+ rp = (AliPHOSEmcRecPoint*)fEMCRecPoints->At(iRP);
+ minDist = 1.e+07;
+
+ for(Int_t iBad=0; iBad<fgCalibData->GetNumOfEmcBadChannels(); iBad++) {
+ rp->GetGlobalPosition(gposRecPoint,gmat);
+ fGeom->RelPosInAlice(badIds[iBad],gposBadChannel);
+ AliDebug(2,Form("BC position:[%.3f,%.3f,%.3f], RP position:[%.3f,%.3f,%.3f]. E=%.3f\n",
+ gposBadChannel.X(),gposBadChannel.Y(),gposBadChannel.Z(),
+ gposRecPoint.X(),gposRecPoint.Y(),gposRecPoint.Z(),rp->GetEnergy()));
+ dR = gposBadChannel-gposRecPoint;
+ dist = dR.Mag();
+ if(dist<minDist) minDist = dist;
+ }
+
+ rp->SetDistanceToBadCrystal(minDist);
+ }
+
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff