FIT/AliFITReconstructor.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 /*********************************************************************
  17  *  FIT reconstruction and filling ESD
  18  *  - reconstruct mean time (interation time)
  19  *  - vertex position
  20  *  -  multiplicity
  21  ********************************************************************/
  22
  23 #include <AliESDEvent.h>
  24 #include "AliLog.h"
  25 #include "AliFITRecPoint.h"
  26 #include "AliFITDigit.h"
  27 #include "AliFITReconstructor.h"
  28 #include "AliESDFIT.h"
  29 #include "AliLog.h"
  30 #include "AliESDRun.h"
  31 #include "AliFITRawReader.h"
  32
  33 #include <TArrayI.h>
  34 #include <TGraph.h>
  35 #include <TMath.h>
  36 #include <Riostream.h>
  37
  38 using std::cout;
  39 using std::endl;
  40
  41 ClassImp(AliFITReconstructor)
  42
  43 AliFITReconstructor:: AliFITReconstructor(): AliReconstructor(),
  44                                              fESD(NULL),
  45                                              fESDFIT(NULL),
  46                                              fDigits(NULL)
  47
  48 {
  49
  50   //constructor
  51
  52
  53    fESDFIT  = new AliESDFIT();
  54
  55
  56 }
  57 //_____________________________________________________________________________
  58 void AliFITReconstructor::Init()
  59 {
  60 // initializer
  61
  62   fESDFIT  = new AliESDFIT;
  63  }
  64
  65 //______________________________________________________________________
  66 void AliFITReconstructor::ConvertDigits(AliRawReader* rawReader, TTree* digitsTree) const
  67 {
  68 // converts RAW to digits
  69   cout<<" AliFITReconstructor::ConvertDigits "<<rawReader<<" "<<digitsTree<<endl;
  70   if (!digitsTree) {
  71     AliError("No digits tree!");
  72     return;
  73   }
  74
  75   fDigits = new TClonesArray ("AliFITDigit", 100);
  76   digitsTree->Branch("FIT", &fDigits);
  77
  78   AliFITRawReader myrawreader(rawReader);
  79   if (myrawreader.Next())
  80     {
  81       Int_t allData[500];
  82       for (Int_t i=0; i<500; i++)  allData[i]=0;
  83       for (Int_t i=0; i<500; i++)
  84         if(myrawreader.GetData(i)>0)  {
  85           allData[i]=myrawreader.GetData(i);
  86         }
  87
  88       Int_t timeCFD, timeLED, timeQT1, timeQT0;
  89       for (Int_t ipmt=0; ipmt<160; ipmt++) {
  90         if(allData[ipmt]>0) {
  91           timeCFD = allData[ipmt];
  92           timeLED = allData[ipmt];
  93           timeQT0= allData[ipmt+160];
  94           timeQT1 = allData[ipmt+320];
  95           //add digit
  96           new((*fDigits)[fDigits->GetEntriesFast()] )
  97             AliFITDigit( ipmt,timeCFD, timeLED, timeQT0, timeQT1, 0);
  98         }
  99       }
 100
 101
 102
 103       digitsTree->Fill();
 104     }
 105
 106
 107 }
 108
 109  //____________________________________________________________
 110 void AliFITReconstructor::FillESD(TTree *digitsTree, TTree * /*clustersTree*/, AliESDEvent *pESD) const
 111 {
 112
 113   /***************************************************
 114   Resonstruct digits to vertex position
 115   ****************************************************/
 116
 117   AliDebug(1,Form("Start FillESD FIT"));
 118   if(!pESD) {
 119     AliError("No ESD Event");
 120     return;
 121   }
 122
 123   Float_t channelWidth = 24.4;
 124   Float_t c = 0.0299792458; // cm/ps
 125   Float_t currentVertex, shift=0;
 126   Int_t ncont=-1;
 127    const AliESDVertex* vertex = pESD->GetPrimaryVertex();
 128   if (!vertex)        vertex = pESD->GetPrimaryVertexSPD();
 129   if (!vertex)        vertex = pESD->GetPrimaryVertexTPC();
 130   if (!vertex)        vertex = pESD->GetVertex();
 131
 132   if (vertex) {
 133     AliDebug(2, Form("Got %s (%s) from ESD: %f",
 134                      vertex->GetName(), vertex->GetTitle(), vertex->GetZ()));
 135     currentVertex = vertex->GetZ();
 136
 137     ncont = vertex->GetNContributors();
 138     if(ncont>0 ) {
 139       shift = currentVertex/c;
 140     }
 141   } //vertex
 142
 143
 144   // FIT digits reconstruction
 145
 146   if (!digitsTree) {
 147     AliError("No digits tree!");
 148     return;
 149   }
 150
 151   fDigits = new TClonesArray("AliFITDigit", 100);
 152   // digitsTree->Print();
 153   TBranch *brDigits=digitsTree->GetBranch("FIT");
 154   if (brDigits)
 155     brDigits->SetAddress(&fDigits);
 156   else {
 157     AliError("no FIT branch");
 158     return;
 159   }
 160   const Float_t max_time = 1e7;
 161   Int_t pmt=-1;
 162   Float_t  time[160],amp[160];
 163   for(Int_t i=0; i<160; i++)   {  time[i]=max_time; amp[i]=0;}
 164
 165   Int_t nEntries = (Int_t)digitsTree->GetEntries();
 166   for (Int_t iev=0; iev<nEntries; iev++) {
 167     digitsTree->GetEvent(iev,1);
 168     brDigits->GetEvent(iev);
 169     Int_t nDigits = fDigits->GetEntriesFast();
 170     for (Int_t dig=0; dig<nDigits; dig++) {
 171       AliFITDigit* digit = (AliFITDigit*) fDigits->At(dig);
 172       pmt = digit->NPMT();
 173       time[pmt] = Float_t (digit->TimeCFD() );
 174       amp[pmt] = Float_t (digit->TimeQT1() - digit->TimeQT0() );
 175     }
 176     fESDFIT->SetFITtime(time);         // best TOF on each PMT
 177     fESDFIT->SetFITamplitude(amp);     // number of particles(MIPs) on each
 178
 179     Float_t firsttime[3] = {max_time,max_time,max_time};
 180
 181     Float_t vertexFIT = 9999999;
 182     for (Int_t ipmt=0; ipmt<80; ipmt++)//timeC
 183           if(time[ipmt]<firsttime[2]) firsttime[2]=time[ipmt];
 184
 185     for ( Int_t ipmt=80; ipmt<160; ipmt++)
 186           if(time[ipmt]<firsttime[1]) firsttime[1]=time[ipmt];
 187
 188
 189     AliDebug(5,Form("1stimeA %f , 1sttimeC %f  ",
 190                     firsttime[1],
 191                     firsttime[2] ) );
 192     if (firsttime[1]<max_time && firsttime[2]<max_time)  {
 193         firsttime[0] =  channelWidth *(firsttime[1] + firsttime[2])/2;
 194         vertexFIT =  c*channelWidth*(firsttime[1] - firsttime[2])/2;
 195     }
 196       if (firsttime[1]<max_time)
 197         firsttime[1] = firsttime[1] * channelWidth + shift;
 198
 199       if (firsttime[2]<max_time)
 200         firsttime[2] = firsttime[2] * channelWidth - shift;
 201
 202
 203
 204       for(Int_t i=0; i<3; i++) {
 205         fESDFIT->SetFITT0(i,firsttime[i]);   // interaction time (ns)
 206       }
 207       fESDFIT->SetFITzVertex(vertexFIT); //vertex Z position
 208
 209
 210       AliDebug(1,Form("FIT: SPDshift %f Vertex %f  FITsignal %f ps FITA %f ps FITC %f ps \n",shift, vertexFIT, firsttime[0], firsttime[1],firsttime[2]));
 211   }
 212
 213    }