Reverting the changes, we need a consistent commit in ESD and T0
[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;
92f6bfd3 120 // printf(" calulated mean %i %f \n",i, fTime0vertex[i]);
612737bb 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
830172c8 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() ;
830172c8 146 Float_t meanVertex = fParam->GetMeanVertex();
147 Float_t c = 0.0299792; // cm/ps
148 Double32_t vertex = 9999999;
adf36b9d 149 Double32_t timeDiff=999999, meanTime=999999, timeclock=999999;
830172c8 150
94c27e4f 151
dc7ca31d 152 AliDebug(1,Form("Start DIGITS reconstruction "));
94c27e4f 153
830172c8 154 Float_t lowAmpThreshold = GetRecoParam()->GetLow(200);
155 Float_t highAmpThreshold = GetRecoParam()->GetHigh(200);
156 Int_t badpmt = GetRecoParam()->GetRefPoint();
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"));
830172c8 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();
830172c8 175
adf36b9d 176 Bool_t tr[5];
177 for (Int_t i=0; i<5; i++) tr[i]=false;
c41ceaac 178
830172c8 179 Double32_t besttimeA=999999;
180 Double32_t besttimeC=999999;
181 Int_t pmtBestA=99999;
182 Int_t pmtBestC=99999;
dc7ca31d 183
7e08771f 184 AliT0RecPoint frecpoints;
185 AliT0RecPoint * pfrecpoints = &frecpoints;
186 clustersTree->Branch( "T0", "AliT0RecPoint" ,&pfrecpoints);
94c27e4f 187
195e1353 188 Float_t time[24], adc[24], adcmip[24];
dc7ca31d 189 for (Int_t ipmt=0; ipmt<24; ipmt++) {
830172c8 190 if(timeCFD->At(ipmt)>0 && ipmt != badpmt) {
191 if(( chargeQT1->At(ipmt) - chargeQT0->At(ipmt))>0)
d0bcd1fb 192 adc[ipmt] = chargeQT1->At(ipmt) - chargeQT0->At(ipmt);
a7027400 193 else
d0bcd1fb 194 adc[ipmt] = 0;
195
830172c8 196 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD->At(ipmt)) ;
197
d0bcd1fb 198 Double_t sl = Double_t(timeLED->At(ipmt) - timeCFD->At(ipmt));
8f620945 199 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD->At(ipmt) ) ;
669dc07f 200 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
830172c8 201 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
202
644a358e 203 Double_t ampMip = 0;
204 TGraph* ampGraph = (TGraph*)fAmpLED.At(ipmt);
205 if (ampGraph) ampMip = ampGraph->Eval(sl);
206 Double_t qtMip = 0;
207 TGraph* qtGraph = (TGraph*)fQTC.At(ipmt);
208 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
669dc07f 209 AliDebug(5,Form(" Amlitude in MIPS LED %f , QTC %f in channels %f\n ",ampMip,qtMip, adc[ipmt]));
830172c8 210
7e08771f 211 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
212 frecpoints.SetAmpLED(ipmt, Float_t( ampMip));
213 frecpoints.SetAmp(ipmt, Float_t(qtMip));
195e1353 214 adcmip[ipmt]=qtMip;
d0bcd1fb 215
dc7ca31d 216 }
217 else {
218 time[ipmt] = 0;
219 adc[ipmt] = 0;
b0ab3f59 220 adcmip[ipmt] = 0;
830172c8 221
dc7ca31d 222 }
223 }
94c27e4f 224
dc7ca31d 225 for (Int_t ipmt=0; ipmt<12; ipmt++){
830172c8 226 if(time[ipmt] > 1 && ipmt != badpmt && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold) {
227 if(time[ipmt]<besttimeC){
228 besttimeC=time[ipmt]; //timeC
229 pmtBestC=ipmt;
dc7ca31d 230 }
830172c8 231 }
dc7ca31d 232 }
830172c8 233 for ( Int_t ipmt=12; ipmt<24; ipmt++){
234 if(time[ipmt] > 1 && ipmt != badpmt && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold) {
235 if(time[ipmt]<besttimeA) {
236 besttimeA=time[ipmt]; //timeA
237 pmtBestA=ipmt;}
dc7ca31d 238 }
830172c8 239 }
240 if(besttimeA < 999999) {
241 frecpoints.SetTimeBestA(Int_t(besttimeA *channelWidth - fdZonA/c));
adf36b9d 242 tr[1]=true;
243 }
830172c8 244 if( besttimeC < 999999 ) {
245 frecpoints.SetTimeBestC(Int_t(besttimeC *channelWidth - fdZonA/c));
adf36b9d 246 tr[2]=true;
247 }
830172c8 248 AliDebug(5,Form(" besttimeA %f ch, besttimeC %f ch",besttimeA, besttimeC));
adf36b9d 249 if(besttimeA <999999 && besttimeC < 999999 ){
9b83615d 250 // timeDiff = (besttimeC - besttimeA)*channelWidth;
251 timeDiff = (besttimeA - besttimeC)*channelWidth;
830172c8 252 meanTime = (besttimeA + besttimeC)/2;// * channelWidth);
253 timeclock = meanTime *channelWidth -fdZonA/c ;
adf36b9d 254 vertex = meanVertex - c*(timeDiff)/2.;// + (fdZonA - fdZonC)/2;
255 tr[0]=true;
256 }
7e08771f 257 frecpoints.SetVertex(vertex);
830172c8 258 frecpoints.SetMeanTime(meanTime);
259 frecpoints.SetT0clock(timeclock);
7e08771f 260 frecpoints.SetT0Trig(tr);
830172c8 261
669dc07f 262 AliDebug(5,Form("T0 triggers %d %d %d %d %d",tr[0],tr[1],tr[2],tr[3],tr[4]));
830172c8 263
adf36b9d 264 //online mean
7e08771f 265 frecpoints.SetOnlineMean(Int_t(onlineMean));
b5a9f753 266 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 267
0518959c 268
830172c8 269
270
271
272 clustersTree->Fill();
273
bd375212 274 delete timeCFD;
275 delete timeLED;
276 delete chargeQT0;
277 delete chargeQT1;
dc7ca31d 278}
279
280
c41ceaac 281//_______________________________________________________________________
282
283void AliT0Reconstructor::Reconstruct(AliRawReader* rawReader, TTree*recTree) const
284{
94c27e4f 285 // T0 raw ->
539b9cb9 286 //
b0e13b29 287
288 Float_t meanOrA = fTime0vertex[0] + 587;
289 Float_t meanOrC = fTime0vertex[0] + 678;
290 Float_t meanTVDC = fTime0vertex[0] + 2564;
82540e64 291 Int_t timeDelayCFD[24];
292
d3e04608 293
612737bb 294 Int_t badpmt[24];
38cbfa7c 295 //Bad channel
58432641 296 for (Int_t i=0; i<24; i++) {
297 badpmt[i] = GetRecoParam() -> GetBadChannels(i);
82540e64 298 timeDelayCFD[i] = Int_t (fParam->GetTimeDelayCFD(i));
58432641 299 }
d3e04608 300 Int_t low[500], high[500];
94b4f7e2 301 Float_t timefull=-9999;;
302 Float_t tvdc = -9999; Float_t ora = -9999; Float_t orc = -9999;
58bd3a16 303
e8ed1cd0 304 Int_t allData[110][5];
2e6a5ee0 305
e8ed1cd0 306 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
612737bb 307 Double32_t timeDiff, meanTime, timeclock;
308 timeDiff = meanTime = timeclock = 9999999;
adf36b9d 309 Float_t c = 29.9792458; // cm/ns
310 Double32_t vertex = 9999999;
776de217 311 Int_t onlineMean=0;
9480f05f 312 Float_t meanVertex = 0;
612737bb 313 for (Int_t i0=0; i0<24; i0++) {
314 low[i0] = Int_t(fTime0vertex[i0]) - 200;
315 high[i0] = Int_t(fTime0vertex[i0]) + 200;
316 }
317
669dc07f 318 for (Int_t i0=0; i0<110; i0++)
b0e13b29 319 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
612737bb 320
321 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
322 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
830172c8 323
324 Double32_t besttimeA=9999999;
325 Double32_t besttimeC=9999999;
326 Int_t pmtBestA=99999;
327 Int_t pmtBestC=99999;
328 Float_t channelWidth = fParam->GetChannelWidth() ;
b0e13b29 329
7e08771f 330 AliT0RecPoint frecpoints;
331 AliT0RecPoint * pfrecpoints = &frecpoints;
bce12dc5 332
7e08771f 333 recTree->Branch( "T0", "AliT0RecPoint" ,&pfrecpoints);
2e6a5ee0 334
bce12dc5 335 AliDebug(10," before read data ");
336 AliT0RawReader myrawreader(rawReader);
776de217 337
338 UInt_t type =rawReader->GetType();
339
bce12dc5 340 if (!myrawreader.Next())
341 AliDebug(1,Form(" no raw data found!!"));
342 else
343 {
38cbfa7c 344 for (Int_t i=0; i<24; i++)
345 {
346 timeCFD[i]=0; timeLED[i]=0; chargeQT0[i]=0; chargeQT1[i]=0;
347 }
348 Int_t fBCID=Int_t (rawReader->GetBCID());
8f620945 349 Int_t trmbunch= myrawreader.GetTRMBunchID();
86fd0587 350 AliDebug(10,Form(" CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
351
38cbfa7c 352 if(type == 7 ) { //only physics
669dc07f 353 for (Int_t i=0; i<107; i++) {
bce12dc5 354 for (Int_t iHit=0; iHit<5; iHit++)
355 {
356 allData[i][iHit] = myrawreader.GetData(i,iHit);
357 }
8f620945 358 }
830172c8 359 Int_t ref=0;
8f620945 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] ;
830172c8 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]));
85f61e3b 403 }
404
8f226345 405 onlineMean = allData[49][0];
406
195e1353 407 Double32_t time[24], adc[24], adcmip[24], noncalibtime[24];
8f620945 408 for (Int_t ipmt=0; ipmt<24; ipmt++) {
36bfca7d 409 if(timeCFD[ipmt] > 0 /* && badpmt[ipmt]==0*/ ){
bce12dc5 410 //for simulated data
d0bcd1fb 411 //for physics data
73df58ab 412 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>0) {
541b42c4 413 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 414 }
d0bcd1fb 415 else
416 adc[ipmt] = 0;
612737bb 417 // time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
bce12dc5 418
612737bb 419 time[ipmt] = fCalib-> WalkCorrection(Int_t (fTime0vertex[ipmt]), ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
d0bcd1fb 420 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 421 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 422 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 423 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
644a358e 424 Double_t ampMip = 0;
425 TGraph * ampGraph = (TGraph*)fAmpLED.At(ipmt);
426 if (ampGraph) ampMip = ampGraph->Eval(sl);
427 Double_t qtMip = 0;
428 TGraph * qtGraph = (TGraph*)fQTC.At(ipmt);
429 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
345f03db 430 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %f\n ",ampMip,qtMip, adc[ipmt]));
830172c8 431 //bad peak removing
7e08771f 432 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
830172c8 433 // frecpoints.SetTime(ipmt,Double32_t(timeCFD[ipmt]));
434 frecpoints.SetAmp(ipmt, Double32_t( qtMip));
435 adcmip[ipmt]=qtMip;
436 frecpoints.SetAmpLED(ipmt, Double32_t(ampMip));
437 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
2e6a5ee0 438 }
439 else {
830172c8 440 time[ipmt] = 0;
2e6a5ee0 441 adc[ipmt] = 0;
b0ab3f59 442 adcmip[ipmt] = 0;
830172c8 443 noncalibtime[ipmt] = 0;
2e6a5ee0 444 }
445 }
830172c8 446 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
2e6a5ee0 447 for (Int_t ipmt=0; ipmt<12; ipmt++){
830172c8 448 if(time[ipmt] !=0 /*&& badpmt[ipmt]==0 */&& adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
73df58ab 449 {
830172c8 450 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC)) {
451 if(time[ipmt]<besttimeC){
452 besttimeC=time[ipmt]; //timeC
453 pmtBestC=ipmt;
454 }
2e6a5ee0 455 }
2e6a5ee0 456 }
73df58ab 457 for ( Int_t ipmt=12; ipmt<24; ipmt++)
458 {
830172c8 459 if(time[ipmt] != 0 /* && badpmt[ipmt]==0*/ && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
73df58ab 460 {
830172c8 461 if(time[ipmt]<besttimeA) {
462 // if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeA)) {
463 besttimeA=time[ipmt]; //timeA
464 pmtBestA=ipmt;
465 }
73df58ab 466 }
2e6a5ee0 467 }
830172c8 468 if(besttimeA < 999999)
469 frecpoints.SetTimeBestA((besttimeA * channelWidth)- 1000.*fLatencyHPTDC + 1000.*fLatencyL1A - 1000.*fGRPdelays - fTimeMeanShift[1] );
470 // frecpoints.SetTimeBestA((besttimeA * channelWidth- fTimeMeanShift[1]));
471
472 if( besttimeC < 999999 )
473 frecpoints.SetTimeBestC((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1C - 1000.*fGRPdelays - fTimeMeanShift[2]);
474 // frecpoints.SetTimeBestC((besttimeC * channelWidth - fTimeMeanShift[2]));
475 AliDebug(5,Form(" pmtA %i besttimeA %f shift A %f ps, pmtC %i besttimeC %f shiftC %f ps",
476 pmtBestA,besttimeA, fTimeMeanShift[1],
477 pmtBestC, besttimeC,fTimeMeanShift[2]));
478 if(besttimeA <999999 && besttimeC < 999999 ){
479 // timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth + fLatencyL1A - fLatencyL1C;
480 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - 1000.*fLatencyHPTDC + 1000.*fLatencyL1 - 1000.*fGRPdelays - fTimeMeanShift[0] ;
481 meanTime = (besttimeA+besttimeC-2.*Float_t(ref))/2.;
612737bb 482 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth ;
830172c8 483 // timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - fTimeMeanShift[0] ;
9480f05f 484 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
830172c8 485 }
776de217 486 } //if phys event
830172c8 487 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));
488 frecpoints.SetT0clock(timeclock);
489 frecpoints.SetVertex(vertex);
490 frecpoints.SetMeanTime(meanTime);
491 frecpoints.SetOnlineMean(Int_t(onlineMean));
492 // Set triggers
adf36b9d 493
494 Bool_t tr[5];
92f6bfd3 495 Int_t trchan[5] = {50,51,52,55,56};
496 Float_t lowtr[5] = {meanTVDC-700, meanOrA-700, meanOrC-700, meanOrC-1000, meanOrC-1000 };
497 Float_t hightr[5] = {meanTVDC+700, meanOrA+700, meanOrC+700, meanOrC+1000, meanOrC+1000};
830172c8 498
82540e64 499 for (Int_t i=0; i<5; i++) tr[i] = false;
adf36b9d 500 for (Int_t itr=0; itr<5; itr++) {
b0e13b29 501 for (Int_t iHit=0; iHit<1; iHit++)
38cbfa7c 502 {
503 Int_t trr=trchan[itr];
92f6bfd3 504 if( allData[trr][iHit] > lowtr[itr] && allData[trr][iHit] < hightr[itr]) tr[itr]=true;
830172c8 505
506 AliDebug(5,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]));
b0e13b29 507 }
adf36b9d 508 }
7e08771f 509 frecpoints.SetT0Trig(tr);
b0e13b29 510
511 for (Int_t iHit=0; iHit<5; iHit++)
512 {
830172c8 513 timefull = -9999;
82540e64 514 tvdc = ora = orc = -9999;
b0e13b29 515 if(allData[50][iHit]>0)
516 tvdc = (Float_t(allData[50][iHit]) - meanTVDC) * channelWidth* 0.001;
517 if(allData[51][iHit]>0)
518 ora = (Float_t(allData[51][iHit]) - meanOrA) * channelWidth* 0.001;
519
520 if(allData[52][iHit]>0)
521 orc = (Float_t(allData[52][iHit]) - meanOrC) * channelWidth* 0.001;
522
7e08771f 523 frecpoints.SetOrC( iHit, orc);
524 frecpoints.SetOrA( iHit, ora);
525 frecpoints.SetTVDC( iHit, tvdc);
b0e13b29 526 for (Int_t i0=0; i0<12; i0++) {
527 timefull = -9999;
528 if(allData[i0+1][iHit]>1)
830172c8 529 timefull = (Float_t(allData[i0+1][iHit])-fTime0vertex[i0] - timeDelayCFD[i0])* channelWidth* 0.001;
7e08771f 530 frecpoints.SetTimeFull(i0, iHit,timefull) ;
82540e64 531 // if(allData[i0+1][iHit]>1) printf("i0 %d iHit %d data %d fTime0vertex %f timefull %f \n",i0, iHit, allData[i0+1][iHit], fTime0vertex[i0], timefull);
b0e13b29 532
533 }
534
535 for (Int_t i0=12; i0<24; i0++) {
536 timefull = -9999;
94b4f7e2 537 if(allData[i0+45][iHit]>1) {
830172c8 538 timefull = (Float_t(allData[i0+45][iHit]) - fTime0vertex[i0] - timeDelayCFD[i0])* channelWidth* 0.001;
94b4f7e2 539 }
82540e64 540 // if(allData[i0+45][iHit]>1) printf("i0 %d iHit %d data %d fTime0vertex %f timefull %f \n",i0, iHit, allData[i0+45][iHit], fTime0vertex[i0], timefull);
7e08771f 541 frecpoints.SetTimeFull(i0, iHit, timefull) ;
b0e13b29 542 }
543 }
544
545
669dc07f 546 //Set MPD
547 if(allData[53][0]>0 && allData[54][0])
7e08771f 548 frecpoints.SetMultA(allData[53][0]-allData[54][0]);
b0e13b29 549 if(allData[105][0]>0 && allData[106][0])
7e08771f 550 frecpoints.SetMultC(allData[105][0]-allData[106][0]);
b0e13b29 551
552
669dc07f 553 } // if (else )raw data
58bd3a16 554 recTree->Fill();
58bd3a16 555}
adf36b9d 556
557
830172c8 558 //____________________________________________________________
adf36b9d 559
560 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
561 {
dc7ca31d 562
563 /***************************************************
564 Resonstruct digits to vertex position
565 ****************************************************/
566
dc7ca31d 567 AliDebug(1,Form("Start FillESD T0"));
b0ab3f59 568 if(!pESD) {
569 AliError("No ESD Event");
570 return;
571 }
4cbe597e 572 pESD ->SetT0spread(fTimeSigmaShift);
573
36cde487 574
58bd3a16 575 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 576 Float_t c = 0.0299792458; // cm/ps
adf36b9d 577 Float_t currentVertex=0, shift=0;
291f31a1 578 Int_t ncont=-1;
adf36b9d 579 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
580 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
581 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
582 if (!vertex) vertex = pESD->GetVertex();
583
584 if (vertex) {
585 AliDebug(2, Form("Got %s (%s) from ESD: %f",
586 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
587 currentVertex = vertex->GetZ();
588
589 ncont = vertex->GetNContributors();
291f31a1 590 if(ncont>0 ) {
adf36b9d 591 shift = currentVertex/c;
adf36b9d 592 }
593 }
d76c31f4 594 TTree *treeR = clustersTree;
dc7ca31d 595
7e08771f 596 AliT0RecPoint frecpoints;
597 AliT0RecPoint * pfrecpoints = &frecpoints;
dc7ca31d 598
599 AliDebug(1,Form("Start FillESD T0"));
600 TBranch *brRec = treeR->GetBranch("T0");
601 if (brRec) {
7e08771f 602 brRec->SetAddress(&pfrecpoints);
dc7ca31d 603 }else{
f16935f7 604 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 605 return;
606 }
73df58ab 607
608 brRec->GetEntry(0);
609 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
610 Double32_t* tcorr;
611 for(Int_t i=0; i<24; i++)
612 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
613
830172c8 614
73df58ab 615 Double32_t timeClock[3];
7e08771f 616 Double32_t zPosition = frecpoints.GetVertex();
830172c8 617 Double32_t timeStart = frecpoints.GetMeanTime();
0518959c 618 timeClock[0] = frecpoints.GetT0clock() ;
830172c8 619 timeClock[1] = frecpoints.GetBestTimeA() + shift;
620 timeClock[2] = frecpoints.GetBestTimeC() - shift;
621
73df58ab 622 for ( Int_t i=0; i<24; i++) {
7e08771f 623 time[i] = frecpoints.GetTime(i); // ps to ns
830172c8 624 // if ( time[i] >1) {
625 if ( time[i] != 0) {
7e08771f 626 ampQTC[i] = frecpoints.GetAmp(i);
627 amp[i] = frecpoints.AmpLED(i);
612737bb 628 AliDebug(1,Form("T0: %i time %f ampQTC %f ampLED %f \n", i, time[i], ampQTC[i], amp[i]));
85f61e3b 629 }
73df58ab 630 }
7e08771f 631 Int_t trig= frecpoints.GetT0Trig();
632 frecpoints.PrintTriggerSignals( trig);
92f6bfd3 633 // printf(" !!!!! FillESD trigger %i \n",trig);
b0e13b29 634 fESDTZERO->SetT0Trig(trig);
830172c8 635 //pESD->SetT0Trig(trig);
636 // pESD->SetT0zVertex(zPosition); //vertex Z position
b0e13b29 637 fESDTZERO->SetT0zVertex(zPosition); //vertex Z position
669dc07f 638
7e08771f 639 Double32_t multA=frecpoints.GetMultA();
640 Double32_t multC=frecpoints.GetMultC();
830172c8 641 // pESD->SetT0(multC); // multiplicity Cside
642 // pESD->SetT0clock(multA); // multiplicity Aside
b0e13b29 643 fESDTZERO->SetMultA(multA); // for backward compatubility
644 fESDTZERO->SetMultC(multC); // for backward compatubility
645
646
647 for (Int_t iHit =0; iHit<5; iHit++ ) {
82540e64 648 AliDebug(10,Form("FillESD ::: iHit %i tvdc %f orA %f orC %f\n", iHit,
7e08771f 649 frecpoints.GetTVDC(iHit),
650 frecpoints.GetOrA(iHit),
651 frecpoints.GetOrC(iHit) ));
652 fESDTZERO->SetTVDC(iHit,frecpoints.GetTVDC(iHit));
653 fESDTZERO->SetOrA(iHit,frecpoints.GetOrA(iHit));
654 fESDTZERO->SetOrC(iHit,frecpoints.GetOrC(iHit));
b0e13b29 655
656 for (Int_t i0=0; i0<24; i0++) {
7e08771f 657 // if(frecpoints.GetTimeFull(i0,iHit)>0){
658 // printf("FillESD ::: iHit %i cfd %i time %f \n", iHit, i0, frecpoints.GetTimeFull(i0,iHit));
659 fESDTZERO->SetTimeFull(i0, iHit,frecpoints.GetTimeFull(i0,iHit));
94b4f7e2 660 // }
b0e13b29 661
662 }
663 }
830172c8 664 for(Int_t i=0; i<3; i++)
665 fESDTZERO->SetT0TOF(i,timeClock[i]); // interaction time (ns)
666 fESDTZERO->SetT0time(time); // best TOF on each PMT
667 fESDTZERO->SetT0amplitude(ampQTC); // number of particles(MIPs) on each PMT
73df58ab 668
830172c8 669 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 670
ea1a8005 671 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
672 // background flags
673 Bool_t background = BackgroundFlag();
674 fESDTZERO->SetBackgroundFlag(background);
675 Bool_t pileup = PileupFlag();
676 fESDTZERO->SetPileupFlag(pileup);
528890e5 677 for (Int_t i=0; i<5; i++) {
7e08771f 678 fESDTZERO->SetPileupTime(i, frecpoints.GetTVDC(i) ) ;
679 // printf("!!!!!! FillESD :: pileup %i %f %f \n", i,fESDTZERO->GetPileupTime(i), frecpoints.GetTVDC(i));
528890e5 680 }
ea1a8005 681 Bool_t sat = SatelliteFlag();
682 fESDTZERO->SetSatelliteFlag(sat);
73df58ab 683
ea1a8005 684 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
73df58ab 685 if (pESD) {
adf36b9d 686
73df58ab 687 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
688 if (fr) {
82540e64 689 AliDebug(10, Form("Writing TZERO friend data to ESD tree"));
73df58ab 690
85f61e3b 691 // if (ncont>2) {
73df58ab 692 tcorr = fESDTZEROfriend->GetT0timeCorr();
693 for ( Int_t i=0; i<24; i++) {
85f61e3b 694 if(i<12 && time[i]>1) timecorr[i] = tcorr[i] - shift/channelWidth;
695 if(i>11 && time[i]>1) timecorr[i] = tcorr[i] + shift/channelWidth;
612737bb 696 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 697 }
73df58ab 698 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
699 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
700 fESDTZEROfriend->SetT0ampQTC(ampQTC);
85f61e3b 701 fr->SetTZEROfriend(fESDTZEROfriend);
702 // }//
58bd3a16 703 }
b0e13b29 704
705 pESD->SetTZEROData(fESDTZERO);
58bd3a16 706 }
73df58ab 707
dc7ca31d 708} // vertex in 3 sigma
709
ea1a8005 710//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
711 //____________________________________________________________
712
713Bool_t AliT0Reconstructor::PileupFlag() const
714{
715 //
716 Bool_t pileup = false;
528890e5 717 Float_t tvdc[5];
ea1a8005 718 for (Int_t ih=0; ih<5; ih++)
719 {
720 tvdc[ih] = fESDTZERO->GetTVDC(ih);
721
528890e5 722 if( tvdc[0] !=0 && tvdc[0]> -10 && tvdc[0]< 10 )
ea1a8005 723 if(ih>0 && tvdc[ih]>20 ) pileup = true;
528890e5 724 if( tvdc[0] >20 || (tvdc[0] < -20 && tvdc[0] > -9000) ) pileup =true;
725 // if (pileup) printf(" !!!!! pile up %i tvdc %f \n",ih, tvdc[ih]);
ea1a8005 726 }
727
728
729 return pileup;
730
731}
732
733 //____________________________________________________________
734
735Bool_t AliT0Reconstructor::BackgroundFlag() const
736{
737 Bool_t background = false;
738
739 Float_t orA = fESDTZERO->GetOrA(0);
740 Float_t orC = fESDTZERO->GetOrC(0);
741 Float_t tvdc = fESDTZERO->GetTVDC(0);
dc7ca31d 742
ea1a8005 743 if ( (orA > -5 && orA <5) && (orC > -5 && orC <5) && (tvdc < -5 || tvdc > 5)) {
744 background = true;
745 // printf(" orA %f orC %f tvdc %f\n", orA, orC, tvdc);
746 }
747 return background;
748
749
750}
751
752
753 //____________________________________________________________
754
755Bool_t AliT0Reconstructor::SatelliteFlag() const
756{
757
758 Bool_t satellite = false;
759 for (Int_t i0=0; i0<24; i0++) {
760 Float_t timefull = fESDTZERO -> GetTimeFull(i0,0);
830172c8 761 if( timefull < -1.5 && timefull > -10) satellite=true;
ea1a8005 762 }
763
764 return satellite;
765
766}