* Adding the material manager class, still no other class use it however

* changing vis.mac file to use OGLI instead of OGLIQt interface

NPSimulation/include/MaterialManager.hh

+#ifndef __MATERIALMANAGER__
+#define __MATERIALMANAGER__
+/*****************************************************************************
+ * Copyright (C) 2009-2014   this file is part of the NPTool Project         *
+ *                                                                           *
+ * For the licensing terms see $NPTOOL/Licence/NPTool_Licence                *
+ * For the list of contributors see $NPTOOL/Licence/Contributors             *
+ *****************************************************************************/
+
+/*****************************************************************************
+ * Original Author: Adrien MATTA  contact address: matta@ipno.in2p3.fr       *
+ *                                                                           *
+ * Creation Date  : October 2014                                             *
+ * Last update    :                                                          *
+ *---------------------------------------------------------------------------*
+ * Decription:                                                               *
+ * This singleton class contain a librairy of material available for         *
+ * the geometry. Instantiate the needed material on demand, and generate the *
+ * associate DEDX table.                                                     *
+ *                                                                           *
+ *---------------------------------------------------------------------------*
+ * Comment:                                                                  *
+ *                                                                           *
+ *****************************************************************************/
+
+// Geant4
+#include"G4Material.hh"
+#include"G4Element.hh"
+#include"G4ParticleDefinition.hh"
+
+// STL
+#include<map>
+using namespace std;
+
+class MaterialManager{
+
+public:
+  // Designed to be a singleton (i.e. only one instance
+  // can exist). A member function called Instance is defined, which allows
+  // the user to get a pointer to the existing instance or to create it if
+  // it does not yet exist:
+  // (see the constructor for an explanation of the arguments)
+  static MaterialManager* getInstance();
+
+  // The analysis class instance can be deleted by calling the Destroy
+  // method (NOTE: The class destructor is protected, and can thus not be
+  // called directly):
+  static void Destroy();
+
+protected:
+    // Constructor and Destructor are not public
+    MaterialManager();
+    ~MaterialManager();
+
+private:
+  // The static instance of MaterialManager:
+  static MaterialManager* instance;
+  
+  // Map of element and material:
+  map<string,G4Element*> m_Element;
+  map<string,G4Material*> m_Material; 
+ 
+public:
+  G4Material* GetMaterialFromLibrary(string Name);
+  G4Element*  GetElementFromLibrary(string Name);
+  void WriteDEDXTable(G4ParticleDefinition* Particle ,G4double Emin,G4double Emax);
+};
+#endif

NPSimulation/include/ParticleStack.hh

@@ -32,38 +32,38 @@
 #include <vector>
 #include <string>
 using namespace std;
-using namespace CLHEP;
+using namespace CLHEP;
 
 // G4
 #include "G4ParticleGun.hh"
 #include "G4Event.hh"
 
 class ParticleStack{
-public:
+  public:
     // Designed to be a singleton (i.e. only one instance
     // can exist). A member function called Instance is defined, which allows
     // the user to get a pointer to the existing instance or to create it if
     // it does not yet exist:
     // (see the constructor for an explanation of the arguments)
     static ParticleStack* getInstance();
-    
+
     // The analysis class instance can be deleted by calling the Destroy
     // method (NOTE: The class destructor is protected, and can thus not be
     // called directly):
     static void Destroy();
-    
-protected:
+
+  protected:
     // Constructor (protected)
     ParticleStack();
-    
+
     // Destructor (protected)
     ~ParticleStack();
-    
+
     // Prevent copying
     ParticleStack(const ParticleStack& only);
     const ParticleStack& operator=(const ParticleStack& only);
-    
-private:
+
+  private:
     // The static instance of the ParticleStack class:
     static ParticleStack* instance;    
     // The particle gun
@@ -73,28 +73,28 @@ private:
     // Needed to now if the IC has to be cleared or not
     // true if the particle to be added is the first one
     bool m_First;
-  
-private: // Private Member
+
+  private: // Private Member
     vector<Particle> m_ParticleStack;
 
-public: // Getter and Setter
+  public: // Getter and Setter
     vector<Particle> GetParticleStack();
     void             SetParticleStack(vector<Particle>);
     G4ParticleGun    GetParticleGun();
 
-public: // Particle management and shooting method
+  public: // Particle management and shooting method
     // EventGenerator use this method to add particle in the stack
     void AddParticleToStack(Particle&);
-  
+
     // Add the Particle to the stack and fill Initial Condition Beam field 
-  void AddBeamParticleToStack(Particle&);
-  
+    void AddBeamParticleToStack(Particle&);
+
     // Search for a specific particle in the stack
     Particle SearchAndRemoveParticle(string);
-    
+
     // Transform the particle name to G4 standard: i.e: 10He -> He10
     string ChangeNameToG4Standard(string);
-    
+
     // Shoot everything in the stack at the end of the event
     void ShootAllParticle(G4Event* anEvent);
 };

NPSimulation/src/MaterialManager.cc

+/*****************************************************************************
+ * Copyright (C) 2009-2014   this file is part of the NPTool Project         *
+ *                                                                           *
+ * For the licensing terms see $NPTOOL/Licence/NPTool_Licence                *
+ * For the list of contributors see $NPTOOL/Licence/Contributors             *
+ *****************************************************************************/
+
+/*****************************************************************************
+ * Original Author: Adrien MATTA  contact address: matta@ipno.in2p3.fr       *
+ *                                                                           *
+ * Creation Date  : October 2014                                             *
+ * Last update    :                                                          *
+ *---------------------------------------------------------------------------*
+ * Decription:                                                               *
+ * This singleton class contain a librairy of material available for         *
+ * the geometry. Instantiate the needed material on demand, and generate the *
+ * associate DEDX table.                                                     *
+ *                                                                           *
+ *---------------------------------------------------------------------------*
+ * Comment:                                                                  *
+ *                                                                           *
+ *****************************************************************************/
+
+// NPS
+#include "MaterialManager.hh"
+
+// Geant4
+#include "G4EmCalculator.hh"
+
+// STL
+#include<iostream>
+using namespace std;
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+MaterialManager* MaterialManager::MaterialManager::instance = 0 ;
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+MaterialManager* MaterialManager::getInstance(){
+  if(instance == 0) instance = new MaterialManager();
+  return instance;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+MaterialManager::~MaterialManager(){
+
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+void MaterialManager::Destroy(){
+
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
+G4Material* MaterialManager::GetMaterialFromLibrary(string Name){
+  // Search if the material is already instantiate  
+  map<string,G4Material*>::iterator it;
+  it = m_Material.find(Name);
+
+  // The element is not found
+  if(it==m_Material.end()){
+
+    //Usual compound
+    if(Name == "Vacuum"){
+      G4Material* material = new G4Material(Name, 0.000000001 * mg / cm3,2);
+      material->AddElement(GetElementFromLibrary("N"),7);
+      material->AddElement(GetElementFromLibrary("O"),3);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else if(Name == "PCB"){
+      // Actually taken value fron Epoxy
+      G4Material* material = new G4Material(Name, 1.85 * g / cm3,3);
+      material->AddElement(GetElementFromLibrary("H"),.475);
+      material->AddElement(GetElementFromLibrary("C"),.45);
+      material->AddElement(GetElementFromLibrary("O"),.075);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Mylar"){
+      G4Material* material = new G4Material(Name, 1.397*g,3);
+      material->AddElement(GetElementFromLibrary("H"),8);
+      material->AddElement(GetElementFromLibrary("C"),10);
+      material->AddElement(GetElementFromLibrary("O"),4);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Kapton"){
+      G4Material* material = new G4Material(Name, 1.39*g,3);
+      material->AddElement(GetElementFromLibrary("H"),4);
+      material->AddElement(GetElementFromLibrary("C"),5);
+      material->AddElement(GetElementFromLibrary("O"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Harvar"){
+      G4Material* material = new G4Material(Name, 8.3*g / cm3,5);
+      material->AddElement(GetElementFromLibrary("Co"),42);
+      material->AddElement(GetElementFromLibrary("Cr"),20);
+      material->AddElement(GetElementFromLibrary("Ni"),13);
+      material->AddElement(GetElementFromLibrary("Fe"),19);
+      material->AddElement(GetElementFromLibrary("W"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+    
+    // Usual Target
+     else if(Name == "CD2"){
+      G4Material* material = new G4Material(Name, 1.06*g/cm3,2);
+      material->AddElement(GetElementFromLibrary("C"),1);
+      material->AddElement(GetElementFromLibrary("D"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "CH2"){
+      G4Material* material = new G4Material(Name, 0.93*g/cm3,2);
+      material->AddElement(GetElementFromLibrary("C"),1);
+      material->AddElement(GetElementFromLibrary("D"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Pb"){
+      G4Material* material = new G4Material(Name, 11.342*g / cm3,1);
+      material->AddElement(GetElementFromLibrary("Pb"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "D2_solid"){ // density should be checked
+      G4Material* material = new G4Material(Name, 0.0715*g/cm3,1);
+      material->AddElement(GetElementFromLibrary("D"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "H2_solid"){
+      G4Material* material = new G4Material(Name, 0.0715*g/cm3,1);
+      material->AddElement(GetElementFromLibrary("H"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+    
+    // Usual detector material
+    else  if(Name == "Si"){
+      G4Material* material = new G4Material(Name,2.321*g/cm3 ,1);
+      material->AddElement(GetElementFromLibrary("Si"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Ge"){
+      G4Material* material = new G4Material(Name,5.323*g/cm3 ,1);
+      material->AddElement(GetElementFromLibrary("Ge"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+   
+    else  if(Name == "BC400"){
+      G4Material* material = new G4Material(Name,1.032*g/cm3 ,2);
+      material->AddElement(GetElementFromLibrary("H"),10);
+      material->AddElement(GetElementFromLibrary("C"),9);
+      m_Material[Name]=material;
+      return material; 
+    }
+    
+    else  if(Name == "NaI"){
+      G4Material* material = new G4Material(Name, 3.67*g/cm3,2);
+      material->AddElement(GetElementFromLibrary("Na"),1);
+      material->AddElement(GetElementFromLibrary("I"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "CsI"){
+      G4Material* material = new G4Material(Name, 4.51*g/cm3,2);
+      material->AddElement(GetElementFromLibrary("Cs"),1);
+      material->AddElement(GetElementFromLibrary("I"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "LaBr3"){
+      G4Material* material = new G4Material(Name, 5.29*g/cm3 , 2);
+      material->AddElement(GetElementFromLibrary("La"),1);
+      material->AddElement(GetElementFromLibrary("Br"),3);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "BaF2"){
+      G4Material* material = new G4Material(Name, 4.89*g/cm3 , 2);
+      material->AddElement(GetElementFromLibrary("Ba"),1);
+      material->AddElement(GetElementFromLibrary("F"),2);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    // Misc
+    else  if(Name == "Al"){
+      G4Material* material = new G4Material(Name, 2.702*g/cm3,1);
+      material->AddElement(GetElementFromLibrary("Al"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else  if(Name == "Fe"){
+      G4Material* material = new G4Material(Name, 7.874*g/cm3,1);
+      material->AddElement(GetElementFromLibrary("Fe"),1);
+      m_Material[Name]=material;
+      return material; 
+    }
+
+    else{
+      cout << "ERROR: Material requested \""<< Name <<"\" is not available in the Material Librairy" << endl;
+      exit(1);
+    }
+
+  }
+
+  else return it->second;
+
+}
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+G4Element* MaterialManager::GetElementFromLibrary(string Name){
+  // Search if the element is already instantiate  
+  map<string,G4Element*>::iterator it;
+  it = m_Element.find(Name);
+
+  // The element is not found
+  if(it==m_Element.end()){
+    if(Name == "H"){
+      G4Element* element = new G4Element(Name,Name,1, 1.008*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "D"){
+      G4Element* element = new G4Element(Name,Name,1, 2.01410178*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "T"){
+      G4Element* element = new G4Element(Name,Name,1, 3.0160492*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "C"){
+      G4Element* element = new G4Element(Name,Name,6, 12.011*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "O"){
+      G4Element* element = new G4Element(Name,Name,8, 15.999*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "N"){
+      G4Element* element = new G4Element(Name,Name,7, 14.007*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Si"){
+      G4Element* element = new G4Element(Name,Name,14, 28.085*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Al"){
+      G4Element* element = new G4Element(Name,Name,13, 26.9815385*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Co"){
+      G4Element* element = new G4Element(Name,Name,27, 58.933194*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Cr"){
+      G4Element* element = new G4Element(Name,Name,24, 51.9961*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Ni"){
+      G4Element* element = new G4Element(Name,Name,28, 58.6934*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Fe"){
+      G4Element* element = new G4Element(Name,Name,26, 55.845*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "W"){
+      G4Element* element = new G4Element(Name,Name,74, 183.84*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Pb"){
+      G4Element* element = new G4Element(Name,Name,82, 207.2*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Ge"){
+      G4Element* element = new G4Element(Name,Name,32, 72.630*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Na"){
+      G4Element* element = new G4Element(Name,Name,11, 22.98976928*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "I"){
+      G4Element* element = new G4Element(Name,Name,53, 126.90447*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "La"){
+      G4Element* element = new G4Element(Name,Name,57, 138.90547*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Br"){
+      G4Element* element = new G4Element(Name,Name,35, 79.904*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Bi"){
+      G4Element* element = new G4Element(Name,Name,83, 208.98040*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "Ba"){
+      G4Element* element = new G4Element(Name,Name,56, 137.327*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else if(Name == "F"){
+      G4Element* element = new G4Element(Name,Name,9, 18.998403163*g/mole);
+      m_Element[Name] = element;
+      return element;
+    }
+
+    else{
+      cout << "ERROR: Element requested \""<< Name <<"\" is not available in the Material Librairy" << endl;
+      exit(1);
+    }
+
+
+  }
+
+  else return it->second;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//   Generate a DEDX file table using the material used in the geometry
+void MaterialManager::WriteDEDXTable(G4ParticleDefinition* Particle ,G4double Emin,G4double Emax){
+  map<string,G4Material*>::iterator it;
+  for(it = m_Material.begin() ; it != m_Material.end() ; it++){
+    //   Opening hte output file
+    G4String GlobalPath = getenv("NPTOOL");
+    G4String Path = GlobalPath + "/Inputs/EnergyLoss/" + Particle->GetParticleName() + "_" + it->second->GetName() + ".G4table";
+  
+    ofstream File ;
+    File.open(Path) ;
+
+    if(!File) return ;
+
+    File   << "Table from Geant4 generate using NPSimulation \t"
+      << "Particle: " << Particle->GetParticleName() << "\tMaterial: " << it->second->GetName() << endl ;
+
+    G4EmCalculator emCalculator;
+    G4double dedx ;
+    // Tipical Range needed, if Emax is larger, then adapted
+    if(Emax < 25000) Emax = 25000;    
+    double step = 1*keV;
+    double before = 0 ;
+    for (G4double E=Emin; E < Emax; E+=step){
+       
+      dedx = emCalculator.ComputeTotalDEDX(E, Particle, it->second)/(MeV/micrometer);
+
+      if(before){
+        if(abs(before-dedx)/abs(before)<0.01) step*=10; 
+      }
+
+      before = dedx;
+      File << E/MeV << "\t" << dedx << endl ;
+    }
+
+    File.close();
+  }
+}
+

NPSimulation/vis.mac

@@ -8,7 +8,7 @@
 /run/verbose 2
 #
 # Use this open statement to create an OpenGL view:
-/vis/open OGLIQt 600x600-0+0
+/vis/open OGLI 600x600-0+0
 #
 # Use this open statement to create a .prim file suitable for
 # viewing in DAWN: