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