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

#include <RE02DetectorConstruction.hh>

Inheritance diagram for RE02DetectorConstruction:
Collaboration diagram for RE02DetectorConstruction:

Public Member Functions

 RE02DetectorConstruction ()
 
virtual ~RE02DetectorConstruction ()
 
virtual G4VPhysicalVolumeConstruct ()
 
virtual void ConstructSDandField ()
 
void SetPhantomSize (G4ThreeVector size)
 
const G4ThreeVectorGetPhantomSize () const
 
void SetNumberOfSegmentsInPhantom (G4int nx, G4int ny, G4int nz)
 
void GetNumberOfSegmentsInPhantom (G4int &nx, G4int &ny, G4int &nz) const
 
void SetLeadSegment (G4bool flag=TRUE)
 
G4bool IsLeadSegment ()
 
- 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
 

Additional Inherited Members

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

Detailed Description

Uer detector construction class

(Description)

Detector construction for example RE02.

[Geometry] The world volume is defined as 200 cm x 200 cm x 200 cm box with Air. Water phantom is defined as 200 mm x 200 mm x 400 mm box with Water. The water phantom is divided into 100 segments in x,y plane using replication, and then divided into 200 segments perpendicular to z axis using nested parameterised volume. These values are defined at constructor, e.g. the size of water phantom (fPhantomSize), and number of segmentation of water phantom (fNx, fNy, fNz).

By default, lead plates are inserted into the position of even order segments. NIST database is used for materials.

[Scorer] Assignment of G4MultiFunctionalDetector and G4PrimitiveScorer

is demonstrated in this example.

The collection names of defined Primitives are 0 PhantomSD/totalEDep 1 PhantomSD/protonEDep 2 PhantomSD/protonNStep 3 PhantomSD/chargedPassCellFlux 4 PhantomSD/chargedCellFlux 5 PhantomSD/chargedSurfFlux 6 PhantomSD/gammaSurfCurr000 7 PhantomSD/gammaSurfCurr001 9 PhantomSD/gammaSurdCurr002

10 PhantomSD/gammaSurdCurr003

Please see README for detail description.

Definition at line 111 of file RE02DetectorConstruction.hh.

Constructor & Destructor Documentation

RE02DetectorConstruction::RE02DetectorConstruction ( )

Definition at line 106 of file RE02DetectorConstruction.cc.

108 {
109  // Default size of water phantom,and segmentation.
110  fPhantomSize.setX(200.*mm);
111  fPhantomSize.setY(200.*mm);
112  fPhantomSize.setZ(400.*mm);
113  fNx = fNy = fNz = 100;
114  fInsertLead = TRUE;
115 }
static constexpr double mm
Definition: G4SIunits.hh:115
void setY(double)
void setZ(double)
void setX(double)
#define TRUE
Definition: globals.hh:55

Here is the call graph for this function:

RE02DetectorConstruction::~RE02DetectorConstruction ( )
virtual

Definition at line 118 of file RE02DetectorConstruction.cc.

119 {;}

Member Function Documentation

G4VPhysicalVolume * RE02DetectorConstruction::Construct ( void  )
virtual

Implements G4VUserDetectorConstruction.

Definition at line 122 of file RE02DetectorConstruction.cc.

123 {
124  //=====================
125  // Material Definitions
126  //=====================
127  //
128  //-------- NIST Materials ----------------------------------------------------
129  // Material Information imported from NIST database.
130  //
132  G4Material* air = NISTman->FindOrBuildMaterial("G4_AIR");
133  G4Material* water = NISTman->FindOrBuildMaterial("G4_WATER");
134  G4Material* lead = NISTman->FindOrBuildMaterial("G4_Pb");
135 
136  //
137  // Print all the materials defined.
138  G4cout << G4endl << "The materials defined are : " << G4endl << G4endl;
139  G4cout << *(G4Material::GetMaterialTable()) << G4endl;
140 
141  //============================================================================
142  // Definitions of Solids, Logical Volumes, Physical Volumes
143  //============================================================================
144 
145  //-------------
146  // World Volume
147  //-------------
148 
149  G4ThreeVector worldSize = G4ThreeVector(200*cm, 200*cm, 200*cm);
150 
151  G4Box * solidWorld
152  = new G4Box("world", worldSize.x()/2., worldSize.y()/2., worldSize.z()/2.);
153  G4LogicalVolume * logicWorld
154  = new G4LogicalVolume(solidWorld, air, "World", 0, 0, 0);
155 
156  //
157  // Must place the World Physical volume unrotated at (0,0,0).
158  G4VPhysicalVolume * physiWorld
159  = new G4PVPlacement(0, // no rotation
160  G4ThreeVector(), // at (0,0,0)
161  logicWorld, // its logical volume
162  "World", // its name
163  0, // its mother volume
164  false, // no boolean operations
165  0); // copy number
166 
167  //---------------
168  // Water Phantom
169  //---------------
170 
171  //................................
172  // Mother Volume of Water Phantom
173  //................................
174 
175  //-- Default size of water phantom is defined at constructor.
176  G4ThreeVector phantomSize = fPhantomSize;
177 
178  G4Box * solidPhantom
179  = new G4Box("phantom",
180  phantomSize.x()/2., phantomSize.y()/2., phantomSize.z()/2.);
181  G4LogicalVolume * logicPhantom
182  = new G4LogicalVolume(solidPhantom, water, "Phantom", 0, 0, 0);
183 
184  G4RotationMatrix* rot = new G4RotationMatrix();
185  //rot->rotateY(30.*deg);
186  G4ThreeVector positionPhantom;
187  //G4VPhysicalVolume * physiPhantom =
188  new G4PVPlacement(rot, // no rotation
189  positionPhantom, // at (x,y,z)
190  logicPhantom, // its logical volume
191  "Phantom", // its name
192  logicWorld, // its mother volume
193  false, // no boolean operations
194  0); // copy number
195 
196  //..............................................
197  // Phantom segmentation using Parameterisation
198  //..............................................
199  //
200  G4cout << "<-- RE02DetectorConstruction::Construct-------" <<G4endl;
201  G4cout << " Water Phantom Size " << fPhantomSize/mm << G4endl;
202  G4cout << " Segmentation ("<< fNx<<","<<fNy<<","<<fNz<<")"<< G4endl;
203  G4cout << " Lead plate at even copy # (0-False,1-True): " << IsLeadSegment()
204  << G4endl;
205  G4cout << "<---------------------------------------------"<< G4endl;
206  // Number of segmentation.
207  // - Default number of segmentation is defined at constructor.
208  G4int nxCells = fNx;
209  G4int nyCells = fNy;
210  G4int nzCells = fNz;
211 
212  G4ThreeVector sensSize;
213  sensSize.setX(phantomSize.x()/(G4double)nxCells);
214  sensSize.setY(phantomSize.y()/(G4double)nyCells);
215  sensSize.setZ(phantomSize.z()/(G4double)nzCells);
216  // i.e Voxel size will be 2.0 x 2.0 x 2.0 mm3 cube by default.
217  //
218 
219  // Replication of Water Phantom Volume.
220  // Y Slice
221  G4String yRepName("RepY");
222  G4VSolid* solYRep =
223  new G4Box(yRepName,phantomSize.x()/2.,sensSize.y()/2.,phantomSize.z()/2.);
224  G4LogicalVolume* logYRep =
225  new G4LogicalVolume(solYRep,water,yRepName);
226  //G4PVReplica* yReplica =
227  new G4PVReplica(yRepName,logYRep,logicPhantom,kYAxis,fNy,sensSize.y());
228  // X Slice
229  G4String xRepName("RepX");
230  G4VSolid* solXRep =
231  new G4Box(xRepName,sensSize.x()/2.,sensSize.y()/2.,phantomSize.z()/2.);
232  G4LogicalVolume* logXRep =
233  new G4LogicalVolume(solXRep,water,xRepName);
234  //G4PVReplica* xReplica =
235  new G4PVReplica(xRepName,logXRep,logYRep,kXAxis,fNx,sensSize.x());
236 
237  //
238  //..................................
239  // Voxel solid and logical volumes
240  //..................................
241  // Z Slice
242  G4String zVoxName("phantomSens");
243  G4VSolid* solVoxel =
244  new G4Box(zVoxName,sensSize.x()/2.,sensSize.y()/2.,sensSize.z()/2.);
245  fLVPhantomSens = new G4LogicalVolume(solVoxel,water,zVoxName);
246  //
247  //
248  std::vector<G4Material*> phantomMat(2,water);
249  if ( IsLeadSegment() ) phantomMat[1]=lead;
250  //
251  // Parameterisation for transformation of voxels.
252  // (voxel size is fixed in this example.
253  // e.g. nested parameterisation handles material and transfomation of voxels.)
255  = new RE02NestedPhantomParameterisation(sensSize/2.,nzCells,phantomMat);
256  //G4VPhysicalVolume * physiPhantomSens =
257  new G4PVParameterised("PhantomSens", // their name
258  fLVPhantomSens, // their logical volume
259  logXRep, // Mother logical volume
260  kUndefined, // Are placed along this axis
261  nzCells, // Number of cells
262  paramPhantom); // Parameterisation.
263  // Optimization flag is avaiable for,
264  // kUndefined, kXAxis, kYAxis, kZAxis.
265  //
266 
267  //===============================
268  // Visualization attributes
269  //===============================
270 
271  G4VisAttributes* boxVisAtt = new G4VisAttributes(G4Colour(1.0,1.0,1.0));
272  logicWorld ->SetVisAttributes(boxVisAtt);
273  //logicWorld->SetVisAttributes(G4VisAttributes::GetInvisible());
274 
275  // Mother volume of WaterPhantom
276  G4VisAttributes* phantomVisAtt = new G4VisAttributes(G4Colour(1.0,1.0,0.0));
277  logicPhantom->SetVisAttributes(phantomVisAtt);
278 
279  // Replica
280  G4VisAttributes* yRepVisAtt = new G4VisAttributes(G4Colour(0.0,1.0,0.0));
281  logYRep->SetVisAttributes(yRepVisAtt);
282  G4VisAttributes* xRepVisAtt = new G4VisAttributes(G4Colour(0.0,1.0,0.0));
283  logXRep->SetVisAttributes(xRepVisAtt);
284 
285  // Skip the visualization for those voxels.
287 
288 
289  return physiWorld;
290 }
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
static constexpr double mm
Definition: G4SIunits.hh:115
CLHEP::Hep3Vector G4ThreeVector
double x() const
CLHEP::HepRotation G4RotationMatrix
Definition: G4Box.hh:64
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
int G4int
Definition: G4Types.hh:78
void setY(double)
static G4NistManager * Instance()
double z() const
void setZ(double)
void setX(double)
G4GLOB_DLL std::ostream G4cout
static constexpr double cm
Definition: G4SIunits.hh:119
double y() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static const G4VisAttributes & GetInvisible()
void SetVisAttributes(const G4VisAttributes *pVA)

Here is the call graph for this function:

void RE02DetectorConstruction::ConstructSDandField ( )
virtual

Reimplemented from G4VUserDetectorConstruction.

Definition at line 292 of file RE02DetectorConstruction.cc.

292  {
293 
294  //================================================
295  // Sensitive detectors : MultiFunctionalDetector
296  //================================================
297  //
298  // Sensitive Detector Manager.
300  //
301  // Sensitive Detector Name
302  G4String phantomSDname = "PhantomSD";
303 
304  //------------------------
305  // MultiFunctionalDetector
306  //------------------------
307  //
308  // Define MultiFunctionalDetector with name.
310  = new G4MultiFunctionalDetector(phantomSDname);
311  pSDman->AddNewDetector( mFDet ); // Register SD to SDManager.
312  fLVPhantomSens->SetSensitiveDetector(mFDet); // Assign SD to the logical volume.
313 
314  //---------------------------------------
315  // SDFilter : Sensitive Detector Filters
316  //---------------------------------------
317  //
318  // Particle Filter for Primitive Scorer with filter name(fltName)
319  // and particle name(particleName),
320  // or particle names are given by add("particle name"); method.
321  //
322  G4String fltName,particleName;
323  //
324  //-- proton filter
325  G4SDParticleFilter* protonFilter =
326  new G4SDParticleFilter(fltName="protonFilter", particleName="proton");
327  //
328  //-- electron filter
329  G4SDParticleFilter* electronFilter =
330  new G4SDParticleFilter(fltName="electronFilter");
331  electronFilter->add(particleName="e+"); // accept electrons.
332  electronFilter->add(particleName="e-"); // accept positorons.
333  //
334  //-- charged particle filter
335  G4SDChargedFilter* chargedFilter =
336  new G4SDChargedFilter(fltName="chargedFilter");
337 
338  //------------------------
339  // PS : Primitive Scorers
340  //------------------------
341  // Primitive Scorers are used with SDFilters according to your purpose.
342  //
343  //
344  //-- Primitive Scorer for Energy Deposit.
345  // Total, by protons, by electrons.
346  G4String psName;
347  G4PSEnergyDeposit3D * scorer0 = new G4PSEnergyDeposit3D(psName="totalEDep",
348  fNx,fNy,fNz);
349  G4PSEnergyDeposit3D * scorer1 = new G4PSEnergyDeposit3D(psName="protonEDep",
350  fNx,fNy,fNz);
351  scorer1->SetFilter(protonFilter);
352 
353  //
354  //-- Number of Steps for protons
355  G4PSNofStep3D * scorer2 =
356  new G4PSNofStep3D(psName="protonNStep",fNx,fNy,fNz);
357  scorer2->SetFilter(protonFilter);
358 
359  //
360  //-- CellFlux for charged particles
361  G4PSPassageCellFlux3D * scorer3 =
362  new G4PSPassageCellFlux3D(psName="chargedPassCellFlux", fNx,fNy,fNz);
363  G4PSCellFlux3D * scorer4 =
364  new G4PSCellFlux3D(psName="chargedCellFlux", fNx,fNy,fNz);
365  G4PSFlatSurfaceFlux3D * scorer5 =
366  new G4PSFlatSurfaceFlux3D(psName="chargedSurfFlux", fFlux_InOut,fNx,fNy,fNz);
367  scorer3->SetFilter(chargedFilter);
368  scorer4->SetFilter(chargedFilter);
369  scorer5->SetFilter(chargedFilter);
370 
371  //
372  //------------------------------------------------------------
373  // Register primitive scorers to MultiFunctionalDetector
374  //------------------------------------------------------------
375  mFDet->RegisterPrimitive(scorer0);
376  mFDet->RegisterPrimitive(scorer1);
377  mFDet->RegisterPrimitive(scorer2);
378  mFDet->RegisterPrimitive(scorer3);
379  mFDet->RegisterPrimitive(scorer4);
380  mFDet->RegisterPrimitive(scorer5);
381 
382  //========================
383  // More additional Primitive Scoreres
384  //========================
385  //
386  //--- Surface Current for gamma with energy bin.
387  // This example creates four primitive scorers.
388  // 4 bins with energy --- Primitive Scorer Name
389  // 1. to 10 KeV, gammaSurfCurr000
390  // 10 keV to 100 KeV, gammaSurfCurr001
391  // 100 keV to 1 MeV, gammaSurfCurr002
392  // 1 MeV to 10 MeV. gammaSurfCurr003
393  //
394  char name[17];
395  for ( G4int i = 0; i < 4; i++){
396  std::sprintf(name,"gammaSurfCurr%03d",i);
397  G4String psgName(name);
398  G4double kmin = std::pow(10.,(G4double)i)*keV;
399  G4double kmax = std::pow(10.,(G4double)(i+1))*keV;
400  //-- Particle with kinetic energy filter.
401  G4SDParticleWithEnergyFilter* pkinEFilter =
402  new G4SDParticleWithEnergyFilter(fltName="gammaE filter",kmin,kmax);
403  pkinEFilter->add("gamma"); // Accept only gamma.
404  pkinEFilter->show(); // Show accepting condition to stdout.
405  //-- Surface Current Scorer which scores number of tracks in unit area.
406  G4PSFlatSurfaceCurrent3D * scorer =
407  new G4PSFlatSurfaceCurrent3D(psgName,fCurrent_InOut,fNx,fNy,fNz);
408  scorer->SetFilter(pkinEFilter); // Assign filter.
409  mFDet->RegisterPrimitive(scorer); // Register it to MultiFunctionalDetector.
410  }
411 
412 }
const XML_Char * name
Definition: expat.h:151
G4bool RegisterPrimitive(G4VPrimitiveScorer *)
void SetFilter(G4VSDFilter *f)
int G4int
Definition: G4Types.hh:78
void AddNewDetector(G4VSensitiveDetector *aSD)
Definition: G4SDManager.cc:71
void add(const G4String &particleName)
static G4SDManager * GetSDMpointer()
Definition: G4SDManager.cc:40
double G4double
Definition: G4Types.hh:76
static constexpr double keV
Definition: G4SIunits.hh:216
void add(const G4String &particleName)
void SetSensitiveDetector(G4VSensitiveDetector *pSDetector)

Here is the call graph for this function:

void RE02DetectorConstruction::GetNumberOfSegmentsInPhantom ( G4int nx,
G4int ny,
G4int nz 
) const
inline

Definition at line 131 of file RE02DetectorConstruction.hh.

132  { nx=fNx; ny = fNy; nz = fNz; }

Here is the caller graph for this function:

const G4ThreeVector& RE02DetectorConstruction::GetPhantomSize ( ) const
inline

Definition at line 127 of file RE02DetectorConstruction.hh.

127 { return fPhantomSize; }
G4bool RE02DetectorConstruction::IsLeadSegment ( )
inline

Definition at line 135 of file RE02DetectorConstruction.hh.

135 { return fInsertLead; }

Here is the caller graph for this function:

void RE02DetectorConstruction::SetLeadSegment ( G4bool  flag = TRUE)
inline

Definition at line 134 of file RE02DetectorConstruction.hh.

134 { fInsertLead = flag; }

Here is the caller graph for this function:

void RE02DetectorConstruction::SetNumberOfSegmentsInPhantom ( G4int  nx,
G4int  ny,
G4int  nz 
)
inline

Definition at line 129 of file RE02DetectorConstruction.hh.

130  { fNx=nx; fNy=ny; fNz=nz; }

Here is the caller graph for this function:

void RE02DetectorConstruction::SetPhantomSize ( G4ThreeVector  size)
inline

Definition at line 126 of file RE02DetectorConstruction.hh.

126 { fPhantomSize=size; }

Here is the caller graph for this function:


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