BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSSamplerCylinder.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#include "BDSExtent.hh"
20#include "BDSSamplerCylinder.hh"
21#include "BDSSDManager.hh"
22#include "BDSSDSamplerCylinder.hh"
23#include "BDSUtilities.hh"
24
25#include "G4CutTubs.hh"
26#include "G4LogicalVolume.hh"
27#include "G4String.hh"
28#include "G4ThreeVector.hh"
29#include "G4Tubs.hh"
30#include "G4Types.hh"
31
32#include "CLHEP/Units/SystemOfUnits.h"
33
34
35BDSSamplerCylinder::BDSSamplerCylinder(const G4String& nameIn,
36 G4double radiusIn,
37 G4double fullLength,
38 G4double startAngle,
39 G4double sweepAngle,
40 G4int filterSetIDIn):
41 BDSSampler(nameIn, filterSetIDIn)
42{
43 G4double thickness = 1e-6 * radiusIn;
44 containerSolid = new G4Tubs(nameIn + "_solid", // name
45 radiusIn, // inner radius
46 radiusIn + thickness, // outer radius
47 fullLength*0.5, // half-length
48 startAngle, // start angle
49 sweepAngle); // sweep angle
50
51 SetExtent(BDSExtent(radiusIn, radiusIn, fullLength*0.5));
53}
54
55BDSSamplerCylinder::BDSSamplerCylinder(const G4String& nameIn,
56 G4double radiusIn,
57 G4double fullLength,
58 const G4ThreeVector& inputFaceNormal,
59 const G4ThreeVector& outputFaceNormal,
60 G4double startAngle,
61 G4double sweepAngle,
62 G4int filterSetIDIn):
63 BDSSampler(nameIn, filterSetIDIn)
64{
65 G4double thickness = 1e-6 * radiusIn;
66 if (!BDS::IsFinite(inputFaceNormal.angle({0,0,-1}),1e-6) && !BDS::IsFinite(outputFaceNormal.angle({0,0,1}),1e-6))
67 {
68 containerSolid = new G4Tubs(nameIn + "_solid", // name
69 radiusIn, // inner radius
70 radiusIn + thickness, // outer radius
71 fullLength*0.5, // half-length
72 startAngle, // start angle
73 sweepAngle); // sweep angle
74 }
75 else
76 {
77 containerSolid = new G4CutTubs(nameIn + "_solid", // name
78 radiusIn, // inner radius
79 radiusIn + thickness, // outer radius
80 fullLength * 0.5, // half-length
81 startAngle, // start angle
82 sweepAngle, // sweep angle
83 inputFaceNormal,
84 outputFaceNormal);
85 }
86 SetExtent(BDSExtent(radiusIn, radiusIn, fullLength*0.5));
88}
89
91{
92 auto sdMan = BDSSDManager::Instance();
93 BDSSDSamplerCylinder* sd = filterSetID > -1 ? sdMan->SamplerCylinderWithFilter(filterSetID) : sdMan->SamplerCylinder();
94 containerLogicalVolume->SetSensitiveDetector(sd);
95}
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
void SetExtent(const BDSExtent &extIn)
Set extent.
The sensitive detector class that provides sensitivity to BDSSamplerCylinder instances.
virtual void SetSensitivity()
Attach sensitive detector to containerLogicalVolume.
Base class and registry of sampler instances.
Definition: BDSSampler.hh:36
virtual void CommonConstruction()
Definition: BDSSampler.cc:36
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())