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>
32 #include "AliRunLoader.h"
33 #include "AliLoader.h"
34 #include "AliRawReader.h"
36 #include "AliPMDdigit.h"
37 #include "AliPMDClusterFinder.h"
38 #include "AliPMDClustering.h"
39 #include "AliPMDClusteringV1.h"
40 #include "AliPMDcluster.h"
41 #include "AliPMDrecpoint1.h"
42 #include "AliPMDrechit.h"
43 #include "AliPMDRawStream.h"
44 #include "AliPMDCalibData.h"
45 #include "AliPMDPedestal.h"
46 #include "AliPMDddldata.h"
47 #include "AliPMDHotData.h"
48 #include "AliPMDNoiseCut.h"
49 #include "AliPMDddlinfoData.h"
50 #include "AliPMDRecoParam.h"
51 #include "AliRecoParam.h"
52 #include "AliPMDReconstructor.h"
55 #include "AliCDBManager.h"
56 #include "AliCDBEntry.h"
60 ClassImp(AliPMDClusterFinder)
62 AliPMDClusterFinder::AliPMDClusterFinder():
65 fCalibGain(GetCalibGain()),
66 fCalibPed(GetCalibPed()),
67 fCalibHot(GetCalibHot()),
68 fNoiseCut(GetNoiseCut()),
69 fDdlinfo(GetDdlinfoData()),
73 fDigits(new TClonesArray("AliPMDdigit", 1000)),
74 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
75 fRechits(new TClonesArray("AliPMDrechit", 1000)),
84 // ------------------------------------------------------------------------- //
85 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
86 fRunLoader(runLoader),
87 fPMDLoader(runLoader->GetLoader("PMDLoader")),
88 fCalibGain(GetCalibGain()),
89 fCalibPed(GetCalibPed()),
90 fCalibHot(GetCalibHot()),
91 fNoiseCut(GetNoiseCut()),
92 fDdlinfo(GetDdlinfoData()),
96 fDigits(new TClonesArray("AliPMDdigit", 1000)),
97 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
98 fRechits(new TClonesArray("AliPMDrechit", 1000)),
107 // ------------------------------------------------------------------------- //
108 AliPMDClusterFinder::AliPMDClusterFinder(const AliPMDClusterFinder & finder):
112 fCalibGain(GetCalibGain()),
113 fCalibPed(GetCalibPed()),
114 fCalibHot(GetCalibHot()),
115 fNoiseCut(GetNoiseCut()),
116 fDdlinfo(GetDdlinfoData()),
128 AliError("Copy constructor not allowed");
130 // ------------------------------------------------------------------------- //
131 AliPMDClusterFinder &AliPMDClusterFinder::operator=(const AliPMDClusterFinder & /*finder*/)
134 AliError("Assignment Operator not allowed");
137 // ------------------------------------------------------------------------- //
138 AliPMDClusterFinder::~AliPMDClusterFinder()
155 // ------------------------------------------------------------------------- //
157 void AliPMDClusterFinder::Digits2RecPoints(TTree *digitsTree,
158 TTree *clustersTree, Int_t gRecoMode)
160 // Converts digits to recpoints after running clustering
161 // algorithm on CPV plane and PREshower plane
163 // This algorithm is called during the reconstruction from digits
165 Int_t det = 0, smn = 0;
166 Int_t xpos = 0, ypos = 0;
170 Float_t clusdata[6] = {0.,0.,0.,0.,0.,0.};
172 AliPMDcluster *pmdcl = 0x0;
174 TObjArray *pmdcont = new TObjArray();
176 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
178 // Fetch the reco param object
180 fRecoParam = AliPMDReconstructor::GetRecoParam();
181 if(fRecoParam == 0x0)
183 AliFatal("No Reco Param found for PMD!!!");
187 AliPMDdigit *pmddigit;
188 TBranch *branch = digitsTree->GetBranch("PMDDigit");
189 branch->SetAddress(&fDigits);
193 Int_t bufsize = 16000;
194 TBranch * branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
195 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
197 Int_t nmodules = (Int_t) digitsTree->GetEntries();
199 for (Int_t imodule = 0; imodule < nmodules; imodule++)
204 digitsTree->GetEntry(imodule);
205 Int_t nentries = fDigits->GetLast();
206 for (Int_t ient = 0; ient < nentries+1; ient++)
208 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
210 det = pmddigit->GetDetector();
211 smn = pmddigit->GetSMNumber();
212 xpos = pmddigit->GetRow();
213 ypos = pmddigit->GetColumn();
214 adc = pmddigit->GetADC();
216 if(det < 0 || det > 1)
218 AliError(Form("*CPV/PRE NUMBER WRONG %d *",det));
221 if(smn == -1 || smn > 23)
223 AliError(Form("*MODULE NUMBER WRONG %d *",smn));
227 if(xpos < 0 || xpos > 47 || ypos < 0 || ypos > 95)
229 AliError(Form("*Row %d and Column NUMBER %d NOT Valid *",
234 // Pedestal Subtraction
235 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,xpos,ypos);
236 Int_t pedrms1 = (Int_t) pedmeanrms%100;
237 Float_t pedrms = (Float_t)pedrms1/10.;
238 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
239 //printf("%f %f\n",pedmean, pedrms);
241 Float_t adc1 = adc - (pedmean + 3.0*pedrms);
243 // Hot cell - set the cell adc = 0
244 Float_t hotflag = fCalibHot->GetHotChannel(det,smn,xpos,ypos);
245 if (hotflag == 1.) adc1 = 0;
248 Float_t gain = fCalibGain->GetGainFact(det,smn,xpos,ypos);
249 // printf("adc = %d gain = %f\n",adc,gain);
253 fCellTrack[xpos][ypos] = pmddigit->GetTrackNumber();
254 fCellPid[xpos][ypos] = pmddigit->GetTrackPid();
255 fCellADC[xpos][ypos] = (Double_t) adc;
257 totADCMod += (Int_t) adc;
264 if (totADCMod <= 0) continue;
266 // Set the minimum noise cut per module before clustering
268 // Int_t imod = idet*24 + ismn;
271 // Int_t cluspar = fRecoParam->GetPbPbParam()->GetClusteringParam();
272 AliPMDRecoParam * par = fRecoParam->GetPPParam();
273 Int_t cluspar = par->GetClusteringParam();
276 // Int_t cluspar = fRecoParam->GetCosmicParam()->GetClusteringParam();
278 //_______________________________________________________//
279 //Added to switch Refine and crude Clustering - satya//
280 // temporary solution - will be sorted out later
282 static AliPMDRecoParam *reconp = NULL;
283 reconp = (AliPMDRecoParam*)AliPMDReconstructor::GetRecoParam();
289 if( reconp->GetClusteringParam() == 1)
291 if( reconp->GetClusteringParam() == 2)
297 //_______________________________________________________//
299 pmdclust->SetClusteringParam(cluspar);
302 pmdclust->SetEdepCut(encut);
303 pmdclust->DoClust(idet,ismn,fCellTrack,fCellPid,fCellADC,pmdcont);
305 Int_t nentries1 = pmdcont->GetEntries();
307 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
309 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
311 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
312 idet = pmdcl->GetDetector();
313 ismn = pmdcl->GetSMN();
314 clusdata[0] = pmdcl->GetClusX();
315 clusdata[1] = pmdcl->GetClusY();
316 clusdata[2] = pmdcl->GetClusADC();
317 clusdata[3] = pmdcl->GetClusCells();
318 clusdata[4] = pmdcl->GetClusSigmaX();
319 clusdata[5] = pmdcl->GetClusSigmaY();
321 AddRecPoint(idet,ismn,clusdata);
323 Int_t ncell = (Int_t) clusdata[3];
324 if (ncell > 19) ncell = 19;
325 for(Int_t ihit = 0; ihit < ncell; ihit++)
327 Int_t celldataX = pmdcl->GetClusCellX(ihit);
328 Int_t celldataY = pmdcl->GetClusCellY(ihit);
329 Int_t celldataTr = pmdcl->GetClusCellTrack(ihit);
330 Int_t celldataPid = pmdcl->GetClusCellPid(ihit);
331 Float_t celldataAdc = pmdcl->GetClusCellAdc(ihit);
332 AddRecHit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
347 // delete the pointers
351 // ------------------------------------------------------------------------- //
353 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
354 TTree *clustersTree, Int_t gRecoMode)
356 // Converts RAW data to recpoints after running clustering
357 // algorithm on CPV and PREshower plane
359 // This method is called at the time of reconstruction from RAW data
362 AliPMDddldata *pmdddl = 0x0;
363 AliPMDcluster *pmdcl = 0x0;
366 TObjArray pmdddlcont;
368 TObjArray *pmdcont = new TObjArray();
370 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
372 // access the ddlinfo database to fetch the no of modules per DDL
374 Int_t moduleddl[6] = {0,0,0,0,0,0};
376 for(Int_t jddl = 0; jddl < 6; jddl++)
378 moduleddl[jddl] = fDdlinfo->GetNoOfModulePerDdl(jddl);
381 // Set the minimum noise cut per module before clustering
383 fRecoParam = AliPMDReconstructor::GetRecoParam();
385 if(fRecoParam == 0x0)
387 AliFatal("No Reco Param found for PMD!!!");
392 Int_t bufsize = 16000;
393 TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
395 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
397 const Int_t kRow = 48;
398 const Int_t kCol = 96;
404 AliPMDRawStream pmdinput(rawReader);
406 while ((indexDDL = pmdinput.DdlData(&pmdddlcont)) >=0)
408 iSMN = moduleddl[indexDDL];
411 precpvADC = new int **[iSMN];
412 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
413 for (Int_t i=0; i<iSMN;i++)
415 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
417 for (Int_t i = 0; i < iSMN; i++)
419 for (Int_t j = 0; j < kRow; j++)
421 for (Int_t k = 0; k < kCol; k++)
423 precpvADC[i][j][k] = 0;
430 Int_t ientries = pmdddlcont.GetEntries();
431 for (Int_t ient = 0; ient < ientries; ient++)
433 pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
435 Int_t det = pmdddl->GetDetector();
436 Int_t smn = pmdddl->GetSMN();
437 //Int_t mcm = pmdddl->GetMCM();
438 //Int_t chno = pmdddl->GetChannel();
439 Int_t row = pmdddl->GetRow();
440 Int_t col = pmdddl->GetColumn();
441 Int_t sig = pmdddl->GetSignal();
444 if(det < 0 || det > 1)
446 AliError(Form("*CPV/PRE NUMBER WRONG %d *",det));
449 if(smn < 0 || smn > 23)
451 AliError(Form("*MODULE NUMBER WRONG %d *",smn));
454 if(row < 0 || row > 47 || col < 0 || col > 95)
456 AliError(Form("*Row %d and Column NUMBER %d NOT Valid *",
462 // Pedestal Subtraction
463 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,row,col);
464 Int_t pedrms1 = (Int_t) pedmeanrms%100;
465 Float_t pedrms = (Float_t)pedrms1/10.;
466 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
468 //printf("%f %f\n",pedmean, pedrms);
470 // Float_t sig1 = (Float_t) sig;
471 Float_t sig1 = (Float_t) sig - (pedmean + 3.0*pedrms);
473 // Hot cell - set the cell adc = 0
474 Float_t hotflag = fCalibHot->GetHotChannel(det,smn,row,col);
475 if (hotflag == 1.) sig1 = 0;
478 Float_t gain = fCalibGain->GetGainFact(det,smn,row,col);
479 //printf("sig = %d gain = %f\n",sig,gain);
480 sig = (Int_t) (sig1*gain);
485 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
495 else if (smn >= 18 && smn < 24)
499 else if (indexDDL >= 1 && indexDDL < 4)
502 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
504 indexsmn = smn - indexDDL * 6;
506 else if (indexDDL == 4)
509 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
515 else if (smn >= 18 && smn < 24)
520 else if (indexDDL == 5)
523 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
525 if (smn >= 6 && smn < 18)
531 precpvADC[indexsmn][row][col] = sig;
539 for (indexsmn = 0; indexsmn < iSMN; indexsmn++)
543 for (Int_t irow = 0; irow < kRow; irow++)
545 for (Int_t icol = 0; icol < kCol; icol++)
547 fCellTrack[irow][icol] = -1;
548 fCellPid[irow][icol] = -1;
550 fCellADC[irow][icol] =
551 (Double_t) precpvADC[indexsmn][irow][icol];
552 totAdcMod += precpvADC[indexsmn][irow][icol];
567 else if (indexsmn >= 6 && indexsmn < 12)
568 ismn = indexsmn + 12;
572 else if (indexDDL >= 1 && indexDDL < 4)
574 ismn = indexsmn + indexDDL * 6;
577 else if (indexDDL == 4)
583 else if (indexsmn >= 6 && indexsmn < 12)
585 ismn = indexsmn + 12;
589 else if (indexDDL == 5)
595 if (totAdcMod <= 0) continue;
597 Int_t imod = idet*24 + ismn;
599 // Int_t cluspar = fRecoParam->GetPbPbParam()->GetClusteringParam();
600 AliPMDRecoParam * par = fRecoParam->GetPPParam();
601 Int_t cluspar = par->GetClusteringParam();
603 // Int_t cluspar = fRecoParam->GetCosmicParam()->GetClusteringParam();
605 //_______________________________________________________//
606 //Added to switch Refine and crude Clustering - satya//
607 // temporary solution - will be sorted out later
609 static AliPMDRecoParam *reconp = NULL;
610 reconp = (AliPMDRecoParam*)AliPMDReconstructor::GetRecoParam();
615 if( reconp->GetClusteringParam() == 1)
617 if( reconp->GetClusteringParam() == 2)
622 cluspar = gRecoMode; // permanent solution
624 //_______________________________________________________//
626 pmdclust->SetClusteringParam(cluspar);
627 Float_t encut = fNoiseCut->GetNoiseCut(imod);
629 pmdclust->SetEdepCut(encut);
630 pmdclust->DoClust(idet,ismn,fCellTrack,fCellPid,fCellADC,pmdcont);
632 Int_t nentries1 = pmdcont->GetEntries();
634 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
636 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
638 pmdcl = (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 if (ncell > 19) ncell = 19;
652 for(Int_t ihit = 0; ihit < ncell; ihit++)
654 Int_t celldataX = pmdcl->GetClusCellX(ihit);
655 Int_t celldataY = pmdcl->GetClusCellY(ihit);
656 Int_t celldataTr = pmdcl->GetClusCellTrack(ihit);
657 Int_t celldataPid = pmdcl->GetClusCellPid(ihit);
658 Float_t celldataAdc = pmdcl->GetClusCellAdc(ihit);
659 AddRecHit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
673 for (Int_t i=0; i<iSMN; i++)
675 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
677 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
685 // delete the pointers
689 // ------------------------------------------------------------------------- //
690 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
692 // Add Reconstructed points
694 TClonesArray &lrecpoints = *fRecpoints;
695 AliPMDrecpoint1 *newrecpoint;
696 newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
697 new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
700 // ------------------------------------------------------------------------- //
701 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY,
702 Int_t celldataTr, Int_t celldataPid,
705 // Add associated cell hits to the Reconstructed points
707 TClonesArray &lrechits = *fRechits;
708 AliPMDrechit *newrechit;
709 newrechit = new AliPMDrechit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
710 new(lrechits[fNhit++]) AliPMDrechit(newrechit);
713 // ------------------------------------------------------------------------- //
714 void AliPMDClusterFinder::ResetCellADC()
716 // Reset the individual cell ADC value to zero
718 for(Int_t irow = 0; irow < fgkRow; irow++)
720 for(Int_t icol = 0; icol < fgkCol; icol++)
722 fCellTrack[irow][icol] = -1;
723 fCellPid[irow][icol] = -1;
724 fCellADC[irow][icol] = 0.;
728 // ------------------------------------------------------------------------- //
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 // ------------------------------------------------------------------------- //
773 AliPMDCalibData* AliPMDClusterFinder::GetCalibGain() const
775 // The run number will be centralized in AliCDBManager,
776 // you don't need to set it here!
778 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
780 if(!entry) AliFatal("Calibration object retrieval failed! ");
782 AliPMDCalibData *calibdata=0;
783 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
785 if (!calibdata) AliFatal("No calibration data from calibration database !");
789 // ------------------------------------------------------------------------- //
790 AliPMDPedestal* AliPMDClusterFinder::GetCalibPed() const
792 // The run number will be centralized in AliCDBManager,
793 // you don't need to set it here!
794 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
796 if(!entry) AliFatal("Pedestal object retrieval failed!");
798 AliPMDPedestal *pedestal = 0;
799 if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
801 if (!pedestal) AliFatal("No pedestal data from pedestal database !");
805 //--------------------------------------------------------------------//
806 AliPMDHotData* AliPMDClusterFinder::GetCalibHot() const
808 // The run number will be centralized in AliCDBManager,
809 // you don't need to set it here!
810 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Hot");
812 if(!entry) AliFatal("HotData object retrieval failed!");
814 AliPMDHotData *hot = 0;
815 if (entry) hot = (AliPMDHotData*) entry->GetObject();
817 if (!hot) AliFatal("No hot data from database !");
821 //--------------------------------------------------------------------//
822 AliPMDNoiseCut* AliPMDClusterFinder::GetNoiseCut() const
824 // The run number will be centralized in AliCDBManager,
825 // you don't need to set it here!
826 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/NoiseCut");
828 if(!entry) AliFatal("Noisecut object retrieval failed!");
830 AliPMDNoiseCut *ncut = 0;
831 if (entry) ncut = (AliPMDNoiseCut*) entry->GetObject();
833 if (!ncut) AliFatal("No noise cut data from database !");
837 //--------------------------------------------------------------------//
838 AliPMDddlinfoData* AliPMDClusterFinder::GetDdlinfoData() const
840 // The run number will be centralized in AliCDBManager,
841 // you don't need to set it here!
842 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ddlinfo");
844 if(!entry) AliFatal("ddlinfo object retrieval failed!");
846 AliPMDddlinfoData *ddlinfo = 0;
847 if (entry) ddlinfo = (AliPMDddlinfoData*) entry->GetObject();
849 if (!ddlinfo) AliFatal("No ddl info data from database !");