CbmMvdDigitizer.cxx

Go to the documentation of this file.

// -------------------------------------------------------------------------
// -----                    CbmMvdDigitizer source file                -----
// -------------------------------------------------------------------------
 
// Includes from MVD
#include "CbmMvdDigitizer.h"
#include "CbmMvdDetector.h"
#include "CbmMvdPileupManager.h"
#include "SensorDataSheets/CbmMvdMimosa26AHR.h"
#include "plugins/tasks/CbmMvdSensorDigitizerTask.h"
#include "tools/CbmMvdGeoHandler.h"
 
#include "CbmMatch.h"
#include "CbmMvdDigi.h"
#include "CbmMvdPoint.h"
 
// Includes from FAIR
#include "FairLogger.h"
#include "FairModule.h"
#include "FairRootManager.h"
 
// Includes from ROOT
#include "TClonesArray.h"
#include "TStopwatch.h"
 
// Includes from C++
#include <cassert>
#include <iomanip>
#include <iostream>
#include <vector>
 
using std::cout;
using std::endl;
using std::vector;
 
// -----   Default constructor   ------------------------------------------
CbmMvdDigitizer::CbmMvdDigitizer()
  : CbmDigitize<CbmMvdDigi>("MvdDigitize")
  , fMode(0)
  , fShowDebugHistos(kFALSE)
  , fNoiseSensors(kFALSE)
  , fDetector(nullptr)
  , fInputPoints(nullptr)
  , fMcPileUp(nullptr)
  , fTmpMatch(nullptr)
  , fTmpDigi(nullptr)
  , fDigiVect()
  , fMatchVect()
  , fPerformanceDigi()
  , fDigiPluginNr(0)
  , fFakeRate(-1.)
  , fNPileup(-1)
  , fNDeltaElect(-1)
  , fDeltaBufferSize(-1)
  , fBgBufferSize(-1)
  , epsilon()
  , fInputBranchName("")
  , fBgFileName("")
  , fDeltaFileName("")
  , fTimer()
  , fPileupManager(nullptr)
  , fDeltaManager(nullptr) {}
// -------------------------------------------------------------------------
 
// -----   Standard constructor   ------------------------------------------
CbmMvdDigitizer::CbmMvdDigitizer(const char* name,
                                 Int_t iMode,
                                 Int_t /*iVerbose*/)
  : CbmDigitize<CbmMvdDigi>(name)
  , fMode(iMode)
  , fShowDebugHistos(kFALSE)
  , fNoiseSensors(kFALSE)
  , fDetector(nullptr)
  , fInputPoints(nullptr)
  , fMcPileUp(nullptr)
  , fTmpMatch(nullptr)
  , fTmpDigi(nullptr)
  , fDigiVect()
  , fMatchVect()
  , fPerformanceDigi()
  , fDigiPluginNr(0)
  , fFakeRate(-1.)
  , fNPileup(0)
  , fNDeltaElect(0)
  , fDeltaBufferSize(-1)
  , fBgBufferSize(-1)
  , epsilon()
  , fInputBranchName("MvdPoint")
  , fBgFileName("")
  , fDeltaFileName("")
  , fTimer()
  , fPileupManager(nullptr)
  , fDeltaManager(nullptr) {}
// -------------------------------------------------------------------------
 
// -----   Destructor   ----------------------------------------------------
CbmMvdDigitizer::~CbmMvdDigitizer() {
 
  if (fMcPileUp) {
    fMcPileUp->Delete();
    delete fMcPileUp;
  }
  delete fPileupManager;
  delete fDeltaManager;
}
// -----------------------------------------------------------------------------
 
// -----   Exec   --------------------------------------------------------------
void CbmMvdDigitizer::Exec(Option_t* /*opt*/) {
  // --- Start timer
  using namespace std;
 
  fTimer.Start();
  GetEventInfo();  // event number and time
 
  BuildEvent();
  Int_t nPoints = fInputPoints->GetEntriesFast();
  Int_t nDigis  = 0;
  if (fInputPoints->GetEntriesFast() > 0) {
    LOG(debug) << "//----------------------------------------//";
    LOG(debug) << fName << ": Send Input";
    fDetector->SendInput(fInputPoints);
    LOG(debug) << fName << ": Execute DigitizerPlugin Nr. " << fDigiPluginNr;
    fDetector->Exec(fDigiPluginNr);
    LOG(debug) << fName << ": End Chain";
 
    // --- Send produced digis to DAQ
    fTmpDigi  = fDetector->GetOutputDigis();
    fTmpMatch = fDetector->GetOutputDigiMatchs();
 
    Int_t nEntries = fTmpDigi->GetEntriesFast();
    for (Int_t index = 0; index < nEntries; index++) {
 
      CbmMvdDigi* digi  = dynamic_cast<CbmMvdDigi*>(fTmpDigi->At(index));
      CbmMvdDigi* digi1 = new CbmMvdDigi(*digi);
      assert(digi1);
      fDigiVect.push_back(digi1);
 
      //     CbmMatch match{*(dynamic_cast<CbmMatch*>(fTmpMatch->At(index)))};
      //     CbmMatch* match1 = new CbmMatch(match);
      CbmMatch* match  = dynamic_cast<CbmMatch*>(fTmpMatch->At(index));
      CbmMatch* match1 = new CbmMatch(*match);
      fMatchVect.push_back(match1);
 
      //digi1->SetMatch(match1);
      SendData(digi1, match1);
      nDigis++;
    }
 
 
    /*
   for (Int_t index = 0; index < fTmpDigi->GetEntriesFast(); index++) {
 
     LOG(info) << "Size: " << fTmpDigi->GetEntriesFast() << ", " 
               << fTmpDigi->GetEntries();
 
     CbmMvdDigi* digi = dynamic_cast<CbmMvdDigi*>(fTmpDigi->Remove(fTmpDigi->At(index)));
     digi->Print();
     fDigiVect.push_back(digi);
        
     CbmMatch* match = dynamic_cast<CbmMatch*>(fTmpMatch->Remove(fTmpMatch->At(index)));
     match->Print();
     fMatchVect.push_back(match);
 
     digi->SetMatch(match);
     SendDigi(digi);
     nDigis++;
   }
*/
    // TODO: (VF) There seem to be no entries in the match array, nor matches
    // attached to the digi object
    LOG(debug) << fName << ": Sent " << nDigis << " digis to DAQ";
 
    //fDigis->AbsorbObjects(fDetector->GetOutputDigis(),0,fDetector->GetOutputArray(fDigiPluginNr)->GetEntriesFast()-1);
    //LOG(debug) << "Total of " << fDigis->GetEntriesFast() << " digis in this Event";
    //if(fVerbose) cout << "Start writing DigiMatchs" << endl;
    //fDigiMatch->AbsorbObjects(fDetector->GetOutputDigiMatchs(),0,fDetector->GetOutputDigiMatchs()->GetEntriesFast()-1);
    //if(fVerbose) cout << "Total of " << fDigiMatch->GetEntriesFast() << " digisMatch in this Event" << endl;
    LOG(debug) << "//----------------------------------------//";
    //LOG(info) << "+ " << setw(20) << GetName() << ": Created: "
    //          << fDigis->GetEntriesFast() << " digis in "
    //          << fixed << setprecision(6) << fTimer.RealTime() << " s";
  }
 
  // --- Event log
  fTimer.Stop();
  LOG(info) << "+ " << setw(15) << GetName() << ": Event " << setw(6) << right
            << fCurrentEvent << " at " << fixed << setprecision(3)
            << fCurrentEventTime << " ns, points: " << nPoints
            << ", digis: " << nDigis << ". Exec time " << setprecision(6)
            << fTimer.RealTime() << " s.";
}
// -----------------------------------------------------------------------------
 
// -----   Init   --------------------------------------------------------------
InitStatus CbmMvdDigitizer::Init() {
  cout << "-I- " << GetName() << ": Initialisation..." << endl;
  cout << endl;
  cout << "---------------------------------------------" << endl;
  cout << "-I- Initialising " << GetName() << " ...." << endl;
 
  // **********  RootManager
  FairRootManager* ioman = FairRootManager::Instance();
  if (!ioman) {
    cout << "-E- " << GetName() << "::Init: No FairRootManager!" << endl;
    return kFATAL;
  }
 
  // **********  Get input arrays
  fInputPoints = (TClonesArray*) ioman->GetObject("MvdPoint");
 
  if (!fInputPoints) {
    LOG(error) << "No MvdPoint branch found. There was no MVD in the "
                  "simulation. Switch this task off";
    return kERROR;
  }
 
 
  // **********  Register output arrays
  RegisterOutput();
  fMcPileUp = new TClonesArray("CbmMvdPoint", 10000);
  ioman->Register("MvdPileUpMC",
                  "Mvd MC Points after Pile Up",
                  fMcPileUp,
                  IsOutputBranchPersistent("MvdPileUpMC"));
 
  fDetector = CbmMvdDetector::Instance();
 
  // **********  Create pileup manager if necessary
  if (fNPileup >= 1 && !(fPileupManager) && fMode == 0) {
    if (fBgFileName == "") {
      cout << "-E- " << GetName() << "::Init: Pileup events needed, but no "
           << " background file name given! " << endl;
      return kERROR;
    }
    fPileupManager =
      new CbmMvdPileupManager(fBgFileName, fInputBranchName, fBgBufferSize);
    if (fPileupManager->GetNEvents() < 2 * fNPileup) {
      cout << "-E- " << GetName()
           << ": The size of your BG-File is insufficient to perform the "
              "requested pileup"
           << endl;
      cout << "    You need at least events > 2* fNPileup." << endl;
      cout << "    Detected: fPileUp = " << fNPileup << ", events in file "
           << fPileupManager->GetNEvents() << endl;
      Fatal(GetName(), "The size of your BG-File is insufficient");
      return kERROR;
    }
  }
 
  // **********   Create delta electron manager if required
  if (fNDeltaElect >= 1 && !(fDeltaManager) && fMode == 0) {
    if (fDeltaFileName == "") {
      cout << "-E- " << GetName() << "::Init: Pileup events needed, but no "
           << " background file name given! " << endl;
      return kERROR;
    }
    fDeltaManager = new CbmMvdPileupManager(
      fDeltaFileName, fInputBranchName, fDeltaBufferSize);
    if (fDeltaManager->GetNEvents() < 2 * fNDeltaElect) {
      cout << "-E- " << GetName()
           << ": The size of your Delta-File is insufficient to perform the "
              "requested pileup"
           << endl;
      cout << "    You need at least events > 2* fNDeltaElect." << endl;
      cout << "    Detected: fNDeltaElect = " << fNDeltaElect
           << ", events in file " << fDeltaManager->GetNEvents() << endl;
      Fatal(GetName(), "The size of your Delta-File is insufficient");
      return kERROR;
    }
  }
  CbmMvdSensorDigitizerTask* digiTask = new CbmMvdSensorDigitizerTask();
 
  fDetector->AddPlugin(digiTask);
  fDigiPluginNr = (UInt_t)(fDetector->GetPluginArraySize());
  if (fShowDebugHistos) fDetector->ShowDebugHistos();
  fDetector->Init();
 
  if (fNoiseSensors) fDetector->SetProduceNoise();
 
  // Screen output
  cout << GetName() << " initialised with parameters: " << endl;
  //PrintParameters();
  cout << "---------------------------------------------" << endl;
 
 
  return kSUCCESS;
}
 
// -----   Virtual public method Reinit   ----------------------------------
InitStatus CbmMvdDigitizer::ReInit() { return kSUCCESS; }
// -------------------------------------------------------------------------
 
 
// -----   Virtual method Finish   -----------------------------------------
void CbmMvdDigitizer::Finish() {
  // cout<< endl << "finishing" << endl;
  fDetector->Finish();
  PrintParameters();
}
// -------------------------------------------------------------------------
 
 
// -----   Private method Reset   ------------------------------------------
void CbmMvdDigitizer::Reset() {}
// -------------------------------------------------------------------------
 
 
// -----   Reset output arrays   -------------------------------------------
void CbmMvdDigitizer::ResetArrays() {
  if (fMcPileUp) fMcPileUp->Delete();
  if (fTmpMatch) fTmpMatch->Delete();
  if (fTmpDigi) fTmpDigi->Delete();
  fDigiVect.clear();
  fMatchVect.clear();
}
// -------------------------------------------------------------------------
 
 
// -----   Private method GetMvdGeometry   ---------------------------------
void CbmMvdDigitizer::GetMvdGeometry() {}
// -------------------------------------------------------------------------
 
 
// -----   Private method PrintParameters   --------------------------------
void CbmMvdDigitizer::PrintParameters() {
 
  using namespace std;
 
  cout.setf(ios_base::fixed, ios_base::floatfield);
  cout << "============================================================"
       << endl;
  cout << "============== Parameters MvdDigitizer ====================="
       << endl;
  cout << "============================================================"
       << endl;
  cout << "=============== End Task ==================================="
       << endl;
}
// -------------------------------------------------------------------------
 
// ---------------  Build Event  ------------------------------------------------------
void CbmMvdDigitizer::BuildEvent() {
 
  // Some frequently used variables
  CbmMvdPoint* point = nullptr;
  Int_t nOrig        = 0;
  Int_t nPile        = 0;
  Int_t nElec        = 0;
 
  // ----- First treat standard input file
  nOrig = fInputPoints->GetEntriesFast();
 
 
  // ----- Then treat event pileup
  if (fNPileup > 0) {
 
    // --- Vector of available background events from pile-up.
    // --- Each event is used only once.
    Int_t nBuffer = fPileupManager->GetNEvents();
    vector<Int_t> freeEvents(nBuffer);
    for (Int_t i = 0; i < nBuffer; i++)
      freeEvents[i] = i;
 
    // --- Loop over pile-up events
    for (Int_t iBg = 0; iBg < fNPileup; iBg++) {
 
      // Select random event from vector and remove it after usage
      Int_t index = gRandom->Integer(freeEvents.size());
 
      Int_t iEvent         = freeEvents[index];
      TClonesArray* points = fPileupManager->GetEvent(iEvent);
      freeEvents.erase(freeEvents.begin() + index);
 
      // Add points from this event to the input arrays
      for (Int_t iPoint = 0; iPoint < points->GetEntriesFast(); iPoint++) {
        point = (CbmMvdPoint*) points->At(iPoint);
        point->SetTrackID(-2);
        nPile++;
        new ((*fInputPoints)[fInputPoints->GetEntriesFast()])
          CbmMvdPoint(*((CbmMvdPoint*) points->At(iPoint)));
      }
 
 
    }  // Pileup event loop
 
  }  // Usage of pile-up
 
 
  // ----- Finally, treat delta electrons
  if (fNDeltaElect > 0) {
 
    // --- Vector of available delta events.
    // --- Each event is used only once.
    Int_t nDeltaBuffer = fDeltaManager->GetNEvents();
    vector<Int_t> freeDeltaEvents(nDeltaBuffer);
    for (Int_t i = 0; i < nDeltaBuffer; i++)
      freeDeltaEvents[i] = i;
 
    // --- Loop over delta events
    for (Int_t it = 0; it < fNDeltaElect; it++) {
 
      // Select random event from vector and remove it after usage
      Int_t indexD          = gRandom->Integer(freeDeltaEvents.size());
      Int_t iEventD         = freeDeltaEvents[indexD];
      TClonesArray* pointsD = fDeltaManager->GetEvent(iEventD);
      freeDeltaEvents.erase(freeDeltaEvents.begin() + indexD);
 
      // Add points from this event to the input arrays
      for (Int_t iPoint = 0; iPoint < pointsD->GetEntriesFast(); iPoint++) {
        point = (CbmMvdPoint*) pointsD->At(iPoint);
        point->SetTrackID(-3);  // Mark the points as delta electron
        new ((*fInputPoints)[fInputPoints->GetEntriesFast()])
          CbmMvdPoint(*((CbmMvdPoint*) pointsD->At(iPoint)));
        nElec++;
      }
 
 
    }  // Delta electron event loop
 
  }  // Usage of delta
 
  for (Int_t i = 0; i < fInputPoints->GetEntriesFast(); i++)
    new ((*fMcPileUp)[i]) CbmMvdPoint(*((CbmMvdPoint*) fInputPoints->At(i)));
 
 
  // ----- At last: Screen output
  LOG(debug) << "+ " << GetName() << "::BuildEvent: original " << nOrig
             << ", pileup " << nPile << ", delta " << nElec << ", total "
             << nOrig + nPile + nElec << " MvdPoints";
}
// -----------------------------------------------------------------------------
 
 
ClassImp(CbmMvdDigitizer);