]> git.uio.no Git - u/mrichter/AliRoot.git/blame - T0/AliT0Reconstructor.cxx
coverity fixed
[u/mrichter/AliRoot.git] / T0 / AliT0Reconstructor.cxx
CommitLineData
dc7ca31d 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
72e48d95 17/*********************************************************************
18 * T0 reconstruction and filling ESD
19 * - reconstruct mean time (interation time)
20 * - vertex position
21 * - multiplicity
22 ********************************************************************/
dc7ca31d 23
af885e0f 24#include <AliESDEvent.h>
dc7ca31d 25#include "AliLog.h"
dc7ca31d 26#include "AliT0RecPoint.h"
27#include "AliRawReader.h"
28#include "AliT0RawReader.h"
dc7ca31d 29#include "AliT0digit.h"
30#include "AliT0Reconstructor.h"
31#include "AliT0Parameters.h"
c41ceaac 32#include "AliT0Calibrator.h"
58bd3a16 33#include "AliESDfriend.h"
b0e13b29 34#include "AliESDTZERO.h"
73df58ab 35#include "AliESDTZEROfriend.h"
8f8d0732 36#include "AliLog.h"
85f61e3b 37#include "AliCDBEntry.h"
38#include "AliCDBManager.h"
39#include "AliCTPTimeParams.h"
40#include "AliLHCClockPhase.h"
669dc07f 41#include "AliT0CalibSeasonTimeShift.h"
4cbe597e 42#include "AliESDRun.h"
dc7ca31d 43
44#include <TArrayI.h>
45#include <TGraph.h>
aad72f45 46#include <TMath.h>
b09247a2 47#include <Riostream.h>
dc7ca31d 48
49ClassImp(AliT0Reconstructor)
50
c41ceaac 51 AliT0Reconstructor:: AliT0Reconstructor(): AliReconstructor(),
f16935f7 52 fdZonA(0),
53 fdZonC(0),
54 fZposition(0),
55 fParam(NULL),
2e6a5ee0 56 fAmpLEDrec(),
c883fdf2 57 fQTC(0),
58 fAmpLED(0),
58bd3a16 59 fCalib(),
60 fLatencyHPTDC(9000),
61 fLatencyL1(0),
62 fLatencyL1A(0),
73df58ab 63 fLatencyL1C(0),
85f61e3b 64 fGRPdelays(0),
669dc07f 65 fTimeMeanShift(0x0),
66 fTimeSigmaShift(0x0),
b0e13b29 67 fESDTZEROfriend(NULL),
68 fESDTZERO(NULL)
58bd3a16 69
e0bba6cc 70{
612737bb 71 for (Int_t i=0; i<24; i++) fTime0vertex[i] =0;
4cbe597e 72
72e48d95 73 //constructor
85f61e3b 74 AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
75 if (!entry) AliFatal("CTP timing parameters are not found in OCDB !");
76 AliCTPTimeParams *ctpParams = (AliCTPTimeParams*)entry->GetObject();
77 Float_t l1Delay = (Float_t)ctpParams->GetDelayL1L0()*25.0;
78
79 AliCDBEntry *entry1 = AliCDBManager::Instance()->Get("GRP/CTP/TimeAlign");
80 if (!entry1) AliFatal("CTP time-alignment is not found in OCDB !");
81 AliCTPTimeParams *ctpTimeAlign = (AliCTPTimeParams*)entry1->GetObject();
82 l1Delay += ((Float_t)ctpTimeAlign->GetDelayL1L0()*25.0);
83
84 AliCDBEntry *entry4 = AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
85 if (!entry4) AliFatal("LHC clock-phase shift is not found in OCDB !");
86 AliLHCClockPhase *phase = (AliLHCClockPhase*)entry4->GetObject();
87
85f61e3b 88 fGRPdelays = l1Delay - phase->GetMeanPhase();
72e48d95 89
669dc07f 90 AliCDBEntry *entry5 = AliCDBManager::Instance()->Get("T0/Calib/TimeAdjust");
91 if (entry5) {
92 AliT0CalibSeasonTimeShift *timeshift = (AliT0CalibSeasonTimeShift*)entry5->GetObject();
93 fTimeMeanShift = timeshift->GetT0Means();
94 fTimeSigmaShift = timeshift->GetT0Sigmas();
4cbe597e 95 }
669dc07f 96 else
97 AliWarning("Time Adjust is not found in OCDB !");
612737bb 98
74adb36a 99 fParam = AliT0Parameters::Instance();
100 fParam->Init();
c883fdf2 101
74adb36a 102 for (Int_t i=0; i<24; i++){
2e6a5ee0 103 TGraph* gr = fParam ->GetAmpLEDRec(i);
29ed1d0e 104 if (gr) fAmpLEDrec.AddAtAndExpand(gr,i) ;
c883fdf2 105 TGraph* gr1 = fParam ->GetAmpLED(i);
106 if (gr1) fAmpLED.AddAtAndExpand(gr1,i) ;
107 TGraph* gr2 = fParam ->GetQTC(i);
539b9cb9 108 if (gr2) fQTC.AddAtAndExpand(gr2,i) ;
612737bb 109 fTime0vertex[i] = fParam->GetCFD(i);
b0e13b29 110 AliDebug(2,Form("OCDB mean CFD time %i %f \n",i, fTime0vertex[i]));
612737bb 111 }
58bd3a16 112 fLatencyL1 = fParam->GetLatencyL1();
612737bb 113 fLatencyL1A = fParam->GetLatencyL1A();
58bd3a16 114 fLatencyL1C = fParam->GetLatencyL1C();
115 fLatencyHPTDC = fParam->GetLatencyHPTDC();
9d026202 116 AliDebug(2,Form(" LatencyL1 %f latencyL1A %f latencyL1C %f latencyHPTDC %f \n",fLatencyL1, fLatencyL1A, fLatencyL1C, fLatencyHPTDC));
612737bb 117
118 for (Int_t i=0; i<24; i++) {
59fe6376 119 if( fTime0vertex[i] < 500 || fTime0vertex[i] > 50000) fTime0vertex[i] =( 1000.*fLatencyHPTDC - 1000.*fLatencyL1 + 1000.*fGRPdelays)/24.4;
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
dc7ca31d 131}
c41ceaac 132
133//_____________________________________________________________________________
dc7ca31d 134void AliT0Reconstructor::Reconstruct(TTree*digitsTree, TTree*clustersTree) const
135{
94c27e4f 136 // T0 digits reconstruction
ce50812a 137 Int_t refAmp = 0 ; /*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
ce50812a 145 Float_t c = 29.9792458; // cm/ns
8955c6b4 146 Float_t channelWidth = fParam->GetChannelWidth() ;
ce50812a 147 Double32_t vertex = 9999999, meanVertex = 9999999 ;
adf36b9d 148 Double32_t timeDiff=999999, meanTime=999999, timeclock=999999;
ce50812a 149
94c27e4f 150
dc7ca31d 151 AliDebug(1,Form("Start DIGITS reconstruction "));
94c27e4f 152
ce50812a 153 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
154 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
155 printf( "AliT0Reconstructor::Reconstruct::: RecoParam amplitude %f %f \n",lowAmpThreshold, highAmpThreshold);
156
157 Double32_t besttimeA=9999999; Double32_t besttimeA_best=9999999;
158 Double32_t besttimeC=9999999; Double32_t besttimeC_best=9999999;
159 Int_t timeDelayCFD[24];
160 Int_t badpmt[24];
161 //Bad channel
162 for (Int_t i=0; i<24; i++) {
163 badpmt[i] = GetRecoParam() -> GetBadChannels(i);
164 timeDelayCFD[i] = Int_t (fParam->GetTimeDelayCFD(i));
165 }
166 fCalib->SetEq(0);
d0bcd1fb 167 TBranch *brDigits=digitsTree->GetBranch("T0");
e0bba6cc 168 AliT0digit *fDigits = new AliT0digit() ;
dc7ca31d 169 if (brDigits) {
170 brDigits->SetAddress(&fDigits);
171 }else{
f16935f7 172 AliError(Form("EXEC Branch T0 digits not found"));
ce50812a 173 return;
dc7ca31d 174 }
e0bba6cc 175
c41ceaac 176 digitsTree->GetEvent(0);
177 digitsTree->GetEntry(0);
178 brDigits->GetEntry(0);
179 fDigits->GetTimeCFD(*timeCFD);
180 fDigits->GetTimeLED(*timeLED);
181 fDigits->GetQT0(*chargeQT0);
182 fDigits->GetQT1(*chargeQT1);
446d6ec4 183 Int_t onlineMean = fDigits->MeanTime();
ce50812a 184
adf36b9d 185 Bool_t tr[5];
186 for (Int_t i=0; i<5; i++) tr[i]=false;
c41ceaac 187
dc7ca31d 188
7e08771f 189 AliT0RecPoint frecpoints;
190 AliT0RecPoint * pfrecpoints = &frecpoints;
191 clustersTree->Branch( "T0", "AliT0RecPoint" ,&pfrecpoints);
94c27e4f 192
195e1353 193 Float_t time[24], adc[24], adcmip[24];
dc7ca31d 194 for (Int_t ipmt=0; ipmt<24; ipmt++) {
ce50812a 195 if(timeCFD->At(ipmt)>0 ) {
196 if(( chargeQT1->At(ipmt) - chargeQT0->At(ipmt))>0)
d0bcd1fb 197 adc[ipmt] = chargeQT1->At(ipmt) - chargeQT0->At(ipmt);
a7027400 198 else
d0bcd1fb 199 adc[ipmt] = 0;
200
ce50812a 201 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD->At(ipmt)) ;
202
d0bcd1fb 203 Double_t sl = Double_t(timeLED->At(ipmt) - timeCFD->At(ipmt));
8f620945 204 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD->At(ipmt) ) ;
669dc07f 205 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
ce50812a 206 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
207
644a358e 208 Double_t ampMip = 0;
209 TGraph* ampGraph = (TGraph*)fAmpLED.At(ipmt);
210 if (ampGraph) ampMip = ampGraph->Eval(sl);
211 Double_t qtMip = 0;
212 TGraph* qtGraph = (TGraph*)fQTC.At(ipmt);
213 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
669dc07f 214 AliDebug(5,Form(" Amlitude in MIPS LED %f , QTC %f in channels %f\n ",ampMip,qtMip, adc[ipmt]));
7e08771f 215 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
216 frecpoints.SetAmpLED(ipmt, Float_t( ampMip));
217 frecpoints.SetAmp(ipmt, Float_t(qtMip));
195e1353 218 adcmip[ipmt]=qtMip;
d0bcd1fb 219
dc7ca31d 220 }
221 else {
222 time[ipmt] = 0;
223 adc[ipmt] = 0;
b0ab3f59 224 adcmip[ipmt] = 0;
ce50812a 225
dc7ca31d 226 }
227 }
94c27e4f 228
dc7ca31d 229 for (Int_t ipmt=0; ipmt<12; ipmt++){
ce50812a 230 if(time[ipmt] !=0 && time[ipmt] > -9000
231 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
232 {
233 if(time[ipmt]<besttimeC) besttimeC=time[ipmt]; //timeC
234 if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC_best))
235 besttimeC_best=time[ipmt]; //timeC
dc7ca31d 236 }
dc7ca31d 237 }
ce50812a 238 for ( Int_t ipmt=12; ipmt<24; ipmt++)
239 {
240 if(time[ipmt] != 0 && time[ipmt] > -9000
241 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
242 {
243 if(time[ipmt]<besttimeA) besttimeA=time[ipmt];
244 if(TMath::Abs(time[ipmt] ) < TMath::Abs(besttimeA_best))
245 besttimeA_best=time[ipmt]; //timeA
246 }
dc7ca31d 247 }
ce50812a 248
249 if( besttimeA < 999999 && besttimeA!=0) {
250 frecpoints.SetTimeBestA((besttimeA_best * channelWidth - fdZonA/c) - 12500.);
251 frecpoints.SetTime1stA((besttimeA * channelWidth - fdZonA/c) - 12500.);
adf36b9d 252 tr[1]=true;
253 }
ce50812a 254
255 if( besttimeC < 999999 && besttimeC!=0) {
256 frecpoints.SetTimeBestC((besttimeC_best * channelWidth - fdZonC/c) - 12500.);
257 frecpoints.SetTime1stC((besttimeC * channelWidth - fdZonC/c) - 12500.);
adf36b9d 258 tr[2]=true;
259 }
ce50812a 260
261 AliDebug(5,Form("1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f ",
262 besttimeA, besttimeA_best,
263 besttimeC, besttimeC_best) );
264
adf36b9d 265 if(besttimeA <999999 && besttimeC < 999999 ){
9b83615d 266 // timeDiff = (besttimeC - besttimeA)*channelWidth;
267 timeDiff = (besttimeA - besttimeC)*channelWidth;
ce50812a 268 meanTime = channelWidth * (besttimeA_best + besttimeC_best)/2. ;
269 timeclock = channelWidth * (besttimeA + besttimeC)/2. ;
adf36b9d 270 vertex = meanVertex - c*(timeDiff)/2.;// + (fdZonA - fdZonC)/2;
271 tr[0]=true;
272 }
7e08771f 273 frecpoints.SetVertex(vertex);
ce50812a 274 frecpoints.SetMeanTime(meanTime - 12500.);
275 frecpoints.SetT0clock(timeclock - 12500.);
7e08771f 276 frecpoints.SetT0Trig(tr);
ce50812a 277
278 AliDebug(5,Form("fRecPoints::: 1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f shiftA %f shiftC %f ",
279 frecpoints.Get1stTimeA(), frecpoints.GetBestTimeA(),
280 frecpoints.Get1stTimeC(), frecpoints.GetBestTimeC(),
281 fTimeMeanShift[1],fTimeMeanShift[2] ) );
282
669dc07f 283 AliDebug(5,Form("T0 triggers %d %d %d %d %d",tr[0],tr[1],tr[2],tr[3],tr[4]));
ce50812a 284
adf36b9d 285 //online mean
7e08771f 286 frecpoints.SetOnlineMean(Int_t(onlineMean));
b5a9f753 287 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 288
830172c8 289 clustersTree->Fill();
ce50812a 290
bd375212 291 delete timeCFD;
292 delete timeLED;
293 delete chargeQT0;
294 delete chargeQT1;
dc7ca31d 295}
296
297
c41ceaac 298//_______________________________________________________________________
299
300void AliT0Reconstructor::Reconstruct(AliRawReader* rawReader, TTree*recTree) const
301{
94c27e4f 302 // T0 raw ->
539b9cb9 303 //
c8d939c8 304
b0e13b29 305 Float_t meanOrA = fTime0vertex[0] + 587;
306 Float_t meanOrC = fTime0vertex[0] + 678;
307 Float_t meanTVDC = fTime0vertex[0] + 2564;
82540e64 308 Int_t timeDelayCFD[24];
309
d3e04608 310
612737bb 311 Int_t badpmt[24];
38cbfa7c 312 //Bad channel
58432641 313 for (Int_t i=0; i<24; i++) {
314 badpmt[i] = GetRecoParam() -> GetBadChannels(i);
82540e64 315 timeDelayCFD[i] = Int_t (fParam->GetTimeDelayCFD(i));
58432641 316 }
ce50812a 317 Int_t equalize = GetRecoParam() -> GetEq();
c8d939c8 318 // printf( "AliT0Reconstructor::Reconstruct::: RecoParam %i \n",equalize);
ce50812a 319 fCalib->SetEq(equalize);
d3e04608 320 Int_t low[500], high[500];
94b4f7e2 321 Float_t timefull=-9999;;
322 Float_t tvdc = -9999; Float_t ora = -9999; Float_t orc = -9999;
58bd3a16 323
e8ed1cd0 324 Int_t allData[110][5];
2e6a5ee0 325
e8ed1cd0 326 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
ce50812a 327 Float_t time2zero[24];
612737bb 328 Double32_t timeDiff, meanTime, timeclock;
329 timeDiff = meanTime = timeclock = 9999999;
adf36b9d 330 Float_t c = 29.9792458; // cm/ns
331 Double32_t vertex = 9999999;
776de217 332 Int_t onlineMean=0;
9480f05f 333 Float_t meanVertex = 0;
612737bb 334 for (Int_t i0=0; i0<24; i0++) {
335 low[i0] = Int_t(fTime0vertex[i0]) - 200;
336 high[i0] = Int_t(fTime0vertex[i0]) + 200;
ce50812a 337 time2zero[i0] = 99999;
612737bb 338 }
339
669dc07f 340 for (Int_t i0=0; i0<110; i0++)
b0e13b29 341 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
612737bb 342
343 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
344 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
ce50812a 345 printf( "AliT0Reconstructor::Reconstruct::: RecoParam amplitude %f %f \n",lowAmpThreshold, highAmpThreshold);
346
347 Double32_t besttimeA=9999999; Double32_t besttimeA_best=9999999;
348 Double32_t besttimeC=9999999; Double32_t besttimeC_best=9999999;
349
350 Float_t channelWidth = fParam->GetChannelWidth() ;
b0e13b29 351
7e08771f 352 AliT0RecPoint frecpoints;
353 AliT0RecPoint * pfrecpoints = &frecpoints;
bce12dc5 354
7e08771f 355 recTree->Branch( "T0", "AliT0RecPoint" ,&pfrecpoints);
2e6a5ee0 356
bce12dc5 357 AliDebug(10," before read data ");
358 AliT0RawReader myrawreader(rawReader);
776de217 359
360 UInt_t type =rawReader->GetType();
361
bce12dc5 362 if (!myrawreader.Next())
363 AliDebug(1,Form(" no raw data found!!"));
364 else
365 {
38cbfa7c 366 for (Int_t i=0; i<24; i++)
367 {
368 timeCFD[i]=0; timeLED[i]=0; chargeQT0[i]=0; chargeQT1[i]=0;
369 }
370 Int_t fBCID=Int_t (rawReader->GetBCID());
8f620945 371 Int_t trmbunch= myrawreader.GetTRMBunchID();
86fd0587 372 AliDebug(10,Form(" CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
373
38cbfa7c 374 if(type == 7 ) { //only physics
669dc07f 375 for (Int_t i=0; i<107; i++) {
bce12dc5 376 for (Int_t iHit=0; iHit<5; iHit++)
377 {
378 allData[i][iHit] = myrawreader.GetData(i,iHit);
379 }
8f620945 380 }
8f620945 381
85f61e3b 382 for (Int_t in=0; in<12; in++)
383 {
384 for (Int_t iHit=0; iHit<5; iHit++)
38cbfa7c 385 {
612737bb 386 if(allData[in+1][iHit] > low[in] &&
387 allData[in+1][iHit] < high[in])
38cbfa7c 388 {
85f61e3b 389 timeCFD[in] = allData[in+1][iHit] ;
390 break;
38cbfa7c 391 }
85f61e3b 392 }
393 for (Int_t iHit=0; iHit<5; iHit++)
394 {
612737bb 395 if(allData[in+1+56][iHit] > low[in] &&
396 allData[in+1+56][iHit] < high[in])
38cbfa7c 397 {
85f61e3b 398 timeCFD[in+12] = allData[in+56+1][iHit] ;
399 break;
38cbfa7c 400 }
38cbfa7c 401 }
612737bb 402 timeLED[in+12] = allData[in+68+1][0] ;
403 timeLED[in] = allData[in+12+1][0] ;
ce50812a 404 AliDebug(50, Form(" readed i %i cfdC %i cfdA %i ledC %i ledA%i ",
85f61e3b 405 in, timeCFD[in],timeCFD[in+12],timeLED[in],
612737bb 406 timeLED[in+12]));
38cbfa7c 407
8f620945 408 }
409
8f620945 410
85f61e3b 411 for (Int_t in=0; in<12; in++)
412 {
413 chargeQT0[in]=allData[2*in+25][0];
414 chargeQT1[in]=allData[2*in+26][0];
612737bb 415 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 416 in, chargeQT0[in],chargeQT1[in]));
417 }
418 for (Int_t in=12; in<24; in++)
419 {
420 chargeQT0[in]=allData[2*in+57][0];
421 chargeQT1[in]=allData[2*in+58][0];
612737bb 422 AliDebug(25, Form(" readed Raw %i %i %i",
85f61e3b 423 in, chargeQT0[in],chargeQT1[in]));
85f61e3b 424 }
425
8f226345 426 onlineMean = allData[49][0];
427
195e1353 428 Double32_t time[24], adc[24], adcmip[24], noncalibtime[24];
8f620945 429 for (Int_t ipmt=0; ipmt<24; ipmt++) {
36bfca7d 430 if(timeCFD[ipmt] > 0 /* && badpmt[ipmt]==0*/ ){
bce12dc5 431 //for simulated data
d0bcd1fb 432 //for physics data
73df58ab 433 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>0) {
541b42c4 434 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 435 }
d0bcd1fb 436 else
437 adc[ipmt] = 0;
612737bb 438 // time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
c8d939c8 439 Int_t refAmp = Int_t (fTime0vertex[ipmt]);
440 time[ipmt] = fCalib-> WalkCorrection( refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
d0bcd1fb 441 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 442 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 443 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 444 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
644a358e 445 Double_t ampMip = 0;
446 TGraph * ampGraph = (TGraph*)fAmpLED.At(ipmt);
447 if (ampGraph) ampMip = ampGraph->Eval(sl);
448 Double_t qtMip = 0;
449 TGraph * qtGraph = (TGraph*)fQTC.At(ipmt);
450 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
345f03db 451 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %f\n ",ampMip,qtMip, adc[ipmt]));
c8d939c8 452 if( qtMip>lowAmpThreshold && qtMip<highAmpThreshold )
453 {
454 if( equalize ==0 )
455 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
ce50812a 456 else
457 frecpoints.SetTime(ipmt, Float_t(time[ipmt] + fTime0vertex[ipmt]) );
c8d939c8 458 // frecpoints.SetTime(ipmt, Float_t(time[ipmt] ) );
459 }
ce50812a 460 frecpoints.SetAmp(ipmt, Double32_t( qtMip));
461 adcmip[ipmt]=qtMip;
462 frecpoints.SetAmpLED(ipmt, Double32_t(ampMip));
463 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
2e6a5ee0 464 }
465 else {
ce50812a 466 time[ipmt] = -9999;
2e6a5ee0 467 adc[ipmt] = 0;
b0ab3f59 468 adcmip[ipmt] = 0;
ce50812a 469 noncalibtime[ipmt] = -9999;
2e6a5ee0 470 }
471 }
ce50812a 472 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
473
2e6a5ee0 474 for (Int_t ipmt=0; ipmt<12; ipmt++){
ce50812a 475 if(time[ipmt] !=0 && time[ipmt] > -9000
476 /*&& badpmt[ipmt]==0 */
477 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
73df58ab 478 {
ce50812a 479 if(time[ipmt]<besttimeC) besttimeC=time[ipmt]; //timeC
480 if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC_best))
481 besttimeC_best=time[ipmt]; //timeC
2e6a5ee0 482 }
2e6a5ee0 483 }
73df58ab 484 for ( Int_t ipmt=12; ipmt<24; ipmt++)
485 {
ce50812a 486 if(time[ipmt] != 0 && time[ipmt] > -9000
487 /* && badpmt[ipmt]==0*/
488 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
73df58ab 489 {
ce50812a 490 if(time[ipmt]<besttimeA) besttimeA=time[ipmt];
491 if(TMath::Abs(time[ipmt] ) < TMath::Abs(besttimeA_best))
492 besttimeA_best=time[ipmt]; //timeA
73df58ab 493 }
2e6a5ee0 494 }
ce50812a 495
496 if(besttimeA < 999999 && besttimeA!=0 ) {
497 if( equalize ==0 )
498 frecpoints.SetTime1stA((besttimeA * channelWidth)- 1000.*fLatencyHPTDC + 1000.*fLatencyL1A - 1000.*fGRPdelays - fTimeMeanShift[1] );
499 else
500 {
501 frecpoints.SetTimeBestA((besttimeA_best * channelWidth ));
502 frecpoints.SetTime1stA((besttimeA * channelWidth - fTimeMeanShift[1]));
503 }
504 }
505 if( besttimeC < 999999 && besttimeC!=0) {
506 if( equalize ==0 )
507 frecpoints.SetTime1stC((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1C - 1000.*fGRPdelays - fTimeMeanShift[2]);
508 else
509 {
510 frecpoints.SetTimeBestC((besttimeC_best * channelWidth ));
511 frecpoints.SetTime1stC((besttimeC * channelWidth - fTimeMeanShift[2]));
512 }
513 }
514 AliDebug(5,Form("1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f ",
515 besttimeA, besttimeA_best,
516 besttimeC, besttimeC_best) );
517 AliDebug(5,Form("fRecPoints::: 1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f shiftA %f shiftC %f ",
518 frecpoints.Get1stTimeA(), frecpoints.GetBestTimeA(),
519 frecpoints.Get1stTimeC(), frecpoints.GetBestTimeC(),
520 fTimeMeanShift[1],fTimeMeanShift[2] ) );
521 if( besttimeC < 999999 && besttimeA < 999999) {
c8d939c8 522 if(equalize ==0 )
523 timeclock = ((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1 - 1000.*fGRPdelays - fTimeMeanShift[0]);
524 else
525 {
526 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - fTimeMeanShift[0];
527 meanTime = channelWidth * Float_t(besttimeA_best + besttimeC_best )/2.;
528 }
612737bb 529 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth ;
9480f05f 530 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
ce50812a 531 }
532
776de217 533 } //if phys event
ce50812a 534 AliDebug(1,Form(" timeDiff %f #channel, meanTime %f #ps, TOFmean%f vertex %f cm meanVertex %f \n",timeDiff, meanTime,timeclock, vertex,meanVertex));
535 frecpoints.SetT0clock(timeclock);
536 frecpoints.SetVertex(vertex);
537 frecpoints.SetMeanTime(meanTime);
538 frecpoints.SetOnlineMean(Int_t(onlineMean));
adf36b9d 539
ce50812a 540 // Set triggers
adf36b9d 541 Bool_t tr[5];
92f6bfd3 542 Int_t trchan[5] = {50,51,52,55,56};
543 Float_t lowtr[5] = {meanTVDC-700, meanOrA-700, meanOrC-700, meanOrC-1000, meanOrC-1000 };
544 Float_t hightr[5] = {meanTVDC+700, meanOrA+700, meanOrC+700, meanOrC+1000, meanOrC+1000};
ce50812a 545
82540e64 546 for (Int_t i=0; i<5; i++) tr[i] = false;
adf36b9d 547 for (Int_t itr=0; itr<5; itr++) {
b0e13b29 548 for (Int_t iHit=0; iHit<1; iHit++)
38cbfa7c 549 {
550 Int_t trr=trchan[itr];
92f6bfd3 551 if( allData[trr][iHit] > lowtr[itr] && allData[trr][iHit] < hightr[itr]) tr[itr]=true;
ce50812a 552
553 AliDebug(15,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 554 }
adf36b9d 555 }
7e08771f 556 frecpoints.SetT0Trig(tr);
b0e13b29 557
ce50812a 558 // all times with amplitude correction
559 Float_t timecent;
b0e13b29 560 for (Int_t iHit=0; iHit<5; iHit++)
561 {
ce50812a 562 timefull = timecent = -9999;
82540e64 563 tvdc = ora = orc = -9999;
b0e13b29 564 if(allData[50][iHit]>0)
565 tvdc = (Float_t(allData[50][iHit]) - meanTVDC) * channelWidth* 0.001;
566 if(allData[51][iHit]>0)
567 ora = (Float_t(allData[51][iHit]) - meanOrA) * channelWidth* 0.001;
568
569 if(allData[52][iHit]>0)
570 orc = (Float_t(allData[52][iHit]) - meanOrC) * channelWidth* 0.001;
571
7e08771f 572 frecpoints.SetOrC( iHit, orc);
573 frecpoints.SetOrA( iHit, ora);
574 frecpoints.SetTVDC( iHit, tvdc);
b0e13b29 575 for (Int_t i0=0; i0<12; i0++) {
ce50812a 576 if (equalize ==0 )
577 timecent = fTime0vertex[i0] + timeDelayCFD[i0];
578 else
579 timecent = fTime0vertex[i0];
b0e13b29 580 timefull = -9999;
581 if(allData[i0+1][iHit]>1)
ce50812a 582 timefull = (Float_t(allData[i0+1][iHit]) - timecent)* channelWidth* 0.001;
7e08771f 583 frecpoints.SetTimeFull(i0, iHit,timefull) ;
82540e64 584 // 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 585
586 }
587
588 for (Int_t i0=12; i0<24; i0++) {
589 timefull = -9999;
ce50812a 590 if (equalize ==0 )
591 timecent = fTime0vertex[i0] + timeDelayCFD[i0];
592 else
593 timecent = fTime0vertex[i0];
94b4f7e2 594 if(allData[i0+45][iHit]>1) {
ce50812a 595 timefull = (Float_t(allData[i0+45][iHit]) - timecent)* channelWidth* 0.001;
94b4f7e2 596 }
82540e64 597 // 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 598 frecpoints.SetTimeFull(i0, iHit, timefull) ;
b0e13b29 599 }
600 }
601
602
669dc07f 603 //Set MPD
604 if(allData[53][0]>0 && allData[54][0])
7e08771f 605 frecpoints.SetMultA(allData[53][0]-allData[54][0]);
b0e13b29 606 if(allData[105][0]>0 && allData[106][0])
7e08771f 607 frecpoints.SetMultC(allData[105][0]-allData[106][0]);
b0e13b29 608
609
669dc07f 610 } // if (else )raw data
58bd3a16 611 recTree->Fill();
58bd3a16 612}
adf36b9d 613
614
ce50812a 615//____________________________________________________________
adf36b9d 616
617 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
618 {
dc7ca31d 619
620 /***************************************************
621 Resonstruct digits to vertex position
622 ****************************************************/
623
dc7ca31d 624 AliDebug(1,Form("Start FillESD T0"));
b0ab3f59 625 if(!pESD) {
626 AliError("No ESD Event");
627 return;
628 }
4cbe597e 629 pESD ->SetT0spread(fTimeSigmaShift);
630
36cde487 631
58bd3a16 632 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 633 Float_t c = 0.0299792458; // cm/ps
adf36b9d 634 Float_t currentVertex=0, shift=0;
291f31a1 635 Int_t ncont=-1;
adf36b9d 636 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
637 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
638 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
639 if (!vertex) vertex = pESD->GetVertex();
640
641 if (vertex) {
642 AliDebug(2, Form("Got %s (%s) from ESD: %f",
643 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
644 currentVertex = vertex->GetZ();
645
646 ncont = vertex->GetNContributors();
291f31a1 647 if(ncont>0 ) {
adf36b9d 648 shift = currentVertex/c;
adf36b9d 649 }
650 }
d76c31f4 651 TTree *treeR = clustersTree;
dc7ca31d 652
7e08771f 653 AliT0RecPoint frecpoints;
654 AliT0RecPoint * pfrecpoints = &frecpoints;
dc7ca31d 655
656 AliDebug(1,Form("Start FillESD T0"));
657 TBranch *brRec = treeR->GetBranch("T0");
658 if (brRec) {
7e08771f 659 brRec->SetAddress(&pfrecpoints);
dc7ca31d 660 }else{
f16935f7 661 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 662 return;
663 }
73df58ab 664
665 brRec->GetEntry(0);
666 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
667 Double32_t* tcorr;
668 for(Int_t i=0; i<24; i++)
669 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
670
ce50812a 671 //1st time
73df58ab 672 Double32_t timeClock[3];
7e08771f 673 Double32_t zPosition = frecpoints.GetVertex();
830172c8 674
ce50812a 675 timeClock[0] = frecpoints.GetT0clock() ;
676 timeClock[1] = frecpoints.Get1stTimeA() + shift;
677 timeClock[2] = frecpoints.Get1stTimeC() - shift;
678 //best time
679 Double32_t timemean[3];
680 timemean[0] = frecpoints.GetMeanTime();
681 timemean[1] = frecpoints.GetBestTimeA() + shift;
682 timemean[2] = frecpoints.GetBestTimeC() - shift;
683
684 for(Int_t i=0; i<3; i++) {
685 fESDTZERO->SetT0TOF(i,timeClock[i]); // interaction time (ns)
686 fESDTZERO->SetT0TOFbest(i,timemean[i]); // interaction time (ns)
687 }
73df58ab 688 for ( Int_t i=0; i<24; i++) {
7e08771f 689 time[i] = frecpoints.GetTime(i); // ps to ns
ce50812a 690 if ( time[i] != 0 && time[i]>-9999) {
7e08771f 691 ampQTC[i] = frecpoints.GetAmp(i);
692 amp[i] = frecpoints.AmpLED(i);
612737bb 693 AliDebug(1,Form("T0: %i time %f ampQTC %f ampLED %f \n", i, time[i], ampQTC[i], amp[i]));
85f61e3b 694 }
73df58ab 695 }
ce50812a 696 fESDTZERO->SetT0time(time); // best TOF on each PMT
697 fESDTZERO->SetT0amplitude(ampQTC); // number of particles(MIPs) on each PMT
7e08771f 698 Int_t trig= frecpoints.GetT0Trig();
699 frecpoints.PrintTriggerSignals( trig);
92f6bfd3 700 // printf(" !!!!! FillESD trigger %i \n",trig);
b0e13b29 701 fESDTZERO->SetT0Trig(trig);
b0e13b29 702 fESDTZERO->SetT0zVertex(zPosition); //vertex Z position
669dc07f 703
7e08771f 704 Double32_t multA=frecpoints.GetMultA();
705 Double32_t multC=frecpoints.GetMultC();
b0e13b29 706 fESDTZERO->SetMultA(multA); // for backward compatubility
707 fESDTZERO->SetMultC(multC); // for backward compatubility
708
709
710 for (Int_t iHit =0; iHit<5; iHit++ ) {
82540e64 711 AliDebug(10,Form("FillESD ::: iHit %i tvdc %f orA %f orC %f\n", iHit,
7e08771f 712 frecpoints.GetTVDC(iHit),
713 frecpoints.GetOrA(iHit),
714 frecpoints.GetOrC(iHit) ));
715 fESDTZERO->SetTVDC(iHit,frecpoints.GetTVDC(iHit));
716 fESDTZERO->SetOrA(iHit,frecpoints.GetOrA(iHit));
717 fESDTZERO->SetOrC(iHit,frecpoints.GetOrC(iHit));
b0e13b29 718
719 for (Int_t i0=0; i0<24; i0++) {
7e08771f 720 // if(frecpoints.GetTimeFull(i0,iHit)>0){
721 // printf("FillESD ::: iHit %i cfd %i time %f \n", iHit, i0, frecpoints.GetTimeFull(i0,iHit));
722 fESDTZERO->SetTimeFull(i0, iHit,frecpoints.GetTimeFull(i0,iHit));
94b4f7e2 723 // }
b0e13b29 724
725 }
726 }
73df58ab 727
ce50812a 728 AliDebug(1,Form("T0: SPDshift %f Vertex %f (T0A+T0C)/2 best %f #ps T0signal %f ps OrA %f ps OrC %f ps T0trig %i\n",shift, zPosition, timemean[0], timeClock[0], timeClock[1], timeClock[2], trig));
b0e13b29 729
ea1a8005 730 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
731 // background flags
732 Bool_t background = BackgroundFlag();
733 fESDTZERO->SetBackgroundFlag(background);
734 Bool_t pileup = PileupFlag();
735 fESDTZERO->SetPileupFlag(pileup);
528890e5 736 for (Int_t i=0; i<5; i++) {
7e08771f 737 fESDTZERO->SetPileupTime(i, frecpoints.GetTVDC(i) ) ;
738 // printf("!!!!!! FillESD :: pileup %i %f %f \n", i,fESDTZERO->GetPileupTime(i), frecpoints.GetTVDC(i));
528890e5 739 }
ea1a8005 740 Bool_t sat = SatelliteFlag();
741 fESDTZERO->SetSatelliteFlag(sat);
73df58ab 742
ea1a8005 743 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
73df58ab 744 if (pESD) {
adf36b9d 745
73df58ab 746 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
747 if (fr) {
82540e64 748 AliDebug(10, Form("Writing TZERO friend data to ESD tree"));
73df58ab 749
85f61e3b 750 // if (ncont>2) {
73df58ab 751 tcorr = fESDTZEROfriend->GetT0timeCorr();
752 for ( Int_t i=0; i<24; i++) {
85f61e3b 753 if(i<12 && time[i]>1) timecorr[i] = tcorr[i] - shift/channelWidth;
754 if(i>11 && time[i]>1) timecorr[i] = tcorr[i] + shift/channelWidth;
612737bb 755 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 756 }
73df58ab 757 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
758 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
759 fESDTZEROfriend->SetT0ampQTC(ampQTC);
85f61e3b 760 fr->SetTZEROfriend(fESDTZEROfriend);
761 // }//
58bd3a16 762 }
b0e13b29 763
764 pESD->SetTZEROData(fESDTZERO);
58bd3a16 765 }
73df58ab 766
dc7ca31d 767} // vertex in 3 sigma
768
ea1a8005 769//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
770 //____________________________________________________________
771
772Bool_t AliT0Reconstructor::PileupFlag() const
773{
774 //
775 Bool_t pileup = false;
528890e5 776 Float_t tvdc[5];
ea1a8005 777 for (Int_t ih=0; ih<5; ih++)
778 {
779 tvdc[ih] = fESDTZERO->GetTVDC(ih);
780
528890e5 781 if( tvdc[0] !=0 && tvdc[0]> -10 && tvdc[0]< 10 )
ea1a8005 782 if(ih>0 && tvdc[ih]>20 ) pileup = true;
528890e5 783 if( tvdc[0] >20 || (tvdc[0] < -20 && tvdc[0] > -9000) ) pileup =true;
784 // if (pileup) printf(" !!!!! pile up %i tvdc %f \n",ih, tvdc[ih]);
ea1a8005 785 }
786
787
788 return pileup;
789
790}
791
792 //____________________________________________________________
793
794Bool_t AliT0Reconstructor::BackgroundFlag() const
795{
796 Bool_t background = false;
797
798 Float_t orA = fESDTZERO->GetOrA(0);
799 Float_t orC = fESDTZERO->GetOrC(0);
800 Float_t tvdc = fESDTZERO->GetTVDC(0);
dc7ca31d 801
ea1a8005 802 if ( (orA > -5 && orA <5) && (orC > -5 && orC <5) && (tvdc < -5 || tvdc > 5)) {
803 background = true;
804 // printf(" orA %f orC %f tvdc %f\n", orA, orC, tvdc);
805 }
806 return background;
807
808
809}
810
811
812 //____________________________________________________________
813
814Bool_t AliT0Reconstructor::SatelliteFlag() const
815{
816
ce50812a 817 Float_t satelliteLow = GetRecoParam() -> GetLowSatelliteThreshold();
818 Float_t satelliteHigh = GetRecoParam() -> GetHighSatelliteThreshold();
ea1a8005 819 Bool_t satellite = false;
820 for (Int_t i0=0; i0<24; i0++) {
821 Float_t timefull = fESDTZERO -> GetTimeFull(i0,0);
ce50812a 822 if( timefull > satelliteLow && timefull < satelliteHigh) satellite=true;
ea1a8005 823 }
824
825 return satellite;
826
827}