#include "AliMFTSegmentation.h"
#include "TTree.h"
#include "TMath.h"
+#include "AliMFTConstants.h"
#include "AliMFTClusterFinder.h"
+#include "TRandom.h"
+
+const Double_t AliMFTClusterFinder::fCutForAvailableDigits = AliMFTConstants::fCutForAvailableDigits;
+const Double_t AliMFTClusterFinder::fCutForAttachingDigits = AliMFTConstants::fCutForAttachingDigits;
+const Double_t AliMFTClusterFinder::fMisalignmentMagnitude = AliMFTConstants::fMisalignmentMagnitude;
ClassImp(AliMFTClusterFinder)
AliMFTClusterFinder::AliMFTClusterFinder() :
TObject(),
fDigitsInCluster(0),
- fCurrentDig(0),
+ fCurrentDigit(0),
+ fCurrentCluster(0),
fSegmentation(0),
- fNPlanes(0)
+ fNPlanes(0),
+ fApplyMisalignment(kFALSE),
+ fStopWatch(0)
{
// Default constructor
for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) fClustersPerPlane[iPlane] = NULL;
fDigitsInCluster = new TClonesArray("AliMFTDigit", fNMaxDigitsPerCluster);
+ fDigitsInCluster -> SetOwner(kTRUE);
+ fStopWatch = new TStopwatch();
+ fStopWatch -> Reset();
}
//====================================================================================================================================================
AliMFTClusterFinder::~AliMFTClusterFinder() {
-
+
AliDebug(1, "Deleting AliMFTClusterFinder...");
+
+ for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
+ if (fClustersPerPlane[iPlane]) fClustersPerPlane[iPlane]->Delete(); delete fClustersPerPlane[iPlane]; fClustersPerPlane[iPlane] = 0x0;
+ }
+ if (fDigitsInCluster) fDigitsInCluster->Delete(); delete fDigitsInCluster; fDigitsInCluster = NULL;
+ delete fSegmentation; fSegmentation=NULL;
- for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) delete fClustersPerPlane[iPlane];
-
+ delete fStopWatch;
+
AliDebug(1, "... done!");
+
+}
+
+//====================================================================================================================================================
+void AliMFTClusterFinder::Clear(const Option_t* /*opt*/) {
+
+ AliDebug(1, "Clearing AliMFTClusterFinder...");
+
+ for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
+ if(fClustersPerPlane[iPlane]) fClustersPerPlane[iPlane]->Clear("C");
+ }
+ if(fDigitsInCluster) fDigitsInCluster->Clear("C");
+ fSegmentation->Clear("C");
+
+ AliDebug(1, "... done!");
+
}
//====================================================================================================================================================
AliDebug(1, "Starting Event...");
for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
- fClustersPerPlane[iPlane]->Clear();
+ fClustersPerPlane[iPlane]->Delete();
}
AliDebug(1, "... done!");
AliInfo("Starting Clusterization for MFT");
AliDebug(1, Form("nPlanes = %d", fNPlanes));
+ if (!fStopWatch) fStopWatch = new TStopwatch();
+ fStopWatch -> Reset();
+
StartEvent();
+ Bool_t isDigAvailableForNewCluster = kTRUE;
+ TClonesArray *myDigitList = 0;
+
for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
AliDebug(1, Form("Plane %02d", iPlane));
- TClonesArray *myDigitList = (TClonesArray*) pDigitList->At(iPlane);
- // myDigitList->Sort();
-
- AliDebug(1, Form("myDigitList->GetEntries() = %d", myDigitList->GetEntries()));
-
- while (myDigitList->GetEntries()) {
-
- fDigitsInCluster->Clear();
-
- Bool_t clusterUpdated=kTRUE;
-
- //---------------------------------------------------------------------------------------------------------
-
- while (clusterUpdated) { // repeat the loop on the digits until no new compatible digit is found
-
- clusterUpdated = kFALSE;
-
- for (Int_t iDig=0; iDig<myDigitList->GetEntries(); iDig++) {
- fCurrentDig = (AliMFTDigit*) myDigitList->At(iDig);
- if (fDigitsInCluster->GetEntries()<fNMaxDigitsPerCluster) {
- if (fDigitsInCluster->GetEntries()==0) {
- new ((*fDigitsInCluster)[fDigitsInCluster->GetEntries()]) AliMFTDigit(*fCurrentDig);
- myDigitList->Remove(fCurrentDig);
- clusterUpdated = kTRUE;
- }
- else if (IsCurrentDigitCompatible()) {
- new ((*fDigitsInCluster)[fDigitsInCluster->GetEntries()]) AliMFTDigit(*fCurrentDig);
- myDigitList->Remove(fCurrentDig);
- clusterUpdated = kTRUE;
- }
- }
- }
-
- if (clusterUpdated) myDigitList->Compress();
-
- }
-
- //---------------------------------------------------------------------------------------------------------
-
- AliDebug(1, "Building new cluster");
+ Int_t nDetElem = fSegmentation->GetPlane(iPlane)->GetNActiveElements();
+ TClonesArray *clustersPerDetElem[fNMaxDetElemPerPlane] = {0};
+ for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) clustersPerDetElem[iDetElem] = new TClonesArray("AliMFTCluster");
- BuildNewCluster(iPlane); // here the new cluster is built
+ myDigitList = (TClonesArray*) pDigitList->At(iPlane);
- }
+ AliDebug(1, Form("myDigitList->GetEntries() = %d", myDigitList->GetEntries()));
- }
+ Int_t cycleOverDigits = 0;
+ Double_t myCutForAvailableDigits = 0;
+
+ Int_t currentDetElem = -1;
+ Int_t currentDetElemLocal = -1;
+ Bool_t areThereSkippedDigits = kFALSE;
-}
+ fStopWatch -> Start();
-//====================================================================================================================================================
-
-Bool_t AliMFTClusterFinder::IsCurrentDigitCompatible() {
+ while (myDigitList->GetEntries()) {
- // where it is decided if the current digit (fCurrentDig) is compatible with the current digits array (fDigitsInCluster)
+ for (Int_t iDig=0; iDig<myDigitList->GetEntries(); iDig++) {
- for (Int_t iDigit=0; iDigit<fDigitsInCluster->GetEntries(); iDigit++) {
- AliMFTDigit *tmpDig = (AliMFTDigit*) fDigitsInCluster->At(iDigit);
- Int_t distX = TMath::Abs(tmpDig->GetPixelX() - fCurrentDig->GetPixelX());
- Int_t distY = TMath::Abs(tmpDig->GetPixelY() - fCurrentDig->GetPixelY());
- if (distX<=1 && distY<=1) return kTRUE;
- }
+ fCurrentDigit = (AliMFTDigit*) myDigitList->At(iDig);
- return kFALSE;
-
-}
+ if (!iDig) {
+ if (fCurrentDigit->GetDetElemID() != currentDetElem) {
+ // first iteration over the digits of a specific detection element
+ currentDetElem = fCurrentDigit->GetDetElemID();
+ currentDetElemLocal = fSegmentation->GetDetElemLocalID(currentDetElem);
+ cycleOverDigits = 0;
+ myCutForAvailableDigits = fCutForAvailableDigits;
+ }
+ else if (fCurrentDigit->GetDetElemID()==currentDetElem && areThereSkippedDigits) {
+ // second (or further) iteration over the digits of a specific detection element
+ cycleOverDigits++;
+ myCutForAvailableDigits -= 0.5;
+ }
+ areThereSkippedDigits = kFALSE;
+ }
+ else {
+ areThereSkippedDigits = kTRUE;
+ if (fCurrentDigit->GetDetElemID() != currentDetElem) break;
+ }
-//====================================================================================================================================================
+ isDigAvailableForNewCluster = kTRUE;
+
+ for (Int_t iCluster=0; iCluster<clustersPerDetElem[currentDetElemLocal]->GetEntries(); iCluster++) {
+ fCurrentCluster = (AliMFTCluster*) clustersPerDetElem[currentDetElemLocal]->At(iCluster);
+ if (fCurrentCluster->GetDistanceFromPixel(fCurrentDigit) < fCutForAttachingDigits) {
+ fCurrentCluster->AddPixel(fCurrentDigit);
+ myDigitList->Remove(fCurrentDigit);
+ myDigitList->Compress();
+ iDig--;
+ isDigAvailableForNewCluster = kFALSE;
+ break;
+ }
+ if (fCurrentCluster->GetDistanceFromPixel(fCurrentDigit) < myCutForAvailableDigits) isDigAvailableForNewCluster=kFALSE;
+ }
-void AliMFTClusterFinder::BuildNewCluster(Int_t plane) {
+ if (isDigAvailableForNewCluster) {
+ AliMFTCluster *newCluster = new AliMFTCluster();
+ newCluster->AddPixel(fCurrentDigit);
+ myDigitList->Remove(fCurrentDigit);
+ myDigitList->Compress();
+ iDig--;
+ new ((*clustersPerDetElem[currentDetElemLocal])[clustersPerDetElem[currentDetElemLocal]->GetEntries()]) AliMFTCluster(*newCluster);
+ delete newCluster;
+ }
+
+ } // end of cycle over the digits
- // where a new cluster is built, starting from the array of compatible digits (fDigitsInCluster)
+ } // no more digits to check in current plane!
- AliDebug(1, Form("Starting cluster building from %d digits", fDigitsInCluster->GetEntries()));
+ fStopWatch -> Print("m");
- AliMFTCluster *newCluster = new AliMFTCluster();
+ printf("Plane %d: clusters found in %f seconds\n",iPlane,fStopWatch->CpuTime());
+ fStopWatch->Start();
- Double_t xCenters[fNMaxDigitsPerCluster] = {0};
- Double_t yCenters[fNMaxDigitsPerCluster] = {0};
- Double_t nElectrons = 0.;
+ // Now we merge the cluster lists coming from each detection element, to build the cluster list of the plane
- for (Int_t iDigit=0; iDigit<fDigitsInCluster->GetEntries(); iDigit++) {
- AliMFTDigit *tmpDig = (AliMFTDigit*) fDigitsInCluster->At(iDigit);
- xCenters[iDigit] = tmpDig->GetPixelCenterX();
- yCenters[iDigit] = tmpDig->GetPixelCenterY();
- nElectrons += tmpDig->GetNElectrons();
- for (Int_t iTrack=0; iTrack<tmpDig->GetNMCTracks(); iTrack++) newCluster->AddMCLabel(tmpDig->GetMCLabel(iTrack));
- }
+ Double_t misalignmentPhi = 2.*TMath::Pi()*gRandom->Rndm();
+ Double_t misalignmentX = fMisalignmentMagnitude*TMath::Cos(misalignmentPhi);
+ Double_t misalignmentY = fMisalignmentMagnitude*TMath::Sin(misalignmentPhi);
- newCluster -> SetX(TMath::Mean(fDigitsInCluster->GetEntries(), xCenters));
- newCluster -> SetY(TMath::Mean(fDigitsInCluster->GetEntries(), yCenters));
- newCluster -> SetZ(((AliMFTDigit*) fDigitsInCluster->At(0))->GetPixelCenterZ());
+ AliMFTCluster *newCluster = NULL;
+ for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) {
+ for (Int_t iCluster=0; iCluster<clustersPerDetElem[iDetElem]->GetEntries(); iCluster++) {
+ newCluster = (AliMFTCluster*) (clustersPerDetElem[iDetElem]->At(iCluster));
+ newCluster -> TerminateCluster();
+ if (fApplyMisalignment) {
+ newCluster -> SetClusterEditable(kTRUE);
+ newCluster -> SetX(newCluster->GetX()+misalignmentX);
+ newCluster -> SetY(newCluster->GetY()+misalignmentY);
+ newCluster -> SetClusterEditable(kFALSE);
+ }
+ new ((*fClustersPerPlane[iPlane])[fClustersPerPlane[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
+ }
+ }
- Double_t minErrX = ((AliMFTDigit*) fDigitsInCluster->At(0))->GetPixelWidthX() / TMath::Sqrt(12.);
- Double_t minErrY = ((AliMFTDigit*) fDigitsInCluster->At(0))->GetPixelWidthY() / TMath::Sqrt(12.);
- Double_t minErrZ = ((AliMFTDigit*) fDigitsInCluster->At(0))->GetPixelWidthZ() / TMath::Sqrt(12.);
- newCluster -> SetErrX( TMath::Max(TMath::RMS(fDigitsInCluster->GetEntries(), xCenters), minErrX) );
- newCluster -> SetErrY( TMath::Max(TMath::RMS(fDigitsInCluster->GetEntries(), yCenters), minErrY) );
- newCluster -> SetErrZ( minErrZ );
-
- newCluster -> SetNElectrons(nElectrons);
- newCluster -> SetPlane(plane);
+ printf("%d Clusters in plane %02d merged in %f seconds\n", fClustersPerPlane[iPlane]->GetEntries(), iPlane, fStopWatch->CpuTime());
- newCluster -> SetSize(fDigitsInCluster->GetEntries());
+ for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) {
+ clustersPerDetElem[iDetElem] -> Delete();
+ delete clustersPerDetElem[iDetElem];
+ }
- AliDebug(1, Form("Adding cluster #%02d to tree (%f, %f, %f)",
- fClustersPerPlane[plane]->GetEntries(), newCluster->GetX(), newCluster->GetY(), newCluster->GetZ()));
+ myDigitList -> Delete();
- new ((*fClustersPerPlane[plane])[fClustersPerPlane[plane]->GetEntries()]) AliMFTCluster(*newCluster);
+ } // end of cycle over the planes
}
if (treeCluster) {
for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
AliDebug(1, Form("Setting Branch Plane_%02d for treeCluster",iPlane));
- TBranch *branch = treeCluster->GetBranch(Form("Plane_%02d",iPlane));
- if (branch) continue;
+ if (treeCluster->GetBranch(Form("Plane_%02d",iPlane))) continue;
AliDebug(1, Form("Branch Plane_%02d does not exist, creating!",iPlane));
- branch = treeCluster->Branch(Form("Plane_%02d",iPlane), &(fClustersPerPlane[iPlane]));
+ treeCluster->Branch(Form("Plane_%02d",iPlane), &(fClustersPerPlane[iPlane]));
}
}
for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
AliDebug(1, Form("plane %02d", iPlane));
fClustersPerPlane[iPlane] = new TClonesArray("AliMFTCluster");
+ fClustersPerPlane[iPlane] -> SetOwner(kTRUE);
+
}
}