00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037 #include "G4ReflectedSolid.hh"
00038
00039 #include <sstream>
00040
00041 #include "G4Point3D.hh"
00042 #include "G4Normal3D.hh"
00043
00044 #include "G4VoxelLimits.hh"
00045
00046 #include "G4VPVParameterisation.hh"
00047
00048 #include "G4VGraphicsScene.hh"
00049 #include "G4Polyhedron.hh"
00050 #include "G4NURBS.hh"
00051
00052
00053
00055
00056
00057
00058 G4ReflectedSolid::G4ReflectedSolid( const G4String& pName,
00059 G4VSolid* pSolid ,
00060 const G4Transform3D& transform )
00061 : G4VSolid(pName), fpPolyhedron(0)
00062 {
00063 fPtrSolid = pSolid ;
00064 G4RotationMatrix rotMatrix ;
00065
00066 fDirectTransform =
00067 new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
00068 fPtrTransform =
00069 new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
00070 fPtrTransform->Invert() ;
00071
00072 fDirectTransform3D = new G4Transform3D(transform) ;
00073 fPtrTransform3D = new G4Transform3D(transform.inverse()) ;
00074 }
00075
00077
00078
00079 G4ReflectedSolid::~G4ReflectedSolid()
00080 {
00081 if(fPtrTransform)
00082 {
00083 delete fPtrTransform; fPtrTransform=0;
00084 delete fDirectTransform; fDirectTransform=0;
00085 }
00086 if(fPtrTransform3D)
00087 {
00088 delete fPtrTransform3D; fPtrTransform3D=0;
00089 delete fDirectTransform3D; fDirectTransform3D=0;
00090 }
00091 delete fpPolyhedron;
00092 }
00093
00095
00096
00097 G4ReflectedSolid::G4ReflectedSolid(const G4ReflectedSolid& rhs)
00098 : G4VSolid(rhs), fPtrSolid(rhs.fPtrSolid), fpPolyhedron(0)
00099 {
00100 fPtrTransform = new G4AffineTransform(*rhs.fPtrTransform);
00101 fDirectTransform = new G4AffineTransform(*rhs.fDirectTransform);
00102 fPtrTransform3D = new G4Transform3D(*rhs.fPtrTransform3D);
00103 fDirectTransform3D = new G4Transform3D(*rhs.fDirectTransform3D);
00104 }
00105
00107
00108
00109 G4ReflectedSolid& G4ReflectedSolid::operator=(const G4ReflectedSolid& rhs)
00110 {
00111
00112
00113 if (this == &rhs) { return *this; }
00114
00115
00116
00117 G4VSolid::operator=(rhs);
00118
00119
00120
00121 fPtrSolid= rhs.fPtrSolid; fpPolyhedron= 0;
00122 delete fPtrTransform;
00123 fPtrTransform= new G4AffineTransform(*rhs.fPtrTransform);
00124 delete fDirectTransform;
00125 fDirectTransform= new G4AffineTransform(*rhs.fDirectTransform);
00126 delete fPtrTransform3D;
00127 fPtrTransform3D= new G4Transform3D(*rhs.fPtrTransform3D);
00128 delete fDirectTransform3D;
00129 fDirectTransform3D= new G4Transform3D(*rhs.fDirectTransform3D);
00130
00131 return *this;
00132 }
00133
00135
00136
00137 G4GeometryType G4ReflectedSolid::GetEntityType() const
00138 {
00139 return G4String("G4ReflectedSolid");
00140 }
00141
00142 const G4ReflectedSolid* G4ReflectedSolid::GetReflectedSolidPtr() const
00143 {
00144 return this;
00145 }
00146
00147 G4ReflectedSolid* G4ReflectedSolid::GetReflectedSolidPtr()
00148 {
00149 return this;
00150 }
00151
00152 G4VSolid* G4ReflectedSolid::GetConstituentMovedSolid() const
00153 {
00154 return fPtrSolid;
00155 }
00156
00158
00159 G4AffineTransform G4ReflectedSolid::GetTransform() const
00160 {
00161 G4AffineTransform aTransform = *fPtrTransform;
00162 return aTransform;
00163 }
00164
00165 void G4ReflectedSolid::SetTransform(G4AffineTransform& transform)
00166 {
00167 fPtrTransform = &transform ;
00168 fpPolyhedron = 0;
00169 }
00170
00172
00173 G4AffineTransform G4ReflectedSolid::GetDirectTransform() const
00174 {
00175 G4AffineTransform aTransform= *fDirectTransform;
00176 return aTransform;
00177 }
00178
00179 void G4ReflectedSolid::SetDirectTransform(G4AffineTransform& transform)
00180 {
00181 fDirectTransform = &transform ;
00182 fpPolyhedron = 0;
00183 }
00184
00186
00187 G4Transform3D G4ReflectedSolid::GetTransform3D() const
00188 {
00189 G4Transform3D aTransform = *fPtrTransform3D;
00190 return aTransform;
00191 }
00192
00193 void G4ReflectedSolid::SetTransform3D(G4Transform3D& transform)
00194 {
00195 fPtrTransform3D = &transform ;
00196 fpPolyhedron = 0;
00197 }
00198
00200
00201 G4Transform3D G4ReflectedSolid::GetDirectTransform3D() const
00202 {
00203 G4Transform3D aTransform= *fDirectTransform3D;
00204 return aTransform;
00205 }
00206
00207 void G4ReflectedSolid::SetDirectTransform3D(G4Transform3D& transform)
00208 {
00209 fDirectTransform3D = &transform ;
00210 fpPolyhedron = 0;
00211 }
00212
00214
00215 G4RotationMatrix G4ReflectedSolid::GetFrameRotation() const
00216 {
00217 G4RotationMatrix InvRotation= fDirectTransform->NetRotation();
00218 return InvRotation;
00219 }
00220
00221 void G4ReflectedSolid::SetFrameRotation(const G4RotationMatrix& matrix)
00222 {
00223 fDirectTransform->SetNetRotation(matrix);
00224 }
00225
00227
00228 G4ThreeVector G4ReflectedSolid::GetFrameTranslation() const
00229 {
00230 return fPtrTransform->NetTranslation();
00231 }
00232
00233 void G4ReflectedSolid::SetFrameTranslation(const G4ThreeVector& vector)
00234 {
00235 fPtrTransform->SetNetTranslation(vector);
00236 }
00237
00239
00240 G4RotationMatrix G4ReflectedSolid::GetObjectRotation() const
00241 {
00242 G4RotationMatrix Rotation= fPtrTransform->NetRotation();
00243 return Rotation;
00244 }
00245
00246 void G4ReflectedSolid::SetObjectRotation(const G4RotationMatrix& matrix)
00247 {
00248 fPtrTransform->SetNetRotation(matrix);
00249 }
00250
00252
00253 G4ThreeVector G4ReflectedSolid::GetObjectTranslation() const
00254 {
00255 return fDirectTransform->NetTranslation();
00256 }
00257
00258 void G4ReflectedSolid::SetObjectTranslation(const G4ThreeVector& vector)
00259 {
00260 fDirectTransform->SetNetTranslation(vector);
00261 }
00262
00264
00265
00266
00267 G4bool
00268 G4ReflectedSolid::CalculateExtent( const EAxis pAxis,
00269 const G4VoxelLimits& pVoxelLimit,
00270 const G4AffineTransform& pTransform,
00271 G4double& pMin,
00272 G4double& pMax ) const
00273 {
00274
00275 G4VoxelLimits unLimit;
00276 G4AffineTransform unTransform;
00277
00278 G4double x1 = -kInfinity, x2 = kInfinity,
00279 y1 = -kInfinity, y2 = kInfinity,
00280 z1 = -kInfinity, z2 = kInfinity;
00281
00282 G4bool existsAfterClip = false ;
00283 existsAfterClip =
00284 fPtrSolid->CalculateExtent(kXAxis,unLimit,unTransform,x1,x2);
00285 existsAfterClip =
00286 fPtrSolid->CalculateExtent(kYAxis,unLimit,unTransform,y1,y2);
00287 existsAfterClip =
00288 fPtrSolid->CalculateExtent(kZAxis,unLimit,unTransform,z1,z2);
00289
00290 existsAfterClip = false;
00291 pMin = +kInfinity ;
00292 pMax = -kInfinity ;
00293
00294 G4Transform3D pTransform3D = G4Transform3D(pTransform.NetRotation().inverse(),
00295 pTransform.NetTranslation());
00296
00297 G4Transform3D transform3D = pTransform3D*(*fDirectTransform3D);
00298
00299 G4Point3D tmpPoint;
00300
00301
00302
00303 G4ThreeVectorList* vertices = new G4ThreeVectorList();
00304
00305 if (vertices)
00306 {
00307 vertices->reserve(8);
00308
00309 G4ThreeVector vertex0(x1,y1,z1) ;
00310 tmpPoint = transform3D*G4Point3D(vertex0);
00311 vertex0 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00312 vertices->push_back(vertex0);
00313
00314 G4ThreeVector vertex1(x2,y1,z1) ;
00315 tmpPoint = transform3D*G4Point3D(vertex1);
00316 vertex1 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00317 vertices->push_back(vertex1);
00318
00319 G4ThreeVector vertex2(x2,y2,z1) ;
00320 tmpPoint = transform3D*G4Point3D(vertex2);
00321 vertex2 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00322 vertices->push_back(vertex2);
00323
00324 G4ThreeVector vertex3(x1,y2,z1) ;
00325 tmpPoint = transform3D*G4Point3D(vertex3);
00326 vertex3 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00327 vertices->push_back(vertex3);
00328
00329 G4ThreeVector vertex4(x1,y1,z2) ;
00330 tmpPoint = transform3D*G4Point3D(vertex4);
00331 vertex4 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00332 vertices->push_back(vertex4);
00333
00334 G4ThreeVector vertex5(x2,y1,z2) ;
00335 tmpPoint = transform3D*G4Point3D(vertex5);
00336 vertex5 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00337 vertices->push_back(vertex5);
00338
00339 G4ThreeVector vertex6(x2,y2,z2) ;
00340 tmpPoint = transform3D*G4Point3D(vertex6);
00341 vertex6 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00342 vertices->push_back(vertex6);
00343
00344 G4ThreeVector vertex7(x1,y2,z2) ;
00345 tmpPoint = transform3D*G4Point3D(vertex7);
00346 vertex7 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
00347 vertices->push_back(vertex7);
00348 }
00349 else
00350 {
00351 DumpInfo();
00352 G4Exception("G4ReflectedSolid::CalculateExtent()",
00353 "GeomMgt0003", FatalException,
00354 "Error in allocation of vertices. Out of memory !");
00355 }
00356
00357 ClipCrossSection(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
00358 ClipCrossSection(vertices,4,pVoxelLimit,pAxis,pMin,pMax) ;
00359 ClipBetweenSections(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
00360
00361 if (pVoxelLimit.IsLimited(pAxis) == false)
00362 {
00363 if ( pMin != kInfinity || pMax != -kInfinity )
00364 {
00365 existsAfterClip = true ;
00366
00367
00368
00369 pMin -= kCarTolerance;
00370 pMax += kCarTolerance;
00371 }
00372 }
00373 else
00374 {
00375 G4ThreeVector clipCentre(
00376 ( pVoxelLimit.GetMinXExtent()+pVoxelLimit.GetMaxXExtent())*0.5,
00377 ( pVoxelLimit.GetMinYExtent()+pVoxelLimit.GetMaxYExtent())*0.5,
00378 ( pVoxelLimit.GetMinZExtent()+pVoxelLimit.GetMaxZExtent())*0.5);
00379
00380 if ( pMin != kInfinity || pMax != -kInfinity )
00381 {
00382 existsAfterClip = true ;
00383
00384
00385
00386
00387 clipCentre(pAxis) = pVoxelLimit.GetMinExtent(pAxis);
00388
00389 if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
00390 {
00391 pMin = pVoxelLimit.GetMinExtent(pAxis);
00392 }
00393 else
00394 {
00395 pMin -= kCarTolerance;
00396 }
00397 clipCentre(pAxis) = pVoxelLimit.GetMaxExtent(pAxis);
00398
00399 if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
00400 {
00401 pMax = pVoxelLimit.GetMaxExtent(pAxis);
00402 }
00403 else
00404 {
00405 pMax += kCarTolerance;
00406 }
00407 }
00408
00409
00410
00411
00412
00413 else if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
00414 {
00415 existsAfterClip = true ;
00416 pMin = pVoxelLimit.GetMinExtent(pAxis) ;
00417 pMax = pVoxelLimit.GetMaxExtent(pAxis) ;
00418 }
00419 }
00420 delete vertices;
00421 return existsAfterClip;
00422 }
00423
00425
00426
00427
00428 EInside G4ReflectedSolid::Inside(const G4ThreeVector& p) const
00429 {
00430
00431 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
00432
00433
00434 return fPtrSolid->Inside(G4ThreeVector(newPoint.x(),
00435 newPoint.y(),
00436 newPoint.z())) ;
00437 }
00438
00440
00441
00442
00443 G4ThreeVector
00444 G4ReflectedSolid::SurfaceNormal( const G4ThreeVector& p ) const
00445 {
00446 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
00447 G4ThreeVector normal =
00448 fPtrSolid->SurfaceNormal(G4ThreeVector(newPoint.x(),
00449 newPoint.y(),
00450 newPoint.z() ) ) ;
00451 G4Point3D newN = (*fDirectTransform3D)*G4Point3D(normal) ;
00452 newN.unit() ;
00453
00454 return G4ThreeVector(newN.x(),newN.y(),newN.z()) ;
00455 }
00456
00458
00459
00460
00461 G4double
00462 G4ReflectedSolid::DistanceToIn( const G4ThreeVector& p,
00463 const G4ThreeVector& v ) const
00464 {
00465 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
00466 G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v) ;
00467 newDirection.unit() ;
00468 return fPtrSolid->DistanceToIn(
00469 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
00470 G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z())) ;
00471 }
00472
00474
00475
00476
00477
00478 G4double
00479 G4ReflectedSolid::DistanceToIn( const G4ThreeVector& p) const
00480 {
00481 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
00482 return fPtrSolid->DistanceToIn(
00483 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z())) ;
00484 }
00485
00487
00488
00489
00490 G4double
00491 G4ReflectedSolid::DistanceToOut( const G4ThreeVector& p,
00492 const G4ThreeVector& v,
00493 const G4bool calcNorm,
00494 G4bool *validNorm,
00495 G4ThreeVector *n ) const
00496 {
00497 G4ThreeVector solNorm ;
00498
00499 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
00500 G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v);
00501 newDirection.unit() ;
00502
00503 G4double dist =
00504 fPtrSolid->DistanceToOut(
00505 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
00506 G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z()),
00507 calcNorm, validNorm, &solNorm) ;
00508 if(calcNorm)
00509 {
00510 G4Point3D newN = (*fDirectTransform3D)*G4Point3D(solNorm);
00511 newN.unit() ;
00512 *n = G4ThreeVector(newN.x(),newN.y(),newN.z());
00513 }
00514 return dist ;
00515 }
00516
00518
00519
00520
00521 G4double
00522 G4ReflectedSolid::DistanceToOut( const G4ThreeVector& p ) const
00523 {
00524 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
00525 return fPtrSolid->DistanceToOut(
00526 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()));
00527 }
00528
00530
00531
00532
00533 void
00534 G4ReflectedSolid::ComputeDimensions( G4VPVParameterisation*,
00535 const G4int,
00536 const G4VPhysicalVolume* )
00537 {
00538 DumpInfo();
00539 G4Exception("G4ReflectedSolid::ComputeDimensions()",
00540 "GeomMgt0001", FatalException,
00541 "Method not applicable in this context!");
00542 }
00543
00545
00546
00547
00548
00549 G4ThreeVector G4ReflectedSolid::GetPointOnSurface() const
00550 {
00551 G4ThreeVector p = fPtrSolid->GetPointOnSurface();
00552 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
00553
00554 return G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z());
00555 }
00556
00558
00559
00560
00561 G4VSolid* G4ReflectedSolid::Clone() const
00562 {
00563 return new G4ReflectedSolid(*this);
00564 }
00565
00566
00568
00569
00570
00571 std::ostream& G4ReflectedSolid::StreamInfo(std::ostream& os) const
00572 {
00573 os << "-----------------------------------------------------------\n"
00574 << " *** Dump for Reflected solid - " << GetName() << " ***\n"
00575 << " ===================================================\n"
00576 << " Solid type: " << GetEntityType() << "\n"
00577 << " Parameters of constituent solid: \n"
00578 << "===========================================================\n";
00579 fPtrSolid->StreamInfo(os);
00580 os << "===========================================================\n"
00581 << " Transformations: \n"
00582 << " Direct transformation - translation : \n"
00583 << " " << fDirectTransform->NetTranslation() << "\n"
00584 << " - rotation : \n"
00585 << " ";
00586 fDirectTransform->NetRotation().print(os);
00587 os << "\n"
00588 << "===========================================================\n";
00589
00590 return os;
00591 }
00592
00594
00595
00596
00597 void
00598 G4ReflectedSolid::DescribeYourselfTo ( G4VGraphicsScene& scene ) const
00599 {
00600 scene.AddSolid (*this);
00601 }
00602
00604
00605
00606
00607 G4Polyhedron*
00608 G4ReflectedSolid::CreatePolyhedron () const
00609 {
00610 G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
00611 if (polyhedron)
00612 {
00613 polyhedron->Transform(*fDirectTransform3D);
00614 return polyhedron;
00615 }
00616 else
00617 {
00618 std::ostringstream message;
00619 message << "Solid - " << GetName()
00620 << " - original solid has no" << G4endl
00621 << "corresponding polyhedron. Returning NULL!";
00622 G4Exception("G4ReflectedSolid::CreatePolyhedron()",
00623 "GeomMgt1001", JustWarning, message);
00624 return 0;
00625 }
00626 }
00627
00629
00630
00631
00632 G4NURBS*
00633 G4ReflectedSolid::CreateNURBS () const
00634 {
00635
00636
00637 return 0;
00638 }
00639
00641
00642
00643
00644 G4Polyhedron*
00645 G4ReflectedSolid::GetPolyhedron () const
00646 {
00647 if (!fpPolyhedron ||
00648 fpPolyhedron->GetNumberOfRotationStepsAtTimeOfCreation() !=
00649 fpPolyhedron->GetNumberOfRotationSteps())
00650 {
00651 delete fpPolyhedron;
00652 fpPolyhedron = CreatePolyhedron ();
00653 }
00654 return fpPolyhedron;
00655 }