Correct treatment of regions of interest in TPC
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 7 Jan 2004 11:03:01 +0000 (11:03 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 7 Jan 2004 11:03:01 +0000 (11:03 +0000)
STEER/AliDigitizer.cxx
STEER/AliDigitizer.h
STEER/AliSimulation.cxx
STEER/AliSimulation.h
TPC/AliTPCDigitizer.cxx

index 1510d4f..1e531d7 100644 (file)
@@ -35,7 +35,8 @@ ClassImp(AliDigitizer)
 //_______________________________________________________________________
 AliDigitizer::AliDigitizer(const Text_t* name, const Text_t* title):
   TTask(name,title),
-  fManager(0)
+  fManager(0),
+  fRegionOfInterest(kTRUE)
 {
   //
   // Default ctor with name and title
@@ -45,7 +46,8 @@ AliDigitizer::AliDigitizer(const Text_t* name, const Text_t* title):
 //_______________________________________________________________________
 AliDigitizer::AliDigitizer(const AliDigitizer &dig):
   TTask(dig.GetName(),dig.GetTitle()),
-  fManager(0)
+  fManager(0),
+  fRegionOfInterest(kTRUE)
 {
   //
   // Copy ctor with
@@ -63,7 +65,8 @@ void AliDigitizer::Copy(TObject &) const
 AliDigitizer::AliDigitizer(AliRunDigitizer *manager, 
                            const Text_t* name, const Text_t* title):
   TTask(name,title),
-  fManager(manager)
+  fManager(manager),
+  fRegionOfInterest(kFALSE)
 {
   //
   // ctor with name and title
index c4a534c..c7f1313 100644 (file)
@@ -34,6 +34,7 @@ class AliDigitizer: public TTask {
       
     virtual ~AliDigitizer();
     virtual Bool_t Init() {return kTRUE;}
+    void SetRegionOfInterest(Bool_t flag) {fRegionOfInterest = flag;};
 //    virtual void Digitize() = 0;
 
 protected:
@@ -41,6 +42,7 @@ protected:
     void Copy(TObject &dig) const;
 
     AliRunDigitizer *fManager;   // Pointer to the Digitizer manager
+    Bool_t fRegionOfInterest;    // Flag for digitization only in region of interest
     
     ClassDef(AliDigitizer,1) // Base class for detector digitizers
 };
index 68ea44b..8b8cddf 100644 (file)
@@ -79,6 +79,7 @@
 #include "AliModule.h"
 #include "AliGenerator.h"
 #include "AliRunDigitizer.h"
+#include "AliDigitizer.h"
 #include <TObjString.h>
 
 
@@ -154,6 +155,7 @@ void AliSimulation::Init()
   fConfigFileName = "Config.C";
   fGAliceFileName = "galice.root";
   fBkgrdFileNames = new TObjArray;
+  fRegionOfInterest = kTRUE;
 
   fRunLoader = NULL;
 }
@@ -217,7 +219,7 @@ Bool_t AliSimulation::Run(Int_t nEvents)
   fRunLoader->LoadgAlice();
   gAlice = fRunLoader->GetAliRun();
   if (!gAlice) {
-    Error("GetLoadersAndDetectors", "no gAlice object found in file %s", 
+    Error("Run", "no gAlice object found in file %s", 
          fGAliceFileName.Data());
     return kFALSE;
   }
@@ -255,6 +257,9 @@ Bool_t AliSimulation::RunSimulation()
 {
 // run the generation and simulation
 
+  TStopwatch stopwatch;
+  stopwatch.Start();
+
   Info("RunSimulation", "initializing gAlice with config file %s",
        fConfigFileName.Data());
   gAlice->Init(fConfigFileName.Data());
@@ -275,9 +280,12 @@ Bool_t AliSimulation::RunSimulation()
     gAlice->Generator()->SetTrackingFlag(0);
   }
 
-  Info("Run", "running gAlice");
+  Info("RunSimulation", "running gAlice");
   gAlice->Run(fNEvents);
 
+  Info("RunSimulation", "execution time:");
+  stopwatch.Print();
+
   return kTRUE;
 }
 
@@ -286,6 +294,9 @@ Bool_t AliSimulation::RunSDigitization(const TString& detectors)
 {
 // run the digitization and produce summable digits
 
+  TStopwatch stopwatch;
+  stopwatch.Start();
+
   TString detStr = detectors;
   TObjArray* detArray = gAlice->Detectors();
   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
@@ -304,6 +315,9 @@ Bool_t AliSimulation::RunSDigitization(const TString& detectors)
     if (fStopOnError) return kFALSE;
   }
 
+  Info("RunSDigitization", "execution time:");
+  stopwatch.Print();
+
   return kTRUE;
 }
 
@@ -314,6 +328,9 @@ Bool_t AliSimulation::RunDigitization(const TString& detectors,
 {
 // run the digitization and produce digits from sdigits
 
+  TStopwatch stopwatch;
+  stopwatch.Start();
+
   Int_t nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
   Int_t signalPerBkgrd = 1;
   if (nStreams > 1) signalPerBkgrd = fBkgrdFileNames->At(0)->GetUniqueID();
@@ -333,9 +350,12 @@ Bool_t AliSimulation::RunDigitization(const TString& detectors,
     if (!det || !det->IsActive()) continue;
     if (IsSelected(det->GetName(), detStr) && 
        !IsSelected(det->GetName(), detExcl)) {
-      if (!det->CreateDigitizer(manager)) {
+      AliDigitizer* digitizer = det->CreateDigitizer(manager);
+      if (!digitizer) {
        Error("RunDigitization", "no digitizer for %s", det->GetName());
        if (fStopOnError) return kFALSE;
+      } else {
+       digitizer->SetRegionOfInterest(fRegionOfInterest);
       }
     }
   }
@@ -352,6 +372,9 @@ Bool_t AliSimulation::RunDigitization(const TString& detectors,
   }
   delete manager;
 
+  Info("RunDigitization", "execution time:");
+  stopwatch.Print();
+
   return kTRUE;
 }
 
@@ -360,6 +383,9 @@ Bool_t AliSimulation::RunHitsDigitization(const TString& detectors)
 {
 // run the digitization and produce digits from hits
 
+  TStopwatch stopwatch;
+  stopwatch.Start();
+
   TString detStr = detectors;
   TObjArray* detArray = gAlice->Detectors();
   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
@@ -378,6 +404,9 @@ Bool_t AliSimulation::RunHitsDigitization(const TString& detectors)
     if (fStopOnError) return kFALSE;
   }
 
+  Info("RunHitsDigitization", "execution time:");
+  stopwatch.Print();
+
   return kTRUE;
 }
 
index 61b12be..d09e93c 100644 (file)
@@ -27,6 +27,7 @@ public:
   void           SetMakeSDigits(const char* detectors) 
                    {fMakeSDigits = detectors;};
   void           MergeWith(const char* fileName, Int_t nSignalPerBkgrd = 1);
+  void           SetRegionOfInterest(Bool_t flag) {fRegionOfInterest = flag;};
   void           SetMakeDigits(const char* detectors)
                    {fMakeDigits = detectors;};
   void           SetMakeDigitsFromHits(const char* detectors)
@@ -56,6 +57,7 @@ private:
   TString        fConfigFileName;     // name of the config file
   TString        fGAliceFileName;     // name of the galice file
   TObjArray*     fBkgrdFileNames;     // names of background files for merging
+  Bool_t         fRegionOfInterest;   // digitization in region of interest
 
   AliRunLoader*  fRunLoader;          //! current run loader object
 
index 71162bc..252cfdb 100644 (file)
@@ -132,11 +132,11 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
     masks[i]= fManager->GetMask(i);
   Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
   Int_t **ptr=  new Int_t*[nInputs];       //pointers to teh expanded tracks array
+  Bool_t *active=  new Bool_t[nInputs];    //flag for active input segments
 
   
   //create digits array for given sectors
   // make indexes
-  Stat_t nentries = 0;//number of entries in TreeS
   
   AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
   for (Int_t i1=0;i1<nInputs; i1++)
@@ -154,7 +154,6 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
             <<" input "<< i1<<endl;
         return;
        }
-      if (i1 == 0) nentries=treear->GetEntries();
     
       if (treear->GetIndex()==0)  
        treear->BuildIndex("fSegmentID","fSegmentID");
@@ -174,7 +173,7 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
      ogime->MakeTree("D");
      tree  = ogime->TreeD();
    }
-  tree->Branch("Segment","AliSimDigits",&digrow);
+  TBranch* branch = tree->Branch("Segment","AliSimDigits",&digrow);
   //
 
   param->SetZeroSup(2);
@@ -183,39 +182,41 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
   //
   //Loop over segments of the TPC
     
-  for (Int_t n=0; n<nentries; n++) 
+  for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++) 
    {
-    rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
-    gime = rl->GetLoader("TPCLoader");
-    gime->TreeS()->GetEvent(n);
-    digarr[0]->ExpandBuffer();
-    digarr[0]->ExpandTrackBuffer();
-           
-    for (Int_t i=1;i<nInputs; i++) 
+    Int_t sec, row;
+    if (!param->AdjustSectorRow(segmentID,sec,row)) 
+     {
+      cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
+      continue;
+     }
+
+    digrow->SetID(segmentID);
+
+    Int_t nrows = 0;
+    Int_t ncols = 0;
+
+    Bool_t digitize = kFALSE;
+    for (Int_t i=0;i<nInputs; i++) 
      { 
 
       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
       gime = rl->GetLoader("TPCLoader");
       
-      gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
-      if ((digarr[0]->GetID()-digarr[i]->GetID())>0) 
-       printf("problem - not corresponding segment in background event\n");
-      
-      digarr[i]->ExpandBuffer();
-      digarr[i]->ExpandTrackBuffer();
-      
+      if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
+       digarr[i]->ExpandBuffer();
+       digarr[i]->ExpandTrackBuffer();
+        nrows = digarr[i]->GetNRows();
+        ncols = digarr[i]->GetNCols();
+       active[i] = kTRUE;
+       if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
+      } else {
+       active[i] = kFALSE;
+      }
+      if (fRegionOfInterest && !digitize) break;
      }   
-    Int_t sec, row;
-    if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) 
-     {
-      cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
-      continue;
-     }
+    if (!digitize) continue;
 
-    digrow->SetID(digarr[0]->GetID());
-
-    Int_t nrows = digarr[0]->GetNRows();
-    Int_t ncols = digarr[0]->GetNCols();
     digrow->Allocate(nrows,ncols);
     digrow->AllocateTrack(3);
 
@@ -226,10 +227,10 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
     Int_t nElems = nrows*ncols;     
  
     for (Int_t i=0;i<nInputs; i++)
-     { 
+     if (active[i]) { 
        pdig[i] = digarr[i]->GetDigits();
        ptr[i]  = digarr[i]->GetTracks();
-     }
+      }
      
     Short_t *pdig1= digrow->GetDigits();
     Int_t   *ptr1= digrow->GetTracks() ;
@@ -242,7 +243,7 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
        q=0;
        labptr=0;
        // looop over digits 
-        for (Int_t i=0;i<nInputs; i++)
+        for (Int_t i=0;i<nInputs; i++) if (active[i]) 
          { 
           //          q  += digarr[i]->GetDigitFast(rows,col);
             q  += *(pdig[i]);
@@ -282,9 +283,10 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
     
     digrow->CompresBuffer(1,zerosup);
     digrow->CompresTrackBuffer(1);
+    branch->SetAddress(&digrow);
     tree->Fill();
     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
-   } //for (Int_t n=0; n<nentries; n++) 
+   } //for (Int_t n=0; n<param->GetNRowsTotal(); n++) 
   
 
   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
@@ -296,7 +298,10 @@ void AliTPCDigitizer::ExecFast(Option_t* option)
   delete digrow;     
   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
   delete []masks;
-  delete digarr;  
+  delete []pdig;
+  delete []ptr;
+  delete []active;
+  delete []digarr;  
 }