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 "AliPMDRecoParam.h"
50 #include "AliPMDReconstructor.h"
53 #include "AliCDBManager.h"
54 #include "AliCDBEntry.h"
58 ClassImp(AliPMDClusterFinder)
60 AliPMDClusterFinder::AliPMDClusterFinder():
63 fCalibGain(GetCalibGain()),
64 fCalibPed(GetCalibPed()),
65 fCalibHot(GetCalibHot()),
66 fNoiseCut(GetNoiseCut()),
70 fDigits(new TClonesArray("AliPMDdigit", 1000)),
71 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
72 fRechits(new TClonesArray("AliPMDrechit", 1000)),
81 // ------------------------------------------------------------------------- //
82 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
83 fRunLoader(runLoader),
84 fPMDLoader(runLoader->GetLoader("PMDLoader")),
85 fCalibGain(GetCalibGain()),
86 fCalibPed(GetCalibPed()),
87 fCalibHot(GetCalibHot()),
88 fNoiseCut(GetNoiseCut()),
92 fDigits(new TClonesArray("AliPMDdigit", 1000)),
93 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
94 fRechits(new TClonesArray("AliPMDrechit", 1000)),
103 // ------------------------------------------------------------------------- //
104 AliPMDClusterFinder::AliPMDClusterFinder(const AliPMDClusterFinder & finder):
108 fCalibGain(GetCalibGain()),
109 fCalibPed(GetCalibPed()),
110 fCalibHot(GetCalibHot()),
111 fNoiseCut(GetNoiseCut()),
123 AliError("Copy constructor not allowed");
125 // ------------------------------------------------------------------------- //
126 AliPMDClusterFinder &AliPMDClusterFinder::operator=(const AliPMDClusterFinder & /*finder*/)
129 AliError("Assignment Operator not allowed");
132 // ------------------------------------------------------------------------- //
133 AliPMDClusterFinder::~AliPMDClusterFinder()
150 // ------------------------------------------------------------------------- //
152 void AliPMDClusterFinder::Digits2RecPoints(TTree *digitsTree,
155 // Converts digits to recpoints after running clustering
156 // algorithm on CPV plane and PREshower plane
158 // This algorithm is called during the reconstruction from digits
160 Int_t det = 0,smn = 0;
167 AliPMDcluster *pmdcl = 0x0;
169 TObjArray *pmdcont = new TObjArray();
171 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
173 // Fetch the reco param object
175 fRecoParam = AliPMDReconstructor::GetRecoParam();
176 if(fRecoParam == 0x0)
178 AliFatal("No Reco Param found for PMD!!!");
182 AliPMDdigit *pmddigit;
183 TBranch *branch = digitsTree->GetBranch("PMDDigit");
184 branch->SetAddress(&fDigits);
188 Int_t bufsize = 16000;
189 TBranch * branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
190 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
192 Int_t nmodules = (Int_t) digitsTree->GetEntries();
194 for (Int_t imodule = 0; imodule < nmodules; imodule++)
199 digitsTree->GetEntry(imodule);
200 Int_t nentries = fDigits->GetLast();
201 for (Int_t ient = 0; ient < nentries+1; ient++)
203 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
205 det = pmddigit->GetDetector();
206 smn = pmddigit->GetSMNumber();
207 xpos = pmddigit->GetRow();
208 ypos = pmddigit->GetColumn();
209 adc = pmddigit->GetADC();
211 if(det < 0 || det > 1)
213 AliError(Form("*CPV/PRE NUMBER WRONG %d *",det));
216 if(smn == -1 || smn > 23)
218 AliError(Form("*MODULE NUMBER WRONG %d *",smn));
222 if(xpos < 0 || xpos > 47 || ypos < 0 || ypos > 95)
224 AliError(Form("*Row %d and Column NUMBER %d NOT Valid *",
229 // Pedestal Subtraction
230 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,xpos,ypos);
231 Int_t pedrms1 = (Int_t) pedmeanrms%100;
232 Float_t pedrms = (Float_t)pedrms1/10.;
233 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
234 //printf("%f %f\n",pedmean, pedrms);
236 Float_t adc1 = adc - (pedmean + 3.0*pedrms);
238 // Hot cell - set the cell adc = 0
239 Float_t hotflag = fCalibHot->GetHotChannel(det,smn,xpos,ypos);
240 if (hotflag == 1.) adc1 = 0;
243 Float_t gain = fCalibGain->GetGainFact(det,smn,xpos,ypos);
244 // printf("adc = %d gain = %f\n",adc,gain);
248 fCellTrack[xpos][ypos] = pmddigit->GetTrackNumber();
249 fCellPid[xpos][ypos] = pmddigit->GetTrackPid();
250 fCellADC[xpos][ypos] = (Double_t) adc;
252 totADCMod += (Int_t) adc;
259 if (totADCMod <= 0) continue;
261 // Set the minimum noise cut per module before clustering
263 // Int_t imod = idet*24 + ismn;
266 // Int_t cluspar = fRecoParam->GetPbPbParam()->GetClusteringParam();
267 Int_t cluspar = fRecoParam->GetPPParam()->GetClusteringParam();
268 // Int_t cluspar = fRecoParam->GetCosmicParam()->GetClusteringParam();
270 pmdclust->SetClusteringParam(cluspar);
273 pmdclust->SetEdepCut(encut);
274 pmdclust->DoClust(idet,ismn,fCellTrack,fCellPid,fCellADC,pmdcont);
276 Int_t nentries1 = pmdcont->GetEntries();
278 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
280 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
282 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
283 idet = pmdcl->GetDetector();
284 ismn = pmdcl->GetSMN();
285 clusdata[0] = pmdcl->GetClusX();
286 clusdata[1] = pmdcl->GetClusY();
287 clusdata[2] = pmdcl->GetClusADC();
288 clusdata[3] = pmdcl->GetClusCells();
289 clusdata[4] = pmdcl->GetClusSigmaX();
290 clusdata[5] = pmdcl->GetClusSigmaY();
292 AddRecPoint(idet,ismn,clusdata);
294 Int_t ncell = (Int_t) clusdata[3];
295 if (ncell > 19) ncell = 19;
296 for(Int_t ihit = 0; ihit < ncell; ihit++)
298 Int_t celldataX = pmdcl->GetClusCellX(ihit);
299 Int_t celldataY = pmdcl->GetClusCellY(ihit);
300 Int_t celldataTr = pmdcl->GetClusCellTrack(ihit);
301 Int_t celldataPid = pmdcl->GetClusCellPid(ihit);
302 Float_t celldataAdc = pmdcl->GetClusCellAdc(ihit);
303 AddRecHit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
318 // delete the pointers
322 // ------------------------------------------------------------------------- //
324 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
327 // Converts RAW data to recpoints after running clustering
328 // algorithm on CPV and PREshower plane
330 // This method is called at the time of reconstruction from RAW data
333 AliPMDddldata *pmdddl = 0x0;
334 AliPMDcluster *pmdcl = 0x0;
337 TObjArray pmdddlcont;
339 TObjArray *pmdcont = new TObjArray();
341 AliPMDClustering *pmdclust = new AliPMDClusteringV1();
343 // open the ddl file info to know the module
344 TString ddlinfofileName(gSystem->Getenv("ALICE_ROOT"));
345 ddlinfofileName += "/PMD/PMD_ddl_info.dat";
348 infileddl.open(ddlinfofileName.Data(), ios::in); // ascii file
349 if(!infileddl) AliError("Could not read the ddl info file");
356 for(Int_t jddl = 0; jddl < 6; jddl++)
358 if (infileddl.eof()) break;
359 infileddl >> ddlno >> modulePerDDL;
360 moduleddl[jddl] = modulePerDDL;
362 if (modulePerDDL == 0) continue;
363 for (Int_t im = 0; im < modulePerDDL; im++)
371 // Set the minimum noise cut per module before clustering
373 fRecoParam = AliPMDReconstructor::GetRecoParam();
375 if(fRecoParam == 0x0)
377 AliFatal("No Reco Param found for PMD!!!");
382 Int_t bufsize = 16000;
383 TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
385 TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
387 const Int_t kRow = 48;
388 const Int_t kCol = 96;
394 AliPMDRawStream pmdinput(rawReader);
396 while ((indexDDL = pmdinput.DdlData(&pmdddlcont)) >=0)
398 iSMN = moduleddl[indexDDL];
401 precpvADC = new int **[iSMN];
402 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
403 for (Int_t i=0; i<iSMN;i++)
405 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
407 for (Int_t i = 0; i < iSMN; i++)
409 for (Int_t j = 0; j < kRow; j++)
411 for (Int_t k = 0; k < kCol; k++)
413 precpvADC[i][j][k] = 0;
420 Int_t ientries = pmdddlcont.GetEntries();
421 for (Int_t ient = 0; ient < ientries; ient++)
423 pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
425 Int_t det = pmdddl->GetDetector();
426 Int_t smn = pmdddl->GetSMN();
427 //Int_t mcm = pmdddl->GetMCM();
428 //Int_t chno = pmdddl->GetChannel();
429 Int_t row = pmdddl->GetRow();
430 Int_t col = pmdddl->GetColumn();
431 Int_t sig = pmdddl->GetSignal();
434 if(det < 0 || det > 1)
436 AliError(Form("*CPV/PRE NUMBER WRONG %d *",det));
439 if(smn < 0 || smn > 23)
441 AliError(Form("*MODULE NUMBER WRONG %d *",smn));
444 if(row < 0 || row > 47 || col < 0 || col > 95)
446 AliError(Form("*Row %d and Column NUMBER %d NOT Valid *",
453 // Pedestal Subtraction
454 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,row,col);
455 Int_t pedrms1 = (Int_t) pedmeanrms%100;
456 Float_t pedrms = (Float_t)pedrms1/10.;
457 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
459 //printf("%f %f\n",pedmean, pedrms);
461 // Float_t sig1 = (Float_t) sig;
462 Float_t sig1 = (Float_t) sig - (pedmean + 3.0*pedrms);
465 // Hot cell - set the cell adc = 0
466 Float_t hotflag = fCalibHot->GetHotChannel(det,smn,row,col);
467 if (hotflag == 1.) sig1 = 0;
470 Float_t gain = fCalibGain->GetGainFact(det,smn,row,col);
471 //printf("sig = %d gain = %f\n",sig,gain);
472 sig = (Int_t) (sig1*gain);
477 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
487 else if (smn >= 18 && smn < 24)
491 else if (indexDDL >= 1 && indexDDL < 4)
494 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
496 indexsmn = smn - indexDDL * 6;
498 else if (indexDDL == 4)
501 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
507 else if (smn >= 18 && smn < 24)
512 else if (indexDDL == 5)
515 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
517 if (smn >= 6 && smn < 18)
523 precpvADC[indexsmn][row][col] = sig;
531 for (indexsmn = 0; indexsmn < iSMN; indexsmn++)
535 for (Int_t irow = 0; irow < kRow; irow++)
537 for (Int_t icol = 0; icol < kCol; icol++)
539 fCellTrack[irow][icol] = -1;
540 fCellPid[irow][icol] = -1;
542 fCellADC[irow][icol] =
543 (Double_t) precpvADC[indexsmn][irow][icol];
544 totAdcMod += precpvADC[indexsmn][irow][icol];
559 else if (indexsmn >= 6 && indexsmn < 12)
560 ismn = indexsmn + 12;
564 else if (indexDDL >= 1 && indexDDL < 4)
566 ismn = indexsmn + indexDDL * 6;
569 else if (indexDDL == 4)
575 else if (indexsmn >= 6 && indexsmn < 12)
577 ismn = indexsmn + 12;
581 else if (indexDDL == 5)
587 if (totAdcMod <= 0) continue;
589 Int_t imod = idet*24 + ismn;
591 // Int_t cluspar = fRecoParam->GetPbPbParam()->GetClusteringParam();
592 Int_t cluspar = fRecoParam->GetPPParam()->GetClusteringParam();
593 // Int_t cluspar = fRecoParam->GetCosmicParam()->GetClusteringParam();
594 pmdclust->SetClusteringParam(cluspar);
595 Float_t encut = fNoiseCut->GetNoiseCut(imod);
597 pmdclust->SetEdepCut(encut);
598 pmdclust->DoClust(idet,ismn,fCellTrack,fCellPid,fCellADC,pmdcont);
600 Int_t nentries1 = pmdcont->GetEntries();
602 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
604 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
606 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
607 idet = pmdcl->GetDetector();
608 ismn = pmdcl->GetSMN();
609 clusdata[0] = pmdcl->GetClusX();
610 clusdata[1] = pmdcl->GetClusY();
611 clusdata[2] = pmdcl->GetClusADC();
612 clusdata[3] = pmdcl->GetClusCells();
613 clusdata[4] = pmdcl->GetClusSigmaX();
614 clusdata[5] = pmdcl->GetClusSigmaY();
616 AddRecPoint(idet,ismn,clusdata);
618 Int_t ncell = (Int_t) clusdata[3];
619 if (ncell > 19) ncell = 19;
620 for(Int_t ihit = 0; ihit < ncell; ihit++)
622 Int_t celldataX = pmdcl->GetClusCellX(ihit);
623 Int_t celldataY = pmdcl->GetClusCellY(ihit);
624 Int_t celldataTr = pmdcl->GetClusCellTrack(ihit);
625 Int_t celldataPid = pmdcl->GetClusCellPid(ihit);
626 Float_t celldataAdc = pmdcl->GetClusCellAdc(ihit);
627 AddRecHit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
641 for (Int_t i=0; i<iSMN; i++)
643 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
645 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
653 // delete the pointers
657 // ------------------------------------------------------------------------- //
658 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
660 // Add Reconstructed points
662 TClonesArray &lrecpoints = *fRecpoints;
663 AliPMDrecpoint1 *newrecpoint;
664 newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
665 new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
668 // ------------------------------------------------------------------------- //
669 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY,
670 Int_t celldataTr, Int_t celldataPid,
673 // Add associated cell hits to the Reconstructed points
675 TClonesArray &lrechits = *fRechits;
676 AliPMDrechit *newrechit;
677 newrechit = new AliPMDrechit(celldataX, celldataY, celldataTr, celldataPid, celldataAdc);
678 new(lrechits[fNhit++]) AliPMDrechit(newrechit);
681 // ------------------------------------------------------------------------- //
682 void AliPMDClusterFinder::ResetCellADC()
684 // Reset the individual cell ADC value to zero
686 for(Int_t irow = 0; irow < fgkRow; irow++)
688 for(Int_t icol = 0; icol < fgkCol; icol++)
690 fCellTrack[irow][icol] = -1;
691 fCellPid[irow][icol] = -1;
692 fCellADC[irow][icol] = 0.;
696 // ------------------------------------------------------------------------- //
697 void AliPMDClusterFinder::ResetRecpoint()
699 // Clear the list of reconstructed points
701 if (fRecpoints) fRecpoints->Clear();
703 // ------------------------------------------------------------------------- //
704 void AliPMDClusterFinder::ResetRechit()
706 // Clear the list of reconstructed points
708 if (fRechits) fRechits->Clear();
710 // ------------------------------------------------------------------------- //
711 void AliPMDClusterFinder::Load()
713 // Load all the *.root files
715 fPMDLoader->LoadDigits("READ");
716 fPMDLoader->LoadRecPoints("recreate");
718 // ------------------------------------------------------------------------- //
719 void AliPMDClusterFinder::LoadClusters()
721 // Load all the *.root files
723 fPMDLoader->LoadRecPoints("recreate");
725 // ------------------------------------------------------------------------- //
726 void AliPMDClusterFinder::UnLoad()
728 // Unload all the *.root files
730 fPMDLoader->UnloadDigits();
731 fPMDLoader->UnloadRecPoints();
733 // ------------------------------------------------------------------------- //
734 void AliPMDClusterFinder::UnLoadClusters()
736 // Unload all the *.root files
738 fPMDLoader->UnloadRecPoints();
740 // ------------------------------------------------------------------------- //
741 AliPMDCalibData* AliPMDClusterFinder::GetCalibGain() const
743 // The run number will be centralized in AliCDBManager,
744 // you don't need to set it here!
746 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
748 if(!entry) AliFatal("Calibration object retrieval failed! ");
750 AliPMDCalibData *calibdata=0;
751 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
753 if (!calibdata) AliFatal("No calibration data from calibration database !");
757 // ------------------------------------------------------------------------- //
758 AliPMDPedestal* AliPMDClusterFinder::GetCalibPed() const
760 // The run number will be centralized in AliCDBManager,
761 // you don't need to set it here!
762 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
764 if(!entry) AliFatal("Pedestal object retrieval failed!");
766 AliPMDPedestal *pedestal = 0;
767 if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
769 if (!pedestal) AliFatal("No pedestal data from pedestal database !");
773 //--------------------------------------------------------------------//
774 AliPMDHotData* AliPMDClusterFinder::GetCalibHot() const
776 // The run number will be centralized in AliCDBManager,
777 // you don't need to set it here!
778 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Hot");
780 if(!entry) AliFatal("HotData object retrieval failed!");
782 AliPMDHotData *hot = 0;
783 if (entry) hot = (AliPMDHotData*) entry->GetObject();
785 if (!hot) AliFatal("No hot data from database !");
789 //--------------------------------------------------------------------//
790 AliPMDNoiseCut* AliPMDClusterFinder::GetNoiseCut() 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/NoiseCut");
796 if(!entry) AliFatal("Noisecut object retrieval failed!");
798 AliPMDNoiseCut *ncut = 0;
799 if (entry) ncut = (AliPMDNoiseCut*) entry->GetObject();
801 if (!ncut) AliFatal("No noise cut data from database !");