T0/AliT0Reconstructor.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 #include <Riostream.h>
  19
  20 #include <TDirectory.h>
  21
  22 #include "AliRunLoader.h"
  23 #include "AliRun.h"
  24 #include <AliESD.h>
  25 #include "AliLog.h"
  26 #include <TClonesArray.h>
  27 #include "AliT0RecPoint.h"
  28 #include "AliT0.h"
  29 #include "AliRawReader.h"
  30 #include "AliT0RawReader.h"
  31 #include "AliT0Loader.h"
  32 #include "AliT0digit.h"
  33 #include "AliT0Reconstructor.h"
  34 #include "AliT0Parameters.h"
  35 #include "AliCDBLocal.h"
  36 #include "AliCDBStorage.h"
  37 #include "AliCDBManager.h"
  38 #include "AliCDBEntry.h"
  39
  40 #include <TArrayI.h>
  41 #include <TGraph.h>
  42
  43 ClassImp(AliT0Reconstructor)
  44
  45 //____________________________________________________________________
  46 void
  47 AliT0Reconstructor::Init(AliRunLoader* runLoader,TTree* digitsTree ) const
  48 {
  49
  50  // Initialize the reconstructor
  51   AliDebug(2, Form("Init called with runloader 0x%x", runLoader));
  52
  53   // Initialize the parameters
  54    AliT0Loader* pStartLoader = (AliT0Loader*) runLoader->GetLoader("T0Loader");
  55
  56   pStartLoader->LoadDigits("READ");
  57
  58   digitsTree = pStartLoader->TreeD();
  59
  60 }
  61
  62 //____________________________________________________________________
  63   void  AliT0Reconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const
  64 {
  65   //T0 raw data-> digits conversion
  66  // reconstruct time information from raw data
  67     AliT0 *baseT0;
  68     baseT0->Raw2Digits(rawReader, digitsTree);
  69
  70
  71 }
  72 void AliT0Reconstructor::Reconstruct(TTree*digitsTree, TTree*clustersTree) const
  73   //void AliT0Reconstructor::Reconstruct(AliRunLoader *runLoader) const
  74 {
  75 // T0 digits reconstruction
  76 // T0RecPoint writing
  77
  78     //Q->T-> coefficients !!!! should be asked!!!
  79   //  Float_t ph2MIP=500;
  80   Float_t gain[24], timeDelayCFD[24], timeDelayLED[24];
  81   Float_t zdetA,zdetC;
  82   TObjArray slewingLED;
  83
  84   TArrayI * fADC = new TArrayI(24);
  85   TArrayI * fTimeCFD = new TArrayI(24);
  86   TArrayI * fADCLED = new TArrayI(24);
  87   TArrayI * fTimeLED = new TArrayI(24);
  88
  89   AliT0Parameters* param = AliT0Parameters::Instance();
  90   param->Init();
  91
  92   Int_t mV2Mip = param->GetmV2Mip();
  93   //mV2Mip = param->GetmV2Mip();
  94   Int_t channelWidth = param->GetChannelWidth() ;
  95
  96   for (Int_t i=0; i<24; i++){
  97     timeDelayCFD[i] = param->GetTimeDelayCFD(i);
  98     timeDelayLED[i] = param->GetTimeDelayLED(i);
  99     gain[i] = param->GetGain(i);
 100     //gain[i] = 1;
 101     slewingLED.AddAtAndExpand(param->GetSlew(i),i);
 102   }
 103   zdetC = param->GetZposition(0);
 104   zdetA  = param->GetZposition(1);
 105
 106   AliDebug(1,Form("Start DIGITS reconstruction "));
 107
 108   TBranch *brDigits=digitsTree->GetBranch("T0");
 109   AliT0digit *fDigits = new AliT0digit() ;
 110   if (brDigits) {
 111     brDigits->SetAddress(&fDigits);
 112   }else{
 113     cerr<<"EXEC Branch T0 digits not found"<<endl;
 114     return;
 115   }
 116
 117    digitsTree->GetEvent(0);
 118    digitsTree->GetEntry(0);
 119    brDigits->GetEntry(0);
 120   fDigits->GetTime(*fTimeCFD);
 121   fDigits->GetADC(*fADC);
 122   fDigits->GetTimeAmp(*fTimeLED);
 123   fDigits->GetADCAmp(*fADCLED);
 124
 125
 126   Float_t besttimeright=999999;
 127   Float_t besttimeleft=999999;
 128   Int_t pmtBestRight=99999;
 129   Int_t pmtBestLeft=99999;
 130   Float_t timeDiff=999999, meanTime=0;
 131
 132
 133
 134    AliT0RecPoint* frecpoints= new AliT0RecPoint ();
 135    clustersTree->Branch( "T0", "AliT0RecPoint" ,&frecpoints, 405,1);
 136
 137   Float_t time[24], adc[24];
 138   for (Int_t ipmt=0; ipmt<24; ipmt++) {
 139     if(fTimeCFD->At(ipmt)>0 ){
 140       time[ipmt] = channelWidth *( fTimeCFD->At(ipmt)) - 1000*timeDelayCFD[ipmt];
 141       //   Float_t adc_digPs = channelWidth * Float_t (fADCLED->At(ipmt)) ;
 142       //      adc[ipmt] = TMath::Exp(adc_digPs/1000) /gain[ipmt];
 143       adc[ipmt]=1;
 144       AliDebug(1,Form(" time %f ps,  adc %f mv in MIP %i\n ",
 145                       time[ipmt], adc[ipmt], Int_t (adc[ipmt]/mV2Mip +0.5)));
 146       frecpoints->SetTime(ipmt,time[ipmt]);
 147       frecpoints->SetAmp(ipmt,adc[ipmt]);
 148     }
 149     else {
 150       time[ipmt] = 0;
 151       adc[ipmt] = 0;
 152     }
 153   }
 154
 155   for (Int_t ipmt=0; ipmt<12; ipmt++){
 156     if(time[ipmt] > 1 ) {
 157       if(time[ipmt]<besttimeleft){
 158         besttimeleft=time[ipmt]; //timeleft
 159         pmtBestLeft=ipmt;
 160       }
 161     }
 162   }
 163   for ( Int_t ipmt=12; ipmt<24; ipmt++){
 164     if(time[ipmt] > 1) {
 165       if(time[ipmt]<besttimeright) {
 166         besttimeright=time[ipmt]; //timeright
 167         pmtBestRight=ipmt;}
 168     }
 169   }
 170   if(besttimeright !=999999)  frecpoints->SetTimeBestRight(Int_t(besttimeright));
 171   if( besttimeleft != 999999 ) frecpoints->SetTimeBestLeft(Int_t(besttimeleft));
 172   AliDebug(1,Form(" besttimeA %f ps,  besttimeC %f ps",besttimeright, besttimeleft));
 173   Float_t c = 0.0299792; // cm/ps
 174   Float_t vertex = 0;
 175   if(besttimeright !=999999 && besttimeleft != 999999 ){
 176     timeDiff = besttimeleft - besttimeright;
 177     meanTime = (besttimeright + besttimeleft)/2.;
 178     vertex = c*(timeDiff)/2.; //-(lenr-lenl))/2;
 179     AliDebug(1,Form("  timeDiff %f ps,  meanTime %f ps, vertex %f cm",timeDiff, meanTime,vertex ));
 180     frecpoints->SetVertex(vertex);
 181     frecpoints->SetMeanTime(Int_t(meanTime));
 182
 183   }
 184   clustersTree->Fill();
 185 }
 186
 187
 188 void AliT0Reconstructor::FillESD(AliRunLoader* runLoader, AliESD *pESD) const
 189 {
 190
 191   /***************************************************
 192   Resonstruct digits to vertex position
 193   ****************************************************/
 194
 195
 196   if (!runLoader) {
 197     AliError("Reconstruct >> No run loader");
 198     return;
 199   }
 200
 201   AliDebug(1,Form("Start FillESD T0"));
 202
 203   AliT0Loader* pStartLoader = (AliT0Loader*) runLoader->GetLoader("T0Loader");
 204
 205   pStartLoader->LoadRecPoints("READ");
 206
 207   TTree *treeR = pStartLoader->TreeR();
 208
 209    AliT0RecPoint* frecpoints= new AliT0RecPoint ();
 210     if (!frecpoints) {
 211     AliError("Reconstruct Fill ESD >> no recpoints found");
 212     return;
 213   }
 214
 215   AliDebug(1,Form("Start FillESD T0"));
 216   TBranch *brRec = treeR->GetBranch("T0");
 217   if (brRec) {
 218     brRec->SetAddress(&frecpoints);
 219   }else{
 220     cerr<<"EXEC Branch T0 rec not found"<<endl;
 221     // exit(111);
 222     return;
 223   }
 224
 225     brRec->GetEntry(0);
 226     Float_t timeStart, Zposition, amp[24], time[24];
 227     Int_t i;
 228     Zposition = frecpoints -> GetVertex();
 229     timeStart = frecpoints -> GetMeanTime();
 230     for ( i=0; i<24; i++) {
 231        time[i] = Float_t (frecpoints -> GetTime(i)) / 1000.; // ps to ns
 232       amp[i] = frecpoints -> GetAmp(i);
 233     }
 234     pESD->SetT0zVertex(Zposition); //vertex Z position
 235     pESD->SetT0(timeStart);        // interaction time
 236     pESD->SetT0time(time);         // best TOF on each PMT
 237     pESD->SetT0amplitude(amp);     // number of particles(MIPs) on each PMT
 238     cout<<" ESD >> "<<Zposition<<" "<<timeStart<<endl;
 239
 240     pStartLoader->UnloadRecPoints();
 241
 242 } // vertex in 3 sigma
 243
 244
 245
 246
 247
 248