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