From ad5b9521711d5576b5c244f08e95067329843ffe Mon Sep 17 00:00:00 2001 From: matta <matta@npt> Date: Wed, 21 Oct 2009 05:00:11 +0000 Subject: [PATCH] * Adding support for Strip Number in ThinSi * Little improvment on MUST2 scorer - Map for X,Y and Pad/cristalNumber are now G4int instead of G4double - StripPitch is now a private member of X and Y scorer and is no longer calculated at every event * Adding some file that should be here like the Dummy detector, the ThinSiData or the DummyDetector.detector input file --- .../DUMMYDetector.detector | 25 + NPSimulation/include/DummyDetector.hh | 128 +++++ NPSimulation/include/Must2Scorers.hh | 8 +- NPSimulation/include/ThinSi.hh | 2 +- NPSimulation/include/ThinSiScorers.hh | 59 +++ NPSimulation/src/DummyDetector.cc | 473 ++++++++++++++++++ NPSimulation/src/MUST2Array.cc | 28 +- NPSimulation/src/Must2Scorers.cc | 43 +- NPSimulation/src/ThinSi.cc | 91 ++-- NPSimulation/src/ThinSiScorers.cc | 100 ++++ 10 files changed, 876 insertions(+), 81 deletions(-) create mode 100755 Inputs/DetectorConfiguration/DUMMYDetector.detector create mode 100755 NPSimulation/include/DummyDetector.hh create mode 100644 NPSimulation/include/ThinSiScorers.hh create mode 100644 NPSimulation/src/DummyDetector.cc create mode 100644 NPSimulation/src/ThinSiScorers.cc diff --git a/Inputs/DetectorConfiguration/DUMMYDetector.detector b/Inputs/DetectorConfiguration/DUMMYDetector.detector new file mode 100755 index 000000000..135a20792 --- /dev/null +++ b/Inputs/DetectorConfiguration/DUMMYDetector.detector @@ -0,0 +1,25 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +TheDUMMYDetector +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +DUMMYDetector + THETA= 0 + PHI= 0 + R= 150 + Thickness= 20 + Radius= 50 + Material= material1 + LeadThickness= 0 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +DUMMYDetector + THETA= 90 + PHI= 0 + R= 175 + Thickness= 100 + Radius= 50 + Material= material2 + LeadThickness= 0 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% Avaible material: +% material1 +% material2 diff --git a/NPSimulation/include/DummyDetector.hh b/NPSimulation/include/DummyDetector.hh new file mode 100755 index 000000000..700cec3a1 --- /dev/null +++ b/NPSimulation/include/DummyDetector.hh @@ -0,0 +1,128 @@ +#ifndef DUMMYDetector_h +#define DUMMYDetector_h 1 +/***************************************************************************** + * Copyright (C) 2009 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 2009 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class describe a simple Dummy Detector : * + * a simple cylinder of predifined material. user can choose to place it * + * where he want using spherical coordinate and choose between two naterial * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * This detector respect all the NPTool convention in naming volume, * + * reading scorers and file. Use it as a basis for your own detector * + *****************************************************************************/ +// C++ header +#include <string> +#include <vector> + +// G4 header defining G4 types +#include "globals.hh" + +// G4 headers +#include "G4ThreeVector.hh" +#include "G4RotationMatrix.hh" +#include "G4LogicalVolume.hh" +#include "G4SDManager.hh" +#include "G4MultiFunctionalDetector.hh" + +// NPTool header +#include "VDetector.hh" +#include "TDUMMYDetectorData.h" + +using namespace std; + +class DUMMYDetector : public VDetector +{ + //////////////////////////////////////////////////// + /////// Default Constructor and Destructor ///////// + //////////////////////////////////////////////////// +public: + DUMMYDetector() ; + virtual ~DUMMYDetector() ; + + //////////////////////////////////////////////////// + //////// Specific Function of this Class /////////// + //////////////////////////////////////////////////// +public: + // By Angle Method + void AddDUMMYDetector( G4double R , + G4double Theta , + G4double Phi , + G4double DUMMYDetectorThickness , + G4double DUMMYDetectorRadius , + G4String Scintillator ); + + void VolumeMaker(G4ThreeVector Det_pos, int DetNumber,G4LogicalVolume* world) ; + //////////////////////////////////////////////////// + ///////// Inherite from VDetector class /////////// + //////////////////////////////////////////////////// +public: + // Read stream at Configfile to pick-up parameters of detector (Position,...) + // Called in DetecorConstruction::ReadDetextorConfiguration Method + void ReadConfiguration(string Path) ; + + // Construct detector and inialise sensitive part. + // Called After DetecorConstruction::AddDetector Method + void ConstructDetector(G4LogicalVolume* world) ; + + // Add Detector branch to the EventTree. + // Called After DetecorConstruction::AddDetector Method + void InitializeRootOutput() ; + + // Read sensitive part and fill the Root tree. + // Called at in the EventAction::EndOfEventAvtion + void ReadSensitive(const G4Event* event) ; + +public: // Material + void InitializeMaterial() ; + // The dummy materials + G4Material* m_MaterialDUMMYDetector_material1 ; //the dummy material 2 + G4Material* m_MaterialDUMMYDetector_material2 ; //the dummy material 2 + + +public: // Scorer + // Initialize all Scorer used by the MUST2Array + void InitializeScorers() ; + + // Silicon Associate Scorer + G4MultiFunctionalDetector* m_DUMMYDetectorScorer ; + + + //////////////////////////////////////////////////// + ///////////Event class to store Data//////////////// + //////////////////////////////////////////////////// +private: + TDUMMYDetectorData* m_Event ; + + //////////////////////////////////////////////////// + ///////////////Private intern Data////////////////// + //////////////////////////////////////////////////// +private: + + vector<double> m_DUMMYDetectorThickness ; + vector<double> m_DUMMYDetectorRadius ; + + // Used for "By Angle Definition" + vector<G4double> m_R ; // | + vector<G4double> m_Theta ; // > Spherical coordinate DUMMYDetector volume center + vector<G4double> m_Phi ; // | + + // Scintillator type + vector<G4String> m_Scintillator ; + + +}; +#endif diff --git a/NPSimulation/include/Must2Scorers.hh b/NPSimulation/include/Must2Scorers.hh index 1ed9b019f..f5867f6ad 100644 --- a/NPSimulation/include/Must2Scorers.hh +++ b/NPSimulation/include/Must2Scorers.hh @@ -48,8 +48,9 @@ namespace MUST2 { private: G4double m_StripPlaneSize; G4int m_NumberOfStrip ; + G4double m_StripPitch ; G4int HCID; - G4THitsMap<G4double>* EvtMap; + G4THitsMap<G4int>* EvtMap; }; @@ -74,8 +75,9 @@ namespace MUST2 { private: G4double m_StripPlaneSize; G4int m_NumberOfStrip ; + G4double m_StripPitch ; G4int HCID; - G4THitsMap<G4double>* EvtMap; + G4THitsMap<G4int>* EvtMap; }; @@ -98,7 +100,7 @@ namespace MUST2 { private: G4int HCID; - G4THitsMap<G4double>* EvtMap; + G4THitsMap<G4int>* EvtMap; }; } diff --git a/NPSimulation/include/ThinSi.hh b/NPSimulation/include/ThinSi.hh index c0a01c935..48febf6c2 100644 --- a/NPSimulation/include/ThinSi.hh +++ b/NPSimulation/include/ThinSi.hh @@ -52,7 +52,7 @@ namespace THINSI const G4double FrameThickness = 3*mm ; const G4double SiliconSize = 50*mm ; const G4double AluThickness = 0.4*micrometer ; - const G4int NumberOfStrip = 32 ; + const G4int NumberOfStrip = 16 ; const G4double AluStripFront_PosZ = -0.5*SiliconThickness - 0.5*AluThickness ; const G4double Si_PosZ = 0 ; diff --git a/NPSimulation/include/ThinSiScorers.hh b/NPSimulation/include/ThinSiScorers.hh new file mode 100644 index 000000000..9449cffa8 --- /dev/null +++ b/NPSimulation/include/ThinSiScorers.hh @@ -0,0 +1,59 @@ +#ifndef ThinSiScorers_h +#define ThinSiScorers_h 1 +/***************************************************************************** + * Copyright (C) 2009 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 2009 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * File old the scorer specific to the ThinSi Detector * + * * + *---------------------------------------------------------------------------* + * Comment: * + * Those Scorer use TrackID as map index. This way ones can rebuild energy * + * deposit, time of flight or position,... particle by particle for each * + * event.Because standard scorer provide by G4 don't work this way but using* + * a global ID for each event you should not use those scorer with some G4 * + * provided ones or being very carefull doing so. * + *****************************************************************************/ +#include "G4VPrimitiveScorer.hh" +#include "G4THitsMap.hh" + +namespace THINSI { + + class PSStripNumber : public G4VPrimitiveScorer + { + + public: // with description + PSStripNumber(G4String name, G4int depth = 0, G4double StripPlaneSize = 50*mm, G4int NumberOfStrip = 16); + virtual ~PSStripNumber(); + + protected: // with description + virtual G4bool ProcessHits(G4Step*, G4TouchableHistory*); + + public: + virtual void Initialize(G4HCofThisEvent*); + virtual void EndOfEvent(G4HCofThisEvent*); + virtual void clear(); + virtual void DrawAll(); + virtual void PrintAll(); + + private: + G4double m_StripPlaneSize; + G4int m_NumberOfStrip ; + G4double m_StripPitch ; + G4int HCID; + G4THitsMap<G4int>* EvtMap; + }; + +} + +#endif diff --git a/NPSimulation/src/DummyDetector.cc b/NPSimulation/src/DummyDetector.cc new file mode 100644 index 000000000..0e4e17d85 --- /dev/null +++ b/NPSimulation/src/DummyDetector.cc @@ -0,0 +1,473 @@ +/***************************************************************************** + * Copyright (C) 2009 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 2009 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * This class describe a simple Dummy Detector : * + * a simple cylinder of predifined material. user can choose to place it * + * where he want using spherical coordinate and choose between two naterial * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * This detector respect all the NPTool convention in naming volume, * + * reading scorers and file. Use it as a basis for your own detector * + *****************************************************************************/ + +// C++ headers +#include <sstream> +#include <cmath> +#include <limits> +//G4 Geometry object +#include "G4Tubs.hh" + +//G4 sensitive +#include "G4SDManager.hh" +#include "G4MultiFunctionalDetector.hh" + +//G4 various object +#include "G4MaterialTable.hh" +#include "G4Element.hh" +#include "G4ElementTable.hh" +#include "G4Transform3D.hh" +#include "G4PVPlacement.hh" +#include "G4VisAttributes.hh" +#include "G4Colour.hh" + +// NPTool header +#include "DummyDetector.hh" +#include "GeneralScorers.hh" +#include "RootOutput.h" +using namespace GENERALSCORERS ; +// CLHEP header +#include "CLHEP/Random/RandGauss.h" + +using namespace std; +using namespace CLHEP; + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +namespace DUMMYDETECTOR +{ + // Energy and time Resolution + const G4double ResoTime = 4.2 ;// = 10ns of Resolution // Unit is MeV/2.35 + const G4double ResoEnergy = 5.0 ;// Resolution in % + +} + +using namespace DUMMYDETECTOR ; +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +// DUMMYDetector Specific Method +DUMMYDetector::DUMMYDetector() +{ + InitializeMaterial(); + m_Event = new TDUMMYDetectorData() ; +} + +DUMMYDetector::~DUMMYDetector() +{ + delete m_MaterialDUMMYDetector_material1 ; + delete m_MaterialDUMMYDetector_material2 ; + delete m_DUMMYDetectorScorer ; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +void DUMMYDetector::AddDUMMYDetector( G4double R , + G4double Theta , + G4double Phi , + G4double DUMMYDetectorThickness , + G4double DUMMYDetectorRadius , + G4String Scintillator ) +{ + + m_R.push_back(R) ; + m_Theta.push_back(Theta) ; + m_Phi.push_back(Phi) ; + m_DUMMYDetectorThickness.push_back(DUMMYDetectorThickness) ; + m_DUMMYDetectorRadius.push_back(DUMMYDetectorRadius) ; + m_Scintillator.push_back(Scintillator) ; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +// Virtual Method of VDetector class + + +// Read stream at Configfile to pick-up parameters of detector (Position,...) +// Called in DetecorConstruction::ReadDetextorConfiguration Method +void DUMMYDetector::ReadConfiguration(string Path) +{ + ifstream ConfigFile ; + ConfigFile.open(Path.c_str()) ; + string LineBuffer ; + string DataBuffer ; + + G4double Theta = 0 , Phi = 0 , R = 0 , Thickness = 0 , Radius = 0 ; + G4String Scintillator ; + + bool check_Theta = false ; + bool check_Phi = false ; + bool check_R = false ; + bool check_Thickness = false ; + bool check_Radius = false ;\ + bool check_Scintillator = false ; + bool ReadingStatus = false ; + + + while (!ConfigFile.eof()) + { + + getline(ConfigFile, LineBuffer); + + // If line is a Start Up DUMMYDetector bloc, Reading toggle to true + if (LineBuffer.compare(0, 13, "DUMMYDetector") == 0) + { + G4cout << "///" << G4endl ; + G4cout << "Dummy Module found: " << G4endl ; + ReadingStatus = true ; + + } + + // Else don't toggle to Reading Block Status + else ReadingStatus = false ; + + // Reading Block + while(ReadingStatus) + { + // Pickup Next Word + ConfigFile >> DataBuffer ; + + // Comment Line + if (DataBuffer.compare(0, 1, "%") == 0) { ConfigFile.ignore ( std::numeric_limits<std::streamsize>::max(), '\n' );} + + // Finding another telescope (safety), toggle out + else if (DataBuffer.compare(0, 13, "DUMMYDetector") == 0) { + cout << "WARNING: Another Detector is find before standard sequence of Token, Error may occured in Telecope definition" << endl ; + ReadingStatus = false ; + } + + //Angle method + else if (DataBuffer.compare(0, 6, "THETA=") == 0) { + check_Theta = true; + ConfigFile >> DataBuffer ; + Theta = atof(DataBuffer.c_str()) ; + Theta = Theta * deg; + cout << "Theta: " << Theta / deg << endl; + } + + else if (DataBuffer.compare(0, 4, "PHI=") == 0) { + check_Phi = true; + ConfigFile >> DataBuffer ; + Phi = atof(DataBuffer.c_str()) ; + Phi = Phi * deg; + cout << "Phi: " << Phi / deg << endl; + } + + else if (DataBuffer.compare(0, 2, "R=") == 0) { + check_R = true; + ConfigFile >> DataBuffer ; + R = atof(DataBuffer.c_str()) ; + R = R * mm; + cout << "R: " << R/mm << endl; + } + + else if (DataBuffer.compare(0, 7, "Radius=") == 0) { + check_Radius = true; + ConfigFile >> DataBuffer ; + Radius = atof(DataBuffer.c_str()) ; + Radius = Radius * mm; + cout << "DUMMYDetector Radius: " << Radius/mm << endl; + } + + else if (DataBuffer.compare(0, 10, "Thickness=") == 0) { + check_Thickness = true; + ConfigFile >> DataBuffer ; + Thickness = atof(DataBuffer.c_str()) ; + Thickness = Thickness * mm; + cout << "DUMMYDetector Thickness: " << Thickness/mm << endl; + } + + else if (DataBuffer.compare(0, 13, "Material=") == 0) { + check_Scintillator = true ; + ConfigFile >> DataBuffer ; + Scintillator = DataBuffer ; + cout << "DUMMYDetector material type: " << Scintillator << endl; + } + + /////////////////////////////////////////////////// + // If no Detector Token and no comment, toggle out + else + {ReadingStatus = false; G4cout << "Wrong Token Sequence: Getting out " << DataBuffer << G4endl ;} + + ///////////////////////////////////////////////// + // If All necessary information there, toggle out + + if ( check_Theta && check_Phi && check_R && check_Thickness && check_Radius && check_Scintillator) + { + AddDUMMYDetector( R , + Theta , + Phi , + Thickness , + Radius , + Scintillator ); + + // Reinitialisation of Check Boolean + + check_Theta = false ; + check_Phi = false ; + check_R = false ; + check_Thickness = false ; + check_Radius = false ; + check_Scintillator = false ; + ReadingStatus = false ; + cout << "///"<< endl ; + } + + } + } + +} + +// Construct detector and inialise sensitive part. +// Called After DetecorConstruction::AddDetector Method +void DUMMYDetector::ConstructDetector(G4LogicalVolume* world) +{ + G4ThreeVector Det_pos = G4ThreeVector(0, 0, 0) ; + + for (unsigned short i = 0 ; i < m_R.size() ; i++) + { + G4double wX = m_R[i] * sin(m_Theta[i] ) * cos(m_Phi[i] ) ; + G4double wY = m_R[i] * sin(m_Theta[i] ) * sin(m_Phi[i] ) ; + G4double wZ = m_R[i] * cos(m_Theta[i] ) ; + + Det_pos = G4ThreeVector(wX, wY, wZ) ; +// G4LogicalVolume* logicDUMMYDetector = NULL ; + + VolumeMaker(Det_pos , i+1, world) ; + } + +} + +void DUMMYDetector::VolumeMaker(G4ThreeVector Det_pos, int DetNumber, G4LogicalVolume* world) + { + //////////////////////////////////////////////////////////////// + ////////////// Starting Volume Definition ////////////////////// + //////////////////////////////////////////////////////////////// + G4PVPlacement* PVPBuffer ; + + // Name of the module + std::ostringstream DetectorNumber ; + DetectorNumber << DetNumber ; + G4String Name = "DUMMYDetector" + DetectorNumber.str() ; + + int i = DetNumber-1; + + G4Material* DUMMYDetectorMaterial ; + + if(m_Scintillator[i] == "material1" ) DUMMYDetectorMaterial = m_MaterialDUMMYDetector_material1 ; + else if(m_Scintillator[i] == "material2" ) DUMMYDetectorMaterial = m_MaterialDUMMYDetector_material2 ; + else { + G4cout << "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX" << endl ; + G4cout << "WARNING: Material Not found, default material set : material1" << endl ; + G4cout << "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX" << endl ; + DUMMYDetectorMaterial = m_MaterialDUMMYDetector_material1; + } + + + // Definition of the volume containing the sensitive detector + if(m_DUMMYDetectorThickness[i]>0 && m_DUMMYDetectorRadius[i]>0) + { + G4Tubs* solidDUMMYDetector = new G4Tubs( Name , + 0 , + m_DUMMYDetectorRadius[i] , + m_DUMMYDetectorThickness[i]/2 , + 0*deg , + 360*deg ); + + G4LogicalVolume* logicDUMMYDetector = new G4LogicalVolume(solidDUMMYDetector, DUMMYDetectorMaterial, Name+ "_Scintillator", 0, 0, 0); + logicDUMMYDetector->SetSensitiveDetector(m_DUMMYDetectorScorer); + + G4VisAttributes* PlastVisAtt = new G4VisAttributes(G4Colour(0.0, 0.0, 0.9)) ; + logicDUMMYDetector->SetVisAttributes(PlastVisAtt) ; + + + + PVPBuffer = new G4PVPlacement( 0 , + Det_pos , + logicDUMMYDetector , + Name + "_Scintillator" , + world , + false , + 0 ); + + + + } + } + +// Add Detector branch to the EventTree. +// Called After DetecorConstruction::AddDetector Method +void DUMMYDetector::InitializeRootOutput() +{ + RootOutput *pAnalysis = RootOutput::getInstance(); + TTree *pTree = pAnalysis->GetTree(); + pTree->Branch("DUMMYDetector", "TDUMMYDetectorData", &m_Event) ; +} + +// Read sensitive part and fill the Root tree. +// Called at in the EventAction::EndOfEventAvtion +void DUMMYDetector::ReadSensitive(const G4Event* event) +{ + G4String DetectorNumber ; + m_Event->Clear() ; + +////////////////////////////////////////////////////////////////////////////////////// +//////////////////////// Used to Read Event Map of detector ////////////////////////// +////////////////////////////////////////////////////////////////////////////////////// + + std::map<G4int, G4int*>::iterator DetectorNumber_itr ; + std::map<G4int, G4double*>::iterator Energy_itr ; + std::map<G4int, G4double*>::iterator Time_itr ; + + G4THitsMap<G4int>* DetectorNumberHitMap ; + G4THitsMap<G4double>* EnergyHitMap ; + G4THitsMap<G4double>* TimeHitMap ; + +////////////////////////////////////////////////////////////////////////////////////// +////////////////////////////////////////////////////////////////////////////////////// + + // Read the Scorer associate to the Silicon Strip + + //Detector Number + G4int StripDetCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("DUMMYDetectorScorer/DUMMYDetectorNumber") ; + DetectorNumberHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(StripDetCollectionID)) ; + DetectorNumber_itr = DetectorNumberHitMap->GetMap()->begin() ; + + //Energy + G4int StripEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("DUMMYDetectorScorer/Energy") ; + EnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripEnergyCollectionID)) ; + Energy_itr = EnergyHitMap->GetMap()->begin() ; + + //Time of Flight + G4int StripTimeCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("DUMMYDetectorScorer/Time") ; + TimeHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripTimeCollectionID)) ; + Time_itr = TimeHitMap->GetMap()->begin() ; + + G4int sizeN = DetectorNumberHitMap->entries() ; + G4int sizeE = EnergyHitMap->entries() ; + G4int sizeT = TimeHitMap->entries() ; + + // Loop on DUMMYDetector Number + for (G4int l = 0 ; l < sizeN ; l++) { + G4int N = *(DetectorNumber_itr->second) ; + G4int NTrackID = DetectorNumber_itr->first - N ; + + + if (N > 0) { + + m_Event->SetDUMMYDetectorNumber(N) ; + + // Energy + Energy_itr = EnergyHitMap->GetMap()->begin(); + for (G4int h = 0 ; h < sizeE ; h++) { + G4int ETrackID = Energy_itr->first - N ; + G4double E = *(Energy_itr->second) ; + + if (ETrackID == NTrackID) { + m_Event->SetEnergy(RandGauss::shoot(E, E*ResoEnergy/100./2.35)) ; + } + + Energy_itr++; + } + + + // Time + Time_itr = TimeHitMap->GetMap()->begin(); + for (G4int h = 0 ; h < sizeT ; h++) { + G4int TTrackID = Time_itr->first - N ; + G4double T = *(Time_itr->second) ; + + if (TTrackID == NTrackID) { + m_Event->SetTime(RandGauss::shoot(T, ResoTime)) ; + } + + Time_itr++; + } + + } + + DetectorNumber_itr++; + } + + // clear map for next event + TimeHitMap ->clear() ; + DetectorNumberHitMap ->clear() ; + EnergyHitMap ->clear() ; + +} + +//////////////////////////////////////////////////////////////// +void DUMMYDetector::InitializeMaterial() + { + + //////////////////////////////////////////////////////////////// + /////////////////Element Definition /////////////////////////// + //////////////////////////////////////////////////////////////// + G4String symbol ; + G4double density = 0. , a = 0, z = 0 ; + G4int ncomponents = 0, natoms = 0, fractionmass = 0 ; + + // for DUMMYDetector + G4Element* H = new G4Element("Hydrogen" , symbol = "H" , z = 1 , a = 1.01 * g / mole); + G4Element* C = new G4Element("Carbon" , symbol = "C" , z = 6 , a = 12.011 * g / mole); + G4Element* Pb = new G4Element("Lead" , symbol = "Pb" , z = 82 , a = 207.2 * g / mole); + //////////////////////////////////////////////////////////////// + /////////////////Material Definition /////////////////////////// + //////////////////////////////////////////////////////////////// + + + // DUMMYDetector BC-400 + density = 1.032 * g / cm3; + m_MaterialDUMMYDetector_material1 = new G4Material("DUMMYDetector_material1", density, ncomponents = 2); + m_MaterialDUMMYDetector_material1->AddElement(H , natoms = 10); + m_MaterialDUMMYDetector_material1->AddElement(C , natoms = 9); + + // DUMMYDetector BC-452 Pb 2% + density = 1.05 * g / cm3; + m_MaterialDUMMYDetector_material2 = new G4Material("DUMMYDetector_material2", density, ncomponents = 3); + m_MaterialDUMMYDetector_material2->AddElement(H , natoms = 10); + m_MaterialDUMMYDetector_material2->AddElement(C , natoms = 9); + m_MaterialDUMMYDetector_material2->AddElement(Pb , fractionmass=2*perCent); + + } + +//////////////////////////////////////////////////////////////// +void DUMMYDetector::InitializeScorers() + { + m_DUMMYDetectorScorer = new G4MultiFunctionalDetector("DUMMYDetectorScorer") ; + G4SDManager::GetSDMpointer()->AddNewDetector(m_DUMMYDetectorScorer); + + G4VPrimitiveScorer* DetNbr = new PSDetectorNumber("DUMMYDetectorNumber","DUMMYDetector", 0) ; + G4VPrimitiveScorer* Energy = new PSEnergy("Energy","DUMMYDetector", 0) ; + G4VPrimitiveScorer* Time = new PSTOF("Time","DUMMYDetector", 0) ; + + //and register it to the multifunctionnal detector + m_DUMMYDetectorScorer->RegisterPrimitive(DetNbr) ; + m_DUMMYDetectorScorer->RegisterPrimitive(Energy) ; + m_DUMMYDetectorScorer->RegisterPrimitive(Time) ; + + + } +//////////////////////////////////////////////////////////////// diff --git a/NPSimulation/src/MUST2Array.cc b/NPSimulation/src/MUST2Array.cc index 6a45a84bc..701f0275b 100644 --- a/NPSimulation/src/MUST2Array.cc +++ b/NPSimulation/src/MUST2Array.cc @@ -925,8 +925,8 @@ void MUST2Array::ReadSensitive(const G4Event* event) std::map<G4int, G4int*>::iterator DetectorNumber_itr ; std::map<G4int, G4double*>::iterator Energy_itr ; std::map<G4int, G4double*>::iterator Time_itr ; - std::map<G4int, G4double*>::iterator X_itr ; - std::map<G4int, G4double*>::iterator Y_itr ; + std::map<G4int, G4int*>::iterator X_itr ; + std::map<G4int, G4int*>::iterator Y_itr ; std::map<G4int, G4double*>::iterator Pos_X_itr ; std::map<G4int, G4double*>::iterator Pos_Y_itr ; std::map<G4int, G4double*>::iterator Pos_Z_itr ; @@ -936,8 +936,8 @@ void MUST2Array::ReadSensitive(const G4Event* event) G4THitsMap<G4int>* DetectorNumberHitMap ; G4THitsMap<G4double>* EnergyHitMap ; G4THitsMap<G4double>* TimeHitMap ; - G4THitsMap<G4double>* XHitMap ; - G4THitsMap<G4double>* YHitMap ; + G4THitsMap<G4int>* XHitMap ; + G4THitsMap<G4int>* YHitMap ; G4THitsMap<G4double>* PosXHitMap ; G4THitsMap<G4double>* PosYHitMap ; G4THitsMap<G4double>* PosZHitMap ; @@ -946,15 +946,15 @@ void MUST2Array::ReadSensitive(const G4Event* event) // Si(Li) std::map<G4int, G4double*>::iterator SiLiEnergy_itr ; - std::map<G4int, G4double*>::iterator SiLiPadNbr_itr ; + std::map<G4int, G4int*>::iterator SiLiPadNbr_itr ; G4THitsMap<G4double>* SiLiEnergyHitMap ; - G4THitsMap<G4double>* SiLiPadNbrHitMap ; + G4THitsMap<G4int>* SiLiPadNbrHitMap ; // CsI std::map<G4int, G4double*>::iterator CsIEnergy_itr ; - std::map<G4int, G4double*>::iterator CsICristalNbr_itr ; + std::map<G4int, G4int*>::iterator CsICristalNbr_itr ; G4THitsMap<G4double>* CsIEnergyHitMap ; - G4THitsMap<G4double>* CsICristalNbrHitMap ; + G4THitsMap<G4int>* CsICristalNbrHitMap ; ////////////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////////// @@ -977,12 +977,12 @@ void MUST2Array::ReadSensitive(const G4Event* event) //Strip Number X G4int StripXCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("MUST2_StripScorer/StripNumberX") ; - XHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripXCollectionID)) ; + XHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(StripXCollectionID)) ; X_itr = XHitMap->GetMap()->begin() ; //Strip Number Y G4int StripYCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("MUST2_StripScorer/StripNumberY") ; - YHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripYCollectionID)) ; + YHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(StripYCollectionID)) ; Y_itr = YHitMap->GetMap()->begin() ; //Interaction Coordinate X @@ -1017,7 +1017,7 @@ void MUST2Array::ReadSensitive(const G4Event* event) SiLiEnergy_itr = SiLiEnergyHitMap->GetMap()->begin() ; G4int SiLiPadNbrCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("MUST2_SiLiScorer/SiLiPadNbr") ; - SiLiPadNbrHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(SiLiPadNbrCollectionID)) ; + SiLiPadNbrHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(SiLiPadNbrCollectionID)) ; SiLiPadNbr_itr = SiLiPadNbrHitMap->GetMap()->begin() ; // Read the Scorer associate to the CsI crystal @@ -1027,7 +1027,7 @@ void MUST2Array::ReadSensitive(const G4Event* event) CsIEnergy_itr = CsIEnergyHitMap->GetMap()->begin() ; G4int CsICristalNbrCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("MUST2_CsIScorer/CsICristalNbr") ; - CsICristalNbrHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(CsICristalNbrCollectionID)) ; + CsICristalNbrHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(CsICristalNbrCollectionID)) ; CsICristalNbr_itr = CsICristalNbrHitMap->GetMap()->begin() ; @@ -1086,7 +1086,7 @@ void MUST2Array::ReadSensitive(const G4Event* event) X_itr = XHitMap->GetMap()->begin(); for (G4int h = 0 ; h < sizeX ; h++) { G4int XTrackID = X_itr->first - N ; - G4double X = *(X_itr->second) ; + G4int X = *(X_itr->second) ; if (XTrackID == NTrackID) { m_Event->SetMMStripXEStripNbr(X) ; m_Event->SetMMStripXTStripNbr(X) ; @@ -1099,7 +1099,7 @@ void MUST2Array::ReadSensitive(const G4Event* event) Y_itr = YHitMap->GetMap()->begin() ; for (G4int h = 0 ; h < sizeY ; h++) { G4int YTrackID = Y_itr->first - N ; - G4double Y = *(Y_itr->second) ; + G4int Y = *(Y_itr->second) ; if (YTrackID == NTrackID) { m_Event->SetMMStripYEStripNbr(Y) ; m_Event->SetMMStripYTStripNbr(Y) ; diff --git a/NPSimulation/src/Must2Scorers.cc b/NPSimulation/src/Must2Scorers.cc index cf93d01d3..2b057adf1 100644 --- a/NPSimulation/src/Must2Scorers.cc +++ b/NPSimulation/src/Must2Scorers.cc @@ -28,14 +28,6 @@ #include <string> using namespace MUST2 ; -//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -// Added by Adrien MATTA: -// Those Scorer use TrackID as map index. This way ones can rebuild energy deposit, -// time of flight or position,... particle by particle for each event. Because standard -// scorer provide by G4 don't work this way but using a global ID for each event you should -// not use those scorer with some G4 provided ones or being very carefull doing so. - //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... //Strip position Scorer @@ -43,8 +35,9 @@ using namespace MUST2 ; PSStripNumberX::PSStripNumberX(G4String name, G4int depth, G4double StripPlaneSize, G4int NumberOfStrip) : G4VPrimitiveScorer(name, depth), HCID(-1) { - m_StripPlaneSize = StripPlaneSize ; - m_NumberOfStrip = NumberOfStrip ; + m_StripPlaneSize = StripPlaneSize ; + m_NumberOfStrip = NumberOfStrip ; + m_StripPitch = m_StripPlaneSize / m_NumberOfStrip; } PSStripNumberX::~PSStripNumberX() @@ -59,12 +52,10 @@ G4bool PSStripNumberX::ProcessHits(G4Step* aStep, G4TouchableHistory*) G4ThreeVector POS = aStep->GetPreStepPoint()->GetPosition(); POS = aStep->GetPreStepPoint()->GetTouchableHandle()->GetHistory()->GetTopTransform().TransformPoint(POS); - G4double StripPitch = m_StripPlaneSize / m_NumberOfStrip; - - G4double temp = (POS(0) + m_StripPlaneSize / 2.) / StripPitch ; - G4double X = int(temp) + 1 ; + G4double temp = (POS(0) + m_StripPlaneSize / 2.) / m_StripPitch ; + G4int X = int(temp) + 1 ; //Rare case where particle is close to edge of silicon plan - if (X == 129) X = 128; + if (X == m_NumberOfStrip+1) X = m_NumberOfStrip; G4double edep = aStep->GetTotalEnergyDeposit(); if (edep < 100*keV) return FALSE; G4int index = aStep->GetTrack()->GetTrackID(); @@ -74,7 +65,7 @@ G4bool PSStripNumberX::ProcessHits(G4Step* aStep, G4TouchableHistory*) void PSStripNumberX::Initialize(G4HCofThisEvent* HCE) { - EvtMap = new G4THitsMap<G4double>(GetMultiFunctionalDetector()->GetName(), GetName()); + EvtMap = new G4THitsMap<G4int>(GetMultiFunctionalDetector()->GetName(), GetName()); if (HCID < 0) { HCID = GetCollectionID(0); } @@ -111,6 +102,7 @@ PSStripNumberY::PSStripNumberY(G4String name, G4int depth, G4double StripPlaneSi { m_StripPlaneSize = StripPlaneSize ; m_NumberOfStrip = NumberOfStrip ; + m_StripPitch = m_StripPlaneSize / m_NumberOfStrip; } PSStripNumberY::~PSStripNumberY() @@ -126,13 +118,11 @@ G4bool PSStripNumberY::ProcessHits(G4Step* aStep, G4TouchableHistory*) G4ThreeVector POS = aStep->GetPreStepPoint()->GetPosition(); POS = aStep->GetPreStepPoint()->GetTouchableHandle()->GetHistory()->GetTopTransform().TransformPoint(POS); - G4double StripPitch = m_StripPlaneSize / m_NumberOfStrip; - - G4double temp = (POS(1) + m_StripPlaneSize / 2.) / StripPitch ; + G4double temp = (POS(1) + m_StripPlaneSize / 2.) / m_StripPitch ; G4int temp2 = temp ; - G4double Y = temp2 + 1 ; + G4int Y = temp2 + 1 ; //Rare case where particle is close to edge of silicon plan - if (Y == 129) Y = 128; + if (Y == m_NumberOfStrip+1) Y = m_NumberOfStrip; G4double edep = aStep->GetTotalEnergyDeposit(); if (edep < 100*keV) return FALSE; @@ -143,7 +133,7 @@ G4bool PSStripNumberY::ProcessHits(G4Step* aStep, G4TouchableHistory*) void PSStripNumberY::Initialize(G4HCofThisEvent* HCE) { - EvtMap = new G4THitsMap<G4double>(GetMultiFunctionalDetector()->GetName(), GetName()); + EvtMap = new G4THitsMap<G4int>(GetMultiFunctionalDetector()->GetName(), GetName()); if (HCID < 0) { HCID = GetCollectionID(0); } @@ -192,10 +182,9 @@ G4bool PSPadOrCristalNumber::ProcessHits(G4Step* aStep, G4TouchableHistory*) { std::string name = aStep->GetTrack()->GetVolume()->GetName(); std::string nbr ; - unsigned int numberOfCharacterInDetectorNumber ; - int temp1,temp2 ; - double VolumeNumber; + G4int temp1,temp2 ; + G4int VolumeNumber; nbr = name[name.length()-1] ; temp1 = atoi( nbr.c_str() ) ; @@ -208,7 +197,7 @@ G4bool PSPadOrCristalNumber::ProcessHits(G4Step* aStep, G4TouchableHistory*) else { nbr= name[name.length()-1] ; VolumeNumber = atoi( nbr.c_str() ) ;} - int DetNbr = GENERALSCORERS::PickUpDetectorNumber(aStep, "MUST2Telescope"); + G4int DetNbr = GENERALSCORERS::PickUpDetectorNumber(aStep, "MUST2Telescope"); G4double edep = aStep->GetTotalEnergyDeposit(); if (edep < 100*keV) return FALSE; @@ -219,7 +208,7 @@ G4bool PSPadOrCristalNumber::ProcessHits(G4Step* aStep, G4TouchableHistory*) void PSPadOrCristalNumber::Initialize(G4HCofThisEvent* HCE) { - EvtMap = new G4THitsMap<G4double>(GetMultiFunctionalDetector()->GetName(), GetName()); + EvtMap = new G4THitsMap<G4int>(GetMultiFunctionalDetector()->GetName(), GetName()); if (HCID < 0) { HCID = GetCollectionID(0); } diff --git a/NPSimulation/src/ThinSi.cc b/NPSimulation/src/ThinSi.cc index 2aa8a5e37..748eadeda 100644 --- a/NPSimulation/src/ThinSi.cc +++ b/NPSimulation/src/ThinSi.cc @@ -46,9 +46,9 @@ // NPTool header #include "ThinSi.hh" #include "GeneralScorers.hh" -#include "Must2Scorers.hh" -#include "MUST2Array.hh" +#include "ThinSiScorers.hh" #include "RootOutput.h" +using namespace THINSI; // CLHEP header #include "CLHEP/Random/RandGauss.h" @@ -561,53 +561,72 @@ void ThinSi::InitializeRootOutput() // Called at in the EventAction::EndOfEventAvtion void ThinSi::ReadSensitive(const G4Event* event) { - G4String DetectorNumber ; m_Event->Clear(); ////////////////////////////////////////////////////////////////////////////////////// //////////////////////// Used to Read Event Map of detector ////////////////////////// ////////////////////////////////////////////////////////////////////////////////////// // Si - G4THitsMap<G4int>* DetNbrHitMap ; - G4THitsMap<G4double>* EnergyHitMap ; - G4THitsMap<G4double>* TimeHitMap ; - - std::map<G4int, G4int*>::iterator DetNbr_itr ; - std::map<G4int, G4double*>::iterator Energy_itr ; - std::map<G4int, G4double*>::iterator Time_itr ; + G4THitsMap<G4int>* DetNbrHitMap ; + G4THitsMap<G4int>* StripNbrHitMap ; + G4THitsMap<G4double>* EnergyHitMap ; + G4THitsMap<G4double>* TimeHitMap ; + + std::map<G4int, G4int*>::iterator DetNbr_itr ; + std::map<G4int, G4int*>::iterator StripNbr_itr ; + std::map<G4int, G4double*>::iterator Energy_itr ; + std::map<G4int, G4double*>::iterator Time_itr ; ////////////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////////// // Read the Scorer associate to the Silicon Strip //DetectorNumber - G4int StripDetNbrCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/DetectorNumber") ; - DetNbrHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(StripDetNbrCollectionID)) ; - DetNbr_itr = DetNbrHitMap->GetMap()->begin() ; + G4int DetNbrCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/DetectorNumber") ; + DetNbrHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(DetNbrCollectionID)) ; + DetNbr_itr = DetNbrHitMap->GetMap()->begin() ; + + //StripNumber + G4int StripNbrCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/StripNumber") ; + StripNbrHitMap = (G4THitsMap<G4int>*)(event->GetHCofThisEvent()->GetHC(StripNbrCollectionID)) ; //Energy - G4int StripEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/StripEnergy") ; - EnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripEnergyCollectionID)) ; - Energy_itr = EnergyHitMap->GetMap()->begin() ; + G4int StripEnergyCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/StripEnergy") ; + EnergyHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripEnergyCollectionID)) ; //Time - G4int StripTimeCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/StripTime") ; - TimeHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripTimeCollectionID)) ; - Time_itr = TimeHitMap->GetMap()->begin() ; + G4int StripTimeCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID("ThinSi_StripScorer/StripTime") ; + TimeHitMap = (G4THitsMap<G4double>*)(event->GetHCofThisEvent()->GetHC(StripTimeCollectionID)) ; - G4int sizeN = DetNbrHitMap->entries() ; - G4int sizeE = EnergyHitMap->entries() ; - G4int sizeT = TimeHitMap->entries() ; + G4int sizeN = DetNbrHitMap ->entries() ; + G4int sizeS = StripNbrHitMap ->entries() ; + G4int sizeE = EnergyHitMap ->entries() ; + G4int sizeT = TimeHitMap ->entries() ; - // Loop on Plastic Number + // Loop on Det Number for (G4int l = 0 ; l < sizeN ; l++) { G4int N = *(DetNbr_itr->second) ; - G4int NTrackID = DetNbr_itr->first - N ; - + G4int NTrackID = DetNbr_itr->first - N ; if (N > 0) { - + m_Event->SetStripEDetectorNbr(N) ; + m_Event->SetStripTDetectorNbr(N) ; + + // Strip Number + StripNbr_itr = StripNbrHitMap->GetMap()->begin(); + for (G4int h = 0 ; h < sizeS ; h++) { + G4int STrackID = StripNbr_itr->first - N ; + G4int S = *(StripNbr_itr->second) ; + + if (STrackID == NTrackID) { + m_Event->SetStripEStripNbr(S) ; + m_Event->SetStripTStripNbr(S) ; + } + + StripNbr_itr++; + } + // Energy Energy_itr = EnergyHitMap->GetMap()->begin(); for (G4int h = 0 ; h < sizeE ; h++) { @@ -615,8 +634,6 @@ void ThinSi::ReadSensitive(const G4Event* event) G4double E = *(Energy_itr->second) ; if (ETrackID == NTrackID) { - m_Event->SetStripEDetectorNbr(1) ; - m_Event->SetStripEStripNbr(1) ; m_Event->SetStripEEnergy( RandGauss::shoot(E, ResoEnergy ) ) ; } @@ -631,8 +648,6 @@ void ThinSi::ReadSensitive(const G4Event* event) G4double T = *(Time_itr->second) ; if (TTrackID == NTrackID) { - m_Event->SetStripTDetectorNbr(1) ; - m_Event->SetStripTStripNbr(1) ; m_Event->SetStripTTime( RandGauss::shoot(T, ResoTime ) ) ; } @@ -645,9 +660,11 @@ void ThinSi::ReadSensitive(const G4Event* event) } // clear map for next event + + DetNbrHitMap ->clear() ; + StripNbrHitMap ->clear() ; + EnergyHitMap ->clear() ; TimeHitMap ->clear() ; - DetNbrHitMap ->clear() ; - EnergyHitMap ->clear() ; } @@ -657,13 +674,15 @@ void ThinSi::InitializeScorers() // Silicon Associate Scorer m_StripScorer = new G4MultiFunctionalDetector("ThinSi_StripScorer"); - G4VPrimitiveScorer* DetNbr = new GENERALSCORERS::PSDetectorNumber("DetectorNumber","ThinSi_", 0) ; - G4VPrimitiveScorer* Energy = new GENERALSCORERS::PSEnergy("StripEnergy","ThinSi_", 0) ; - G4VPrimitiveScorer* TOF = new GENERALSCORERS::PSTOF("StripTime","ThinSi_", 0) ; - + G4VPrimitiveScorer* DetNbr = new GENERALSCORERS::PSDetectorNumber("DetectorNumber","ThinSi_", 0) ; + G4VPrimitiveScorer* StripNbr = new PSStripNumber("StripNumber",0,SiliconSize, NumberOfStrip) ; + G4VPrimitiveScorer* Energy = new GENERALSCORERS::PSEnergy("StripEnergy","ThinSi_", 0) ; + G4VPrimitiveScorer* TOF = new GENERALSCORERS::PSTOF("StripTime","ThinSi_", 0) ; + //and register it to the multifunctionnal detector m_StripScorer->RegisterPrimitive(DetNbr) ; + m_StripScorer->RegisterPrimitive(StripNbr) ; m_StripScorer->RegisterPrimitive(Energy) ; m_StripScorer->RegisterPrimitive(TOF) ; diff --git a/NPSimulation/src/ThinSiScorers.cc b/NPSimulation/src/ThinSiScorers.cc new file mode 100644 index 000000000..53faa61b0 --- /dev/null +++ b/NPSimulation/src/ThinSiScorers.cc @@ -0,0 +1,100 @@ +/***************************************************************************** + * Copyright (C) 2009 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 2009 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * File holding the scorer specific to the ThinSi Detector * + * * + *---------------------------------------------------------------------------* + * Comment: * + * Those Scorer use TrackID as map index. This way ones can rebuild energy * + * deposit, time of flight or position,... particle by particle for each * + * event.Because standard scorer provide by G4 don't work this way but using* + * a global ID for each event you should not use those scorer with some G4 * + * provided ones or being very carefull doing so. * + *****************************************************************************/ +#include "ThinSiScorers.hh" +#include "G4UnitsTable.hh" +#include "GeneralScorers.hh" +#include <string> +using namespace THINSI ; + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +//Strip position Scorer +//X +PSStripNumber::PSStripNumber(G4String name, G4int depth, G4double StripPlaneSize, G4int NumberOfStrip) + : G4VPrimitiveScorer(name, depth), HCID(-1) +{ + m_StripPlaneSize = StripPlaneSize ; + m_NumberOfStrip = NumberOfStrip ; + m_StripPitch = m_StripPlaneSize / m_NumberOfStrip ; +} + +PSStripNumber::~PSStripNumber() +{ + ; +} + +G4bool PSStripNumber::ProcessHits(G4Step* aStep, G4TouchableHistory*) +{ + int DetNbr = GENERALSCORERS::PickUpDetectorNumber(aStep, "ThinSi"); + + G4ThreeVector POS = aStep->GetPreStepPoint()->GetPosition(); + POS = aStep->GetPreStepPoint()->GetTouchableHandle()->GetHistory()->GetTopTransform().TransformPoint(POS); + + G4double temp = (POS(0) + m_StripPlaneSize / 2.) / m_StripPitch ; + G4double X = int(temp) + 1 ; + //Rare case where particle is close to edge of silicon plan + if (X == m_NumberOfStrip+1) X = m_NumberOfStrip; + G4double edep = aStep->GetTotalEnergyDeposit(); + if (edep < 100*keV) return FALSE; + G4int index = aStep->GetTrack()->GetTrackID(); + G4int S = X ; + EvtMap->set(index+DetNbr, S); + return TRUE; +} + +void PSStripNumber::Initialize(G4HCofThisEvent* HCE) +{ + EvtMap = new G4THitsMap<G4int>(GetMultiFunctionalDetector()->GetName(), GetName()); + if (HCID < 0) { + HCID = GetCollectionID(0); + } + HCE->AddHitsCollection(HCID, (G4VHitsCollection*)EvtMap); +} + +void PSStripNumber::EndOfEvent(G4HCofThisEvent*) +{ + ; +} + +void PSStripNumber::clear() +{ + EvtMap->clear(); +} + +void PSStripNumber::DrawAll() +{ + ; +} + +void PSStripNumber::PrintAll() +{ + G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl ; + G4cout << " PrimitiveScorer " << GetName() << G4endl ; + G4cout << " Number of entries " << EvtMap->entries() << G4endl ; +} + + + + -- GitLab