source: trunk/source/geometry/solids/specific/test/testG4Polyhedra.cc

Last change on this file was 1350, checked in by garnier, 15 years ago

update to last version 4.9.4
File size: 9.9 KB
Line 
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//
27// testG4Polyhedra.cc is derived from testG4Polycone.cc 13/10/2010
28//
29// Cvs version: $ Id $
30// Cvs tag : $ Name $
31
32#include "G4Timer.hh"
33#include <assert.h>
34#include <cmath>
35#include <fstream>
36#include <stdlib.h>
37#include "G4ios.hh"
38#include "G4Polyhedra.hh"
39
40#include "globals.hh"
41#include "G4Vector3D.hh"
42#include "G4Point3D.hh"
43#include "G4GeometryTolerance.hh"
44G4double kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
45int main(int, char **)
46{
47
48 double RMINVec[8];
49 RMINVec[0] = 30;
50 RMINVec[1] = 30;
51 RMINVec[2] = 0;
52 RMINVec[3] = 0;
53 RMINVec[4] = 0;
54 RMINVec[5] = 0;
55 RMINVec[6] = 40;
56 RMINVec[7] = 40;
57
58 double RMAXVec[8];
59 RMAXVec[0] = 70;
60 RMAXVec[1] = 70;
61 RMAXVec[2] = 70;
62 RMAXVec[3] = 40;
63 RMAXVec[4] = 40;
64 RMAXVec[5] = 80;
65 RMAXVec[6] = 80;
66 RMAXVec[7] = 60;
67
68 double Z_Values[8];
69 Z_Values[0] =-30;
70 Z_Values[1] =-20;
71 Z_Values[2] =-10;
72 Z_Values[3] = 0;
73 Z_Values[4] = 10;
74 Z_Values[5] = 20;
75 Z_Values[6] = 30;
76 Z_Values[7] = 40;
77
78 double Phi_Values[1];
79 Phi_Values[0]=0*deg;
80
81 Phi_Values[1]=45.*deg;
82 // Phi_Values[0]=0;
83 // Phi_Values[1]=2*pi;
84
85 G4cout << "\n======= Polyhedra test ========";
86
87 G4Polyhedra *MyPCone = new G4Polyhedra ("MyPCone",
88 Phi_Values[0],
89 Phi_Values[1],
90 5 ,
91 8 ,
92 Z_Values ,
93 RMINVec ,
94 RMAXVec );
95
96 G4cout << "\n\nPolyhedra is created ! "<<G4endl;
97 // -> Check methods :
98 // - Inside
99 // - DistanceToIn
100 // - DistanceToOut
101
102 EInside in;
103
104 G4cout<<"\n\n==================================================";
105 G4ThreeVector pt(0, -100, 24);
106 G4int y;
107 for (y = -100; y<=100; y+=10)
108 {
109 pt.setY(y);
110 in = MyPCone->Inside(pt);
111
112 G4cout << "\nx=" << pt.x() << " y=" << pt.y() << " z=" << pt.z();
113
114 if( in == kInside )
115 G4cout <<" is inside";
116 else
117 if( in == kOutside )
118 G4cout <<" is outside";
119 else
120 G4cout <<" is on the surface";
121 }
122
123 G4cout<<"\n\n==================================================";
124 G4ThreeVector start( 0, 0, -30);
125 G4ThreeVector dir(1./std::sqrt(2.), 1./std::sqrt(2.), 0);
126 G4double d;
127 G4int z;
128
129 G4cout<<"\nPdep is (0, 0, z)";
130 G4cout<<"\nDir is (1, 1, 0)\n";
131
132 for(z=-30; z<=50; z+=5)
133 {
134 start.setZ(z);
135
136 in = MyPCone->Inside(start);
137 G4cout<< "x=" << start.x() << " y=" << start.y() << " z=" << start.z();
138
139 if( in == kInside )
140 {
141 G4cout <<" is inside";
142
143 d = MyPCone->DistanceToOut(start, dir);
144 G4cout<<" distance to out="<<d;
145 d = MyPCone->DistanceToOut(start);
146 G4cout<<" closest distance to out="<<d<<G4endl;
147 }
148 else if( in == kOutside )
149 {
150 G4cout <<" is outside";
151
152 d = MyPCone->DistanceToIn(start, dir);
153 G4cout<<" distance to in="<<d;
154 d = MyPCone->DistanceToIn(start);
155 G4cout<<" closest distance to in="<<d<<G4endl;
156 }
157 else
158 G4cout <<" is on the surface"<<G4endl;
159
160 }
161
162 G4cout<<"\n\n==================================================";
163 G4ThreeVector start2( 0, -100, -30);
164 G4ThreeVector dir2(0, 1, 0);
165 G4double d2;
166
167 G4cout<<"\nPdep is (0, -100, z)";
168 G4cout<<"\nDir is (0, 1, 0)\n";
169
170 for(z=-30; z<=50; z+=5)
171 {
172 G4cout<<" z="<<z;
173 start2.setZ(z);
174 d2 = MyPCone->DistanceToIn(start2, dir2);
175 G4cout<<" distance to in="<<d2;
176 d2 = MyPCone->DistanceToIn(start2);
177 G4cout<<" distance to in="<<d2<<G4endl;
178 }
179
180 G4cout<<"\n\n==================================================";
181 G4ThreeVector start3( 0, 0, -50);
182 G4ThreeVector dir3(0, 0, 1);
183 G4double d3;
184
185 G4cout<<"\nPdep is (0, y, -50)";
186 G4cout<<"\nDir is (0, 0, 1)\n";
187
188 for(y=-0; y<=90; y+=5)
189 {
190 G4cout<<" y="<<y;
191 start3.setY(y);
192 d3 = MyPCone->DistanceToIn(start3, dir3);
193 G4cout<<" distance to in="<<d3<<G4endl;
194 }
195 //
196 // Add checks in Phi direction
197 // Point move in Phi direction for differents Z
198 //
199 G4cout<<"\n\n==================================================";
200
201 for(z=-10; z<=50; z+=5)
202{
203 G4cout<<"\n\n===================Z="<<z<<"==============================";
204
205 G4ThreeVector start4( 0, 0, z);
206 G4double phi=45.*pi/180.*rad;
207 G4ThreeVector dir4(std::cos(phi), std::sin(phi), 0);
208 G4double d4;
209
210 G4cout<<"\nPdep is (0<<R<<50, phi, z)";
211 G4cout<<"\nDir is (std::cos(phi), std::sin(phi), 0)\n";
212 G4cout<<"Ndirection is="<<dir4 <<G4endl;
213
214 for(y=-0; y<=50; y+=5)
215 {
216
217 start4.setX(y*std::cos(phi));
218 start4.setY(y*std::sin(phi));
219 G4cout<<" R="<<y<<" with Start"<<start4;
220 in = MyPCone->Inside(start4);
221 if( in == kInside )
222 {
223 G4cout <<" is inside";
224 d4 = MyPCone->DistanceToOut(start4, dir4);
225 G4cout<<" distance to out="<<d4;
226 d4 = MyPCone->DistanceToOut(start4);
227 G4cout<<" closest distance to out="<<d4<<G4endl;
228 }
229 else
230 if( in == kOutside )
231 {
232 G4cout <<" is outside";
233 d4 = MyPCone->DistanceToIn(start4, dir4);
234 G4cout<<" distance to in="<<d4;
235 d4 = MyPCone->DistanceToIn(start4);
236 G4cout<<" closest distance to in="<<d4<<G4endl;
237 }
238 else
239 {G4cout <<" is on the surface";
240 d4 = MyPCone->DistanceToIn(start4, dir4);
241 G4cout<<" distance to in="<<d4;
242 d4 = MyPCone->DistanceToIn(start4);
243 G4cout<<" closest distance to in="<<d4<<G4endl;
244 }
245
246 }
247}
248 //
249 // Add checks in Phi direction
250 // Point move in X direction for differents Z
251 // and 'schoot' on rhi edge
252 G4cout<<"\n\n==================================================";
253
254 for(z=-10; z<=50; z+=5)
255 {
256 G4cout<<"\n\n===================Z="<<z<<"==============================";
257
258 G4ThreeVector start5( 0., 1, z);
259 G4ThreeVector dir5(0,-1, 0);
260 G4double d5;
261
262 G4cout<<"\nPdep is (0<<X<<50, 1, z)";
263 G4cout<<"\nDir is (0, -1, 0)\n";
264 G4cout<<"Ndirection is="<<dir5 <<G4endl;
265
266 for(y=-0; y<=50; y+=5)
267 {
268
269 start5.setX(y);
270 G4cout<<" Start"<<start5;
271 in = MyPCone->Inside(start5);
272 if( in == kInside )
273 {
274 G4cout <<" is inside";
275 d5 = MyPCone->DistanceToOut(start5, dir5);
276 G4cout<<" distance to out="<<d5;
277 d5 = MyPCone->DistanceToOut(start5);
278 G4cout<<" closest distance to out="<<d5<<G4endl;
279 }
280 else
281 if( in == kOutside )
282 {
283 G4cout <<" is outside";
284 d5 = MyPCone->DistanceToIn(start5, dir5);
285 G4cout<<" distance to in="<<d5;
286 d5 = MyPCone->DistanceToIn(start5);
287 G4cout<<" closest distance to in="<<d5<<G4endl;
288 }
289 else
290 {
291 G4cout <<" is on the surface";
292 d5 = MyPCone->DistanceToIn(start5, dir5);
293 G4cout<<" distance to in="<<d5;
294 d5 = MyPCone->DistanceToIn(start5);
295 G4cout<<" closest distance to in="<<d5<<G4endl;
296 }
297
298 }
299 }
300
301
302 // Asserts
303 G4ThreeVector p1,p2,p3,p4,p5,p6,dirx,diry,dirz;
304 p1=G4ThreeVector(0,0,-5);
305 p2=G4ThreeVector(50,0,40);
306 p3=G4ThreeVector(5,1,20 );
307 p4=G4ThreeVector(45,5,30);
308 p5=G4ThreeVector(0,0,30);
309 p6=G4ThreeVector(41,0,10);
310
311 dirx=G4ThreeVector(1,0,0);
312 diry=G4ThreeVector(0,1,0);
313 dirz=G4ThreeVector(0,0,1);
314
315 // Inside
316 assert(MyPCone->Inside(p1) == kSurface);
317 assert(MyPCone->Inside(p2) == kSurface);
318 assert(MyPCone->Inside(p3) == kInside);
319 assert(MyPCone->Inside(p4) == kInside);
320 assert(MyPCone->Inside(p5) == kOutside);
321 assert(MyPCone->Inside(p6) == kOutside);
322 // DistanceToIn
323
324 assert(std::fabs((MyPCone->DistanceToIn(p1,dirx))) < kCarTolerance);
325 assert(std::fabs((MyPCone->DistanceToIn(p1,-diry)))< kCarTolerance);
326 assert(std::fabs((MyPCone->DistanceToIn(p2,diry))) < kCarTolerance);
327 assert(std::fabs((MyPCone->DistanceToIn(p5,dirx) -40.12368793931)) < kCarTolerance);
328 assert(std::fabs((MyPCone->DistanceToIn(p6,-dirx) -0.87631206069)) < kCarTolerance);
329 assert(std::fabs((MyPCone->DistanceToIn(p6,dirz) -0.218402670765))< kCarTolerance);
330 // DistanceToOut
331 assert(std::fabs((MyPCone->DistanceToOut (p1,-dirx))) < kCarTolerance);
332 assert(std::fabs((MyPCone->DistanceToOut (p3,-diry) -1.) ) < kCarTolerance);
333 assert(std::fabs((MyPCone->DistanceToOut (p3,dirz) -1.27382374146))< kCarTolerance);
334 assert(std::fabs((MyPCone->DistanceToOut(p4,dirz) -10.)) < kCarTolerance);
335 assert(std::fabs((MyPCone->DistanceToOut(p4,dirx) -34.8538673445)) < kCarTolerance);
336 assert(std::fabs((MyPCone->DistanceToOut(p4,diry) -40.)) < kCarTolerance);
337
338
339
340 return EXIT_SUCCESS;
341}
342
Note: See TracBrowser for help on using the repository browser.