correct warnings
[u/mrichter/AliRoot.git] / T0 / AliT0Reconstructor.cxx
CommitLineData
73df58ab 1
dc7ca31d 2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17/* $Id$ */
72e48d95 18/*********************************************************************
19 * T0 reconstruction and filling ESD
20 * - reconstruct mean time (interation time)
21 * - vertex position
22 * - multiplicity
23 ********************************************************************/
dc7ca31d 24
af885e0f 25#include <AliESDEvent.h>
dc7ca31d 26#include "AliLog.h"
dc7ca31d 27#include "AliT0RecPoint.h"
28#include "AliRawReader.h"
29#include "AliT0RawReader.h"
dc7ca31d 30#include "AliT0digit.h"
31#include "AliT0Reconstructor.h"
32#include "AliT0Parameters.h"
c41ceaac 33#include "AliT0Calibrator.h"
58bd3a16 34#include "AliESDfriend.h"
73df58ab 35#include "AliESDTZEROfriend.h"
8f8d0732 36#include "AliLog.h"
dc7ca31d 37
38#include <TArrayI.h>
39#include <TGraph.h>
aad72f45 40#include <TMath.h>
b09247a2 41#include <Riostream.h>
dc7ca31d 42
43ClassImp(AliT0Reconstructor)
44
c41ceaac 45 AliT0Reconstructor:: AliT0Reconstructor(): AliReconstructor(),
f16935f7 46 fdZonA(0),
47 fdZonC(0),
48 fZposition(0),
49 fParam(NULL),
2e6a5ee0 50 fAmpLEDrec(),
c883fdf2 51 fQTC(0),
52 fAmpLED(0),
58bd3a16 53 fCalib(),
54 fLatencyHPTDC(9000),
55 fLatencyL1(0),
56 fLatencyL1A(0),
73df58ab 57 fLatencyL1C(0),
58 fESDTZEROfriend(NULL)
58bd3a16 59
e0bba6cc 60{
72e48d95 61 //constructor
62
74adb36a 63 fParam = AliT0Parameters::Instance();
64 fParam->Init();
c883fdf2 65
74adb36a 66 for (Int_t i=0; i<24; i++){
2e6a5ee0 67 TGraph* gr = fParam ->GetAmpLEDRec(i);
29ed1d0e 68 if (gr) fAmpLEDrec.AddAtAndExpand(gr,i) ;
c883fdf2 69 TGraph* gr1 = fParam ->GetAmpLED(i);
70 if (gr1) fAmpLED.AddAtAndExpand(gr1,i) ;
71 TGraph* gr2 = fParam ->GetQTC(i);
539b9cb9 72 if (gr2) fQTC.AddAtAndExpand(gr2,i) ;
c883fdf2 73 }
539b9cb9 74
58bd3a16 75 fLatencyL1 = fParam->GetLatencyL1();
76 fLatencyL1A = fParam->GetLatencyL1A();
77 fLatencyL1C = fParam->GetLatencyL1C();
78 fLatencyHPTDC = fParam->GetLatencyHPTDC();
79 AliDebug(10,Form(" LatencyL1 %f latencyL1A %f latencyL1C %f latencyHPTDC %f \n",fLatencyL1, fLatencyL1A, fLatencyL1C, fLatencyHPTDC));
29ed1d0e 80
adf36b9d 81 // fdZonC = TMath::Abs(fParam->GetZPositionShift("T0/C/PMT1"));
82 //fdZonA = TMath::Abs(fParam->GetZPositionShift("T0/A/PMT15"));
8f620945 83 //here real Z position
84 fdZonC = TMath::Abs(fParam->GetZPosition("T0/C/PMT1"));
85 fdZonA = TMath::Abs(fParam->GetZPosition("T0/A/PMT15"));
539b9cb9 86
12e9daf9 87 fCalib = new AliT0Calibrator();
73df58ab 88 fESDTZEROfriend = new AliESDTZEROfriend();
12e9daf9 89
dc7ca31d 90}
c41ceaac 91
92//_____________________________________________________________________________
dc7ca31d 93void AliT0Reconstructor::Reconstruct(TTree*digitsTree, TTree*clustersTree) const
94c27e4f 94
dc7ca31d 95{
94c27e4f 96 // T0 digits reconstruction
c555f418 97 Int_t refAmp = GetRecoParam()->GetRefAmp();
776de217 98
c41ceaac 99 TArrayI * timeCFD = new TArrayI(24);
100 TArrayI * timeLED = new TArrayI(24);
101 TArrayI * chargeQT0 = new TArrayI(24);
102 TArrayI * chargeQT1 = new TArrayI(24);
74adb36a 103
d0bcd1fb 104
8955c6b4 105 Float_t channelWidth = fParam->GetChannelWidth() ;
b95e8d87 106 Float_t meanVertex = fParam->GetMeanVertex();
776de217 107 Float_t c = 0.0299792; // cm/ps
adf36b9d 108 Double32_t vertex = 9999999;
109 Double32_t timeDiff=999999, meanTime=999999, timeclock=999999;
776de217 110
94c27e4f 111
dc7ca31d 112 AliDebug(1,Form("Start DIGITS reconstruction "));
94c27e4f 113
2e6a5ee0 114
d0bcd1fb 115 TBranch *brDigits=digitsTree->GetBranch("T0");
e0bba6cc 116 AliT0digit *fDigits = new AliT0digit() ;
dc7ca31d 117 if (brDigits) {
118 brDigits->SetAddress(&fDigits);
119 }else{
f16935f7 120 AliError(Form("EXEC Branch T0 digits not found"));
121 return;
dc7ca31d 122 }
e0bba6cc 123
c41ceaac 124 digitsTree->GetEvent(0);
125 digitsTree->GetEntry(0);
126 brDigits->GetEntry(0);
127 fDigits->GetTimeCFD(*timeCFD);
128 fDigits->GetTimeLED(*timeLED);
129 fDigits->GetQT0(*chargeQT0);
130 fDigits->GetQT1(*chargeQT1);
446d6ec4 131 Int_t onlineMean = fDigits->MeanTime();
c883fdf2 132
adf36b9d 133 Bool_t tr[5];
134 for (Int_t i=0; i<5; i++) tr[i]=false;
c41ceaac 135
adf36b9d 136 Double32_t besttimeA=999999;
137 Double32_t besttimeC=999999;
c41ceaac 138 Int_t pmtBestA=99999;
139 Int_t pmtBestC=99999;
dc7ca31d 140
94c27e4f 141 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
142 clustersTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints, 405,1);
143
b95e8d87 144 Float_t time[24], adc[24];
dc7ca31d 145 for (Int_t ipmt=0; ipmt<24; ipmt++) {
c41ceaac 146 if(timeCFD->At(ipmt)>0 ){
d0bcd1fb 147 if(( chargeQT1->At(ipmt) - chargeQT0->At(ipmt))>0)
148 adc[ipmt] = chargeQT1->At(ipmt) - chargeQT0->At(ipmt);
a7027400 149 else
d0bcd1fb 150 adc[ipmt] = 0;
151
8f620945 152 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, adc[ipmt], timeCFD->At(ipmt)) ;
d0bcd1fb 153
154 Double_t sl = Double_t(timeLED->At(ipmt) - timeCFD->At(ipmt));
8f620945 155 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD->At(ipmt) ) ;
d0bcd1fb 156 AliDebug(10,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
157 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
776de217 158
d0bcd1fb 159 Double_t ampMip =((TGraph*)fAmpLED.At(ipmt))->Eval(sl);
160 Double_t qtMip = ((TGraph*)fQTC.At(ipmt))->Eval(adc[ipmt]);
b5a9f753 161 AliDebug(10,Form(" Amlitude in MIPS LED %f , QTC %f in channels %f\n ",ampMip,qtMip, adc[ipmt]));
d0bcd1fb 162
163 frecpoints->SetTime(ipmt, Float_t(time[ipmt]) );
164 frecpoints->SetAmp(ipmt, Float_t( ampMip)); //for cosmic &pp beam
165 frecpoints->SetAmpLED(ipmt, Float_t(qtMip));
166
dc7ca31d 167 }
168 else {
169 time[ipmt] = 0;
170 adc[ipmt] = 0;
171 }
172 }
94c27e4f 173
dc7ca31d 174 for (Int_t ipmt=0; ipmt<12; ipmt++){
175 if(time[ipmt] > 1 ) {
c41ceaac 176 if(time[ipmt]<besttimeC){
177 besttimeC=time[ipmt]; //timeC
178 pmtBestC=ipmt;
dc7ca31d 179 }
180 }
181 }
182 for ( Int_t ipmt=12; ipmt<24; ipmt++){
183 if(time[ipmt] > 1) {
c41ceaac 184 if(time[ipmt]<besttimeA) {
185 besttimeA=time[ipmt]; //timeA
186 pmtBestA=ipmt;}
dc7ca31d 187 }
188 }
adf36b9d 189 if(besttimeA < 999999) {
8f620945 190 frecpoints->SetTimeBestA(Int_t(besttimeA *channelWidth - fdZonA/c));
adf36b9d 191 tr[1]=true;
192 }
193 if( besttimeC < 999999 ) {
8f620945 194 frecpoints->SetTimeBestC(Int_t(besttimeC *channelWidth - fdZonA/c));
adf36b9d 195 tr[2]=true;
196 }
d0bcd1fb 197 AliDebug(10,Form(" besttimeA %f ch, besttimeC %f ch",besttimeA, besttimeC));
adf36b9d 198 if(besttimeA <999999 && besttimeC < 999999 ){
9b83615d 199 // timeDiff = (besttimeC - besttimeA)*channelWidth;
200 timeDiff = (besttimeA - besttimeC)*channelWidth;
adf36b9d 201 meanTime = (besttimeA + besttimeC)/2;// * channelWidth);
8f620945 202 timeclock = meanTime *channelWidth -fdZonA/c ;
adf36b9d 203 vertex = meanVertex - c*(timeDiff)/2.;// + (fdZonA - fdZonC)/2;
204 tr[0]=true;
205 }
206 frecpoints->SetVertex(vertex);
207 frecpoints->SetMeanTime(meanTime);
208 frecpoints->SetT0clock(timeclock);
209 frecpoints->SetT0Trig(tr);
210
8f620945 211 AliDebug(10,Form("T0 triggers %d %d %d %d %d",tr[0],tr[1],tr[2],tr[3],tr[4]));
adf36b9d 212
213 //online mean
214 frecpoints->SetOnlineMean(Int_t(onlineMean));
b5a9f753 215 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 216
217
218
219
b95e8d87 220
dc7ca31d 221 clustersTree->Fill();
bd375212 222
223 delete timeCFD;
224 delete timeLED;
225 delete chargeQT0;
226 delete chargeQT1;
dc7ca31d 227}
228
229
c41ceaac 230//_______________________________________________________________________
231
232void AliT0Reconstructor::Reconstruct(AliRawReader* rawReader, TTree*recTree) const
233{
94c27e4f 234 // T0 raw ->
539b9cb9 235 //
236 // reference amplitude and time ref. point from reco param
237
c555f418 238 Int_t refAmp = GetRecoParam()->GetRefAmp();
239 Int_t refPoint = GetRecoParam()->GetRefPoint();
58bd3a16 240
e8ed1cd0 241 Int_t allData[110][5];
2e6a5ee0 242
e8ed1cd0 243 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
adf36b9d 244 Double32_t timeDiff=999999, meanTime=999999, timeclock=999999;
245 Float_t c = 29.9792458; // cm/ns
246 Double32_t vertex = 9999999;
776de217 247 Int_t onlineMean=0;
9480f05f 248 // Float_t meanVertex = fParam->GetMeanVertex();
249 Float_t meanVertex = 0;
bce12dc5 250 for (Int_t i0=0; i0<105; i0++)
251 {
252 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
253 }
2e6a5ee0 254
adf36b9d 255 Double32_t besttimeA=9999999;
256 Double32_t besttimeC=9999999;
bce12dc5 257 Int_t pmtBestA=99999;
258 Int_t pmtBestC=99999;
29a60970 259
bce12dc5 260 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
261
262 recTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints, 405,1);
2e6a5ee0 263
bce12dc5 264 AliDebug(10," before read data ");
265 AliT0RawReader myrawreader(rawReader);
776de217 266
267 UInt_t type =rawReader->GetType();
268
bce12dc5 269 if (!myrawreader.Next())
270 AliDebug(1,Form(" no raw data found!!"));
271 else
272 {
8f620945 273 Int_t fBCID=Int_t (rawReader->GetBCID());
274 Int_t trmbunch= myrawreader.GetTRMBunchID();
275 if(type == 7 &&( trmbunch -fBCID )==37 ) { //only physics
276 for (Int_t i=0; i<105; i++) {
bce12dc5 277 for (Int_t iHit=0; iHit<5; iHit++)
278 {
279 allData[i][iHit] = myrawreader.GetData(i,iHit);
280 }
8f620945 281 }
282 Int_t ref=0;
283 if (refPoint>0)
284 ref = allData[refPoint][0]-5000;
285
286 Float_t channelWidth = fParam->GetChannelWidth() ;
287
288 // Int_t meanT0 = fParam->GetMeanT0();
289
290
291 for (Int_t in=0; in<12; in++)
292 {
293 timeCFD[in] = allData[in+1][0] ;
294 timeCFD[in+12] = allData[in+56+1][0] ;
295 timeLED[in] = allData[in+12+1][0] ;
296 timeLED[in+12] = allData[in+68+1][0] ;
297 AliDebug(10, Form(" readed i %i cfdC %i cfdA %i ledC %i ledA%i ",
298 in, timeCFD[in],timeCFD[in+12],timeLED[in],
299 timeLED[in+12]));
300 }
301
302 for (Int_t in=0; in<12; in++)
303 {
304 chargeQT0[in]=allData[2*in+25][0];
305 chargeQT1[in]=allData[2*in+26][0];
306 }
307
308 for (Int_t in=12; in<24; in++)
309 {
310 chargeQT0[in]=allData[2*in+57][0];
311 chargeQT1[in]=allData[2*in+58][0];
312 }
313
541b42c4 314
8f620945 315 for (Int_t in=0; in<24; in++)
316 AliDebug(10, Form(" readed Raw %i %i %i %i %i",
317 in, timeLED[in],timeCFD[in],chargeQT0[in],chargeQT1[in]));
776de217 318 onlineMean = allData[49][0];
8f620945 319
73df58ab 320 Double32_t time[24], adc[24], noncalibtime[24];
8f620945 321 for (Int_t ipmt=0; ipmt<24; ipmt++) {
2e6a5ee0 322 if(timeCFD[ipmt]>0 && timeLED[ipmt]>0){
bce12dc5 323 //for simulated data
d0bcd1fb 324 //for physics data
73df58ab 325 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>0) {
541b42c4 326 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 327 }
d0bcd1fb 328 else
329 adc[ipmt] = 0;
bce12dc5 330
539b9cb9 331
8f620945 332 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, adc[ipmt], timeCFD[ipmt] ) ;
bce12dc5 333
d0bcd1fb 334 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 335 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 336 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 337 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
338 Double_t ampMip =( (TGraph*)fAmpLED.At(ipmt))->Eval(sl);
339 Double_t qtMip = ((TGraph*)fQTC.At(ipmt))->Eval(adc[ipmt]);
f7c2c2fc 340 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %i\n ",ampMip,qtMip, adc[ipmt]));
73df58ab 341 //bad peak removing
8f620945 342 if(sl<550) {
73df58ab 343 frecpoints->SetTime(ipmt, Float_t(time[ipmt]) );
344 // frecpoints->SetTime(ipmt,Double32_t(timeCFD[ipmt]));
345 frecpoints->SetAmpLED(ipmt, Double32_t( qtMip)); //for cosmic &pp beam
346 frecpoints->SetAmp(ipmt, Double32_t(ampMip));
347 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
348 }
2e6a5ee0 349 }
350 else {
351 time[ipmt] = 0;
352 adc[ipmt] = 0;
73df58ab 353 noncalibtime[ipmt] = 0;
2e6a5ee0 354 }
355 }
73df58ab 356 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
2e6a5ee0 357 for (Int_t ipmt=0; ipmt<12; ipmt++){
73df58ab 358 if(time[ipmt] > 1
f7c2c2fc 359 && (timeLED[ipmt] - timeCFD[ipmt])<540 )
73df58ab 360 {
361 if(time[ipmt]<besttimeC){
362 besttimeC=time[ipmt]; //timeC
363 pmtBestC=ipmt;
364 }
2e6a5ee0 365 }
2e6a5ee0 366 }
73df58ab 367 for ( Int_t ipmt=12; ipmt<24; ipmt++)
368 {
369 if(time[ipmt] > 1
f7c2c2fc 370 && (timeLED[ipmt] - timeCFD[ipmt])<540 )
73df58ab 371 {
372 if(time[ipmt]<besttimeA) {
373 besttimeA=time[ipmt]; //timeA
374 pmtBestA=ipmt;
375 }
376 }
2e6a5ee0 377 }
adf36b9d 378 if(besttimeA < 999999)
73df58ab 379 frecpoints->SetTimeBestA(besttimeA * channelWidth - 1000.*fLatencyHPTDC + 1000.*fLatencyL1A);
adf36b9d 380 if( besttimeC < 999999 )
73df58ab 381 frecpoints->SetTimeBestC(besttimeC * channelWidth - 1000.*fLatencyHPTDC +1000.*fLatencyL1C);
adf36b9d 382 AliDebug(10,Form(" pmtA %i besttimeA %f ps, pmtC %i besttimeC %f ps",
d0bcd1fb 383 pmtBestA,besttimeA, pmtBestC, besttimeC));
adf36b9d 384 if(besttimeA <999999 && besttimeC < 999999 ){
73df58ab 385 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth + fLatencyL1A - fLatencyL1C;
386 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - 1000.*fLatencyHPTDC + 1000.*fLatencyL1;
adf36b9d 387 meanTime = (besttimeA+besttimeC-2.*Float_t(ref))/2.;
9480f05f 388 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
adf36b9d 389 }
776de217 390 } //if phys event
adf36b9d 391 AliDebug(5,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 392 frecpoints->SetT0clock(timeclock);
393 frecpoints->SetVertex(vertex);
adf36b9d 394 frecpoints->SetMeanTime(meanTime);
776de217 395 frecpoints->SetOnlineMean(Int_t(onlineMean));
adf36b9d 396 // Set triggers
397
398 Bool_t tr[5];
399 Int_t trchan[5]= {50,51,52,55,56};
400 for (Int_t i=0; i<5; i++) tr[i]=false;
401 for (Int_t itr=0; itr<5; itr++) {
402 if(allData[trchan[itr]][0]>0) tr[itr]=true;
403 frecpoints->SetT0Trig(tr);
404 }
58bd3a16 405 } // if (else )raw data
406 recTree->Fill();
407 if(frecpoints) delete frecpoints;
408}
adf36b9d 409
410
411 //____________________________________________________________
412
413 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
414 {
dc7ca31d 415
416 /***************************************************
417 Resonstruct digits to vertex position
418 ****************************************************/
419
dc7ca31d 420 AliDebug(1,Form("Start FillESD T0"));
58bd3a16 421 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 422 Float_t c = 0.0299792458; // cm/ps
adf36b9d 423 Float_t currentVertex=0, shift=0;
424 Int_t ncont=0;
425 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
426 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
427 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
428 if (!vertex) vertex = pESD->GetVertex();
429
430 if (vertex) {
431 AliDebug(2, Form("Got %s (%s) from ESD: %f",
432 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
433 currentVertex = vertex->GetZ();
434
435 ncont = vertex->GetNContributors();
436 // cout<<" spdver "<<spdver<<" ncont "<<ncont<<endl;
58bd3a16 437 if(ncont>2 ) {
adf36b9d 438 shift = currentVertex/c;
439 // cout<<" vertex shif "<<shift<<" vertex "<<spdver<<" IsFromVertexer3D "<<fverSPD->IsFromVertexer3D()<<endl;
440 }
441 }
d76c31f4 442 TTree *treeR = clustersTree;
dc7ca31d 443
73df58ab 444 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
445 if (!frecpoints) {
dc7ca31d 446 AliError("Reconstruct Fill ESD >> no recpoints found");
447 return;
448 }
449
450 AliDebug(1,Form("Start FillESD T0"));
451 TBranch *brRec = treeR->GetBranch("T0");
452 if (brRec) {
453 brRec->SetAddress(&frecpoints);
454 }else{
f16935f7 455 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 456 return;
457 }
73df58ab 458
459 brRec->GetEntry(0);
460 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
461 Double32_t* tcorr;
462 for(Int_t i=0; i<24; i++)
463 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
464
465 Double32_t timeClock[3];
466 Double32_t zPosition = frecpoints -> GetVertex();
467 Double32_t timeStart = frecpoints -> GetMeanTime();
468 timeClock[0] = frecpoints -> GetT0clock() ;
469 timeClock[1] = frecpoints -> GetBestTimeA() + shift;
470 timeClock[2] = frecpoints -> GetBestTimeC() - shift;
471 for ( Int_t i=0; i<24; i++) {
472 time[i] = frecpoints -> GetTime(i); // ps to ns
473 if ( time[i] >1) {
dc7ca31d 474 amp[i] = frecpoints -> GetAmp(i);
58bd3a16 475 ampQTC[i] = frecpoints -> AmpLED(i);
dc7ca31d 476 }
73df58ab 477 }
478 Int_t trig= frecpoints ->GetT0Trig();
479 pESD->SetT0Trig(trig);
480
481 pESD->SetT0zVertex(zPosition); //vertex Z position
482 pESD->SetT0(timeStart); // interaction time
483 for(Int_t i=0; i<3; i++)
484 pESD->SetT0TOF(i,timeClock[i]); // interaction time (ns)
485 pESD->SetT0time(time); // best TOF on each PMT
486 pESD->SetT0amplitude(amp); // number of particles(MIPs) on each PMT
487
488 AliDebug(1,Form("T0: Vertex %f (T0A+T0C)/2 %f #channels T0signal %f ns OrA %f ns OrC %f T0trig %i\n",zPosition, timeStart, timeClock[0], timeClock[1], timeClock[2], trig));
489
490 if (pESD) {
adf36b9d 491
73df58ab 492 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
493 if (fr) {
494 AliDebug(1, Form("Writing TZERO friend data to ESD tree"));
495
496 if (ncont>2) {
497 tcorr = fESDTZEROfriend->GetT0timeCorr();
498 for ( Int_t i=0; i<24; i++) {
499 timecorr[i]=tcorr[i];
f7c2c2fc 500 if(i<12 && time[i]>1) timecorr[i] -= shift/channelWidth;
501 if(i>11 && time[i]>1) timecorr[i] += shift/channelWidth;
58bd3a16 502 }
73df58ab 503 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
504 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
505 fESDTZEROfriend->SetT0ampQTC(ampQTC);
506 fr->SetTZEROfriend(fESDTZEROfriend);
507 }//
58bd3a16 508 }
509 }
73df58ab 510
511
512
dc7ca31d 513} // vertex in 3 sigma
514
515
516
517
518
519