Bug fix
[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];
58bd3a16 294
e8ed1cd0 295 Int_t allData[110][5];
2e6a5ee0 296
e8ed1cd0 297 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
612737bb 298 Double32_t timeDiff, meanTime, timeclock;
299 timeDiff = meanTime = timeclock = 9999999;
adf36b9d 300 Float_t c = 29.9792458; // cm/ns
301 Double32_t vertex = 9999999;
776de217 302 Int_t onlineMean=0;
9480f05f 303 Float_t meanVertex = 0;
612737bb 304 for (Int_t i0=0; i0<24; i0++) {
305 low[i0] = Int_t(fTime0vertex[i0]) - 200;
306 high[i0] = Int_t(fTime0vertex[i0]) + 200;
307 }
308
669dc07f 309 for (Int_t i0=0; i0<110; i0++)
b0e13b29 310 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
612737bb 311
312 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
313 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
d3e04608 314
adf36b9d 315 Double32_t besttimeA=9999999;
316 Double32_t besttimeC=9999999;
bce12dc5 317 Int_t pmtBestA=99999;
318 Int_t pmtBestC=99999;
b0e13b29 319 Float_t channelWidth = fParam->GetChannelWidth() ;
320
bce12dc5 321 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
322
1b544d5a 323 recTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints);
2e6a5ee0 324
bce12dc5 325 AliDebug(10," before read data ");
326 AliT0RawReader myrawreader(rawReader);
776de217 327
328 UInt_t type =rawReader->GetType();
329
bce12dc5 330 if (!myrawreader.Next())
331 AliDebug(1,Form(" no raw data found!!"));
332 else
333 {
38cbfa7c 334 for (Int_t i=0; i<24; i++)
335 {
336 timeCFD[i]=0; timeLED[i]=0; chargeQT0[i]=0; chargeQT1[i]=0;
337 }
338 Int_t fBCID=Int_t (rawReader->GetBCID());
8f620945 339 Int_t trmbunch= myrawreader.GetTRMBunchID();
86fd0587 340 AliDebug(10,Form(" CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
341
38cbfa7c 342 if(type == 7 ) { //only physics
669dc07f 343 for (Int_t i=0; i<107; i++) {
bce12dc5 344 for (Int_t iHit=0; iHit<5; iHit++)
345 {
346 allData[i][iHit] = myrawreader.GetData(i,iHit);
347 }
8f620945 348 }
349 Int_t ref=0;
8f620945 350
85f61e3b 351 for (Int_t in=0; in<12; in++)
352 {
353 for (Int_t iHit=0; iHit<5; iHit++)
38cbfa7c 354 {
612737bb 355 if(allData[in+1][iHit] > low[in] &&
356 allData[in+1][iHit] < high[in])
38cbfa7c 357 {
85f61e3b 358 timeCFD[in] = allData[in+1][iHit] ;
359 break;
38cbfa7c 360 }
85f61e3b 361 }
362 for (Int_t iHit=0; iHit<5; iHit++)
363 {
612737bb 364 if(allData[in+1+56][iHit] > low[in] &&
365 allData[in+1+56][iHit] < high[in])
38cbfa7c 366 {
85f61e3b 367 timeCFD[in+12] = allData[in+56+1][iHit] ;
368 break;
38cbfa7c 369 }
38cbfa7c 370 }
612737bb 371 timeLED[in+12] = allData[in+68+1][0] ;
372 timeLED[in] = allData[in+12+1][0] ;
669dc07f 373 AliDebug(5, Form(" readed i %i cfdC %i cfdA %i ledC %i ledA%i ",
85f61e3b 374 in, timeCFD[in],timeCFD[in+12],timeLED[in],
612737bb 375 timeLED[in+12]));
38cbfa7c 376
8f620945 377 }
378
8f620945 379
85f61e3b 380 for (Int_t in=0; in<12; in++)
381 {
382 chargeQT0[in]=allData[2*in+25][0];
383 chargeQT1[in]=allData[2*in+26][0];
612737bb 384 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 385 in, chargeQT0[in],chargeQT1[in]));
386 }
387 for (Int_t in=12; in<24; in++)
388 {
389 chargeQT0[in]=allData[2*in+57][0];
390 chargeQT1[in]=allData[2*in+58][0];
612737bb 391 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 392 in, chargeQT0[in],chargeQT1[in]));
393
394 }
395
8f226345 396 onlineMean = allData[49][0];
397
195e1353 398 Double32_t time[24], adc[24], adcmip[24], noncalibtime[24];
8f620945 399 for (Int_t ipmt=0; ipmt<24; ipmt++) {
36bfca7d 400 if(timeCFD[ipmt] > 0 /* && badpmt[ipmt]==0*/ ){
bce12dc5 401 //for simulated data
d0bcd1fb 402 //for physics data
73df58ab 403 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>0) {
541b42c4 404 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 405 }
d0bcd1fb 406 else
407 adc[ipmt] = 0;
612737bb 408 // time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
bce12dc5 409
612737bb 410 time[ipmt] = fCalib-> WalkCorrection(Int_t (fTime0vertex[ipmt]), ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
d0bcd1fb 411 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 412 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 413 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 414 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
415 Double_t ampMip =( (TGraph*)fAmpLED.At(ipmt))->Eval(sl);
416 Double_t qtMip = ((TGraph*)fQTC.At(ipmt))->Eval(adc[ipmt]);
345f03db 417 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %f\n ",ampMip,qtMip, adc[ipmt]));
73df58ab 418 //bad peak removing
73df58ab 419 frecpoints->SetTime(ipmt, Float_t(time[ipmt]) );
420 // frecpoints->SetTime(ipmt,Double32_t(timeCFD[ipmt]));
195e1353 421 frecpoints->SetAmp(ipmt, Double32_t( qtMip));
422 adcmip[ipmt]=qtMip;
345f03db 423 frecpoints->SetAmpLED(ipmt, Double32_t(ampMip));
73df58ab 424 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
2e6a5ee0 425 }
426 else {
427 time[ipmt] = 0;
428 adc[ipmt] = 0;
b0ab3f59 429 adcmip[ipmt] = 0;
73df58ab 430 noncalibtime[ipmt] = 0;
2e6a5ee0 431 }
432 }
73df58ab 433 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
2e6a5ee0 434 for (Int_t ipmt=0; ipmt<12; ipmt++){
36bfca7d 435 if(time[ipmt] !=0 /*&& badpmt[ipmt]==0 */&& adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
73df58ab 436 {
b0e13b29 437 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC)) {
58432641 438 if(time[ipmt]<besttimeC){
612737bb 439 besttimeC=time[ipmt]; //timeC
669dc07f 440 pmtBestC=ipmt;
73df58ab 441 }
2e6a5ee0 442 }
2e6a5ee0 443 }
73df58ab 444 for ( Int_t ipmt=12; ipmt<24; ipmt++)
445 {
36bfca7d 446 if(time[ipmt] != 0 /* && badpmt[ipmt]==0*/ && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
73df58ab 447 {
58432641 448 if(time[ipmt]<besttimeA) {
449 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeA)) {
73df58ab 450 besttimeA=time[ipmt]; //timeA
451 pmtBestA=ipmt;
452 }
453 }
2e6a5ee0 454 }
adf36b9d 455 if(besttimeA < 999999)
58432641 456 frecpoints->SetTimeBestA((besttimeA * channelWidth)- 1000.*fLatencyHPTDC + 1000.*fLatencyL1A - 1000.*fGRPdelays - fTimeMeanShift[1] );
457 // frecpoints->SetTimeBestA((besttimeA * channelWidth- fTimeMeanShift[1]));
669dc07f 458
adf36b9d 459 if( besttimeC < 999999 )
58432641 460 frecpoints->SetTimeBestC((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1C - 1000.*fGRPdelays - fTimeMeanShift[2]);
461 // frecpoints->SetTimeBestC((besttimeC * channelWidth - fTimeMeanShift[2]));
669dc07f 462 AliDebug(5,Form(" pmtA %i besttimeA %f shift A %f ps, pmtC %i besttimeC %f shiftC %f ps",
463 pmtBestA,besttimeA, fTimeMeanShift[1],
464 pmtBestC, besttimeC,fTimeMeanShift[2]));
adf36b9d 465 if(besttimeA <999999 && besttimeC < 999999 ){
612737bb 466 // timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth + fLatencyL1A - fLatencyL1C;
58432641 467 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - 1000.*fLatencyHPTDC + 1000.*fLatencyL1 - 1000.*fGRPdelays - fTimeMeanShift[0] ;
adf36b9d 468 meanTime = (besttimeA+besttimeC-2.*Float_t(ref))/2.;
612737bb 469 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth ;
58432641 470 // timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - fTimeMeanShift[0] ;
9480f05f 471 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
adf36b9d 472 }
776de217 473 } //if phys event
291f31a1 474 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 475 frecpoints->SetT0clock(timeclock);
476 frecpoints->SetVertex(vertex);
adf36b9d 477 frecpoints->SetMeanTime(meanTime);
776de217 478 frecpoints->SetOnlineMean(Int_t(onlineMean));
adf36b9d 479 // Set triggers
480
481 Bool_t tr[5];
482 Int_t trchan[5]= {50,51,52,55,56};
483 for (Int_t i=0; i<5; i++) tr[i]=false;
484 for (Int_t itr=0; itr<5; itr++) {
b0e13b29 485 for (Int_t iHit=0; iHit<1; iHit++)
38cbfa7c 486 {
487 Int_t trr=trchan[itr];
b0e13b29 488 if( allData[trr][iHit] > 0) tr[itr]=true;
489
490 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]));
491 }
adf36b9d 492 }
38cbfa7c 493 frecpoints->SetT0Trig(tr);
b0e13b29 494
495 for (Int_t iHit=0; iHit<5; iHit++)
496 {
497 Float_t tvdc = -9999; Float_t ora = -9999; Float_t orc = -9999;
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);
509 Float_t timefull;
510 for (Int_t i0=0; i0<12; i0++) {
511 timefull = -9999;
512 if(allData[i0+1][iHit]>1)
513 timefull = (Float_t(allData[i0+1][iHit])-fTime0vertex[i0])* channelWidth* 0.001;
514 frecpoints->SetTimeFull(i0, iHit,timefull) ;
515
516 }
517
518 for (Int_t i0=12; i0<24; i0++) {
519 timefull = -9999;
520 if(allData[i0+45][iHit]>1)
521 timefull = (Float_t(allData[i0+45][iHit])-fTime0vertex[i0])* channelWidth* 0.001;
522 frecpoints->SetTimeFull(i0, iHit, timefull) ;
523 }
524 }
525
526
669dc07f 527 //Set MPD
528 if(allData[53][0]>0 && allData[54][0])
529 frecpoints->SetMultA(allData[53][0]-allData[54][0]);
b0e13b29 530 if(allData[105][0]>0 && allData[106][0])
531 frecpoints->SetMultC(allData[105][0]-allData[106][0]);
532
533
669dc07f 534 } // if (else )raw data
58bd3a16 535 recTree->Fill();
536 if(frecpoints) delete frecpoints;
537}
adf36b9d 538
539
540 //____________________________________________________________
541
542 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
543 {
dc7ca31d 544
545 /***************************************************
546 Resonstruct digits to vertex position
547 ****************************************************/
548
dc7ca31d 549 AliDebug(1,Form("Start FillESD T0"));
b0ab3f59 550 if(!pESD) {
551 AliError("No ESD Event");
552 return;
553 }
4cbe597e 554 pESD ->SetT0spread(fTimeSigmaShift);
555
36cde487 556
58bd3a16 557 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 558 Float_t c = 0.0299792458; // cm/ps
adf36b9d 559 Float_t currentVertex=0, shift=0;
291f31a1 560 Int_t ncont=-1;
adf36b9d 561 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
562 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
563 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
564 if (!vertex) vertex = pESD->GetVertex();
565
566 if (vertex) {
567 AliDebug(2, Form("Got %s (%s) from ESD: %f",
568 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
569 currentVertex = vertex->GetZ();
570
571 ncont = vertex->GetNContributors();
291f31a1 572 if(ncont>0 ) {
adf36b9d 573 shift = currentVertex/c;
adf36b9d 574 }
575 }
d76c31f4 576 TTree *treeR = clustersTree;
dc7ca31d 577
73df58ab 578 AliT0RecPoint* frecpoints= new AliT0RecPoint ();
579 if (!frecpoints) {
dc7ca31d 580 AliError("Reconstruct Fill ESD >> no recpoints found");
581 return;
582 }
583
584 AliDebug(1,Form("Start FillESD T0"));
585 TBranch *brRec = treeR->GetBranch("T0");
586 if (brRec) {
587 brRec->SetAddress(&frecpoints);
588 }else{
f16935f7 589 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 590 return;
591 }
73df58ab 592
593 brRec->GetEntry(0);
594 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
595 Double32_t* tcorr;
596 for(Int_t i=0; i<24; i++)
597 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
598
669dc07f 599
73df58ab 600 Double32_t timeClock[3];
601 Double32_t zPosition = frecpoints -> GetVertex();
602 Double32_t timeStart = frecpoints -> GetMeanTime();
603 timeClock[0] = frecpoints -> GetT0clock() ;
604 timeClock[1] = frecpoints -> GetBestTimeA() + shift;
605 timeClock[2] = frecpoints -> GetBestTimeC() - shift;
612737bb 606
73df58ab 607 for ( Int_t i=0; i<24; i++) {
608 time[i] = frecpoints -> GetTime(i); // ps to ns
612737bb 609 // if ( time[i] >1) {
610 if ( time[i] != 0) {
345f03db 611 ampQTC[i] = frecpoints -> GetAmp(i);
612 amp[i] = frecpoints -> AmpLED(i);
612737bb 613 AliDebug(1,Form("T0: %i time %f ampQTC %f ampLED %f \n", i, time[i], ampQTC[i], amp[i]));
85f61e3b 614 }
73df58ab 615 }
616 Int_t trig= frecpoints ->GetT0Trig();
b0e13b29 617 frecpoints->PrintTriggerSignals( trig);
618 printf(" FillESD trigger %i \n",trig);
619 fESDTZERO->SetT0Trig(trig);
620 //pESD->SetT0Trig(trig);
621 // pESD->SetT0zVertex(zPosition); //vertex Z position
622 fESDTZERO->SetT0zVertex(zPosition); //vertex Z position
669dc07f 623
624 Double32_t multA=frecpoints ->GetMultA();
625 Double32_t multC=frecpoints ->GetMultC();
b0e13b29 626 // pESD->SetT0(multC); // multiplicity Cside
627 // pESD->SetT0clock(multA); // multiplicity Aside
628 fESDTZERO->SetMultA(multA); // for backward compatubility
629 fESDTZERO->SetMultC(multC); // for backward compatubility
630
631
632 for (Int_t iHit =0; iHit<5; iHit++ ) {
633 AliDebug(1,Form("FillESD ::: iHit %i tvdc %f orA %f orC %f\n", iHit,
634 frecpoints->GetTVDC(iHit),
635 frecpoints->GetOrA(iHit),
636 frecpoints->GetOrC(iHit) ));
637 fESDTZERO->SetTVDC(iHit,frecpoints->GetTVDC(iHit));
638 fESDTZERO->SetOrA(iHit,frecpoints->GetOrA(iHit));
639 fESDTZERO->SetOrC(iHit,frecpoints->GetOrC(iHit));
640
641 for (Int_t i0=0; i0<24; i0++) {
642 if(frecpoints->GetTimeFull(i0,iHit)>0){
643 // printf("FillESD ::: iHit %i cfd %i time %f \n", iHit, i0,
644 // frecpoints->GetTimeFull(i0,iHit));
645 fESDTZERO->SetTimeFull(i0, iHit,frecpoints->GetTimeFull(i0,iHit));
646 }
647
648 }
649 }
73df58ab 650 for(Int_t i=0; i<3; i++)
b0e13b29 651 fESDTZERO->SetT0TOF(i,timeClock[i]); // interaction time (ns)
652 fESDTZERO->SetT0time(time); // best TOF on each PMT
653 fESDTZERO->SetT0amplitude(ampQTC); // number of particles(MIPs) on each PMT
73df58ab 654
612737bb 655 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 656
73df58ab 657
b0e13b29 658
73df58ab 659 if (pESD) {
adf36b9d 660
73df58ab 661 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
662 if (fr) {
36cde487 663 AliDebug(1, Form("Writing TZERO friend data to ESD tree"));
73df58ab 664
85f61e3b 665 // if (ncont>2) {
73df58ab 666 tcorr = fESDTZEROfriend->GetT0timeCorr();
667 for ( Int_t i=0; i<24; i++) {
85f61e3b 668 if(i<12 && time[i]>1) timecorr[i] = tcorr[i] - shift/channelWidth;
669 if(i>11 && time[i]>1) timecorr[i] = tcorr[i] + shift/channelWidth;
612737bb 670 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 671 }
73df58ab 672 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
673 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
674 fESDTZEROfriend->SetT0ampQTC(ampQTC);
85f61e3b 675 fr->SetTZEROfriend(fESDTZEROfriend);
676 // }//
58bd3a16 677 }
b0e13b29 678
679 pESD->SetTZEROData(fESDTZERO);
58bd3a16 680 }
73df58ab 681
682
dc7ca31d 683} // vertex in 3 sigma
684
685
686
687
688
689