diff --git a/README.md b/README.md index 1842173..89aeef0 100644 --- a/README.md +++ b/README.md @@ -1,16 +1,16 @@ # MASK: Monte cArlo Simulation of Kinematics MASK is a Monte Carlo simulation of reaction kinematics for use detector systems at Florida State University. -MASK is capable of simulating multi-step kinematic reaction-decay sequences, storing data in a lightweight binary format, after which the kinematic data can be fed to a detector geometry for efficiency testing. Currently geometries for ANASEN and SABRE are included in the code. Mask depends on the CERN ROOT analysis framework. All data is stored in the ROOT format using a custom dictionary. +MASK is capable of simulating multi-step kinematic reaction-decay sequences, storing data in ROOT Trees, after which the kinematic data can be fed to a detector geometry for efficiency testing. Currently geometries for ANASEN and SABRE are included in the code. ## Building MASK -To clone MASK and all submodules use `git clone --recursive https://github.com/gwm17/Kinematics.git`. -Mask is built using the CMake build tool. On most systems this can be done using the following commands inside the Mask repository: - +Dowload the repository from github. CMake is use to build the project; in most environments you can build Mask using the following methods: - `mkdir build` - `cd build` - `cmake ..` - `make` +Executables will be installed to the repostiory's `bin` directory. Libraries will be installed to the `lib` directory. + By default Mask builds for release. To build for debug replace `cmake ..` with `cmake -DCMAKE_BUILD_TYPE=Debug ..`. Mask uses CMake to find the installed ROOT libraries and headers. ## Using the kinematics simulation @@ -21,7 +21,7 @@ By default MASK is capable of simulating reactions of up to three steps. Here is 2. A two step reaction is a single step reaction where the residual subsequently decays into two daughters (called breakups). Again, all sampling is allowed. 3. A three step reaction is a two step reaction where one of the breakup particles subsequently decays into two more breakup particles. Again, all sampling is allowed -For decays, a specific angular distribution can be given as input as a text file with values of coefficiencts of a Legendre polynomial series. Examples can be found in the `./etc` directory, including an isotropic case. It is assumed that the decays in the center-of-mass frame are isotropic in phi (i.e. m=0). Decay1 corresponds to the first decay, if there are multiple steps, Decay2 to the second. If there are no decays, these parameters are not used (or if only one decay, Decay2_AngularMomentum is not used). +For decays, a specific angular distribution can be given as input as a text file with values of coefficiencts of a Legendre polynomial series. Examples can be found in the `./etc` directory, including an isotropic case. It is assumed that the decays in the center-of-mass frame are isotropic in phi (i.e. m=0). Decay1 corresponds to the first decay, if there are multiple steps, Decay2 to the second. If there are no decays, these parameters are not used (or if only one decay, Decay2_AngularMomentum is not used). The input file requires that the user include target information, which will be used to calculate energy loss for all of the reactants and reaction products. The energy loss through materials is calculated using the `catima` library (found in `src/vendor`), which is a C/C++ interface to the `atima` library (the same energy loss methods used by LISE). The target can contain layers, and each layer can be composed of a compound of elements with a given stoichiometry. If the user wishes to not include energy loss in the kinematics, simply give all target layers a thickness of 0. Note that more layers and more thickness = more time spent calculating energy loss. These energy loss methods are only applicable for solid targets, and should not be applied to gas or liquid targets. Energy loss calculations have a stated uncertainty of approximately five percent. The input file requires that the user include target information, which will be used to calculate energy loss for all of the reactants and reaction products. The target can contain layers, and each layer can be composed of a compound of elements with a given stoichiometry. If the user wishes to not include energy loss in the kinematics, simply give all target layers a thickness of 0. Note that more layers and more thickness = more time spent calculating energy loss. These energy loss methods are only applicable for solid targets, and should not be applied to gas or liquid targets. Energy loss calculations have a stated uncertainty of approximately five percent. The energy loss library used is called `catima` and can be found [here](https://github.com/gwm17/catima). @@ -42,7 +42,7 @@ To run DetEff use the format `./bin/DetEff ` -where the detection datafile contains all of the kinematics data as well as information about which particles are detected (this is in the mask file format) and the statsfile is a text file containing efficiency statistics. +where the detection datafile contains all of the kinematics data as well as information about which particles are detected and the statsfile is a text file containing efficiency statistics. ## Data visualization All data is saved as ROOT trees of std::vectors of Mask::Nucleus classes. To enable this, a ROOT dictionary is generated and linked into a shared library found in the `lib` directory of the repository. This allows the user to link to the shared library for accessing and analyzing the data generated by MASK. @@ -51,8 +51,9 @@ Mask also provides a default visualization tool called RootPlot. RootPlot is run `./bin/RootPlot ` -where the datafile can be either the datafile from Mask or the datafile from DetEff.s +where the datafile can be either the datafile from Mask or the datafile from DetEff. The outputfile is saved in the ROOT file format. ## Requirements -- Requires CMake > 3.16 -- Requires ROOT > 6.16 \ No newline at end of file +ROOT version 6.22 or greater is required +CMake version 3.0 or greater is required + diff --git a/src/Mask/MassLookup.cpp b/src/Mask/MassLookup.cpp index 4bd1e7c..ad26256 100644 --- a/src/Mask/MassLookup.cpp +++ b/src/Mask/MassLookup.cpp @@ -41,7 +41,7 @@ namespace Mask { MassLookup::~MassLookup() {} //Returns nuclear mass in MeV - double MassLookup::FindMass(int Z, int A) + double MassLookup::FindMass(uint32_t Z, uint32_t A) { KeyPair key({Z, A}); auto data = massTable.find(key.GetID()); @@ -52,7 +52,7 @@ namespace Mask { } //returns element symbol - std::string MassLookup::FindSymbol(int Z, int A) + std::string MassLookup::FindSymbol(uint32_t Z, uint32_t A) { KeyPair key({Z, A}); auto data = elementTable.find(key.GetID()); diff --git a/src/Mask/MassLookup.h b/src/Mask/MassLookup.h index 61220cf..fac5e9e 100644 --- a/src/Mask/MassLookup.h +++ b/src/Mask/MassLookup.h @@ -35,9 +35,9 @@ namespace Mask { }; ~MassLookup(); - double FindMass(int Z, int A); - double FindMassU(int Z, int A) { return FindMass(Z, A)/u_to_mev; } - std::string FindSymbol(int Z, int A); + double FindMass(uint32_t Z, uint32_t A); + double FindMassU(uint32_t Z, uint32_t A) { return FindMass(Z, A)/u_to_mev; } + std::string FindSymbol(uint32_t Z, uint32_t A); static MassLookup& GetInstance() { return *s_instance; } diff --git a/src/Mask/Nucleus.cpp b/src/Mask/Nucleus.cpp index d375fb3..c82b7f5 100644 --- a/src/Mask/Nucleus.cpp +++ b/src/Mask/Nucleus.cpp @@ -2,7 +2,7 @@ namespace Mask { - Nucleus CreateNucleus(int z, int a) + Nucleus CreateNucleus(uint32_t z, uint32_t a) { Nucleus nuc; nuc.Z = z; diff --git a/src/Mask/Nucleus.h b/src/Mask/Nucleus.h index 7baef5c..2769cf6 100644 --- a/src/Mask/Nucleus.h +++ b/src/Mask/Nucleus.h @@ -50,7 +50,7 @@ namespace Mask { double detectedPhi = 0.0; }; - Nucleus CreateNucleus(int z, int a); + Nucleus CreateNucleus(uint32_t z, uint32_t a); bool EnforceDictionaryLinked(); diff --git a/src/Plotters/RootPlotter.cpp b/src/Plotters/RootPlotter.cpp index 476f934..34b5068 100644 --- a/src/Plotters/RootPlotter.cpp +++ b/src/Plotters/RootPlotter.cpp @@ -5,6 +5,11 @@ #include +static double FullPhi(double phi) +{ + return phi < 0.0 ? 2.0*M_PI + phi : phi; +} + RootPlotter::RootPlotter() { TH1::AddDirectory(kFALSE); @@ -53,7 +58,7 @@ void RootPlotter::Run(const std::string& inputname, const std::string& outputnam output->cd(); for(auto& obj : m_map) - obj.second->Write(obj.second->GetName(), TObject::kOverwrite); + obj.second->Write(); output->Close(); } @@ -76,15 +81,15 @@ void RootPlotter::FillData(const Mask::Nucleus& nuc) std::string angdist_title = angdist_name+";cos#right(#theta_{CM}#left);counts"; MyFill(ke_vs_th_name.c_str(), ke_vs_th_title.c_str(), nuc.vec4.Theta()*s_rad2deg, nuc.GetKE(), 2); - MyFill(ke_vs_ph_name.c_str(), ke_vs_ph_title.c_str(), nuc.vec4.Phi()*s_rad2deg, nuc.GetKE(), 4); - MyFill(th_vs_ph_name.c_str(), th_vs_ph_title.c_str(), nuc.vec4.Theta()*s_rad2deg, nuc.vec4.Phi()*s_rad2deg, 2); + MyFill(ke_vs_ph_name.c_str(), ke_vs_ph_title.c_str(), FullPhi(nuc.vec4.Phi())*s_rad2deg, nuc.GetKE(), 4); + MyFill(th_vs_ph_name.c_str(), th_vs_ph_title.c_str(), nuc.vec4.Theta()*s_rad2deg, FullPhi(nuc.vec4.Phi())*s_rad2deg, 2); MyFill(ex_name.c_str(),ex_title.c_str(),260,-1.0,25,nuc.GetExcitationEnergy()); MyFill(angdist_name.c_str(), angdist_title.c_str(),20,-1.0,1.0,std::cos(nuc.thetaCM)); if(nuc.isDetected) { MyFill(ke_vs_th_name.c_str(), ke_vs_th_title.c_str(), nuc.vec4.Theta()*s_rad2deg, nuc.detectedKE, 2); - MyFill(ke_vs_ph_name.c_str(), ke_vs_ph_title.c_str(), nuc.vec4.Phi()*s_rad2deg, nuc.detectedKE, 4); - MyFill(th_vs_ph_name.c_str(), th_vs_ph_title.c_str(), nuc.vec4.Theta()*s_rad2deg, nuc.vec4.Phi()*s_rad2deg, 2); + MyFill(ke_vs_ph_name.c_str(), ke_vs_ph_title.c_str(), FullPhi(nuc.vec4.Phi())*s_rad2deg, nuc.detectedKE, 4); + MyFill(th_vs_ph_name.c_str(), th_vs_ph_title.c_str(), nuc.vec4.Theta()*s_rad2deg, FullPhi(nuc.vec4.Phi())*s_rad2deg, 2); MyFill(ex_name.c_str(),ex_title.c_str(),260,-1.0,25,nuc.GetExcitationEnergy()); MyFill(angdist_name.c_str(), angdist_title.c_str(),20,-1.0,1.0,std::cos(nuc.thetaCM)); }