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// ********************************************************************
26// G4AffineTransformation Inline implementation
28// --------------------------------------------------------------------
30inline G4AffineTransform::G4AffineTransform()
31 : rxx(1),rxy(0),rxz(0),
37inline G4AffineTransform::G4AffineTransform(const G4ThreeVector& tlate)
38 : rxx(1),rxy(0),rxz(0),
41 tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
44inline G4AffineTransform::G4AffineTransform(const G4RotationMatrix& rot)
45 : rxx(rot.xx()),rxy(rot.xy()),rxz(rot.xz()),
46 ryx(rot.yx()),ryy(rot.yy()),ryz(rot.yz()),
47 rzx(rot.zx()),rzy(rot.zy()),rzz(rot.zz()),
51inline G4AffineTransform::G4AffineTransform( const G4RotationMatrix& rot,
52 const G4ThreeVector& tlate )
53 : rxx(rot.xx()),rxy(rot.xy()),rxz(rot.xz()),
54 ryx(rot.yx()),ryy(rot.yy()),ryz(rot.yz()),
55 rzx(rot.zx()),rzy(rot.zy()),rzz(rot.zz()),
56 tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
59inline G4AffineTransform::G4AffineTransform( const G4RotationMatrix* rot,
60 const G4ThreeVector& tlate )
61 : tx(tlate.x()),ty(tlate.y()),tz(tlate.z())
65 rxx=rot->xx();rxy=rot->xy();rxz=rot->xz();
66 ryx=rot->yx();ryy=rot->yy();ryz=rot->yz();
67 rzx=rot->zx();rzy=rot->zy();rzz=rot->zz();
79G4AffineTransform( const G4double prxx,const G4double prxy,const G4double prxz,
80 const G4double pryx,const G4double pryy,const G4double pryz,
81 const G4double przx,const G4double przy,const G4double przz,
82 const G4double ptx,const G4double pty,const G4double ptz )
83 : rxx(prxx),rxy(prxy),rxz(prxz),
84 ryx(pryx),ryy(pryy),ryz(pryz),
85 rzx(przx),rzy(przy),rzz(przz),
86 tx(ptx),ty(pty),tz(ptz)
89inline G4AffineTransform&
90G4AffineTransform::operator = (const G4AffineTransform& rhs)
92 if (this == &rhs) { return *this; }
93 rxx = rhs.rxx; rxy = rhs.rxy; rxz = rhs.rxz;
94 ryx = rhs.ryx; ryy = rhs.ryy; ryz = rhs.ryz;
95 rzx = rhs.rzx; rzy = rhs.rzy; rzz = rhs.rzz;
96 tx = rhs.tx; ty = rhs.ty; tz = rhs.tz;
100inline G4AffineTransform
101G4AffineTransform::operator * (const G4AffineTransform& tf) const
103 return G4AffineTransform(rxx*tf.rxx+rxy*tf.ryx+rxz*tf.rzx,
104 rxx*tf.rxy+rxy*tf.ryy+rxz*tf.rzy,
105 rxx*tf.rxz+rxy*tf.ryz+rxz*tf.rzz,
107 ryx*tf.rxx+ryy*tf.ryx+ryz*tf.rzx,
108 ryx*tf.rxy+ryy*tf.ryy+ryz*tf.rzy,
109 ryx*tf.rxz+ryy*tf.ryz+ryz*tf.rzz,
111 rzx*tf.rxx+rzy*tf.ryx+rzz*tf.rzx,
112 rzx*tf.rxy+rzy*tf.ryy+rzz*tf.rzy,
113 rzx*tf.rxz+rzy*tf.ryz+rzz*tf.rzz,
115 tx*tf.rxx+ty*tf.ryx+tz*tf.rzx+tf.tx,
116 tx*tf.rxy+ty*tf.ryy+tz*tf.rzy+tf.ty,
117 tx*tf.rxz+ty*tf.ryz+tz*tf.rzz+tf.tz);
120inline G4AffineTransform&
121G4AffineTransform::operator *= (const G4AffineTransform& tf)
123 // Use temporaries for `in place' compound transform computation
125 G4double nrxx=rxx*tf.rxx+rxy*tf.ryx+rxz*tf.rzx;
126 G4double nrxy=rxx*tf.rxy+rxy*tf.ryy+rxz*tf.rzy;
127 G4double nrxz=rxx*tf.rxz+rxy*tf.ryz+rxz*tf.rzz;
129 G4double nryx=ryx*tf.rxx+ryy*tf.ryx+ryz*tf.rzx;
130 G4double nryy=ryx*tf.rxy+ryy*tf.ryy+ryz*tf.rzy;
131 G4double nryz=ryx*tf.rxz+ryy*tf.ryz+ryz*tf.rzz;
133 G4double nrzx=rzx*tf.rxx+rzy*tf.ryx+rzz*tf.rzx;
134 G4double nrzy=rzx*tf.rxy+rzy*tf.ryy+rzz*tf.rzy;
135 G4double nrzz=rzx*tf.rxz+rzy*tf.ryz+rzz*tf.rzz;
137 G4double ntx=tx*tf.rxx+ty*tf.ryx+tz*tf.rzx+tf.tx;
138 G4double nty=tx*tf.rxy+ty*tf.ryy+tz*tf.rzy+tf.ty;
139 G4double ntz=tx*tf.rxz+ty*tf.ryz+tz*tf.rzz+tf.tz;
141 tx=ntx; ty=nty; tz=ntz;
142 rxx=nrxx; rxy=nrxy; rxz=nrxz;
143 ryx=nryx; ryy=nryy; ryz=nryz;
144 rzx=nrzx; rzy=nrzy; rzz=nrzz;
149inline G4AffineTransform&
150G4AffineTransform::Product( const G4AffineTransform& tf1,
151 const G4AffineTransform& tf2 )
153 rxx = tf1.rxx*tf2.rxx + tf1.rxy*tf2.ryx + tf1.rxz*tf2.rzx;
154 rxy = tf1.rxx*tf2.rxy + tf1.rxy*tf2.ryy + tf1.rxz*tf2.rzy;
155 rxz = tf1.rxx*tf2.rxz + tf1.rxy*tf2.ryz + tf1.rxz*tf2.rzz;
157 ryx = tf1.ryx*tf2.rxx + tf1.ryy*tf2.ryx + tf1.ryz*tf2.rzx;
158 ryy = tf1.ryx*tf2.rxy + tf1.ryy*tf2.ryy + tf1.ryz*tf2.rzy;
159 ryz = tf1.ryx*tf2.rxz + tf1.ryy*tf2.ryz + tf1.ryz*tf2.rzz;
161 rzx = tf1.rzx*tf2.rxx + tf1.rzy*tf2.ryx + tf1.rzz*tf2.rzx;
162 rzy = tf1.rzx*tf2.rxy + tf1.rzy*tf2.ryy + tf1.rzz*tf2.rzy;
163 rzz = tf1.rzx*tf2.rxz + tf1.rzy*tf2.ryz + tf1.rzz*tf2.rzz;
165 tx = tf1.tx*tf2.rxx + tf1.ty*tf2.ryx + tf1.tz*tf2.rzx + tf2.tx;
166 ty = tf1.tx*tf2.rxy + tf1.ty*tf2.ryy + tf1.tz*tf2.rzy + tf2.ty;
167 tz = tf1.tx*tf2.rxz + tf1.ty*tf2.ryz + tf1.tz*tf2.rzz + tf2.tz;
172inline G4AffineTransform&
173G4AffineTransform::InverseProduct( const G4AffineTransform& tf1,
174 const G4AffineTransform& tf2 )
176 if ((tf2.rxx + tf2.ryy + tf2.rzz) == 3.)
190 tx = tf1.tx - tf2.tx;
191 ty = tf1.ty - tf2.ty;
192 tz = tf1.tz - tf2.tz;
196 G4double tf1rxx = tf1.rxx, tf1rxy = tf1.rxy, tf1rxz = tf1.rxz;
197 rxx = tf1rxx*tf2.rxx + tf1rxy*tf2.rxy + tf1rxz*tf2.rxz;
198 rxy = tf1rxx*tf2.ryx + tf1rxy*tf2.ryy + tf1rxz*tf2.ryz;
199 rxz = tf1rxx*tf2.rzx + tf1rxy*tf2.rzy + tf1rxz*tf2.rzz;
201 G4double tf1ryx = tf1.ryx, tf1ryy = tf1.ryy, tf1ryz = tf1.ryz;
202 ryx = tf1ryx*tf2.rxx + tf1ryy*tf2.rxy + tf1ryz*tf2.rxz;
203 ryy = tf1ryx*tf2.ryx + tf1ryy*tf2.ryy + tf1ryz*tf2.ryz;
204 ryz = tf1ryx*tf2.rzx + tf1ryy*tf2.rzy + tf1ryz*tf2.rzz;
206 G4double tf1rzx = tf1.rzx, tf1rzy = tf1.rzy, tf1rzz = tf1.rzz;
207 rzx = tf1rzx*tf2.rxx + tf1rzy*tf2.rxy + tf1rzz*tf2.rxz;
208 rzy = tf1rzx*tf2.ryx + tf1rzy*tf2.ryy + tf1rzz*tf2.ryz;
209 rzz = tf1rzx*tf2.rzx + tf1rzy*tf2.rzy + tf1rzz*tf2.rzz;
211 G4double tf1_2tx = tf1.tx - tf2.tx;
212 G4double tf1_2ty = tf1.ty - tf2.ty;
213 G4double tf1_2tz = tf1.tz - tf2.tz;
214 tx = tf1_2tx*tf2.rxx + tf1_2ty*tf2.rxy + tf1_2tz*tf2.rxz;
215 ty = tf1_2tx*tf2.ryx + tf1_2ty*tf2.ryy + tf1_2tz*tf2.ryz;
216 tz = tf1_2tx*tf2.rzx + tf1_2ty*tf2.rzy + tf1_2tz*tf2.rzz;
222G4AffineTransform::TransformPoint(const G4ThreeVector& vec) const
224 G4double vecx = vec.x(), vecy = vec.y(), vecz = vec.z();
225 return G4ThreeVector( vecx*rxx + vecy*ryx + vecz*rzx + tx,
226 vecx*rxy + vecy*ryy + vecz*rzy + ty,
227 vecx*rxz + vecy*ryz + vecz*rzz + tz );
231G4AffineTransform::InverseTransformPoint(const G4ThreeVector& vec) const
233 G4double vecx = vec.x()-tx, vecy = vec.y()-ty, vecz = vec.z()-tz;
234 return G4ThreeVector( vecx*rxx + vecy*rxy + vecz*rxz,
235 vecx*ryx + vecy*ryy + vecz*ryz,
236 vecx*rzx + vecy*rzy + vecz*rzz );
240G4AffineTransform::TransformAxis(const G4ThreeVector& axis) const
242 G4double axisx = axis.x(), axisy = axis.y(), axisz = axis.z();
243 return G4ThreeVector( axisx*rxx + axisy*ryx + axisz*rzx,
244 axisx*rxy + axisy*ryy + axisz*rzy,
245 axisx*rxz + axisy*ryz + axisz*rzz );
249G4AffineTransform::InverseTransformAxis(const G4ThreeVector& axis) const
251 G4double axisx = axis.x(), axisy = axis.y(), axisz = axis.z();
252 return G4ThreeVector( axisx*rxx + axisy*rxy + axisz*rxz,
253 axisx*ryx + axisy*ryy + axisz*ryz,
254 axisx*rzx + axisy*rzy + axisz*rzz );
258void G4AffineTransform::ApplyPointTransform(G4ThreeVector& vec) const
260 G4double vecx = vec.x(), vecy = vec.y(), vecz = vec.z();
261 vec.setX( vecx*rxx + vecy*ryx + vecz*rzx + tx );
262 vec.setY( vecx*rxy + vecy*ryy + vecz*rzy + ty );
263 vec.setZ( vecx*rxz + vecy*ryz + vecz*rzz + tz );
267void G4AffineTransform::ApplyAxisTransform(G4ThreeVector& axis) const
269 G4double axisx = axis.x(), axisy = axis.y(), axisz = axis.z();
270 axis.setX( axisx*rxx + axisy*ryx + axisz*rzx );
271 axis.setY( axisx*rxy + axisy*ryy + axisz*rzy );
272 axis.setZ( axisx*rxz + axisy*ryz + axisz*rzz );
276G4AffineTransform G4AffineTransform::Inverse() const
278 G4double ttx = -tx, tty = -ty, ttz = -tz;
279 return G4AffineTransform( rxx, ryx, rzx,
282 ttx*rxx + tty*rxy + ttz*rxz,
283 ttx*ryx + tty*ryy + ttz*ryz,
284 ttx*rzx + tty*rzy + ttz*rzz );
288G4AffineTransform& G4AffineTransform::Invert()
290 G4double ttx = -tx, tty = -ty, ttz = -tz;
291 tx = ttx*rxx + tty*rxy + ttz*rxz;
292 ty = ttx*ryx + tty*ryy + ttz*ryz;
293 tz = ttx*rzx + tty*rzy + ttz*rzz;
295 G4double tmp1=ryx; ryx=rxy; rxy=tmp1;
296 G4double tmp2=rzx; rzx=rxz; rxz=tmp2;
297 G4double tmp3=rzy; rzy=ryz; ryz=tmp3;
303G4AffineTransform& G4AffineTransform::operator +=(const G4ThreeVector& tlate)
313G4AffineTransform& G4AffineTransform::operator -=(const G4ThreeVector& tlate)
323G4bool G4AffineTransform::operator == (const G4AffineTransform& tf) const
325 return (tx==tf.tx&&ty==tf.ty&&tz==tf.tz&&
326 rxx==tf.rxx&&rxy==tf.rxy&&rxz==tf.rxz&&
327 ryx==tf.ryx&&ryy==tf.ryy&&ryz==tf.ryz&&
328 rzx==tf.rzx&&rzy==tf.rzy&&rzz==tf.rzz) ? true : false;
332G4bool G4AffineTransform::operator != (const G4AffineTransform& tf) const
334 return (tx!=tf.tx||ty!=tf.ty||tz!=tf.tz||
335 rxx!=tf.rxx||rxy!=tf.rxy||rxz!=tf.rxz||
336 ryx!=tf.ryx||ryy!=tf.ryy||ryz!=tf.ryz||
337 rzx!=tf.rzx||rzy!=tf.rzy||rzz!=tf.rzz) ? true : false;
341G4double G4AffineTransform::operator [] (const G4int n) const
394G4bool G4AffineTransform::IsRotated() const
396 return (rxx==1.0 && ryy==1.0 && rzz==1.0) ? false : true;
400G4bool G4AffineTransform::IsTranslated() const
402 return (tx || ty || tz) ? true:false;
405inline G4RotationMatrix G4AffineTransform::NetRotation() const
407 return G4Rep3x3(rxx,rxy,rxz,
412inline G4RotationMatrix G4AffineTransform::InverseNetRotation() const
414 return G4Rep3x3(rxx,ryx,rzx,
420G4ThreeVector G4AffineTransform::NetTranslation() const
422 return G4ThreeVector(tx,ty,tz);
426G4ThreeVector G4AffineTransform::InverseNetTranslation() const
428 G4double ttx = -tx, tty = -ty, ttz = -tz;
429 G4double invtx = ttx*rxx + tty*rxy + ttz*rxz;
430 G4double invty = ttx*ryx + tty*ryy + ttz*ryz;
431 G4double invtz = ttx*rzx + tty*rzy + ttz*rzz;
432 return G4ThreeVector(invtx,invty,invtz);
436void G4AffineTransform::SetNetRotation(const G4RotationMatrix& rot)
450void G4AffineTransform::SetNetTranslation(const G4ThreeVector& tlate)
458G4AffineTransform::operator G4Transform3D () const
460 return G4Transform3D(NetRotation().inverse(),NetTranslation());
464std::ostream& operator << (std::ostream& os, const G4AffineTransform& transf)
466 std::streamsize oldPrec = os.precision(6);
467 G4double DeviationTolerance = 1.0e-05;
469 G4double diagDeviation = 0.0;
470 G4double offdDeviationUL = 0.0;
471 G4double offdDeviationDR = 0.0;
472 G4double offdDeviation = 0.0;
474 os << " Transformation: " << G4endl;
476 // double a = std::max ( 1, 2, 3 ) ;
478 G4bool UnitTr = ! transf.IsRotated();
479 diagDeviation = std::max( std::fabs( transf[0] - 1.0 ) , // |rxx - 1|
480 std::fabs( transf[5] - 1.0 ) ); // |ryy - 1|
481 diagDeviation = std::max( diagDeviation,
482 std::fabs( transf[10] - 1.0 ) ); // |rzz - 1|
484 offdDeviationUL = std::max( std::fabs( transf[1] ) , // |rxy|
485 std::fabs( transf[2] ) ); // |rxz|
486 offdDeviationUL = std::max( offdDeviationUL,
487 std::fabs( transf[4] ) ); // |ryx|
489 offdDeviationDR = std::max( std::fabs( transf[6] ) , // |ryz|
490 std::fabs( transf[8] ) ); // |rzx|
491 offdDeviationDR = std::max( offdDeviationDR,
492 std::fabs( transf[9] ) ); // |rzy|
493 offdDeviation = std::max( offdDeviationUL, offdDeviationDR );
495 if( UnitTr || std::max(diagDeviation, offdDeviation) < DeviationTolerance )
497 os << " UNIT Rotation " << G4endl;
502 << transf[0] << " " << transf[1] << " " << transf[2] << G4endl
504 << transf[4] << " " << transf[5] << " " << transf[6] << G4endl
506 << transf[8] << " " << transf[9] << " " << transf[10] << G4endl;
509 os << "tr/x,y,z: " << transf[12] << " " << transf[13] << " " << transf[14]
512 os.precision(oldPrec);