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 // Date : August 05 2003 //
19 // This reads the file PMD.digits.root(TreeD), //
20 // calls the Clustering algorithm and stores the //
21 // clustering output in PMD.RecPoints.root(TreeR) //
23 //-----------------------------------------------------//
25 #include <Riostream.h>
27 #include <TObjArray.h>
28 #include <TClonesArray.h>
31 #include "AliRunLoader.h"
32 #include "AliLoader.h"
33 #include "AliRawReader.h"
35 #include "AliPMDdigit.h"
36 #include "AliPMDClusterFinder.h"
37 #include "AliPMDClustering.h"
38 #include "AliPMDClusteringV1.h"
39 #include "AliPMDcluster.h"
40 #include "AliPMDrecpoint1.h"
41 #include "AliPMDrechit.h"
42 #include "AliPMDRawStream.h"
43 #include "AliPMDCalibData.h"
44 #include "AliPMDddldata.h"
46 #include "AliCDBManager.h"
47 #include "AliCDBEntry.h"
51 ClassImp(AliPMDClusterFinder)
53 AliPMDClusterFinder::AliPMDClusterFinder():
58 fDigits(new TClonesArray("AliPMDdigit", 1000)),
59 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
60 fRechits(new TClonesArray("AliPMDrechit", 1000)),
68 fCalibData = GetCalibData();
70 // ------------------------------------------------------------------------- //
71 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
72 fRunLoader(runLoader),
73 fPMDLoader(runLoader->GetLoader("PMDLoader")),
76 fDigits(new TClonesArray("AliPMDdigit", 1000)),
77 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
78 fRechits(new TClonesArray("AliPMDrechit", 1000)),
86 fCalibData = GetCalibData();
88 // ------------------------------------------------------------------------- //
89 AliPMDClusterFinder::~AliPMDClusterFinder()
100 fRecpoints->Delete();
111 // ------------------------------------------------------------------------- //
113 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
115 // Converts digits to recpoints after running clustering
116 // algorithm on CPV plane and PREshower plane
118 Int_t det = 0,smn = 0;
125 TObjArray *pmdcont = new TObjArray();
126 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
128 pmdclust->SetEdepCut(fEcut);
130 fRunLoader->GetEvent(ievt);
133 fTreeD = fPMDLoader->TreeD();
136 AliFatal("AliPMDClusterFinder: Can not get TreeD");
139 AliPMDdigit *pmddigit;
140 TBranch *branch = fTreeD->GetBranch("PMDDigit");
141 branch->SetAddress(&fDigits);
145 fTreeR = fPMDLoader->TreeR();
148 fPMDLoader->MakeTree("R");
149 fTreeR = fPMDLoader->TreeR();
152 Int_t bufsize = 16000;
153 TBranch * branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
154 TBranch * branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
156 Int_t nmodules = (Int_t) fTreeD->GetEntries();
158 for (Int_t imodule = 0; imodule < nmodules; imodule++)
161 fTreeD->GetEntry(imodule);
162 Int_t nentries = fDigits->GetLast();
163 for (Int_t ient = 0; ient < nentries+1; ient++)
165 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
167 det = pmddigit->GetDetector();
168 smn = pmddigit->GetSMNumber();
169 xpos = pmddigit->GetRow();
170 ypos = pmddigit->GetColumn();
171 adc = pmddigit->GetADC();
174 Float_t gain = fCalibData->GetGainFact(det,smn,xpos,ypos);
176 // printf("adc = %d gain = %f\n",adc,gain);
181 //Int_t trno = pmddigit->GetTrackNumber();
182 fCellADC[xpos][ypos] = (Double_t) adc;
187 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
189 Int_t nentries1 = pmdcont->GetEntries();
191 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
193 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
195 AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
196 idet = pmdcl->GetDetector();
197 ismn = pmdcl->GetSMN();
198 clusdata[0] = pmdcl->GetClusX();
199 clusdata[1] = pmdcl->GetClusY();
200 clusdata[2] = pmdcl->GetClusADC();
201 clusdata[3] = pmdcl->GetClusCells();
202 clusdata[4] = pmdcl->GetClusSigmaX();
203 clusdata[5] = pmdcl->GetClusSigmaY();
205 AddRecPoint(idet,ismn,clusdata);
207 Int_t ncell = (Int_t) clusdata[3];
208 for(Int_t ihit = 0; ihit < ncell; ihit++)
210 Int_t celldataX = pmdcl->GetClusCellX(ihit);
211 Int_t celldataY = pmdcl->GetClusCellY(ihit);
212 AddRecHit(celldataX, celldataY);
225 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
226 fPMDLoader->WriteRecPoints("OVERWRITE");
228 // delete the pointers
233 // ------------------------------------------------------------------------- //
235 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
238 // Converts RAW data to recpoints after running clustering
239 // algorithm on CPV and PREshower plane
243 TObjArray pmdddlcont;
245 TObjArray *pmdcont = new TObjArray();
246 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
248 pmdclust->SetEdepCut(fEcut);
252 Int_t bufsize = 16000;
253 TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
255 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
257 const Int_t kDDL = 6;
258 const Int_t kRow = 48;
259 const Int_t kCol = 96;
264 for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
270 else if (indexDDL >= 4)
275 precpvADC = new int **[iSMN];
276 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
277 for (Int_t i=0; i<iSMN;i++)
279 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
281 for (Int_t i = 0; i < iSMN; i++)
283 for (Int_t j = 0; j < kRow; j++)
285 for (Int_t k = 0; k < kCol; k++)
287 precpvADC[i][j][k] = 0;
293 AliPMDRawStream pmdinput(rawReader);
294 rawReader->Select("PMD", indexDDL, indexDDL);
295 pmdinput.DdlData(&pmdddlcont);
297 Int_t ientries = pmdddlcont.GetEntries();
298 for (Int_t ient = 0; ient < ientries; ient++)
300 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
302 Int_t det = pmdddl->GetDetector();
303 Int_t smn = pmdddl->GetSMN();
304 //Int_t mcm = pmdddl->GetMCM();
305 //Int_t chno = pmdddl->GetChannel();
306 Int_t row = pmdddl->GetRow();
307 Int_t col = pmdddl->GetColumn();
308 Int_t sig = pmdddl->GetSignal();
310 Float_t sig1 = (Float_t) sig;
312 Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
313 //printf("sig = %d gain = %f\n",sig,gain);
314 sig = (Int_t) (sig1*gain);
319 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
321 indexsmn = smn - indexDDL * 6;
323 else if (indexDDL == 4)
326 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
332 else if (smn >= 18 && smn < 24)
337 else if (indexDDL == 5)
340 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
342 if (smn >= 6 && smn < 12)
346 else if (smn >= 12 && smn < 18)
351 precpvADC[indexsmn][row][col] = sig;
357 for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
360 for (Int_t irow = 0; irow < kRow; irow++)
362 for (Int_t icol = 0; icol < kCol; icol++)
364 fCellADC[irow][icol] =
365 (Double_t) precpvADC[indexsmn][irow][icol];
371 ismn = indexsmn + indexDDL * 6;
374 else if (indexDDL == 4)
380 else if (indexsmn >= 6 && indexsmn < 12)
382 ismn = indexsmn + 12;
386 else if (indexDDL == 5)
392 else if (indexsmn >= 6 && indexsmn < 12)
399 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
400 Int_t nentries1 = pmdcont->GetEntries();
402 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
404 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
406 AliPMDcluster *pmdcl =
407 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
408 idet = pmdcl->GetDetector();
409 ismn = pmdcl->GetSMN();
410 clusdata[0] = pmdcl->GetClusX();
411 clusdata[1] = pmdcl->GetClusY();
412 clusdata[2] = pmdcl->GetClusADC();
413 clusdata[3] = pmdcl->GetClusCells();
414 clusdata[4] = pmdcl->GetClusSigmaX();
415 clusdata[5] = pmdcl->GetClusSigmaY();
417 AddRecPoint(idet,ismn,clusdata);
419 Int_t ncell = (Int_t) clusdata[3];
420 for(Int_t ihit = 0; ihit < ncell; ihit++)
422 Int_t celldataX = pmdcl->GetClusCellX(ihit);
423 Int_t celldataY = pmdcl->GetClusCellY(ihit);
424 AddRecHit(celldataX, celldataY);
438 for (Int_t i=0; i<iSMN; i++)
440 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
442 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
448 // delete the pointers
453 // ------------------------------------------------------------------------- //
455 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
457 // Converts RAW data to recpoints after running clustering
458 // algorithm on CPV and PREshower plane
462 TObjArray pmdddlcont;
463 TObjArray *pmdcont = new TObjArray();
465 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
467 pmdclust->SetEdepCut(fEcut);
469 fRunLoader->GetEvent(ievt);
473 fTreeR = fPMDLoader->TreeR();
476 fPMDLoader->MakeTree("R");
477 fTreeR = fPMDLoader->TreeR();
479 Int_t bufsize = 16000;
480 TBranch *branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
481 TBranch *branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
483 const Int_t kDDL = 6;
484 const Int_t kRow = 48;
485 const Int_t kCol = 96;
490 for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
496 else if (indexDDL >= 4)
501 precpvADC = new int **[iSMN];
502 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
503 for (Int_t i=0; i<iSMN;i++)
505 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
507 for (Int_t i = 0; i < iSMN; i++)
509 for (Int_t j = 0; j < kRow; j++)
511 for (Int_t k = 0; k < kCol; k++)
513 precpvADC[i][j][k] = 0;
519 AliPMDRawStream pmdinput(rawReader);
520 rawReader->Select("PMD", indexDDL, indexDDL);
522 pmdinput.DdlData(&pmdddlcont);
525 Int_t ientries = pmdddlcont.GetEntries();
526 for (Int_t ient = 0; ient < ientries; ient++)
528 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
530 Int_t det = pmdddl->GetDetector();
531 Int_t smn = pmdddl->GetSMN();
532 //Int_t mcm = pmdddl->GetMCM();
533 //Int_t chno = pmdddl->GetChannel();
534 Int_t row = pmdddl->GetRow();
535 Int_t col = pmdddl->GetColumn();
536 Int_t sig = pmdddl->GetSignal();
538 Float_t sig1 = (Float_t) sig;
540 Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
541 //printf("sig = %d gain = %f\n",sig,gain);
542 sig = (Int_t) (sig1*gain);
548 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
550 indexsmn = smn - indexDDL * 6;
552 else if (indexDDL == 4)
555 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
561 else if (smn >= 18 && smn < 24)
566 else if (indexDDL == 5)
569 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
571 if (smn >= 6 && smn < 12)
575 else if (smn >= 12 && smn < 18)
580 precpvADC[indexsmn][row][col] = sig;
587 for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
590 for (Int_t irow = 0; irow < kRow; irow++)
592 for (Int_t icol = 0; icol < kCol; icol++)
594 fCellADC[irow][icol] =
595 (Double_t) precpvADC[indexsmn][irow][icol];
602 ismn = indexsmn + indexDDL * 6;
605 else if (indexDDL == 4)
611 else if (indexsmn >= 6 && indexsmn < 12)
613 ismn = indexsmn + 12;
617 else if (indexDDL == 5)
623 else if (indexsmn >= 6 && indexsmn < 12)
630 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
631 Int_t nentries1 = pmdcont->GetEntries();
633 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
635 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
637 AliPMDcluster *pmdcl =
638 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
639 idet = pmdcl->GetDetector();
640 ismn = pmdcl->GetSMN();
641 clusdata[0] = pmdcl->GetClusX();
642 clusdata[1] = pmdcl->GetClusY();
643 clusdata[2] = pmdcl->GetClusADC();
644 clusdata[3] = pmdcl->GetClusCells();
645 clusdata[4] = pmdcl->GetClusSigmaX();
646 clusdata[5] = pmdcl->GetClusSigmaY();
648 AddRecPoint(idet,ismn,clusdata);
650 Int_t ncell = (Int_t) clusdata[3];
651 for(Int_t ihit = 0; ihit < ncell; ihit++)
653 Int_t celldataX = pmdcl->GetClusCellX(ihit);
654 Int_t celldataY = pmdcl->GetClusCellY(ihit);
655 AddRecHit(celldataX, celldataY);
669 for (Int_t i=0; i<iSMN; i++)
671 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
673 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
679 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
680 fPMDLoader->WriteRecPoints("OVERWRITE");
682 // delete the pointers
687 // ------------------------------------------------------------------------- //
688 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
692 // ------------------------------------------------------------------------- //
693 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
695 // Add Reconstructed points
697 TClonesArray &lrecpoints = *fRecpoints;
698 AliPMDrecpoint1 *newrecpoint;
699 newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
700 new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
703 // ------------------------------------------------------------------------- //
704 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY)
706 // Add associated cell hits to the Reconstructed points
708 TClonesArray &lrechits = *fRechits;
709 AliPMDrechit *newrechit;
710 newrechit = new AliPMDrechit(celldataX, celldataY);
711 new(lrechits[fNhit++]) AliPMDrechit(newrechit);
714 // ------------------------------------------------------------------------- //
715 void AliPMDClusterFinder::ResetCellADC()
717 // Reset the individual cell ADC value to zero
719 for(Int_t irow = 0; irow < fgkRow; irow++)
721 for(Int_t icol = 0; icol < fgkCol; icol++)
723 fCellADC[irow][icol] = 0.;
727 // ------------------------------------------------------------------------- //
729 void AliPMDClusterFinder::ResetRecpoint()
731 // Clear the list of reconstructed points
733 if (fRecpoints) fRecpoints->Clear();
735 // ------------------------------------------------------------------------- //
736 void AliPMDClusterFinder::ResetRechit()
738 // Clear the list of reconstructed points
740 if (fRechits) fRechits->Clear();
742 // ------------------------------------------------------------------------- //
743 void AliPMDClusterFinder::Load()
745 // Load all the *.root files
747 fPMDLoader->LoadDigits("READ");
748 fPMDLoader->LoadRecPoints("recreate");
750 // ------------------------------------------------------------------------- //
751 void AliPMDClusterFinder::LoadClusters()
753 // Load all the *.root files
755 fPMDLoader->LoadRecPoints("recreate");
757 // ------------------------------------------------------------------------- //
758 void AliPMDClusterFinder::UnLoad()
760 // Unload all the *.root files
762 fPMDLoader->UnloadDigits();
763 fPMDLoader->UnloadRecPoints();
765 // ------------------------------------------------------------------------- //
766 void AliPMDClusterFinder::UnLoadClusters()
768 // Unload all the *.root files
770 fPMDLoader->UnloadRecPoints();
772 // ------------------------------------------------------------------------- //
774 AliPMDCalibData* AliPMDClusterFinder::GetCalibData() const
776 // The run number will be centralized in AliCDBManager,
777 // you don't need to set it here!
779 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Data");
782 AliWarning("Calibration object retrieval failed! Dummy calibration will be used.");
784 // this just remembers the actual default storage. No problem if it is null.
785 AliCDBStorage *origStorage = AliCDBManager::Instance()->GetDefaultStorage();
786 AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
788 entry = AliCDBManager::Instance()->Get("PMD/Calib/Data");
790 // now reset the original default storage to AliCDBManager...
791 AliCDBManager::Instance()->SetDefaultStorage(origStorage);
794 AliPMDCalibData *calibdata=0;
795 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
797 if (!calibdata) AliError("No calibration data from calibration database !");