G4QuadrangularFacet.hh

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration and of QinetiQ Ltd,   *
00020 // * subject to DEFCON 705 IPR conditions.                            *
00021 // * By using,  copying,  modifying or  distributing the software (or *
00022 // * any work based  on the software)  you  agree  to acknowledge its *
00023 // * use  in  resulting  scientific  publications,  and indicate your *
00024 // * acceptance of all terms of the Geant4 Software license.          *
00025 // ********************************************************************
00026 //
00027 // $Id: G4QuadrangularFacet.hh 67011 2013-01-29 16:17:41Z gcosmo $
00028 //
00029 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
00030 //
00031 // Class G4QuadrangularFacet
00032 //
00033 // Class description:
00034 //
00035 //   The G4QuadrangularFacet class is used for the contruction of
00036 //   G4TessellatedSolid.
00037 //   It is defined by four fVertices, which shall be in the same plane and be
00038 //   supplied in anti-clockwise order looking from the outsider of the solid
00039 //   where it belongs. Its constructor
00040 //   
00041 //     G4QuadrangularFacet (const G4ThreeVector Pt0, const G4ThreeVector vt1,
00042 //                          const G4ThreeVector vt2, const G4ThreeVector vt3,
00043 //                          G4FacetVertexType);
00044 //
00045 //   takes 5 parameters to define the four fVertices:
00046 //     1) G4FacetvertexType = "ABSOLUTE": in this case Pt0, vt1, vt2 and vt3
00047 //        are the four fVertices required in anti-clockwise order when looking
00048 //        from the outsider.
00049 //     2) G4FacetvertexType = "RELATIVE": in this case the first vertex is Pt0,
00050 //        the second vertex is Pt0+vt, the third vertex is Pt0+vt2 and 
00051 //        the fourth vertex is Pt0+vt3, in anti-clockwise order when looking 
00052 //        from the outsider.
00053 
00054 // CHANGE HISTORY
00055 // --------------
00056 //
00057 // 31 October 2004, P R Truscott, QinetiQ Ltd, UK - Created.
00058 // 12 October 2012, M Gayer, CERN, - Reviewed optimized implementation.
00059 //
00061 #ifndef G4QuadrangularFacet_HH
00062 #define G4QuadrangularFacet_HH 1
00063 
00064 #include "G4VFacet.hh"
00065 #include "G4Types.hh"
00066 #include "G4ThreeVector.hh"
00067 #include "G4TriangularFacet.hh"
00068 
00069 class G4QuadrangularFacet : public G4VFacet
00070 {
00071   public:  // with description
00072 
00073     G4QuadrangularFacet (const G4ThreeVector &Pt0, const G4ThreeVector &vt1,
00074                          const G4ThreeVector &vt2, const G4ThreeVector &vt3,
00075                                G4FacetVertexType);
00076     G4QuadrangularFacet (const G4QuadrangularFacet &right);
00077    ~G4QuadrangularFacet ();
00078 
00079     G4QuadrangularFacet &operator=(const G4QuadrangularFacet &right);    
00080 
00081     G4VFacet *GetClone ();
00082 
00083     G4ThreeVector Distance (const G4ThreeVector &p);
00084     G4double Distance (const G4ThreeVector &p, G4double minDist);
00085     G4double Distance (const G4ThreeVector &p, G4double minDist,
00086                        const G4bool outgoing);
00087     G4double Extent   (const G4ThreeVector axis);
00088     G4bool Intersect  (const G4ThreeVector &p, const G4ThreeVector &v,
00089                        const G4bool outgoing, G4double &distance,
00090                              G4double &distFromSurface, G4ThreeVector &normal);
00091     G4ThreeVector GetSurfaceNormal () const;
00092 
00093     G4double GetArea ();
00094     G4ThreeVector GetPointOnFace () const;
00095 
00096     G4GeometryType GetEntityType () const;
00097 
00098     inline G4bool IsDefined () const;
00099     inline G4int GetNumberOfVertices () const;
00100     inline G4ThreeVector GetVertex (G4int i) const;
00101     inline void SetVertex (G4int i, const G4ThreeVector &val);
00102     inline void SetVertices(std::vector<G4ThreeVector> *v);
00103 
00104     inline G4double GetRadius () const;
00105     inline G4ThreeVector GetCircumcentre () const;
00106 
00107   private:
00108 
00109     inline G4int GetVertexIndex (G4int i) const;
00110     inline void SetVertexIndex (G4int i, G4int val);
00111 
00112     inline G4int AllocatedMemory();
00113 
00114   private:
00115 
00116     G4double fRadius;
00117     G4ThreeVector fCircumcentre;
00118 
00119     G4TriangularFacet fFacet1, fFacet2;
00120 };
00121 
00123 // Inlined Methods
00125 
00126 inline G4int G4QuadrangularFacet::GetNumberOfVertices () const
00127 {
00128   return 4;
00129 }
00130 
00131 inline G4ThreeVector G4QuadrangularFacet::GetVertex (G4int i) const
00132 {
00133   return i == 3 ? fFacet2.GetVertex(2) : fFacet1.GetVertex(i);
00134 }
00135 
00136 
00137 inline G4double G4QuadrangularFacet::GetRadius () const
00138 {
00139   return fRadius;
00140 }
00141 
00142 inline G4ThreeVector G4QuadrangularFacet::GetCircumcentre () const
00143 {
00144   return fCircumcentre;
00145 }
00146 
00147 inline void G4QuadrangularFacet::SetVertex (G4int i, const G4ThreeVector &val)
00148 {
00149   switch (i)
00150   {
00151     case 0:
00152       fFacet1.SetVertex(0, val);
00153       fFacet2.SetVertex(0, val);
00154       break;
00155     case 1:
00156       fFacet1.SetVertex(1, val);
00157       break;
00158     case 2:
00159       fFacet1.SetVertex(2, val);
00160       fFacet2.SetVertex(1, val);
00161       break;
00162     case 3:
00163       fFacet2.SetVertex(2, val);
00164       break;
00165     }
00166 }
00167 
00168 inline void G4QuadrangularFacet::SetVertices(std::vector<G4ThreeVector> *v)
00169 {
00170   fFacet1.SetVertices(v);
00171   fFacet2.SetVertices(v);
00172 }
00173 
00174 inline G4bool G4QuadrangularFacet::IsDefined () const
00175 {
00176   return fFacet1.IsDefined();
00177 }
00178 
00179 inline G4int G4QuadrangularFacet::GetVertexIndex (G4int i) const
00180 {
00181   return i == 3 ? fFacet2.GetVertexIndex(2) : fFacet1.GetVertexIndex(i);
00182 }
00183 
00184 
00185 inline void G4QuadrangularFacet::SetVertexIndex (G4int i, G4int val)
00186 {
00187   switch (i)
00188   {
00189     case 0:
00190       fFacet1.SetVertexIndex(0, val);
00191       fFacet2.SetVertexIndex(0, val);
00192       break;
00193     case 1:
00194       fFacet1.SetVertexIndex(1, val);
00195       break;
00196     case 2:
00197       fFacet1.SetVertexIndex(2, val);
00198       fFacet2.SetVertexIndex(1, val);
00199       break;
00200     case 3:
00201       fFacet2.SetVertexIndex(2, val);
00202       break;
00203   }
00204 }
00205 
00206 inline G4int G4QuadrangularFacet::AllocatedMemory()
00207 {
00208   return sizeof(*this) + fFacet1.AllocatedMemory() + fFacet2.AllocatedMemory();
00209 }
00210 
00211 #endif

Generated on Mon May 27 17:49:40 2013 for Geant4 by  doxygen 1.4.7