cgma
geometry.hpp
Go to the documentation of this file.
```00001 #ifndef MCNP2CAD_GEOMETRY_H
00003
00004 #include <vector>
00005 #include <iosfwd>
00006
00007 #include "dataref.hpp"
00008
00009 // under visual studio, the following macro must be defined
00010 // in order to get M_PI and friends from the cmath header.
00011 #ifdef _MSC_VER
00012 #define _USE_MATH_DEFINES
00013 #endif
00014
00015 #include <cmath>
00016
00017 class Vector3d{
00018
00019 public:
00020
00021   double v[3];
00022   Vector3d(){
00023     v[2] = v[1] = v[0] = 0;
00024   }
00025
00026   Vector3d( const double p[3] ){
00027     v[0] = p[0];
00028     v[1] = p[1];
00029     v[2] = p[2];
00030   }
00031
00032   Vector3d( double x, double y, double z ){
00033     v[0] = x;
00034     v[1] = y;
00035     v[2] = z;
00036   }
00037
00038   Vector3d( const std::vector<double>& p, int idx = 0 ){
00039     v[0] = p.at(idx+0);
00040     v[1] = p.at(idx+1);
00041     v[2] = p.at(idx+2);
00042   }
00043
00044   double length() const{
00045     return sqrt( v[0]*v[0] + v[1]*v[1] + v[2]*v[2] );
00046   }
00047
00048   Vector3d normalize() const {
00049     double length = this->length();
00050     return Vector3d( v[0]/length, v[1]/length, v[2]/length );
00051   }
00052
00053   Vector3d operator-() const {
00054     return Vector3d(-v[0], -v[1], -v[2]);
00055   }
00056
00057   Vector3d reverse() const {
00058     return -(*this);
00059   }
00060
00061   Vector3d scale( double d ) const {
00062     return Vector3d(v[0]*d, v[1]*d, v[2]*d);
00063   }
00064
00065   Vector3d operator*( double d ) const {
00066     return scale(d);
00067   }
00068
00069   Vector3d add( const Vector3d& alt ) const {
00070     return Vector3d( v[0]+alt.v[0], v[1]+alt.v[1], v[2]+alt.v[2] );
00071   }
00072
00073   Vector3d operator+( const Vector3d& alt ) const {
00075   }
00076
00077   double dot( const Vector3d& alt ) const {
00078     return v[0]*alt.v[0] + v[1]*alt.v[1] + v[2]*alt.v[2];
00079   }
00080
00081   Vector3d cross( const Vector3d& alt ) const {
00082     Vector3d c;
00083     c.v[0] = v[1]*alt.v[2] - v[2]*alt.v[1];
00084     c.v[1] = v[2]*alt.v[0] - v[0]*alt.v[2];
00085     c.v[2] = v[0]*alt.v[1] - v[1]*alt.v[0];
00086     return c;
00087   }
00088
00089   Vector3d rotate_about( const Vector3d& v_p, double theta_p, bool degrees = true ) const {
00090     Vector3d v = v_p.normalize();
00091     double theta = theta_p;
00092     if( degrees ) { theta *= M_PI / 180.0; }
00093     double cos_t = cos( theta );
00094     Vector3d ret = scale(cos_t) + v.scale((1.0-cos_t)*(dot(v))) + cross(v).scale( sin(theta) );
00095     return ret;
00096   }
00097
00099   Vector3d projection( const Vector3d& v ) const {
00100     return scale( this->dot(v) / this->dot(*this) );
00101   }
00102
00103 };
00104
00105 std::ostream& operator<<(std::ostream& str, const Vector3d& v );
00106
00107
00108 // determinant of 3x3 matrix (C-style matrix ordering)
00109 double matrix_det( double mat[9] );
00110
00111 class Transform{
00112
00113 public:
00114   enum mat_format{ C_STYLE, FORTRAN_STYLE };
00115
00116 protected:
00117   Vector3d translation;
00118   bool has_rot;
00119   double theta; Vector3d axis;
00120   bool invert;
00121
00122   void set_rots_from_matrix( double raw_matrix[9], enum mat_format  );
00123
00124 public:
00125   Transform():translation(),has_rot(false),invert(false){}
00126   Transform( const Vector3d& v ):translation(v),has_rot(false),invert(false){}
00127   Transform( const std::vector< double >& inputs, bool degree_format_p = false, enum mat_format = FORTRAN_STYLE );
00128   Transform( double rot[9], const Vector3d& trans,  enum mat_format = C_STYLE );
00129
00130   const Vector3d& getTranslation() const { return translation; }
00131   bool hasRot() const{ return has_rot; }
00132   bool hasInversion() const{ return invert; }
00133   double getTheta() const { return theta; }
00134   const Vector3d& getAxis() const { return axis; }
00135   void print( std::ostream& str ) const;
00136
00137   Transform reverse() const;
00138
00139 };
00140
00141 std::ostream& operator<<(std::ostream& str, const Transform& t );
00142
00144 class FillNode {
00145
00146 protected:
00147   int universe;
00148   DataRef<Transform>* tr;
00149
00150 public:
00151   FillNode():
00152     universe(0), tr(new NullRef<Transform>())
00153   {}
00154
00155   FillNode( int universe_p ):
00156     universe(universe_p), tr(new NullRef<Transform>())
00157   {}
00158
00159   FillNode( int universe_p, DataRef<Transform>* tr_p ):
00160     universe(universe_p), tr(tr_p)
00161   {}
00162
00163   FillNode( const FillNode& node_p ):
00164     universe(node_p.universe), tr(node_p.tr->clone())
00165   {}
00166
00167   FillNode& operator=( const FillNode& node_p ){
00168     if( this != &node_p ){
00169       universe = node_p.universe;
00170       tr = node_p.tr->clone();
00171     }
00172     return *this;
00173   }
00174
00175   ~FillNode(){
00176     delete tr;
00177   }
00178
00179   int getFillingUniverse() const { return universe; }
00180
00181   bool hasTransform() const{ return tr->hasData();}
00182   const Transform& getTransform() const { return tr->getData(); }
00183
00184   void setTransform( DataRef<Transform>* tr_p ){
00185     delete tr;
00186     tr = tr_p;
00187   }
00188 };
00189
00190 typedef std::pair<int,int> irange;
00191
00193 class Fill{
00194
00195 protected:
00196   std::vector<FillNode> nodes;
00197   bool has_grid;
00198   irange xrange, yrange, zrange;
00199
00200
00201  size_t indicesToSerialIndex( int x, int y, int z ) const ;
00202
00203 public:
00204   Fill():
00205     nodes(1,FillNode()),has_grid(false)
00206   {}
00207
00208   Fill( const FillNode& origin_p ):
00209     nodes(1,origin_p), has_grid(false)
00210   {}
00211
00212   Fill( irange x, irange y, irange z, std::vector<FillNode> nodes_p ):
00213     nodes(nodes_p), has_grid(true), xrange(x), yrange(y), zrange(z)
00214   {}
00215
00216   const FillNode& getOriginNode() const;
00217   const FillNode& getNode( int x, int y, int z ) const;
00218
00219   friend class Lattice;
00220 };
00221
00222
00223 class Lattice{
00224
00225 protected:
00226   int num_finite_dims;
00227   Vector3d v1, v2, v3;
00228
00229   DataRef<Fill> *fill;
00230
00231 public:
00232   Lattice() : fill(new NullRef<Fill>()){}
00233   Lattice( int dims, const Vector3d& v1_p, const Vector3d& v2_p, const Vector3d& v3_p, const FillNode& singleton_fill );
00234   Lattice( int dims, const Vector3d& v1_p, const Vector3d& v2_p, const Vector3d& v3_p, const Fill& full_fill );
00235
00236   Lattice( const Lattice& l );
00237   Lattice& operator=( const Lattice& l );
00238   ~Lattice(){ delete fill; }
00239
00240
00241   int numFiniteDirections() const { return num_finite_dims; }
00242   Transform getTxForNode( int x, int y, int z ) const ;
00243
00244   const FillNode& getFillForNode( int x, int y, int z ) const ;
00245
00246   bool isFixedSize() const { return fill->getData().has_grid; }
00247   irange getXRange() const { return fill->getData().xrange; }
00248   irange getYRange() const { return fill->getData().yrange; }
00249   irange getZRange() const { return fill->getData().zrange; }
00250
00251
00252 };
00253
00254
00255 #endif /* MCNP2CAD_GEOMETRY_H */
```