| Geant4 Cross Reference |
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // G4TwistedTrd
27 //
28 // Author: 18/03/2005 - O.Link (Oliver.Link@cern.ch)
29 // --------------------------------------------------------------------
30
31 #include "G4TwistedTrd.hh"
32 #include "G4SystemOfUnits.hh"
33 #include "G4Polyhedron.hh"
34
35 //=====================================================================
36 //* Constructor -------------------------------------------------------
37
38 G4TwistedTrd::G4TwistedTrd( const G4String& pName,
39 G4double pDx1,
40 G4double pDx2,
41 G4double pDy1,
42 G4double pDy2,
43 G4double pDz,
44 G4double pPhiTwist )
45 : G4VTwistedFaceted( pName, pPhiTwist,pDz,0.,0.,
46 pDy1, pDx1, pDx1, pDy2, pDx2, pDx2,0.)
47 {
48 }
49
50 //=====================================================================
51 // Fake default constructor - sets only member data and allocates memory
52 // for usage restricted to object persistency.
53
54 G4TwistedTrd::G4TwistedTrd( __void__& a )
55 : G4VTwistedFaceted(a)
56 {
57 }
58
59 //=====================================================================
60 //* Destructor --------------------------------------------------------
61
62 G4TwistedTrd::~G4TwistedTrd() = default;
63
64 //=====================================================================
65 //* Copy constructor --------------------------------------------------
66
67 G4TwistedTrd::G4TwistedTrd(const G4TwistedTrd& rhs)
68 : G4VTwistedFaceted(rhs)
69 {
70 fpPolyhedron = GetPolyhedron();
71 }
72
73 //=====================================================================
74 //* Assignment operator -----------------------------------------------
75
76 G4TwistedTrd& G4TwistedTrd::operator = (const G4TwistedTrd& rhs)
77 {
78 // Check assignment to self
79 //
80 if (this == &rhs) { return *this; }
81
82 // Copy base class data
83 //
84 G4VTwistedFaceted::operator=(rhs);
85 fpPolyhedron = GetPolyhedron();
86
87 return *this;
88 }
89
90 //=====================================================================
91 //* StreamInfo --------------------------------------------------------
92
93 std::ostream& G4TwistedTrd::StreamInfo(std::ostream& os) const
94 {
95 //
96 // Stream object contents to an output stream
97 //
98 os << "-----------------------------------------------------------\n"
99 << " *** Dump for solid - " << GetName() << " ***\n"
100 << " ===================================================\n"
101 << " Solid type: G4TwistedTrd\n"
102 << " Parameters: \n"
103 << " pDx1 = " << GetX1HalfLength()/cm << " cm" << G4endl
104 << " pDx2 = " << GetX2HalfLength()/cm << " cm" << G4endl
105 << " pDy1 = " << GetY1HalfLength()/cm << " cm" << G4endl
106 << " pDy2 = " << GetY2HalfLength()/cm << " cm" << G4endl
107 << " pDz = " << GetZHalfLength()/cm << " cm" << G4endl
108 << " pPhiTwist = " << GetPhiTwist()/degree << " deg" << G4endl
109 << "-----------------------------------------------------------\n";
110
111 return os;
112 }
113
114 //=====================================================================
115 //* GetEntityType -----------------------------------------------------
116
117 G4GeometryType G4TwistedTrd::GetEntityType() const
118 {
119 return {"G4TwistedTrd"};
120 }
121
122 //=====================================================================
123 //* Clone -------------------------------------------------------------
124
125 G4VSolid* G4TwistedTrd::Clone() const
126 {
127 return new G4TwistedTrd(*this);
128 }
129
130 //=====================================================================
131 //* GetCubicVolume ----------------------------------------------------
132
133 double G4TwistedTrd::GetCubicVolume()
134 {
135 if (fCubicVolume == 0.)
136 {
137 G4double x1 = GetX1HalfLength();
138 G4double x2 = GetX2HalfLength();
139 G4double y1 = GetY1HalfLength();
140 G4double y2 = GetY2HalfLength();
141 G4double h = 2.*GetZHalfLength();
142 fCubicVolume = h*((x1 + x2)*(y1 + y2) + (x2 - x1)*(y2 - y1)/3.);
143 }
144 return fCubicVolume;
145 }
146
147 //=====================================================================
148 //* GetSurfaceArea ----------------------------------------------------
149
150 double G4TwistedTrd::GetSurfaceArea()
151 {
152 if (fSurfaceArea == 0.)
153 {
154 G4double ang = GetPhiTwist();
155 G4double x1 = GetX1HalfLength();
156 G4double x2 = GetX2HalfLength();
157 G4double y1 = GetY1HalfLength();
158 G4double y2 = GetY2HalfLength();
159 G4double h = 2.*GetZHalfLength();
160 G4double hh = h*h;
161 G4double delX = x2 - x1;
162 G4double delY = y2 - y1;
163 if (ang == 0.)
164 {
165 G4double hx = std::sqrt(delY*delY + hh);
166 G4double hy = std::sqrt(delX*delX + hh);
167 return fSurfaceArea =
168 2.*(x1 + x2)*hx + 2.*(y1 + y2)*hy + 4.*(x1*y1 + x2*y2);
169 }
170
171 // compute area of x-faces
172 G4double U1, U2, V1, V2;
173 G4double areaX = 0.;
174 U1 = delY + x1*ang;
175 U2 = delY + x2*ang;
176 V1 = delY - x1*ang;
177 V2 = delY - x2*ang;
178 if (std::abs(delX) < kCarTolerance) // case x1 == x2
179 {
180 areaX = (U1*std::sqrt(hh + U1*U1) + hh*std::asinh(U1/h) -
181 V1*std::sqrt(hh + V1*V1) - hh*std::asinh(V1/h))/ang;
182 }
183 else
184 {
185 // U contribution
186 areaX += ((hh + U2*U2)*std::sqrt(hh + U2*U2) -
187 (hh + U1*U1)*std::sqrt(hh + U1*U1))/3.
188 + hh*(U2*std::asinh(U2/h) - U1*std::asinh(U1/h))
189 - hh*(std::sqrt(hh + U2*U2) - std::sqrt(hh + U1*U1));
190 // V contribution
191 areaX += ((hh + V2*V2)*std::sqrt(hh + V2*V2) -
192 (hh + V1*V1)*std::sqrt(hh + V1*V1))/3.
193 + hh*(V2*std::asinh(V2/h) - V1*std::asinh(V1/h))
194 - hh*(std::sqrt(hh + V2*V2) - std::sqrt(hh + V1*V1));
195 areaX /= delX*ang*ang;
196 }
197
198 // compute area of y-faces
199 G4double areaY = 0.;
200 U1 = delX + y1*ang;
201 U2 = delX + y2*ang;
202 V1 = delX - y1*ang;
203 V2 = delX - y2*ang;
204 if (std::abs(delY) < kCarTolerance) // case y1 == y2
205 {
206 areaY = (U1*std::sqrt(hh + U1*U1) + hh*std::asinh(U1/h) -
207 V1*std::sqrt(hh + V1*V1) - hh*std::asinh(V1/h))/ang;
208 }
209 else
210 {
211 // U contribution
212 areaY += ((hh + U2*U2)*std::sqrt(hh + U2*U2) -
213 (hh + U1*U1)*std::sqrt(hh + U1*U1))/3.
214 + hh*(U2*std::asinh(U2/h) - U1*std::asinh(U1/h))
215 - hh*(std::sqrt(hh + U2*U2) - std::sqrt(hh + U1*U1));
216 // V contribution
217 areaY += ((hh + V2*V2)*std::sqrt(hh + V2*V2) -
218 (hh + V1*V1)*std::sqrt(hh + V1*V1))/3.
219 + hh*(V2*std::asinh(V2/h) - V1*std::asinh(V1/h))
220 - hh*(std::sqrt(hh + V2*V2) - std::sqrt(hh + V1*V1));
221 areaY /= delY*ang*ang;
222 }
223 fSurfaceArea = areaX + areaY + 4.*(x1*y1 + x2*y2);
224 }
225 return fSurfaceArea;
226 }
227