Adding a bit that specifies if the online (beam-beam,beam-gas) flags are filled ...
[u/mrichter/AliRoot.git] / VZERO / AliVZEROReconstructor.cxx
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$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 ///                                                                          //
20 /// class for VZERO reconstruction                                           //
21 ///                                                                          //
22 ///////////////////////////////////////////////////////////////////////////////
23
24 #include <TH1F.h>
25 #include <TF1.h>
26 #include <TParameter.h>
27
28 #include "AliRunLoader.h"
29 #include "AliRawReader.h"
30 #include "AliGRPObject.h"
31 #include "AliCDBManager.h"
32 #include "AliCDBStorage.h"
33 #include "AliCDBEntry.h"
34 #include "AliVZEROReconstructor.h"
35 #include "AliVZERORawStream.h"
36 #include "AliVZEROConst.h"
37 #include "AliESDEvent.h"
38 #include "AliVZEROTriggerMask.h"
39 #include "AliESDfriend.h"
40 #include "AliESDVZEROfriend.h"
41 #include "AliVZEROdigit.h"
42 #include "AliVZEROCalibData.h"
43 #include "AliRunInfo.h"
44 #include "AliCTPTimeParams.h"
45
46 ClassImp(AliVZEROReconstructor)
47
48 //_____________________________________________________________________________
49 AliVZEROReconstructor:: AliVZEROReconstructor(): AliReconstructor(),
50                         fESDVZERO(0x0),
51                         fESD(0x0),
52                         fESDVZEROfriend(0x0),
53                         fCalibData(NULL),
54                         fTimeSlewing(NULL),
55                         fCollisionMode(0),
56                         fBeamEnergy(0.),
57                         fDigitsArray(0)
58 {
59   // Default constructor  
60   // Get calibration data
61   
62   fCalibData = GetCalibData();
63
64   AliCDBEntry *entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
65   if (!entry) AliFatal("CTP timing parameters are not found in OCDB !");
66   AliCTPTimeParams *ctpParams = (AliCTPTimeParams*)entry->GetObject();
67   Float_t l1Delay = (Float_t)ctpParams->GetDelayL1L0()*25.0;
68
69   AliCDBEntry *entry2 = AliCDBManager::Instance()->Get("VZERO/Calib/TimeDelays");
70   if (!entry2) AliFatal("VZERO time delays are not found in OCDB !");
71   TH1F *delays = (TH1F*)entry2->GetObject();
72
73   AliCDBEntry *entry3 = AliCDBManager::Instance()->Get("VZERO/Calib/TimeSlewing");
74   if (!entry3) AliFatal("VZERO time slewing function is not found in OCDB !");
75   fTimeSlewing = (TF1*)entry3->GetObject();
76
77   for(Int_t i = 0 ; i < 64; ++i) {
78     Int_t board = AliVZEROCalibData::GetBoardNumber(i);
79     fTimeOffset[i] = (((Float_t)fCalibData->GetTriggerCountOffset(board)-
80                         (Float_t)fCalibData->GetRollOver(board))*25.0+
81                        fCalibData->GetTimeOffset(i)+
82                        l1Delay+
83                        delays->GetBinContent(i+1)+
84                        kV0Offset);
85   }
86 }
87
88
89 //_____________________________________________________________________________
90 AliVZEROReconstructor& AliVZEROReconstructor::operator = 
91   (const AliVZEROReconstructor& /*reconstructor*/)
92 {
93 // assignment operator
94
95   Fatal("operator =", "assignment operator not implemented");
96   return *this;
97 }
98
99 //_____________________________________________________________________________
100 AliVZEROReconstructor::~AliVZEROReconstructor()
101 {
102 // destructor
103
104    delete fESDVZERO;
105    delete fESDVZEROfriend;
106    delete fDigitsArray;
107 }
108
109 //_____________________________________________________________________________
110 void AliVZEROReconstructor::Init()
111 {
112 // initializer
113
114   fESDVZERO  = new AliESDVZERO;
115   fESDVZEROfriend = new AliESDVZEROfriend;
116   
117   GetCollisionMode();  // fCollisionMode =1 for Pb-Pb simulated data
118 }
119
120 //______________________________________________________________________
121 void AliVZEROReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const
122 {
123 // converts RAW to digits 
124
125   if (!digitsTree) {
126     AliError("No digits tree!");
127     return;
128   }
129
130   if (!fDigitsArray)
131     fDigitsArray = new TClonesArray("AliVZEROdigit", 64);
132   digitsTree->Branch("VZERODigit", &fDigitsArray);
133
134   fESDVZEROfriend->Reset();
135
136   rawReader->Reset();
137   AliVZERORawStream rawStream(rawReader);
138   if (rawStream.Next()) { 
139
140     Int_t aBBflagsV0A = 0;
141     Int_t aBBflagsV0C = 0;
142     Int_t aBGflagsV0A = 0;
143     Int_t aBGflagsV0C = 0;
144
145     for(Int_t iChannel=0; iChannel < 64; ++iChannel) {
146       Int_t offlineCh = rawStream.GetOfflineChannel(iChannel);
147       // ADC charge samples
148       Short_t chargeADC[AliVZEROdigit::kNClocks];
149       for(Int_t iClock=0; iClock < AliVZEROdigit::kNClocks; ++iClock) {
150         chargeADC[iClock] = rawStream.GetPedestal(iChannel,iClock);
151       }
152       // Integrator flag
153       Bool_t integrator = rawStream.GetIntegratorFlag(iChannel,AliVZEROdigit::kNClocks/2);
154       // Beam-beam and beam-gas flags
155       if(offlineCh<32) {
156         if (rawStream.GetBBFlag(iChannel,AliVZEROdigit::kNClocks/2)) aBBflagsV0C |= (1 << offlineCh);
157         if (rawStream.GetBGFlag(iChannel,AliVZEROdigit::kNClocks/2)) aBGflagsV0C |= (1 << offlineCh);
158       } else {
159         if (rawStream.GetBBFlag(iChannel,AliVZEROdigit::kNClocks/2)) aBBflagsV0A |= (1 << (offlineCh-32));
160         if (rawStream.GetBGFlag(iChannel,AliVZEROdigit::kNClocks/2)) aBGflagsV0A |= (1 << (offlineCh-32));
161       }
162       // HPTDC data (leading time and width)
163       Int_t board = AliVZEROCalibData::GetBoardNumber(offlineCh);
164       Float_t time = rawStream.GetTime(iChannel)*fCalibData->GetTimeResolution(board);
165       Float_t width = rawStream.GetWidth(iChannel)*fCalibData->GetWidthResolution(board);
166       // Add a digit
167       if(!fCalibData->IsChannelDead(iChannel)){
168           new ((*fDigitsArray)[fDigitsArray->GetEntriesFast()])
169             AliVZEROdigit(offlineCh, 0.0, time,
170                           width,integrator,
171                           chargeADC);
172       }
173
174       // Filling the part of esd friend object that is available only for raw data
175       fESDVZEROfriend->SetBBScalers(offlineCh,rawStream.GetBBScalers(iChannel));
176       fESDVZEROfriend->SetBGScalers(offlineCh,rawStream.GetBGScalers(iChannel));
177       for (Int_t iBunch = 0; iBunch < AliESDVZEROfriend::kNBunches; iBunch++) {
178         fESDVZEROfriend->SetChargeMB(offlineCh,iBunch,rawStream.GetChargeMB(iChannel,iBunch));
179         fESDVZEROfriend->SetIntMBFlag(offlineCh,iBunch,rawStream.GetIntMBFlag(iChannel,iBunch));
180         fESDVZEROfriend->SetBBMBFlag(offlineCh,iBunch,rawStream.GetBBMBFlag(iChannel,iBunch));
181         fESDVZEROfriend->SetBGMBFlag(offlineCh,iBunch,rawStream.GetBGMBFlag(iChannel,iBunch));
182       }
183
184     }  
185
186     // Filling the global part of esd friend object that is available only for raw data
187     fESDVZEROfriend->SetTriggerInputs(rawStream.GetTriggerInputs());
188     fESDVZEROfriend->SetTriggerInputsMask(rawStream.GetTriggerInputsMask());
189
190     for(Int_t iScaler = 0; iScaler < AliESDVZEROfriend::kNScalers; iScaler++)
191       fESDVZEROfriend->SetTriggerScalers(iScaler,rawStream.GetTriggerScalers(iScaler));
192
193     for(Int_t iBunch = 0; iBunch < AliESDVZEROfriend::kNBunches; iBunch++)
194       fESDVZEROfriend->SetBunchNumbersMB(iBunch,rawStream.GetBunchNumbersMB(iBunch));
195      
196     // Store the BB and BG flags in the digits tree (user info)
197     digitsTree->GetUserInfo()->Add(new TParameter<int>("BBflagsV0A",aBBflagsV0A));
198     digitsTree->GetUserInfo()->Add(new TParameter<int>("BBflagsV0C",aBBflagsV0C));
199     digitsTree->GetUserInfo()->Add(new TParameter<int>("BGflagsV0A",aBGflagsV0A));
200     digitsTree->GetUserInfo()->Add(new TParameter<int>("BGflagsV0C",aBGflagsV0C));
201
202     digitsTree->Fill();
203   }
204
205   fDigitsArray->Clear();
206 }      
207
208 //______________________________________________________________________
209 void AliVZEROReconstructor::FillESD(TTree* digitsTree, TTree* /*clustersTree*/,
210                                     AliESDEvent* esd) const
211 {
212 // fills multiplicities to the ESD - pedestal is now subtracted
213     
214   if (!digitsTree) {
215       AliError("No digits tree!");
216       return;
217   }
218
219   TBranch* digitBranch = digitsTree->GetBranch("VZERODigit");
220   digitBranch->SetAddress(&fDigitsArray);
221
222   Float_t   mult[64];  
223   Float_t    adc[64]; 
224   Float_t   time[64]; 
225   Float_t  width[64];
226   Bool_t aBBflag[64];
227   Bool_t aBGflag[64];
228    
229   for (Int_t i=0; i<64; i++){
230        adc[i]    = 0.0;
231        mult[i]   = 0.0;
232        time[i]   = kInvalidTime;
233        width[i]  = 0.0;
234        aBBflag[i] = kFALSE;
235        aBGflag[i] = kFALSE;
236   }
237      
238   Int_t aBBflagsV0A = 0;
239   Int_t aBBflagsV0C = 0;
240   Int_t aBGflagsV0A = 0;
241   Int_t aBGflagsV0C = 0;
242
243   if (digitsTree->GetUserInfo()->FindObject("BBflagsV0A")) {
244     aBBflagsV0A = ((TParameter<int>*)digitsTree->GetUserInfo()->FindObject("BBflagsV0A"))->GetVal();
245   }
246   else
247     AliWarning("V0A beam-beam flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
248
249   if (digitsTree->GetUserInfo()->FindObject("BBflagsV0C")) {
250     aBBflagsV0C = ((TParameter<int>*)digitsTree->GetUserInfo()->FindObject("BBflagsV0C"))->GetVal();
251   }
252   else
253     AliWarning("V0C beam-beam flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
254
255   if (digitsTree->GetUserInfo()->FindObject("BGflagsV0A")) {
256     aBGflagsV0A = ((TParameter<int>*)digitsTree->GetUserInfo()->FindObject("BGflagsV0A"))->GetVal();
257   }
258   else
259     AliWarning("V0A beam-gas flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
260
261   if (digitsTree->GetUserInfo()->FindObject("BGflagsV0C")) {
262     aBGflagsV0C = ((TParameter<int>*)digitsTree->GetUserInfo()->FindObject("BGflagsV0C"))->GetVal();
263   }
264   else
265     AliWarning("V0C beam-gas flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
266   
267   // Beam-beam and beam-gas flags (hardware)
268   for (Int_t iChannel = 0; iChannel < 64; ++iChannel) {
269     if(iChannel < 32) {
270       aBBflag[iChannel] = (aBBflagsV0C >> iChannel) & 0x1;
271       aBGflag[iChannel] = (aBGflagsV0C >> iChannel) & 0x1;
272     }
273     else {
274       aBBflag[iChannel] = (aBBflagsV0A >> (iChannel-32)) & 0x1;
275       aBGflag[iChannel] = (aBGflagsV0A >> (iChannel-32)) & 0x1;
276     }
277   }
278
279   Int_t nEntries = (Int_t)digitsTree->GetEntries();
280   for (Int_t e=0; e<nEntries; e++) {
281     digitsTree->GetEvent(e);
282
283     Int_t nDigits = fDigitsArray->GetEntriesFast();
284     
285     for (Int_t d=0; d<nDigits; d++) {    
286         AliVZEROdigit* digit = (AliVZEROdigit*) fDigitsArray->At(d);      
287         Int_t  pmNumber = digit->PMNumber();
288
289         // Pedestal retrieval and suppression
290         Bool_t integrator = digit->Integrator();
291         Float_t maxadc = 0;
292         Int_t imax = -1;
293         Float_t adcPedSub[AliVZEROdigit::kNClocks];
294         for(Int_t iClock=0; iClock < AliVZEROdigit::kNClocks; ++iClock) {
295           Short_t charge = digit->ChargeADC(iClock);
296           Bool_t iIntegrator = (iClock%2 == 0) ? integrator : !integrator;
297           Int_t k = pmNumber + 64*iIntegrator;
298           adcPedSub[iClock] = (Float_t)charge - fCalibData->GetPedestal(k);
299           if(adcPedSub[iClock] <= GetRecoParam()->GetNSigmaPed()*fCalibData->GetSigma(k)) {
300             adcPedSub[iClock] = 0;
301             continue;
302           }
303           if(iClock < GetRecoParam()->GetStartClock() || iClock > GetRecoParam()->GetEndClock()) continue;
304           if(adcPedSub[iClock] > maxadc) {
305             maxadc = adcPedSub[iClock];
306             imax   = iClock;
307           }
308         }
309
310         if (imax != -1) {
311           Int_t start = imax - GetRecoParam()->GetNPreClocks();
312           if (start < 0) start = 0;
313           Int_t end = imax + GetRecoParam()->GetNPostClocks();
314           if (end > 20) end = 20;
315           for(Int_t iClock = start; iClock <= end; iClock++) {
316             adc[pmNumber] += adcPedSub[iClock];
317           }
318         }
319
320         // HPTDC leading time and width
321         // Correction for slewing and various time delays
322         time[pmNumber]  =  CorrectLeadingTime(pmNumber,digit->Time(),adc[pmNumber]);
323         width[pmNumber] =  digit->Width();
324
325         if (adc[pmNumber] > 0) AliDebug(1,Form("PM = %d ADC = %f TDC %f (%f)   Int %d (%d %d %d %d %d)    %f %f   %f %f    %d %d",pmNumber, adc[pmNumber],
326                                                digit->Time(),time[pmNumber],
327                                                integrator,
328                                                digit->ChargeADC(8),digit->ChargeADC(9),digit->ChargeADC(10),
329                                                digit->ChargeADC(11),digit->ChargeADC(12),
330                                                fCalibData->GetPedestal(pmNumber),fCalibData->GetSigma(pmNumber),
331                                                fCalibData->GetPedestal(pmNumber+64),fCalibData->GetSigma(pmNumber+64),
332                                                aBBflag[pmNumber],aBGflag[pmNumber]));
333
334         mult[pmNumber] = adc[pmNumber]*fCalibData->GetMIPperADC(pmNumber);
335
336         // Fill ESD friend object
337         for (Int_t iEv = 0; iEv < AliESDVZEROfriend::kNEvOfInt; iEv++) {
338           fESDVZEROfriend->SetPedestal(pmNumber,iEv,(Float_t)digit->ChargeADC(iEv));
339           fESDVZEROfriend->SetIntegratorFlag(pmNumber,iEv,(iEv%2 == 0) ? integrator : !integrator);
340           fESDVZEROfriend->SetBBFlag(pmNumber,iEv,aBBflag[pmNumber]);
341           fESDVZEROfriend->SetBGFlag(pmNumber,iEv,aBGflag[pmNumber]);
342         }
343         fESDVZEROfriend->SetTime(pmNumber,digit->Time());
344         fESDVZEROfriend->SetWidth(pmNumber,digit->Width());
345
346     } // end of loop over digits
347   } // end of loop over events in digits tree
348          
349   fESDVZERO->SetBit(AliESDVZERO::kCorrectedLeadingTime,kTRUE);
350   fESDVZERO->SetMultiplicity(mult);
351   fESDVZERO->SetADC(adc);
352   fESDVZERO->SetTime(time);
353   fESDVZERO->SetWidth(width);
354   fESDVZERO->SetBit(AliESDVZERO::kOnlineBitsFilled,kTRUE);
355   fESDVZERO->SetBBFlag(aBBflag);
356   fESDVZERO->SetBGFlag(aBGflag);
357
358   // now fill the V0 decision and channel flags
359   {
360     AliVZEROTriggerMask triggerMask;
361     triggerMask.FillMasks(fESDVZERO, fCalibData, fTimeSlewing);
362   }
363
364   if (esd) { 
365      AliDebug(1, Form("Writing VZERO data to ESD tree"));
366      esd->SetVZEROData(fESDVZERO);
367   }
368
369   if (esd) {
370      AliESDfriend *fr = (AliESDfriend*)esd->FindListObject("AliESDfriend");
371      if (fr) {
372         AliDebug(1, Form("Writing VZERO friend data to ESD tree"));
373         fr->SetVZEROfriend(fESDVZEROfriend);
374     }
375   }
376
377   fDigitsArray->Clear();
378 }
379
380 //_____________________________________________________________________________
381 AliCDBStorage* AliVZEROReconstructor::SetStorage(const char *uri) 
382 {
383 // Sets the storage  
384
385   Bool_t deleteManager = kFALSE;
386   
387   AliCDBManager *manager = AliCDBManager::Instance();
388   AliCDBStorage *defstorage = manager->GetDefaultStorage();
389   
390   if(!defstorage || !(defstorage->Contains("VZERO"))){ 
391      AliWarning("No default storage set or default storage doesn't contain VZERO!");
392      manager->SetDefaultStorage(uri);
393      deleteManager = kTRUE;
394   }
395  
396   AliCDBStorage *storage = manager->GetDefaultStorage();
397
398   if(deleteManager){
399      AliCDBManager::Instance()->UnsetDefaultStorage();
400      defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
401   }
402
403   return storage; 
404 }
405
406 //____________________________________________________________________________
407 void AliVZEROReconstructor::GetCollisionMode()
408 {
409   // Retrieval of collision mode 
410
411   TString beamType = GetRunInfo()->GetBeamType();
412   if(beamType==AliGRPObject::GetInvalidString()){
413      AliError("VZERO cannot retrieve beam type");
414      return;
415   }
416
417   if( (beamType.CompareTo("P-P") ==0)  || (beamType.CompareTo("p-p") ==0) ){
418     fCollisionMode=0;
419   }
420   else if( (beamType.CompareTo("Pb-Pb") ==0)  || (beamType.CompareTo("A-A") ==0) ){
421     fCollisionMode=1;
422   }
423     
424   fBeamEnergy = GetRunInfo()->GetBeamEnergy();
425   if(fBeamEnergy==AliGRPObject::GetInvalidFloat()) {
426      AliError("Missing value for the beam energy ! Using 0");
427      fBeamEnergy = 0.;
428   }
429   
430   AliDebug(1,Form("\n ++++++ Beam type and collision mode retrieved as %s %d @ %1.3f GeV ++++++\n\n",beamType.Data(), fCollisionMode, fBeamEnergy));
431
432 }
433
434 //_____________________________________________________________________________
435 AliVZEROCalibData* AliVZEROReconstructor::GetCalibData() const
436 {
437   // Gets calibration object for VZERO set
438
439   AliCDBManager *man = AliCDBManager::Instance();
440
441   AliCDBEntry *entry=0;
442
443   entry = man->Get("VZERO/Calib/Data");
444
445   AliVZEROCalibData *calibdata = 0;
446
447   if (entry) calibdata = (AliVZEROCalibData*) entry->GetObject();
448   if (!calibdata)  AliFatal("No calibration data from calibration database !");
449
450   return calibdata;
451 }
452
453 Float_t AliVZEROReconstructor::CorrectLeadingTime(Int_t i, Float_t time, Float_t adc) const
454 {
455   // Correct the leading time
456   // for slewing effect and
457   // misalignment of the channels
458   if (time < 1e-6) return kInvalidTime;
459
460   // Channel alignment and general offset subtraction
461   if (i < 32) time -= kV0CDelayCables;
462   time -= fTimeOffset[i];
463
464   // In case of pathological signals
465   if (adc < 1e-6) return time;
466
467   // Slewing correction
468   Float_t thr = fCalibData->GetDiscriThr(i);
469   time -= fTimeSlewing->Eval(adc/thr);
470
471   return time;
472 }