1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 //====================================================================================================================================================
18 // Class for finding and building the clusters of the ALICE Muon Forward Tracker
20 // Contact author: antonio.uras@cern.ch
22 //====================================================================================================================================================
25 #include "TObjArray.h"
26 #include "TClonesArray.h"
27 #include "AliMFTDigit.h"
28 #include "AliMFTCluster.h"
29 #include "AliMFTSegmentation.h"
32 #include "AliMFTConstants.h"
33 #include "AliMFTClusterFinder.h"
35 const Double_t AliMFTClusterFinder::fCutForAvailableDigits = AliMFTConstants::fCutForAvailableDigits;
36 const Double_t AliMFTClusterFinder::fCutForAttachingDigits = AliMFTConstants::fCutForAttachingDigits;
39 ClassImp(AliMFTClusterFinder)
41 //====================================================================================================================================================
43 AliMFTClusterFinder::AliMFTClusterFinder() :
53 // Default constructor
55 for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) fClustersPerPlane[iPlane] = NULL;
56 fDigitsInCluster = new TClonesArray("AliMFTDigit", fNMaxDigitsPerCluster);
57 fDigitsInCluster -> SetOwner(kTRUE);
58 sw = new TStopwatch();
63 //====================================================================================================================================================
65 AliMFTClusterFinder::~AliMFTClusterFinder() {
67 AliDebug(1, "Deleting AliMFTClusterFinder...");
69 for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
70 if (fClustersPerPlane[iPlane]) fClustersPerPlane[iPlane]->Delete(); delete fClustersPerPlane[iPlane]; fClustersPerPlane[iPlane] = 0x0;
72 if (fDigitsInCluster) fDigitsInCluster->Delete(); delete fDigitsInCluster; fDigitsInCluster = NULL;
73 delete fSegmentation; fSegmentation=NULL;
77 AliDebug(1, "... done!");
81 //====================================================================================================================================================
83 void AliMFTClusterFinder::Clear(const Option_t* /*opt*/) {
85 AliDebug(1, "Clearing AliMFTClusterFinder...");
87 for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
88 if(fClustersPerPlane[iPlane]) fClustersPerPlane[iPlane]->Clear("C");
90 if(fDigitsInCluster) fDigitsInCluster->Clear("C");
91 fSegmentation->Clear("C");
93 AliDebug(1, "... done!");
97 //====================================================================================================================================================
99 void AliMFTClusterFinder::Init(const Char_t *nameGeomFile) {
101 fSegmentation = new AliMFTSegmentation(nameGeomFile);
102 fNPlanes = fSegmentation -> GetNPlanes();
106 //====================================================================================================================================================
108 void AliMFTClusterFinder::StartEvent() {
110 // Cleaning up and preparation for the clustering procedure
112 AliDebug(1, "Starting Event...");
114 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
115 fClustersPerPlane[iPlane]->Delete();
118 AliDebug(1, "... done!");
122 //====================================================================================================================================================
124 void AliMFTClusterFinder::DigitsToClusters(const TObjArray *pDigitList) {
126 // where the clusterization is performed
128 AliInfo("Starting Clusterization for MFT");
129 AliDebug(1, Form("nPlanes = %d", fNPlanes));
131 if (!sw) sw = new TStopwatch();
135 Bool_t isDigAvailableForNewCluster = kTRUE;
137 TClonesArray *myDigitList = 0;
139 for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
141 AliDebug(1, Form("Plane %02d", iPlane));
143 Int_t nDetElem = fSegmentation->GetPlane(iPlane)->GetNActiveElements();
144 TClonesArray *clustersPerDetElem[fNMaxDetElemPerPlane] = {0};
145 for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) clustersPerDetElem[iDetElem] = new TClonesArray("AliMFTCluster");
147 myDigitList = (TClonesArray*) pDigitList->At(iPlane);
149 AliDebug(1, Form("myDigitList->GetEntries() = %d", myDigitList->GetEntries()));
151 Int_t cycleOverDigits = 0;
152 Double_t myCutForAvailableDigits = 0;
154 Int_t currentDetElem = -1;
155 Int_t currentDetElemLocal = -1;
156 Bool_t areThereSkippedDigits = kFALSE;
160 while (myDigitList->GetEntries()) {
162 for (Int_t iDig=0; iDig<myDigitList->GetEntries(); iDig++) {
164 fCurrentDigit = (AliMFTDigit*) myDigitList->At(iDig);
167 if (fCurrentDigit->GetDetElemID() != currentDetElem) {
168 // first iteration over the digits of a specific detection element
169 currentDetElem = fCurrentDigit->GetDetElemID();
170 currentDetElemLocal = fSegmentation->GetDetElemLocalID(currentDetElem);
172 myCutForAvailableDigits = fCutForAvailableDigits;
174 else if (fCurrentDigit->GetDetElemID()==currentDetElem && areThereSkippedDigits) {
175 // second (or further) iteration over the digits of a specific detection element
177 myCutForAvailableDigits -= 0.5;
179 areThereSkippedDigits = kFALSE;
182 areThereSkippedDigits = kTRUE;
183 if (fCurrentDigit->GetDetElemID() != currentDetElem) break;
186 isDigAvailableForNewCluster = kTRUE;
188 for (Int_t iCluster=0; iCluster<clustersPerDetElem[currentDetElemLocal]->GetEntries(); iCluster++) {
189 fCurrentCluster = (AliMFTCluster*) clustersPerDetElem[currentDetElemLocal]->At(iCluster);
190 if (fCurrentCluster->GetDistanceFromPixel(fCurrentDigit) < fCutForAttachingDigits) {
191 fCurrentCluster->AddPixel(fCurrentDigit);
192 myDigitList->Remove(fCurrentDigit);
193 myDigitList->Compress();
195 isDigAvailableForNewCluster = kFALSE;
198 if (fCurrentCluster->GetDistanceFromPixel(fCurrentDigit) < myCutForAvailableDigits) isDigAvailableForNewCluster=kFALSE;
201 if (isDigAvailableForNewCluster) {
202 AliMFTCluster *newCluster = new AliMFTCluster();
203 newCluster->AddPixel(fCurrentDigit);
204 myDigitList->Remove(fCurrentDigit);
205 myDigitList->Compress();
207 new ((*clustersPerDetElem[currentDetElemLocal])[clustersPerDetElem[currentDetElemLocal]->GetEntries()]) AliMFTCluster(*newCluster);
211 } // end of cycle over the digits
213 } // no more digits to check in current plane!
217 printf("Plane %d: clusters found in %f seconds\n",iPlane,sw->CpuTime());
220 // Now we merge the cluster lists coming from each detection element, to build the cluster list of the plane
222 AliMFTCluster *newCluster = NULL;
223 for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) {
224 for (Int_t iCluster=0; iCluster<clustersPerDetElem[iDetElem]->GetEntries(); iCluster++) {
225 newCluster = (AliMFTCluster*) (clustersPerDetElem[iDetElem]->At(iCluster));
226 newCluster -> TerminateCluster();
227 new ((*fClustersPerPlane[iPlane])[fClustersPerPlane[iPlane]->GetEntries()]) AliMFTCluster(*newCluster);
231 printf("%d Clusters in plane %02d merged in %f seconds\n", fClustersPerPlane[iPlane]->GetEntries(), iPlane, sw->CpuTime());
233 for (Int_t iDetElem=0; iDetElem<nDetElem; iDetElem++) {
234 clustersPerDetElem[iDetElem] -> Delete();
235 delete clustersPerDetElem[iDetElem];
238 myDigitList -> Delete();
240 } // end of cycle over the planes
244 //====================================================================================================================================================
246 void AliMFTClusterFinder::MakeClusterBranch(TTree *treeCluster) {
248 // Creating the cluster branches, one for each plane (see AliMFTReconstructor::Reconstruct)
250 AliDebug(1, "Making Cluster Branch");
255 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
256 AliDebug(1, Form("Setting Branch Plane_%02d for treeCluster",iPlane));
257 if (treeCluster->GetBranch(Form("Plane_%02d",iPlane))) continue;
258 AliDebug(1, Form("Branch Plane_%02d does not exist, creating!",iPlane));
259 treeCluster->Branch(Form("Plane_%02d",iPlane), &(fClustersPerPlane[iPlane]));
265 //====================================================================================================================================================
267 void AliMFTClusterFinder::SetClusterTreeAddress(TTree *treeCluster) {
269 // Addressing the cluster branches, one for each plane (see AliMFTReconstructor::Reconstruct)
271 if (treeCluster && treeCluster->GetBranch("Plane_00")) {
273 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
274 if (treeCluster->GetBranch(Form("Plane_%02d",iPlane))) {
275 treeCluster->SetBranchAddress(Form("Plane_%02d",iPlane), &(fClustersPerPlane[iPlane]));
277 else AliError(Form("No branch available with name Plane_%02d", iPlane));
283 //====================================================================================================================================================
285 void AliMFTClusterFinder::CreateClusters() {
287 // create cluster list
289 AliDebug(1, Form("Creating clusters list: nPlanes = %d",fNPlanes));
291 if (fClustersPerPlane[0]) return;
293 for(Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
294 AliDebug(1, Form("plane %02d", iPlane));
295 fClustersPerPlane[iPlane] = new TClonesArray("AliMFTCluster");
296 fClustersPerPlane[iPlane] -> SetOwner(kTRUE);
302 //====================================================================================================================================================