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