GGSHitsAction Class Reference

An action which saves hits in sensitive volumes on ROOT file. More...

#include <GGSHitsAction.h>

Inheritance diagram for GGSHitsAction:

Collaboration diagram for GGSHitsAction:

Public Member Functions
	GGSHitsAction ()
	Constructor.
	~GGSHitsAction ()
	Destructor.
void	BeginOfRunAction (const G4Run *run)
	Actions executed at beginning of run. More...
void	EndOfRunAction (const G4Run *run)
	Actions executed at end of run. More...
void	BeginOfEventAction (const G4Event *event)
	Actions executed at beginning of event. More...
void	EndOfEventAction (const G4Event *event)
	Actions executed at end of event. More...
void	SetOutputFileBase (const std::string &outFileBase)
	Sets the output file base name. More...
void	SetOutputTreeName (const std::string &outTreeName)
	Sets the output tree name. More...
void	SetOutputIntHitClass (std::string detectorName, const std::string &className)
	Sets the name of the integrated hit class to be used for output of a given detector. More...
void	SetOutputPartHitClass (std::string detectorName, const std::string &className)
	Sets the name of the particle hit class to be used for output. More...
void	SetOutputPosHitClass (std::string detectorName, const std::string &className)
	Sets the name of the position hit class to be used for output. More...
void	SetHitThreshold (const std::string &detectorName, const std::string &hitType, const std::string &valueStr, const std::string &unit)
	Sets a lower energy deposit threshold for hits. More...
Public Member Functions inherited from GGSUserAction
	GGSUserAction ()
	Constructor. More...
virtual	~GGSUserAction ()
	Destructor.
G4ClassificationOfNewTrack	ClassifyNewTrack (const G4Track *)
	Override of the ClassifyNewTrack method. More...

Friends
class	GGSHitsActionPrivateMessenger

Detailed Description

An action which saves hits in sensitive volumes on ROOT file.

This action reads the hits in sensitive volumes, as defined by GGS sensitive logical volumes naming convention. It then converts the hits in ROOT format and saves them on a ROOT file. For each sensitive logical volume, a branch will be created whose name is equal to that of logical volume. This branch is made by a TClonesArray containing all the hits in physical sensitive volumes associated to the logical. This action can handle different kind of output hit classes derived from standard GGS ones (GGSTIntHitBase, GGSTPartHitBase and GGSTPosHitBase).The actual classes to be used is specified by calling the related setter methods (e.g. SetOutputIntHitClass); the hit classes must have a Root dictionary. If no class is specified then the standard GGS hit classes (GGSTIntHit, GGSTPartHit, GGSTPosHit) will be used.

See Also

GGSIntHit

GGSPartHit

GGSPosHit

Definition at line 64 of file GGSHitsAction.h.

Member Function Documentation

void GGSHitsAction::BeginOfEventAction

(

const G4Event *

event

)

Actions executed at beginning of event.

Method executed at the beginning of each event. Hits buffers are cleared.

Parameters

event

Pointer to a G4Event

Definition at line 254 of file GGSHitsAction.cpp.

                                                      {

  // Give back arrays to the TCLonesArray services.
  _GiveBackAllArrays();

  for (auto &handleElem : _outputHandles) {
    OutputHandle &handle = handleElem.second;
    if (handle.partHitArray) {

      // Clear the TObjarray persistence structures. These contains TClonesArrays which
      // must not be deleted since they are owned and managed by the GGSTClonesArrayService...

      // ... so set the TObjArray to not own them...
      handle.partHitArray->SetOwner(false);
      // ... then clean to make size=0 without deleting the contained objects...
      handle.partHitArray->Clear();
      // ... and then restore ownership to avoid memory leak on readout of ROOT file.
      // (see http://root.cern.ch/root/roottalk/roottalk02/0444.html)
      handle.partHitArray->SetOwner(true);
    }

    if (handle.posHitArray) {
      handle.posHitArray->SetOwner(false);
      handle.posHitArray->Clear();
      handle.posHitArray->SetOwner(true);
    }
  }
}

void GGSHitsAction::BeginOfRunAction

(

const G4Run *

run

)

Actions executed at beginning of run.

Method executed at the beginning of each run. A new ROOT file is created for each run.

Parameters

run	Pointer to the G4Run

Definition at line 74 of file GGSHitsAction.cpp.

                                                     {
  const std::string routineName("GGSHitsAction::BeginOfRunAction");

  // -----------------
  // Open the run file
  // -----------------
  _outFile = GGSRootFileService::GetInstance().GetFileForThisRun(_outFileBase, run);
  if (!_outFile || (_outFile && _outFile->IsZombie())) {
    GGSCOUT(ERROR) << "Cannot open ROOT file " << _outFileBase << "for run" << run->GetRunID() << GGSENDL;
    throw -1;
  }

  // -----------------
  // Create a new tree
  // -----------------
  if (_outTreeName == "") {
    _outTree = GGSRootFileService::GetInstance().GetDefaultTree(_outFile);
    _outTree->SetTitle(TString(_outTree->GetTitle()) + "Hits ");
  } else
    _outTree = new TTree(_outTreeName, "GGS hits");

  // Loop over logical volumes...
  G4LogicalVolumeStore *logicalVS = G4LogicalVolumeStore::GetInstance();
  if (!logicalVS)
    return;

  std::vector<G4VSensitiveDetector *> processedSDs;
  for (unsigned int i = 0; i < logicalVS->size(); i++) {
    G4LogicalVolume *logicalV = logicalVS->at(i);

    // ...and select sensitive volumes with GGS hits SD
    std::vector<G4VSensitiveDetector *> sDets(0);
    G4VSensitiveDetector *sd = logicalV->GetSensitiveDetector();
    if (sd) {
      GGSMultiSensitiveDetector *msd = dynamic_cast<GGSMultiSensitiveDetector *>(sd);
      if (msd) {
        sDets = msd->GetListOfSensitiveDetectors();
      } else {
        sDets.push_back(sd);
      }
    }

    for (auto &sdPtr : sDets) {
      if (std::find(processedSDs.begin(), processedSDs.end(), sdPtr) != processedSDs.end()) {
        continue;
      }
      auto sdName = sdPtr->GetName();
      auto ggssdNamePos = sdName.find(".GGSIntHitSD", 0);
      if (ggssdNamePos != std::string::npos) {
        GGSIntHitSD *intHitSD = dynamic_cast<GGSIntHitSD *>(sdPtr);
        if (intHitSD && intHitSD->isActive()) {
          std::string sdNameReduced = sdName;
          sdNameReduced.erase(ggssdNamePos, 12); // Remove ".GGSIntHitSD"

          std::unordered_map<std::string, OutputHandle>::iterator handleIter = _outputHandles.find(sdName);
          if (handleIter == _outputHandles.end()) {
            // The handle for this detector does not exist yet, so create a standard one

            auto insertResult = _outputHandles.insert(std::pair<std::string, OutputHandle>{sdName, OutputHandle()});
            if (insertResult.second) {
              handleIter = insertResult.first;
            }
          }

          // Create arrays et al.
          auto &newHandle = handleIter->second;
          newHandle.intHitArray = std::unique_ptr<TClonesArray>(new TClonesArray(
              newHandle.intHitClassName.c_str())); // Create it here and Clear("C") it at beginning of event
          // Create a branch to store hits in the output file
          _outTree->Branch((sdNameReduced + "-GGSIntHits").data(), "TClonesArray", newHandle.intHitArray.get(), 64000);

          // Set output threshold
          {
            auto thresholdIter = std::find_if(
                _thresholds.begin(), _thresholds.end(), [&sdNameReduced](const ThresholdStruct &threshold) {
                  return (threshold.hitType == 0) && (sdNameReduced == threshold.detectorName);
                });
            if (thresholdIter != _thresholds.end()) {
              newHandle.intHitThreshold = thresholdIter->value;
            }
          }
          if (intHitSD->GetPartHitsStorage()) {
            newHandle.partHitArray = std::unique_ptr<TObjArray>(new TObjArray());
            // Set ownership to avoid memory leaks when reading back the TObjArray.
            // http://root.cern.ch/root/roottalk/roottalk02/0444.html
            newHandle.partHitArray->SetOwner(true);
            _outTree->Branch((sdNameReduced + "-GGSPartHits").data(), "TObjArray", newHandle.partHitArray.get(), 64000);
            newHandle.partHitCAService =
                std::unique_ptr<GGSTClonesArrayService>(new GGSTClonesArrayService(newHandle.partHitClassName.c_str()));
            // Set output threshold
            auto thresholdIter = std::find_if(
                _thresholds.begin(), _thresholds.end(), [&sdNameReduced](const ThresholdStruct &threshold) {
                  return (threshold.hitType == 1) && (sdNameReduced == threshold.detectorName);
                });
            if (thresholdIter != _thresholds.end()) {
              newHandle.partHitThreshold = thresholdIter->value;
            }
          }
          if (intHitSD->GetPosHitsStorage()) {
            newHandle.posHitArray = std::unique_ptr<TObjArray>(new TObjArray());
            newHandle.posHitArray->SetOwner(true);
            _outTree->Branch((sdNameReduced + "-GGSPosHits").data(), "TObjArray", newHandle.posHitArray.get(), 64000);
            newHandle.posHitCAService =
                std::unique_ptr<GGSTClonesArrayService>(new GGSTClonesArrayService(newHandle.posHitClassName.c_str()));
            // Set output threshold
            auto thresholdIter = std::find_if(
                _thresholds.begin(), _thresholds.end(), [&sdNameReduced](const ThresholdStruct &threshold) {
                  return (threshold.hitType == 2) && (sdNameReduced == threshold.detectorName);
                });
            if (thresholdIter != _thresholds.end()) {
              newHandle.posHitThreshold = thresholdIter->value;
            }
          }

          // Save the time bins as run products
          if (GGSRootFileService::GetInstance().GetOutputMode() == GGSRootFileService::OutputMode::MULTIFILE ||
              G4Threading::G4GetThreadId() <= 0) {
            unsigned int nBins = intHitSD->GetTimeBins().size();
            TArrayF timeBins(nBins);
            for (unsigned int iBin = 0; iBin < nBins; iBin++)
              timeBins[iBin] = intHitSD->GetTimeBins()[iBin];
            _outFile->WriteObject(&timeBins, (sdNameReduced + "-GGSIntHits-TimeBins").data());

            // Save the hit thresholds
            std::vector<float> detThresholds{(float)(newHandle.intHitThreshold / GeV),
                                             (float)(newHandle.partHitThreshold / GeV),
                                             (float)(newHandle.posHitThreshold / GeV)};
            _outFile->WriteObject(&detThresholds, (sdNameReduced + "-GGSIntHits-Thresholds").data());
          }

        } // Check on active volume

      } // Check on sensitive volume name
      processedSDs.push_back(sdPtr);
    } // End loop on sensitive detectors

  } // End loop on logical volumes
}

void GGSHitsAction::EndOfEventAction

(

const G4Event *

event

)

Actions executed at end of event.

Method executed at the end of each event. Hits are read, converted to ROOT format and filled into output tree.

Parameters

event

Pointer to a G4Event

Definition at line 285 of file GGSHitsAction.cpp.

                                                         {

  //-------------------------------
  // Set hits info
  //-------------------------------

  // Retrieve SD manager
  G4SDManager *SDman = G4SDManager::GetSDMpointer();

  // Retrieve event hits
  G4HCofThisEvent *evHitColl = event->GetHCofThisEvent();

  // Loop on detectors
  for (auto &handleElem : _outputHandles) {

    const std::string &detectorName = handleElem.first;
    OutputHandle &handle = handleElem.second;
    // Get hits of current detector
    G4int hitsCollID = SDman->GetCollectionID(detectorName);

    // Convert hits to persistent data format
    handle.intHitArray->Clear(); // Set entries to 0 but do not clear elements (will be cleared in _Convert)
    _Convert(*((GGSIntHitsCollection *)(evHitColl->GetHC(hitsCollID))), *(handle.intHitArray), handle, detectorName);

    // Link particle and position hit collections to persistence array
    for (int iIntHit = 0; iIntHit < handle.intHitArray->GetEntries(); ++iIntHit) {
      GGSTIntHitBase *tIntHit = (GGSTIntHitBase *)(handle.intHitArray->At(iIntHit));
      if (tIntHit->_partHits) {
        handle.partHitArray->AddAtAndExpand(tIntHit->_partHits, handle.partHitArray->GetEntries());
        if (tIntHit->_partHits->GetEntries() > 0 && ((GGSTPartHitBase *)(tIntHit->_partHits->At(0)))->_posHits) {
          static TClonesArray *partHits = nullptr;
          static Int_t nPartHits = 0;
          partHits = tIntHit->_partHits;
          nPartHits = partHits->GetEntries();
          for (Int_t iPartHit = 0; iPartHit < nPartHits; iPartHit++) {
            GGSTPartHitBase *partHit = (GGSTPartHitBase *)(tIntHit->_partHits->At(iPartHit));
            partHit->_posHitIndex = handle.posHitArray->GetEntries(); // Save reference to pos hit array
            handle.posHitArray->AddAtAndExpand(partHit->_posHits, handle.posHitArray->GetEntries());
          }
        }
      }
    }
  }

  //-----------
  // Fill tree
  //-----------
  if (_outTreeName != "") {
    // Call fill only if not using the default tree
    _outFile->cd();
    _outTree->Fill();
  }
}

void GGSHitsAction::EndOfRunAction

(

const G4Run *

run

)

Actions executed at end of run.

Method executed at the end of each run. The ROOT tree is written to file and the file is closed.

Parameters

run	Pointer to a G4Run

Definition at line 215 of file GGSHitsAction.cpp.

                                                   {

  // -----------------
  // close output file
  // -----------------

  if (_outFile) {
    _outFile->cd();
    if (_outTreeName != "") {
      // Write the tree if we are not using the default one
      _outTree->Write();
    }

    GGSRootFileService::GetInstance().CloseFileForThisRun(_outFileBase,
                                                          run); // Will automatically save detector informations
    _outFile = nullptr;
    _outTree = nullptr; // The tree is delete by the file's destructor.
  } else {
    throw -1; // TODO handle errors
  }

  // -----------------
  // delete detectors
  // -----------------

  // Prepare for an eventual next run
  _GiveBackAllArrays();
  for (auto &handleElem : _outputHandles) {
    if (handleElem.second.posHitCAService) {
      handleElem.second.posHitCAService->DeleteAll();
    }
    if (handleElem.second.partHitCAService) {
      handleElem.second.partHitCAService->DeleteAll();
    }
  }
}

void GGSHitsAction::SetHitThreshold	(	const std::string &	detectorName,
		const std::string &	hitType,
		const std::string &	valueStr,
		const std::string &	unit
	)

Sets a lower energy deposit threshold for hits.

Hits with a total energy deposit lower than the corresponding threshold will not be saved in the output file. If a finer-grained hit is above its threshold then the parent hit will be saved regardless of its value. For example, if particle hits have threshold 1 and position hits 0.1, then a particle hit with a total release of 0.2 containing a single position hit with release 0.2 will be saved (since position hits are contained inside particle hits, so the particle hit is necessary to store the above-threshold position hit). The energy deposit of parent hits will be the total one computed by the MC simulation, irrespective of the total energy deposit of above-threshold children hits. Extending the example above, if the particle hit have a release of 1.05 due to five position hits with release 0.2 and one position hit with release 0.05 then only the five position hits with release 0.2 will be saved; but the energy deposit of the particle hit will still be 1.05 although the sum of stored position hits is 1. The above consideration apply also to integrated-particle hits hierarchy. In case custom hits are used, which eventually customize the meaning of energy deposit, then the thresholds will be relative to the customized deposits.

Parameters

detectorName	The name of the detector.
hitType	The type of the hit (integrated, particle or position).
valueStr	The threshold value.
unit	The unit of the threshold value.

Definition at line 516 of file GGSHitsAction.cpp.

                                                                                      {

  const std::string routineName("GGSHitsAction::SetHitThreshold");

  int hitTypeNo = -1;
  if (hitType == "integrated") {
    hitTypeNo = 0;
  } else if (hitType == "particle") {
    hitTypeNo = 1;
  } else if (hitType == "position") {
    hitTypeNo = 2;
  } else {
    throw std::runtime_error(routineName + ": Unknown hit type " + hitType + " for detector " + detectorName);
  }

  float value = 0.;
  try {
    value = std::stof(valueStr);
  } catch (std::invalid_argument &exc) {
    GGSCOUT(ERROR) << "Invalid threshold for " << hitType << " hits of detector " << detectorName << GGSENDL;
    throw exc;
  } catch (std::out_of_range &exc) {
    GGSCOUT(ERROR) << "Invalid value for threshold of " << hitType << " hits of detector " << detectorName << GGSENDL;
    throw exc;
  }

  auto &unitsTable = G4UnitDefinition::GetUnitsTable();
  auto unitCatIter = std::find_if(unitsTable.begin(), unitsTable.end(),
                                  [](const G4UnitsCategory *unitCat) { return (unitCat->GetName() == "Energy"); });
  if (unitCatIter == unitsTable.end()) {
    throw std::runtime_error(routineName + ": energy unit category not found");
  }

  auto &unitsList = (*unitCatIter)->GetUnitsList();
  auto unitIter = std::find_if(unitsList.begin(), unitsList.end(), [&unit](const G4UnitDefinition *unitDef) {
    return unit == unitDef->GetSymbol().c_str();
  });
  if (unitIter == unitsList.end()) {
    throw std::runtime_error(routineName + ": unit " + unit + " not found");
  }

  auto thresholdIter = std::find_if(
      _thresholds.begin(), _thresholds.end(), [&detectorName, hitTypeNo](const ThresholdStruct &threshold) {
        return (threshold.hitType == hitTypeNo) && (threshold.detectorName == detectorName);
      });
  if (thresholdIter == _thresholds.end()) {
    _thresholds.push_back(ThresholdStruct{detectorName, hitTypeNo, -1.});
    thresholdIter = _thresholds.end() - 1;
  }

  auto *sd =
      dynamic_cast<GGSIntHitSD *>(G4SDManager::GetSDMpointer()->FindSensitiveDetector(detectorName + ".GGSIntHitSD"));
  if (!sd) {
    throw std::invalid_argument(routineName + ": detector " + detectorName + " not found");
  }

  switch (hitTypeNo) {
  case 0:
    thresholdIter->value = value * (*unitIter)->GetValue();
    sd->SetIntHitThreshold(thresholdIter->value);
    break;
  case 1:
    thresholdIter->value = value * (*unitIter)->GetValue();
    sd->SetPartHitThreshold(thresholdIter->value);
    break;
  case 2:
    thresholdIter->value = value * (*unitIter)->GetValue();
    sd->SetPosHitThreshold(thresholdIter->value);
    break;
  }
}

void GGSHitsAction::SetOutputFileBase ( const std::string &  outFileBase )

inline

Sets the output file base name.

The file base name can be with or without extension (.root will be automatically used as extension). For each run, the run number will be appended to the base name before the .root extension. If no value is provided the file base name will fallback to the default value set in GSRootFileservice.

Parameters

outFileBase

The output file base name.

See Also

GGSRootFileservice

Definition at line 117 of file GGSHitsAction.h.

{ _outFileBase = outFileBase; }

void GGSHitsAction::SetOutputIntHitClass	(	std::string	detectorName,
		const std::string &	className
	)

Sets the name of the integrated hit class to be used for output of a given detector.

Parameters

detectorName	name of the detector
className	name of the output class for integrated hits

Definition at line 450 of file GGSHitsAction.cpp.

                                                                                           {
  const std::string routineName("[GGSHitsAction::SetOutputIntHitClass] ");

  // Check if the given class has a dictionary
  TClass *intHitClass = TClass::GetClass(className.c_str());
  if (!intHitClass) {
    throw std::runtime_error(routineName + "Dictionary for class " + className + " not found");
  }

  // Check if the given detector already has a specified output class
  detectorName.insert(detectorName.find_first_of('.'), ".GGSIntHitSD");
  auto handleIter = _outputHandles.find(detectorName);
  if (handleIter == _outputHandles.end()) {
    handleIter = _outputHandles.insert(std::pair<std::string, OutputHandle>{detectorName, OutputHandle()}).first;
  }

  // Set output facilities
  handleIter->second.intHitClassName = className;
}

void GGSHitsAction::SetOutputPartHitClass	(	std::string	detectorName,
		const std::string &	className
	)

Sets the name of the particle hit class to be used for output.

Parameters

detectorName	name of the detector
className	name of the output class for particle hits

Definition at line 472 of file GGSHitsAction.cpp.

                                                                                            {
  const std::string routineName("[GGSHitsAction::SetOutputPartHitClass] ");

  // Check if the given class has a dictionary
  TClass *partHitClass = TClass::GetClass(className.c_str());
  if (!partHitClass) {
    throw std::runtime_error(routineName + "Dictionary for class " + className + " not found");
  }

  // Check if the given detector already has a specified output class
  detectorName.insert(detectorName.find_first_of('.'), ".GGSIntHitSD");
  auto handleIter = _outputHandles.find(detectorName);
  if (handleIter == _outputHandles.end()) {
    handleIter = _outputHandles.insert(std::pair<std::string, OutputHandle>{detectorName, OutputHandle()}).first;
  }

  // Set output facilities
  handleIter->second.partHitClassName = className;
}

void GGSHitsAction::SetOutputPosHitClass	(	std::string	detectorName,
		const std::string &	className
	)

Sets the name of the position hit class to be used for output.

Parameters

detectorName	name of the detector
className	name of the output class for position hits

Definition at line 494 of file GGSHitsAction.cpp.

                                                                                           {
  const std::string routineName("[GGSHitsAction::SetOutputPosHitClass] ");

  // Check if the given class has a dictionary
  TClass *posHitClass = TClass::GetClass(className.c_str());
  if (!posHitClass) {
    throw std::runtime_error(routineName + "Dictionary for class " + className + " not found");
  }

  // Check if the given detector already has a specified output class
  detectorName.insert(detectorName.find_first_of('.'), ".GGSIntHitSD");
  auto handleIter = _outputHandles.find(detectorName);
  if (handleIter == _outputHandles.end()) {
    handleIter = _outputHandles.insert(std::pair<std::string, OutputHandle>{detectorName, OutputHandle()}).first;
  }

  // Set output facilities
  handleIter->second.posHitClassName = className;
}

void GGSHitsAction::SetOutputTreeName ( const std::string &  outTreeName )

inline

Sets the output tree name.

This name will be used for the TTree object where th hits for each event will be stored. By default, this value is set to "Hits" in constructor.

Parameters

outTreeName

The output tree name.

Definition at line 126 of file GGSHitsAction.h.

{ _outTreeName = outTreeName; }

---

The documentation for this class was generated from the following files:

• GGSHitsAction.h
• GGSHitsAction.cpp