memory leak fixed
[u/mrichter/AliRoot.git] / T0 / AliT0Reconstructor.cxx
CommitLineData
dc7ca31d 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
15
16/* $Id$ */
72e48d95 17/*********************************************************************
18 * T0 reconstruction and filling ESD
19 * - reconstruct mean time (interation time)
20 * - vertex position
21 * - multiplicity
22 ********************************************************************/
dc7ca31d 23
af885e0f 24#include <AliESDEvent.h>
dc7ca31d 25#include "AliLog.h"
dc7ca31d 26#include "AliT0RecPoint.h"
27#include "AliRawReader.h"
28#include "AliT0RawReader.h"
dc7ca31d 29#include "AliT0digit.h"
30#include "AliT0Reconstructor.h"
31#include "AliT0Parameters.h"
c41ceaac 32#include "AliT0Calibrator.h"
58bd3a16 33#include "AliESDfriend.h"
b0e13b29 34#include "AliESDTZERO.h"
73df58ab 35#include "AliESDTZEROfriend.h"
8f8d0732 36#include "AliLog.h"
85f61e3b 37#include "AliCDBEntry.h"
38#include "AliCDBManager.h"
39#include "AliCTPTimeParams.h"
40#include "AliLHCClockPhase.h"
669dc07f 41#include "AliT0CalibSeasonTimeShift.h"
4cbe597e 42#include "AliESDRun.h"
dc7ca31d 43
44#include <TArrayI.h>
45#include <TGraph.h>
aad72f45 46#include <TMath.h>
b09247a2 47#include <Riostream.h>
dc7ca31d 48
49ClassImp(AliT0Reconstructor)
50
c41ceaac 51 AliT0Reconstructor:: AliT0Reconstructor(): AliReconstructor(),
f16935f7 52 fdZonA(0),
53 fdZonC(0),
54 fZposition(0),
55 fParam(NULL),
2e6a5ee0 56 fAmpLEDrec(),
c883fdf2 57 fQTC(0),
58 fAmpLED(0),
58bd3a16 59 fCalib(),
60 fLatencyHPTDC(9000),
61 fLatencyL1(0),
62 fLatencyL1A(0),
73df58ab 63 fLatencyL1C(0),
85f61e3b 64 fGRPdelays(0),
669dc07f 65 fTimeMeanShift(0x0),
66 fTimeSigmaShift(0x0),
b0e13b29 67 fESDTZEROfriend(NULL),
68 fESDTZERO(NULL)
58bd3a16 69
e0bba6cc 70{
612737bb 71 for (Int_t i=0; i<24; i++) fTime0vertex[i] =0;
4cbe597e 72
72e48d95 73 //constructor
85f61e3b 74 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
75 if (!entry) AliFatal("CTP timing parameters are not found in OCDB !");
76 AliCTPTimeParams *ctpParams = (AliCTPTimeParams*)entry->GetObject();
77 Float_t l1Delay = (Float_t)ctpParams->GetDelayL1L0()*25.0;
78
79 AliCDBEntry *entry1 = AliCDBManager::Instance()->Get("GRP/CTP/TimeAlign");
80 if (!entry1) AliFatal("CTP time-alignment is not found in OCDB !");
81 AliCTPTimeParams *ctpTimeAlign = (AliCTPTimeParams*)entry1->GetObject();
82 l1Delay += ((Float_t)ctpTimeAlign->GetDelayL1L0()*25.0);
83
84 AliCDBEntry *entry4 = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
85 if (!entry4) AliFatal("LHC clock-phase shift is not found in OCDB !");
86 AliLHCClockPhase *phase = (AliLHCClockPhase*)entry4->GetObject();
87
85f61e3b 88 fGRPdelays = l1Delay - phase->GetMeanPhase();
72e48d95 89
669dc07f 90 AliCDBEntry *entry5 = AliCDBManager::Instance()->Get("T0/Calib/TimeAdjust");
91 if (entry5) {
92 AliT0CalibSeasonTimeShift *timeshift = (AliT0CalibSeasonTimeShift*)entry5->GetObject();
93 fTimeMeanShift = timeshift->GetT0Means();
94 fTimeSigmaShift = timeshift->GetT0Sigmas();
4cbe597e 95 }
669dc07f 96 else
97 AliWarning("Time Adjust is not found in OCDB !");
612737bb 98
74adb36a 99 fParam = AliT0Parameters::Instance();
100 fParam->Init();
c883fdf2 101
74adb36a 102 for (Int_t i=0; i<24; i++){
2e6a5ee0 103 TGraph* gr = fParam ->GetAmpLEDRec(i);
29ed1d0e 104 if (gr) fAmpLEDrec.AddAtAndExpand(gr,i) ;
c883fdf2 105 TGraph* gr1 = fParam ->GetAmpLED(i);
106 if (gr1) fAmpLED.AddAtAndExpand(gr1,i) ;
107 TGraph* gr2 = fParam ->GetQTC(i);
539b9cb9 108 if (gr2) fQTC.AddAtAndExpand(gr2,i) ;
612737bb 109 fTime0vertex[i] = fParam->GetCFD(i);
b0e13b29 110 AliDebug(2,Form("OCDB mean CFD time %i %f \n",i, fTime0vertex[i]));
612737bb 111 }
58bd3a16 112 fLatencyL1 = fParam->GetLatencyL1();
612737bb 113 fLatencyL1A = fParam->GetLatencyL1A();
58bd3a16 114 fLatencyL1C = fParam->GetLatencyL1C();
115 fLatencyHPTDC = fParam->GetLatencyHPTDC();
9d026202 116 AliDebug(2,Form(" LatencyL1 %f latencyL1A %f latencyL1C %f latencyHPTDC %f \n",fLatencyL1, fLatencyL1A, fLatencyL1C, fLatencyHPTDC));
612737bb 117
118 for (Int_t i=0; i<24; i++) {
59fe6376 119 if( fTime0vertex[i] < 500 || fTime0vertex[i] > 50000) fTime0vertex[i] =( 1000.*fLatencyHPTDC - 1000.*fLatencyL1 + 1000.*fGRPdelays)/24.4;
612737bb 120
121 }
8f620945 122 //here real Z position
123 fdZonC = TMath::Abs(fParam->GetZPosition("T0/C/PMT1"));
124 fdZonA = TMath::Abs(fParam->GetZPosition("T0/A/PMT15"));
539b9cb9 125
12e9daf9 126 fCalib = new AliT0Calibrator();
73df58ab 127 fESDTZEROfriend = new AliESDTZEROfriend();
b0e13b29 128 fESDTZERO = new AliESDTZERO();
12e9daf9 129
b0e13b29 130
c41ceaac 131}
c41ceaac 132
133//_____________________________________________________________________________
dc7ca31d 134void AliT0Reconstructor::Reconstruct(TTree*digitsTree, TTree*clustersTree) const
135{
94c27e4f 136 // T0 digits reconstruction
38cbfa7c 137 Int_t refAmp = Int_t (GetRecoParam()->GetRefAmp());
776de217 138
c41ceaac 139 TArrayI * timeCFD = new TArrayI(24);
140 TArrayI * timeLED = new TArrayI(24);
141 TArrayI * chargeQT0 = new TArrayI(24);
142 TArrayI * chargeQT1 = new TArrayI(24);
74adb36a 143
d0bcd1fb 144
8955c6b4 145 Float_t channelWidth = fParam->GetChannelWidth() ;
b95e8d87 146 Float_t meanVertex = fParam->GetMeanVertex();
776de217 147 Float_t c = 0.0299792; // cm/ps
adf36b9d 148 Double32_t vertex = 9999999;
149 Double32_t timeDiff=999999, meanTime=999999, timeclock=999999;
776de217 150
94c27e4f 151
dc7ca31d 152 AliDebug(1,Form("Start DIGITS reconstruction "));
94c27e4f 153
d3e04608 154 Float_t lowAmpThreshold = GetRecoParam()->GetLow(200);
155 Float_t highAmpThreshold = GetRecoParam()->GetHigh(200);
156 Int_t badpmt = GetRecoParam()->GetRefPoint();
2e6a5ee0 157
d0bcd1fb 158 TBranch *brDigits=digitsTree->GetBranch("T0");
e0bba6cc 159 AliT0digit *fDigits = new AliT0digit() ;
dc7ca31d 160 if (brDigits) {
161 brDigits->SetAddress(&fDigits);
162 }else{
f16935f7 163 AliError(Form("EXEC Branch T0 digits not found"));
164 return;
dc7ca31d 165 }
e0bba6cc 166
c41ceaac 167 digitsTree->GetEvent(0);
168 digitsTree->GetEntry(0);
169 brDigits->GetEntry(0);
170 fDigits->GetTimeCFD(*timeCFD);
171 fDigits->GetTimeLED(*timeLED);
172 fDigits->GetQT0(*chargeQT0);
173 fDigits->GetQT1(*chargeQT1);
446d6ec4 174 Int_t onlineMean = fDigits->MeanTime();
c883fdf2 175
adf36b9d 176 Bool_t tr[5];
177 for (Int_t i=0; i<5; i++) tr[i]=false;
c41ceaac 178
adf36b9d 179 Double32_t besttimeA=999999;
180 Double32_t besttimeC=999999;
c41ceaac 181 Int_t pmtBestA=99999;
182 Int_t pmtBestC=99999;
dc7ca31d 183
94c27e4f 184 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
1b544d5a 185 clustersTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints);
94c27e4f 186
195e1353 187 Float_t time[24], adc[24], adcmip[24];
dc7ca31d 188 for (Int_t ipmt=0; ipmt<24; ipmt++) {
d3e04608 189 if(timeCFD->At(ipmt)>0 && ipmt != badpmt) {
d0bcd1fb 190 if(( chargeQT1->At(ipmt) - chargeQT0->At(ipmt))>0)
191 adc[ipmt] = chargeQT1->At(ipmt) - chargeQT0->At(ipmt);
a7027400 192 else
d0bcd1fb 193 adc[ipmt] = 0;
194
38cbfa7c 195 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD->At(ipmt)) ;
d0bcd1fb 196
197 Double_t sl = Double_t(timeLED->At(ipmt) - timeCFD->At(ipmt));
8f620945 198 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD->At(ipmt) ) ;
669dc07f 199 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 200 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
776de217 201
d0bcd1fb 202 Double_t ampMip =((TGraph*)fAmpLED.At(ipmt))->Eval(sl);
203 Double_t qtMip = ((TGraph*)fQTC.At(ipmt))->Eval(adc[ipmt]);
669dc07f 204 AliDebug(5,Form(" Amlitude in MIPS LED %f , QTC %f in channels %f\n ",ampMip,qtMip, adc[ipmt]));
d0bcd1fb 205
206 frecpoints->SetTime(ipmt, Float_t(time[ipmt]) );
612737bb 207 frecpoints->SetAmpLED(ipmt, Float_t( ampMip));
345f03db 208 frecpoints->SetAmp(ipmt, Float_t(qtMip));
195e1353 209 adcmip[ipmt]=qtMip;
d0bcd1fb 210
dc7ca31d 211 }
212 else {
213 time[ipmt] = 0;
214 adc[ipmt] = 0;
b0ab3f59 215 adcmip[ipmt] = 0;
216
dc7ca31d 217 }
218 }
94c27e4f 219
dc7ca31d 220 for (Int_t ipmt=0; ipmt<12; ipmt++){
195e1353 221 if(time[ipmt] > 1 && ipmt != badpmt && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold) {
c41ceaac 222 if(time[ipmt]<besttimeC){
223 besttimeC=time[ipmt]; //timeC
224 pmtBestC=ipmt;
dc7ca31d 225 }
226 }
227 }
228 for ( Int_t ipmt=12; ipmt<24; ipmt++){
195e1353 229 if(time[ipmt] > 1 && ipmt != badpmt && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold) {
c41ceaac 230 if(time[ipmt]<besttimeA) {
231 besttimeA=time[ipmt]; //timeA
232 pmtBestA=ipmt;}
dc7ca31d 233 }
234 }
adf36b9d 235 if(besttimeA < 999999) {
8f620945 236 frecpoints->SetTimeBestA(Int_t(besttimeA *channelWidth - fdZonA/c));
adf36b9d 237 tr[1]=true;
238 }
239 if( besttimeC < 999999 ) {
8f620945 240 frecpoints->SetTimeBestC(Int_t(besttimeC *channelWidth - fdZonA/c));
adf36b9d 241 tr[2]=true;
242 }
669dc07f 243 AliDebug(5,Form(" besttimeA %f ch, besttimeC %f ch",besttimeA, besttimeC));
adf36b9d 244 if(besttimeA <999999 && besttimeC < 999999 ){
9b83615d 245 // timeDiff = (besttimeC - besttimeA)*channelWidth;
246 timeDiff = (besttimeA - besttimeC)*channelWidth;
adf36b9d 247 meanTime = (besttimeA + besttimeC)/2;// * channelWidth);
8f620945 248 timeclock = meanTime *channelWidth -fdZonA/c ;
adf36b9d 249 vertex = meanVertex - c*(timeDiff)/2.;// + (fdZonA - fdZonC)/2;
250 tr[0]=true;
251 }
252 frecpoints->SetVertex(vertex);
253 frecpoints->SetMeanTime(meanTime);
254 frecpoints->SetT0clock(timeclock);
255 frecpoints->SetT0Trig(tr);
256
669dc07f 257 AliDebug(5,Form("T0 triggers %d %d %d %d %d",tr[0],tr[1],tr[2],tr[3],tr[4]));
adf36b9d 258
259 //online mean
260 frecpoints->SetOnlineMean(Int_t(onlineMean));
b5a9f753 261 AliDebug(10,Form(" timeDiff %f #channel, meanTime %f #channel, vertex %f cm online mean %i timeclock %f ps",timeDiff, meanTime,vertex, Int_t(onlineMean), timeclock));
adf36b9d 262
263
264
265
b95e8d87 266
dc7ca31d 267 clustersTree->Fill();
bd375212 268
269 delete timeCFD;
270 delete timeLED;
271 delete chargeQT0;
272 delete chargeQT1;
dc7ca31d 273}
274
275
c41ceaac 276//_______________________________________________________________________
277
278void AliT0Reconstructor::Reconstruct(AliRawReader* rawReader, TTree*recTree) const
279{
94c27e4f 280 // T0 raw ->
539b9cb9 281 //
b0e13b29 282
283 Float_t meanOrA = fTime0vertex[0] + 587;
284 Float_t meanOrC = fTime0vertex[0] + 678;
285 Float_t meanTVDC = fTime0vertex[0] + 2564;
286
d3e04608 287
612737bb 288 Int_t badpmt[24];
38cbfa7c 289 //Bad channel
58432641 290 for (Int_t i=0; i<24; i++) {
291 badpmt[i] = GetRecoParam() -> GetBadChannels(i);
292 }
d3e04608 293 Int_t low[500], high[500];
94b4f7e2 294 Float_t timefull=-9999;;
295 Float_t tvdc = -9999; Float_t ora = -9999; Float_t orc = -9999;
58bd3a16 296
e8ed1cd0 297 Int_t allData[110][5];
2e6a5ee0 298
e8ed1cd0 299 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
612737bb 300 Double32_t timeDiff, meanTime, timeclock;
301 timeDiff = meanTime = timeclock = 9999999;
adf36b9d 302 Float_t c = 29.9792458; // cm/ns
303 Double32_t vertex = 9999999;
776de217 304 Int_t onlineMean=0;
9480f05f 305 Float_t meanVertex = 0;
612737bb 306 for (Int_t i0=0; i0<24; i0++) {
307 low[i0] = Int_t(fTime0vertex[i0]) - 200;
308 high[i0] = Int_t(fTime0vertex[i0]) + 200;
309 }
310
669dc07f 311 for (Int_t i0=0; i0<110; i0++)
b0e13b29 312 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
612737bb 313
314 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
315 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
d3e04608 316
adf36b9d 317 Double32_t besttimeA=9999999;
318 Double32_t besttimeC=9999999;
bce12dc5 319 Int_t pmtBestA=99999;
320 Int_t pmtBestC=99999;
b0e13b29 321 Float_t channelWidth = fParam->GetChannelWidth() ;
322
bce12dc5 323 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
324
1b544d5a 325 recTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints);
2e6a5ee0 326
bce12dc5 327 AliDebug(10," before read data ");
328 AliT0RawReader myrawreader(rawReader);
776de217 329
330 UInt_t type =rawReader->GetType();
331
bce12dc5 332 if (!myrawreader.Next())
333 AliDebug(1,Form(" no raw data found!!"));
334 else
335 {
38cbfa7c 336 for (Int_t i=0; i<24; i++)
337 {
338 timeCFD[i]=0; timeLED[i]=0; chargeQT0[i]=0; chargeQT1[i]=0;
339 }
340 Int_t fBCID=Int_t (rawReader->GetBCID());
8f620945 341 Int_t trmbunch= myrawreader.GetTRMBunchID();
86fd0587 342 AliDebug(10,Form(" CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
343
38cbfa7c 344 if(type == 7 ) { //only physics
669dc07f 345 for (Int_t i=0; i<107; i++) {
bce12dc5 346 for (Int_t iHit=0; iHit<5; iHit++)
347 {
348 allData[i][iHit] = myrawreader.GetData(i,iHit);
349 }
8f620945 350 }
351 Int_t ref=0;
8f620945 352
85f61e3b 353 for (Int_t in=0; in<12; in++)
354 {
355 for (Int_t iHit=0; iHit<5; iHit++)
38cbfa7c 356 {
612737bb 357 if(allData[in+1][iHit] > low[in] &&
358 allData[in+1][iHit] < high[in])
38cbfa7c 359 {
85f61e3b 360 timeCFD[in] = allData[in+1][iHit] ;
361 break;
38cbfa7c 362 }
85f61e3b 363 }
364 for (Int_t iHit=0; iHit<5; iHit++)
365 {
612737bb 366 if(allData[in+1+56][iHit] > low[in] &&
367 allData[in+1+56][iHit] < high[in])
38cbfa7c 368 {
85f61e3b 369 timeCFD[in+12] = allData[in+56+1][iHit] ;
370 break;
38cbfa7c 371 }
38cbfa7c 372 }
612737bb 373 timeLED[in+12] = allData[in+68+1][0] ;
374 timeLED[in] = allData[in+12+1][0] ;
669dc07f 375 AliDebug(5, Form(" readed i %i cfdC %i cfdA %i ledC %i ledA%i ",
85f61e3b 376 in, timeCFD[in],timeCFD[in+12],timeLED[in],
612737bb 377 timeLED[in+12]));
38cbfa7c 378
8f620945 379 }
380
8f620945 381
85f61e3b 382 for (Int_t in=0; in<12; in++)
383 {
384 chargeQT0[in]=allData[2*in+25][0];
385 chargeQT1[in]=allData[2*in+26][0];
612737bb 386 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 387 in, chargeQT0[in],chargeQT1[in]));
388 }
389 for (Int_t in=12; in<24; in++)
390 {
391 chargeQT0[in]=allData[2*in+57][0];
392 chargeQT1[in]=allData[2*in+58][0];
612737bb 393 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 394 in, chargeQT0[in],chargeQT1[in]));
85f61e3b 395 }
396
8f226345 397 onlineMean = allData[49][0];
398
195e1353 399 Double32_t time[24], adc[24], adcmip[24], noncalibtime[24];
8f620945 400 for (Int_t ipmt=0; ipmt<24; ipmt++) {
36bfca7d 401 if(timeCFD[ipmt] > 0 /* && badpmt[ipmt]==0*/ ){
bce12dc5 402 //for simulated data
d0bcd1fb 403 //for physics data
73df58ab 404 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>0) {
541b42c4 405 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 406 }
d0bcd1fb 407 else
408 adc[ipmt] = 0;
612737bb 409 // time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
bce12dc5 410
612737bb 411 time[ipmt] = fCalib-> WalkCorrection(Int_t (fTime0vertex[ipmt]), ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
d0bcd1fb 412 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 413 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 414 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 415 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
416 Double_t ampMip =( (TGraph*)fAmpLED.At(ipmt))->Eval(sl);
417 Double_t qtMip = ((TGraph*)fQTC.At(ipmt))->Eval(adc[ipmt]);
345f03db 418 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %f\n ",ampMip,qtMip, adc[ipmt]));
73df58ab 419 //bad peak removing
73df58ab 420 frecpoints->SetTime(ipmt, Float_t(time[ipmt]) );
421 // frecpoints->SetTime(ipmt,Double32_t(timeCFD[ipmt]));
195e1353 422 frecpoints->SetAmp(ipmt, Double32_t( qtMip));
423 adcmip[ipmt]=qtMip;
345f03db 424 frecpoints->SetAmpLED(ipmt, Double32_t(ampMip));
73df58ab 425 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
2e6a5ee0 426 }
427 else {
428 time[ipmt] = 0;
429 adc[ipmt] = 0;
b0ab3f59 430 adcmip[ipmt] = 0;
73df58ab 431 noncalibtime[ipmt] = 0;
2e6a5ee0 432 }
433 }
73df58ab 434 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
2e6a5ee0 435 for (Int_t ipmt=0; ipmt<12; ipmt++){
36bfca7d 436 if(time[ipmt] !=0 /*&& badpmt[ipmt]==0 */&& adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
73df58ab 437 {
b0e13b29 438 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC)) {
58432641 439 if(time[ipmt]<besttimeC){
612737bb 440 besttimeC=time[ipmt]; //timeC
669dc07f 441 pmtBestC=ipmt;
73df58ab 442 }
2e6a5ee0 443 }
2e6a5ee0 444 }
73df58ab 445 for ( Int_t ipmt=12; ipmt<24; ipmt++)
446 {
36bfca7d 447 if(time[ipmt] != 0 /* && badpmt[ipmt]==0*/ && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
73df58ab 448 {
58432641 449 if(time[ipmt]<besttimeA) {
450 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeA)) {
73df58ab 451 besttimeA=time[ipmt]; //timeA
452 pmtBestA=ipmt;
453 }
454 }
2e6a5ee0 455 }
adf36b9d 456 if(besttimeA < 999999)
58432641 457 frecpoints->SetTimeBestA((besttimeA * channelWidth)- 1000.*fLatencyHPTDC + 1000.*fLatencyL1A - 1000.*fGRPdelays - fTimeMeanShift[1] );
458 // frecpoints->SetTimeBestA((besttimeA * channelWidth- fTimeMeanShift[1]));
669dc07f 459
adf36b9d 460 if( besttimeC < 999999 )
58432641 461 frecpoints->SetTimeBestC((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1C - 1000.*fGRPdelays - fTimeMeanShift[2]);
462 // frecpoints->SetTimeBestC((besttimeC * channelWidth - fTimeMeanShift[2]));
669dc07f 463 AliDebug(5,Form(" pmtA %i besttimeA %f shift A %f ps, pmtC %i besttimeC %f shiftC %f ps",
464 pmtBestA,besttimeA, fTimeMeanShift[1],
465 pmtBestC, besttimeC,fTimeMeanShift[2]));
adf36b9d 466 if(besttimeA <999999 && besttimeC < 999999 ){
612737bb 467 // timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth + fLatencyL1A - fLatencyL1C;
58432641 468 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - 1000.*fLatencyHPTDC + 1000.*fLatencyL1 - 1000.*fGRPdelays - fTimeMeanShift[0] ;
adf36b9d 469 meanTime = (besttimeA+besttimeC-2.*Float_t(ref))/2.;
612737bb 470 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth ;
58432641 471 // timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - fTimeMeanShift[0] ;
9480f05f 472 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
adf36b9d 473 }
776de217 474 } //if phys event
291f31a1 475 AliDebug(10,Form(" timeDiff %f #channel, meanTime %f #channel, TOFmean%f vertex %f cm meanVertex %f online mean %i \n",timeDiff, meanTime,timeclock, vertex,meanVertex, onlineMean));
776de217 476 frecpoints->SetT0clock(timeclock);
477 frecpoints->SetVertex(vertex);
adf36b9d 478 frecpoints->SetMeanTime(meanTime);
776de217 479 frecpoints->SetOnlineMean(Int_t(onlineMean));
adf36b9d 480 // Set triggers
481
482 Bool_t tr[5];
483 Int_t trchan[5]= {50,51,52,55,56};
484 for (Int_t i=0; i<5; i++) tr[i]=false;
485 for (Int_t itr=0; itr<5; itr++) {
b0e13b29 486 for (Int_t iHit=0; iHit<1; iHit++)
38cbfa7c 487 {
488 Int_t trr=trchan[itr];
b0e13b29 489 if( allData[trr][iHit] > 0) tr[itr]=true;
490
491 AliDebug(1,Form("Reconstruct ::: T0 triggers iHit %i tvdc %d orA %d orC %d centr %d semicentral %d",iHit, tr[0],tr[1],tr[2],tr[3],tr[4]));
492 }
adf36b9d 493 }
38cbfa7c 494 frecpoints->SetT0Trig(tr);
b0e13b29 495
496 for (Int_t iHit=0; iHit<5; iHit++)
497 {
b0e13b29 498 if(allData[50][iHit]>0)
499 tvdc = (Float_t(allData[50][iHit]) - meanTVDC) * channelWidth* 0.001;
500 if(allData[51][iHit]>0)
501 ora = (Float_t(allData[51][iHit]) - meanOrA) * channelWidth* 0.001;
502
503 if(allData[52][iHit]>0)
504 orc = (Float_t(allData[52][iHit]) - meanOrC) * channelWidth* 0.001;
505
506 frecpoints->SetOrC( iHit, orc);
507 frecpoints->SetOrA( iHit, ora);
508 frecpoints->SetTVDC( iHit, tvdc);
b0e13b29 509 for (Int_t i0=0; i0<12; i0++) {
510 timefull = -9999;
511 if(allData[i0+1][iHit]>1)
512 timefull = (Float_t(allData[i0+1][iHit])-fTime0vertex[i0])* channelWidth* 0.001;
513 frecpoints->SetTimeFull(i0, iHit,timefull) ;
94b4f7e2 514 // printf("i0 %d iHit %d data %d fTime0vertex %f timefull %f \n",i0, iHit, allData[i0+1][iHit], fTime0vertex[i0], timefull);
b0e13b29 515
516 }
517
518 for (Int_t i0=12; i0<24; i0++) {
519 timefull = -9999;
94b4f7e2 520 if(allData[i0+45][iHit]>1) {
b0e13b29 521 timefull = (Float_t(allData[i0+45][iHit])-fTime0vertex[i0])* channelWidth* 0.001;
94b4f7e2 522 }
523 // printf("i0 %d iHit %d data %d fTime0vertex %f timefull %f \n",i0, iHit, allData[i0+45][iHit], fTime0vertex[i0], timefull);
b0e13b29 524 frecpoints->SetTimeFull(i0, iHit, timefull) ;
525 }
526 }
527
528
669dc07f 529 //Set MPD
530 if(allData[53][0]>0 && allData[54][0])
531 frecpoints->SetMultA(allData[53][0]-allData[54][0]);
b0e13b29 532 if(allData[105][0]>0 && allData[106][0])
533 frecpoints->SetMultC(allData[105][0]-allData[106][0]);
534
535
669dc07f 536 } // if (else )raw data
58bd3a16 537 recTree->Fill();
538 if(frecpoints) delete frecpoints;
539}
adf36b9d 540
541
542 //____________________________________________________________
543
544 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
545 {
dc7ca31d 546
547 /***************************************************
548 Resonstruct digits to vertex position
549 ****************************************************/
550
dc7ca31d 551 AliDebug(1,Form("Start FillESD T0"));
b0ab3f59 552 if(!pESD) {
553 AliError("No ESD Event");
554 return;
555 }
4cbe597e 556 pESD ->SetT0spread(fTimeSigmaShift);
557
36cde487 558
58bd3a16 559 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 560 Float_t c = 0.0299792458; // cm/ps
adf36b9d 561 Float_t currentVertex=0, shift=0;
291f31a1 562 Int_t ncont=-1;
adf36b9d 563 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
564 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
565 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
566 if (!vertex) vertex = pESD->GetVertex();
567
568 if (vertex) {
569 AliDebug(2, Form("Got %s (%s) from ESD: %f",
570 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
571 currentVertex = vertex->GetZ();
572
573 ncont = vertex->GetNContributors();
291f31a1 574 if(ncont>0 ) {
adf36b9d 575 shift = currentVertex/c;
adf36b9d 576 }
577 }
d76c31f4 578 TTree *treeR = clustersTree;
dc7ca31d 579
73df58ab 580 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
581 if (!frecpoints) {
dc7ca31d 582 AliError("Reconstruct Fill ESD >> no recpoints found");
583 return;
584 }
585
586 AliDebug(1,Form("Start FillESD T0"));
587 TBranch *brRec = treeR->GetBranch("T0");
588 if (brRec) {
589 brRec->SetAddress(&frecpoints);
590 }else{
f16935f7 591 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 592 return;
593 }
73df58ab 594
595 brRec->GetEntry(0);
596 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
597 Double32_t* tcorr;
598 for(Int_t i=0; i<24; i++)
599 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
600
669dc07f 601
73df58ab 602 Double32_t timeClock[3];
603 Double32_t zPosition = frecpoints -> GetVertex();
604 Double32_t timeStart = frecpoints -> GetMeanTime();
605 timeClock[0] = frecpoints -> GetT0clock() ;
606 timeClock[1] = frecpoints -> GetBestTimeA() + shift;
607 timeClock[2] = frecpoints -> GetBestTimeC() - shift;
612737bb 608
73df58ab 609 for ( Int_t i=0; i<24; i++) {
610 time[i] = frecpoints -> GetTime(i); // ps to ns
612737bb 611 // if ( time[i] >1) {
612 if ( time[i] != 0) {
345f03db 613 ampQTC[i] = frecpoints -> GetAmp(i);
614 amp[i] = frecpoints -> AmpLED(i);
612737bb 615 AliDebug(1,Form("T0: %i time %f ampQTC %f ampLED %f \n", i, time[i], ampQTC[i], amp[i]));
85f61e3b 616 }
73df58ab 617 }
618 Int_t trig= frecpoints ->GetT0Trig();
b0e13b29 619 frecpoints->PrintTriggerSignals( trig);
620 printf(" FillESD trigger %i \n",trig);
621 fESDTZERO->SetT0Trig(trig);
622 //pESD->SetT0Trig(trig);
623 // pESD->SetT0zVertex(zPosition); //vertex Z position
624 fESDTZERO->SetT0zVertex(zPosition); //vertex Z position
669dc07f 625
626 Double32_t multA=frecpoints ->GetMultA();
627 Double32_t multC=frecpoints ->GetMultC();
b0e13b29 628 // pESD->SetT0(multC); // multiplicity Cside
629 // pESD->SetT0clock(multA); // multiplicity Aside
630 fESDTZERO->SetMultA(multA); // for backward compatubility
631 fESDTZERO->SetMultC(multC); // for backward compatubility
632
633
634 for (Int_t iHit =0; iHit<5; iHit++ ) {
635 AliDebug(1,Form("FillESD ::: iHit %i tvdc %f orA %f orC %f\n", iHit,
636 frecpoints->GetTVDC(iHit),
637 frecpoints->GetOrA(iHit),
638 frecpoints->GetOrC(iHit) ));
639 fESDTZERO->SetTVDC(iHit,frecpoints->GetTVDC(iHit));
640 fESDTZERO->SetOrA(iHit,frecpoints->GetOrA(iHit));
641 fESDTZERO->SetOrC(iHit,frecpoints->GetOrC(iHit));
642
643 for (Int_t i0=0; i0<24; i0++) {
94b4f7e2 644 // if(frecpoints->GetTimeFull(i0,iHit)>0){
645 // printf("FillESD ::: iHit %i cfd %i time %f \n", iHit, i0, frecpoints->GetTimeFull(i0,iHit));
b0e13b29 646 fESDTZERO->SetTimeFull(i0, iHit,frecpoints->GetTimeFull(i0,iHit));
94b4f7e2 647 // }
b0e13b29 648
649 }
650 }
73df58ab 651 for(Int_t i=0; i<3; i++)
b0e13b29 652 fESDTZERO->SetT0TOF(i,timeClock[i]); // interaction time (ns)
653 fESDTZERO->SetT0time(time); // best TOF on each PMT
654 fESDTZERO->SetT0amplitude(ampQTC); // number of particles(MIPs) on each PMT
73df58ab 655
612737bb 656 AliDebug(1,Form("T0: SPDshift %f Vertex %f (T0A+T0C)/2 %f #channels T0signal %f ns OrA %f ns OrC %f T0trig %i\n",shift, zPosition, timeStart, timeClock[0], timeClock[1], timeClock[2], trig));
b0e13b29 657
73df58ab 658
b0e13b29 659
73df58ab 660 if (pESD) {
adf36b9d 661
73df58ab 662 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
663 if (fr) {
36cde487 664 AliDebug(1, Form("Writing TZERO friend data to ESD tree"));
73df58ab 665
85f61e3b 666 // if (ncont>2) {
73df58ab 667 tcorr = fESDTZEROfriend->GetT0timeCorr();
668 for ( Int_t i=0; i<24; i++) {
85f61e3b 669 if(i<12 && time[i]>1) timecorr[i] = tcorr[i] - shift/channelWidth;
670 if(i>11 && time[i]>1) timecorr[i] = tcorr[i] + shift/channelWidth;
612737bb 671 if(time[i]>1) AliDebug(10,Form("T0 friend : %i time %f ampQTC %f ampLED %f \n", i, timecorr[i], ampQTC[i], amp[i]));
58bd3a16 672 }
73df58ab 673 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
674 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
675 fESDTZEROfriend->SetT0ampQTC(ampQTC);
85f61e3b 676 fr->SetTZEROfriend(fESDTZEROfriend);
677 // }//
58bd3a16 678 }
b0e13b29 679
680 pESD->SetTZEROData(fESDTZERO);
58bd3a16 681 }
73df58ab 682
683
dc7ca31d 684} // vertex in 3 sigma
685
686
687
688
689
690