GGSRootFileService.cpp

/*
 * GGSRootFileService.cpp
 *
 *  Created on: 03 Jun 2011
 *      Author: Nicola Mori
 */

#include "montecarlo/services/GGSRootFileService.h"
#include "geometry/pluginmanagers/GGSGeoPluginManager.h"
#include "montecarlo/GGSRunManager.h"
#include "montecarlo/dataobjs/GGSTGeoParams.h"
#include "montecarlo/dataobjs/GGSTHitDetInfo.h"
#include "montecarlo/dataobjs/GGSTHitVolInfo.h"
#include "montecarlo/dataobjs/GGSTParameters.h"
#include "utils/GGSSmartLog.h"

#include "G4LogicalVolume.hh"
#include "G4SystemOfUnits.hh"

#include "TTree.h"

#include <sstream>

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

GGSRootFileService::FileInfo::FileInfo()
    : absBaseName(""), filePtr(NULL), defaultEventsTree(NULL), nRequests(0), detectorsMap(NULL), detectorsArray(NULL) {}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

GGSRootFileService &GGSRootFileService::GetInstance() {
  static GGSRootFileService instance;
  return instance;
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

GGSRootFileService::GGSRootFileService() : _suffix(""), _defaultOutBase("GGSRootOutput"), _currVolStorageFile(NULL) {}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

GGSRootFileService::~GGSRootFileService() {
  for (FileInfoContainer::iterator iter = _files.begin(); iter != _files.end(); iter++)
    if (iter->filePtr != NULL) {
      iter->filePtr->Close();
      delete iter->filePtr;
    }
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "G4Run.hh"
TFile *GGSRootFileService::GetFileForThisRun(const path &baseName, const G4Run *run) {

  path absBaseName;

  // Set default if no base name has been provided
  if (baseName.string() == "") {
    absBaseName = _defaultOutBase;
  } else
    absBaseName = baseName;

  // Retrieve absolute path
  absBaseName = _GetAbsolutePath(absBaseName);

  // check if the file is already opened
  for (FileInfoContainer::iterator iter = _files.begin(); iter != _files.end(); iter++)
    if (absBaseName == iter->absBaseName) {
      iter->nRequests++;
      return iter->filePtr;
    }

  // Create new file info struct
  FileInfo fInfo;
  fInfo.nRequests = 1;
  fInfo.absBaseName = absBaseName;

  // Strip extension and append suffix and extension
  TString absFileName = _AppendSuffixAndExt(absBaseName, run);

  // Open file and store it
  fInfo.filePtr = new TFile(absFileName, "RECREATE");
  _files.push_back(fInfo);

  // Write simulation informations
  _simInfo.Write();

  // Write geometry informations
  GGSVGeometryConstruction *geoCons = GGSGeoPluginManager::GetInstance().GetGeoConstruction();
  GGSTGeoParams geoParams;
  geoParams.SetIntGeoParams(geoCons->GetIntParameters());
  geoParams.SetBoolGeoParams(geoCons->GetBoolParameters());
  geoParams.SetRealGeoParams(geoCons->GetRealParameters());
  geoParams.SetStringGeoParams(geoCons->GetStringParameters());
  geoParams.SetVectIntGeoParams(geoCons->GetVectIntParameters());
  geoParams.SetVectBoolGeoParams(geoCons->GetVectBoolParameters());
  geoParams.SetVectRealGeoParams(geoCons->GetVectRealParameters());
  geoParams.SetVectStringGeoParams(geoCons->GetVectStringParameters());
  geoParams.Write("GGSGeoParams");

  // Write generator informations
  // These are written only for GGS generators that export a parameter object containing the "generator" name setting.
  auto *genAction =
      dynamic_cast<const GGSGeneratorAction *>(GGSRunManager::GetRunManager()->GetUserPrimaryGeneratorAction());
  if (genAction) {
    GGSTParameters generatorParams(genAction->GetParameters());
    try {
      generatorParams.GetParam<std::string>("generator"); // Throws if "generator" parameter is not present
      generatorParams.Write("GGSGenParams");
    } catch (...) {
    }
  }
  // Set storage of detector to null values to force proper initialization in StoreVolume
  _currVolStorageFile = NULL;
  _currDetVolName = "";

  return fInfo.filePtr;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void GGSRootFileService::CloseFileForThisRun(const path &baseName) {

  path absBaseName;

  // Set default if no base name has been provided
  if (baseName.string() == "") {
    absBaseName = _defaultOutBase;
  } else
    absBaseName = baseName;
  absBaseName = _GetAbsolutePath(absBaseName);

  for (FileInfoContainer::iterator iter = _files.begin(); iter != _files.end(); iter++)
    if (iter->absBaseName == absBaseName) {
      if (iter->nRequests > 1) {
        iter->nRequests--;
        return;
      } else {
        iter->filePtr->cd();
        if (iter->defaultEventsTree)
          iter->defaultEventsTree->Write();
        if (iter->detectorsArray)
          iter->detectorsArray->Write(
              "GGSHitDetInfo",
              TObject::kSingleKey); //  Automatically set iter->detectorsArray ownership to iter->filePtr
        iter->filePtr->Close();
        delete iter->filePtr;
        delete iter->detectorsMap;
        iter = _files.erase(iter);
      }
    }
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

TTree *GGSRootFileService::GetDefaultTree(TFile *file) {

  for (FileInfoContainer::iterator iter = _files.begin(); iter != _files.end(); iter++) {
    if (file == iter->filePtr) {
      if (iter->defaultEventsTree == NULL) {
        iter->filePtr->cd();
        iter->defaultEventsTree = new TTree("GGSEventsTree", "GGS events tree. Info: ");
      }
      return iter->defaultEventsTree;
    }
  }
  return NULL;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

path GGSRootFileService::_GetAbsolutePath(const path &baseName) {

  path absBaseName(baseName);

  absBaseName = absolute(baseName);

  return absBaseName;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

TString GGSRootFileService::_AppendSuffixAndExt(const path &baseName, const G4Run *run) {

  int runID = run->GetRunID();

  path fileName(baseName);
  fileName.replace_extension(); // Remove the extension
  TString newFileName = fileName.string();

  // if runID > 0 append _RunID
  if (runID > 0) {
    newFileName += "_Run";
    newFileName.Form(newFileName + "%i", runID);
  }

  // Append suffix
  if (_suffix != "")
    newFileName = newFileName + "_" + _suffix;

  // Append extension
  newFileName += ".root";

  return newFileName;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void GGSRootFileService::_FillDefaultEventsTrees() {
  for (FileInfoContainer::iterator iter = _files.begin(); iter != _files.end(); iter++) {
    if (iter->filePtr && iter->defaultEventsTree) {
      iter->filePtr->cd();
      iter->defaultEventsTree->Fill();
    }
  }
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

int GGSRootFileService::StoreVolume(TFile *filePtr, const std::string &detector, const G4VPhysicalVolume *volume,
                                    const G4ThreeVector &position, int id) {
  static const std::string routineName("GGSRootFileService::StoreVolume");

  std::string detVolName = detector.substr(0, detector.find_first_of('.'));

  // 1. Retrieve file info structure
  bool fileChanged = false;
  if (filePtr != _currVolStorageFile) {
    fileChanged = true;
    _currVolStorageFile = filePtr;
    for (_currVolStorageFileInfo = _files.begin(); _currVolStorageFileInfo != _files.end(); _currVolStorageFileInfo++) {
      if (_currVolStorageFileInfo->filePtr == filePtr)
        break;
    }
    if (_currVolStorageFileInfo == _files.end()) {
      COUT(WARNING) << "The requested file is not managed by GGSRootFileService. Volume informations will not be saved."
                    << ENDL;
      return -1;
    }
  }

  // 2. Find the detector in the index
  static HitDetMap::iterator currDetector;
  static GGSTHitDetInfo *currPersDetector = NULL;
  if (fileChanged || _currDetVolName != detVolName) {
    if (_currVolStorageFileInfo->detectorsMap == NULL)
      _currVolStorageFileInfo->detectorsMap = new HitDetMap;
    _currDetVolName = detVolName;
    int currPersDetectorIndex = 0;
    for (currDetector = _currVolStorageFileInfo->detectorsMap->begin();
         currDetector != _currVolStorageFileInfo->detectorsMap->end(); currDetector++, currPersDetectorIndex++) {
      if (currDetector->first == _currDetVolName)
        break;
    }
    if (currDetector == _currVolStorageFileInfo->detectorsMap->end()) {
      // 2.1 Detector not found. Create it both in index and persistence array
      _currVolStorageFileInfo->detectorsMap->push_back(std::pair<std::string, HitVolMap>(_currDetVolName, HitVolMap()));
      currDetector = _currVolStorageFileInfo->detectorsMap->end(); // Past-the-end iterator
      currDetector--;                                              // Set the iterator to the last element
      if (_currVolStorageFileInfo->detectorsArray == NULL) {
        _currVolStorageFileInfo->detectorsArray = new TClonesArray("GGSTHitDetInfo");
      }
      currPersDetectorIndex = _currVolStorageFileInfo->detectorsArray->GetEntries(); // Redundant
      currPersDetector = new ((*(_currVolStorageFileInfo->detectorsArray))[currPersDetectorIndex]) GGSTHitDetInfo;
      currPersDetector->detectorName = _currDetVolName.data();
    } else {
      // 2.2 Detector found. Set the current persistent detector
      currPersDetector = (GGSTHitDetInfo *)(_currVolStorageFileInfo->detectorsArray->At(currPersDetectorIndex));
    }
  }

  // 3. Insert the volume in the index and in the persistency structure

  // Append volume ID to volume name to avoid placing all the replicated volumes at the same place in the index
  static std::stringstream ss;
  ss.str("");
  ss << id;
  std::string volAndId = volume->GetName() + ss.str();

  std::pair<VolumeKey, G4int> insertValue(VolumeKey(volAndId, volume, position),
                                          currPersDetector->volumes.GetEntries());
  std::pair<HitVolMap::iterator, bool> insertResult;
  insertResult = currDetector->second.insert(insertValue); // TODO: using emplace instead of insert would improve speed?
  if (insertResult.second) {
    // 3.1 No hash collision: this is a new volume. Add it also to the volumes array in detector object...
    GGSTHitVolInfo *volInfo = new ((currPersDetector->volumes)[currPersDetector->volumes.GetEntries()]) GGSTHitVolInfo;
    volInfo->volumeName = volume->GetName().data();
    for (int i = 0; i < 3; i++) {
      volInfo->volumePos[i] = position[i] / cm;
    }
    volInfo->id = id;
    // ... then return its position in the array
    return currPersDetector->volumes.GetEntries() - 1;
  } else {
    // 3.2 A hash collision has happened: the volume is already present. So take its position from the index
    return insertResult.first->second;
  }
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......