2 // ********************************************************************
 
    3 // * License and Disclaimer                                           *
 
    5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
 
    6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
 
    7 // * conditions of the Geant4 Software License,  included in the file *
 
    8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
 
    9 // * include a list of copyright holders.                             *
 
   11 // * Neither the authors of this software system, nor their employing *
 
   12 // * institutes,nor the agencies providing financial support for this *
 
   13 // * work  make  any representation or  warranty, express or implied, *
 
   14 // * regarding  this  software system or assume any liability for its *
 
   15 // * use.  Please see the license in the file  LICENSE  and URL above *
 
   16 // * for the full disclaimer and the limitation of liability.         *
 
   18 // * This  code  implementation is the result of  the  scientific and *
 
   19 // * technical work of the GEANT4 collaboration.                      *
 
   20 // * By using,  copying,  modifying or  distributing the software (or *
 
   21 // * any work based  on the software)  you  agree  to acknowledge its *
 
   22 // * use  in  resulting  scientific  publications,  and indicate your *
 
   23 // * acceptance of all terms of the Geant4 Software license.          *
 
   24 // ********************************************************************
 
   27 // $Id: G4AffineTransform.icc 93152 2015-10-08 11:53:57Z gcosmo $
 
   30 // G4AffineTransformation Inline implementation
 
   32 // --------------------------------------------------------------------
 
   34 inline G4AffineTransform::G4AffineTransform()
 
   35  : rxx(1),rxy(0),rxz(0),
 
   41 inline G4AffineTransform::G4AffineTransform(const G4ThreeVector& tlate)
 
   42  : rxx(1),rxy(0),rxz(0),
 
   45    tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
 
   48 inline G4AffineTransform::G4AffineTransform(const G4RotationMatrix& rot)
 
   49  : rxx(rot.xx()),rxy(rot.xy()),rxz(rot.xz()),
 
   50    ryx(rot.yx()),ryy(rot.yy()),ryz(rot.yz()),
 
   51    rzx(rot.zx()),rzy(rot.zy()),rzz(rot.zz()),
 
   55 inline G4AffineTransform::G4AffineTransform( const G4RotationMatrix& rot,
 
   56                                              const G4ThreeVector& tlate   )
 
   57  : rxx(rot.xx()),rxy(rot.xy()),rxz(rot.xz()),
 
   58    ryx(rot.yx()),ryy(rot.yy()),ryz(rot.yz()),
 
   59    rzx(rot.zx()),rzy(rot.zy()),rzz(rot.zz()),
 
   60    tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
 
   63 inline G4AffineTransform::G4AffineTransform( const G4RotationMatrix* rot,
 
   64                                              const G4ThreeVector& tlate)
 
   65  : tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
 
   69       rxx=rot->xx();rxy=rot->xy();rxz=rot->xz();
 
   70       ryx=rot->yx();ryy=rot->yy();ryz=rot->yz();
 
   71       rzx=rot->zx();rzy=rot->zy();rzz=rot->zz();
 
   83 G4AffineTransform( const G4double prxx,const G4double prxy,const G4double prxz,
 
   84                    const G4double pryx,const G4double pryy,const G4double pryz,
 
   85                    const G4double przx,const G4double przy,const G4double przz,
 
   86                    const G4double ptx,const G4double pty,const G4double ptz)
 
   87  : rxx(prxx),rxy(prxy),rxz(prxz),
 
   88    ryx(pryx),ryy(pryy),ryz(pryz),
 
   89    rzx(przx),rzy(przy),rzz(przz),
 
   90    tx(ptx),ty(pty),tz(ptz)
 
   93 inline G4AffineTransform
 
   94 G4AffineTransform::operator * (const G4AffineTransform& tf) const
 
   96         return G4AffineTransform(
 
   97         rxx*tf.rxx+rxy*tf.ryx+rxz*tf.rzx,
 
   98         rxx*tf.rxy+rxy*tf.ryy+rxz*tf.rzy,
 
   99         rxx*tf.rxz+rxy*tf.ryz+rxz*tf.rzz,
 
  101         ryx*tf.rxx+ryy*tf.ryx+ryz*tf.rzx,
 
  102         ryx*tf.rxy+ryy*tf.ryy+ryz*tf.rzy,
 
  103         ryx*tf.rxz+ryy*tf.ryz+ryz*tf.rzz,
 
  105         rzx*tf.rxx+rzy*tf.ryx+rzz*tf.rzx,
 
  106         rzx*tf.rxy+rzy*tf.ryy+rzz*tf.rzy,
 
  107         rzx*tf.rxz+rzy*tf.ryz+rzz*tf.rzz,
 
  109         tx*tf.rxx+ty*tf.ryx+tz*tf.rzx+tf.tx,
 
  110         tx*tf.rxy+ty*tf.ryy+tz*tf.rzy+tf.ty,
 
  111         tx*tf.rxz+ty*tf.ryz+tz*tf.rzz+tf.tz);
 
  114 inline G4AffineTransform&
 
  115 G4AffineTransform::operator *= (const G4AffineTransform& tf)
 
  117          // Use temporaries for `in place' compound transform computation
 
  119         G4double nrxx=rxx*tf.rxx+rxy*tf.ryx+rxz*tf.rzx;
 
  120         G4double nrxy=rxx*tf.rxy+rxy*tf.ryy+rxz*tf.rzy;
 
  121         G4double nrxz=rxx*tf.rxz+rxy*tf.ryz+rxz*tf.rzz;
 
  123         G4double nryx=ryx*tf.rxx+ryy*tf.ryx+ryz*tf.rzx;
 
  124         G4double nryy=ryx*tf.rxy+ryy*tf.ryy+ryz*tf.rzy;
 
  125         G4double nryz=ryx*tf.rxz+ryy*tf.ryz+ryz*tf.rzz;
 
  127         G4double nrzx=rzx*tf.rxx+rzy*tf.ryx+rzz*tf.rzx;
 
  128         G4double nrzy=rzx*tf.rxy+rzy*tf.ryy+rzz*tf.rzy;
 
  129         G4double nrzz=rzx*tf.rxz+rzy*tf.ryz+rzz*tf.rzz;
 
  131         G4double ntx=tx*tf.rxx+ty*tf.ryx+tz*tf.rzx+tf.tx;
 
  132         G4double nty=tx*tf.rxy+ty*tf.ryy+tz*tf.rzy+tf.ty;
 
  133         G4double ntz=tx*tf.rxz+ty*tf.ryz+tz*tf.rzz+tf.tz;
 
  135         tx=ntx; ty=nty; tz=ntz;
 
  136         rxx=nrxx; rxy=nrxy; rxz=nrxz;
 
  137         ryx=nryx; ryy=nryy; ryz=nryz;
 
  138         rzx=nrzx; rzy=nrzy; rzz=nrzz;
 
  143 inline G4AffineTransform&
 
  144 G4AffineTransform::Product(const G4AffineTransform& tf1,
 
  145                            const G4AffineTransform& tf2)
 
  147         rxx=tf1.rxx*tf2.rxx + tf1.rxy*tf2.ryx + tf1.rxz*tf2.rzx;
 
  148         rxy=tf1.rxx*tf2.rxy + tf1.rxy*tf2.ryy + tf1.rxz*tf2.rzy;
 
  149         rxz=tf1.rxx*tf2.rxz + tf1.rxy*tf2.ryz + tf1.rxz*tf2.rzz;
 
  151         ryx=tf1.ryx*tf2.rxx + tf1.ryy*tf2.ryx + tf1.ryz*tf2.rzx;
 
  152         ryy=tf1.ryx*tf2.rxy + tf1.ryy*tf2.ryy + tf1.ryz*tf2.rzy;
 
  153         ryz=tf1.ryx*tf2.rxz + tf1.ryy*tf2.ryz + tf1.ryz*tf2.rzz;
 
  155         rzx=tf1.rzx*tf2.rxx + tf1.rzy*tf2.ryx + tf1.rzz*tf2.rzx;
 
  156         rzy=tf1.rzx*tf2.rxy + tf1.rzy*tf2.ryy + tf1.rzz*tf2.rzy;
 
  157         rzz=tf1.rzx*tf2.rxz + tf1.rzy*tf2.ryz + tf1.rzz*tf2.rzz;
 
  159         tx=tf1.tx*tf2.rxx + tf1.ty*tf2.ryx + tf1.tz*tf2.rzx   + tf2.tx;
 
  160         ty=tf1.tx*tf2.rxy + tf1.ty*tf2.ryy + tf1.tz*tf2.rzy   + tf2.ty;
 
  161         tz=tf1.tx*tf2.rxz + tf1.ty*tf2.ryz + tf1.tz*tf2.rzz   + tf2.tz; 
 
  166 inline G4AffineTransform&
 
  167 G4AffineTransform::InverseProduct( const G4AffineTransform& tf1,
 
  168                                    const G4AffineTransform& tf2)
 
  170         G4double itf2tx = - tf2.tx*tf2.rxx - tf2.ty*tf2.rxy - tf2.tz*tf2.rxz;
 
  171         G4double itf2ty = - tf2.tx*tf2.ryx - tf2.ty*tf2.ryy - tf2.tz*tf2.ryz;
 
  172         G4double itf2tz = - tf2.tx*tf2.rzx - tf2.ty*tf2.rzy - tf2.tz*tf2.rzz;
 
  174         rxx=tf1.rxx*tf2.rxx+tf1.rxy*tf2.rxy+tf1.rxz*tf2.rxz;
 
  175         rxy=tf1.rxx*tf2.ryx+tf1.rxy*tf2.ryy+tf1.rxz*tf2.ryz;
 
  176         rxz=tf1.rxx*tf2.rzx+tf1.rxy*tf2.rzy+tf1.rxz*tf2.rzz;
 
  178         ryx=tf1.ryx*tf2.rxx+tf1.ryy*tf2.rxy+tf1.ryz*tf2.rxz;
 
  179         ryy=tf1.ryx*tf2.ryx+tf1.ryy*tf2.ryy+tf1.ryz*tf2.ryz;
 
  180         ryz=tf1.ryx*tf2.rzx+tf1.ryy*tf2.rzy+tf1.ryz*tf2.rzz;
 
  182         rzx=tf1.rzx*tf2.rxx+tf1.rzy*tf2.rxy+tf1.rzz*tf2.rxz;
 
  183         rzy=tf1.rzx*tf2.ryx+tf1.rzy*tf2.ryy+tf1.rzz*tf2.ryz;
 
  184         rzz=tf1.rzx*tf2.rzx+tf1.rzy*tf2.rzy+tf1.rzz*tf2.rzz;
 
  186         tx=tf1.tx*tf2.rxx+tf1.ty*tf2.rxy+tf1.tz*tf2.rxz+itf2tx;
 
  187         ty=tf1.tx*tf2.ryx+tf1.ty*tf2.ryy+tf1.tz*tf2.ryz+itf2ty;
 
  188         tz=tf1.tx*tf2.rzx+tf1.ty*tf2.rzy+tf1.tz*tf2.rzz+itf2tz;
 
  194 G4ThreeVector G4AffineTransform::TransformPoint(const G4ThreeVector& vec) const
 
  196         return G4ThreeVector( vec.x()*rxx + vec.y()*ryx + vec.z()*rzx   + tx,
 
  197                               vec.x()*rxy + vec.y()*ryy + vec.z()*rzy   + ty,
 
  198                               vec.x()*rxz + vec.y()*ryz + vec.z()*rzz   + tz  );
 
  202 G4ThreeVector G4AffineTransform::TransformAxis(const G4ThreeVector& axis) const
 
  204         return G4ThreeVector( axis.x()*rxx + axis.y()*ryx + axis.z()*rzx,
 
  205                               axis.x()*rxy + axis.y()*ryy + axis.z()*rzy,
 
  206                               axis.x()*rxz + axis.y()*ryz + axis.z()*rzz  );
 
  210 void G4AffineTransform::ApplyPointTransform(G4ThreeVector& vec) const
 
  212         G4double x = vec.x()*rxx + vec.y()*ryx + vec.z()*rzx    + tx;
 
  213         G4double y = vec.x()*rxy + vec.y()*ryy + vec.z()*rzy    + ty;
 
  214         G4double z = vec.x()*rxz + vec.y()*ryz + vec.z()*rzz    + tz;
 
  222 void G4AffineTransform::ApplyAxisTransform(G4ThreeVector& axis) const
 
  224         G4double x = axis.x()*rxx + axis.y()*ryx + axis.z()*rzx;
 
  225         G4double y = axis.x()*rxy + axis.y()*ryy + axis.z()*rzy;
 
  226         G4double z = axis.x()*rxz + axis.y()*ryz + axis.z()*rzz;
 
  234 G4AffineTransform G4AffineTransform::Inverse() const
 
  236         return G4AffineTransform( rxx, ryx, rzx,
 
  240                                  -tx*rxx - ty*rxy - tz*rxz,
 
  241                                  -tx*ryx - ty*ryy - tz*ryz,
 
  242                                  -tx*rzx - ty*rzy - tz*rzz  );
 
  246 G4AffineTransform& G4AffineTransform::Invert()
 
  248         G4double v1 = -tx*rxx - ty*rxy - tz*rxz;
 
  249         G4double v2 = -tx*ryx - ty*ryy - tz*ryz;
 
  250         G4double v3 = -tx*rzx - ty*rzy - tz*rzz;
 
  254         G4double tmp1=ryx; ryx=rxy; rxy=tmp1;
 
  255         G4double tmp2=rzx; rzx=rxz; rxz=tmp2;
 
  256         G4double tmp3=rzy; rzy=ryz; ryz=tmp3;
 
  263 G4AffineTransform& G4AffineTransform::operator +=(const G4ThreeVector& tlate)
 
  273 G4AffineTransform& G4AffineTransform::operator -=(const G4ThreeVector& tlate)
 
  283 G4bool G4AffineTransform::operator == (const G4AffineTransform& tf) const
 
  285         return (tx==tf.tx&&ty==tf.ty&&tz==tf.tz&&
 
  286                 rxx==tf.rxx&&rxy==tf.rxy&&rxz==tf.rxz&&
 
  287                 ryx==tf.ryx&&ryy==tf.ryy&&ryz==tf.ryz&&
 
  288                 rzx==tf.rzx&&rzy==tf.rzy&&rzz==tf.rzz) ? true : false;
 
  291 G4bool G4AffineTransform::operator != (const G4AffineTransform& tf) const
 
  293         return (tx!=tf.tx||ty!=tf.ty||tz!=tf.tz||
 
  294                 rxx!=tf.rxx||rxy!=tf.rxy||rxz!=tf.rxz||
 
  295                 ryx!=tf.ryx||ryy!=tf.ryy||ryz!=tf.ryz||
 
  296                 rzx!=tf.rzx||rzy!=tf.rzy||rzz!=tf.rzz) ? true : false;
 
  300 G4double G4AffineTransform::operator [] (const G4int n) const
 
  353 G4bool G4AffineTransform::IsRotated() const
 
  355         return (rxx==1.0 && ryy==1.0 && rzz==1.0) ? false : true;
 
  359 G4bool G4AffineTransform::IsTranslated() const
 
  361         return (tx || ty || tz) ? true:false;
 
  364 inline G4RotationMatrix G4AffineTransform::NetRotation() const
 
  366   G4RotationMatrix mat;
 
  367   return mat.rotateAxes(G4ThreeVector(rxx,ryx,rzx),
 
  368                         G4ThreeVector(rxy,ryy,rzy),
 
  369                         G4ThreeVector(rxz,ryz,rzz));
 
  373 G4ThreeVector G4AffineTransform::NetTranslation() const
 
  375         return G4ThreeVector(tx,ty,tz);
 
  379 void G4AffineTransform::SetNetRotation(const G4RotationMatrix& rot)
 
  393 void G4AffineTransform::SetNetTranslation(const G4ThreeVector& tlate)
 
  401 std::ostream& operator << (std::ostream& os, const G4AffineTransform& transf)
 
  403   std::streamsize oldPrec = os.precision(6); 
 
  404   G4double DeviationTolerance = 1.0e-05;
 
  406   G4double diagDeviation = 0.0; 
 
  407   G4double offdDeviationUL = 0.0;
 
  408   G4double offdDeviationDR = 0.0;
 
  409   G4double offdDeviation = 0.0; 
 
  411   os << "  Transformation: " << G4endl; 
 
  413   // double  a = std::max ( 1, 2, 3 ) ; 
 
  415   G4bool UnitTr = ! transf.IsRotated(); 
 
  416   diagDeviation = std::max( std::fabs( transf[0] - 1.0 ) ,   //  |rxx - 1|
 
  417                             std::fabs( transf[5] - 1.0 ) );  //  |ryy - 1|
 
  418   diagDeviation = std::max( diagDeviation, 
 
  419                             std::fabs( transf[10] - 1.0 ) ); //  |rzz - 1|
 
  421   offdDeviationUL = std::max( std::fabs( transf[1] ) ,       //  |rxy|
 
  422                               std::fabs( transf[2] ) );      //  |rxz|
 
  423   offdDeviationUL = std::max( offdDeviationUL, 
 
  424                               std::fabs( transf[4] ) );      //  |ryx|
 
  426   offdDeviationDR = std::max( std::fabs( transf[6] ) ,       //  |ryz|
 
  427                               std::fabs( transf[8] ) );      //  |rzx|
 
  428   offdDeviationDR = std::max( offdDeviationDR, 
 
  429                               std::fabs( transf[9] ) );      //  |rzy|
 
  430   offdDeviation = std::max( offdDeviationUL, offdDeviationDR ); 
 
  432   if( UnitTr || std::max(diagDeviation, offdDeviation) < DeviationTolerance )
 
  434      os << "    UNIT  Rotation " << G4endl;
 
  439         << transf[0]  << " " << transf[1] << " " << transf[2] << G4endl
 
  441         << transf[4]  << " " << transf[5] << " " << transf[6] << G4endl
 
  443         << transf[8]  << " " << transf[9] << " " << transf[10] << G4endl;
 
  446   os << "tr/x,y,z: " << transf[12]  << " " << transf[13] << " " << transf[14]
 
  449   os.precision(oldPrec);