]> 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
trigger class added (M.Ivanov)
[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"));
5773ee77 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() ;
eadc5bd4 147 Double32_t vertex = 9999999, meanVertex = 0 ;
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 ) {
5773ee77 196 Float_t timefull = 0.001*( timeCFD->At(ipmt) - 511 - timeDelayCFD[ipmt]) * channelWidth;
eadc5bd4 197 frecpoints.SetTimeFull(ipmt, 0 ,timefull) ;
ce50812a 198 if(( chargeQT1->At(ipmt) - chargeQT0->At(ipmt))>0)
d0bcd1fb 199 adc[ipmt] = chargeQT1->At(ipmt) - chargeQT0->At(ipmt);
a7027400 200 else
d0bcd1fb 201 adc[ipmt] = 0;
202
ce50812a 203 time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD->At(ipmt)) ;
5773ee77 204 time[ipmt] = time[ipmt] - 511;
d0bcd1fb 205 Double_t sl = Double_t(timeLED->At(ipmt) - timeCFD->At(ipmt));
8f620945 206 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD->At(ipmt) ) ;
eadc5bd4 207 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
ce50812a 208 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
209
644a358e 210 Double_t ampMip = 0;
211 TGraph* ampGraph = (TGraph*)fAmpLED.At(ipmt);
212 if (ampGraph) ampMip = ampGraph->Eval(sl);
213 Double_t qtMip = 0;
214 TGraph* qtGraph = (TGraph*)fQTC.At(ipmt);
215 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
669dc07f 216 AliDebug(5,Form(" Amlitude in MIPS LED %f , QTC %f in channels %f\n ",ampMip,qtMip, adc[ipmt]));
7e08771f 217 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
218 frecpoints.SetAmpLED(ipmt, Float_t( ampMip));
219 frecpoints.SetAmp(ipmt, Float_t(qtMip));
195e1353 220 adcmip[ipmt]=qtMip;
d0bcd1fb 221
dc7ca31d 222 }
223 else {
eadc5bd4 224 time[ipmt] = -99999;
dc7ca31d 225 adc[ipmt] = 0;
b0ab3f59 226 adcmip[ipmt] = 0;
ce50812a 227
dc7ca31d 228 }
229 }
94c27e4f 230
dc7ca31d 231 for (Int_t ipmt=0; ipmt<12; ipmt++){
ce50812a 232 if(time[ipmt] !=0 && time[ipmt] > -9000
233 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
234 {
235 if(time[ipmt]<besttimeC) besttimeC=time[ipmt]; //timeC
236 if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC_best))
237 besttimeC_best=time[ipmt]; //timeC
dc7ca31d 238 }
dc7ca31d 239 }
ce50812a 240 for ( Int_t ipmt=12; ipmt<24; ipmt++)
241 {
242 if(time[ipmt] != 0 && time[ipmt] > -9000
243 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
244 {
245 if(time[ipmt]<besttimeA) besttimeA=time[ipmt];
246 if(TMath::Abs(time[ipmt] ) < TMath::Abs(besttimeA_best))
247 besttimeA_best=time[ipmt]; //timeA
248 }
dc7ca31d 249 }
ce50812a 250
251 if( besttimeA < 999999 && besttimeA!=0) {
eadc5bd4 252 frecpoints.SetTimeBestA((besttimeA_best * channelWidth - fdZonA/c) );
5773ee77 253 frecpoints.SetTime1stA((besttimeA * channelWidth - fdZonA/c) );
adf36b9d 254 tr[1]=true;
255 }
ce50812a 256
257 if( besttimeC < 999999 && besttimeC!=0) {
eadc5bd4 258 frecpoints.SetTimeBestC((besttimeC_best * channelWidth - fdZonC/c) );
5773ee77 259 frecpoints.SetTime1stC((besttimeC * channelWidth - fdZonC/c) );
adf36b9d 260 tr[2]=true;
261 }
ce50812a 262
263 AliDebug(5,Form("1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f ",
264 besttimeA, besttimeA_best,
265 besttimeC, besttimeC_best) );
266
adf36b9d 267 if(besttimeA <999999 && besttimeC < 999999 ){
9b83615d 268 // timeDiff = (besttimeC - besttimeA)*channelWidth;
269 timeDiff = (besttimeA - besttimeC)*channelWidth;
ce50812a 270 meanTime = channelWidth * (besttimeA_best + besttimeC_best)/2. ;
271 timeclock = channelWidth * (besttimeA + besttimeC)/2. ;
eadc5bd4 272 vertex = meanVertex - 0.001* c*(timeDiff)/2.;// + (fdZonA - fdZonC)/2;
adf36b9d 273 tr[0]=true;
274 }
7e08771f 275 frecpoints.SetVertex(vertex);
eadc5bd4 276 frecpoints.SetMeanTime(meanTime );
277 frecpoints.SetT0clock(timeclock );
7e08771f 278 frecpoints.SetT0Trig(tr);
ce50812a 279
eadc5bd4 280 AliDebug(5,Form("fRecPoints::: 1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f vertex %f",
ce50812a 281 frecpoints.Get1stTimeA(), frecpoints.GetBestTimeA(),
282 frecpoints.Get1stTimeC(), frecpoints.GetBestTimeC(),
eadc5bd4 283 vertex ) );
ce50812a 284
669dc07f 285 AliDebug(5,Form("T0 triggers %d %d %d %d %d",tr[0],tr[1],tr[2],tr[3],tr[4]));
ce50812a 286
adf36b9d 287 //online mean
7e08771f 288 frecpoints.SetOnlineMean(Int_t(onlineMean));
b5a9f753 289 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 290
830172c8 291 clustersTree->Fill();
ce50812a 292
bd375212 293 delete timeCFD;
294 delete timeLED;
295 delete chargeQT0;
296 delete chargeQT1;
dc7ca31d 297}
298
299
c41ceaac 300//_______________________________________________________________________
301
302void AliT0Reconstructor::Reconstruct(AliRawReader* rawReader, TTree*recTree) const
303{
94c27e4f 304 // T0 raw ->
539b9cb9 305 //
c8d939c8 306
b0e13b29 307 Float_t meanOrA = fTime0vertex[0] + 587;
308 Float_t meanOrC = fTime0vertex[0] + 678;
309 Float_t meanTVDC = fTime0vertex[0] + 2564;
82540e64 310 Int_t timeDelayCFD[24];
5773ee77 311 Int_t corridor = GetRecoParam() -> GetCorridor();
794ab872 312 printf("!!!! corrior %i \n",corridor);
612737bb 313 Int_t badpmt[24];
38cbfa7c 314 //Bad channel
58432641 315 for (Int_t i=0; i<24; i++) {
316 badpmt[i] = GetRecoParam() -> GetBadChannels(i);
82540e64 317 timeDelayCFD[i] = Int_t (fParam->GetTimeDelayCFD(i));
58432641 318 }
ce50812a 319 Int_t equalize = GetRecoParam() -> GetEq();
c8d939c8 320 // printf( "AliT0Reconstructor::Reconstruct::: RecoParam %i \n",equalize);
ce50812a 321 fCalib->SetEq(equalize);
d3e04608 322 Int_t low[500], high[500];
94b4f7e2 323 Float_t timefull=-9999;;
324 Float_t tvdc = -9999; Float_t ora = -9999; Float_t orc = -9999;
58bd3a16 325
e8ed1cd0 326 Int_t allData[110][5];
2e6a5ee0 327
e8ed1cd0 328 Int_t timeCFD[24], timeLED[24], chargeQT0[24], chargeQT1[24];
ce50812a 329 Float_t time2zero[24];
612737bb 330 Double32_t timeDiff, meanTime, timeclock;
331 timeDiff = meanTime = timeclock = 9999999;
adf36b9d 332 Float_t c = 29.9792458; // cm/ns
333 Double32_t vertex = 9999999;
776de217 334 Int_t onlineMean=0;
9480f05f 335 Float_t meanVertex = 0;
b7fd8e80 336 Int_t pedestal[24];
612737bb 337 for (Int_t i0=0; i0<24; i0++) {
5773ee77 338 low[i0] = Int_t(fTime0vertex[i0]) - corridor;
339 high[i0] = Int_t(fTime0vertex[i0]) + corridor;
340 time2zero[i0] = 99999;
341 pedestal[i0]=Int_t (GetRecoParam()->GetLow(100+i0) );
342 // printf("pmt %i pedestal %f\n", i0,pedestal[i0]);
343
612737bb 344 }
345
669dc07f 346 for (Int_t i0=0; i0<110; i0++)
b0e13b29 347 for (Int_t j0=0; j0<5; j0++) allData[i0][j0]=0;
612737bb 348
349 Float_t lowAmpThreshold = GetRecoParam()->GetAmpLowThreshold();
350 Float_t highAmpThreshold = GetRecoParam()->GetAmpHighThreshold();
5773ee77 351 printf( "AliT0Reconstructor::Reconstruct::: RecoParam amplitude %f %f \n",lowAmpThreshold, highAmpThreshold);
ce50812a 352
353 Double32_t besttimeA=9999999; Double32_t besttimeA_best=9999999;
354 Double32_t besttimeC=9999999; Double32_t besttimeC_best=9999999;
355
356 Float_t channelWidth = fParam->GetChannelWidth() ;
b0e13b29 357
7e08771f 358 AliT0RecPoint frecpoints;
359 AliT0RecPoint * pfrecpoints = &frecpoints;
bce12dc5 360
7e08771f 361 recTree->Branch( "T0", "AliT0RecPoint" ,&pfrecpoints);
2e6a5ee0 362
bce12dc5 363 AliDebug(10," before read data ");
364 AliT0RawReader myrawreader(rawReader);
776de217 365
366 UInt_t type =rawReader->GetType();
367
bce12dc5 368 if (!myrawreader.Next())
369 AliDebug(1,Form(" no raw data found!!"));
370 else
371 {
38cbfa7c 372 for (Int_t i=0; i<24; i++)
373 {
374 timeCFD[i]=0; timeLED[i]=0; chargeQT0[i]=0; chargeQT1[i]=0;
375 }
9120829d 376
38cbfa7c 377 if(type == 7 ) { //only physics
669dc07f 378 for (Int_t i=0; i<107; i++) {
bce12dc5 379 for (Int_t iHit=0; iHit<5; iHit++)
380 {
381 allData[i][iHit] = myrawreader.GetData(i,iHit);
382 }
8f620945 383 }
8f620945 384
9120829d 385 Int_t fBCID=Int_t (rawReader->GetBCID());
386 Int_t trmbunch= myrawreader.GetTRMBunchID();
387 AliDebug(10,Form(" CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
388 if( (trmbunch-fBCID)!=37 ) {
389 AliDebug(0,Form("wrong :::: CDH BC ID %i, TRM BC ID %i \n", fBCID, trmbunch ));
390 // type = -1;
391 }
85f61e3b 392 for (Int_t in=0; in<12; in++)
393 {
394 for (Int_t iHit=0; iHit<5; iHit++)
38cbfa7c 395 {
612737bb 396 if(allData[in+1][iHit] > low[in] &&
397 allData[in+1][iHit] < high[in])
38cbfa7c 398 {
85f61e3b 399 timeCFD[in] = allData[in+1][iHit] ;
794ab872 400 break;
38cbfa7c 401 }
85f61e3b 402 }
403 for (Int_t iHit=0; iHit<5; iHit++)
404 {
794ab872 405 if(allData[in+1+56][iHit] > low[in+12] &&
406 allData[in+1+56][iHit] < high[in+12])
38cbfa7c 407 {
85f61e3b 408 timeCFD[in+12] = allData[in+56+1][iHit] ;
409 break;
38cbfa7c 410 }
38cbfa7c 411 }
612737bb 412 timeLED[in+12] = allData[in+68+1][0] ;
413 timeLED[in] = allData[in+12+1][0] ;
ce50812a 414 AliDebug(50, Form(" readed i %i cfdC %i cfdA %i ledC %i ledA%i ",
85f61e3b 415 in, timeCFD[in],timeCFD[in+12],timeLED[in],
612737bb 416 timeLED[in+12]));
38cbfa7c 417
8f620945 418 }
419
8f620945 420
85f61e3b 421 for (Int_t in=0; in<12; in++)
422 {
794ab872 423 if(allData[2*in+26][0] > 4700 &&
424 allData[2*in+26][0] < 5700)
425 {
426
427 chargeQT0[in]=allData[2*in+25][0];
428 chargeQT1[in]=allData[2*in+26][0];
429 AliDebug(25, Form(" readed Raw %i %i %i",
430 in, chargeQT0[in],chargeQT1[in]));
431 }
85f61e3b 432 }
433 for (Int_t in=12; in<24; in++)
434 {
794ab872 435 if(allData[2*in+58][0] > 4700 &&
436 allData[2*in+58][0] < 5700)
437 {
438 chargeQT0[in]=allData[2*in+57][0];
439 chargeQT1[in]=allData[2*in+58][0];
440 AliDebug(25, Form(" readed Raw %i %i %i",
441 in, chargeQT0[in],chargeQT1[in]));
442 }
85f61e3b 443 }
444
8f226345 445 onlineMean = allData[49][0];
446
195e1353 447 Double32_t time[24], adc[24], adcmip[24], noncalibtime[24];
8f620945 448 for (Int_t ipmt=0; ipmt<24; ipmt++) {
5773ee77 449 // if(timeCFD[ipmt] > 0 && (chargeQT0[ipmt] - chargeQT1[ipmt])>pedestal[ipmt] ){
450 if(timeCFD[ipmt] > 0 && (chargeQT0[ipmt] - chargeQT1[ipmt])> 0 ){
451 //for simulated data
d0bcd1fb 452 //for physics data
5773ee77 453 if(( chargeQT0[ipmt] - chargeQT1[ipmt])>pedestal[ipmt]) {
541b42c4 454 adc[ipmt] = chargeQT0[ipmt] - chargeQT1[ipmt];
73df58ab 455 }
d0bcd1fb 456 else
457 adc[ipmt] = 0;
612737bb 458 // time[ipmt] = fCalib-> WalkCorrection(refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
c8d939c8 459 Int_t refAmp = Int_t (fTime0vertex[ipmt]);
460 time[ipmt] = fCalib-> WalkCorrection( refAmp, ipmt, Int_t(adc[ipmt]), timeCFD[ipmt] ) ;
d0bcd1fb 461 Double_t sl = timeLED[ipmt] - timeCFD[ipmt];
8f620945 462 // time[ipmt] = fCalib-> WalkCorrection( refAmp,ipmt, Int_t(sl), timeCFD[ipmt] ) ;
73df58ab 463 AliDebug(5,Form(" ipmt %i QTC %i , time in chann %i (led-cfd) %i ",
d0bcd1fb 464 ipmt, Int_t(adc[ipmt]) ,Int_t(time[ipmt]),Int_t( sl)));
644a358e 465 Double_t ampMip = 0;
466 TGraph * ampGraph = (TGraph*)fAmpLED.At(ipmt);
467 if (ampGraph) ampMip = ampGraph->Eval(sl);
468 Double_t qtMip = 0;
469 TGraph * qtGraph = (TGraph*)fQTC.At(ipmt);
470 if (qtGraph) qtMip = qtGraph->Eval(adc[ipmt]);
345f03db 471 AliDebug(10,Form(" Amlitude in MIPS LED %f ; QTC %f; in channels %f\n ",ampMip,qtMip, adc[ipmt]));
5773ee77 472 // if( qtMip>lowAmpThreshold && qtMip<highAmpThreshold )
c8d939c8 473 {
474 if( equalize ==0 )
475 frecpoints.SetTime(ipmt, Float_t(time[ipmt]) );
5773ee77 476 else
477 frecpoints.SetTime(ipmt, Float_t(time[ipmt] + fTime0vertex[ipmt]) );
c8d939c8 478 // frecpoints.SetTime(ipmt, Float_t(time[ipmt] ) );
5773ee77 479 }
ce50812a 480 frecpoints.SetAmp(ipmt, Double32_t( qtMip));
481 adcmip[ipmt]=qtMip;
5773ee77 482 frecpoints.SetAmpLED(ipmt, Double32_t(ampMip));
ce50812a 483 noncalibtime[ipmt]= Double32_t (timeCFD[ipmt]);
2e6a5ee0 484 }
485 else {
ce50812a 486 time[ipmt] = -9999;
2e6a5ee0 487 adc[ipmt] = 0;
b0ab3f59 488 adcmip[ipmt] = 0;
ce50812a 489 noncalibtime[ipmt] = -9999;
2e6a5ee0 490 }
491 }
ce50812a 492 fESDTZEROfriend->SetT0timeCorr(noncalibtime) ;
493
2e6a5ee0 494 for (Int_t ipmt=0; ipmt<12; ipmt++){
ce50812a 495 if(time[ipmt] !=0 && time[ipmt] > -9000
496 /*&& badpmt[ipmt]==0 */
497 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold )
73df58ab 498 {
ce50812a 499 if(time[ipmt]<besttimeC) besttimeC=time[ipmt]; //timeC
500 if(TMath::Abs(time[ipmt])<TMath::Abs(besttimeC_best))
501 besttimeC_best=time[ipmt]; //timeC
2e6a5ee0 502 }
2e6a5ee0 503 }
73df58ab 504 for ( Int_t ipmt=12; ipmt<24; ipmt++)
505 {
ce50812a 506 if(time[ipmt] != 0 && time[ipmt] > -9000
507 /* && badpmt[ipmt]==0*/
5773ee77 508 && adcmip[ipmt]>lowAmpThreshold && adcmip[ipmt]<highAmpThreshold)
73df58ab 509 {
ce50812a 510 if(time[ipmt]<besttimeA) besttimeA=time[ipmt];
511 if(TMath::Abs(time[ipmt] ) < TMath::Abs(besttimeA_best))
512 besttimeA_best=time[ipmt]; //timeA
73df58ab 513 }
2e6a5ee0 514 }
ce50812a 515
516 if(besttimeA < 999999 && besttimeA!=0 ) {
517 if( equalize ==0 )
518 frecpoints.SetTime1stA((besttimeA * channelWidth)- 1000.*fLatencyHPTDC + 1000.*fLatencyL1A - 1000.*fGRPdelays - fTimeMeanShift[1] );
519 else
520 {
521 frecpoints.SetTimeBestA((besttimeA_best * channelWidth ));
522 frecpoints.SetTime1stA((besttimeA * channelWidth - fTimeMeanShift[1]));
523 }
524 }
525 if( besttimeC < 999999 && besttimeC!=0) {
526 if( equalize ==0 )
527 frecpoints.SetTime1stC((besttimeC * channelWidth)- 1000.*fLatencyHPTDC +1000.*fLatencyL1C - 1000.*fGRPdelays - fTimeMeanShift[2]);
528 else
529 {
530 frecpoints.SetTimeBestC((besttimeC_best * channelWidth ));
531 frecpoints.SetTime1stC((besttimeC * channelWidth - fTimeMeanShift[2]));
532 }
533 }
534 AliDebug(5,Form("1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f ",
535 besttimeA, besttimeA_best,
536 besttimeC, besttimeC_best) );
537 AliDebug(5,Form("fRecPoints::: 1stimeA %f , besttimeA %f 1sttimeC %f besttimeC %f shiftA %f shiftC %f ",
538 frecpoints.Get1stTimeA(), frecpoints.GetBestTimeA(),
539 frecpoints.Get1stTimeC(), frecpoints.GetBestTimeC(),
540 fTimeMeanShift[1],fTimeMeanShift[2] ) );
541 if( besttimeC < 999999 && besttimeA < 999999) {
c8d939c8 542 if(equalize ==0 )
597a8434 543 timeclock = (channelWidth*(besttimeC + besttimeA)/2.- 1000.*fLatencyHPTDC +1000.*fLatencyL1 - 1000.*fGRPdelays - fTimeMeanShift[0]);
c8d939c8 544 else
545 {
546 timeclock = channelWidth * Float_t( besttimeA+besttimeC)/2. - fTimeMeanShift[0];
547 meanTime = channelWidth * Float_t(besttimeA_best + besttimeC_best )/2.;
548 }
612737bb 549 timeDiff = ( besttimeA - besttimeC)* 0.001* channelWidth ;
9480f05f 550 vertex = meanVertex - c*(timeDiff)/2. ; //+ (fdZonA - fdZonC)/2;
ce50812a 551 }
552
776de217 553 } //if phys event
ce50812a 554 AliDebug(1,Form(" timeDiff %f #channel, meanTime %f #ps, TOFmean%f vertex %f cm meanVertex %f \n",timeDiff, meanTime,timeclock, vertex,meanVertex));
555 frecpoints.SetT0clock(timeclock);
556 frecpoints.SetVertex(vertex);
557 frecpoints.SetMeanTime(meanTime);
558 frecpoints.SetOnlineMean(Int_t(onlineMean));
adf36b9d 559
ce50812a 560 // Set triggers
adf36b9d 561 Bool_t tr[5];
92f6bfd3 562 Int_t trchan[5] = {50,51,52,55,56};
563 Float_t lowtr[5] = {meanTVDC-700, meanOrA-700, meanOrC-700, meanOrC-1000, meanOrC-1000 };
564 Float_t hightr[5] = {meanTVDC+700, meanOrA+700, meanOrC+700, meanOrC+1000, meanOrC+1000};
ce50812a 565
82540e64 566 for (Int_t i=0; i<5; i++) tr[i] = false;
adf36b9d 567 for (Int_t itr=0; itr<5; itr++) {
b0e13b29 568 for (Int_t iHit=0; iHit<1; iHit++)
38cbfa7c 569 {
570 Int_t trr=trchan[itr];
92f6bfd3 571 if( allData[trr][iHit] > lowtr[itr] && allData[trr][iHit] < hightr[itr]) tr[itr]=true;
ce50812a 572
573 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 574 }
adf36b9d 575 }
7e08771f 576 frecpoints.SetT0Trig(tr);
b0e13b29 577
ce50812a 578 // all times with amplitude correction
579 Float_t timecent;
b0e13b29 580 for (Int_t iHit=0; iHit<5; iHit++)
581 {
ce50812a 582 timefull = timecent = -9999;
82540e64 583 tvdc = ora = orc = -9999;
b0e13b29 584 if(allData[50][iHit]>0)
585 tvdc = (Float_t(allData[50][iHit]) - meanTVDC) * channelWidth* 0.001;
586 if(allData[51][iHit]>0)
587 ora = (Float_t(allData[51][iHit]) - meanOrA) * channelWidth* 0.001;
588
589 if(allData[52][iHit]>0)
590 orc = (Float_t(allData[52][iHit]) - meanOrC) * channelWidth* 0.001;
591
7e08771f 592 frecpoints.SetOrC( iHit, orc);
593 frecpoints.SetOrA( iHit, ora);
594 frecpoints.SetTVDC( iHit, tvdc);
b0e13b29 595 for (Int_t i0=0; i0<12; i0++) {
ce50812a 596 if (equalize ==0 )
597 timecent = fTime0vertex[i0] + timeDelayCFD[i0];
598 else
599 timecent = fTime0vertex[i0];
b0e13b29 600 timefull = -9999;
601 if(allData[i0+1][iHit]>1)
ce50812a 602 timefull = (Float_t(allData[i0+1][iHit]) - timecent)* channelWidth* 0.001;
7e08771f 603 frecpoints.SetTimeFull(i0, iHit,timefull) ;
82540e64 604 // 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 605
606 }
607
608 for (Int_t i0=12; i0<24; i0++) {
609 timefull = -9999;
ce50812a 610 if (equalize ==0 )
611 timecent = fTime0vertex[i0] + timeDelayCFD[i0];
612 else
613 timecent = fTime0vertex[i0];
94b4f7e2 614 if(allData[i0+45][iHit]>1) {
ce50812a 615 timefull = (Float_t(allData[i0+45][iHit]) - timecent)* channelWidth* 0.001;
94b4f7e2 616 }
82540e64 617 // 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 618 frecpoints.SetTimeFull(i0, iHit, timefull) ;
b0e13b29 619 }
620 }
621
622
669dc07f 623 //Set MPD
624 if(allData[53][0]>0 && allData[54][0])
7e08771f 625 frecpoints.SetMultA(allData[53][0]-allData[54][0]);
b0e13b29 626 if(allData[105][0]>0 && allData[106][0])
7e08771f 627 frecpoints.SetMultC(allData[105][0]-allData[106][0]);
b0e13b29 628
629
669dc07f 630 } // if (else )raw data
58bd3a16 631 recTree->Fill();
58bd3a16 632}
adf36b9d 633
634
ce50812a 635//____________________________________________________________
adf36b9d 636
637 void AliT0Reconstructor::FillESD(TTree */*digitsTree*/, TTree *clustersTree, AliESDEvent *pESD) const
638 {
dc7ca31d 639
640 /***************************************************
641 Resonstruct digits to vertex position
642 ****************************************************/
643
dc7ca31d 644 AliDebug(1,Form("Start FillESD T0"));
b0ab3f59 645 if(!pESD) {
646 AliError("No ESD Event");
647 return;
648 }
4cbe597e 649 pESD ->SetT0spread(fTimeSigmaShift);
650
36cde487 651
58bd3a16 652 Float_t channelWidth = fParam->GetChannelWidth() ;
f7c2c2fc 653 Float_t c = 0.0299792458; // cm/ps
adf36b9d 654 Float_t currentVertex=0, shift=0;
291f31a1 655 Int_t ncont=-1;
adf36b9d 656 const AliESDVertex* vertex = pESD->GetPrimaryVertex();
657 if (!vertex) vertex = pESD->GetPrimaryVertexSPD();
658 if (!vertex) vertex = pESD->GetPrimaryVertexTPC();
659 if (!vertex) vertex = pESD->GetVertex();
660
661 if (vertex) {
662 AliDebug(2, Form("Got %s (%s) from ESD: %f",
663 vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
664 currentVertex = vertex->GetZ();
665
666 ncont = vertex->GetNContributors();
291f31a1 667 if(ncont>0 ) {
adf36b9d 668 shift = currentVertex/c;
adf36b9d 669 }
670 }
d76c31f4 671 TTree *treeR = clustersTree;
dc7ca31d 672
7e08771f 673 AliT0RecPoint frecpoints;
674 AliT0RecPoint * pfrecpoints = &frecpoints;
dc7ca31d 675
676 AliDebug(1,Form("Start FillESD T0"));
677 TBranch *brRec = treeR->GetBranch("T0");
678 if (brRec) {
7e08771f 679 brRec->SetAddress(&pfrecpoints);
dc7ca31d 680 }else{
f16935f7 681 AliError(Form("EXEC Branch T0 rec not found"));
dc7ca31d 682 return;
683 }
73df58ab 684
685 brRec->GetEntry(0);
686 Double32_t amp[24], time[24], ampQTC[24], timecorr[24];
687 Double32_t* tcorr;
688 for(Int_t i=0; i<24; i++)
689 amp[i]=time[i]=ampQTC[i]=timecorr[i]=0;
690
ce50812a 691 //1st time
73df58ab 692 Double32_t timeClock[3];
7e08771f 693 Double32_t zPosition = frecpoints.GetVertex();
830172c8 694
ce50812a 695 timeClock[0] = frecpoints.GetT0clock() ;
696 timeClock[1] = frecpoints.Get1stTimeA() + shift;
697 timeClock[2] = frecpoints.Get1stTimeC() - shift;
698 //best time
699 Double32_t timemean[3];
700 timemean[0] = frecpoints.GetMeanTime();
701 timemean[1] = frecpoints.GetBestTimeA() + shift;
702 timemean[2] = frecpoints.GetBestTimeC() - shift;
703
704 for(Int_t i=0; i<3; i++) {
705 fESDTZERO->SetT0TOF(i,timeClock[i]); // interaction time (ns)
706 fESDTZERO->SetT0TOFbest(i,timemean[i]); // interaction time (ns)
707 }
73df58ab 708 for ( Int_t i=0; i<24; i++) {
7e08771f 709 time[i] = frecpoints.GetTime(i); // ps to ns
ce50812a 710 if ( time[i] != 0 && time[i]>-9999) {
7e08771f 711 ampQTC[i] = frecpoints.GetAmp(i);
712 amp[i] = frecpoints.AmpLED(i);
612737bb 713 AliDebug(1,Form("T0: %i time %f ampQTC %f ampLED %f \n", i, time[i], ampQTC[i], amp[i]));
85f61e3b 714 }
73df58ab 715 }
ce50812a 716 fESDTZERO->SetT0time(time); // best TOF on each PMT
717 fESDTZERO->SetT0amplitude(ampQTC); // number of particles(MIPs) on each PMT
7e08771f 718 Int_t trig= frecpoints.GetT0Trig();
719 frecpoints.PrintTriggerSignals( trig);
92f6bfd3 720 // printf(" !!!!! FillESD trigger %i \n",trig);
b0e13b29 721 fESDTZERO->SetT0Trig(trig);
b0e13b29 722 fESDTZERO->SetT0zVertex(zPosition); //vertex Z position
669dc07f 723
7e08771f 724 Double32_t multA=frecpoints.GetMultA();
725 Double32_t multC=frecpoints.GetMultC();
b0e13b29 726 fESDTZERO->SetMultA(multA); // for backward compatubility
727 fESDTZERO->SetMultC(multC); // for backward compatubility
728
729
730 for (Int_t iHit =0; iHit<5; iHit++ ) {
82540e64 731 AliDebug(10,Form("FillESD ::: iHit %i tvdc %f orA %f orC %f\n", iHit,
7e08771f 732 frecpoints.GetTVDC(iHit),
733 frecpoints.GetOrA(iHit),
734 frecpoints.GetOrC(iHit) ));
735 fESDTZERO->SetTVDC(iHit,frecpoints.GetTVDC(iHit));
736 fESDTZERO->SetOrA(iHit,frecpoints.GetOrA(iHit));
737 fESDTZERO->SetOrC(iHit,frecpoints.GetOrC(iHit));
b0e13b29 738
739 for (Int_t i0=0; i0<24; i0++) {
7e08771f 740 // if(frecpoints.GetTimeFull(i0,iHit)>0){
741 // printf("FillESD ::: iHit %i cfd %i time %f \n", iHit, i0, frecpoints.GetTimeFull(i0,iHit));
742 fESDTZERO->SetTimeFull(i0, iHit,frecpoints.GetTimeFull(i0,iHit));
94b4f7e2 743 // }
b0e13b29 744
745 }
746 }
73df58ab 747
ce50812a 748 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 749
ea1a8005 750 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
751 // background flags
752 Bool_t background = BackgroundFlag();
753 fESDTZERO->SetBackgroundFlag(background);
754 Bool_t pileup = PileupFlag();
755 fESDTZERO->SetPileupFlag(pileup);
528890e5 756 for (Int_t i=0; i<5; i++) {
7e08771f 757 fESDTZERO->SetPileupTime(i, frecpoints.GetTVDC(i) ) ;
758 // printf("!!!!!! FillESD :: pileup %i %f %f \n", i,fESDTZERO->GetPileupTime(i), frecpoints.GetTVDC(i));
528890e5 759 }
ea1a8005 760 Bool_t sat = SatelliteFlag();
761 fESDTZERO->SetSatelliteFlag(sat);
73df58ab 762
ea1a8005 763 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
73df58ab 764 if (pESD) {
adf36b9d 765
73df58ab 766 AliESDfriend *fr = (AliESDfriend*)pESD->FindListObject("AliESDfriend");
767 if (fr) {
82540e64 768 AliDebug(10, Form("Writing TZERO friend data to ESD tree"));
73df58ab 769
85f61e3b 770 // if (ncont>2) {
73df58ab 771 tcorr = fESDTZEROfriend->GetT0timeCorr();
772 for ( Int_t i=0; i<24; i++) {
85f61e3b 773 if(i<12 && time[i]>1) timecorr[i] = tcorr[i] - shift/channelWidth;
774 if(i>11 && time[i]>1) timecorr[i] = tcorr[i] + shift/channelWidth;
612737bb 775 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 776 }
73df58ab 777 fESDTZEROfriend->SetT0timeCorr( timecorr) ;
778 fESDTZEROfriend->SetT0ampLEDminCFD(amp);
779 fESDTZEROfriend->SetT0ampQTC(ampQTC);
85f61e3b 780 fr->SetTZEROfriend(fESDTZEROfriend);
781 // }//
58bd3a16 782 }
b0e13b29 783
784 pESD->SetTZEROData(fESDTZERO);
58bd3a16 785 }
73df58ab 786
dc7ca31d 787} // vertex in 3 sigma
788
ea1a8005 789//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
790 //____________________________________________________________
791
792Bool_t AliT0Reconstructor::PileupFlag() const
793{
794 //
795 Bool_t pileup = false;
528890e5 796 Float_t tvdc[5];
ea1a8005 797 for (Int_t ih=0; ih<5; ih++)
798 {
799 tvdc[ih] = fESDTZERO->GetTVDC(ih);
800
528890e5 801 if( tvdc[0] !=0 && tvdc[0]> -10 && tvdc[0]< 10 )
ea1a8005 802 if(ih>0 && tvdc[ih]>20 ) pileup = true;
528890e5 803 if( tvdc[0] >20 || (tvdc[0] < -20 && tvdc[0] > -9000) ) pileup =true;
804 // if (pileup) printf(" !!!!! pile up %i tvdc %f \n",ih, tvdc[ih]);
ea1a8005 805 }
806
807
808 return pileup;
809
810}
811
812 //____________________________________________________________
813
814Bool_t AliT0Reconstructor::BackgroundFlag() const
815{
9120829d 816
ea1a8005 817 Bool_t background = false;
9120829d 818 /*
819 Float_t orA = fESDTZERO->GetOrA(0);
820 Float_t orC = fESDTZERO->GetOrC(0);
821 Float_t tvdc = fESDTZERO->GetTVDC(ih);
822
823 if ( (orA > -5 && orA <5) && (orC > -5 && orC <5) && (tvdc < -5 || tvdc > 5)) {
824 background = true;
825 // printf(" orA %f orC %f tvdc %f\n", orA, orC, tvdc);
826 }
827 */
ea1a8005 828 return background;
829
830
831}
832
833
834 //____________________________________________________________
835
836Bool_t AliT0Reconstructor::SatelliteFlag() const
837{
838
ce50812a 839 Float_t satelliteLow = GetRecoParam() -> GetLowSatelliteThreshold();
840 Float_t satelliteHigh = GetRecoParam() -> GetHighSatelliteThreshold();
ea1a8005 841 Bool_t satellite = false;
842 for (Int_t i0=0; i0<24; i0++) {
843 Float_t timefull = fESDTZERO -> GetTimeFull(i0,0);
ce50812a 844 if( timefull > satelliteLow && timefull < satelliteHigh) satellite=true;
ea1a8005 845 }
846
847 return satellite;
848
849}