Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
TSDetectorConstruction Class Reference

#include <TSDetectorConstruction.hh>

Inheritance diagram for TSDetectorConstruction:
Collaboration diagram for TSDetectorConstruction:

Public Types

typedef std::map< G4String,
G4Material * > 
MaterialCollection_t
 
typedef std::set
< G4LogicalVolume * > 
ScoringVolumes_t
 

Public Member Functions

 TSDetectorConstruction ()
 
virtual ~TSDetectorConstruction ()
 
G4VPhysicalVolumeConstruct ()
 
const G4ThreeVectorGetWorldDimensions () const
 
const ScoringVolumes_tGetScoringVolumes () const
 
const G4StringGetMFDName () const
 
- Public Member Functions inherited from G4VUserDetectorConstruction
 G4VUserDetectorConstruction ()
 
virtual ~G4VUserDetectorConstruction ()
 
virtual void CloneSD ()
 
virtual void CloneF ()
 
void RegisterParallelWorld (G4VUserParallelWorld *)
 
G4int ConstructParallelGeometries ()
 
void ConstructParallelSD ()
 
G4int GetNumberOfParallelWorld () const
 
G4VUserParallelWorldGetParallelWorld (G4int i) const
 

Static Public Member Functions

static TSDetectorConstructionInstance ()
 

Protected Member Functions

virtual MaterialCollection_t ConstructMaterials ()
 
virtual G4VPhysicalVolumeConstructWorld (const MaterialCollection_t &)
 
virtual void ConstructSDandField ()
 
- Protected Member Functions inherited from G4VUserDetectorConstruction
void SetSensitiveDetector (const G4String &logVolName, G4VSensitiveDetector *aSD, G4bool multi=false)
 
void SetSensitiveDetector (G4LogicalVolume *logVol, G4VSensitiveDetector *aSD)
 

Detailed Description

Definition at line 74 of file TSDetectorConstruction.hh.

Member Typedef Documentation

Constructor & Destructor Documentation

TSDetectorConstruction::TSDetectorConstruction ( )

Definition at line 92 of file TSDetectorConstruction.cc.

93 : fWorldPhys(0),
94  fWorldMaterialName("G4_Galactic"),
95  fTargetMaterialName("G4_B"),
96  fCasingMaterialName("G4_WATER"),
97  fWorldDim(G4ThreeVector(0.5*m, 0.5*m, 0.5*m)),
98  fTargetDim(G4ThreeVector(0.5*m, 0.5*m, 0.5*m)),
99  fTargetSections(G4ThreeVector(5, 5, 5)),
100  fMfdName("Target_MFD")
101 {
102  fgInstance = this;
103 }
CLHEP::Hep3Vector G4ThreeVector
static constexpr double m
Definition: G4SIunits.hh:129
TSDetectorConstruction::~TSDetectorConstruction ( )
virtual

Definition at line 107 of file TSDetectorConstruction.cc.

108 {
109  fgInstance = 0;
110 }

Member Function Documentation

G4VPhysicalVolume * TSDetectorConstruction::Construct ( void  )
virtual

Implements G4VUserDetectorConstruction.

Definition at line 114 of file TSDetectorConstruction.cc.

115 {
117 }
virtual MaterialCollection_t ConstructMaterials()
virtual G4VPhysicalVolume * ConstructWorld(const MaterialCollection_t &)

Here is the call graph for this function:

TSDetectorConstruction::MaterialCollection_t TSDetectorConstruction::ConstructMaterials ( )
protectedvirtual

Definition at line 122 of file TSDetectorConstruction.cc.

123 {
124  MaterialCollection_t materials;
126 
127  materials["World"] = nist->FindOrBuildMaterial(fWorldMaterialName);
128  materials["Target"] = nist->FindOrBuildMaterial(fTargetMaterialName);
129  materials["Casing"] = nist->FindOrBuildMaterial(fCasingMaterialName);
130 
131  return materials;
132 }
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
static G4NistManager * Instance()
std::map< G4String, G4Material * > MaterialCollection_t

Here is the call graph for this function:

Here is the caller graph for this function:

void TSDetectorConstruction::ConstructSDandField ( )
protectedvirtual

Reimplemented from G4VUserDetectorConstruction.

Definition at line 277 of file TSDetectorConstruction.cc.

278 {
279  //------------------------------------------------//
280  // MultiFunctionalDetector //
281  //------------------------------------------------//
282  // Define MultiFunctionalDetector with name.
285  G4VPrimitiveScorer* edep = new G4PSEnergyDeposit("EnergyDeposit");
286  MFDet->RegisterPrimitive(edep);
287  G4VPrimitiveScorer* nstep = new G4PSNofStep("NumberOfSteps");
288  MFDet->RegisterPrimitive(nstep);
289 
290  // add scoring volumes
291  for(auto ite : fScoringVolumes)
292  {
293  SetSensitiveDetector(ite, MFDet);
294  }
295 }
G4bool RegisterPrimitive(G4VPrimitiveScorer *)
void SetSensitiveDetector(const G4String &logVolName, G4VSensitiveDetector *aSD, G4bool multi=false)
void AddNewDetector(G4VSensitiveDetector *aSD)
Definition: G4SDManager.cc:71
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40

Here is the call graph for this function:

G4VPhysicalVolume * TSDetectorConstruction::ConstructWorld ( const MaterialCollection_t materials)
protectedvirtual

Definition at line 137 of file TSDetectorConstruction.cc.

138 {
139  G4UserLimits* steplimit = new G4UserLimits(0.1*(fTargetDim.z()
140  /fTargetSections.z()));
141  G4bool check_overlap = true;
142 
143  G4Box* world_solid = new G4Box("World",
144  0.5*fWorldDim.x(),
145  0.5*fWorldDim.y(),
146  0.5*fWorldDim.z());
147  G4LogicalVolume* world_log = new G4LogicalVolume(world_solid,
148  materials.find("World")
149  ->second,
150  "World");
151  fWorldPhys = new G4PVPlacement(0,
152  G4ThreeVector(0.),
153  "World",
154  world_log,
155  0, false, 0, check_overlap);
156 
157  G4int nz = fTargetSections.z();
158  G4int ny = fTargetSections.y();
159  G4int nx = fTargetSections.x();
160 
161  // spacing between sections
162  G4double sx = fTargetDim.x()/fTargetSections.x();
163  G4double sy = fTargetDim.y()/fTargetSections.y();
164  G4double sz = fTargetDim.z()/fTargetSections.z();
165 
166  //G4cout << "World has dimensions : "
167  //<< G4BestUnit(fWorldDim, "Length") << G4endl;
168 
169  //------------------------------------------------------------------------//
170  // Set Visual Attributes
171  //------------------------------------------------------------------------//
172  G4VisAttributes* red = new G4VisAttributes(G4Color(1., 0., 0., 1.0));
173  G4VisAttributes* green = new G4VisAttributes(G4Color(0., 1., 0., 0.25));
174  G4VisAttributes* blue = new G4VisAttributes(G4Color(0., 0., 1., 0.1));
175  G4VisAttributes* white = new G4VisAttributes(G4Color(1., 1., 1., 1.));
176 
177  white->SetVisibility(true);
178  red->SetVisibility(true);
179  green->SetVisibility(true);
180  blue->SetVisibility(true);
181 
182  white->SetForceWireframe(true);
183  red->SetForceSolid(true);
184  green->SetForceSolid(true);
185  blue->SetForceSolid(true);
186 
187  world_log->SetVisAttributes(white);
188 
189  for(G4int k = 0; k < nz; ++k)
190  for(G4int j = 0; j < ny; ++j)
191  for(G4int i = 0; i < nx; ++i)
192  {
193  // displacement of section
194  G4double dx
195  = 0.5*sx + static_cast<G4double>(i)*sx - 0.5*fWorldDim.x();
196  G4double dy
197  = 0.5*sy + static_cast<G4double>(j)*sy - 0.5*fWorldDim.y();
198  G4double dz
199  = 0.5*sz + static_cast<G4double>(k)*sz - 0.5*fWorldDim.z();
200  G4ThreeVector td = G4ThreeVector(dx, dy, -dz);
201  // make unique name
202  std::stringstream ss_name;
203  ss_name << "Target_" << i << "_" << j << "_" << k;
204 
205  G4Box* target_solid = new G4Box(ss_name.str(),
206  0.5*sx,
207  0.5*sy,
208  0.5*sz);
209 
210 
211  G4Material* target_material = 0;
212  G4bool is_casing = true;
213 
214  if(j == 0 || j+1 == ny || i == 0 || i+1 == nx ||
215  (nz > 1 && (k == 0 || k+1 == nz)))
216  target_material = materials.find("Casing")->second;
217  else {
218  target_material = materials.find("Target")->second;
219  is_casing = false;
220  }
221 
222  G4LogicalVolume* target_log =
223  new G4LogicalVolume(target_solid,
224  target_material,
225  ss_name.str());
226 
227  target_log->SetUserLimits(steplimit);
228 
229  new G4PVPlacement(0,
230  td,
231  ss_name.str(),
232  target_log,
233  fWorldPhys,
234  true,
235  k*nx*ny + j*nx + i,
236  check_overlap);
237 
238  fScoringVolumes.insert(target_log);
239 
240  if(is_casing)
241  target_log->SetVisAttributes(blue);
242  else
243  {
244  // making a checkerboard for kicks...
245  G4bool even_z = (k%2 == 0) ? true : false;
246  G4bool even_y = (j%2 == 0) ? true : false;
247  G4bool even_x = (i%2 == 0) ? true : false;
248 
249  G4VisAttributes* theColor = nullptr;
250 
251  if((even_z))
252  {
253  if((even_y && even_x) || (!even_y && !even_x))
254  theColor = red;
255  else
256  theColor = green;
257  } else // ! even_z
258  {
259  if((!even_y && even_x) || (even_y && !even_x))
260  theColor = red;
261  else
262  theColor = green;
263  }
264 
265  target_log->SetVisAttributes(theColor);
266  }
267 
268  }
269 
270 
271 
272  return fWorldPhys;
273 }
Definition: test07.cc:36
CLHEP::Hep3Vector G4ThreeVector
double x() const
Definition: G4Box.hh:64
void SetUserLimits(G4UserLimits *pULimits)
Definition: test07.cc:36
int G4int
Definition: G4Types.hh:78
double z() const
void SetForceSolid(G4bool=true)
void SetVisibility(G4bool=true)
bool G4bool
Definition: G4Types.hh:79
G4Colour G4Color
Definition: G4Color.hh:42
double y() const
double G4double
Definition: G4Types.hh:76
void SetForceWireframe(G4bool=true)
void SetVisAttributes(const G4VisAttributes *pVA)

Here is the call graph for this function:

Here is the caller graph for this function:

const G4String& TSDetectorConstruction::GetMFDName ( ) const
inline

Definition at line 91 of file TSDetectorConstruction.hh.

91 { return fMfdName; }
const ScoringVolumes_t& TSDetectorConstruction::GetScoringVolumes ( ) const
inline

Definition at line 89 of file TSDetectorConstruction.hh.

90  { return fScoringVolumes; }
const G4ThreeVector& TSDetectorConstruction::GetWorldDimensions ( ) const
inline

Definition at line 88 of file TSDetectorConstruction.hh.

88 { return fWorldDim; }

Here is the caller graph for this function:

TSDetectorConstruction * TSDetectorConstruction::Instance ( void  )
static

Definition at line 85 of file TSDetectorConstruction.cc.

86 {
87  return fgInstance;
88 }

Here is the caller graph for this function:


The documentation for this class was generated from the following files: