Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:59:14

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 // G4TwistBoxSide
0027 //
0028 // Class description:
0029 //
0030 // G4TwistBoxSide describes a twisted boundary surface for a trapezoid.
0031 
0032 // Author: 27-Oct-2004 - O.Link (Oliver.Link@cern.ch)
0033 // --------------------------------------------------------------------
0034 #ifndef G4TWISTBOXSIDE_HH
0035 #define G4TWISTBOXSIDE_HH
0036 
0037 #include "G4VTwistSurface.hh"
0038 
0039 #include <vector>
0040 
0041 class G4TwistBoxSide : public G4VTwistSurface
0042 {
0043   public:
0044    
0045     G4TwistBoxSide(const G4String& name,
0046                          G4double  PhiTwist, // twist angle
0047                          G4double  pDz,      // half z lenght
0048                          G4double  pTheta,   // direction between end planes
0049                          G4double  pPhi,     // by polar and azimutal angles
0050                          G4double  pDy1,     // half y length at -pDz
0051                          G4double  pDx1,     // half x length at -pDz,-pDy
0052                          G4double  pDx2,     // half x length at -pDz,+pDy
0053                          G4double  pDy2,     // half y length at +pDz
0054                          G4double  pDx3,     // half x length at +pDz,-pDy
0055                          G4double  pDx4,     // half x length at +pDz,+pDy
0056                          G4double  pAlph,    // tilt angle at +pDz
0057                          G4double  AngleSide // parity
0058                    );
0059   
0060     ~G4TwistBoxSide() override;
0061    
0062     G4ThreeVector  GetNormal(const G4ThreeVector& xx,
0063                                    G4bool isGlobal = false) override ;   
0064    
0065     G4int DistanceToSurface(const G4ThreeVector& gp,
0066                             const G4ThreeVector& gv,
0067                                   G4ThreeVector  gxx[],
0068                                   G4double  distance[],
0069                                   G4int     areacode[],
0070                                   G4bool    isvalid[],
0071                             EValidate validate = kValidateWithTol) override;
0072                                                   
0073     G4int DistanceToSurface(const G4ThreeVector& gp,
0074                                   G4ThreeVector  gxx[],
0075                                   G4double       distance[],
0076                                   G4int          areacode[]) override;
0077 
0078     G4TwistBoxSide(__void__&);
0079       // Fake default constructor for usage restricted to direct object
0080       // persistency for clients requiring preallocation of memory for
0081       // persistifiable objects.
0082 
0083   private:
0084 
0085     G4int GetAreaCode(const G4ThreeVector& xx, 
0086                             G4bool withTol = true) override;
0087     void SetCorners() override;
0088     void SetBoundaries() override;
0089 
0090     void GetPhiUAtX(const G4ThreeVector& p, G4double& phi, G4double& u);
0091     G4ThreeVector ProjectPoint(const G4ThreeVector& p,
0092                                      G4bool isglobal = false);
0093 
0094     inline G4ThreeVector SurfacePoint(G4double phi, G4double u,
0095                                       G4bool isGlobal = false) override;
0096     inline G4double GetBoundaryMin(G4double phi) override;
0097     inline G4double GetBoundaryMax(G4double phi) override;
0098     inline G4double GetSurfaceArea() override;
0099     void GetFacets( G4int m, G4int n, G4double xyz[][3],
0100                     G4int faces[][4], G4int iside ) override;
0101 
0102     inline G4double GetValueA(G4double phi);
0103     inline G4double GetValueB(G4double phi);
0104     inline G4ThreeVector NormAng(G4double phi, G4double u);
0105     inline G4double Xcoef(G4double u, G4double phi);
0106       // To calculate the w(u) function
0107 
0108   private:
0109 
0110     G4double fTheta;   
0111     G4double fPhi ;
0112 
0113     G4double fDy1;   
0114     G4double fDx1;     
0115     G4double fDx2;     
0116 
0117     G4double fDy2;   
0118     G4double fDx3;     
0119     G4double fDx4;     
0120 
0121     G4double fDz;         // Half-length along the z axis
0122 
0123     G4double fAlph;
0124     G4double fTAlph;      // std::tan(fAlph)
0125     
0126     G4double fPhiTwist;   // twist angle ( dphi in surface equation)
0127 
0128     G4double fAngleSide;
0129 
0130     G4double fdeltaX;
0131     G4double fdeltaY;
0132 
0133     G4double fDx4plus2;   // fDx4 + fDx2  == a2/2 + a1/2
0134     G4double fDx4minus2;  // fDx4 - fDx2          -
0135     G4double fDx3plus1;   // fDx3 + fDx1  == d2/2 + d1/2
0136     G4double fDx3minus1;  // fDx3 - fDx1          -
0137     G4double fDy2plus1;   // fDy2 + fDy1  == b2/2 + b1/2
0138     G4double fDy2minus1;  // fDy2 - fDy1          -
0139     G4double fa1md1;      // 2 fDx2 - 2 fDx1  == a1 - d1
0140     G4double fa2md2;      // 2 fDx4 - 2 fDx3
0141 };   
0142 
0143 //========================================================
0144 // inline functions
0145 //========================================================
0146 
0147 inline
0148 G4double G4TwistBoxSide::GetValueA(G4double phi)
0149 {
0150   return ( fDx4plus2 + fDx4minus2 * ( 2 * phi ) / fPhiTwist  ) ;
0151 }
0152 
0153 
0154 inline 
0155 G4double G4TwistBoxSide::GetValueB(G4double phi) 
0156 {
0157   return ( fDy2plus1 + fDy2minus1 * ( 2 * phi ) / fPhiTwist ) ;
0158 }
0159 
0160 inline
0161 G4double G4TwistBoxSide::Xcoef(G4double u, G4double phi)
0162 {
0163   
0164   return GetValueA(phi)/2. + u*fTAlph    ;
0165 
0166 }
0167 
0168 inline G4ThreeVector
0169 G4TwistBoxSide::SurfacePoint( G4double phi, G4double u, G4bool isGlobal ) 
0170 {
0171   // function to calculate a point on the surface, given by parameters phi,u
0172 
0173   G4ThreeVector SurfPoint ( Xcoef(u,phi) * std::cos(phi)
0174                           - u * std::sin(phi) + fdeltaX*phi/fPhiTwist,
0175                             Xcoef(u,phi) * std::sin(phi)
0176                           + u * std::cos(phi) + fdeltaY*phi/fPhiTwist,
0177                             2*fDz*phi/fPhiTwist  );
0178 
0179   if (isGlobal) { return (fRot * SurfPoint + fTrans); }
0180   return SurfPoint;
0181 }
0182 
0183 inline
0184 G4double G4TwistBoxSide::GetBoundaryMin(G4double phi)
0185 {
0186   return -0.5*GetValueB(phi) ;
0187 }
0188 
0189 inline
0190 G4double G4TwistBoxSide::GetBoundaryMax(G4double phi)
0191 {
0192   return 0.5*GetValueB(phi) ;
0193 }
0194 
0195 inline
0196 G4double G4TwistBoxSide::GetSurfaceArea()
0197 {
0198   return (fDz*(std::sqrt(16*fDy1*fDy1
0199          + (fa1md1 + 4*fDy1*fTAlph)*(fa1md1 + 4*fDy1*fTAlph))
0200          + std::sqrt(16*fDy1*fDy1 + (fa2md2 + 4*fDy1*fTAlph)
0201          * (fa2md2 + 4*fDy1*fTAlph))))/2. ;
0202 }
0203 
0204 inline
0205 G4ThreeVector G4TwistBoxSide::NormAng( G4double phi, G4double u ) 
0206 {
0207   // function to calculate the norm at a given point on the surface
0208   // replace a1-d1
0209 
0210   G4ThreeVector nvec( 4*fDz*(std::cos(phi) + fTAlph*std::sin(phi)) ,
0211                       4*fDz*(-(fTAlph*std::cos(phi)) + std::sin(phi)),
0212                       (fDx2 + fDx4)*fPhiTwist*fTAlph
0213                                    + 2*fDx4minus2*(-1 + fTAlph*phi)
0214                                    + 2*fPhiTwist*(1 + fTAlph*fTAlph)*u
0215                                 - 2*(fdeltaX - fdeltaY*fTAlph)*std::cos(phi)
0216                                 - 2*(fdeltaY + fdeltaX*fTAlph)*std::sin(phi) );
0217   return nvec.unit();
0218 }
0219 
0220 #endif