Updated version of TOF digitization, N^2 problem solved (J.Chudoba)
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 22 Nov 2001 11:22:51 +0000 (11:22 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 22 Nov 2001 11:22:51 +0000 (11:22 +0000)
19 files changed:
TOF/AliTOF.cxx
TOF/AliTOF.h
TOF/AliTOFConstants.cxx [new file with mode: 0644]
TOF/AliTOFConstants.h [new file with mode: 0644]
TOF/AliTOFHitMap.cxx [new file with mode: 0644]
TOF/AliTOFHitMap.h [new file with mode: 0644]
TOF/AliTOFSDigit.cxx [new file with mode: 0644]
TOF/AliTOFSDigit.h [new file with mode: 0644]
TOF/AliTOFSDigitizer.cxx
TOF/AliTOFSDigitizer.h
TOF/AliTOFv0.cxx
TOF/AliTOFv1.cxx
TOF/AliTOFv2.cxx
TOF/AliTOFv2FHoles.cxx
TOF/AliTOFv3.cxx
TOF/AliTOFv4.cxx
TOF/Makefile
TOF/TOFLinkDef.h
TOF/libTOF.pkg

index c76333d..d6808dd 100644 (file)
@@ -15,6 +15,9 @@
 
 /*
 $Log$
+Revision 1.30  2001/10/21 18:30:39  hristov
+Several pointers were set to zero in the default constructors to avoid memory management problems
+
 Revision 1.29  2001/10/17 14:19:24  hristov
 delete replaced by delete []
 
@@ -117,6 +120,7 @@ Introduction of the Copyright and cvs Log
 #include "AliTOF.h"
 #include "AliTOFhit.h"
 #include "AliTOFdigit.h"
+#include "AliTOFSDigit.h"
 #include "AliTOFRawSector.h"
 #include "AliTOFRoc.h"
 #include "AliTOFRawDigit.h"
@@ -174,7 +178,7 @@ AliTOF::AliTOF(const char *name, const char *title)
   fHits   = new TClonesArray("AliTOFhit",  1000);
   gAlice->AddHitList(fHits);
   fIshunt  = 0;
-  fSDigits       = new TClonesArray("AliTOFdigit",  1000);
+  fSDigits       = new TClonesArray("AliTOFSDigit",  1000);
   fDigits        = new TClonesArray("AliTOFdigit",  1000);
   //
   // Digitization parameters
@@ -337,7 +341,7 @@ void AliTOF::AddDigit(Int_t *tracks, Int_t *vol, Float_t *digits)
 }
 
 //___________________________________________
-void AliTOF::AddSDigit(Int_t *tracks, Int_t *vol, Float_t *digits)
+void AliTOF::AddSDigit(Int_t tracknum, Int_t *vol, Float_t *digits)
 {
      
 //
@@ -345,7 +349,7 @@ void AliTOF::AddSDigit(Int_t *tracks, Int_t *vol, Float_t *digits)
 //
         
   TClonesArray &lSDigits = *fSDigits;   
-  new(lSDigits[fNSDigits++]) AliTOFdigit(tracks, vol, digits);
+  new(lSDigits[fNSDigits++]) AliTOFSDigit(tracknum, vol, digits);
 }
 
 //_____________________________________________________________________________
@@ -369,7 +373,7 @@ void AliTOF::SetTreeAddress ()
        }
 
     }
-  if (fSDigits)
+//  if (fSDigits)
     //  fSDigits->Clear ();
 
   if (gAlice->TreeS () && fSDigits)
@@ -730,118 +734,6 @@ AliTOFMerger*  AliTOF::Merger()
     return fMerger;
 }
 
-//____________________________________________________________________________
-void AliTOF::TOFHits2SDigits()
-{
-//
-// Starting from the Hits Tree (TreeH), this
-// function writes the TOF SDigits Branch in the TreeS storing 
-// the digits informations.
-// It has to be called just at the end of an event or 
-// at the end of a whole run.
-// It could  also be called by AliTOF::Finish Event()
-// Just for MC events. 
-//
-// Called by the ROOT script Hits2SDigits.C
-//
-// Simulation of detector response.
-
-   Int_t ver = this->IsVersion();
-   if(ver==0) return; // no sdigits for AliTOFv0
-
-  AliTOF *TOF = (AliTOF *) gAlice->GetDetector ("TOF");
-
-  if (fNevents == 0)
-    fNevents = (Int_t) gAlice->TreeE ()->GetEntries ();
-
-       cout << "nevents found on file " << fNevents << endl;
-      // Start Event ------------------------- LOOP
-
-  for (Int_t ievent = 0; ievent < fNevents; ievent++)
-    {
-      gAlice->GetEvent (ievent);
-      if (gAlice->TreeH () == 0)
-       return; // no hits stored 
-      if (gAlice->TreeS () == 0)
-       gAlice->MakeTree ("S");
-
-
-      Int_t nSdigits = 0;
-      
-            //Make branches
-      char branchname[20];
-       sprintf (branchname, "%s", TOF->GetName ());
-      //Make branch for TOF sdigits
-      TOF->MakeBranch ("S");
-
-      Int_t    tracks[3];    // track info
-      Int_t    vol[5];       // location for a digit
-      Float_t  digit[2];     // TOF digit variables
-      Int_t hit, nbytes;
-      TParticle *particle;
-      AliTOFhit *tofHit;
-      TClonesArray *TOFhits = TOF->Hits();
-
-
-      if (TOF)
-       {
-         TOFhits = TOF->Hits();             // pointer to the TClonesArray of TOF hits
-         TTree *TH = gAlice->TreeH();       // pointer to the current TreeH
-         Stat_t ntracks = TH->GetEntries(); // number of tracks for the current event
-         for (Int_t track = 0; track < ntracks; track++)
-           {
-             gAlice->ResetHits ();
-             nbytes += TH->GetEvent(track);
-             particle = gAlice->Particle(track);
-             Int_t nhits = TOFhits->GetEntriesFast(); // number of hits for the current track
-
-             for (hit = 0; hit < nhits; hit++)
-               {
-               tofHit = (AliTOFhit*) TOFhits->UncheckedAt(hit);
-               vol[0] = tofHit->GetSector();
-               vol[1] = tofHit->GetPlate();
-               vol[2] = tofHit->GetStrip();
-               vol[3] = tofHit->GetPadx();
-               vol[4] = tofHit->GetPadz();
-
-               // 95% of efficiency to be inserted here
-               // edge effect to be inserted here
-               // cross talk  to be inserted here
-
-               Float_t idealtime = tofHit->GetTof(); // unit s
-               idealtime *= 1.E+12;  // conversion from s to ps
-                              // fTimeRes is given usually in ps
-               Float_t tdctime   = gRandom->Gaus(idealtime, fTimeRes); 
-               digit[0] = tdctime;
-
-               // typical Landau Distribution to be inserted here
-               // instead of Gaussian Distribution
-               Float_t idealcharge = tofHit->GetEdep();
-               Float_t adccharge = gRandom->Gaus(idealcharge, fChrgRes);
-               digit[1] = adccharge;
-               Int_t tracknum = tofHit->GetTrack();
-               tracks[0] = tracknum;
-               tracks[1] = 0;
-               tracks[2] = 0;
-
-               // check if two digit are on the same pad; in that case we sum
-               // the two or more digits
-               Bool_t overlap = CheckOverlap(vol, digit, tracknum);
-               if(!overlap) 
-                       new ((*fSDigits)[nSdigits++]) AliTOFdigit(tracks, vol, digit);
-               //cout << "nSdigits" << endl; 
-               } // end loop on hits for the current track
-
-            } // end loop on ntracks
-
-         } // close if TOF switched ON
-
-      gAlice->TreeS()->Fill();
-      gAlice->TreeS()->Print();
-
-    }                          //event loop
-
-}
 
 //---------------------------------------------------------------------
 
index 1228d8e..1a12791 100644 (file)
@@ -50,7 +50,7 @@ public:
   virtual void    SetTreeAddress();
   virtual void    AddHit(Int_t track, Int_t* vol, Float_t* hits);
   virtual void    AddDigit(Int_t* tracks, Int_t* vol, Float_t* digits);
-  virtual void    AddSDigit(Int_t* tracks, Int_t* vol, Float_t* digits);
+  virtual void    AddSDigit(Int_t tracknum, Int_t* vol, Float_t* digits);
   virtual void    CreateGeometry();
   virtual void    CreateMaterials();
   virtual void    Init();
@@ -74,7 +74,6 @@ public:
   virtual void    SetMerger(AliTOFMerger* merger);
   virtual AliTOFMerger* Merger();
 //  virtual void    Hits2Digits();   
-  virtual void    TOFHits2SDigits();
   virtual void    Hits2SDigits();
   virtual void    Digits2Reco() {cout << "AliTOF::Digits2Reco()  dummy function called" << endl;}
           void    Digits2Raw (Int_t evNumber=0);
diff --git a/TOF/AliTOFConstants.cxx b/TOF/AliTOFConstants.cxx
new file mode 100644 (file)
index 0000000..1585162
--- /dev/null
@@ -0,0 +1,46 @@
+ /**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*$Log$*/
+
+////////////////////////////////////////////////////////////////////////
+//
+// AliTOFConstants class
+//
+// This class serves to group constants needed by TOF detector in 1
+// easily accessible place. All constants are public const static data 
+// members. The class is never instatiated.
+//
+// Note: only a few constants are in the first version of this class,
+//       more should be added by TOF developpers
+//
+// Author: Jiri Chudoba (CERN)
+//
+////////////////////////////////////////////////////////////////////////
+
+#include "AliTOFConstants.h"
+
+const Int_t AliTOFConstants::fgkNStripA;
+const Int_t AliTOFConstants::fgkNStripB;
+const Int_t AliTOFConstants::fgkNStripC;
+const Int_t AliTOFConstants::fgkNpadX;
+const Int_t AliTOFConstants::fgkNpadZ;
+const Int_t AliTOFConstants::fgkPadXSector;
+const Int_t AliTOFConstants::fgkNSectors;
+const Int_t AliTOFConstants::fgkNPlates;
+
+const Int_t AliTOFConstants::fgkTimeDiff;
+
+ClassImp(AliTOFConstants)
diff --git a/TOF/AliTOFConstants.h b/TOF/AliTOFConstants.h
new file mode 100644 (file)
index 0000000..d10b4a9
--- /dev/null
@@ -0,0 +1,48 @@
+#ifndef ALITOFCONSTANTS_H
+#define ALITOFCONSTANTS_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/*$Id$*/
+
+////////////////////////////////////////////////////////////////////////
+//
+// AliTOFConstants class
+//
+// This class serves to group constants needed by TOF detector in 1
+// easily accessible place. All constants are public const static data 
+// members. The class is never instatiated.
+//
+//
+// Author: Jiri Chudoba (CERN)
+//
+////////////////////////////////////////////////////////////////////////
+
+#include <TObject.h>
+
+class AliTOFConstants {
+ public:
+    // return number of chambers
+    static const Int_t fgkNStripA = 15; // number of strips in A type module 
+    static const Int_t fgkNStripB = 19; // number of strips in B type module 
+    static const Int_t fgkNStripC = 20; // number of strips in C type module 
+    static const Int_t fgkNpadX   = 48; // Number of pads in a strip along the X direction
+    static const Int_t fgkNpadZ   =  2; // Number of pads in a strip along the Z direction
+    static const Int_t fgkPadXSector =
+      (fgkNStripA + 2*fgkNStripB + 2*fgkNStripC)*fgkNpadX*fgkNpadZ;
+    static const Int_t fgkNSectors   =  18;
+    static const Int_t fgkNPlates    =  5;
+
+// if two signals ar eseparated less than fgkTimeDiff, they are merged
+// and considered as one     
+    static const Int_t fgkTimeDiff   =  25000; // time in ps, 
+    
+    
+ private:
+    AliTOFConstants(){}
+    virtual ~AliTOFConstants(){}
+
+    ClassDef(AliTOFConstants, 0)             // TOF global constants 
+};
+       
+#endif
diff --git a/TOF/AliTOFHitMap.cxx b/TOF/AliTOFHitMap.cxx
new file mode 100644 (file)
index 0000000..26d7e90
--- /dev/null
@@ -0,0 +1,196 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Log$ */
+
+////////////////////////////////////////////////////////////////////////
+//
+// AliTOFHitMap class
+//
+// hitmap enables fast check if the pad was already hit
+// The index of a AliTOFSDigit is saved in the each hitmap "cell"
+// (there is an offset +1, because the index can be zero and 
+// zero means empty cell. 
+// In TOF, number of strips varies according plate type, the highest
+// number is in plate C. For all plates is used this number, so
+// the size of the hitmap is a little bit greater than necessary, but
+// it simplifies the access algorithm. 
+// 
+//
+// Author: Jiri Chudoba (CERN), based on AliMUONHitMap
+//
+////////////////////////////////////////////////////////////////////////
+
+#include <iostream.h>
+#include <TMath.h>
+
+#include "AliTOFHitMap.h"
+#include "AliTOFSDigit.h"
+#include "AliTOFConstants.h"
+
+
+#include <TClonesArray.h>
+
+ClassImp(AliTOFHitMap)
+
+AliTOFHitMap::AliTOFHitMap()
+{
+//
+// Default ctor
+//
+  fHitMap = 0;
+  fSDigits = 0;
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFHitMap::AliTOFHitMap(TClonesArray *dig)
+{
+//
+// ctor
+//
+
+// of course, these constants must not be hardwired
+// change later
+
+  fNSector = AliTOFConstants::fgkNSectors;
+  fNplate = AliTOFConstants::fgkNPlates;
+  fNstrip = AliTOFConstants::fgkNStripC;
+  fNpx  = AliTOFConstants::fgkNpadX;
+  fNpy  = AliTOFConstants::fgkNpadZ;
+  fMaxIndex=fNSector*fNplate*fNstrip*fNpx*fNpy;
+  fHitMap = new Int_t[fMaxIndex];
+  fSDigits =  dig;
+  Clear();
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFHitMap::AliTOFHitMap(const AliTOFHitMap & hitMap)
+{
+//
+// Dummy copy constructor
+//
+    ;
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFHitMap::~AliTOFHitMap()
+{
+//
+// Destructor
+//
+    if (fHitMap) delete[] fHitMap;
+}
+
+////////////////////////////////////////////////////////////////////////
+void AliTOFHitMap::Clear(const char *)
+{
+//
+// Clear hitmap
+//
+    memset(fHitMap,0,sizeof(int)*fMaxIndex);
+}
+
+////////////////////////////////////////////////////////////////////////
+Int_t AliTOFHitMap::CheckedIndex(Int_t *vol) const
+{
+//
+// Return checked indices for vol
+//
+  Int_t index=
+    (vol[0]-1)*fNplate*fNstrip*fNpx*fNpy+             // sector
+    (vol[1]-1)*fNstrip*fNpx*fNpy+                     // plate
+    (vol[2]-1)*fNpx*fNpy+                             // strip
+    (vol[3]-1)*fNpy+                                  // padx
+    (vol[4]-1);                                        // pady (=padz)
+
+    if (index >= fMaxIndex) {
+      Error("AliTOFHitMap","CheckedIndex - input outside bounds");
+       return -1;
+    } else {
+       return index;
+    }
+}
+
+////////////////////////////////////////////////////////////////////////
+void  AliTOFHitMap::SetHit(Int_t *vol, Int_t idigit)
+{
+//
+// Assign digit to pad vol
+//
+
+// 0 means empty pad, we need to shift indeces by 1
+    fHitMap[CheckedIndex(vol)]=idigit+1;
+}
+
+////////////////////////////////////////////////////////////////////////
+void  AliTOFHitMap::SetHit(Int_t *vol)
+{
+//
+// Assign last digit to pad vol 
+//
+
+// 0 means empty pad, we need to shift indeces by 1
+    fHitMap[CheckedIndex(vol)]=fSDigits->GetLast()+1;
+}
+
+////////////////////////////////////////////////////////////////////////
+Int_t AliTOFHitMap::GetHitIndex(Int_t *vol) const
+{
+//
+// Get contents of pad vol
+//
+
+// 0 means empty pad, we need to shift indeces by 1
+    return fHitMap[CheckedIndex(vol)]-1;
+}
+
+////////////////////////////////////////////////////////////////////////
+TObject* AliTOFHitMap::GetHit(Int_t *vol) const
+{
+//
+// Get pointer to object at vol
+// return 0 if vol out of bounds
+    Int_t index=GetHitIndex(vol)-1;
+    return (index <0) ? 0 : fSDigits->UncheckedAt(index);
+}
+
+////////////////////////////////////////////////////////////////////////
+FlagType AliTOFHitMap::TestHit(Int_t *vol) const
+{
+//
+// Check if hit cell is empty, used or unused
+//
+    Int_t inf=fHitMap[CheckedIndex(vol)];
+    if (inf < 0) {
+       return kUsed;
+    } else if (inf == 0) {
+       return kEmpty;
+    } else {
+       return kUnused;
+    }
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFHitMap & AliTOFHitMap::operator = (const AliTOFHitMap & rhs) 
+{
+// Dummy assignment operator
+    return *this;
+}
+
+
+
+
+
diff --git a/TOF/AliTOFHitMap.h b/TOF/AliTOFHitMap.h
new file mode 100644 (file)
index 0000000..64222ff
--- /dev/null
@@ -0,0 +1,61 @@
+#ifndef ALITOFHITMAP_H
+#define ALITOFHITMAP_H
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+////////////////////////////////////////////////////////////////////////
+//
+// AliTOFHitMap class
+//
+// hitmap enables fast check if the pad was already hit
+//
+// Author: Jiri Chudoba (CERN)
+//
+////////////////////////////////////////////////////////////////////////
+
+#include "AliHitMap.h"
+#include "TObject.h"
+class TClonesArray;
+
+class AliTOFHitMap : public TObject
+{
+ public:
+    AliTOFHitMap();
+    AliTOFHitMap(TClonesArray *sdig);
+    AliTOFHitMap(const AliTOFHitMap & hitMap);
+    
+    virtual ~AliTOFHitMap();
+    // Clear the hit map
+    virtual  void  Clear(const char *opt = "");
+    // Set a single hit
+    virtual  void  SetHit(Int_t *vol, Int_t idigit);
+    virtual  void  SetHit(Int_t *vol);
+    // Get index of hit in the list of digits
+    virtual Int_t  GetHitIndex(Int_t *vol) const;
+    // Get pointer to digit
+    virtual TObject*  GetHit(Int_t *vol) const;
+    // Test hit status
+    virtual FlagType TestHit(Int_t *vol) const;
+    // Assignment operator
+    AliTOFHitMap& operator = (const AliTOFHitMap& rhs);
+    
+ private:
+    // Check index
+    Int_t CheckedIndex(Int_t *vol) const;
+ private:
+    Int_t fNSector;                       // Number of sectors
+    Int_t fNplate;                        // Number of plates
+    Int_t fNstrip;                        // Maximum number of strips
+    Int_t fNpx;                           // Number of pads in x
+    Int_t fNpy;                           // Number of pads in y
+
+    TClonesArray *fSDigits;               // Pointer to sdigits
+    Int_t fMaxIndex;                      // maximum index in hit map
+    Int_t *fHitMap;                       // ! [fMaxIndex]         
+
+    ClassDef(AliTOFHitMap,0) // Implements HitMap as a 1-dim array
+};
+#endif 
diff --git a/TOF/AliTOFSDigit.cxx b/TOF/AliTOFSDigit.cxx
new file mode 100644 (file)
index 0000000..665d351
--- /dev/null
@@ -0,0 +1,272 @@
+//_________________________________________________________________________
+//  TOF digit: member variables 
+//  fSector  : TOF sector
+//  fPlate   : TOF plate
+//  fStrip   : strips number
+//  fPadx    : pad number along x
+//  fPadz    : pad number along z
+//  fTdc     : TArrayF of TDC values
+//  fAdc     : TArrayF of ADC values
+//              
+//  Getters, setters and member functions  defined here
+//
+//*-- Authors: F. Pierella, A. Seganti, D. Vicinanza
+
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *    
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  * 
+ **************************************************************************/
+
+#include <iostream.h>
+
+#include "TArrayF.h"
+#include "TArrayI.h"
+
+#include "AliTOF.h"
+#include "AliTOFSDigit.h"
+#include "AliTOFConstants.h"
+#include "AliRun.h"
+#include "AliMC.h"
+
+ClassImp(AliTOFSDigit)
+
+////////////////////////////////////////////////////////////////////////
+AliTOFSDigit::AliTOFSDigit()
+{
+//
+// default ctor
+//
+  fNDigits = 0;
+  fTdc = 0;
+  fAdc = 0;
+  fTracks = 0;
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFSDigit::AliTOFSDigit(Int_t tracknum, Int_t *vol,Float_t *digit)
+{
+//
+// Constructor of digit object
+//
+  fSector = vol[0];
+  fPlate  = vol[1];
+  fStrip  = vol[2];
+  fPadx   = vol[3];
+  fPadz   = vol[4];
+  fNDigits = 1;
+  fTdc = new TArrayF(fNDigits);
+  (*fTdc)[0] = digit[0];
+  fAdc = new TArrayF(fNDigits);
+  (*fAdc)[0] = digit[1];
+  fTracks = new TArrayI(kMAXDIGITS*fNDigits);
+  (*fTracks)[0] = tracknum;
+  for (Int_t i = 1; i <kMAXDIGITS*fNDigits; i++) {
+    (*fTracks)[i] = -1;
+  }
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFSDigit::AliTOFSDigit(const AliTOFSDigit & digit)
+{
+  // 
+  // copy ctor for AliTOFSDigit object
+  //
+  fSector = digit.fSector;
+  fPlate  = digit.fPlate;
+  fStrip  = digit.fStrip;
+  fPadx   = digit.fPadx;
+  fPadz   = digit.fPadz;
+  fNDigits = digit.fNDigits;
+  fTdc = new TArrayF(*digit.fTdc);  
+  fAdc = new TArrayF(*digit.fAdc);
+  fTracks = new TArrayI(*digit.fTracks);
+}
+
+////////////////////////////////////////////////////////////////////////
+AliTOFSDigit::AliTOFSDigit(Int_t sector, Int_t plate, Int_t strip, Int_t padx,
+Int_t padz, Float_t tdc, Float_t adc)
+{
+//
+// Constructor for sdigit
+//
+  fSector = sector;
+  fPlate  = plate;
+  fStrip  = strip;
+  fPadx   = padx;
+  fPadz   = padz;  
+  fNDigits = 1;
+  fTdc = new TArrayF(fNDigits);
+  (*fTdc)[0] = tdc;   
+  fAdc = new TArrayF(fNDigits);
+  (*fAdc)[0] = tdc;   
+// no tracks were specified, set them to -1
+  fTracks = new TArrayI(kMAXDIGITS*fNDigits);
+  for (Int_t i = 0; i <kMAXDIGITS*fNDigits; i++) {
+    (*fTracks)[i] = -1;
+  }
+}
+   
+////////////////////////////////////////////////////////////////////////
+void AliTOFSDigit::GetLocation(Int_t *Loc) const
+{
+//
+// Get the coordinates of the digit
+// in terms of Sector - Plate - Strip - Pad
+//
+
+   Loc[0]=fSector;
+   Loc[1]=fPlate;
+   Loc[2]=fStrip;
+   Loc[3]=fPadx;
+   Loc[4]=fPadz;
+}
+
+////////////////////////////////////////////////////////////////////////
+void AliTOFSDigit::Update(Int_t tdc, Int_t adc, Int_t track)
+{
+//
+// Add charge and track
+//
+
+  Int_t sameTime = -1;
+
+  for (Int_t i = 0; i < fNDigits; i++) {
+    if (TMath::Abs(tdc-fTdc->At(i)) < AliTOFConstants::fgkTimeDiff) {
+      sameTime = i;
+      break;
+    }
+  }
+
+  if (sameTime >= 0) {
+    (*fAdc)[sameTime] += static_cast<Float_t>(adc);
+// update track - find the first -1  value and replace it by the
+// track number
+    for (Int_t iTrack=0; iTrack<kMAXDIGITS; iTrack++) {
+      if ((*fTracks)[sameTime*kMAXDIGITS+iTrack] == -1) {
+       (*fTracks)[sameTime*kMAXDIGITS+iTrack] = track;
+       break;
+      }
+// write warning about many tracks going to this pad
+      if (iTrack == kMAXDIGITS) {
+       cerr<<"WARNING: AliTOFSDigit::Update  Many hits in the padhit"<<endl;
+       cerr<<"         ";
+//     PrintPad();
+      }
+    }
+  } else {
+// add new time slot
+    fNDigits++;
+    fTdc->Set(fNDigits);
+    (*fTdc)[fNDigits-1] = tdc;
+    fAdc->Set(fNDigits);
+    (*fAdc)[fNDigits-1] = adc;
+    fTracks->Set(fNDigits*kMAXDIGITS);
+    (*fTracks)[(fNDigits-1)*kMAXDIGITS] = track;
+    for (Int_t i = 1; i <kMAXDIGITS; i++) {
+      (*fTracks)[(fNDigits-1)*kMAXDIGITS+i] = -1;
+    }
+  }
+
+}
+////////////////////////////////////////////////////////////////////////
+AliTOFSDigit::~AliTOFSDigit()
+{
+//
+// dtor
+//
+  delete fTdc;
+  delete fAdc;
+  delete fTracks;
+}
+
+////////////////////////////////////////////////////////////////////////
+
+Int_t AliTOFSDigit::GetTotPad() const
+{
+//
+// Get the "total" index of the pad inside a Sector
+// starting from the digits data.
+//
+
+  AliTOF* tof;
+  
+  if(gAlice){
+     tof =(AliTOF*) gAlice->GetDetector("TOF");
+  }else{
+     printf("AliTOFSDigit::GetTotPad - No AliRun object present, exiting");
+     return 0;
+  }
+  
+  Int_t pad = fPadx+tof->GetNpadX()*(fPadz-1);
+  Int_t before=0;
+
+  switch(fPlate){ 
+  case 1: before = 0;
+          break;
+  case 2: before = tof->GetNStripC();
+          break;
+  case 3: before = tof->GetNStripB() + tof->GetNStripC();
+          break;
+  case 4: before = tof->GetNStripA() + tof->GetNStripB() + tof->GetNStripC();
+          break;
+  case 5: before = tof->GetNStripA() + 2*tof->GetNStripB() + tof->GetNStripC();
+          break;
+  }
+  
+  Int_t strip = fStrip+before;
+  Int_t padTot = tof->GetPadXStr()*(strip-1)+pad;
+  return padTot;
+}
+
+////////////////////////////////////////////////////////////////////////
+//void AliTOFSDigit::AddTrack(Int_t track)
+//{
+//
+// Add a new and different track to the digit -- but to which digit??
+// do not implemet this function
+//
+////////////////////////////////////////////////////////////////////////
+
+// Overloading of Streaming, Sum and Comparison operators
+
+////////////////////////////////////////////////////////////////////////
+/*
+Bool_t AliTOFSDigit::operator==(AliTOFSDigit const &digit) const
+{
+//
+// Overloading of Comparison operator
+//   
+ if (fSector==digit.fSector &&
+     fPlate==digit.fPlate &&
+     fStrip==digit.fStrip &&
+     fPadx==digit.fPadx &&
+     fPadz==digit.fPadz &&
+     fTdc==digit.fTdc &&
+     fAdc==digit.fAdc) return kTRUE;
+     else return kFALSE;
+}
+*/
+////////////////////////////////////////////////////////////////////////
+/*
+ostream& operator << (ostream& out, const AliTOFSDigit &digit)
+{
+//
+// Output streamer: output of the digit data
+//
+out << "Sector " << digit.fSector << ", Plate " << digit.fPlate << ", Strip " << digit.fStrip << endl;
+out << "Padx" << digit.fPadx << ", Padz " << digit.fPadz << endl;
+out << "TDC " << digit.fTdc->At(0) << ", ADC "<< digit.fAdc->At(0) << endl;
+return out;
+}
+*/
diff --git a/TOF/AliTOFSDigit.h b/TOF/AliTOFSDigit.h
new file mode 100644 (file)
index 0000000..d92bfdc
--- /dev/null
@@ -0,0 +1,75 @@
+#ifndef ALITOFSDIGIT_H
+#define ALITOFSDIGIT_H
+
+////////////////////////////////////////////////
+//                                            //
+//  Class for TOF SDigits                     //
+//                                            //
+////////////////////////////////////////////////
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#include "TArrayF.h"
+#include "TArrayI.h"
+#include "AliDigit.h"
+
+//class TArrayF;
+//class TArrayI;
+
+// number 3 is a legacy from AliDigit object
+const Int_t kMAXDIGITS = 3;
+
+class AliTOFSDigit : public TObject {
+
+  //overloading of the streamer << operator
+//friend ostream& operator << ( ostream& , const AliTOFSDigit&) ;
+
+ public:
+ AliTOFSDigit();
+  AliTOFSDigit(Int_t tracknum, Int_t* vol, Float_t* digit);
+// new ctor for sdigits
+  AliTOFSDigit(Int_t sector, Int_t plate, Int_t strip, Int_t padx, Int_t padz, Float_t tdc, Float_t adc);
+// copy ctor
+  AliTOFSDigit(const AliTOFSDigit & digit) ;
+  virtual ~AliTOFSDigit();
+  void            GetLocation(Int_t* Loc) const;
+  Int_t           GetTotPad() const;
+
+  void Update(Int_t tdc, Int_t adc, Int_t track);
+
+// getters for AliTOFSDigit object 
+  Int_t   GetNDigits() const    {return fNDigits;}
+  Float_t GetTdc(Int_t i) const {return fTdc->At(i);}
+  Float_t GetAdc(Int_t i) const {return fAdc->At(i);}
+//  Int_t   GetNTracks(Int_t i) const {return fTracks[i]->GetSize();}
+  Int_t   GetTrack(Int_t i, Int_t j) const {return fTracks->At(i*kMAXDIGITS+j);}
+  Int_t   GetSector() const     {return fSector;}
+  Int_t   GetPlate()  const     {return fPlate;}
+  Int_t   GetStrip()  const     {return fStrip;}
+  Int_t   GetPadx()   const     {return fPadx;}
+  Int_t   GetPadz()   const     {return fPadz;}
+
+protected:
+  Int_t   fSector;  // number of sector
+  Int_t   fPlate;   // number of plate
+  Int_t   fStrip;   // number of strip
+  Int_t   fPadx;    // number of pad along x
+  Int_t   fPadz;    // number of pad along z
+  Int_t   fNDigits;  // dimension of fTdc array
+  TArrayF *fTdc;     // tdc values for sdigit
+  TArrayF *fAdc;     // adc values for sdigit
+  TArrayI *fTracks;  // contributing tracks, kMAXDIGITS entries per
+                     // 1 tdc value
+
+//  Float_t *fTdc;    //[fNDigits] tdc values for sdigit
+//  Float_t *fAdc;    //[fNDigits] adc values for sdigit
+//  Int_t **fTracks;  //[fNDigits] contributing tracks, pointers to
+                    //  arrays with track indices
+  ClassDef(AliTOFSDigit,1)  // SDigit for Time Of Flight
+};
+
+#endif /* ALITOFSDIGIT_H */
index 7a0aae6..a87eef6 100644 (file)
@@ -27,7 +27,8 @@
 #include "TSystem.h"
 #include "TFile.h"
 
-#include "AliTOFdigit.h"
+#include "AliTOFHitMap.h"
+#include "AliTOFSDigit.h"
 #include "AliTOFhit.h"
 #include "AliTOF.h"
 #include "AliTOFv1.h"
@@ -56,7 +57,7 @@ ClassImp(AliTOFSDigitizer)
 {
   // ctor
   fNevents = 0 ;     
-  fSDigits = 0 ;
+//  fSDigits = 0 ;
   fHits = 0 ;
 
 }
@@ -80,105 +81,94 @@ ClassImp(AliTOFSDigitizer)
 void AliTOFSDigitizer::Exec(Option_t *option) { 
 
 
-  // Initialise Hit array
-  fHits = new TClonesArray ("AliTOFhit", 1000);
-  fSDigits = new TClonesArray ("AliTOFdigit", 1000);
-
   AliTOF *TOF = (AliTOF *) gAlice->GetDetector ("TOF");
 
+  if (!TOF) {
+    Error("AliTOFSDigitizer","TOF not found");
+    return;
+  }
+
   if (fNevents == 0)
-    fNevents = (Int_t) gAlice->TreeE ()->GetEntries ();
+    fNevents = (Int_t) gAlice->TreeE()->GetEntries();
 
-  for (Int_t ievent = 0; ievent < fNevents; ievent++)
-    {
-      gAlice->GetEvent (ievent);
-      if (gAlice->TreeH () == 0)
-       return;
-      if (gAlice->TreeS () == 0)
-       gAlice->MakeTree ("S");
+  for (Int_t ievent = 0; ievent < fNevents; ievent++) {
+    gAlice->GetEvent(ievent);
+    TTree *TH = gAlice->TreeH ();
+    if (!TH)
+      return;
+    if (gAlice->TreeS () == 0)
+      gAlice->MakeTree ("S");
 
       
-            //Make branches
-      char branchname[20];
-       sprintf (branchname, "%s", TOF->GetName ());
-      //Make branch for digits
-        TOF->MakeBranch ("S");
+    //Make branches
+    char branchname[20];
+    sprintf (branchname, "%s", TOF->GetName ());
+    //Make branch for digits
+    TOF->MakeBranch ("S");
     
-       //Now made SDigits from hits
-
-      Int_t    tracks[3];    // track info
-      Int_t    vol[5];       // location for a digit
-      Float_t  digit[2];     // TOF digit variables
-      Int_t hit, nbytes;
-      TParticle *particle;
-      AliTOFhit *tofHit;
-      TClonesArray *TOFhits = TOF->Hits();
-
-
-      // Event ------------------------- LOOP  
-
-
-      if (TOF)
-       {
-         TOFhits = TOF->Hits ();
-         TTree *TH = gAlice->TreeH ();
-         Stat_t ntracks = TH->GetEntries ();
-         for (Int_t track = 0; track < ntracks; track++)
-           {
-             gAlice->ResetHits ();
-             nbytes += TH->GetEvent (track);
-             particle = gAlice->Particle (track);
-             Int_t nhits = TOFhits->GetEntriesFast ();
-
-             for (hit = 0; hit < nhits; hit++)
-               {
-                 tofHit = (AliTOFhit *) TOFhits->UncheckedAt(hit);
-               vol[0] = tofHit->GetSector();
-               vol[1] = tofHit->GetPlate();
-               vol[2] = tofHit->GetStrip();
-               vol[3] = tofHit->GetPadx();
-               vol[4] = tofHit->GetPadz();
-
-               // 95% of efficiency to be inserted here
-               // edge effect to be inserted here
-               // cross talk  to be inserted here
-
-               Float_t idealtime = tofHit->GetTof(); // unit s
-               idealtime *= 1.E+12;  // conversion from s to ps
-                              // fTimeRes is given usually in ps
-//             Float_t tdctime   = gRandom->Gaus(idealtime, fTimeRes); 
-               Float_t tdctime   = gRandom->Gaus(idealtime, TOF->GetTimeRes());        
-               digit[0] = tdctime;
-
-               // typical Landau Distribution to be inserted here
-               // instead of Gaussian Distribution
-               Float_t idealcharge = tofHit->GetEdep();
-//             Float_t adccharge = gRandom->Gaus(idealcharge, fChrgRes);
-               Float_t adccharge = gRandom->Gaus(idealcharge, TOF->GetChrgRes());
-               digit[1] = adccharge;
-               Int_t tracknum = tofHit->GetTrack();
-               tracks[0] = tracknum;
-               tracks[1] = 0;
-               tracks[2] = 0;
-
-               // check if two digit are on the same pad; in that case we sum
-               // the two or more digits
-//             Bool_t overlap = AliTOF::CheckOverlap(vol, digit, tracknum);
-               Bool_t overlap = TOF->CheckOverlap(vol, digit, tracknum);
-               if(!overlap) 
-                 //                    new ((*fSDigits)[nSdigits++]) AliTOFdigit(tracks, vol, digit);
-                 TOF->AddSDigit(tracks, vol, digit);
-//             printf("Sect. %d, Plate %d, PadX %d, PadZ %d, Strip %d, Tdc %f",vol[0],vol[1],vol[2],vol[3],vol[5],digit[0]);
-               } // end loop on hits for the current track
-
-            } // end loop on ntracks
-
-         } // close if TOF switched ON
+    //Now made SDigits from hits
+
+    Int_t    vol[5];       // location for a digit
+    Float_t  digit[2];     // TOF digit variables
+    TParticle *particle;
+    AliTOFhit *tofHit;
+    TClonesArray *TOFhits = TOF->Hits();
+
+// create hit map
+    AliTOFHitMap *hitMap = new AliTOFHitMap(TOF->SDigits());
+
+    Int_t ntracks = static_cast<Int_t>(TH->GetEntries());
+    for (Int_t track = 0; track < ntracks; track++)
+    {
+      gAlice->ResetHits();
+      TH->GetEvent(track);
+      particle = gAlice->Particle(track);
+      Int_t nhits = TOFhits->GetEntriesFast();
+
+      for (Int_t hit = 0; hit < nhits; hit++)
+      {
+       tofHit = (AliTOFhit *) TOFhits->UncheckedAt(hit);
+       vol[0] = tofHit->GetSector();
+       vol[1] = tofHit->GetPlate();
+       vol[2] = tofHit->GetStrip();
+       vol[3] = tofHit->GetPadx();
+       vol[4] = tofHit->GetPadz();
+
+       // 95% of efficiency to be inserted here
+       // edge effect to be inserted here
+       // cross talk  to be inserted here
+
+       Float_t idealtime = tofHit->GetTof(); // unit s
+       idealtime *= 1.E+12;  // conversion from s to ps
+       // fTimeRes is given usually in ps
+       Float_t tdctime   = gRandom->Gaus(idealtime, TOF->GetTimeRes());
+       digit[0] = tdctime;
+
+       // typical Landau Distribution to be inserted here
+       // instead of Gaussian Distribution
+       Float_t idealcharge = tofHit->GetEdep();
+       Float_t adccharge = gRandom->Gaus(idealcharge, TOF->GetChrgRes());
+       digit[1] = adccharge;
+       Int_t tracknum = tofHit->GetTrack();
+
+       // check if two digit are on the same pad; in that case we sum
+       // the two or more digits
+       if (hitMap->TestHit(vol) != kEmpty) {
+         AliTOFSDigit *sdig = static_cast<AliTOFSDigit*>(hitMap->GetHit(vol));
+         sdig->Update(tdctime,adccharge,tracknum);
+       } else {
+         TOF->AddSDigit(tracknum, vol, digit);
+         hitMap->SetHit(vol);
+       }
+      } // end loop on hits for the current track
+    } // end loop on ntracks
+
+    delete hitMap;
       
-      gAlice->TreeS()->Reset();
-      gAlice->TreeS()->Fill();
-      gAlice->TreeS()->Write(0,TObject::kOverwrite) ;
-    }                          //event loop
+    gAlice->TreeS()->Reset();
+    gAlice->TreeS()->Fill();
+    gAlice->TreeS()->Write(0,TObject::kOverwrite) ;
+  }                            //event loop
 
 
 }
@@ -199,55 +189,3 @@ void AliTOFSDigitizer::Print(Option_t* option)const
     cout << "    Writing SDigitis to file  " << (char*) fSDigitsFile.Data() << endl ;
 
 }
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
index 3e49bba..9ec6e56 100644 (file)
@@ -19,15 +19,16 @@ class AliTOFSDigitizer: public TTask {
 
 public:
   AliTOFSDigitizer() ;          // ctor
-  AliTOFSDigitizer(char* HeaderFile,char *SdigitsFile = 0) ; 
+  AliTOFSDigitizer(char* HeaderFile, char *SdigitsFile = 0) ; 
 
   virtual ~AliTOFSDigitizer() ; // dtor
   // Int_t    Digitize(Float_t Energy);
 
-  char *GetSDigitsFile()const{return (char*) fSDigitsFile.Data();}  
+//  char *GetSDigitsFile() const {return const_cast<char*>(fSDigitsFile.Data());}  
+  const char *GetSDigitsFile() const {return fSDigitsFile.Data();}  
   virtual void  Exec(Option_t *option); 
-  void SetNEvents(Int_t Nevents){fNevents = Nevents;}
-  Stat_t GetNEvents(){return fNevents;}
+  void  SetNEvents(Int_t Nevents) {fNevents = Nevents;}
+  Int_t GetNEvents() const {return fNevents;}
   void SetSDigitsFile(char * file ) ;
   virtual void Print(Option_t* option) const ;
   TClonesArray *SDigits() const {return fSDigits;}
@@ -36,16 +37,16 @@ public:
 
 
 private:
-  Int_t   fNevents ;        // Number of events to digitize
-  TString fSDigitsFile ;    //output file 
-  TClonesArray *fSDigits      ; // List of summable digits
-  TClonesArray *fHits      ; // List of summable digits
-  TString fHeadersFile ;    //input file
+  Int_t   fNevents;         // Number of events to digitize
+  TString fSDigitsFile;     // output file 
+  TClonesArray *fSDigits;   // array of summable digits
+  TClonesArray *fHits;      // array of summable digits
+  TString fHeadersFile;     // input file
 
  protected:
 
 
-  ClassDef(AliTOFSDigitizer,1)  // description 
+  ClassDef(AliTOFSDigitizer,1)  // creates TOF SDigits
 
 };
 
index f36e355..194d21f 100644 (file)
@@ -1033,8 +1033,7 @@ void AliTOFv0::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
-    if(ppad[1]>1.) ppad[1]=1.;
-    if(ppad[1]<-1.) ppad[1]=-1.;
+
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index 1fd4682..cd4e0a7 100644 (file)
 
 /*
 $Log$
+
+Revision 1.26  2001/11/13 14:36:40  vicinanz
+Updated check for ppad[1] range
+
 Revision 1.24  2001/09/27 10:39:20  vicinanz
 SDigitizer and Merger added
 
@@ -1046,8 +1050,10 @@ void AliTOFv1::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
+
     if(ppad[1]>1.) ppad[1]=1.;
     if(ppad[1]<-1.) ppad[1]=-1.;
+
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index 55ec8fa..dea8bb3 100644 (file)
@@ -1060,8 +1060,9 @@ void AliTOFv2::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
-    if(ppad[1]>1.) ppad[1]=1.;
-    if(ppad[1]<-1.) ppad[1]=-1.;
+
+    if (ppad[1] > 1.) ppad[1]=1;
+    if (ppad[1] < -1.) ppad[1]=-1;
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index 00a93a2..460d926 100644 (file)
@@ -1023,8 +1023,7 @@ void AliTOFv2FHoles::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
-    if(ppad[1]>1.) ppad[1]=1.;
-    if(ppad[1]<-1.) ppad[1]=-1.;
+
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index ce731bb..7a7b87a 100644 (file)
@@ -1047,8 +1047,7 @@ void AliTOFv3::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
-    if(ppad[1]>1.) ppad[1]=1.;
-    if(ppad[1]<-1.) ppad[1]=-1.;
+
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index 20491c7..cfe5ff7 100644 (file)
@@ -1045,8 +1045,7 @@ void AliTOFv4::StepManager()
  
     gMC->Gmtod(xm,xpad,1);
     gMC->Gmtod(pm,ppad,2);
-    if(ppad[1]>1.) ppad[1]=1.;
-    if(ppad[1]<-1.) ppad[1]=-1.;
+
     incidenceAngle = TMath::ACos(ppad[1])*kRaddeg;
 
     z = pos[2];
index c47e2a4..0fd61ed 100644 (file)
@@ -9,7 +9,7 @@ PACKAGE = TOF
 
 # C++ sources
 
-SRCS  = AliTOF.cxx  AliTOFv0.cxx  AliTOFv1.cxx  AliTOFv2.cxx  AliTOFv3.cxx  AliTOFv4.cxx AliTOFhit.cxx  AliTOFdigit.cxx AliTOFRawSector.cxx  AliTOFRoc.cxx  AliTOFRawDigit.cxx AliTOFDigitizer.cxx AliTOFSDigitizer.cxx AliTOFMerger.cxx  AliTOFv2FHoles.cxx
+SRCS  = AliTOF.cxx  AliTOFv0.cxx  AliTOFv1.cxx  AliTOFv2.cxx  AliTOFv3.cxx  AliTOFv4.cxx AliTOFhit.cxx  AliTOFdigit.cxx AliTOFRawSector.cxx  AliTOFRoc.cxx  AliTOFRawDigit.cxx AliTOFDigitizer.cxx AliTOFSDigitizer.cxx AliTOFMerger.cxx  AliTOFv2FHoles.cxx AliTOFSDigit.cxx AliTOFHitMap.cxx AliTOFConstants.cxx
 
 # C++ Headers
 
index d54ce3d..48ba3ac 100644 (file)
@@ -14,7 +14,7 @@
 #pragma link C++ class  AliTOFv2+;
 #pragma link C++ class  AliTOFv3+;
 #pragma link C++ class  AliTOFv4+;
-#pragma link C++ class  AliTOFv2FHoles+;
+//#pragma link C++ class  AliTOFv2FHoles+;
 
 #pragma link C++ class  AliTOFhit+;
 
@@ -25,5 +25,8 @@
 #pragma link C++ class  AliTOFDigitizer+;
 #pragma link C++ class  AliTOFSDigitizer+;
 #pragma link C++ class  AliTOFMerger+;
+#pragma link C++ class  AliTOFSDigit+;
+#pragma link C++ class  AliTOFHitMap+;
+#pragma link C++ class  AliTOFConstants+;
 
 #endif
index b14c6c1..a3ed209 100644 (file)
@@ -1,6 +1,6 @@
-SRCS  = AliTOF.cxx  AliTOFv0.cxx  AliTOFv1.cxx  AliTOFv2.cxx  AliTOFv3.cxx  AliTOFv4.cxx AliTOFhit.cxx  AliTOFdigit.cxx AliTOFRawSector.cxx  AliTOFRoc.cxx  AliTOFRawDigit.cxx AliTOFDigitizer.cxx AliTOFSDigitizer.cxx AliTOFMerger.cxx  AliTOFv2FHoles.cxx
+SRCS  = AliTOF.cxx  AliTOFv0.cxx  AliTOFv1.cxx  AliTOFv2.cxx  AliTOFv3.cxx  AliTOFv4.cxx AliTOFhit.cxx  AliTOFdigit.cxx AliTOFRawSector.cxx  AliTOFRoc.cxx  AliTOFRawDigit.cxx AliTOFDigitizer.cxx AliTOFSDigitizer.cxx AliTOFMerger.cxx  AliTOFv2FHoles.cxx AliTOFSDigit.cxx AliTOFHitMap.cxx AliTOFConstants.cxx
 
 
 HDRS:= $(SRCS:.cxx=.h)
 
-DHDR=TOFLinkDef.h
\ No newline at end of file
+DHDR=TOFLinkDef.h