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