AD/ADrec/AliADReconstructor.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: AliADReconstructor.cxx 20956 2007-09-26 14:22:18Z mrodrigu $ */
  17 //////////////////////////////////////////////////////////////////////////////
  18 //                                                                          //
  19 //  Class for AD reconstruction                                         //
  20 //////////////////////////////////////////////////////////////////////////////
  21 #include <TParameter.h>
  22
  23 #include "AliRawReader.h"
  24 #include "AliGRPObject.h"
  25 #include "AliCDBManager.h"
  26 #include "AliCDBStorage.h"
  27 #include "AliCDBEntry.h"
  28 #include "AliESDEvent.h"
  29 #include "AliRunInfo.h"
  30 #include "AliCTPTimeParams.h"
  31 #include "AliLHCClockPhase.h"
  32
  33 #include "AliADReconstructor.h"
  34 #include "AliADdigit.h"
  35 #include "AliESDAD.h"
  36 #include "AliESDADfriend.h"
  37 #include "AliADConst.h"
  38 #include "AliADCalibData.h"
  39 #include "AliADRawStream.h"
  40
  41 ClassImp(AliADReconstructor)
  42 //_____________________________________________________________________________
  43 AliADReconstructor:: AliADReconstructor():
  44   AliReconstructor(),
  45   fESDAD(0x0),
  46   fESD(0x0),
  47   fESDADfriend(0x0),
  48   fCalibData(NULL),
  49   fDigitsArray(0),
  50   fCollisionMode(0),
  51   fBeamEnergy(0.)
  52
  53 {
  54   fCalibData = GetCalibData();
  55   // Default constructor
  56   // Get calibration data
  57
  58 }
  59
  60 //_____________________________________________________________________________
  61 AliADReconstructor& AliADReconstructor::operator =
  62   (const AliADReconstructor& /*reconstructor*/)
  63 {
  64 // assignment operator
  65
  66   Fatal("operator =", "assignment operator not implemented");
  67   return *this;
  68 }
  69
  70 //_____________________________________________________________________________
  71 AliADReconstructor::~AliADReconstructor()
  72 {
  73 // destructor
  74   delete fESDAD;
  75   delete fESDADfriend;
  76   delete fDigitsArray;
  77 }
  78
  79 //_____________________________________________________________________________
  80 void AliADReconstructor::Init()
  81 {
  82 // initializer
  83     fESDAD  = new AliESDAD;
  84     fESDADfriend  = new AliESDADfriend;
  85 }
  86
  87 //_____________________________________________________________________________
  88 void AliADReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const
  89 {
  90 // converts RAW to digits
  91
  92   if (!digitsTree) {
  93     AliError("No digits tree!");
  94     return;
  95   }
  96
  97   if (!fDigitsArray) fDigitsArray = new TClonesArray("AliADdigit", 16);
  98   digitsTree->Branch("ADDigit", &fDigitsArray);
  99
 100   rawReader->Reset();
 101   AliADRawStream rawStream(rawReader);
 102   if (rawStream.Next()) {
 103
 104     for(Int_t iChannel=0; iChannel < 16; ++iChannel) {
 105       Int_t offlineCh = kOfflineChannel[iChannel];
 106       // ADC charge samples
 107       Short_t chargeADC[kNClocks];
 108       for(Int_t iClock=0; iClock < kNClocks; ++iClock) {
 109         chargeADC[iClock] = rawStream.GetPedestal(iChannel,iClock);
 110       }
 111       // Integrator flag
 112       Bool_t integrator = rawStream.GetIntegratorFlag(iChannel,kNClocks/2);
 113       Bool_t BBflag = rawStream.GetBBFlag(iChannel,kNClocks/2);
 114       Bool_t BGflag = rawStream.GetBGFlag(iChannel,kNClocks/2);
 115
 116       // HPTDC data (leading time and width)
 117       Int_t board = AliADCalibData::GetBoardNumber(offlineCh);
 118       Float_t time = rawStream.GetTime(iChannel)*fCalibData->GetTimeResolution(board);
 119       Float_t width = rawStream.GetWidth(iChannel)*fCalibData->GetWidthResolution(board);
 120       // Add a digit
 121       if(!fCalibData->IsChannelDead(iChannel)){
 122           new ((*fDigitsArray)[fDigitsArray->GetEntriesFast()]) AliADdigit(offlineCh, time, width,integrator, chargeADC, BBflag, BGflag);
 123       }
 124
 125       fESDADfriend->SetTriggerInputs(rawStream.GetTriggerInputs());
 126       fESDADfriend->SetTriggerInputsMask(rawStream.GetTriggerInputsMask());
 127
 128       fESDADfriend->SetBBScalers(offlineCh,rawStream.GetBBScalers(iChannel));
 129       fESDADfriend->SetBGScalers(offlineCh,rawStream.GetBGScalers(iChannel));
 130       for (Int_t iEv = 0; iEv < kNClocks; iEv++) {
 131           fESDADfriend->SetBBFlag(offlineCh,iEv,rawStream.GetBBFlag(iChannel,iEv));
 132           fESDADfriend->SetBGFlag(offlineCh,iEv,rawStream.GetBGFlag(iChannel,iEv));
 133       }
 134     }
 135     for(Int_t iScaler = 0; iScaler < AliESDADfriend::kNScalers; iScaler++)
 136       fESDADfriend->SetTriggerScalers(iScaler,rawStream.GetTriggerScalers(iScaler));
 137
 138     digitsTree->Fill();
 139   }
 140
 141   fDigitsArray->Clear();
 142
 143 }
 144
 145 //_____________________________________________________________________________
 146 void AliADReconstructor::FillESD(TTree* digitsTree, TTree* /*clustersTree*/,AliESDEvent* esd) const
 147 {
 148
 149   printf("Running AD Reconstruction \n");
 150
 151   if (!digitsTree) {
 152       AliError("No digits tree!");
 153       return;
 154   }
 155
 156   TBranch* digitBranch = digitsTree->GetBranch("ADDigit");
 157   digitBranch->SetAddress(&fDigitsArray);
 158
 159   Float_t   mult[16];
 160   Float_t    adc[16];
 161   Float_t   time[16];
 162   Float_t  width[16];
 163   Bool_t aBBflag[16];
 164   Bool_t aBGflag[16];
 165
 166   for (Int_t i=0; i<16; i++){
 167        adc[i]    = 0.0;
 168        mult[i]   = 0.0;
 169        time[i]   = kInvalidTime;
 170        width[i]  = 0.0;
 171        aBBflag[i] = kFALSE;
 172        aBGflag[i] = kFALSE;
 173   }
 174
 175   Int_t nEntries = (Int_t)digitsTree->GetEntries();
 176   for (Int_t e=0; e<nEntries; e++) {
 177     digitsTree->GetEvent(e);
 178
 179     Int_t nDigits = fDigitsArray->GetEntriesFast();
 180
 181     for (Int_t d=0; d<nDigits; d++) {
 182         AliADdigit* digit = (AliADdigit*) fDigitsArray->At(d);
 183         Int_t  pmNumber = digit->PMNumber();
 184
 185         // Pedestal retrieval and suppression
 186         Bool_t integrator = digit->Integrator();
 187         Float_t maxadc = 0;
 188         Int_t imax = -1;
 189         Float_t adcPedSub[kNClocks];
 190         for(Int_t iClock=0; iClock < kNClocks; ++iClock) {
 191           Short_t charge = digit->ChargeADC(iClock);
 192           Bool_t iIntegrator = (iClock%2 == 0) ? integrator : !integrator;
 193           Int_t k = pmNumber + 16*iIntegrator;
 194           adcPedSub[iClock] = (Float_t)charge - fCalibData->GetPedestal(k);
 195           //if(adcPedSub[iClock] <= GetRecoParam()->GetNSigmaPed()*fCalibData->GetSigma(k)) {
 196           if(adcPedSub[iClock] <= 2*fCalibData->GetSigma(k)) {
 197             adcPedSub[iClock] = 0;
 198             continue;
 199           }
 200           //if(iClock < GetRecoParam()->GetStartClock() || iClock > GetRecoParam()->GetEndClock()) continue;
 201           if(iClock < 8 || iClock > 12) continue;
 202           if(adcPedSub[iClock] > maxadc) {
 203             maxadc = adcPedSub[iClock];
 204             imax   = iClock;
 205           }
 206         }
 207
 208         if (imax != -1) {
 209           //Int_t start = imax - GetRecoParam()->GetNPreClocks();
 210           Int_t start = imax - 2;
 211           if (start < 0) start = 0;
 212           //Int_t end = imax + GetRecoParam()->GetNPostClocks();
 213           Int_t end = imax + 1;
 214           if (end > 20) end = 20;
 215           for(Int_t iClock = start; iClock <= end; iClock++) {
 216             adc[pmNumber] += adcPedSub[iClock];
 217           }
 218         }
 219
 220         // HPTDC leading time and width
 221         // Correction for slewing and various time delays
 222         //time[pmNumber]  =  CorrectLeadingTime(pmNumber,digit->Time(),adc[pmNumber]);
 223         time[pmNumber]  =  digit->Time();
 224         width[pmNumber] =  digit->Width();
 225         aBBflag[pmNumber] = digit->GetBBflag();
 226         aBGflag[pmNumber] = digit->GetBGflag();
 227
 228         if (adc[pmNumber] > 0) {
 229           AliDebug(1,Form("PM = %d ADC = %.2f (%.2f) TDC %.2f (%.2f)   Int %d (%d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d)    %.2f %.2f   %.2f %.2f    %d %d",pmNumber, adc[pmNumber],
 230                        digit->ChargeADC(11)+digit->ChargeADC(10)+digit->ChargeADC(9)+digit->ChargeADC(8)+
 231                        digit->ChargeADC(7)+digit->ChargeADC(6)+digit->ChargeADC(5)+digit->ChargeADC(4)-
 232                        4.*fCalibData->GetPedestal(pmNumber)-4.*fCalibData->GetPedestal(pmNumber+16),
 233                           digit->Time(),time[pmNumber],
 234                           integrator,
 235                           digit->ChargeADC(0),digit->ChargeADC(1),digit->ChargeADC(2),digit->ChargeADC(3),digit->ChargeADC(4),digit->ChargeADC(5),digit->ChargeADC(6),digit->ChargeADC(7),
 236                           digit->ChargeADC(8),digit->ChargeADC(9),digit->ChargeADC(10),
 237                           digit->ChargeADC(11),digit->ChargeADC(12),
 238                           digit->ChargeADC(13),digit->ChargeADC(14),digit->ChargeADC(15),digit->ChargeADC(16),digit->ChargeADC(17),digit->ChargeADC(18),digit->ChargeADC(19),digit->ChargeADC(20),
 239                           fCalibData->GetPedestal(pmNumber),fCalibData->GetSigma(pmNumber),
 240                           fCalibData->GetPedestal(pmNumber+16),fCalibData->GetSigma(pmNumber+16),
 241                           aBBflag[pmNumber],aBGflag[pmNumber]));
 242             };
 243
 244         // Fill ESD friend object
 245         for (Int_t iEv = 0; iEv < kNClocks; iEv++) {
 246           fESDADfriend->SetPedestal(pmNumber,iEv,(Float_t)digit->ChargeADC(iEv));
 247           fESDADfriend->SetIntegratorFlag(pmNumber,iEv,(iEv%2 == 0) ? integrator : !integrator);
 248         }
 249         fESDADfriend->SetTime(pmNumber,digit->Time());
 250         fESDADfriend->SetWidth(pmNumber,digit->Width());
 251
 252     } // end of loop over digits
 253   } // end of loop over events in digits tree
 254
 255   //fESDAD->SetBit(AliESDAD::kCorrectedLeadingTime,kTRUE);
 256   //fESDAD->SetMultiplicity(mult);
 257   fESDAD->SetADC(adc);
 258   fESDAD->SetTime(time);
 259   fESDAD->SetWidth(width);
 260   fESDAD->SetBit(AliESDAD::kOnlineBitsFilled,kTRUE);
 261   fESDAD->SetBBFlag(aBBflag);
 262   fESDAD->SetBGFlag(aBGflag);
 263   //fESDAD->SetBit(AliESDAD::kCorrectedForSaturation,kTRUE);
 264
 265   /*/ now fill the AD decision and channel flags
 266   {
 267     AliADDecision offlineDecision;
 268     offlineDecision.SetRecoParam(GetRecoParam());
 269     offlineDecision.FillDecisions(fESDAD, fCalibData, fTimeSlewing);
 270   }/*/
 271
 272   if (esd) {
 273      AliDebug(1, Form("Writing AD data to ESD tree"));
 274      esd->SetADData(fESDAD);
 275
 276      AliESDfriend *fr = (AliESDfriend*)esd->FindListObject("AliESDfriend");
 277      if (fr) {
 278         AliDebug(1, Form("Writing AD friend data to ESD tree"));
 279         fr->SetADfriend(fESDADfriend);
 280     }
 281   }
 282
 283   fDigitsArray->Clear();
 284
 285 }
 286
 287 //_____________________________________________________________________________
 288 AliADCalibData* AliADReconstructor::GetCalibData() const
 289 {
 290   AliCDBManager *man = AliCDBManager::Instance();
 291
 292   AliCDBEntry *entry=0;
 293
 294   entry = man->Get("AD/Calib/Data");
 295   if(!entry){
 296     AliWarning("Load of calibration data from default storage failed!");
 297     AliWarning("Calibration data will be loaded from local storage ($ALICE_ROOT)");
 298
 299     man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
 300     man->SetRun(1);
 301     entry = man->Get("AD/Calib/Data");
 302   }
 303   // Retrieval of data in directory AD/Calib/Data:
 304
 305   AliADCalibData *calibdata = 0;
 306
 307   if (entry) calibdata = (AliADCalibData*) entry->GetObject();
 308   if (!calibdata)  AliFatal("No calibration data from calibration database !");
 309
 310   return calibdata;
 311 }
 312
 313 //_____________________________________________________________________________
 314 AliCDBStorage* AliADReconstructor::SetStorage(const char *uri)
 315 {
 316 // Sets the storage
 317
 318   Bool_t deleteManager = kFALSE;
 319
 320   AliCDBManager *manager = AliCDBManager::Instance();
 321   AliCDBStorage *defstorage = manager->GetDefaultStorage();
 322
 323   if(!defstorage || !(defstorage->Contains("AD"))){
 324      AliWarning("No default storage set or default storage doesn't contain AD!");
 325      manager->SetDefaultStorage(uri);
 326      deleteManager = kTRUE;
 327   }
 328
 329   AliCDBStorage *storage = manager->GetDefaultStorage();
 330
 331   if(deleteManager){
 332      AliCDBManager::Instance()->UnsetDefaultStorage();
 333      defstorage = 0;   // the storage is killed by AliCDBManager::Instance()->Destroy()
 334   }
 335
 336   return storage;
 337 }
 338
 339 //____________________________________________________________________________
 340 void AliADReconstructor::GetCollisionMode()
 341 {
 342   // Retrieval of collision mode
 343
 344   TString beamType = GetRunInfo()->GetBeamType();
 345   if(beamType==AliGRPObject::GetInvalidString()){
 346      AliError("AD cannot retrieve beam type");
 347      return;
 348   }
 349
 350   if( (beamType.CompareTo("P-P") ==0)  || (beamType.CompareTo("p-p") ==0) ){
 351     fCollisionMode=0;
 352   }
 353   else if( (beamType.CompareTo("Pb-Pb") ==0)  || (beamType.CompareTo("A-A") ==0) ){
 354     fCollisionMode=1;
 355   }
 356
 357   fBeamEnergy = GetRunInfo()->GetBeamEnergy();
 358   if(fBeamEnergy==AliGRPObject::GetInvalidFloat()) {
 359      AliError("Missing value for the beam energy ! Using 0");
 360      fBeamEnergy = 0.;
 361   }
 362
 363   AliDebug(1,Form("\n ++++++ Beam type and collision mode retrieved as %s %d @ %1.3f GeV ++++++\n\n",beamType.Data(), fCollisionMode, fBeamEnergy));
 364
 365 }