START/AliSTARTDigitizer.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-2000, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17
  18 #include <TTree.h>
  19 #include <TFile.h>
  20 #include <TDirectory.h>
  21 #include <TRandom.h>
  22 #include <TArrayI.h>
  23 #include <TError.h>
  24
  25
  26 #include "AliSTARTDigitizer.h"
  27 #include "AliSTART.h"
  28 #include "AliSTARThit.h"
  29 #include "AliSTARThitPhoton.h"
  30 #include "AliSTARTdigit.h"
  31 #include "AliRunDigitizer.h"
  32 #include <AliDetector.h>
  33 #include "AliRun.h"
  34 #include <AliLoader.h>
  35 #include <AliRunLoader.h>
  36 #include <stdlib.h>
  37 #include <Riostream.h>
  38 #include <Riostream.h>
  39
  40 ClassImp(AliSTARTDigitizer)
  41
  42 //___________________________________________
  43   AliSTARTDigitizer::AliSTARTDigitizer()  :AliDigitizer()
  44 {
  45 // Default ctor - don't use it
  46   ;
  47 }
  48
  49 //___________________________________________
  50 AliSTARTDigitizer::AliSTARTDigitizer(AliRunDigitizer* manager)
  51     :AliDigitizer(manager)
  52 {
  53   //    cout<<"AliSTARTDigitizer::AliSTARTDigitizer"<<endl;
  54 // ctor which should be used
  55 //  fDebug =0;
  56   if (GetDebug())
  57     Info("(AliRunDigitizer* manager)" ,"processed");
  58
  59   ftimeRightTDC = new TArrayI(12);
  60   ftimeLeftTDC = new TArrayI(12);
  61   fRightADC = new TArrayI(12);
  62   fLeftADC = new TArrayI(12);
  63 }
  64
  65 //------------------------------------------------------------------------
  66 AliSTARTDigitizer::~AliSTARTDigitizer()
  67 {
  68 // Destructor
  69   if(GetDebug()) Info("dtor","START");
  70   delete ftimeRightTDC;
  71   delete ftimeLeftTDC;
  72   delete fRightADC;
  73   delete fLeftADC;
  74 }
  75
  76  //------------------------------------------------------------------------
  77 Bool_t AliSTARTDigitizer::Init()
  78 {
  79 // Initialization
  80 // cout<<"AliSTARTDigitizer::Init"<<endl;
  81  return kTRUE;
  82 }
  83
  84
  85 //---------------------------------------------------------------------
  86
  87 void AliSTARTDigitizer::Exec(Option_t* /*option*/)
  88 {
  89
  90   /*
  91     Produde digits from hits
  92         digits is TObject and includes
  93         We are writing array if left & right  TDC
  94         left & right  ADC (will need for slow simulation)
  95         TOF first particle left & right
  96         mean time and time difference (vertex position)
  97
  98   */
  99
 100   AliRunLoader *inRL, *outRL;//in and out Run Loaders
 101   AliLoader *pInStartLoader, *pOutStartLoader;// in and out STARTLoaders
 102
 103   outRL = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 104   pOutStartLoader = outRL->GetLoader("STARTLoader");
 105
 106 #ifdef DEBUG
 107   cout<<"AliSTARTDigitizer::>Hits2Digits start...\n";
 108 #endif
 109   //
 110   // From hits to digits
 111   //
 112   Int_t hit, nhits;
 113   Float_t meanTime;
 114   Int_t countEr[13],countEl[13];
 115   Int_t volume,pmt,tr,tl,sumRight;
 116   Int_t  bestRightADC,bestLeftADC;
 117   Float_t besttimeleftGaus, besttimerightGaus;
 118   Float_t timeright[13]={13*0};
 119   Float_t timeleft[13]={13*0};
 120   Float_t channelWidth=2.5; //ps
 121   Int_t channelWidthADC=1; //ps
 122   //  Int_t thresholdAmpl=10;
 123
 124
 125   AliSTARThit  *startHit;
 126   TBranch *brHits=0;
 127   TBranch *brHitPhoton=0;
 128   pOutStartLoader->LoadDigits("UPDATE");//probably it is necessary to load them before
 129   fdigits= new AliSTARTdigit();
 130   pOutStartLoader->GetDigitsDataLoader()->GetBaseLoader(0)->Post(fdigits);
 131
 132   Int_t nFiles=fManager->GetNinputs();
 133   for (Int_t inputFile=0; inputFile<nFiles;  inputFile++) {
 134     if (inputFile < nFiles-1) {
 135       Warning("Exec", "ignoring input stream %d", inputFile);
 136       continue;
 137     }
 138
 139     Float_t besttimeright=9999.;
 140     Float_t besttimeleft=9999.;
 141     Int_t iTimeDiff=0;
 142     Int_t iTimeAv=0;
 143     Float_t timeDiff,timeAv;
 144     sumRight=0;
 145     for (Int_t i0=0; i0<13; i0++)
 146       {
 147         timeright[i0]=0; timeleft[i0]=0;
 148         countEr[i0]=0;   countEl[i0]=0;
 149       }
 150
 151     inRL = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(inputFile));
 152     if (!inRL->GetAliRun()) inRL->LoadgAlice();
 153     AliSTART *fSTART  = (AliSTART*) inRL->GetAliRun()->GetDetector("START");
 154     pInStartLoader = inRL->GetLoader("STARTLoader");
 155     pInStartLoader->LoadHits("READ");//probably it is necessary to load them before
 156     TClonesArray *fHits = fSTART->Hits ();
 157
 158     TTree *th = pInStartLoader->TreeH();
 159     brHits = th->GetBranch("START");
 160     brHitPhoton = th->GetBranch("STARThitPhoton");
 161     if (brHits) {
 162       fSTART->SetHitsAddressBranch(brHits,brHitPhoton);
 163     }else{
 164       cerr<<"EXEC Branch START hit not found"<<endl;
 165       exit(111);
 166     }
 167     Int_t ntracks    = (Int_t) th->GetEntries();
 168 #ifdef DEBUG
 169    Info("Digitizer",ntracks);
 170 #endif
 171      if (ntracks<=0) return;
 172     // Start loop on tracks in the hits containers
 173     for (Int_t track=0; track<ntracks;track++) {
 174       brHits->GetEntry(track);
 175       nhits = fHits->GetEntriesFast();
 176       for (hit=0;hit<nhits;hit++) {
 177         startHit   = (AliSTARThit*) fHits->UncheckedAt(hit);
 178         if (!startHit) {
 179           ::Error("Exec","The unchecked hit doesn't exist");
 180           break;
 181         }
 182         pmt=startHit->Pmt();
 183         volume = startHit->Volume();
 184         if(volume==1){
 185           timeright[pmt] = startHit->Time();
 186           if(timeright[pmt]<besttimeright)
 187             {
 188               besttimeright=timeright[pmt];
 189           } //timeright
 190         }//time for right shoulder
 191         if(volume==2){
 192           timeleft[pmt] = startHit->Time();
 193           if(timeleft[pmt]<besttimeleft)
 194             {
 195               besttimeleft=timeleft[pmt];
 196
 197           } //timeleftbest
 198         }//time for left shoulder
 199       } //hit loop
 200     } //track loop
 201
 202     // z position
 203
 204     //folding with experimental time distribution
 205
 206     Float_t koef=69.7/350.;
 207     besttimeright=koef*besttimeright;
 208     besttimeleftGaus=gRandom->Gaus(besttimeleft,0.05);
 209     bestLeftADC=Int_t (besttimeleftGaus*1000/channelWidth);
 210     besttimerightGaus=gRandom->Gaus(besttimeright,0.05);
 211     bestRightADC=Int_t (besttimerightGaus*1000/channelWidth);
 212     timeDiff=besttimerightGaus-besttimeleftGaus;
 213 #ifdef DEBUG
 214     cout<<" timediff in ns "<<timeDiff<<" z= "<<timeDiff*30<<endl;
 215 #endif
 216     meanTime=(besttimerightGaus+besttimeleftGaus)/2.;
 217     if ( TMath::Abs(timeDiff)<TMath::Abs(0.3) )
 218       {
 219         Float_t t1=1000.*besttimeleftGaus;
 220         Float_t t2=1000.*besttimerightGaus;
 221         t1=t1/channelWidth;   //time in ps to channelWidth
 222         t2=t2/channelWidth;   //time in ps to channelWidth
 223         timeAv=(t1+t2)/2.;// time  channel numbres
 224
 225         // Time to TDC signal
 226         // 256 channels for timediff, range 1ns
 227         iTimeAv=(Int_t)timeAv;
 228         timeDiff= 512+1000*timeDiff/channelWidth; // time  channel numbres
 229         iTimeDiff=(Int_t)timeDiff;
 230         //       fill digits
 231         fdigits->SetTimeBestLeft(bestLeftADC);
 232         fdigits->SetTimeBestRight(bestRightADC);
 233         fdigits->SetMeanTime(iTimeAv);
 234         fdigits->SetTimeDiff(iTimeDiff);
 235         for (Int_t i=0; i<12; i++)
 236           {
 237             //  fill TDC
 238             timeright[i+1]=gRandom->Gaus(timeright[i+1],0.05);
 239             timeleft[i+1]=gRandom->Gaus(timeleft[i+1],0.05);
 240             tr= Int_t (timeright[i+1]*1000/channelWidth);
 241             if(tr<200) tr=0;
 242             tl= Int_t (timeleft[i+1]*1000/channelWidth);
 243             if(tl<1000) tl=0;
 244
 245             ftimeRightTDC->AddAt(tr,i);
 246             ftimeLeftTDC->AddAt(tl,i);
 247             //fill ADC
 248             Int_t al=( Int_t ) countEl[i+1]/ channelWidthADC;
 249             Int_t ar=( Int_t ) countEr[i+1]/ channelWidthADC;
 250             fRightADC->AddAt(ar,i);
 251             fLeftADC ->AddAt(al,i);
 252             sumRight+=countEr[i+1];
 253           }
 254         fdigits->SetTimeRight(*ftimeRightTDC);
 255         fdigits->SetTimeLeft(*ftimeLeftTDC);
 256         fdigits->SetADCRight(*fRightADC);
 257         fdigits->SetADCLeft(*fLeftADC);
 258         fdigits->SetSumADCRight(sumRight);
 259       }
 260     else
 261       {timeAv=999999; timeDiff=99999;}
 262
 263     pInStartLoader->UnloadHits();
 264   } //input streams loop
 265
 266   pOutStartLoader->WriteDigits("OVERWRITE");
 267   pOutStartLoader->UnloadDigits();
 268 }
 269
 270
 271 //------------------------------------------------------------------------
 272 Bool_t AliSTARTDigitizer::RegisterPhotoE(/*AliSTARThitPhoton *hit*/)
 273 {
 274     Double_t    pP = 0.2;
 275     Double_t    p;
 276
 277     p = gRandom->Rndm();
 278     if (p > pP)
 279       return kFALSE;
 280
 281     return kTRUE;
 282 }
 283 //----------------------------------------------------------------------------