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"
47 #include "AliCDBManager.h"
48 #include "AliCDBEntry.h"
52 ClassImp(AliPMDClusterFinder)
54 AliPMDClusterFinder::AliPMDClusterFinder():
57 fCalibData(GetCalibData()),
60 fDigits(new TClonesArray("AliPMDdigit", 1000)),
61 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
62 fRechits(new TClonesArray("AliPMDrechit", 1000)),
72 // ------------------------------------------------------------------------- //
73 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
74 fRunLoader(runLoader),
75 fPMDLoader(runLoader->GetLoader("PMDLoader")),
76 fCalibData(GetCalibData()),
79 fDigits(new TClonesArray("AliPMDdigit", 1000)),
80 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
81 fRechits(new TClonesArray("AliPMDrechit", 1000)),
91 // ------------------------------------------------------------------------- //
92 AliPMDClusterFinder::AliPMDClusterFinder(const AliPMDClusterFinder & finder):
96 fCalibData(GetCalibData()),
108 AliError("Copy constructor not allowed");
110 // ------------------------------------------------------------------------- //
111 AliPMDClusterFinder &AliPMDClusterFinder::operator=(const AliPMDClusterFinder & /*finder*/)
114 AliError("Assignment Operator not allowed");
117 // ------------------------------------------------------------------------- //
118 AliPMDClusterFinder::~AliPMDClusterFinder()
129 fRecpoints->Delete();
140 // ------------------------------------------------------------------------- //
142 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
144 // Converts digits to recpoints after running clustering
145 // algorithm on CPV plane and PREshower plane
148 Int_t det = 0,smn = 0;
155 TObjArray *pmdcont = new TObjArray();
156 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
158 pmdclust->SetEdepCut(fEcut);
160 fRunLoader->GetEvent(ievt);
163 fTreeD = fPMDLoader->TreeD();
166 AliFatal("AliPMDClusterFinder: Can not get TreeD");
169 AliPMDdigit *pmddigit;
170 TBranch *branch = fTreeD->GetBranch("PMDDigit");
171 branch->SetAddress(&fDigits);
175 fTreeR = fPMDLoader->TreeR();
178 fPMDLoader->MakeTree("R");
179 fTreeR = fPMDLoader->TreeR();
182 Int_t bufsize = 16000;
183 TBranch * branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
184 TBranch * branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
186 Int_t nmodules = (Int_t) fTreeD->GetEntries();
188 for (Int_t imodule = 0; imodule < nmodules; imodule++)
191 fTreeD->GetEntry(imodule);
192 Int_t nentries = fDigits->GetLast();
193 for (Int_t ient = 0; ient < nentries+1; ient++)
195 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
197 det = pmddigit->GetDetector();
198 smn = pmddigit->GetSMNumber();
199 xpos = pmddigit->GetRow();
200 ypos = pmddigit->GetColumn();
201 adc = pmddigit->GetADC();
204 Float_t gain = fCalibData->GetGainFact(det,smn,xpos,ypos);
205 // printf("adc = %d gain = %f\n",adc,gain);
209 //Int_t trno = pmddigit->GetTrackNumber();
210 fCellADC[xpos][ypos] = (Double_t) adc;
215 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
217 Int_t nentries1 = pmdcont->GetEntries();
219 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
221 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
223 AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
224 idet = pmdcl->GetDetector();
225 ismn = pmdcl->GetSMN();
226 clusdata[0] = pmdcl->GetClusX();
227 clusdata[1] = pmdcl->GetClusY();
228 clusdata[2] = pmdcl->GetClusADC();
229 clusdata[3] = pmdcl->GetClusCells();
230 clusdata[4] = pmdcl->GetClusSigmaX();
231 clusdata[5] = pmdcl->GetClusSigmaY();
233 AddRecPoint(idet,ismn,clusdata);
235 Int_t ncell = (Int_t) clusdata[3];
236 for(Int_t ihit = 0; ihit < ncell; ihit++)
238 Int_t celldataX = pmdcl->GetClusCellX(ihit);
239 Int_t celldataY = pmdcl->GetClusCellY(ihit);
240 AddRecHit(celldataX, celldataY);
253 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
254 fPMDLoader->WriteRecPoints("OVERWRITE");
256 // delete the pointers
261 // ------------------------------------------------------------------------- //
263 void AliPMDClusterFinder::Digits2RecPoints(TTree *digitsTree,
266 // Converts digits to recpoints after running clustering
267 // algorithm on CPV plane and PREshower plane
270 Int_t det = 0,smn = 0;
277 TObjArray *pmdcont = new TObjArray();
278 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
280 pmdclust->SetEdepCut(fEcut);
282 AliPMDdigit *pmddigit;
283 TBranch *branch = digitsTree->GetBranch("PMDDigit");
284 branch->SetAddress(&fDigits);
288 Int_t bufsize = 16000;
289 TBranch * branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
290 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
292 Int_t nmodules = (Int_t) digitsTree->GetEntries();
294 for (Int_t imodule = 0; imodule < nmodules; imodule++)
297 digitsTree->GetEntry(imodule);
298 Int_t nentries = fDigits->GetLast();
299 for (Int_t ient = 0; ient < nentries+1; ient++)
301 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
303 det = pmddigit->GetDetector();
304 smn = pmddigit->GetSMNumber();
305 xpos = pmddigit->GetRow();
306 ypos = pmddigit->GetColumn();
307 adc = pmddigit->GetADC();
310 Float_t gain = fCalibData->GetGainFact(det,smn,xpos,ypos);
311 // printf("adc = %d gain = %f\n",adc,gain);
315 //Int_t trno = pmddigit->GetTrackNumber();
316 fCellADC[xpos][ypos] = (Double_t) adc;
321 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
323 Int_t nentries1 = pmdcont->GetEntries();
325 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
327 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
329 AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
330 idet = pmdcl->GetDetector();
331 ismn = pmdcl->GetSMN();
332 clusdata[0] = pmdcl->GetClusX();
333 clusdata[1] = pmdcl->GetClusY();
334 clusdata[2] = pmdcl->GetClusADC();
335 clusdata[3] = pmdcl->GetClusCells();
336 clusdata[4] = pmdcl->GetClusSigmaX();
337 clusdata[5] = pmdcl->GetClusSigmaY();
339 AddRecPoint(idet,ismn,clusdata);
341 Int_t ncell = (Int_t) clusdata[3];
342 for(Int_t ihit = 0; ihit < ncell; ihit++)
344 Int_t celldataX = pmdcl->GetClusCellX(ihit);
345 Int_t celldataY = pmdcl->GetClusCellY(ihit);
346 AddRecHit(celldataX, celldataY);
360 // delete the pointers
365 // ------------------------------------------------------------------------- //
367 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
370 // Converts RAW data to recpoints after running clustering
371 // algorithm on CPV and PREshower plane
373 // This method is called at the time of reconstruction
377 TObjArray pmdddlcont;
379 TObjArray *pmdcont = new TObjArray();
380 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
382 pmdclust->SetEdepCut(fEcut);
386 Int_t bufsize = 16000;
387 TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
389 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
391 const Int_t kDDL = AliDAQ::NumberOfDdls("PMD");
392 const Int_t kRow = 48;
393 const Int_t kCol = 96;
399 for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
405 else if (indexDDL >= 4)
410 precpvADC = new int **[iSMN];
411 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
412 for (Int_t i=0; i<iSMN;i++)
414 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
416 for (Int_t i = 0; i < iSMN; i++)
418 for (Int_t j = 0; j < kRow; j++)
420 for (Int_t k = 0; k < kCol; k++)
422 precpvADC[i][j][k] = 0;
428 AliPMDRawStream pmdinput(rawReader);
430 rawReader->Select("PMD", indexDDL, indexDDL);
432 pmdinput.DdlData(indexDDL,&pmdddlcont);
435 Int_t ientries = pmdddlcont.GetEntries();
436 for (Int_t ient = 0; ient < ientries; ient++)
438 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
440 Int_t det = pmdddl->GetDetector();
441 Int_t smn = pmdddl->GetSMN();
442 //Int_t mcm = pmdddl->GetMCM();
443 //Int_t chno = pmdddl->GetChannel();
444 Int_t row = pmdddl->GetRow();
445 Int_t col = pmdddl->GetColumn();
446 Int_t sig = pmdddl->GetSignal();
448 Float_t sig1 = (Float_t) sig;
450 Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
451 //printf("sig = %d gain = %f\n",sig,gain);
452 sig = (Int_t) (sig1*gain);
457 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
459 indexsmn = smn - indexDDL * 6;
461 else if (indexDDL == 4)
464 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
470 else if (smn >= 18 && smn < 24)
475 else if (indexDDL == 5)
478 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
480 if (smn >= 6 && smn < 12)
484 else if (smn >= 12 && smn < 18)
489 precpvADC[indexsmn][row][col] = sig;
495 for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
498 for (Int_t irow = 0; irow < kRow; irow++)
500 for (Int_t icol = 0; icol < kCol; icol++)
502 fCellADC[irow][icol] =
503 (Double_t) precpvADC[indexsmn][irow][icol];
509 ismn = indexsmn + indexDDL * 6;
512 else if (indexDDL == 4)
518 else if (indexsmn >= 6 && indexsmn < 12)
520 ismn = indexsmn + 12;
524 else if (indexDDL == 5)
530 else if (indexsmn >= 6 && indexsmn < 12)
537 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
538 Int_t nentries1 = pmdcont->GetEntries();
540 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
542 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
544 AliPMDcluster *pmdcl =
545 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
546 idet = pmdcl->GetDetector();
547 ismn = pmdcl->GetSMN();
548 clusdata[0] = pmdcl->GetClusX();
549 clusdata[1] = pmdcl->GetClusY();
550 clusdata[2] = pmdcl->GetClusADC();
551 clusdata[3] = pmdcl->GetClusCells();
552 clusdata[4] = pmdcl->GetClusSigmaX();
553 clusdata[5] = pmdcl->GetClusSigmaY();
555 AddRecPoint(idet,ismn,clusdata);
557 Int_t ncell = (Int_t) clusdata[3];
558 for(Int_t ihit = 0; ihit < ncell; ihit++)
560 Int_t celldataX = pmdcl->GetClusCellX(ihit);
561 Int_t celldataY = pmdcl->GetClusCellY(ihit);
562 AddRecHit(celldataX, celldataY);
576 for (Int_t i=0; i<iSMN; i++)
578 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
580 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
586 // delete the pointers
591 // ------------------------------------------------------------------------- //
593 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
595 // Converts RAW data to recpoints after running clustering
596 // algorithm on CPV and PREshower plane
600 TObjArray pmdddlcont;
601 TObjArray *pmdcont = new TObjArray();
603 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
605 pmdclust->SetEdepCut(fEcut);
607 fRunLoader->GetEvent(ievt);
611 fTreeR = fPMDLoader->TreeR();
614 fPMDLoader->MakeTree("R");
615 fTreeR = fPMDLoader->TreeR();
617 Int_t bufsize = 16000;
618 TBranch *branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
619 TBranch *branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
621 const Int_t kDDL = AliDAQ::NumberOfDdls("PMD");
622 const Int_t kRow = 48;
623 const Int_t kCol = 96;
628 for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
634 else if (indexDDL >= 4)
639 precpvADC = new int **[iSMN];
640 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
641 for (Int_t i=0; i<iSMN;i++)
643 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
645 for (Int_t i = 0; i < iSMN; i++)
647 for (Int_t j = 0; j < kRow; j++)
649 for (Int_t k = 0; k < kCol; k++)
651 precpvADC[i][j][k] = 0;
657 rawReader->Select("PMD", indexDDL, indexDDL);
659 AliPMDRawStream pmdinput(rawReader);
660 pmdinput.DdlData(indexDDL,&pmdddlcont);
663 Int_t ientries = pmdddlcont.GetEntries();
664 for (Int_t ient = 0; ient < ientries; ient++)
666 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
668 Int_t det = pmdddl->GetDetector();
669 Int_t smn = pmdddl->GetSMN();
670 //Int_t mcm = pmdddl->GetMCM();
671 //Int_t chno = pmdddl->GetChannel();
672 Int_t row = pmdddl->GetRow();
673 Int_t col = pmdddl->GetColumn();
674 Int_t sig = pmdddl->GetSignal();
676 Float_t sig1 = (Float_t) sig;
678 Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
680 //printf("sig = %d gain = %f\n",sig,gain);
681 sig = (Int_t) (sig1*gain);
687 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
689 indexsmn = smn - indexDDL * 6;
691 else if (indexDDL == 4)
694 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
700 else if (smn >= 18 && smn < 24)
705 else if (indexDDL == 5)
708 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
710 if (smn >= 6 && smn < 12)
714 else if (smn >= 12 && smn < 18)
719 precpvADC[indexsmn][row][col] = sig;
726 for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
729 for (Int_t irow = 0; irow < kRow; irow++)
731 for (Int_t icol = 0; icol < kCol; icol++)
733 fCellADC[irow][icol] =
734 (Double_t) precpvADC[indexsmn][irow][icol];
741 ismn = indexsmn + indexDDL * 6;
744 else if (indexDDL == 4)
750 else if (indexsmn >= 6 && indexsmn < 12)
752 ismn = indexsmn + 12;
756 else if (indexDDL == 5)
762 else if (indexsmn >= 6 && indexsmn < 12)
769 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
770 Int_t nentries1 = pmdcont->GetEntries();
772 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
774 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
776 AliPMDcluster *pmdcl =
777 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
778 idet = pmdcl->GetDetector();
779 ismn = pmdcl->GetSMN();
780 clusdata[0] = pmdcl->GetClusX();
781 clusdata[1] = pmdcl->GetClusY();
782 clusdata[2] = pmdcl->GetClusADC();
783 clusdata[3] = pmdcl->GetClusCells();
784 clusdata[4] = pmdcl->GetClusSigmaX();
785 clusdata[5] = pmdcl->GetClusSigmaY();
787 AddRecPoint(idet,ismn,clusdata);
789 Int_t ncell = (Int_t) clusdata[3];
790 for(Int_t ihit = 0; ihit < ncell; ihit++)
792 Int_t celldataX = pmdcl->GetClusCellX(ihit);
793 Int_t celldataY = pmdcl->GetClusCellY(ihit);
794 AddRecHit(celldataX, celldataY);
808 for (Int_t i=0; i<iSMN; i++)
810 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
812 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
818 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
819 fPMDLoader->WriteRecPoints("OVERWRITE");
821 // delete the pointers
826 // ------------------------------------------------------------------------- //
827 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
831 // ------------------------------------------------------------------------- //
832 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
834 // Add Reconstructed points
836 TClonesArray &lrecpoints = *fRecpoints;
837 AliPMDrecpoint1 *newrecpoint;
838 newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
839 new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
842 // ------------------------------------------------------------------------- //
843 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY)
845 // Add associated cell hits to the Reconstructed points
847 TClonesArray &lrechits = *fRechits;
848 AliPMDrechit *newrechit;
849 newrechit = new AliPMDrechit(celldataX, celldataY);
850 new(lrechits[fNhit++]) AliPMDrechit(newrechit);
853 // ------------------------------------------------------------------------- //
854 void AliPMDClusterFinder::ResetCellADC()
856 // Reset the individual cell ADC value to zero
858 for(Int_t irow = 0; irow < fgkRow; irow++)
860 for(Int_t icol = 0; icol < fgkCol; icol++)
862 fCellADC[irow][icol] = 0.;
866 // ------------------------------------------------------------------------- //
868 void AliPMDClusterFinder::ResetRecpoint()
870 // Clear the list of reconstructed points
872 if (fRecpoints) fRecpoints->Clear();
874 // ------------------------------------------------------------------------- //
875 void AliPMDClusterFinder::ResetRechit()
877 // Clear the list of reconstructed points
879 if (fRechits) fRechits->Clear();
881 // ------------------------------------------------------------------------- //
882 void AliPMDClusterFinder::Load()
884 // Load all the *.root files
886 fPMDLoader->LoadDigits("READ");
887 fPMDLoader->LoadRecPoints("recreate");
889 // ------------------------------------------------------------------------- //
890 void AliPMDClusterFinder::LoadClusters()
892 // Load all the *.root files
894 fPMDLoader->LoadRecPoints("recreate");
896 // ------------------------------------------------------------------------- //
897 void AliPMDClusterFinder::UnLoad()
899 // Unload all the *.root files
901 fPMDLoader->UnloadDigits();
902 fPMDLoader->UnloadRecPoints();
904 // ------------------------------------------------------------------------- //
905 void AliPMDClusterFinder::UnLoadClusters()
907 // Unload all the *.root files
909 fPMDLoader->UnloadRecPoints();
911 // ------------------------------------------------------------------------- //
913 AliPMDCalibData* AliPMDClusterFinder::GetCalibData() const
915 // The run number will be centralized in AliCDBManager,
916 // you don't need to set it here!
918 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Data");
920 if(!entry) AliFatal("Calibration object retrieval failed! ");
922 AliPMDCalibData *calibdata=0;
923 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
925 if (!calibdata) AliFatal("No calibration data from calibration database !");