Skip to content
Snippets Groups Projects
Commit fe69bdaf authored by GIRARD ALCINDOR Valérian's avatar GIRARD ALCINDOR Valérian
Browse files

Modifiying Exogam simu

parent 59bf833d
No related branches found
No related tags found
1 merge request!27Draft: [Epic] Preparation of the environement for the new GaseousDetectorScorers...
Pipeline #355010 passed
Stage: build-NPLib
Stage: build-NPSimulation
Stage: test
Hide whitespace changes
Inline Side-by-side
NPSimulation/Detectors/Exogam/Exogam.cc
+ 125
− 23
View file @ fe69bdaf
......@@ -70,7 +70,8 @@ using namespace CLHEP;
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
namespace Exogam_NS {
// Energy and time Resolution
const double EnergyThreshold = 10 * keV;
// const double EnergyThreshold = 10 * keV;
const double EnergyThreshold = -100 * keV;
// const double ResoTime = 4.5*ns ; //not used
const double ResoEnergy = 2. * keV;
} // namespace Exogam_NS
......@@ -79,7 +80,7 @@ namespace Exogam_NS {
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
// Exogam Specific Method
Exogam::Exogam() {
m_Event = new TExogamData();
m_Event = new TExogamCalData();
m_ExogamScorer = 0;
InitializeMaterials();
......@@ -127,12 +128,65 @@ void Exogam::BuildClover(int i_clo, G4LogicalVolume* world) {
G4LogicalVolume* logicSupClover = new G4LogicalVolume(solidSupClover, m_Vacuum, "SupClover");
Offset = dzEnv; //-distCollimatorToGeCan;
// Offset = 0;
G4RotationMatrix rm;
rm.rotateX(m_ThetaX[i_clo] / rad).rotateY(m_ThetaY[i_clo] / rad).rotateZ(m_ThetaZ[i_clo] / rad);
new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(m_X[i_clo] * mm, m_Y[i_clo] * mm, m_Z[i_clo] * mm + Offset)),
logicSupClover, "Clover", world, false, i_clo + 1, false); // this void overlaps the whole setup
if (m_X[i_clo] >= 0) {
G4RotationMatrix rm;
rm.rotateX(m_ThetaX[i_clo] / rad).rotateY(m_ThetaY[i_clo] / rad).rotateZ(m_ThetaZ[i_clo] / rad);
new G4PVPlacement(G4Transform3D(rm, G4ThreeVector(m_X[i_clo] * mm, m_Y[i_clo] * mm, m_Z[i_clo] * mm + Offset)),
logicSupClover, "Clover", world, false, i_clo + 1, false); // this void overlaps the whole setup
}
else if (m_R[i_clo] >= 0) {
G4RotationMatrix* MMrot = NULL;
G4ThreeVector MMpos = G4ThreeVector(0, 0, 0);
G4ThreeVector MMu = G4ThreeVector(0, 0, 0);
G4ThreeVector MMv = G4ThreeVector(0, 0, 0);
G4ThreeVector MMw = G4ThreeVector(0, 0, 0);
G4ThreeVector MMCenter = G4ThreeVector(0, 0, 0);
G4double Theta = m_Theta[i_clo];
G4double Phi = m_Phi[i_clo];
// (u,v,w) unitary vector associated to telescope referencial
// (u,v) // to silicon plan
// w perpendicular to (u,v) plan and pointing ThirdStage
// Phi is angle between X axis and projection in (X,Y) plan
// Theta is angle between position vector and z axis
G4double wX = (m_R[i_clo]+Offset) * sin(Theta / rad) * cos(Phi / rad);
G4double wY = (m_R[i_clo]+Offset) * sin(Theta / rad) * sin(Phi / rad);
G4double wZ = (m_R[i_clo]+Offset) * cos(Theta / rad);
MMw = G4ThreeVector(wX, wY, wZ);
// vector corresponding to the center of the module
G4ThreeVector CT = MMw;
// vector parallel to one axis of silicon plane
G4double ii = cos(Theta / rad) * cos(Phi / rad);
G4double jj = cos(Theta / rad) * sin(Phi / rad);
G4double kk = -sin(Theta / rad);
G4ThreeVector Y = G4ThreeVector(ii, jj, kk);
MMw = MMw.unit();
MMu = MMw.cross(Y);
MMv = MMw.cross(MMu);
MMv = MMv.unit();
MMu = MMu.unit();
// Passage Matrix from Lab Referential to Telescope Referential
// MUST2
MMrot = new G4RotationMatrix(MMu, MMv, MMw);
// Telescope is rotate of Beta angle around MMv axis.
// MMrot->rotate(m_beta_u[i_clo], MMu);
// MMrot->rotate(m_beta_v[i_clo], MMv);
// MMrot->rotate(m_beta_w[i_clo], MMw);
// translation to place Telescope
double Length = 0 * cm;
MMpos = MMw * Length * 0.5 + CT;
std::cout << i_clo << std::endl;
std::string clover_name = "Clover_" + std::to_string(i_clo + 1);
new G4PVPlacement(G4Transform3D(*MMrot, MMpos), logicSupClover, clover_name, world, false, i_clo + 1,
false); // this void overlaps the whole setup
}
// The Cryostat
////////////////
......@@ -157,7 +211,7 @@ void Exogam::BuildClover(int i_clo, G4LogicalVolume* world) {
rOuterCan[1] = 6.2 * cm;
rOuterCan[2] = 6.2 * cm;
rInnerCan[0] = rInnerCan[1] = rInnerCan[2] = 0. * cm;
rInnerCan[0] = rInnerCan[1] = rInnerCan[2] = 0.1 * cm;
G4Polyhedra* solidCloverCan =
new G4Polyhedra("CloverCan", PhiStartCan, PhiTotCan, 4, 3, zPlaneCan, rInnerCan, rOuterCan);
......@@ -414,15 +468,24 @@ void Exogam::BuildClover(int i_clo, G4LogicalVolume* world) {
// Now the individual diode is built; create logical volumes for each of
// the four individual diodes A, B, C and D:
G4LogicalVolume* logicGeA = new G4LogicalVolume(solidGe, m_Germanium, "GeA");
G4LogicalVolume* logicGeB = new G4LogicalVolume(solidGe, m_Germanium, "GeB");
G4LogicalVolume* logicGeC = new G4LogicalVolume(solidGe, m_Germanium, "GeC");
G4LogicalVolume* logicGeD = new G4LogicalVolume(solidGe, m_Germanium, "GeD");
std::string GeA_name = "GeA_" + std::to_string(i_clo + 1);
G4LogicalVolume* logicGeA = new G4LogicalVolume(solidGe, m_Germanium, GeA_name);
std::string GeB_name = "GeB_" + std::to_string(i_clo + 1);
G4LogicalVolume* logicGeB = new G4LogicalVolume(solidGe, m_Germanium, GeB_name);
std::string GeC_name = "GeC_" + std::to_string(i_clo + 1);
G4LogicalVolume* logicGeC = new G4LogicalVolume(solidGe, m_Germanium, GeC_name);
std::string GeD_name = "GeD_" + std::to_string(i_clo + 1);
G4LogicalVolume* logicGeD = new G4LogicalVolume(solidGe, m_Germanium, GeD_name);
std::cout << "SetSensitiveDetector: " << i_clo << std::endl;
logicGeA->SetSensitiveDetector(m_ExogamScorer);
std::cout << "GeA SensitiveDetector: " << logicGeA->GetSensitiveDetector()->GetName() << std::endl;
logicGeB->SetSensitiveDetector(m_ExogamScorer);
std::cout << "GeB SensitiveDetector: " << logicGeB->GetSensitiveDetector()->GetName() << std::endl;
logicGeC->SetSensitiveDetector(m_ExogamScorer);
std::cout << "GeC SensitiveDetector: " << logicGeC->GetSensitiveDetector()->GetName() << std::endl;
logicGeD->SetSensitiveDetector(m_ExogamScorer);
std::cout << "GeD SensitiveDetector: " << logicGeD->GetSensitiveDetector()->GetName() << std::endl;
// positioning the tapered partial diodes (A to D)
// into the real vacuum of the can
......@@ -443,10 +506,15 @@ void Exogam::BuildClover(int i_clo, G4LogicalVolume* world) {
G4ThreeVector posDumVacD(-xDumVac, yDumVac, zDumVac);
G4Transform3D positionVacD(rm90, posDumVacD);
new G4PVPlacement(0, positionVacA, logicGeA, "GeA", logicVac, false, 1, true); // There is an overlap with vacumm Vac
new G4PVPlacement(positionVacB, logicGeB, "GeB", logicVac, false, 2, true);
new G4PVPlacement(positionVacC, logicGeC, "GeC", logicVac, false, 3, true);
new G4PVPlacement(positionVacD, logicGeD, "GeD", logicVac, false, 4, true);
std::string GeA_name2 = "GeA_" + std::to_string(i_clo + 1);
new G4PVPlacement(0, positionVacA, logicGeA, GeA_name, logicVac, false, 1,
true); // There is an overlap with vacumm Vac
std::string GeB_name2 = "GeB_" + std::to_string(i_clo + 1);
new G4PVPlacement(positionVacB, logicGeB, GeB_name, logicVac, false, 2, true);
std::string GeC_name2 = "GeC_" + std::to_string(i_clo + 1);
new G4PVPlacement(positionVacC, logicGeC, GeC_name, logicVac, false, 3, true);
std::string GeD_name2 = "GeD_" + std::to_string(i_clo + 1);
new G4PVPlacement(positionVacD, logicGeD, GeD_name, logicVac, false, 4, true);
//
// some material between the diodes to reproduce the experimental addback factor ...
......@@ -530,6 +598,24 @@ void Exogam::AddDetector(double X, double Y, double Z, double ThetaX, double The
m_ThetaX.push_back(ThetaX);
m_ThetaY.push_back(ThetaY);
m_ThetaZ.push_back(ThetaZ);
m_R.push_back(-1000);
m_Theta.push_back(-1000);
m_Phi.push_back(-1000);
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void Exogam::AddDetector(double R, double Theta, double Phi) {
m_R.push_back(R);
m_Theta.push_back(Theta);
m_Phi.push_back(Phi);
m_X.push_back(-1000);
m_Y.push_back(-1000);
m_Z.push_back(-1000);
m_ThetaX.push_back(-1000);
m_ThetaY.push_back(-1000);
m_ThetaZ.push_back(-1000);
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
......@@ -546,6 +632,7 @@ void Exogam::ReadConfiguration(NPL::InputParser parser) {
cout << "//// " << blocks.size() << " detectors found " << endl;
vector<string> coord = {"X", "Y", "Z", "ThetaX", "ThetaY", "ThetaZ"};
vector<string> sphe = {"R", "Theta", "Phi"};
for (unsigned int i = 0; i < blocks.size(); i++) {
if (blocks[i]->HasTokenList(coord)) {
......@@ -559,6 +646,14 @@ void Exogam::ReadConfiguration(NPL::InputParser parser) {
double ThetaZ = blocks[i]->GetDouble("ThetaZ", "deg");
AddDetector(X, Y, Z, ThetaX, ThetaY, ThetaZ);
}
else if (blocks[i]->HasTokenList(sphe)) {
if (NPOptionManager::getInstance()->GetVerboseLevel())
cout << endl << "//// Exogam " << i + 1 << endl;
double R = blocks[i]->GetDouble("R", "mm");
double Theta = blocks[i]->GetDouble("Theta", "deg");
double Phi = blocks[i]->GetDouble("Phi", "deg");
AddDetector(R, Theta, Phi);
}
else {
cout << "ERROR: check your input file formatting " << endl;
exit(1);
......@@ -596,7 +691,7 @@ void Exogam::InitializeRootOutput() {
RootOutput* pAnalysis = RootOutput::getInstance();
TTree* pTree = pAnalysis->GetTree();
if (!pTree->FindBranch("Exogam")) {
pTree->Branch("Exogam", "TExogamData", &m_Event);
pTree->Branch("Exogam", "TExogamCalData", &m_Event);
}
pTree->SetBranchAddress("Exogam", &m_Event);
}
......@@ -612,17 +707,24 @@ void Exogam::ReadSensitive(const G4Event*) {
CalorimeterScorers::PS_Calorimeter* Scorer = (CalorimeterScorers::PS_Calorimeter*)m_ExogamScorer->GetPrimitive(0);
unsigned int size = Scorer->GetMult();
// if (size > 0)
// std::cout << "/////////: " << size << " "<< Scorer->GetName()<< std::endl;
for (unsigned int i = 0; i < size; i++) {
double Energy = RandGauss::shoot(Scorer->GetEnergy(i), Exogam_NS::ResoEnergy);
// double Energy = RandGauss::shoot(Scorer->GetEnergy(i) * MeV, Exogam_NS::ResoEnergy);
double Energy = Scorer->GetEnergy(i) * MeV;
if (Energy > Exogam_NS::EnergyThreshold) {
double Time = Scorer->GetTime(i);
int CristalNbr = Scorer->GetLevel(i)[0];
int CloverNbr = Scorer->GetLevel(i)[1];
//FIXME
//m_Event->SetCristal(CristalNbr);
//m_Event->SetClover(CloverNbr);
//m_Event->SetEnergy(Energy);
//m_Event->SetTime(Time);
// FIXME
// m_Event->SetCristal(CristalNbr);
// m_Event->SetClover(CloverNbr);
// m_Event->SetEnergy(Energy);
// m_Event->SetTime(Time);
int CristalID = CristalNbr + (CloverNbr - 1) * 4;
m_Event->SetExo(CristalID, Energy, -1000, Time, -1000, -1000, -1000, -1000, -1000, -1000, -1000);
}
}
}
......
NPSimulation/Detectors/Exogam/Exogam.hh
+ 10
− 3
View file @ fe69bdaf
......@@ -34,7 +34,7 @@ using namespace std;
// NPTool header
#include "NPSVDetector.hh"
#include "TExogamData.h"
#include "TExogamCalData.h"
#include "NPInputParser.h"
class Exogam : public NPS::VDetector{
......@@ -49,8 +49,10 @@ class Exogam : public NPS::VDetector{
/////// Specific Function of this Class ///////////
////////////////////////////////////////////////////
public:
// Spherical
// Cartesian
void AddDetector(double X,double Y, double Z, double ThetaX, double ThetaY, double ThetaZ);
// Spherical
void AddDetector(double R,double Theta, double Phi);
G4int InitializeMaterials();
void BuildClover(int i_clo, G4LogicalVolume* world);
......@@ -91,7 +93,7 @@ class Exogam : public NPS::VDetector{
///////////Event class to store Data////////////////
////////////////////////////////////////////////////
private:
TExogamData* m_Event ;
TExogamCalData* m_Event ;
////////////////////////////////////////////////////
///////////////Private intern Data//////////////////
......@@ -105,6 +107,11 @@ class Exogam : public NPS::VDetector{
vector<double> m_ThetaX; //rotation angles to X, Y, Z axis
vector<double> m_ThetaY;
vector<double> m_ThetaZ;
// Detector Coordinate in spherical
vector<double> m_R;
vector<double> m_Theta;
vector<double> m_Phi;
private:
//materials
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment