Ignore:
Timestamp:
Jun 18, 2010, 11:42:07 AM (15 years ago)
Author:
garnier
Message:

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08

Location:
trunk/source/geometry/solids/specific/src
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/source/geometry/solids/specific/src/G4ExtrudedSolid.cc

    r1228 r1315  
    2525//
    2626//
    27 // $Id: G4ExtrudedSolid.cc,v 1.18 2008/10/30 11:47:45 ivana Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4ExtrudedSolid.cc,v 1.19 2010/04/15 10:23:34 ivana Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030//
    … …  
    4040#include <algorithm>
    4141#include <cmath>
     42#include <iomanip>
    4243
    4344#include "G4ExtrudedSolid.hh"
    … …  
    9495        "Z-sections with the same z position are not supported.");
    9596    }
    96   }     
    97 
     97  } 
     98 
     99  // Check if polygon vertices are defined clockwise
     100  // (the area is positive if polygon vertices are defined anti-clockwise)
     101  //
     102  G4double area = 0.;
     103  for ( G4int i=0; i<fNv; ++i ) {
     104    G4int j = i+1;
     105    if ( j == fNv ) j = 0;
     106    area += 0.5 * ( polygon[i].x()*polygon[j].y() - polygon[j].x()*polygon[i].y());
     107  }
     108 
    98109  // Copy polygon
    99110  //
    100   for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[i]); }
     111  if  ( area < 0. ) {   
     112    // Polygon vertices are defined clockwise, we just copy the polygon       
     113    for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[i]); }
     114  }
     115  else {
     116    // Polygon vertices are defined anti-clockwise, we revert them
     117    //G4Exception("G4ExtrudedSolid::G4ExtrudedSolid()", errorDescription,
     118    //            JustWarning,
     119    //            "Polygon vertices defined anti-clockwise, reverting polygon");     
     120    for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[fNv-i-1]); }
     121  }
     122   
    101123 
    102124  // Copy z-sections
    … …  
    148170  }
    149171     
     172  // Check if polygon vertices are defined clockwise
     173  // (the area is positive if polygon vertices are defined anti-clockwise)
     174 
     175  G4double area = 0.;
     176  for ( G4int i=0; i<fNv; ++i ) {
     177    G4int j = i+1;
     178    if ( j == fNv ) j = 0;
     179    area += 0.5 * ( polygon[i].x()*polygon[j].y() - polygon[j].x()*polygon[i].y());
     180  }
     181 
    150182  // Copy polygon
    151183  //
    152   for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[i]); }
     184  if  ( area < 0. ) {   
     185    // Polygon vertices are defined clockwise, we just copy the polygon       
     186    for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[i]); }
     187  }
     188  else {
     189    // Polygon vertices are defined anti-clockwise, we revert them
     190    //G4Exception("G4ExtrudedSolid::G4ExtrudedSolid()", errorDescription,
     191    //            JustWarning,
     192    //            "Polygon vertices defined anti-clockwise, reverting polygon");     
     193    for ( G4int i=0; i<fNv; ++i ) { fPolygon.push_back(polygon[fNv-i-1]); }
     194  }
    153195 
    154196  // Copy z-sections
    … …  
    453495    //
    454496    G4double angle = GetAngle(c2->first, c3->first, c1->first);
    455 
    456     if ( angle > pi )
     497    //G4cout << "angle " << angle  << G4endl;
     498
     499    G4int counter = 0;
     500    while ( angle > pi )
    457501    {
    458502      // G4cout << "Skipping concave vertex " << c2->second << G4endl;
    … …  
    467511      // G4cout << "Looking at triangle : "
    468512      //        << c1->second << "  " << c2->second
    469       //        << "  " << c3->second << G4endl; 
    470 
     513      //        << "  " << c3->second << G4endl;
     514     
     515      angle = GetAngle(c2->first, c3->first, c1->first);
     516      //G4cout << "angle " << angle  << G4endl;
     517     
     518      counter++;
     519     
     520      if ( counter > fNv) {
     521        G4Exception("G4ExtrudedSolid::AddGeneralPolygonFacets", "InvalidSetup" ,
     522                    FatalException, "Triangularisation has failed.");
     523        break;
     524      } 
    471525    }
    472526
    … …  
    734788}
    735789
    736 
    737790//_____________________________________________________________________________
    738791
    … …  
    751804  for ( G4int i=0; i<fNv; ++i )
    752805  {
    753     os << "   vx = " << fPolygon[i].x()/mm << " mm"
     806    os << std::setw(5) << "#" << i
     807       << "   vx = " << fPolygon[i].x()/mm << " mm"
    754808       << "   vy = " << fPolygon[i].y()/mm << " mm" << G4endl;
    755809  }
    … …  
    764818  }     
    765819
     820/*
     821  // Triangles (for debogging)
     822  os << G4endl;
     823  os << " Triangles:" << G4endl;
     824  os << " Triangle #   vertex1   vertex2   vertex3" << G4endl;
     825
     826  G4int counter = 0;
     827  std::vector< std::vector<G4int> >::const_iterator it;
     828  for ( it = fTriangles.begin(); it != fTriangles.end(); it++ ) {
     829     std::vector<G4int> triangle = *it;
     830     os << std::setw(10) << counter++
     831        << std::setw(10) << triangle[0] << std::setw(10)  << triangle[1]  << std::setw(10)  << triangle[2]
     832        << G4endl;
     833  }         
     834*/
    766835  return os;
    767836} 

  • trunk/source/geometry/solids/specific/src/G4PolyconeSide.cc

    r1228 r1315  
    2525//
    2626//
    27 // $Id: G4PolyconeSide.cc,v 1.22 2009/11/11 12:23:37 gcosmo Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4PolyconeSide.cc,v 1.24 2010/02/24 11:31:49 gcosmo Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030//
    … …  
    6767  kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
    6868  fSurfaceArea = 0.0;
     69  fPhi.first = G4ThreeVector(0,0,0);
     70  fPhi.second= 0.0;
    6971
    7072  //
    … …  
    492494  if (phiIsOpen)
    493495  {
    494     G4double phi = axis.phi();
     496    G4double phi = GetPhi(axis);
    495497    while( phi < startPhi ) phi += twopi;
    496498   
    … …  
    853855}
    854856
     857//
     858// GetPhi
     859//
     860// Calculate Phi for a given 3-vector (point), if not already cached for the
     861// same point, in the attempt to avoid consecutive computation of the same
     862// quantity
     863//
     864G4double G4PolyconeSide::GetPhi( const G4ThreeVector& p )
     865{
     866  G4double val=0.;
     867
     868  if (fPhi.first != p)
     869  {
     870    val = p.phi();
     871    fPhi.first = p;
     872    fPhi.second = val;
     873  }
     874  else
     875  {
     876    val = fPhi.second;
     877  }
     878  return val;
     879}
    855880
    856881//
    … …  
    922947    // Finally, check phi
    923948    //
    924     G4double phi = p.phi();
     949    G4double phi = GetPhi(p);
    925950    while( phi < startPhi ) phi += twopi;
    926951   
    … …  
    9781003    // PolyPhiFace. See PolyPhiFace::InsideEdgesExact
    9791004    //
    980     G4double phi = hit.phi();
     1005    G4double phi = GetPhi(hit);
    9811006    while( phi < startPhi-phiTolerant ) phi += twopi;
    9821007   

  • trunk/source/geometry/solids/specific/src/G4PolyhedraSide.cc

    r1228 r1315  
    2525//
    2626//
    27 // $Id: G4PolyhedraSide.cc,v 1.15 2008/05/15 11:41:59 gcosmo Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4PolyhedraSide.cc,v 1.17 2010/02/24 11:31:49 gcosmo Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030//
    … …  
    6464                                        G4bool isAllBehind )
    6565{
    66 
    6766  kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
    6867  fSurfaceArea=0.;
     68  fPhi.first = G4ThreeVector(0,0,0);
     69  fPhi.second= 0.0;
     70
    6971  //
    7072  // Record values
    … …  
    561563  // Try the closest phi segment first
    562564  //
    563   G4int iPhi = ClosestPhiSegment( p.phi() );
     565  G4int iPhi = ClosestPhiSegment( GetPhi(p) );
    564566 
    565567  G4ThreeVector pdotc = p - vecs[iPhi].center;
    … …  
    594596  // Which phi segment is closest to this point?
    595597  //
    596   G4int iPhi = ClosestPhiSegment( p.phi() );
     598  G4int iPhi = ClosestPhiSegment( GetPhi(p) );
    597599 
    598600  G4double norm;
    … …  
    624626  // Which phi segment is closest to this point?
    625627  //
    626   G4int iPhi = ClosestPhiSegment( p.phi() );
     628  G4int iPhi = ClosestPhiSegment( GetPhi(p) );
    627629
    628630  //
    … …  
    656658  // Which phi segment, if any, does the axis belong to
    657659  //
    658   iPhi = PhiSegment( axis.phi() );
     660  iPhi = PhiSegment( GetPhi(axis) );
    659661 
    660662  if (iPhi < 0)
    … …  
    971973
    972974//
     975// GetPhi
     976//
     977// Calculate Phi for a given 3-vector (point), if not already cached for the
     978// same point, in the attempt to avoid consecutive computation of the same
     979// quantity
     980//
     981G4double G4PolyhedraSide::GetPhi( const G4ThreeVector& p )
     982{
     983  G4double val=0.;
     984
     985  if (fPhi.first != p)
     986  {
     987    val = p.phi();
     988    fPhi.first = p;
     989    fPhi.second = val;
     990  }
     991  else
     992  {
     993    val = fPhi.second;
     994  }
     995  return val;
     996}
     997
     998
     999//
    9731000// DistanceToOneSide
    9741001//

  • trunk/source/geometry/solids/specific/src/G4SolidExtentList.cc

    r1228 r1315  
    2525//
    2626//
    27 // $Id: G4SolidExtentList.cc,v 1.5 2007/05/11 13:54:29 gcosmo Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4SolidExtentList.cc,v 1.6 2010/02/10 16:38:37 gcosmo Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030//
    … …  
    5151  axis = kZAxis;
    5252  limited = false;
    53   minLimit = -DBL_MAX;
    54   maxLimit = +DBL_MAX;
     53  minLimit = -kInfinity;
     54  maxLimit = +kInfinity;
    5555}
    5656
    … …  
    7272  else
    7373  {
    74     minLimit = -DBL_MAX;
    75     maxLimit = +DBL_MAX;
     74    minLimit = -kInfinity;
     75    maxLimit = +kInfinity;
    7676  }
    7777}

  • trunk/source/geometry/solids/specific/src/G4TessellatedSolid.cc

    r1228 r1315  
    2525// ********************************************************************
    2626//
    27 // $Id: G4TessellatedSolid.cc,v 1.19 2009/04/27 08:06:27 gcosmo Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4TessellatedSolid.cc,v 1.20 2010/04/28 16:21:21 flei Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    … …  
    4242// CHANGE HISTORY
    4343// --------------
     44//
     45// 12 April 2010      P R Truscott, QinetiQ, bug fixes to treat optical
     46//                    photon transport, in particular internal reflection
     47//                    at surface.
    4448//
    4549// 14 November 2007   P R Truscott, QinetiQ & Stan Seibert, U Texas
    … …  
    670674    if ((*f)->Intersect(p,v,false,dist,distFromSurface,normal))
    671675    {
     676//
     677//
     678// Set minDist to the new distance to current facet if distFromSurface is in
     679// positive direction and point is not at surface.  If the point is within
     680// 0.5*kCarTolerance of the surface, then force distance to be zero and
     681// leave member function immediately (for efficiency), as proposed by & credit
     682// to Akira Okumura.
     683//
    672684      if (distFromSurface > 0.5*kCarTolerance && dist >= 0.0 && dist < minDist)
    673685      {
    674686        minDist  = dist;
     687      }
     688      else if (-0.5*kCarTolerance <= dist && dist <= 0.5*kCarTolerance)
     689      {
     690        return 0.0;
    675691      }
    676692    }

  • trunk/source/geometry/solids/specific/src/G4TriangularFacet.cc

    r1228 r1315  
    2525// ********************************************************************
    2626//
    27 // $Id: G4TriangularFacet.cc,v 1.12 2008/11/13 08:25:07 gcosmo Exp $
    28 // GEANT4 tag $Name: geant4-09-03 $
     27// $Id: G4TriangularFacet.cc,v 1.13 2010/04/28 16:21:21 flei Exp $
     28// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
    2929//
    3030// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    … …  
    5656//                  plane of the triangle.
    5757//
     58// 12 April 2010    P R Truscott, QinetiQ, bug fixes to treat optical
     59//                  photon transport, in particular internal reflection
     60//                  at surface.
    5861// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    5962
    … …  
    350353//
    351354//
    352 // Do a heck for rounding errors in the distance-squared.
     355// Do a check for rounding errors in the distance-squared.  It appears that
     356// the conventional methods for calculating sqrDist breaks down when very near
     357// to or at the surface (as required by transport).  We'll therefore also use
     358// the magnitude-squared of the vector displacement.  (Note that I've also
     359// tried to get around this problem by using the existing equations for
     360//
     361//    sqrDist = function(a,b,c,d,s,t)
     362//
     363// and use a more accurate addition process which minimises errors and
     364// breakdown of cummutitivity [where (A+B)+C != A+(B+C)] but this still
     365// doesn't work.
     366// Calculation from u = D + s*E[0] + t*E[1] is less efficient, but appears
     367// more robust.
    353368//
    354369  if (sqrDist < 0.0) { sqrDist = 0.0; }
    355 
    356   return D + s*E[0] + t*E[1];
     370  G4ThreeVector u = D + s*E[0] + t*E[1];
     371  G4double u2     = u.mag2();
     372//
     373//
     374// The following (part of the roundoff error check) is from Oliver Merle's
     375// updates.
     376//
     377  if ( sqrDist > u2 ) sqrDist = u2;
     378
     379  return u;
    357380}
    358381
    … …  
    542565// we intersect.
    543566//
    544       distance = -0.5*kCarTolerance;
     567//      distance = -0.5*kCarTolerance;
     568      distance = 0.0;
    545569      normal   = surfaceNormal;
    546570      return true;
    … …  
    728752  return area;
    729753}
     754////////////////////////////////////////////////////////////////////////
     755//
Note: See TracChangeset for help on using the changeset viewer.