| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // G4ReduciblePolygon implementation; a utilit << 26 //
27 // test, reduce, and/or otherwise manipulate a << 27 // $Id: G4ReduciblePolygon.cc 92024 2015-08-13 14:16:00Z gcosmo $
>> 28 //
>> 29 //
>> 30 // --------------------------------------------------------------------
>> 31 // GEANT 4 class source file
>> 32 //
>> 33 //
>> 34 // G4ReduciblePolygon.cc
>> 35 //
>> 36 // Implementation of a utility class used to specify, test, reduce,
>> 37 // and/or otherwise manipulate a 2D polygon.
>> 38 //
28 // See G4ReduciblePolygon.hh for more info. 39 // See G4ReduciblePolygon.hh for more info.
29 // 40 //
30 // Author: David C. Williams (davidw@scipp.ucs <<
31 // ------------------------------------------- 41 // --------------------------------------------------------------------
32 42
33 #include "G4ReduciblePolygon.hh" 43 #include "G4ReduciblePolygon.hh"
34 #include "globals.hh" 44 #include "globals.hh"
35 45
>> 46 //
36 // Constructor: with simple arrays 47 // Constructor: with simple arrays
37 // 48 //
38 G4ReduciblePolygon::G4ReduciblePolygon( const 49 G4ReduciblePolygon::G4ReduciblePolygon( const G4double a[],
39 const 50 const G4double b[],
40 51 G4int n )
41 : aMin(0.), aMax(0.), bMin(0.), bMax(0.) << 52 : aMin(0.), aMax(0.), bMin(0.), bMax(0.),
>> 53 vertexHead(0)
42 { 54 {
43 // 55 //
44 // Do all of the real work in Create 56 // Do all of the real work in Create
45 // 57 //
46 Create( a, b, n ); 58 Create( a, b, n );
47 } 59 }
48 60
>> 61
>> 62 //
49 // Constructor: special PGON/PCON case 63 // Constructor: special PGON/PCON case
50 // 64 //
51 G4ReduciblePolygon::G4ReduciblePolygon( const 65 G4ReduciblePolygon::G4ReduciblePolygon( const G4double rmin[],
52 const 66 const G4double rmax[],
53 const 67 const G4double z[], G4int n )
54 : aMin(0.), aMax(0.), bMin(0.), bMax(0.) << 68 : aMin(0.), aMax(0.), bMin(0.), bMax(0.),
>> 69 vertexHead(0)
55 { 70 {
56 // 71 //
57 // Translate 72 // Translate
58 // 73 //
59 auto a = new G4double[n*2]; << 74 G4double *a = new G4double[n*2];
60 auto b = new G4double[n*2]; << 75 G4double *b = new G4double[n*2];
61 76
62 G4double *rOut = a + n, 77 G4double *rOut = a + n,
63 *zOut = b + n, 78 *zOut = b + n,
64 *rIn = rOut-1, 79 *rIn = rOut-1,
65 *zIn = zOut-1; 80 *zIn = zOut-1;
66 81
67 for( G4int i=0; i < n; ++i, ++rOut, ++zOut, << 82 G4int i;
>> 83 for( i=0; i < n; i++, rOut++, zOut++, rIn--, zIn-- )
68 { 84 {
69 *rOut = rmax[i]; 85 *rOut = rmax[i];
70 *rIn = rmin[i]; 86 *rIn = rmin[i];
71 *zOut = *zIn = z[i]; 87 *zOut = *zIn = z[i];
72 } 88 }
73 89
74 Create( a, b, n*2 ); 90 Create( a, b, n*2 );
75 91
76 delete [] a; 92 delete [] a;
77 delete [] b; 93 delete [] b;
78 } 94 }
79 95
>> 96
>> 97 //
80 // Create 98 // Create
81 // 99 //
82 // To be called by constructors, fill in the l 100 // To be called by constructors, fill in the list and statistics for a new
83 // polygon 101 // polygon
84 // 102 //
85 void G4ReduciblePolygon::Create( const G4doubl 103 void G4ReduciblePolygon::Create( const G4double a[],
86 const G4doubl 104 const G4double b[], G4int n )
87 { 105 {
88 if (n<3) 106 if (n<3)
89 G4Exception("G4ReduciblePolygon::Create()", 107 G4Exception("G4ReduciblePolygon::Create()", "GeomSolids0002",
90 FatalErrorInArgument, "Less tha 108 FatalErrorInArgument, "Less than 3 vertices specified.");
91 109
92 const G4double *anext = a, *bnext = b; 110 const G4double *anext = a, *bnext = b;
93 ABVertex* prev = nullptr; << 111 ABVertex *prev = 0;
94 do // Loop checking, 13.08.2015, G.Cosmo 112 do // Loop checking, 13.08.2015, G.Cosmo
95 { 113 {
96 auto newVertex = new ABVertex; << 114 ABVertex *newVertex = new ABVertex;
97 newVertex->a = *anext; 115 newVertex->a = *anext;
98 newVertex->b = *bnext; 116 newVertex->b = *bnext;
99 newVertex->next = nullptr; << 117 newVertex->next = 0;
100 if (prev==nullptr) << 118 if (prev==0)
101 { 119 {
102 vertexHead = newVertex; 120 vertexHead = newVertex;
103 } 121 }
104 else 122 else
105 { 123 {
106 prev->next = newVertex; 124 prev->next = newVertex;
107 } 125 }
108 126
109 prev = newVertex; 127 prev = newVertex;
110 } while( ++anext, ++bnext < b+n ); 128 } while( ++anext, ++bnext < b+n );
111 129
112 numVertices = n; 130 numVertices = n;
113 131
114 CalculateMaxMin(); 132 CalculateMaxMin();
115 } 133 }
116 134
>> 135
>> 136 //
117 // Fake default constructor - sets only member 137 // Fake default constructor - sets only member data and allocates memory
118 // for usage restri 138 // for usage restricted to object persistency.
119 // 139 //
120 G4ReduciblePolygon::G4ReduciblePolygon( __void 140 G4ReduciblePolygon::G4ReduciblePolygon( __void__& )
121 : aMin(0.), aMax(0.), bMin(0.), bMax(0.) << 141 : aMin(0.), aMax(0.), bMin(0.), bMax(0.), numVertices(0), vertexHead(0)
122 { 142 {
123 } 143 }
124 144
125 145
126 // 146 //
127 // Destructor 147 // Destructor
128 // 148 //
129 G4ReduciblePolygon::~G4ReduciblePolygon() 149 G4ReduciblePolygon::~G4ReduciblePolygon()
130 { 150 {
131 ABVertex* curr = vertexHead; << 151 ABVertex *curr = vertexHead;
132 while( curr != nullptr ) // Loop checking << 152 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
133 { 153 {
134 ABVertex* toDelete = curr; << 154 ABVertex *toDelete = curr;
135 curr = curr->next; 155 curr = curr->next;
136 delete toDelete; 156 delete toDelete;
137 } 157 }
138 } 158 }
139 159
>> 160
>> 161 //
140 // CopyVertices 162 // CopyVertices
141 // 163 //
142 // Copy contents into simple linear arrays. 164 // Copy contents into simple linear arrays.
143 // ***** CAUTION ***** Be care to declare the 165 // ***** CAUTION ***** Be care to declare the arrays to a large
144 // enough size! 166 // enough size!
145 // 167 //
146 void G4ReduciblePolygon::CopyVertices( G4doubl 168 void G4ReduciblePolygon::CopyVertices( G4double a[], G4double b[] ) const
147 { 169 {
148 G4double *anext = a, *bnext = b; 170 G4double *anext = a, *bnext = b;
149 ABVertex *curr = vertexHead; 171 ABVertex *curr = vertexHead;
150 while( curr != nullptr ) // Loop checking << 172 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
151 { 173 {
152 *anext++ = curr->a; 174 *anext++ = curr->a;
153 *bnext++ = curr->b; 175 *bnext++ = curr->b;
154 curr = curr->next; 176 curr = curr->next;
155 } 177 }
156 } 178 }
157 179
>> 180
>> 181 //
158 // ScaleA 182 // ScaleA
159 // 183 //
160 // Multiply all a values by a common scale 184 // Multiply all a values by a common scale
161 // 185 //
162 void G4ReduciblePolygon::ScaleA( G4double scal 186 void G4ReduciblePolygon::ScaleA( G4double scale )
163 { 187 {
164 ABVertex* curr = vertexHead; << 188 ABVertex *curr = vertexHead;
165 while( curr != nullptr ) // Loop checking << 189 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
166 { 190 {
167 curr->a *= scale; 191 curr->a *= scale;
168 curr = curr->next; 192 curr = curr->next;
169 } 193 }
170 } 194 }
171 195
>> 196
>> 197 //
172 // ScaleB 198 // ScaleB
173 // 199 //
174 // Multiply all b values by a common scale 200 // Multiply all b values by a common scale
175 // 201 //
176 void G4ReduciblePolygon::ScaleB( G4double scal 202 void G4ReduciblePolygon::ScaleB( G4double scale )
177 { 203 {
178 ABVertex* curr = vertexHead; << 204 ABVertex *curr = vertexHead;
179 while( curr != nullptr ) // Loop checking << 205 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
180 { 206 {
181 curr->b *= scale; 207 curr->b *= scale;
182 curr = curr->next; 208 curr = curr->next;
183 } 209 }
184 } 210 }
185 211
>> 212
>> 213 //
186 // RemoveDuplicateVertices 214 // RemoveDuplicateVertices
187 // 215 //
188 // Remove adjacent vertices that are equal. Re 216 // Remove adjacent vertices that are equal. Returns "false" if there
189 // is a problem (too few vertices remaining). 217 // is a problem (too few vertices remaining).
190 // 218 //
191 G4bool G4ReduciblePolygon::RemoveDuplicateVert 219 G4bool G4ReduciblePolygon::RemoveDuplicateVertices( G4double tolerance )
192 { 220 {
193 ABVertex *curr = vertexHead, 221 ABVertex *curr = vertexHead,
194 *prev = nullptr, *next = nullptr; << 222 *prev = 0, *next = 0;
195 while( curr != nullptr ) // Loop checking << 223 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
196 { 224 {
197 next = curr->next; 225 next = curr->next;
198 if (next == nullptr) next = vertexHead; << 226 if (next == 0) next = vertexHead;
199 227
200 if (std::fabs(curr->a-next->a) < tolerance 228 if (std::fabs(curr->a-next->a) < tolerance &&
201 std::fabs(curr->b-next->b) < tolerance 229 std::fabs(curr->b-next->b) < tolerance )
202 { 230 {
203 // 231 //
204 // Duplicate found: do we have > 3 verti 232 // Duplicate found: do we have > 3 vertices?
205 // 233 //
206 if (numVertices <= 3) 234 if (numVertices <= 3)
207 { 235 {
208 CalculateMaxMin(); 236 CalculateMaxMin();
209 return false; 237 return false;
210 } 238 }
211 239
212 // 240 //
213 // Delete 241 // Delete
214 // 242 //
215 ABVertex* toDelete = curr; << 243 ABVertex *toDelete = curr;
216 curr = curr->next; 244 curr = curr->next;
217 delete toDelete; 245 delete toDelete;
218 246
219 numVertices--; 247 numVertices--;
220 248
221 if (prev != nullptr) << 249 if (prev) prev->next = curr; else vertexHead = curr;
222 prev->next = curr; <<
223 else <<
224 vertexHead = curr; <<
225 } 250 }
226 else 251 else
227 { 252 {
228 prev = curr; 253 prev = curr;
229 curr = curr->next; 254 curr = curr->next;
230 } 255 }
231 } 256 }
232 257
233 // 258 //
234 // In principle, this is not needed, but why 259 // In principle, this is not needed, but why not just play it safe?
235 // 260 //
236 CalculateMaxMin(); 261 CalculateMaxMin();
237 262
238 return true; 263 return true;
239 } 264 }
240 265
>> 266
>> 267 //
241 // RemoveRedundantVertices 268 // RemoveRedundantVertices
242 // 269 //
243 // Remove any unneeded vertices, i.e. those ve 270 // Remove any unneeded vertices, i.e. those vertices which
244 // are on the line connecting the previous and 271 // are on the line connecting the previous and next vertices.
245 // 272 //
246 G4bool G4ReduciblePolygon::RemoveRedundantVert 273 G4bool G4ReduciblePolygon::RemoveRedundantVertices( G4double tolerance )
247 { 274 {
248 // 275 //
249 // Under these circumstances, we can quit no 276 // Under these circumstances, we can quit now!
250 // 277 //
251 if (numVertices <= 2) return false; 278 if (numVertices <= 2) return false;
252 279
253 G4double tolerance2 = tolerance*tolerance; 280 G4double tolerance2 = tolerance*tolerance;
254 281
255 // 282 //
256 // Loop over all vertices 283 // Loop over all vertices
257 // 284 //
258 ABVertex *curr = vertexHead, *next = nullptr << 285 ABVertex *curr = vertexHead, *next = 0;
259 while( curr != nullptr ) // Loop checking << 286 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
260 { 287 {
261 next = curr->next; 288 next = curr->next;
262 if (next == nullptr) next = vertexHead; << 289 if (next == 0) next = vertexHead;
263 290
264 G4double da = next->a - curr->a, 291 G4double da = next->a - curr->a,
265 db = next->b - curr->b; 292 db = next->b - curr->b;
266 293
267 // 294 //
268 // Loop over all subsequent vertices, up t 295 // Loop over all subsequent vertices, up to curr
269 // 296 //
270 for(;;) 297 for(;;)
271 { 298 {
272 // 299 //
273 // Get vertex after next 300 // Get vertex after next
274 // 301 //
275 ABVertex* test = next->next; << 302 ABVertex *test = next->next;
276 if (test == nullptr) test = vertexHead; << 303 if (test == 0) test = vertexHead;
277 304
278 // 305 //
279 // If we are back to the original vertex 306 // If we are back to the original vertex, stop
280 // 307 //
281 if (test==curr) break; 308 if (test==curr) break;
282 309
283 // 310 //
284 // Test for parallel line segments 311 // Test for parallel line segments
285 // 312 //
286 G4double dat = test->a - curr->a, 313 G4double dat = test->a - curr->a,
287 dbt = test->b - curr->b; 314 dbt = test->b - curr->b;
288 315
289 if (std::fabs(dat*db-dbt*da)>tolerance2) 316 if (std::fabs(dat*db-dbt*da)>tolerance2) break;
290 317
291 // 318 //
292 // Redundant vertex found: do we have > 319 // Redundant vertex found: do we have > 3 vertices?
293 // 320 //
294 if (numVertices <= 3) 321 if (numVertices <= 3)
295 { 322 {
296 CalculateMaxMin(); 323 CalculateMaxMin();
297 return false; 324 return false;
298 } 325 }
299 326
300 // 327 //
301 // Delete vertex pointed to by next. Car 328 // Delete vertex pointed to by next. Carefully!
302 // 329 //
303 if (curr->next != nullptr) << 330 if (curr->next)
304 { // next is not head 331 { // next is not head
305 if (next->next != nullptr) << 332 if (next->next)
306 curr->next = test; // next is not t 333 curr->next = test; // next is not tail
307 else 334 else
308 curr->next = nullptr; // New tail << 335 curr->next = 0; // New tail
309 } 336 }
310 else 337 else
311 vertexHead = test; // New head 338 vertexHead = test; // New head
312 339
313 if ((curr != next) && (next != test)) de 340 if ((curr != next) && (next != test)) delete next;
314 341
315 --numVertices; << 342 numVertices--;
316 343
317 // 344 //
318 // Replace next by the vertex we just te 345 // Replace next by the vertex we just tested,
319 // and keep on going... 346 // and keep on going...
320 // 347 //
321 next = test; 348 next = test;
322 da = dat; db = dbt; 349 da = dat; db = dbt;
323 } 350 }
324 curr = curr->next; 351 curr = curr->next;
325 } 352 }
326 353
327 // 354 //
328 // In principle, this is not needed, but why 355 // In principle, this is not needed, but why not just play it safe?
329 // 356 //
330 CalculateMaxMin(); 357 CalculateMaxMin();
331 358
332 return true; 359 return true;
333 } 360 }
334 361
>> 362
>> 363 //
335 // ReverseOrder 364 // ReverseOrder
336 // 365 //
337 // Reverse the order of the vertices 366 // Reverse the order of the vertices
338 // 367 //
339 void G4ReduciblePolygon::ReverseOrder() 368 void G4ReduciblePolygon::ReverseOrder()
340 { 369 {
341 // 370 //
342 // Loop over all vertices 371 // Loop over all vertices
343 // 372 //
344 ABVertex* prev = vertexHead; << 373 ABVertex *prev = vertexHead;
345 if (prev==nullptr) return; // No vertices << 374 if (prev==0) return; // No vertices
346 375
347 ABVertex* curr = prev->next; << 376 ABVertex *curr = prev->next;
348 if (curr==nullptr) return; // Just one ve << 377 if (curr==0) return; // Just one vertex
349 378
350 // 379 //
351 // Our new tail 380 // Our new tail
352 // 381 //
353 vertexHead->next = nullptr; << 382 vertexHead->next = 0;
354 383
355 for(;;) 384 for(;;)
356 { 385 {
357 // 386 //
358 // Save pointer to next vertex (in origina 387 // Save pointer to next vertex (in original order)
359 // 388 //
360 ABVertex *save = curr->next; 389 ABVertex *save = curr->next;
361 390
362 // 391 //
363 // Replace it with a pointer to the previo 392 // Replace it with a pointer to the previous one
364 // (in original order) 393 // (in original order)
365 // 394 //
366 curr->next = prev; 395 curr->next = prev;
367 396
368 // 397 //
369 // Last vertex? 398 // Last vertex?
370 // 399 //
371 if (save == nullptr) break; << 400 if (save == 0) break;
372 401
373 // 402 //
374 // Next vertex 403 // Next vertex
375 // 404 //
376 prev = curr; 405 prev = curr;
377 curr = save; 406 curr = save;
378 } 407 }
379 408
380 // 409 //
381 // Our new head 410 // Our new head
382 // 411 //
383 vertexHead = curr; 412 vertexHead = curr;
384 } 413 }
385 414
386 415
387 // StartWithZMin 416 // StartWithZMin
388 // 417 //
389 // Starting alway with Zmin=bMin 418 // Starting alway with Zmin=bMin
390 // This method is used for GenericPolycone 419 // This method is used for GenericPolycone
391 // 420 //
392 void G4ReduciblePolygon::StartWithZMin() 421 void G4ReduciblePolygon::StartWithZMin()
393 { 422 {
394 ABVertex* curr = vertexHead; << 423 ABVertex *curr = vertexHead;
395 G4double bcurr = curr->b; 424 G4double bcurr = curr->b;
396 ABVertex* prev = curr; << 425 ABVertex *prev = curr;
397 while( curr != nullptr) // Loop checking, << 426 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
398 { 427 {
399 if(curr->b < bcurr) 428 if(curr->b < bcurr)
400 { 429 {
401 bcurr = curr->b; 430 bcurr = curr->b;
402 ABVertex* curr1 = curr; << 431 ABVertex *curr1 = curr;
403 while( curr1 != nullptr ) // Loop che << 432 while( curr1 ) // Loop checking, 13.08.2015, G.Cosmo
404 { 433 {
405 if(curr1->next == nullptr) { curr1->ne << 434 if(curr1->next == 0) { curr1->next = vertexHead; break; }
406 curr1 = curr1->next; 435 curr1 = curr1->next;
407 } 436 }
408 vertexHead = curr; 437 vertexHead = curr;
409 prev->next = nullptr; << 438 prev->next = 0;
410 } 439 }
411 prev = curr; 440 prev = curr;
412 curr = curr->next; 441 curr = curr->next;
413 } 442 }
414 } 443 }
415 444
>> 445
>> 446 //
416 // CrossesItself 447 // CrossesItself
417 // 448 //
418 // Return "true" if the polygon crosses itself 449 // Return "true" if the polygon crosses itself
419 // 450 //
420 // Warning: this routine is not very fast (run 451 // Warning: this routine is not very fast (runs as N**2)
421 // 452 //
422 G4bool G4ReduciblePolygon::CrossesItself( G4do 453 G4bool G4ReduciblePolygon::CrossesItself( G4double tolerance )
423 { 454 {
424 G4double tolerance2 = tolerance*tolerance; 455 G4double tolerance2 = tolerance*tolerance;
425 G4double one = 1.0-tolerance, 456 G4double one = 1.0-tolerance,
426 zero = tolerance; 457 zero = tolerance;
427 // 458 //
428 // Top loop over line segments. By the time 459 // Top loop over line segments. By the time we finish
429 // with the second to last segment, we're do 460 // with the second to last segment, we're done.
430 // 461 //
431 ABVertex *curr1 = vertexHead, *next1 = nullp << 462 ABVertex *curr1 = vertexHead, *next1=0;
432 while (curr1->next != nullptr) // Loop ch << 463 while (curr1->next) // Loop checking, 13.08.2015, G.Cosmo
433 { 464 {
434 next1 = curr1->next; 465 next1 = curr1->next;
435 G4double da1 = next1->a-curr1->a, 466 G4double da1 = next1->a-curr1->a,
436 db1 = next1->b-curr1->b; 467 db1 = next1->b-curr1->b;
437 468
438 // 469 //
439 // Inner loop over subsequent line segment 470 // Inner loop over subsequent line segments
440 // 471 //
441 ABVertex* curr2 = next1->next; << 472 ABVertex *curr2 = next1->next;
442 while( curr2 != nullptr ) // Loop check << 473 while( curr2 ) // Loop checking, 13.08.2015, G.Cosmo
443 { 474 {
444 ABVertex* next2 = curr2->next; << 475 ABVertex *next2 = curr2->next;
445 if (next2==nullptr) next2 = vertexHead; << 476 if (next2==0) next2 = vertexHead;
446 G4double da2 = next2->a-curr2->a, 477 G4double da2 = next2->a-curr2->a,
447 db2 = next2->b-curr2->b; 478 db2 = next2->b-curr2->b;
448 G4double a12 = curr2->a-curr1->a, 479 G4double a12 = curr2->a-curr1->a,
449 b12 = curr2->b-curr1->b; 480 b12 = curr2->b-curr1->b;
450 481
451 // 482 //
452 // Calculate intersection of the two lin 483 // Calculate intersection of the two lines
453 // 484 //
454 G4double deter = da1*db2 - db1*da2; 485 G4double deter = da1*db2 - db1*da2;
455 if (std::fabs(deter) > tolerance2) 486 if (std::fabs(deter) > tolerance2)
456 { 487 {
457 G4double s1, s2; 488 G4double s1, s2;
458 s1 = (a12*db2-b12*da2)/deter; 489 s1 = (a12*db2-b12*da2)/deter;
459 490
460 if (s1 >= zero && s1 < one) 491 if (s1 >= zero && s1 < one)
461 { 492 {
462 s2 = -(da1*b12-db1*a12)/deter; 493 s2 = -(da1*b12-db1*a12)/deter;
463 if (s2 >= zero && s2 < one) return t 494 if (s2 >= zero && s2 < one) return true;
464 } 495 }
465 } 496 }
466 curr2 = curr2->next; 497 curr2 = curr2->next;
467 } 498 }
468 curr1 = next1; 499 curr1 = next1;
469 } 500 }
470 return false; 501 return false;
471 } 502 }
472 503
>> 504
>> 505
>> 506 //
473 // BisectedBy 507 // BisectedBy
474 // 508 //
475 // Decide if a line through two points crosses 509 // Decide if a line through two points crosses the polygon, within tolerance
476 // 510 //
477 G4bool G4ReduciblePolygon::BisectedBy( G4doubl 511 G4bool G4ReduciblePolygon::BisectedBy( G4double a1, G4double b1,
478 G4doubl 512 G4double a2, G4double b2,
479 G4doubl 513 G4double tolerance )
480 { 514 {
481 G4int nNeg = 0, nPos = 0; 515 G4int nNeg = 0, nPos = 0;
482 516
483 G4double a12 = a2-a1, b12 = b2-b1; 517 G4double a12 = a2-a1, b12 = b2-b1;
484 G4double len12 = std::sqrt( a12*a12 + b12*b1 518 G4double len12 = std::sqrt( a12*a12 + b12*b12 );
485 a12 /= len12; b12 /= len12; 519 a12 /= len12; b12 /= len12;
486 520
487 ABVertex* curr = vertexHead; << 521 ABVertex *curr = vertexHead;
488 do // Loop checking, 13.08.2015, G.Cosmo 522 do // Loop checking, 13.08.2015, G.Cosmo
489 { 523 {
490 G4double av = curr->a - a1, 524 G4double av = curr->a - a1,
491 bv = curr->b - b1; << 525 bv = curr->b - b1;
492 526
493 G4double cross = av*b12 - bv*a12; 527 G4double cross = av*b12 - bv*a12;
494 528
495 if (cross < -tolerance) 529 if (cross < -tolerance)
496 { 530 {
497 if (nPos != 0) return true; << 531 if (nPos) return true;
498 ++nNeg; << 532 nNeg++;
499 } 533 }
500 else if (cross > tolerance) 534 else if (cross > tolerance)
501 { 535 {
502 if (nNeg != 0) return true; << 536 if (nNeg) return true;
503 ++nPos; << 537 nPos++;
504 } 538 }
505 curr = curr->next; 539 curr = curr->next;
506 } while( curr != nullptr ); << 540 } while( curr );
507 541
508 return false; 542 return false;
509 } 543 }
510 544
>> 545
>> 546
>> 547 //
511 // Area 548 // Area
512 // 549 //
513 // Calculated signed polygon area, where polyg 550 // Calculated signed polygon area, where polygons specified in a
514 // clockwise manner (where x==a, y==b) have ne 551 // clockwise manner (where x==a, y==b) have negative area
515 // 552 //
516 // References: [O' Rourke (C)] pp. 18-27; [ 553 // References: [O' Rourke (C)] pp. 18-27; [Gems II] pp. 5-6:
517 // "The Area of a Simple Polygon", Jon Rokn 554 // "The Area of a Simple Polygon", Jon Rokne.
518 // 555 //
519 G4double G4ReduciblePolygon::Area() 556 G4double G4ReduciblePolygon::Area()
520 { 557 {
521 G4double answer = 0; 558 G4double answer = 0;
522 559
523 ABVertex *curr = vertexHead, *next = nullptr << 560 ABVertex *curr = vertexHead, *next;
524 do // Loop checking, 13.08.2015, G.Cosmo 561 do // Loop checking, 13.08.2015, G.Cosmo
525 { 562 {
526 next = curr->next; 563 next = curr->next;
527 if (next==nullptr) next = vertexHead; << 564 if (next==0) next = vertexHead;
528 565
529 answer += curr->a*next->b - curr->b*next-> 566 answer += curr->a*next->b - curr->b*next->a;
530 curr = curr->next; 567 curr = curr->next;
531 } while( curr != nullptr ); << 568 } while( curr );
532 569
533 return 0.5*answer; 570 return 0.5*answer;
534 } 571 }
535 572
>> 573
>> 574 //
536 // Print 575 // Print
537 // 576 //
538 void G4ReduciblePolygon::Print() 577 void G4ReduciblePolygon::Print()
539 { 578 {
540 ABVertex* curr = vertexHead; << 579 ABVertex *curr = vertexHead;
541 do // Loop checking, 13.08.2015, G.Cosmo 580 do // Loop checking, 13.08.2015, G.Cosmo
542 { 581 {
543 G4cerr << curr->a << " " << curr->b << G4e 582 G4cerr << curr->a << " " << curr->b << G4endl;
544 curr = curr->next; 583 curr = curr->next;
545 } while( curr != nullptr ); << 584 } while( curr );
546 } 585 }
547 586
>> 587
>> 588 //
548 // CalculateMaxMin 589 // CalculateMaxMin
549 // 590 //
550 // To be called when the vertices are changed, 591 // To be called when the vertices are changed, this
551 // routine re-calculates global values 592 // routine re-calculates global values
552 // 593 //
553 void G4ReduciblePolygon::CalculateMaxMin() 594 void G4ReduciblePolygon::CalculateMaxMin()
554 { 595 {
555 ABVertex* curr = vertexHead; << 596 ABVertex *curr = vertexHead;
556 aMin = aMax = curr->a; 597 aMin = aMax = curr->a;
557 bMin = bMax = curr->b; 598 bMin = bMax = curr->b;
558 curr = curr->next; 599 curr = curr->next;
559 while( curr != nullptr ) // Loop checking << 600 while( curr ) // Loop checking, 13.08.2015, G.Cosmo
560 { 601 {
561 if (curr->a < aMin) 602 if (curr->a < aMin)
562 aMin = curr->a; 603 aMin = curr->a;
563 else if (curr->a > aMax) 604 else if (curr->a > aMax)
564 aMax = curr->a; 605 aMax = curr->a;
565 606
566 if (curr->b < bMin) 607 if (curr->b < bMin)
567 bMin = curr->b; 608 bMin = curr->b;
568 else if (curr->b > bMax) 609 else if (curr->b > bMax)
569 bMax = curr->b; 610 bMax = curr->b;
570 611
571 curr = curr->next; 612 curr = curr->next;
572 } 613 }
573 } 614 }
574 615