G4GeomTestVolume Class Reference

#include <G4GeomTestVolume.hh>


Public Member Functions

 G4GeomTestVolume (const G4VPhysicalVolume *theTarget, G4GeomTestLogger *theLogger, G4double theTolerance=1E-4)
 ~G4GeomTestVolume ()
G4double GetTolerance () const
void SetTolerance (G4double tolerance)
void TestCartGridXYZ (G4int nx=100, G4int ny=100, G4int nz=100)
void TestCartGridX (G4int ny=100, G4int nz=100)
void TestCartGridY (G4int nz=100, G4int nx=100)
void TestCartGridZ (G4int nx=100, G4int ny=100)
void TestCartGrid (const G4ThreeVector &g1, const G4ThreeVector &g2, const G4ThreeVector &v, G4int n1, G4int n2)
void TestRecursiveCartGrid (G4int nx=100, G4int ny=100, G4int nz=100, G4int sLevel=0, G4int depth=-1)
void TestCylinder (G4int nPhi=90, G4int nZ=50, G4int nRho=50, G4double fracZ=0.8, G4double fracRho=0.8, G4bool usePhi=false)
void TestRecursiveCylinder (G4int nPhi=90, G4int nZ=50, G4int nRho=50, G4double fracZ=0.8, G4double fracRho=0.8, G4bool usePhi=false, G4int sLevel=0, G4int depth=-1)
void TestOneLine (const G4ThreeVector &p, const G4ThreeVector &v)
void TestRecursiveLine (const G4ThreeVector &p, const G4ThreeVector &v, G4int sLevel=0, G4int depth=-1)
void ReportErrors ()
void ClearErrors ()

Protected Attributes

const G4VPhysicalVolumetarget
G4GeomTestLoggerlogger
G4double tolerance
G4VisExtent extent
std::map< G4long, G4GeomTestOverlapListoverlaps
std::map< G4long, G4GeomTestOvershootListovershoots


Detailed Description

Definition at line 54 of file G4GeomTestVolume.hh.


Constructor & Destructor Documentation

G4GeomTestVolume::G4GeomTestVolume ( const G4VPhysicalVolume theTarget,
G4GeomTestLogger theLogger,
G4double  theTolerance = 1E-4 
)

Definition at line 56 of file G4GeomTestVolume.cc.

00059   : target(theTarget),
00060     logger(theLogger),
00061     tolerance(theTolerance),
00062     extent(theTarget->GetLogicalVolume()->GetSolid()->GetExtent())
00063 {;}

G4GeomTestVolume::~G4GeomTestVolume (  ) 

Definition at line 68 of file G4GeomTestVolume.cc.

00068 {;}


Member Function Documentation

void G4GeomTestVolume::ClearErrors (  ) 

Definition at line 687 of file G4GeomTestVolume.cc.

References overlaps, and overshoots.

00688 {
00689   overshoots.clear();
00690   overlaps.clear();
00691 }

G4double G4GeomTestVolume::GetTolerance (  )  const

Definition at line 73 of file G4GeomTestVolume.cc.

References tolerance.

00074 {
00075   return tolerance;
00076 }

void G4GeomTestVolume::ReportErrors (  ) 

Definition at line 649 of file G4GeomTestVolume.cc.

References logger, G4GeomTestLogger::NoProblem(), G4GeomTestLogger::OverlappingDaughters(), overlaps, G4GeomTestLogger::OvershootingDaughter(), and overshoots.

Referenced by TestRecursiveCartGrid(), TestRecursiveCylinder(), and TestRecursiveLine().

00650 {
00651   //
00652   // Report overshoots
00653   //
00654   if (overshoots.empty())
00655     logger->NoProblem("GeomTest: no daughter volume extending outside mother detected.");
00656   else
00657   {
00658     std::map<G4long,G4GeomTestOvershootList>::iterator overshoot =
00659       overshoots.begin();
00660     while( overshoot != overshoots.end() )
00661     {
00662       logger->OvershootingDaughter( &(*overshoot).second );
00663       ++overshoot;
00664     }
00665   }
00666 
00667   //
00668   // Report overlaps
00669   //
00670   if (overlaps.empty())
00671     logger->NoProblem("GeomTest: no overlapping daughters detected.");
00672   else
00673   {
00674     std::map<G4long,G4GeomTestOverlapList>::iterator overlap =
00675       overlaps.begin();
00676     while( overlap != overlaps.end() )
00677     {
00678       logger->OverlappingDaughters( &(*overlap).second );
00679       ++overlap;
00680     }
00681   }
00682 }

void G4GeomTestVolume::SetTolerance ( G4double  tolerance  ) 

Definition at line 81 of file G4GeomTestVolume.cc.

References tolerance.

00082 {
00083   tolerance = val;
00084 }

void G4GeomTestVolume::TestCartGrid ( const G4ThreeVector g1,
const G4ThreeVector g2,
const G4ThreeVector v,
G4int  n1,
G4int  n2 
)

Definition at line 333 of file G4GeomTestVolume.cc.

References extent, FatalException, G4Exception(), G4VisExtent::GetXmax(), G4VisExtent::GetXmin(), G4VisExtent::GetYmax(), G4VisExtent::GetYmin(), G4VisExtent::GetZmax(), G4VisExtent::GetZmin(), and TestOneLine().

Referenced by TestCartGridX(), TestCartGridY(), and TestCartGridZ().

00337 {
00338   if (n1 <= 0 || n2 <= 0) 
00339     G4Exception( "G4GeomTestVolume::TestCartGrid()", "GeomNav0002",
00340                  FatalException, "Arguments n1 and n2 must be >= 1" );
00341     
00342   G4ThreeVector xMin( extent.GetXmin(), extent.GetYmin(),
00343                       extent.GetZmin() );
00344   G4ThreeVector xMax( extent.GetXmax(), extent.GetYmax(),
00345                       extent.GetZmax() );
00346   
00347   G4ThreeVector g1(theG1.unit());
00348   G4ThreeVector g2(theG2.unit());
00349   G4ThreeVector v(theV.unit());
00350   
00351   G4double gMin1 = g1.dot(xMin);
00352   G4double gMax1 = g1.dot(xMax);
00353   
00354   G4double gMin2 = g2.dot(xMin);
00355   G4double gMax2 = g2.dot(xMax);
00356   
00357   G4double delta1 = (gMax1-gMin1)/G4double(n1);
00358   G4double delta2 = (gMax2-gMin2)/G4double(n2);
00359     
00360   G4int i1, i2;
00361   
00362   for(i1=0;i1<=n1;++i1)
00363   {
00364     G4ThreeVector p1 = (gMin1 + G4double(i1)*delta1)*g1;
00365     
00366     for(i2=0;i2<=n2;++i2)
00367     {
00368       G4ThreeVector p2 = (gMin2 + G4double(i2)*delta2)*g2;
00369       
00370       TestOneLine( p1+p2, v );
00371     }
00372   }
00373 }  

void G4GeomTestVolume::TestCartGridX ( G4int  ny = 100,
G4int  nz = 100 
)

Definition at line 99 of file G4GeomTestVolume.cc.

References TestCartGrid().

Referenced by TestCartGridXYZ().

00100 {
00101   TestCartGrid( G4ThreeVector(0,1,0), G4ThreeVector(0,0,1),
00102                 G4ThreeVector(1,0,0), ny, nz );
00103 }

void G4GeomTestVolume::TestCartGridXYZ ( G4int  nx = 100,
G4int  ny = 100,
G4int  nz = 100 
)

Definition at line 89 of file G4GeomTestVolume.cc.

References TestCartGridX(), TestCartGridY(), and TestCartGridZ().

Referenced by TestRecursiveCartGrid().

00090 {
00091   TestCartGridX( ny, nz );
00092   TestCartGridY( nz, nx );
00093   TestCartGridZ( nx, ny );
00094 }

void G4GeomTestVolume::TestCartGridY ( G4int  nz = 100,
G4int  nx = 100 
)

Definition at line 108 of file G4GeomTestVolume.cc.

References TestCartGrid().

Referenced by TestCartGridXYZ().

00109 {
00110   TestCartGrid( G4ThreeVector(0,0,1), G4ThreeVector(1,0,0),
00111                 G4ThreeVector(0,1,0), nz, nx );
00112 }

void G4GeomTestVolume::TestCartGridZ ( G4int  nx = 100,
G4int  ny = 100 
)

Definition at line 117 of file G4GeomTestVolume.cc.

References TestCartGrid().

Referenced by TestCartGridXYZ().

00118 {
00119   TestCartGrid( G4ThreeVector(1,0,0), G4ThreeVector(0,1,0),
00120                 G4ThreeVector(0,0,1), nx, ny );
00121 }

void G4GeomTestVolume::TestCylinder ( G4int  nPhi = 90,
G4int  nZ = 50,
G4int  nRho = 50,
G4double  fracZ = 0.8,
G4double  fracRho = 0.8,
G4bool  usePhi = false 
)

Definition at line 249 of file G4GeomTestVolume.cc.

References extent, G4VisExtent::GetXmax(), G4VisExtent::GetXmin(), G4VisExtent::GetYmax(), G4VisExtent::GetYmin(), G4VisExtent::GetZmax(), G4VisExtent::GetZmin(), G4INCL::Math::pi, and TestOneLine().

Referenced by TestRecursiveCylinder().

00252 {
00253   //
00254   // Get size of our volume
00255   //
00256   G4double xMax = std::max(extent.GetXmax(),-extent.GetXmin());
00257   G4double yMax = std::max(extent.GetYmax(),-extent.GetYmin());
00258   G4double rhoMax = std::sqrt(xMax*xMax + yMax*yMax);
00259   
00260   G4double zMax = extent.GetZmax();
00261   G4double zMin = extent.GetZmin();
00262   G4double zHalfLength = 0.5*(zMax-zMin);
00263   G4double z0 = 0.5*(zMax+zMin);
00264   
00265   //
00266   // Loop over phi
00267   //
00268   G4double deltaPhi = 2*pi/G4double(nPhi);
00269   
00270   G4double phi = 0;
00271   G4int iPhi = nPhi;
00272   if ((iPhi&1) == 0) iPhi++;  // Also use odd number phi slices
00273   do
00274   {
00275     G4double cosPhi = std::cos(phi);
00276     G4double sinPhi = std::sin(phi);
00277     G4ThreeVector vPhi1(sinPhi,-cosPhi,0);
00278 
00279     //
00280     // Loop over rho
00281     //
00282     G4double rho = rhoMax;
00283     G4int iRho = nRho;
00284     do
00285     {
00286       G4ThreeVector p(rho*cosPhi,rho*sinPhi,0);
00287       static G4ThreeVector v(0,0,1);
00288       
00289       TestOneLine( p, v );
00290       
00291       if (usePhi)
00292       {
00293         //
00294         // Loop over z
00295         //
00296         G4double zScale = 1.0;
00297         G4int iZ=nZ;
00298         do
00299         {
00300           p.setZ( z0 + zScale*zHalfLength );
00301           TestOneLine(p,vPhi1);
00302           p.setZ( z0 - zScale*zHalfLength );
00303           TestOneLine(p,vPhi1);
00304         } while( zScale *= fracZ, --iZ );
00305       }
00306     } while( rho *= fracRho, --iRho );
00307 
00308     //
00309     // Loop over z
00310     //
00311     G4ThreeVector p(0,0,0);
00312     G4ThreeVector vPhi2(cosPhi,sinPhi,0);
00313     
00314     G4double zScale = 1.0;
00315     G4int iZ=nZ;
00316     do
00317     {
00318       p.setZ( z0 + zScale*zHalfLength );
00319       
00320       TestOneLine(p,vPhi2);
00321       
00322       p.setZ( z0 - zScale*zHalfLength );
00323       
00324       TestOneLine(p,vPhi2);
00325     } while( zScale *= fracZ, --iZ );
00326     
00327   } while( phi += deltaPhi, --iPhi );
00328 }

void G4GeomTestVolume::TestOneLine ( const G4ThreeVector p,
const G4ThreeVector v 
)

Definition at line 441 of file G4GeomTestVolume.cc.

References G4LogicalVolume::GetDaughter(), G4VPhysicalVolume::GetLogicalVolume(), G4LogicalVolume::GetNoDaughters(), G4LogicalVolume::GetSolid(), logger, CLHEP::detail::n, overlaps, overshoots, G4GeomTestLogger::SolidProblem(), target, and tolerance.

Referenced by TestCartGrid(), TestCylinder(), and TestRecursiveLine().

00443 {
00444   //
00445   // Keep track of intersection points
00446   //
00447   std::vector<G4GeomTestVolPoint> points;
00448   
00449   //
00450   // Calculate intersections with the mother volume
00451   //
00452   G4GeomTestSegment targetSegment( target->GetLogicalVolume()->GetSolid(),
00453                                    p, v, logger );
00454   
00455   //
00456   // Add them to our list
00457   // 
00458   G4int n = targetSegment.GetNumberPoints();
00459   G4int i;
00460   for(i=0;i<n;++i)
00461   {
00462     points.push_back( G4GeomTestVolPoint( targetSegment.GetPoint(i), -1 ) );
00463   } 
00464 
00465   //
00466   // Loop over daughter volumes
00467   //
00468   const G4LogicalVolume *logical = target->GetLogicalVolume();
00469   G4int nDaughter = logical->GetNoDaughters();
00470   G4int iDaughter;
00471   for( iDaughter=0; iDaughter<nDaughter; ++iDaughter)
00472   {
00473     const G4VPhysicalVolume *daughter = 
00474           logical->GetDaughter(iDaughter);
00475     
00476     //
00477     // Convert coordinates to daughter local coordinates
00478     //
00479     const G4RotationMatrix *rotation = daughter->GetFrameRotation();
00480     const G4ThreeVector &translation = daughter->GetFrameTranslation();
00481     
00482     G4ThreeVector pLocal = translation + p;
00483     G4ThreeVector vLocal = v;
00484     
00485     if (rotation)
00486     {
00487       pLocal = (*rotation)*pLocal;
00488       vLocal = (*rotation)*vLocal;
00489     }
00490     
00491     //
00492     // Find intersections
00493     //
00494     G4GeomTestSegment
00495        daughterSegment( daughter->GetLogicalVolume()->GetSolid(), 
00496                         pLocal, vLocal, logger );
00497     
00498     //
00499     // Add them to the list
00500     //
00501     n = daughterSegment.GetNumberPoints();
00502     for(i=0;i<n;++i)
00503     {
00504       points.push_back( G4GeomTestVolPoint( daughterSegment.GetPoint(i),
00505             iDaughter, translation, rotation ) );
00506     } 
00507   }
00508   
00509   //
00510   // Now sort the list of intersection points
00511   //
00512   std::sort( points.begin(), points.end() );
00513   
00514   //
00515   // Search for problems:
00516   //    1. Overlap daughters will be indicated by intersecting
00517   //       points that lie within another daughter
00518   //    2. Daughter extending outside the mother will be
00519   //       indicated by intersecting points outside the mother
00520   //
00521   // The search method is always looking forward when
00522   // resolving discrepencies due to roundoff error. Once
00523   // one instance of a daughter intersection is analyzed,
00524   // it is removed from further consideration
00525   //
00526   n = points.size();
00527   
00528   //
00529   // Set true if this point has been analyzed
00530   //
00531   std::vector<G4bool> checked( n, false );
00532   
00533   for(i=0;i<n;++i)
00534   {
00535     if (checked[i]) continue;
00536   
00537     G4int iDaug = points[i].GetDaughterIndex();
00538     if (iDaug < 0) continue;
00539   
00540     //
00541     // Intersection found. Find matching pair.
00542     //
00543     G4double iS = points[i].GetDistance();
00544     G4int j = i;
00545     while(++j<n)
00546     {
00547       if (iDaug == points[j].GetDaughterIndex()) break;
00548     }
00549     if (j>=n)
00550     {
00551       //
00552       // Unmatched? This shouldn't happen
00553       //
00554       logger->SolidProblem( logical->GetDaughter(iDaug)->
00555                             GetLogicalVolume()->GetSolid(),
00556                             "Unmatched intersection point",
00557           points[i].GetPosition() );
00558       continue;
00559     }
00560     
00561     checked[j] = true;
00562     
00563     G4double jS = points[j].GetDistance();
00564 
00565     //
00566     // Otherwise, we could have a problem
00567     //
00568     G4int k = i;
00569     while(++k<j)
00570     {
00571       if (checked[k]) continue;
00572       
00573       G4bool kEntering = points[k].Entering();
00574       G4double kS = points[k].GetDistance();
00575       
00576       
00577       //
00578       // Problem found: catagorize
00579       //
00580       G4int kDaug = points[k].GetDaughterIndex();
00581       if (kDaug < 0)
00582       {
00583         //
00584         // Ignore small overshoots if they are within tolerance
00585         //
00586         if (kS-iS < tolerance &&   kEntering ) continue;
00587         if (jS-kS < tolerance && (!kEntering)) continue;
00588 
00589         //
00590         // We appear to extend outside the mother volume
00591         //
00592         std::map<G4long,G4GeomTestOvershootList>::iterator overshoot =
00593           overshoots.find(iDaug);
00594         if (overshoot == overshoots.end())
00595         {
00596           std::pair<std::map<G4long,G4GeomTestOvershootList>::iterator,G4bool>
00597             result =
00598               overshoots.insert( std::pair<const G4long,G4GeomTestOvershootList>
00599                                (iDaug,G4GeomTestOvershootList(target,iDaug)) );
00600           overshoot = result.first;
00601         }
00602 
00603         if (kEntering)
00604           (*overshoot).second.AddError( points[i].GetPosition(),
00605                                         points[k].GetPosition() );
00606         else
00607           (*overshoot).second.AddError( points[k].GetPosition(),
00608                                         points[j].GetPosition() );
00609       }
00610       else
00611       {
00612         //
00613         // Ignore small overlaps if they are within tolerance
00614         //
00615         if (kS-iS < tolerance && (!kEntering)) continue;
00616         if (jS-kS < tolerance &&   kEntering ) continue;
00617         
00618         //
00619         // We appear to overlap with another daughter
00620         //
00621         G4long key = iDaug < kDaug ?
00622                (iDaug*nDaughter + kDaug) : (kDaug*nDaughter + iDaug);
00623         
00624         std::map<G4long,G4GeomTestOverlapList>::iterator overlap =
00625           overlaps.find(key);
00626         if (overlap == overlaps.end())
00627         {
00628           std::pair<std::map<G4long,G4GeomTestOverlapList>::iterator,G4bool>
00629             result =
00630             overlaps.insert( std::pair<const G4long,G4GeomTestOverlapList>
00631                            (key,G4GeomTestOverlapList(target,iDaug,kDaug)) );
00632           overlap = result.first;
00633         }
00634 
00635         if (points[k].Entering())
00636           (*overlap).second.AddError( points[k].GetPosition(),
00637                                       points[j].GetPosition() );
00638         else
00639           (*overlap).second.AddError( points[i].GetPosition(),
00640                                       points[k].GetPosition() );
00641       }
00642     }
00643   }
00644 }

void G4GeomTestVolume::TestRecursiveCartGrid ( G4int  nx = 100,
G4int  ny = 100,
G4int  nz = 100,
G4int  sLevel = 0,
G4int  depth = -1 
)

Definition at line 126 of file G4GeomTestVolume.cc.

References G4LogicalVolume::GetDaughter(), G4VPhysicalVolume::GetLogicalVolume(), G4LogicalVolume::GetNoDaughters(), G4VPhysicalVolume::IsReplicated(), logger, ReportErrors(), target, TestCartGridXYZ(), TestRecursiveCartGrid(), and tolerance.

Referenced by TestRecursiveCartGrid().

00128 {
00129   // If reached requested level of depth (i.e. set to 0), exit.
00130   // If not depth specified (i.e. set to -1), visit the whole tree.
00131   // If requested initial level of depth is not zero, visit from beginning
00132   //
00133   if (depth == 0) return;
00134   if (depth != -1) depth--;
00135   if (slevel != 0) slevel--;
00136 
00137   //
00138   // As long as we aren't a replica and we reached the requested
00139   // initial level of depth, test ourselves
00140   //
00141   if ( (!target->IsReplicated()) && (slevel==0) )
00142   {
00143     TestCartGridXYZ( nx, ny, nz );
00144     ReportErrors();
00145   }
00146 
00147   //
00148   // Loop over unique daughters
00149   //
00150   std::set<const G4LogicalVolume *> tested;
00151 
00152   const G4LogicalVolume *logical = target->GetLogicalVolume();
00153   G4int nDaughter = logical->GetNoDaughters();
00154   G4int iDaughter;
00155   for( iDaughter=0; iDaughter<nDaughter; ++iDaughter )
00156   {
00157     const G4VPhysicalVolume *daughter =
00158           logical->GetDaughter(iDaughter);
00159     const G4LogicalVolume *daughterLogical =
00160           daughter->GetLogicalVolume();
00161     
00162     //
00163     // Skip empty daughters
00164     //
00165     if (daughterLogical->GetNoDaughters() == 0) continue;
00166     
00167     //
00168     // Tested already?
00169     //
00170     std::pair<std::set<const G4LogicalVolume *>::iterator,G4bool>
00171            there = tested.insert(daughterLogical);
00172     if (!there.second) continue;
00173 
00174     //
00175     // Recurse
00176     //
00177     G4GeomTestVolume vTest( daughter, logger, tolerance );
00178     vTest.TestRecursiveCartGrid( nx,ny,nz,slevel,depth );
00179   }
00180 }

void G4GeomTestVolume::TestRecursiveCylinder ( G4int  nPhi = 90,
G4int  nZ = 50,
G4int  nRho = 50,
G4double  fracZ = 0.8,
G4double  fracRho = 0.8,
G4bool  usePhi = false,
G4int  sLevel = 0,
G4int  depth = -1 
)

Definition at line 186 of file G4GeomTestVolume.cc.

References G4LogicalVolume::GetDaughter(), G4VPhysicalVolume::GetLogicalVolume(), G4LogicalVolume::GetNoDaughters(), G4VPhysicalVolume::IsReplicated(), logger, ReportErrors(), target, TestCylinder(), TestRecursiveCylinder(), and tolerance.

Referenced by TestRecursiveCylinder().

00190 {
00191   // If reached requested level of depth (i.e. set to 0), exit.
00192   // If not depth specified (i.e. set to -1), visit the whole tree.
00193   // If requested initial level of depth is not zero, visit from beginning
00194   //
00195   if (depth == 0) return;
00196   if (depth != -1) depth--;
00197   if (slevel != 0) slevel--;
00198 
00199   //
00200   // As long as we aren't a replica and we reached the requested
00201   // initial level of depth, test ourselves
00202   //
00203   if ( (!target->IsReplicated()) && (slevel==0) )
00204   {
00205     TestCylinder( nPhi, nZ, nRho, fracZ, fracRho, usePhi );
00206     ReportErrors();
00207   }
00208 
00209   //
00210   // Loop over unique daughters
00211   //
00212   std::set<const G4LogicalVolume *> tested;
00213 
00214   const G4LogicalVolume *logical = target->GetLogicalVolume();
00215   G4int nDaughter = logical->GetNoDaughters();
00216   G4int iDaughter;
00217   for( iDaughter=0; iDaughter<nDaughter; ++iDaughter )
00218   {
00219     const G4VPhysicalVolume *daughter =
00220           logical->GetDaughter(iDaughter);
00221     const G4LogicalVolume *daughterLogical =
00222           daughter->GetLogicalVolume();
00223     
00224     //
00225     // Skip empty daughters
00226     //
00227     if (daughterLogical->GetNoDaughters() == 0) continue;
00228     
00229     //
00230     // Tested already?
00231     //
00232     std::pair<std::set<const G4LogicalVolume *>::iterator,G4bool>
00233            there = tested.insert(daughterLogical);
00234     if (!there.second) continue;
00235 
00236     //
00237     // Recurse
00238     //
00239     G4GeomTestVolume vTest( daughter, logger, tolerance );
00240     vTest.TestRecursiveCylinder( nPhi, nZ, nRho,
00241                                  fracZ, fracRho, usePhi,
00242                                  slevel, depth );
00243   }
00244 }

void G4GeomTestVolume::TestRecursiveLine ( const G4ThreeVector p,
const G4ThreeVector v,
G4int  sLevel = 0,
G4int  depth = -1 
)

Definition at line 379 of file G4GeomTestVolume.cc.

References G4LogicalVolume::GetDaughter(), G4VPhysicalVolume::GetLogicalVolume(), G4LogicalVolume::GetNoDaughters(), G4VPhysicalVolume::IsReplicated(), logger, ReportErrors(), target, TestOneLine(), TestRecursiveLine(), and tolerance.

Referenced by TestRecursiveLine().

00382 {
00383   // If reached requested level of depth (i.e. set to 0), exit.
00384   // If not depth specified (i.e. set to -1), visit the whole tree.
00385   // If requested initial level of depth is not zero, visit from beginning
00386   //
00387   if (depth == 0) return;
00388   if (depth != -1) depth--;
00389   if (slevel != 0) slevel--;
00390 
00391   //
00392   // As long as we aren't a replica and we reached the requested
00393   // initial level of depth, test ourselves
00394   //
00395   if ( (!target->IsReplicated()) && (slevel==0) )
00396   {
00397     TestOneLine( p, v );
00398     ReportErrors();
00399   }
00400 
00401   //
00402   // Loop over unique daughters
00403   //
00404   std::set<const G4LogicalVolume *> tested;
00405 
00406   const G4LogicalVolume *logical = target->GetLogicalVolume();
00407   G4int nDaughter = logical->GetNoDaughters();
00408   G4int iDaughter;
00409   for( iDaughter=0; iDaughter<nDaughter; ++iDaughter )
00410   {
00411     const G4VPhysicalVolume *daughter =
00412           logical->GetDaughter(iDaughter);
00413     const G4LogicalVolume *daughterLogical =
00414           daughter->GetLogicalVolume();
00415     
00416     //
00417     // Skip empty daughters
00418     //
00419     if (daughterLogical->GetNoDaughters() == 0) continue;
00420     
00421     //
00422     // Tested already?
00423     //
00424     std::pair<std::set<const G4LogicalVolume *>::iterator,G4bool>
00425            there = tested.insert(daughterLogical);
00426     if (!there.second) continue;
00427 
00428     //
00429     // Recurse
00430     //
00431     G4GeomTestVolume vTest( daughter, logger, tolerance );
00432     vTest.TestRecursiveLine( p, v, slevel, depth );
00433   }
00434 }


Field Documentation

G4VisExtent G4GeomTestVolume::extent [protected]

Definition at line 165 of file G4GeomTestVolume.hh.

Referenced by TestCartGrid(), and TestCylinder().

G4GeomTestLogger* G4GeomTestVolume::logger [protected]

Definition at line 163 of file G4GeomTestVolume.hh.

Referenced by ReportErrors(), TestOneLine(), TestRecursiveCartGrid(), TestRecursiveCylinder(), and TestRecursiveLine().

std::map<G4long,G4GeomTestOverlapList> G4GeomTestVolume::overlaps [protected]

Definition at line 167 of file G4GeomTestVolume.hh.

Referenced by ClearErrors(), ReportErrors(), and TestOneLine().

std::map<G4long,G4GeomTestOvershootList> G4GeomTestVolume::overshoots [protected]

Definition at line 172 of file G4GeomTestVolume.hh.

Referenced by ClearErrors(), ReportErrors(), and TestOneLine().

const G4VPhysicalVolume* G4GeomTestVolume::target [protected]

Definition at line 162 of file G4GeomTestVolume.hh.

Referenced by TestOneLine(), TestRecursiveCartGrid(), TestRecursiveCylinder(), and TestRecursiveLine().

G4double G4GeomTestVolume::tolerance [protected]

Definition at line 164 of file G4GeomTestVolume.hh.

Referenced by GetTolerance(), SetTolerance(), TestOneLine(), TestRecursiveCartGrid(), TestRecursiveCylinder(), and TestRecursiveLine().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:05 2013 for Geant4 by  doxygen 1.4.7