include/Geant4/G4VTwistSurface.icc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // G4VTwistSurface class inline methods
0027 //
0028 // 01-Aug-2002 - Kotoyo Hoshina (hoshina@hepburn.s.chiba-u.ac.jp), created.
0029 // 13-Nov-2003 - O.Link (Oliver.Link@cern.ch), Integration in Geant4
0030 //               from original version in Jupiter-2.5.02 application.
0031 // --------------------------------------------------------------------
0032
0033 //=====================================================================
0034 //* DistanceToPlaneWithV ----------------------------------------------
0035
0036 inline
0037 G4double G4VTwistSurface::DistanceToPlaneWithV(const G4ThreeVector& p,
0038                                           const G4ThreeVector& v,
0039                                           const G4ThreeVector& x0,
0040                                           const G4ThreeVector& n0,
0041                                                 G4ThreeVector& xx)
0042 {
0043    G4double q = n0 * v;
0044    G4double t = kInfinity;
0045    if (q != 0.0) { t = (n0 * (x0 - p)) / q; }
0046    xx = p + t * v;
0047    return t;
0048 }
0049
0050 //=====================================================================
0051 //* DistanceToPlane ---------------------------------------------------
0052
0053 inline
0054 G4double G4VTwistSurface::DistanceToPlane(const G4ThreeVector& p,
0055                                      const G4ThreeVector& x0,
0056                                      const G4ThreeVector& n0,
0057                                            G4ThreeVector& xx)
0058 {
0059    // DistanceToPlane :
0060    // Calculate distance to plane in local coordinate,
0061    // then return distance and global intersection points.
0062    //
0063    // p          - location of flying particle
0064    // x0         - reference point of surface
0065    // xx         - a foot of perpendicular line from p to the plane
0066    // t          - distance from xx to p
0067    // n          - a unit normal of this plane from plane to p.
0068    //
0069    // equation of plane:
0070    //      n*(x - x0) = 0;
0071    //
0072    // vector to xx:
0073    //      xx = p - t*n
0074    //
0075    //         where
0076    //         t = n * (p - x0) / std::fabs(n)
0077    //
0078    G4double t;
0079    G4ThreeVector n = n0.unit();
0080    t = n * (p - x0);
0081    xx = p - t * n;
0082    return t;
0083 }
0084
0085 //=====================================================================
0086 //* DistanceToPlane ---------------------------------------------------
0087
0088 inline
0089 G4double G4VTwistSurface::DistanceToPlane(const G4ThreeVector& p,
0090                                           const G4ThreeVector& x0,
0091                                           const G4ThreeVector& t1,
0092                                           const G4ThreeVector& t2,
0093                                                 G4ThreeVector& xx,
0094                                                 G4ThreeVector& n)
0095 {
0096    // DistanceToPlane :
0097    // Calculate distance to plane in local coordinate,
0098    // then return distance and global intersection points.
0099    // t1         - 1st. vector lying on the plane
0100    // t2         - 2nd. vector lying on the plane
0101
0102    n = (t1.cross(t2)).unit();
0103    return DistanceToPlane(p, x0, n, xx);
0104 }
0105
0106 //=====================================================================
0107 //* DistanceToLine ----------------------------------------------------
0108
0109 inline
0110 G4double G4VTwistSurface::DistanceToLine(const G4ThreeVector& p,
0111                                          const G4ThreeVector& x0,
0112                                          const G4ThreeVector& d,
0113                                                G4ThreeVector& xx)
0114 {
0115    // DistanceToLine :
0116    // Calculate distance to line,
0117    // then return distance and global intersection points.
0118    //
0119    // p          - location of flying particle
0120    // x0         - reference point of line
0121    // d          - direction vector of line
0122    // xx         - a foot of perpendicular line from p to the plane
0123    // t          - distance from xx to p
0124    //
0125    // Equation
0126    //
0127    //    distance^2 = |(xx - p)|^2
0128    //    with
0129    //       xx = x0 + t*d
0130    //
0131    //   (d/dt)distance^2 = (d/dt)|((x0 + t*d) - p)|^2
0132    //                    = 2*t*|d|^2 + 2*d*(x0 - p)
0133    //                    = 0  // smallest distance
0134    //   then
0135    //      t = - d*(x0 - p) / |d|^2
0136    //
0137
0138    G4double t;
0139    G4ThreeVector dir = d.unit();
0140    t  = - dir * (x0 - p);      // |dir|^2 = 1.
0141    xx = x0 + t * dir;
0142
0143    G4ThreeVector dist = xx - p;
0144    return dist.mag();
0145 }
0146
0147 //=====================================================================
0148 //* IsAxis0 -----------------------------------------------------------
0149
0150 inline
0151 G4bool G4VTwistSurface::IsAxis0(G4int areacode) const
0152 {
0153    return (areacode & sAxis0) != 0;
0154 }
0155
0156 //=====================================================================
0157 //* IsAxis1 -----------------------------------------------------------
0158
0159 inline
0160 G4bool G4VTwistSurface::IsAxis1(G4int areacode) const
0161 {
0162    return (areacode & sAxis1) != 0;
0163 }
0164
0165 //=====================================================================
0166 //* IsOutside ---------------------------------------------------------
0167
0168 inline
0169 G4bool G4VTwistSurface::IsOutside(G4int areacode) const
0170 {
0171    return (areacode & sInside) == 0;
0172 }
0173
0174 //=====================================================================
0175 //* IsInside ----------------------------------------------------------
0176
0177 inline
0178 G4bool G4VTwistSurface::IsInside(G4int areacode, G4bool testbitmode) const
0179 {
0180    if ((areacode & sInside) != 0) {
0181       if (testbitmode) {
0182          return true;
0183       } else {
0184        if (((areacode & sBoundary) == 0) && ((areacode & sCorner) == 0)) return true;
0185       }
0186    }
0187    return false;
0188 }
0189
0190 //=====================================================================
0191 //* IsBoundary --------------------------------------------------------
0192
0193 inline
0194 G4bool G4VTwistSurface::IsBoundary(G4int areacode, G4bool testbitmode) const
0195 {
0196    if ((areacode & sBoundary) == sBoundary) {
0197       if (testbitmode) {
0198          return true;
0199       } else {
0200          if ((areacode & sInside) == sInside) return true;
0201       }
0202    }
0203    return false;
0204 }
0205
0206 //=====================================================================
0207 //* IsCorner ----------------------------------------------------------
0208
0209 inline
0210 G4bool G4VTwistSurface::IsCorner(G4int areacode, G4bool testbitmode) const
0211 {
0212    if ((areacode & sCorner) == sCorner) {
0213       if (testbitmode) {
0214          return true;
0215       } else {
0216          if ((areacode & sInside) == sInside) return true;
0217       }
0218    }
0219    return false;
0220 }
0221
0222 //=====================================================================
0223 //* GetAxisType -------------------------------------------------------
0224
0225 inline
0226 G4int G4VTwistSurface::GetAxisType(G4int areacode, G4int whichaxis) const
0227 {
0228    G4int axiscode = areacode & sAxisMask & whichaxis;
0229
0230    if (axiscode == (sAxisX & sAxis0) ||
0231        axiscode == (sAxisX & sAxis1)) {
0232       return sAxisX;
0233    } else if (axiscode == (sAxisY & sAxis0) ||
0234               axiscode == (sAxisY & sAxis1)) {
0235       return sAxisY;
0236    } else if (axiscode == (sAxisZ & sAxis0) ||
0237               axiscode == (sAxisZ & sAxis1)) {
0238       return sAxisZ;
0239    } else if (axiscode == (sAxisRho & sAxis0) ||
0240               axiscode == (sAxisRho & sAxis1)) {
0241       return sAxisRho;
0242    } else if (axiscode == (sAxisPhi & sAxis0) ||
0243               axiscode == (sAxisPhi & sAxis1)) {
0244       return sAxisPhi;
0245    } else {
0246       std::ostringstream message;
0247       message << "Configuration not supported." << G4endl
0248               << "        areacode = " << areacode;
0249       G4Exception("G4VTwistSurface::GetAxisType()","GeomSolids0001",
0250                   FatalException, message);
0251    }
0252    return 1;
0253 }
0254
0255 //=====================================================================
0256 //* ComputeGlobalPoint ------------------------------------------------
0257
0258 inline
0259 G4ThreeVector G4VTwistSurface::ComputeGlobalPoint(const G4ThreeVector& lp) const
0260 {
0261    return fRot * G4ThreeVector(lp) + fTrans;
0262 }
0263
0264 //=====================================================================
0265 //* ComputeGlobalPoint ------------------------------------------------
0266
0267 inline
0268 G4ThreeVector G4VTwistSurface::ComputeLocalPoint(const G4ThreeVector& gp) const
0269 {
0270    return fRot.inverse() * ( G4ThreeVector(gp) - fTrans ) ;
0271 }
0272
0273 //=====================================================================
0274 //* ComputeGlobalDirection --------------------------------------------
0275
0276 inline G4ThreeVector
0277 G4VTwistSurface::ComputeGlobalDirection(const G4ThreeVector& lp) const
0278 {
0279    return fRot * G4ThreeVector(lp);
0280 }
0281
0282 //=====================================================================
0283 //* ComputeLocalDirection ---------------------------------------------
0284
0285 inline G4ThreeVector
0286 G4VTwistSurface::ComputeLocalDirection(const G4ThreeVector& gp) const
0287 {
0288    return fRot.inverse() * G4ThreeVector(gp);
0289 }
0290
0291 //=====================================================================
0292 //* SetNeighbours -----------------------------------------------------
0293
0294 inline void
0295 G4VTwistSurface::SetNeighbours(G4VTwistSurface* ax0min, G4VTwistSurface* ax1min,
0296                                G4VTwistSurface* ax0max, G4VTwistSurface* ax1max)
0297 {
0298    fNeighbours[0] = ax0min;
0299    fNeighbours[1] = ax1min;
0300    fNeighbours[2] = ax0max;
0301    fNeighbours[3] = ax1max;
0302 }
0303
0304 //=====================================================================
0305 //* GetNeighbours -----------------------------------------------------
0306
0307 inline G4int
0308 G4VTwistSurface::GetNeighbours(G4int areacode, G4VTwistSurface** surfaces)
0309 {
0310
0311   G4int sAxis0Min = sAxis0 & sAxisMin ;
0312   G4int sAxis1Min = sAxis1 & sAxisMin ;
0313   G4int sAxis0Max = sAxis0 & sAxisMax ;
0314   G4int sAxis1Max = sAxis1 & sAxisMax ;
0315
0316   G4int i = 0;
0317
0318   if ( (areacode & sAxis0Min ) == sAxis0Min )
0319   {
0320     surfaces[i] = fNeighbours[0] ;
0321     ++i ;
0322   }
0323
0324   if ( ( areacode & sAxis1Min ) == sAxis1Min )
0325   {
0326      surfaces[i] = fNeighbours[1] ;
0327     ++i ;
0328     if ( i == 2 ) return i ;
0329   }
0330
0331   if ( ( areacode & sAxis0Max ) == sAxis0Max )
0332   {
0333     surfaces[i] = fNeighbours[2] ;
0334     ++i ;
0335     if ( i == 2 ) return i ;
0336   }
0337
0338   if ( ( areacode & sAxis1Max ) == sAxis1Max )
0339   {
0340     surfaces[i] = fNeighbours[3] ;
0341     ++i ;
0342     if ( i == 2 ) return i ;
0343   }
0344
0345   return i ;
0346 }
0347
0348 //=====================================================================
0349 //* GetCorner ---------------------------------------------------------
0350
0351 inline
0352 G4ThreeVector G4VTwistSurface::GetCorner(G4int areacode) const
0353 {
0354    if ((areacode & sCorner) == 0)
0355    {
0356       std::ostringstream message;
0357       message << "Area code must represent corner." << G4endl
0358               << "        areacode = " << areacode;
0359       G4Exception("G4VTwistSurface::GetCorner()","GeomSolids0002",
0360                   FatalException, message);
0361    }
0362
0363    if ((areacode & sC0Min1Min) == sC0Min1Min) {
0364       return fCorners[0];
0365    } else if ((areacode & sC0Max1Min) == sC0Max1Min) {
0366       return fCorners[1];
0367    } else if ((areacode & sC0Max1Max) == sC0Max1Max) {
0368       return fCorners[2];
0369    } else if ((areacode & sC0Min1Max) == sC0Min1Max) {
0370       return fCorners[3];
0371    } else {
0372       std::ostringstream message;
0373       message << "Configuration not supported." << G4endl
0374               << "        areacode = " << areacode;
0375       G4Exception("G4VTwistSurface::GetCorner()", "GeomSolids0001",
0376                   FatalException, message);
0377    }
0378    return fCorners[0];
0379 }