Commit 9be914d9 authored by j.a.galan.81@gmail.com's avatar j.a.galan.81@gmail.com
Browse files

First commit

parents
set(contents)
if (NOT DEFINED RESTLIB_DETECTOR OR NOT ${RESTLIB_DETECTOR} MATCHES "ON")
set(excludes TRestGeant4ToHitsProcess)
endif (NOT DEFINED RESTLIB_DETECTOR OR NOT ${RESTLIB_DETECTOR} MATCHES "ON")
COMPILEDIR(RestGeant4)
# The REST Geant4 Library
TOBE written
This file should include basic information to understand the main use trend of the library. It should include a brief description of this class. What is intended for, features, functionalities, connection to other REST libraries, etc. And a description of the main use of this library.
NOTE: The use of this library should be documented through examples and doxygen class description.
## Getting Started
These instructions will get you a copy of the project up and running on your local machine.
TODO: Add a description and explain that the recommended way is to install it as a git submodule in RestFramework project.
Additionaly you may download a copy of this library cloning it using the corresponding git command.
```
git clone git@lfna.unizar.es:rest/RestGeant4Lib.git
```
### Prerequisites
The only mandatory prerequisite is that REST and ROOT6 are installed in your system. Details on the installation of ROOT will be found at the [ROOT's official site](root.cern.ch).
To avoid problems during the installation it is recommended that REST libraries will be compiled using the same version of g++ compiler used to compile ROOT and REST.
### Installing
After ROOT6 and REST haven been installed in the system, the compilation of this REST library should be straight forward.
Before starting the installation process make sure you are running the desired ROOT version and binary,
```
root-config --version
which root
```
and the corresponding REST version and binary,
```
rest-config --version
which restRoot
```
The cmake compilation environment will recognize those and use them to link to your ROOT installation.
Then, go to the root directory of your local REST repository, lets name it here `REST_SOURCE_PATH` and execute the following commands.
```
cd THIS_LIB_SRC_PATH
mkdir build
cd build
cmake ../
make -j4 install
```
After all the compilation and installation process ends, you will end up with a new installed REST library inside `REST_PATH/lib/`.
Three other directories, data, macros and examples - common to any REST library - will be installed.
### Basic tests of this library
You can test now that REST libraries are loaded without errors inside the ROOT environment using the `restRoot` command. After launching `restRoot` you should see a message on screen similar to the following one, where the recently compiled library will be found in the list of loaded libraries.
```
------------------------------------------------------------
| Welcome to ROOT 6.14/06 http://root.cern.ch |
| (c) 1995-2018, The ROOT Team |
| Built for macosx64 |
| From tags/v6-14-06@v6-14-06, Dec 11 2018, 12:58:00 |
| Try '.help', '.demo', '.license', '.credits', '.quit'/'.q' |
------------------------------------------------------------
Loading library : libRestGeant4.dylib
Loading library : libRestFramework.dylib
Loading library : libRestDetector.dylib
```
### A brief description of the examples provided
### Compilation options
TO define any library compilation options if necessary.
Different options can be passed to the `cmake` command to personalize the REST installation. The following options are available.
* **INSTALL_PREFIX** (Default: REST_PATH) : Allows to define the destination of the final REST install directory.
* **REST_DUMMY** (Default: ON) : A dummy option that might be defined in the CMakeLists.
## License
This project is licensed under the GNU License - see the [LICENSE](https://lfna.unizar.es/rest-development/REST_v2/blob/master/LICENCE) file for details
## Acknowledgments
`TOBE written`
* Hat tip to anyone whose code was used
* Inspiration
* etc
/*************************************************************************
* This file is part of the REST software framework. *
* *
* Copyright (C) 2016 GIFNA/TREX (University of Zaragoza) *
* For more information see http://gifna.unizar.es/trex *
* *
* REST is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* REST is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have a copy of the GNU General Public License along with *
* REST in $REST_PATH/LICENSE. *
* If not, see http://www.gnu.org/licenses/. *
* For the list of contributors see $REST_PATH/CREDITS. *
*************************************************************************/
#ifndef RestCore_TRestGeant4AnalysisProcess
#define RestCore_TRestGeant4AnalysisProcess
#include <TRestGeant4Event.h>
#include <TRestGeant4Metadata.h>
#include "TRestEventProcess.h"
//! A pure analysis process to extract information from a TRestGeant4Event
class TRestGeant4AnalysisProcess : public TRestEventProcess {
private:
/// A pointer to the specific TRestGeant4Event input
TRestGeant4Event* fInputG4Event; //!
/// A pointer to the specific TRestGeant4Event output
TRestGeant4Event* fOutputG4Event; //!
/// A pointer to the simulation metadata information accessible to TRestRun
TRestGeant4Metadata* fG4Metadata; //!
/// It stores the name of observables `xxxVolumeEDep` related to the energy deposition in volume `xxx`.
std::vector<std::string> fEnergyInObservables; //!
/// A vector storing the active volume ids of observables `xxxVolumeEDep`.
std::vector<Int_t> fVolumeID; //!
/// It stores the name of observables (xxxMeanPosX,Y,Z) related to mean hits position in volume `xxx`.
std::vector<std::string> fMeanPosObservables; //!
/// A vector storing the active volume ids corresponding mean position observable `xxxMeanPosX,Y,Z`.
std::vector<Int_t> fVolumeID2; //!
/// A vector storing the direction X,Y or Z from corresponding mean position observable `xxxMeanPosX,Y,Z`.
std::vector<std::string> fDirID; //!
/// A vector storing the name of observables related to processes in a particular active volume.
std::vector<std::string> fProcessObservables; //!
/// A vector storing the active volume ids corresponding process observable .
std::vector<Int_t> fVolumeID3; //!
/// A vector storing the name of processes.
std::vector<std::string> fProcessName; //!
/// A vector storing the observable name `xxxTracksCounter` for a given `xxx` particle.
std::vector<std::string> fTrackCounterObservables; //!
/// A vector storing the `xxx` particle name extracted from `xxxTracksCounter`.
std::vector<std::string> fParticleTrackCounter; //!
/// A vector storing the observable name `xxxTracksEDep` for a given `xxx` particle.
std::vector<std::string> fTracksEDepObservables; //!
/// A vector storing the `xxx` particle name extracted from `xxxTracksEDep`.
std::vector<std::string> fParticleTrackEdep; //!
// TODO these two thresholds should become OBSOLETE. We must use now the more generic
// way ApplyCut() with the <cut ...> definition at the RML process definition.
/// A low energy cut. Events below the threshold will be not further processed.
Double_t fLowEnergyCut;
/// A high energy cut. Events above the threshold will be not further processed.
Double_t fHighEnergyCut;
Bool_t fPerProcessSensitiveEnergy = false;
Bool_t fPerProcessSensitiveEnergyNorm = false;
void InitFromConfigFile();
void Initialize();
void LoadDefaultConfig();
protected:
// add here the members of your event process
public:
any GetInputEvent() { return fInputG4Event; }
any GetOutputEvent() { return fOutputG4Event; }
void InitProcess();
TRestEvent* ProcessEvent(TRestEvent* eventInput);
void EndProcess();
void LoadConfig(std::string cfgFilename, std::string name = "");
/// It prints out the process parameters stored in the metadata structure
void PrintMetadata() {
BeginPrintProcess();
metadata << "Low energy cut : " << fLowEnergyCut << " keV" << endl;
metadata << "High energy cut : " << fHighEnergyCut << " keV" << endl;
EndPrintProcess();
}
/// Returns a new instance of this class
TRestEventProcess* Maker() { return new TRestGeant4AnalysisProcess; }
/// Returns the name of this process
TString GetProcessName() { return (TString) "geant4Analysis"; }
TRestGeant4AnalysisProcess();
TRestGeant4AnalysisProcess(char* cfgFileName);
~TRestGeant4AnalysisProcess();
ClassDef(TRestGeant4AnalysisProcess, 2);
};
#endif
///______________________________________________________________________________
///______________________________________________________________________________
///______________________________________________________________________________
///
///
/// RESTSoft : Software for Rare Event Searches with TPCs
///
/// TRestGeant4BiasingVolume.h
///
/// A template class for copy/paste
///
/// jul 2015: First concept
/// J. Galan
///_______________________________________________________________________________
#ifndef RestCore_TRestGeant4BiasingVolume
#define RestCore_TRestGeant4BiasingVolume
#include <iostream>
#include "TObject.h"
#include "TVector3.h"
class TRestGeant4BiasingVolume : public TObject {
protected:
TVector3 fVolumePosition;
Double_t fVolumeSize;
TString fBiasingVolumeType;
Double_t fBiasingFactor;
TVector2 fEnergyRange;
TString fVolumeType;
public:
Double_t GetBiasingFactor() { return fBiasingFactor; }
Double_t GetBiasingVolumeSize() { return fVolumeSize; }
TString GetBiasingVolumeType() { return fVolumeType; }
TVector3 GetBiasingVolumePosition() { return fVolumePosition; }
TVector2 GetEnergyRange() { return fEnergyRange; }
Double_t GetMaxEnergy() { return fEnergyRange.Y(); }
Double_t GetMinEnergy() { return fEnergyRange.X(); }
void SetBiasingVolumeSize(Double_t size) { fVolumeSize = size; }
void SetBiasingVolumeType(TString type) { fVolumeType = type; }
void SetBiasingVolumePosition(TVector3 pos) { fVolumePosition = pos; }
void SetBiasingFactor(Double_t factor) { fBiasingFactor = factor; }
void SetEnergyRange(TVector2 eRange) { fEnergyRange = eRange; }
// Check if it is inside the sphere
Int_t isInside(Double_t x, Double_t y, Double_t z) {
if (fVolumeType == "virtualBox") {
if (x < fVolumeSize / 2. && x > -fVolumeSize / 2.)
if (y < fVolumeSize / 2. && y > -fVolumeSize / 2.)
if (z < fVolumeSize / 2. && z > -fVolumeSize / 2.) return 1;
}
if (fVolumeType == "virtualSphere") {
Double_t r2 = x * x + y * y + z * z;
if (r2 < fVolumeSize * fVolumeSize) return 1;
}
return 0;
}
void PrintBiasingVolume();
// Construtor
TRestGeant4BiasingVolume();
// Destructor
virtual ~TRestGeant4BiasingVolume();
ClassDef(TRestGeant4BiasingVolume, 2); // REST event superclass
};
#endif
///______________________________________________________________________________
///______________________________________________________________________________
///______________________________________________________________________________
///
///
/// RESTSoft : Software for Rare Event Searches with TPCs
///
/// TRestGeant4BlobAnalysisProcess.h
///
///_______________________________________________________________________________
#ifndef RestCore_TRestGeant4BlobAnalysisProcess
#define RestCore_TRestGeant4BlobAnalysisProcess
#include <TRestGeant4Event.h>
#include <TRestGeant4Metadata.h>
#include "TRestEventProcess.h"
class TRestGeant4BlobAnalysisProcess : public TRestEventProcess {
private:
#ifndef __CINT__
TRestGeant4Event* fG4Event; //!
TRestGeant4Metadata* fG4Metadata; //!
std::vector<std::string> fQ1_Observables; //!
std::vector<double> fQ1_Radius; //!
std::vector<std::string> fQ2_Observables; //!
std::vector<double> fQ2_Radius; //!
#endif
void InitFromConfigFile();
void Initialize();
void LoadDefaultConfig();
protected:
// add here the members of your event process
public:
any GetInputEvent() { return fG4Event; }
any GetOutputEvent() { return fG4Event; }
void InitProcess();
TRestEvent* ProcessEvent(TRestEvent* eventInput);
void EndProcess();
void LoadConfig(std::string cfgFilename, std::string name = "");
void PrintMetadata() {
BeginPrintProcess();
EndPrintProcess();
}
TString GetProcessName() { return (TString) "findG4BlobAnalysis"; }
// Constructor
TRestGeant4BlobAnalysisProcess();
TRestGeant4BlobAnalysisProcess(char* cfgFileName);
// Destructor
~TRestGeant4BlobAnalysisProcess();
ClassDef(TRestGeant4BlobAnalysisProcess, 1); // Template for a REST "event process" class inherited from
// TRestEventProcess
};
#endif
/*************************************************************************
* This file is part of the REST software framework. *
* *
* Copyright (C) 2016 GIFNA/TREX (University of Zaragoza) *
* For more information see http://gifna.unizar.es/trex *
* *
* REST is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or *
* (at your option) any later version. *
* *
* REST is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have a copy of the GNU General Public License along with *
* REST in $REST_PATH/LICENSE. *
* If not, see http://www.gnu.org/licenses/. *
* For the list of contributors see $REST_PATH/CREDITS. *
*************************************************************************/
#ifndef RestCore_TRestGeant4Event
#define RestCore_TRestGeant4Event
#include <TGraph.h>
#include <TGraph2D.h>
#include <TH1D.h>
#include <TH2F.h>
#include <TLegend.h>
#include <TMultiGraph.h>
#include <TObject.h>
#include <TRestEvent.h>
#include <TRestGeant4Track.h>
#include <TVector3.h>
#include <iostream>
#include <map>
/// An event class to store geant4 generated event information
class TRestGeant4Event : public TRestEvent {
private:
#ifndef __CINT__
Double_t fMinX, fMaxX; //!
Double_t fMinY, fMaxY; //!
Double_t fMinZ, fMaxZ; //!
Double_t fMinEnergy, fMaxEnergy; //!
#endif
void AddEnergyDepositToVolume(Int_t volID, Double_t eDep);
Bool_t PerProcessEnergyInitFlag = false;
std::map<string, Double_t> PerProcessEnergyInSensitive;
void inline InitializePerProcessEnergyInSensitive() {
PerProcessEnergyInitFlag = true;
PerProcessEnergyInSensitive["photoelectric"] = 0;
PerProcessEnergyInSensitive["compton"] = 0;
PerProcessEnergyInSensitive["electron_ionization"] = 0;
PerProcessEnergyInSensitive["ion_ionization"] = 0;
PerProcessEnergyInSensitive["alpha_ionization"] = 0;
PerProcessEnergyInSensitive["msc"] = 0;
PerProcessEnergyInSensitive["hadronic_ionization"] = 0;
PerProcessEnergyInSensitive["proton_ionization"] = 0;
PerProcessEnergyInSensitive["hadronic_elastic"] = 0;
PerProcessEnergyInSensitive["neutron_elastic"] = 0;
string volume_name;
string process_name;
TRestGeant4Track* track;
TRestGeant4Hits* hits;
Double_t energy;
for (Int_t track_id = 0; track_id < GetNumberOfTracks(); track_id++) {
track = GetTrack(track_id);
if (track->GetEnergyInVolume(0) == 0) {
continue;
}
hits = track->GetHits();
for (Int_t hit_id = 0; hit_id < hits->GetNumberOfHits(); hit_id++) {
if (hits->GetVolumeId(hit_id) != 0) {
continue;
}
process_name = (string)track->GetProcessName(hits->GetHitProcess(hit_id));
energy = hits->GetEnergy(hit_id);
if (process_name == "phot") {
PerProcessEnergyInSensitive["photoelectric"] += energy;
} else if (process_name == "compt") {
PerProcessEnergyInSensitive["compton"] += energy;
} else if (process_name == "eIoni" || process_name == "e-Step" || process_name == "e+Step") {
PerProcessEnergyInSensitive["electron_ionization"] += energy;
} else if (process_name == "ionIoni") {
PerProcessEnergyInSensitive["ion_ionization"] += energy;
if (track->GetParticleName() == "alpha") {
PerProcessEnergyInSensitive["alpha_ionization"] += energy;
}
} else if (process_name == "msc") {
PerProcessEnergyInSensitive["msc"] += energy;
} else if (process_name == "hIoni") {
PerProcessEnergyInSensitive["hadronic_ionization"] += energy;
if (track->GetParticleName() == "proton") {
PerProcessEnergyInSensitive["proton_ionization"] += energy;
}
} else if (process_name == "hadElastic") {
PerProcessEnergyInSensitive["hadronic_elastic"] += energy;
if (track->GetParticleName() == "neutron") {
PerProcessEnergyInSensitive["neutron_elastic"] += energy;
}
} else if (process_name == "Transportation") {
if (track->GetParticleName() == "proton") {
PerProcessEnergyInSensitive["hadronic_ionization"] += energy;
PerProcessEnergyInSensitive["proton_ionization"] += energy;
} else if (track->GetParticleName() == "e-" || track->GetParticleName() == "e+") {
PerProcessEnergyInSensitive["electron_ionization"] += energy;
}
}
}
}
}
protected:
#ifndef __CINT__
// TODO These graphs should be placed in TRestTrack?
// (following GetGraph implementation in TRestDetectorSignal)
TGraph* fXZHitGraph; //!
TGraph* fYZHitGraph; //!
TGraph* fXYHitGraph; //!
// TGraph2D *fXYZHitGraph; //! (TODO to implement XYZ visualization)
std::vector<TGraph*> fXYPcsMarker; //!
std::vector<TGraph*> fYZPcsMarker; //!
std::vector<TGraph*> fXZPcsMarker; //!
TMultiGraph* fXZMultiGraph; //!
TMultiGraph* fYZMultiGraph; //!
TMultiGraph* fXYMultiGraph; //!
// TMultiGraph *fXYZMultiGraph; //! (TODO to implement XYZ visualization)
TH2F* fXYHisto; //!
TH2F* fXZHisto; //!
TH2F* fYZHisto; //!
TH1D* fXHisto; //!
TH1D* fYHisto; //!
TH1D* fZHisto; //!
TLegend* fLegend_XY; //!
TLegend* fLegend_XZ; //!
TLegend* fLegend_YZ; //!
std::vector<Int_t> legendAdded; //!
Int_t fTotalHits; //!
TMultiGraph* GetXZMultiGraph(Int_t gridElement, std::vector<TString> pcsList, Double_t minPointSize = 0.4,
Double_t maxPointSize = 4);
TMultiGraph* GetYZMultiGraph(Int_t gridElement, std::vector<TString> pcsList, Double_t minPointSize = 0.4,
Double_t maxPointSize = 4);
TMultiGraph* GetXYMultiGraph(Int_t gridElement, std::vector<TString> pcsList, Double_t minPointSize = 0.4,
Double_t maxPointSize = 4);
TH2F* GetXYHistogram(Int_t gridElement, std::vector<TString> optList);
TH2F* GetXZHistogram(Int_t gridElement, std::vector<TString> optList);
TH2F* GetYZHistogram(Int_t gridElement, std::vector<TString> optList);
TH1D* GetXHistogram(Int_t gridElement, std::vector<TString> optList);
TH1D* GetYHistogram(Int_t gridElement, std::vector<TString> optList);
TH1D* GetZHistogram(Int_t gridElement, std::vector<TString> optList);
#endif
TVector3 fPrimaryEventOrigin;
std::vector<TString> fPrimaryParticleName;
std::vector<TVector3> fPrimaryEventDirection;
std::vector<Double_t> fPrimaryEventEnergy;
Double_t fTotalDepositedEnergy;
Double_t fSensitiveVolumeEnergy;
Int_t fNVolumes;
std::vector<Int_t> fVolumeStored;
std::vector<string> fVolumeStoredNames;
std::vector<Double_t> fVolumeDepositedEnergy;
Int_t fNTracks;
std::vector<TRestGeant4Track> fTrack;
Int_t fMaxSubEventID;
public:
void SetBoundaries();
void SetBoundaries(Double_t xMin, Double_t xMax, Double_t yMin, Double_t yMax, Double_t zMin,
Double_t zMax);
TString GetPrimaryEventParticleName(int n) {
if (fPrimaryParticleName.size() > n) return fPrimaryParticleName[n];
return "Not defined";
}
TVector3 GetPrimaryEventDirection(int n) { return fPrimaryEventDirection[n]; }
TVector3 GetPrimaryEventOrigin() { return fPrimaryEventOrigin; }
Double_t GetPrimaryEventEnergy(int n) { return fPrimaryEventEnergy[n]; }
Int_t GetNumberOfHits();
Int_t GetNumberOfTracks() { return fNTracks; }
Int_t GetNumberOfPrimaries() { return fPrimaryEventDirection.size(); }
Int_t GetNumberOfActiveVolumes() { return fNVolumes; }
Int_t isVolumeStored(int n) { return fVolumeStored[n]; }
TRestGeant4Track* GetTrack(int n) { return &fTrack[n]; }
TRestGeant4Track* GetTrackByID(int id);
Int_t GetNumberOfSubEventIDTracks() { return fMaxSubEventID + 1; }
Double_t GetTotalDepositedEnergy() { return fTotalDepositedEnergy; }
Double_t GetTotalDepositedEnergyFromTracks();
Double_t GetEnergyDepositedInVolume(Int_t volID) { return fVolumeDepositedEnergy[volID]; }
Double_t GetSensitiveVolumeEnergy() { return fSensitiveVolumeEnergy; }
TVector3 GetMeanPositionInVolume(Int_t volID);
TVector3 GetFirstPositionInVolume(Int_t volID);
TVector3 GetLastPositionInVolume(Int_t volID);
TRestHits GetHits();
Int_t GetNumberOfTracksForParticle(TString parName);
Int_t GetEnergyDepositedByParticle(TString parName);