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 **************************************************************************/
18 Revision 1.26 2007/04/30 19:02:24 arcelli
19 hopefully the last refinements for correct type conversion in calibration
21 Revision 1.25 2007/04/30 15:22:17 arcelli
22 Change TOF digit Time, Tot etc to int type
24 Revision 1.24 2007/04/27 11:19:31 arcelli
25 updates for the new decoder
27 Revision 1.23 2007/04/23 16:51:39 decaro
28 Digits-to-raw_data conversion: correction for a more real description (A.De Caro, R.Preghenella)
30 Revision 1.22 2007/04/19 17:26:32 arcelli
31 Fix a bug (add some debug printout
33 Revision 1.21 2007/04/18 17:28:12 arcelli
34 Set the ToT bin width to the one actually used...
36 Revision 1.20 2007/03/09 09:57:23 arcelli
37 Remove a forgotten include of Riostrem
39 Revision 1.19 2007/03/08 15:41:20 arcelli
40 set uncorrected times when filling RecPoints
42 Revision 1.18 2007/03/06 16:31:20 arcelli
43 Add Uncorrected TOF Time signal
45 Revision 1.17 2007/02/28 18:09:11 arcelli
46 Add protection against failed retrieval of the CDB cal object, now Reconstruction exits with AliFatal
48 Revision 1.16 2007/02/20 15:57:00 decaro
49 Raw data update: to read the TOF raw data defined in UNPACKED mode
52 Revision 0.03 2005/07/28 A. De Caro
53 Implement public method
54 Raw2Digits(Int_t, AliRawReader *)
55 to convert digits from raw data in MC digits
58 Revision 0.02 2005/07/27 A. De Caro
59 Implement public method
60 Digits2RecPoint(Int_t)
61 to convert digits in clusters
63 Revision 0.02 2005/07/26 A. De Caro
64 Implement private methods
65 InsertCluster(AliTOFcluster *)
66 FindClusterIndex(Double_t)
67 originally implemented in AliTOFtracker
68 by S. Arcelli and C. Zampolli
70 Revision 0.01 2005/07/25 A. De Caro
71 Implement public methods
72 Digits2RecPoint(AliRawReader *, TTree *)
73 Digits2RecPoint(Int_t, AliRawReader *)
74 to convert raw data in clusters
77 ////////////////////////////////////////////////////////////////
79 // Class for TOF cluster finder //
81 // Starting from Raw Data, create rec points, //
82 // fill TreeR for TOF, //
83 // write TOF.RecPoints.root file //
85 ////////////////////////////////////////////////////////////////
88 #include "TClonesArray.h"
90 #include "TStopwatch.h"
94 #include "AliLoader.h"
96 #include "AliRawReader.h"
97 #include "AliRunLoader.h"
99 #include "AliTOFCal.h"
100 #include "AliTOFcalib.h"
101 #include "AliTOFChannel.h"
102 #include "AliTOFClusterFinder.h"
103 #include "AliTOFcluster.h"
104 #include "AliTOFdigit.h"
105 #include "AliTOFGeometry.h"
106 #include "AliTOFGeometryV5.h"
107 #include "AliTOFrawData.h"
108 #include "AliTOFRawStream.h"
109 #include "Riostream.h"
111 //extern TFile *gFile;
113 ClassImp(AliTOFClusterFinder)
115 AliTOFClusterFinder::AliTOFClusterFinder(AliTOFcalib *calib):
120 fTOFGeometry(new AliTOFGeometryV5()),
121 fDigits(new TClonesArray("AliTOFdigit", 4000)),
122 fRecPoints(new TClonesArray("AliTOFcluster", 4000)),
123 fNumberOfTofClusters(0),
132 AliInfo("V5 TOF Geometry is taken as the default");
135 //______________________________________________________________________________
137 AliTOFClusterFinder::AliTOFClusterFinder(AliRunLoader* runLoader, AliTOFcalib *calib):
138 fRunLoader(runLoader),
139 fTOFLoader(runLoader->GetLoader("TOFLoader")),
142 fTOFGeometry(new AliTOFGeometryV5()),
143 fDigits(new TClonesArray("AliTOFdigit", 4000)),
144 fRecPoints(new TClonesArray("AliTOFcluster", 4000)),
145 fNumberOfTofClusters(0),
156 //------------------------------------------------------------------------
157 AliTOFClusterFinder::AliTOFClusterFinder(const AliTOFClusterFinder &source)
163 fTOFGeometry(new AliTOFGeometryV5()),
164 fDigits(new TClonesArray("AliTOFdigit", 4000)),
165 fRecPoints(new TClonesArray("AliTOFcluster", 4000)),
166 fNumberOfTofClusters(0),
172 this->fDigits=source.fDigits;
173 this->fRecPoints=source.fRecPoints;
174 this->fTOFGeometry=source.fTOFGeometry;
175 this->fDecoderVersion=source.fDecoderVersion;
176 this->fTOFcalib=source.fTOFcalib;
180 //------------------------------------------------------------------------
181 AliTOFClusterFinder& AliTOFClusterFinder::operator=(const AliTOFClusterFinder &source)
184 this->fDigits=source.fDigits;
185 this->fRecPoints=source.fRecPoints;
186 this->fTOFGeometry=source.fTOFGeometry;
187 this->fVerbose=source.fVerbose;
188 this->fDecoderVersion=source.fDecoderVersion;
189 this->fTOFcalib=source.fTOFcalib;
193 //______________________________________________________________________________
195 AliTOFClusterFinder::~AliTOFClusterFinder()
210 fRecPoints->Delete();
218 //______________________________________________________________________________
220 void AliTOFClusterFinder::Digits2RecPoints(Int_t iEvent)
223 // Converts digits to recpoints for TOF
226 TStopwatch stopwatch;
229 fRunLoader->GetEvent(iEvent);
231 fTreeD = fTOFLoader->TreeD();
234 AliFatal("AliTOFClusterFinder: Can not get TreeD");
237 TBranch *branch = fTreeD->GetBranch("TOF");
239 AliError("can't get the branch with the TOF digits !");
243 TClonesArray *digits = new TClonesArray("AliTOFdigit",10000);
244 branch->SetAddress(&digits);
248 fTreeR = fTOFLoader->TreeR();
251 fTOFLoader->MakeTree("R");
252 fTreeR = fTOFLoader->TreeR();
255 Int_t bufsize = 32000;
256 fTreeR->Branch("TOF", &fRecPoints, bufsize);
259 Int_t nDigits = digits->GetEntriesFast();
260 AliDebug(2,Form("Number of TOF digits: %d",nDigits));
263 Int_t dig[5]; //cluster detector indeces
264 Float_t g[3]; //cluster cartesian coord
265 Double_t h[3]; // the cluster spatial cyl. coordinates
266 Int_t parTOF[5]; //The TOF signal parameters
267 Bool_t status=kTRUE; // assume all sim channels ok in the beginning...
268 for (ii=0; ii<nDigits; ii++) {
269 AliTOFdigit *d = (AliTOFdigit*)digits->UncheckedAt(ii);
270 dig[0]=d->GetSector();
271 dig[1]=d->GetPlate();
272 dig[2]=d->GetStrip();
276 // AliDebug(2,Form(" %2i %1i %2i %1i %2i ",dig[0],dig[1],dig[2],dig[3],dig[4]));
278 for (jj=0; jj<3; jj++) g[jj] = 0.;
279 fTOFGeometry->GetPos(dig,g);
281 h[0] = TMath::Sqrt(g[0]*g[0]+g[1]*g[1]);
282 h[1] = TMath::ATan2(g[1],g[0]);
285 parTOF[0] = d->GetTdc(); //the TDC signal
286 parTOF[1] = d->GetToT(); //the ToT signal
287 parTOF[2] = d->GetAdc(); // the adc charge
288 parTOF[3] = d->GetTdcND(); // non decalibrated sim time
289 parTOF[4] = d->GetTdc(); // raw time, == Tdc time for the moment
290 AliTOFcluster *tofCluster = new AliTOFcluster(h,dig,parTOF,status,d->GetTracks(),ii);
291 InsertCluster(tofCluster);
295 AliInfo(Form("Number of found clusters: %i for event: %i", fNumberOfTofClusters, iEvent));
303 fTOFLoader = fRunLoader->GetLoader("TOFLoader");
304 fTOFLoader->WriteRecPoints("OVERWRITE");
306 AliInfo(Form("Execution time to read TOF digits and to write TOF clusters : R:%.4fs C:%.4fs",
307 stopwatch.RealTime(),stopwatch.CpuTime()));
313 //______________________________________________________________________________
315 void AliTOFClusterFinder::Digits2RecPoints(TTree* digitsTree, TTree* clusterTree)
318 // Converts digits to recpoints for TOF
321 TStopwatch stopwatch;
324 /// fRunLoader->GetEvent(iEvent);
326 if (digitsTree == 0x0)
328 AliFatal("AliTOFClusterFinder: Can not get TreeD");
331 TBranch *branch = digitsTree->GetBranch("TOF");
333 AliError("can't get the branch with the TOF digits !");
337 TClonesArray *digits = new TClonesArray("AliTOFdigit",10000);
338 branch->SetAddress(&digits);
343 Int_t bufsize = 32000;
344 fTreeR->Branch("TOF", &fRecPoints, bufsize);
346 digitsTree->GetEvent(0);
347 Int_t nDigits = digits->GetEntriesFast();
348 AliDebug(2,Form("Number of TOF digits: %d",nDigits));
351 Int_t dig[5]; //cluster detector indeces
352 Float_t g[3]; //cluster cartesian coord
353 Double_t h[3]; // the cluster spatial cyl. coordinates
354 Int_t parTOF[5]; //The TOF signal parameters
355 Bool_t status=kTRUE; // assume all sim channels ok in the beginning...
356 for (ii=0; ii<nDigits; ii++) {
357 AliTOFdigit *d = (AliTOFdigit*)digits->UncheckedAt(ii);
358 dig[0]=d->GetSector();
359 dig[1]=d->GetPlate();
360 dig[2]=d->GetStrip();
364 // AliDebug(2,Form(" %2i %1i %2i %1i %2i ",dig[0],dig[1],dig[2],dig[3],dig[4]));
366 for (jj=0; jj<3; jj++) g[jj] = 0.;
367 fTOFGeometry->GetPos(dig,g);
369 h[0] = TMath::Sqrt(g[0]*g[0]+g[1]*g[1]);
370 h[1] = TMath::ATan2(g[1],g[0]);
373 parTOF[0] = d->GetTdc(); //the TDC signal
374 parTOF[1] = d->GetToT(); //the ToT signal
375 parTOF[2] = d->GetAdc(); // the adc charge
376 parTOF[3] = d->GetTdcND(); // non decalibrated sim time
377 parTOF[4] = d->GetTdc(); // raw time, == Tdc time for the moment
378 AliTOFcluster *tofCluster = new AliTOFcluster(h,dig,parTOF,status,d->GetTracks(),ii);
379 InsertCluster(tofCluster);
383 AliInfo(Form("Number of found clusters: %i", fNumberOfTofClusters));
391 AliInfo(Form("Execution time to read TOF digits and to write TOF clusters : R:%.4fs C:%.4fs",
392 stopwatch.RealTime(),stopwatch.CpuTime()));
397 //______________________________________________________________________________
399 void AliTOFClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
403 // Converts RAW data to recpoints for TOF
406 TStopwatch stopwatch;
409 //const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
410 const Int_t kDDL = AliDAQ::NumberOfDdls("TOF");
414 Int_t bufsize = 32000;
415 clustersTree->Branch("TOF", &fRecPoints, bufsize);
417 TClonesArray * clonesRawData;
422 Int_t detectorIndex[5];
424 Double_t cylindricalPosition[3];
428 if (fVerbose==2) ftxt.open("TOFdigitsRead.txt",ios::app);
431 for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
434 AliTOFRawStream tofInput(rawReader);
435 if (fDecoderVersion) {
436 AliInfo("Using New Decoder \n");
437 tofInput.LoadRawDataBuffers(indexDDL,fVerbose);
439 else tofInput.LoadRawData(indexDDL);
441 clonesRawData = (TClonesArray*)tofInput.GetRawData();
443 for (Int_t iRawData = 0; iRawData<clonesRawData->GetEntriesFast(); iRawData++) {
445 AliTOFrawData *tofRawDatum = (AliTOFrawData*)clonesRawData->UncheckedAt(iRawData);
447 if (tofRawDatum->GetTOT()==-1 || tofRawDatum->GetTOF()==-1) continue;
450 if (indexDDL<10) ftxt << " " << indexDDL;
451 else ftxt << " " << indexDDL;
452 if (tofRawDatum->GetTRM()<10) ftxt << " " << tofRawDatum->GetTRM();
453 else ftxt << " " << tofRawDatum->GetTRM();
454 ftxt << " " << tofRawDatum->GetTRMchain();
455 if (tofRawDatum->GetTDC()<10) ftxt << " " << tofRawDatum->GetTDC();
456 else ftxt << " " << tofRawDatum->GetTDC();
457 ftxt << " " << tofRawDatum->GetTDCchannel();
460 tofInput.EquipmentId2VolumeId(indexDDL, tofRawDatum->GetTRM(), tofRawDatum->GetTRMchain(),
461 tofRawDatum->GetTDC(), tofRawDatum->GetTDCchannel(), detectorIndex);
462 dummy = detectorIndex[3];
463 detectorIndex[3] = detectorIndex[4];
464 detectorIndex[4] = dummy;
467 if (detectorIndex[0]<10) ftxt << " -> " << detectorIndex[0];
468 else ftxt << " -> " << detectorIndex[0];
469 ftxt << " " << detectorIndex[1];
470 if (detectorIndex[2]<10) ftxt << " " << detectorIndex[2];
471 else ftxt << " " << detectorIndex[2];
472 ftxt << " " << detectorIndex[3];
473 if (detectorIndex[4]<10) ftxt << " " << detectorIndex[4];
474 else ftxt << " " << detectorIndex[4];
477 for (ii=0; ii<3; ii++) position[ii] = 0.;
478 fTOFGeometry->GetPos(detectorIndex, position);
480 cylindricalPosition[0] = TMath::Sqrt(position[0]*position[0] + position[1]*position[1]);
481 cylindricalPosition[1] = TMath::ATan2(position[1], position[0]);
482 cylindricalPosition[2] = position[2];
484 parTOF[0] = tofRawDatum->GetTOF(); //TDC
485 parTOF[1] = tofRawDatum->GetTOT(); // TOT
486 parTOF[2] = tofRawDatum->GetTOT(); //ADC==TOF
487 parTOF[3] = -1;//raw data: no track of undecalib sim time
488 parTOF[4] = tofRawDatum->GetTOF(); // RAW time
489 AliTOFcluster *tofCluster = new AliTOFcluster(cylindricalPosition, detectorIndex, parTOF);
490 InsertCluster(tofCluster);
493 if (parTOF[1]<10)ftxt << " " << parTOF[1];
494 else if (parTOF[1]>=10 && parTOF[1]<100) ftxt << " " << parTOF[1];
495 else ftxt << " " << parTOF[1];
496 if (parTOF[0]<10) ftxt << " " << parTOF[0] << endl;
497 else if (parTOF[0]>=10 && parTOF[0]<100) ftxt << " " << parTOF[0] << endl;
498 else if (parTOF[0]>=100 && parTOF[0]<1000) ftxt << " " << parTOF[0] << endl;
499 else ftxt << " " << parTOF[3] << endl;
502 } // closed loop on TOF raw data per current DDL file
504 clonesRawData->Clear();
506 } // closed loop on DDL index
508 if (fVerbose==2) ftxt.close();
510 AliInfo(Form("Number of found clusters: %i", fNumberOfTofClusters));
515 clustersTree->Fill();
519 AliDebug(1, Form("Execution time to read TOF raw data and to write TOF clusters : R:%.4fs C:%.4fs",
520 stopwatch.RealTime(),stopwatch.CpuTime()));
523 //______________________________________________________________________________
525 void AliTOFClusterFinder::Digits2RecPoints(Int_t iEvent, AliRawReader *rawReader)
528 // Converts RAW data to recpoints for TOF
531 TStopwatch stopwatch;
534 //const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
535 const Int_t kDDL = AliDAQ::NumberOfDdls("TOF");
537 fRunLoader->GetEvent(iEvent);
539 AliDebug(2,Form(" Event number %2i ", iEvent));
541 fTreeR = fTOFLoader->TreeR();
544 fTOFLoader->MakeTree("R");
545 fTreeR = fTOFLoader->TreeR();
548 Int_t bufsize = 32000;
549 fTreeR->Branch("TOF", &fRecPoints, bufsize);
551 TClonesArray * clonesRawData;
556 Int_t detectorIndex[5] = {-1, -1, -1, -1, -1};
558 Double_t cylindricalPosition[5];
561 if (fVerbose==2) ftxt.open("TOFdigitsRead.txt",ios::app);
564 for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
567 AliTOFRawStream tofInput(rawReader);
568 if (fDecoderVersion) {
569 AliInfo("Using New Decoder \n");
570 tofInput.LoadRawDataBuffers(indexDDL,fVerbose);
572 else tofInput.LoadRawData(indexDDL);
574 clonesRawData = (TClonesArray*)tofInput.GetRawData();
576 for (Int_t iRawData = 0; iRawData<clonesRawData->GetEntriesFast(); iRawData++) {
578 AliTOFrawData *tofRawDatum = (AliTOFrawData*)clonesRawData->UncheckedAt(iRawData);
580 if (tofRawDatum->GetTOT()==-1 || tofRawDatum->GetTOF()==-1) continue;
583 if (indexDDL<10) ftxt << " " << indexDDL;
584 else ftxt << " " << indexDDL;
585 if (tofRawDatum->GetTRM()<10) ftxt << " " << tofRawDatum->GetTRM();
586 else ftxt << " " << tofRawDatum->GetTRM();
587 ftxt << " " << tofRawDatum->GetTRMchain();
588 if (tofRawDatum->GetTDC()<10) ftxt << " " << tofRawDatum->GetTDC();
589 else ftxt << " " << tofRawDatum->GetTDC();
590 ftxt << " " << tofRawDatum->GetTDCchannel();
593 tofInput.EquipmentId2VolumeId(indexDDL, tofRawDatum->GetTRM(), tofRawDatum->GetTRMchain(),
594 tofRawDatum->GetTDC(), tofRawDatum->GetTDCchannel(), detectorIndex);
595 dummy = detectorIndex[3];
596 detectorIndex[3] = detectorIndex[4];
597 detectorIndex[4] = dummy;
600 if (detectorIndex[0]<10) ftxt << " -> " << detectorIndex[0];
601 else ftxt << " -> " << detectorIndex[0];
602 ftxt << " " << detectorIndex[1];
603 if (detectorIndex[2]<10) ftxt << " " << detectorIndex[2];
604 else ftxt << " " << detectorIndex[2];
605 ftxt << " " << detectorIndex[3];
606 if (detectorIndex[4]<10) ftxt << " " << detectorIndex[4];
607 else ftxt << " " << detectorIndex[4];
610 for (ii=0; ii<3; ii++) position[ii] = 0.;
611 fTOFGeometry->GetPos(detectorIndex, position);
613 cylindricalPosition[0] = TMath::Sqrt(position[0]*position[0] + position[1]*position[1]);
614 cylindricalPosition[1] = TMath::ATan2(position[1], position[0]);
615 cylindricalPosition[2] = position[2];
616 parTOF[0] = tofRawDatum->GetTOF(); // TDC
617 parTOF[1] = tofRawDatum->GetTOT(); // TOT
618 parTOF[2] = tofRawDatum->GetTOT(); // raw data have ADC=TOT
619 parTOF[3] = -1; //raw data: no track of the undecalib sim time
620 parTOF[4] = tofRawDatum->GetTOF(); // Raw time == TDC
621 AliTOFcluster *tofCluster = new AliTOFcluster(cylindricalPosition, detectorIndex, parTOF);
622 InsertCluster(tofCluster);
625 if (parTOF[1]<10)ftxt << " " << parTOF[1];
626 else if (parTOF[1]>=10 && parTOF[1]<100) ftxt << " " << parTOF[1];
627 else ftxt << " " << parTOF[1];
628 if (parTOF[0]<10) ftxt << " " << parTOF[0] << endl;
629 else if (parTOF[0]>=10 && parTOF[0]<100) ftxt << " " << parTOF[0] << endl;
630 else if (parTOF[0]>=100 && parTOF[0]<1000) ftxt << " " << parTOF[0] << endl;
631 else ftxt << " " << parTOF[3] << endl;
634 } // closed loop on TOF raw data per current DDL file
636 clonesRawData->Clear();
638 } // closed loop on DDL index
640 if (fVerbose==2) ftxt.close();
642 AliInfo(Form("Number of found clusters: %i for event: %i", fNumberOfTofClusters, iEvent));
650 fTOFLoader = fRunLoader->GetLoader("TOFLoader");
651 fTOFLoader->WriteRecPoints("OVERWRITE");
653 AliDebug(1, Form("Execution time to read TOF raw data and to write TOF clusters : R:%.4fs C:%.4fs",
654 stopwatch.RealTime(),stopwatch.CpuTime()));
657 //______________________________________________________________________________
659 void AliTOFClusterFinder::Raw2Digits(Int_t iEvent, AliRawReader *rawReader)
662 // Converts RAW data to MC digits for TOF
664 // (temporary solution)
667 TStopwatch stopwatch;
670 //const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
671 const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
673 fRunLoader->GetEvent(iEvent);
675 fTreeD = fTOFLoader->TreeD();
678 AliInfo("TreeD re-creation");
680 fTOFLoader->MakeTree("D");
681 fTreeD = fTOFLoader->TreeD();
684 TClonesArray *tofDigits = new TClonesArray("AliTOFdigit",10000);
685 Int_t bufsize = 32000;
686 fTreeD->Branch("TOF", &tofDigits, bufsize);
688 fRunLoader->GetEvent(iEvent);
690 AliDebug(2,Form(" Event number %2i ", iEvent));
692 TClonesArray * clonesRawData;
696 Int_t detectorIndex[5];
700 for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
703 AliTOFRawStream tofInput(rawReader);
704 if (fDecoderVersion) {
705 AliInfo("Using New Decoder \n");
706 tofInput.LoadRawDataBuffers(indexDDL,fVerbose);
708 else tofInput.LoadRawData(indexDDL);
710 clonesRawData = (TClonesArray*)tofInput.GetRawData();
712 for (Int_t iRawData = 0; iRawData<clonesRawData->GetEntriesFast(); iRawData++) {
714 AliTOFrawData *tofRawDatum = (AliTOFrawData*)clonesRawData->UncheckedAt(iRawData);
716 if (!tofRawDatum->GetTOT() || !tofRawDatum->GetTOF()) continue;
718 tofInput.EquipmentId2VolumeId(indexDDL, tofRawDatum->GetTRM(), tofRawDatum->GetTRMchain(),
719 tofRawDatum->GetTDC(), tofRawDatum->GetTDCchannel(), detectorIndex);
720 dummy = detectorIndex[3];
721 detectorIndex[3] = detectorIndex[4];
722 detectorIndex[4] = dummy;
724 digit[0] = tofInput.GetTofBin();
725 digit[1] = tofInput.GetToTbin();
726 digit[2] = tofInput.GetToTbin();
729 Int_t tracknum[3]={-1,-1,-1};
731 TClonesArray &aDigits = *tofDigits;
732 Int_t last=tofDigits->GetEntriesFast();
733 new (aDigits[last]) AliTOFdigit(tracknum, detectorIndex, digit);
737 clonesRawData->Clear();
743 fTOFLoader = fRunLoader->GetLoader("TOFLoader");
744 fTOFLoader->WriteDigits("OVERWRITE");
746 AliDebug(1, Form("Execution time to read TOF raw data and to write TOF clusters : R:%.2fs C:%.2fs",
747 stopwatch.RealTime(),stopwatch.CpuTime()));
751 //______________________________________________________________________________
753 void AliTOFClusterFinder::Raw2Digits(AliRawReader *rawReader, TTree* digitsTree)
756 // Converts RAW data to MC digits for TOF for the current event
760 TStopwatch stopwatch;
763 const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
767 AliError("No input digits Tree");
771 TClonesArray *tofDigits = new TClonesArray("AliTOFdigit",10000);
772 Int_t bufsize = 32000;
773 digitsTree->Branch("TOF", &tofDigits, bufsize);
775 /// fRunLoader->GetEvent(iEvent);
777 /// AliDebug(2,Form(" Event number %2i ", iEvent));
779 TClonesArray * clonesRawData;
783 Int_t detectorIndex[5];
787 for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
790 AliTOFRawStream tofInput(rawReader);
791 if (fDecoderVersion) {
792 AliInfo("Using New Decoder \n");
793 tofInput.LoadRawDataBuffers(indexDDL,fVerbose);
795 else tofInput.LoadRawData(indexDDL);
797 clonesRawData = (TClonesArray*)tofInput.GetRawData();
799 for (Int_t iRawData = 0; iRawData<clonesRawData->GetEntriesFast(); iRawData++) {
801 AliTOFrawData *tofRawDatum = (AliTOFrawData*)clonesRawData->UncheckedAt(iRawData);
803 if (!tofRawDatum->GetTOT() || !tofRawDatum->GetTOF()) continue;
805 tofInput.EquipmentId2VolumeId(indexDDL, tofRawDatum->GetTRM(), tofRawDatum->GetTRMchain(),
806 tofRawDatum->GetTDC(), tofRawDatum->GetTDCchannel(), detectorIndex);
807 dummy = detectorIndex[3];
808 detectorIndex[3] = detectorIndex[4];
809 detectorIndex[4] = dummy;
811 digit[0] = tofInput.GetTofBin();
812 digit[1] = tofInput.GetToTbin();
813 digit[2] = tofInput.GetToTbin();
816 Int_t tracknum[3]={-1,-1,-1};
818 TClonesArray &aDigits = *tofDigits;
819 Int_t last=tofDigits->GetEntriesFast();
820 new (aDigits[last]) AliTOFdigit(tracknum, detectorIndex, digit);
824 clonesRawData->Clear();
830 AliDebug(1, Form("Got %d digits: ", tofDigits->GetEntries()));
831 AliDebug(1, Form("Execution time to read TOF raw data and fill TOF digit tree : R:%.2fs C:%.2fs",
832 stopwatch.RealTime(),stopwatch.CpuTime()));
835 //______________________________________________________________________________
837 Int_t AliTOFClusterFinder::InsertCluster(AliTOFcluster *tofCluster) {
838 //---------------------------------------------------------------------------//
839 // This function adds a TOF cluster to the array of TOF clusters sorted in Z //
840 //---------------------------------------------------------------------------//
841 if (fNumberOfTofClusters==kTofMaxCluster) {
842 AliError("Too many clusters !");
846 if (fNumberOfTofClusters==0) {
847 fTofClusters[fNumberOfTofClusters++] = tofCluster;
851 Int_t ii = FindClusterIndex(tofCluster->GetZ());
852 memmove(fTofClusters+ii+1 ,fTofClusters+ii,(fNumberOfTofClusters-ii)*sizeof(AliTOFcluster*));
853 fTofClusters[ii] = tofCluster;
854 fNumberOfTofClusters++;
859 //_________________________________________________________________________
861 Int_t AliTOFClusterFinder::FindClusterIndex(Double_t z) const {
862 //--------------------------------------------------------------------
863 // This function returns the index of the nearest cluster
864 //--------------------------------------------------------------------
865 if (fNumberOfTofClusters==0) return 0;
866 if (z <= fTofClusters[0]->GetZ()) return 0;
867 if (z > fTofClusters[fNumberOfTofClusters-1]->GetZ()) return fNumberOfTofClusters;
868 Int_t b = 0, e = fNumberOfTofClusters-1, m = (b+e)/2;
869 for (; b<e; m=(b+e)/2) {
870 if (z > fTofClusters[m]->GetZ()) b=m+1;
877 //_________________________________________________________________________
879 void AliTOFClusterFinder::FillRecPoint()
882 // Copy the global array of AliTOFcluster, i.e. fTofClusters (sorted
883 // in Z) in the global TClonesArray of AliTOFcluster,
889 Int_t detectorIndex[5];
890 Double_t cylindricalPosition[3];
892 Int_t trackLabels[3];
893 Int_t digitIndex = -1;
896 TClonesArray &lRecPoints = *fRecPoints;
898 for (ii=0; ii<fNumberOfTofClusters; ii++) {
900 digitIndex = fTofClusters[ii]->GetIndex();
901 for(jj=0; jj<5; jj++) detectorIndex[jj] = fTofClusters[ii]->GetDetInd(jj);
902 for(jj=0; jj<3; jj++) trackLabels[jj] = fTofClusters[ii]->GetLabel(jj);
903 cylindricalPosition[0] = fTofClusters[ii]->GetR();
904 cylindricalPosition[1] = fTofClusters[ii]->GetPhi();
905 cylindricalPosition[2] = fTofClusters[ii]->GetZ();
906 parTOF[0] = fTofClusters[ii]->GetTDC(); // TDC
907 parTOF[1] = fTofClusters[ii]->GetToT(); // TOT
908 parTOF[2] = fTofClusters[ii]->GetADC(); // ADC=TOT
909 parTOF[3] = fTofClusters[ii]->GetTDCND(); // TDCND
910 parTOF[4] = fTofClusters[ii]->GetTDCRAW();//RAW
911 status=fTofClusters[ii]->GetStatus();
912 new(lRecPoints[ii]) AliTOFcluster(cylindricalPosition, detectorIndex, parTOF,status,trackLabels,digitIndex);
914 } // loop on clusters
918 //_________________________________________________________________________
919 void AliTOFClusterFinder::CalibrateRecPoint()
922 // Copy the global array of AliTOFcluster, i.e. fTofClusters (sorted
923 // in Z) in the global TClonesArray of AliTOFcluster,
929 Int_t detectorIndex[5];
930 Int_t digitIndex = -1;
934 AliInfo(" Calibrating TOF Clusters: ")
936 AliTOFCal *calTOFArray = fTOFcalib->GetTOFCalArray();
938 for (ii=0; ii<fNumberOfTofClusters; ii++) {
939 digitIndex = fTofClusters[ii]->GetIndex();
940 for(jj=0; jj<5; jj++) detectorIndex[jj] = fTofClusters[ii]->GetDetInd(jj);
942 Int_t index = fTOFcalib->GetIndex(detectorIndex);
944 AliTOFChannel * calChannel = calTOFArray->GetChannel(index);
946 // Get channel status
947 Bool_t status=calChannel->GetStatus();
948 if(status)fTofClusters[ii]->SetStatus(!status); //odd convention, to avoid conflict with calibration objects currently in the db (temporary solution).
949 // Get Rough channel online equalization
950 Float_t roughDelay=calChannel->GetDelay();
951 AliDebug(2,Form(" channel delay = %f", roughDelay));
952 // Get Refined channel offline calibration parameters
954 for (Int_t j = 0; j<6; j++){
955 par[j]=calChannel->GetSlewPar(j);
956 AliDebug(2,Form(" Calib Pars = %f, %f, %f, %f, %f, %f ",par[0],par[1],par[2],par[3],par[4],par[5]));
958 AliDebug(2,Form(" The ToT and Time, uncorr (counts) = %i , %i", fTofClusters[ii]->GetToT(),fTofClusters[ii]->GetTDC()));
959 tToT = (Double_t)(fTofClusters[ii]->GetToT()*AliTOFGeometry::ToTBinWidth()); tToT*=1.E-3; //ToT in ns
960 AliDebug(2,Form(" The ToT and Time, uncorr (ns)= %e, %e",fTofClusters[ii]->GetTDC()*AliTOFGeometry::TdcBinWidth()*1.E-3,tToT));
961 timeCorr=par[0]+par[1]*tToT+par[2]*tToT*tToT+par[3]*tToT*tToT*tToT+par[4]*tToT*tToT*tToT*tToT+par[5]*tToT*tToT*tToT*tToT*tToT+roughDelay; // the time correction
962 AliDebug(2,Form(" The time correction (ns) = %f", timeCorr));
963 timeCorr=(Double_t)(fTofClusters[ii]->GetTDC())*AliTOFGeometry::TdcBinWidth()*1.E-3-timeCorr;//redefine the time
965 AliDebug(2,Form(" The channel time, corr (ps)= %e",timeCorr ));
966 tdcCorr=(Int_t)(timeCorr/AliTOFGeometry::TdcBinWidth()); //the corrected time (tdc counts)
967 fTofClusters[ii]->SetTDC(tdcCorr);
968 AliDebug(2,Form(" The channel time, corr (counts) = %i",fTofClusters[ii]->GetTDC()));
970 } // loop on clusters
973 //______________________________________________________________________________
975 void AliTOFClusterFinder::ResetRecpoint()
978 // Clear the list of reconstructed points
981 fNumberOfTofClusters = 0;
982 if (fRecPoints) fRecPoints->Clear();
985 //______________________________________________________________________________
987 void AliTOFClusterFinder::Load()
990 // Load TOF.Digits.root and TOF.RecPoints.root files
993 fTOFLoader->LoadDigits("READ");
994 fTOFLoader->LoadRecPoints("recreate");
997 //______________________________________________________________________________
999 void AliTOFClusterFinder::LoadClusters()
1002 // Load TOF.RecPoints.root file
1005 fTOFLoader->LoadRecPoints("recreate");
1008 //______________________________________________________________________________
1010 void AliTOFClusterFinder::UnLoad()
1013 // Unload TOF.Digits.root and TOF.RecPoints.root files
1016 fTOFLoader->UnloadDigits();
1017 fTOFLoader->UnloadRecPoints();
1020 //______________________________________________________________________________
1022 void AliTOFClusterFinder::UnLoadClusters()
1025 // Unload TOF.RecPoints.root file
1028 fTOFLoader->UnloadRecPoints();
1031 //______________________________________________________________________________