ESyS-Particle
4.0.1
|
00001 00002 // // 00003 // Copyright (c) 2003-2011 by The University of Queensland // 00004 // Earth Systems Science Computational Centre (ESSCC) // 00005 // http://www.uq.edu.au/esscc // 00006 // // 00007 // Primary Business: Brisbane, Queensland, Australia // 00008 // Licensed under the Open Software License version 3.0 // 00009 // http://www.opensource.org/licenses/osl-3.0.php // 00010 // // 00012 00013 00014 #ifndef __VEC3_HPP 00015 #define __VEC3_HPP 00016 00017 #include "Foundation/Matrix3.h" 00018 00019 //the error... 00020 VEC3_INLINE VecErr::VecErr(const string& m):MError(m) 00021 { 00022 message.insert(0,"Vec3 "); 00023 } 00024 00025 // constructors 00026 VEC3_INLINE Vec3::Vec3() 00027 { 00028 data[0]=0; 00029 data[1]=0; 00030 data[2]=0; 00031 } 00032 00033 VEC3_INLINE Vec3::Vec3(double s) 00034 { 00035 data[0]=s; 00036 data[1]=s; 00037 data[2]=s; 00038 } 00039 00040 VEC3_INLINE Vec3::Vec3(double a,double b,double c) 00041 { 00042 data[0]=a; 00043 data[1]=b; 00044 data[2]=c; 00045 } 00046 00047 VEC3_INLINE Vec3::Vec3(const Vec3& rhs) 00048 { 00049 data[0]=rhs.data[0]; 00050 data[1]=rhs.data[1]; 00051 data[2]=rhs.data[2]; 00052 } 00053 00054 // operators 00055 00056 VEC3_INLINE Vec3& Vec3::operator=(const Vec3& rhs) 00057 { 00058 data[0]=rhs.data[0]; 00059 data[1]=rhs.data[1]; 00060 data[2]=rhs.data[2]; 00061 return *this; 00062 } 00063 00064 VEC3_INLINE Vec3& Vec3::operator=(double s) 00065 { 00066 data[0]=s; 00067 data[1]=s; 00068 data[2]=s; 00069 return *this; 00070 } 00071 00072 VEC3_INLINE Vec3& Vec3::operator-=(const Vec3& rhs) 00073 { 00074 data[0]-=rhs.data[0]; 00075 data[1]-=rhs.data[1]; 00076 data[2]-=rhs.data[2]; 00077 return *this; 00078 } 00079 00080 VEC3_INLINE Vec3& Vec3::operator+=(const Vec3& rhs) 00081 { 00082 data[0]+=rhs.data[0]; 00083 data[1]+=rhs.data[1]; 00084 data[2]+=rhs.data[2]; 00085 return *this; 00086 } 00087 00088 VEC3_INLINE Vec3 Vec3::operator+(const Vec3& rhs) const 00089 { 00090 return Vec3(data[0]+rhs.data[0], data[1]+rhs.data[1], data[2]+rhs.data[2]); 00091 } 00092 00093 VEC3_INLINE Vec3 Vec3::operator-(const Vec3& rhs) const 00094 { 00095 return Vec3(data[0]-rhs.data[0], data[1]-rhs.data[1], data[2]-rhs.data[2]); 00096 } 00097 00098 VEC3_INLINE Vec3 Vec3::operator-() const 00099 { 00100 return Vec3( -data[0],-data[1],-data[2] ); 00101 } 00102 00103 VEC3_INLINE Vec3 Vec3::operator*(const Matrix3 &m) const 00104 { 00105 const double x = m(0,0)*data[0] + m(1,0)*data[1] + m(2,0)*data[2]; 00106 const double y = m(0,1)*data[0] + m(1,1)*data[1] + m(2,1)*data[2]; 00107 const double z = m(0,2)*data[0] + m(1,2)*data[1] + m(2,2)*data[2]; 00108 00109 return Vec3(x,y,z); 00110 } 00111 00112 VEC3_INLINE double Vec3::operator*(const Vec3& rhs) const 00113 { 00114 return data[0]*rhs.data[0]+data[1]*rhs.data[1]+data[2]*rhs.data[2]; 00115 } 00116 00117 VEC3_INLINE Vec3 Vec3::operator*(double s) const 00118 { 00119 return Vec3(data[0]*s,data[1]*s,data[2]*s) ; 00120 } 00121 00122 VEC3_INLINE Vec3& Vec3::operator*=(double rhs) 00123 { 00124 data[0]*=rhs; 00125 data[1]*=rhs; 00126 data[2]*=rhs; 00127 return *this; 00128 } 00129 00130 VEC3_INLINE Vec3& Vec3::operator/=(double c) 00131 { 00132 data[0] /= c; 00133 data[1] /= c; 00134 data[2] /= c; 00135 00136 return *this; 00137 } 00138 00139 VEC3_INLINE Vec3 Vec3::operator/(double s) const 00140 { 00141 return Vec3(data[0]/s,data[1]/s,data[2]/s) ; 00142 } 00143 00144 VEC3_INLINE Vec3 Vec3::operator+(double s) const 00145 { 00146 return Vec3(data[0]+s, data[1]+s, data[2]+s) ; 00147 } 00148 00149 VEC3_INLINE Vec3 Vec3::operator-(double s) const 00150 { 00151 return Vec3(data[0]-s, data[1]-s, data[2]-s); 00152 } 00153 00154 VEC3_INLINE Vec3 Vec3::rotate(const Vec3 &axis, const Vec3 &axisPt) const 00155 { 00156 const double phi = axis.norm(); 00157 if (phi > 0.0) 00158 { 00159 const Vec3 r = *this - axisPt; 00160 const Vec3 n = axis/phi; 00161 const double cosPhi = cos(phi); 00162 const Vec3 rotatedR = 00163 r*cosPhi + n*((dot(n, r))*(1-cosPhi)) + cross(r, n)*sin(phi); 00164 return rotatedR + axisPt; 00165 } 00166 return *this; 00167 } 00168 00169 VEC3_INLINE Vec3 &Vec3::operator+=(double s) 00170 { 00171 data[0] += s; 00172 data[1] += s; 00173 data[2] += s; 00174 return *this; 00175 } 00176 00177 VEC3_INLINE Vec3 &Vec3::operator-=(double s) 00178 { 00179 data[0] -= s; 00180 data[1] -= s; 00181 data[2] -= s; 00182 return *this; 00183 } 00184 00185 // vector product 00186 // 9 Flops ( 6 mult, 3 sub ) 00187 VEC3_INLINE Vec3 cross(const Vec3& lhs,const Vec3& rhs) 00188 { 00189 return Vec3(lhs.data[1]*rhs.data[2]-lhs.data[2]*rhs.data[1], 00190 lhs.data[2]*rhs.data[0]-lhs.data[0]*rhs.data[2], 00191 lhs.data[0]*rhs.data[1]-lhs.data[1]*rhs.data[0]); 00192 } 00193 00194 // dot product 00195 00196 VEC3_INLINE double dot(const Vec3& v1, const Vec3& v2) 00197 { 00198 return v1.data[0] * v2.data[0] + 00199 v1.data[1] * v2.data[1] + 00200 v1.data[2] * v2.data[2]; 00201 } 00202 00203 VEC3_INLINE Vec3 operator*(double f,const Vec3& rhs) 00204 { 00205 return Vec3(f*rhs.data[0], f*rhs.data[1], f*rhs.data[2]); 00206 } 00207 00208 00209 // euclidian norm 00210 // 6 Flops ( 3 mult, 2 add, 1 sqrt ) 00211 VEC3_INLINE double Vec3::norm() const 00212 { 00213 return sqrt(data[0]*data[0]+data[1]*data[1]+data[2]*data[2]); 00214 } 00215 00216 // square of the euclidian norm 00217 // 5 Flops ( 3 mult, 2 add) 00218 VEC3_INLINE double Vec3::norm2() const 00219 { 00220 return data[0]*data[0]+data[1]*data[1]+data[2]*data[2]; 00221 } 00222 00223 // returns unit vector in direction of the original vector 00224 // 9 Flops ( 3 mult, 2 add, 3 div, 1 sqrt ) 00225 VEC3_INLINE Vec3 Vec3::unit() const 00226 { 00227 return (*this)/norm(); 00228 } 00229 00230 // per element min/max 00231 VEC3_INLINE Vec3 cmax(const Vec3& v1,const Vec3& v2) 00232 { 00233 Vec3 res; 00234 res.data[0]=v1.data[0]>v2.data[0] ? v1.data[0] : v2.data[0]; 00235 res.data[1]=v1.data[1]>v2.data[1] ? v1.data[1] : v2.data[1]; 00236 res.data[2]=v1.data[2]>v2.data[2] ? v1.data[2] : v2.data[2]; 00237 return res; 00238 } 00239 00240 VEC3_INLINE Vec3 cmin(const Vec3& v1,const Vec3& v2) 00241 { 00242 Vec3 res; 00243 res.data[0]=v1.data[0]<v2.data[0] ? v1.data[0] : v2.data[0]; 00244 res.data[1]=v1.data[1]<v2.data[1] ? v1.data[1] : v2.data[1]; 00245 res.data[2]=v1.data[2]<v2.data[2] ? v1.data[2] : v2.data[2]; 00246 return res; 00247 } 00248 00249 // save version, throws exception if norm()==0 00250 VEC3_INLINE Vec3 Vec3::unit_s() const 00251 { 00252 double n=norm(); 00253 if(n==0) throw VecErr("norm() of data[2]ero-vector"); 00254 Vec3 res(data[0],data[1],data[2]); 00255 return res/n; 00256 } 00257 00258 VEC3_INLINE double Vec3::max() const 00259 { 00260 double m = ( data[0]>data[1] ? data[0] : data[1] ); 00261 00262 return ( m>data[2] ? m : data[2] ); 00263 00264 } 00265 00266 VEC3_INLINE double Vec3::min() const 00267 { 00268 double m = ( data[0]<data[1] ? data[0] : data[1] ); 00269 00270 return ( m<data[2] ? m : data[2] ); 00271 } 00272 00273 VEC3_INLINE bool Vec3::operator==(const Vec3& V) const 00274 { 00275 return((data[0]==V.data[0])&&(data[1]==V.data[1])&&(data[2]==V.data[2])); 00276 } 00277 00278 VEC3_INLINE bool Vec3::operator!=(const Vec3& V) const 00279 { 00280 return((data[0]!=V.data[0])||(data[1]!=V.data[1])||(data[2]!=V.data[2])); 00281 } 00282 00283 // per component min/max 00284 00285 VEC3_INLINE Vec3 comp_max(const Vec3& V1,const Vec3& V2) 00286 { 00287 double x=(V1.X() > V2.X()) ? V1.X() : V2.X(); 00288 double y=(V1.Y() > V2.Y()) ? V1.Y() : V2.Y(); 00289 double z=(V1.Z() > V2.Z()) ? V1.Z() : V2.Z(); 00290 00291 return Vec3(x,y,z); 00292 } 00293 00294 VEC3_INLINE Vec3 comp_min(const Vec3& V1,const Vec3& V2) 00295 { 00296 double x=(V1.X() < V2.X()) ? V1.X() : V2.X(); 00297 double y=(V1.Y() < V2.Y()) ? V1.Y() : V2.Y(); 00298 double z=(V1.Z() < V2.Z()) ? V1.Z() : V2.Z(); 00299 00300 return Vec3(x,y,z); 00301 } 00302 00303 // in/output 00304 00305 VEC3_INLINE std::ostream& operator << (std::ostream& ostr,const Vec3& V) 00306 { 00307 const char delimiter = ' '; 00308 ostr 00309 << V.data[0] << delimiter 00310 << V.data[1] << delimiter 00311 << V.data[2]; 00312 00313 return ostr; 00314 } 00315 00316 VEC3_INLINE std::istream& operator >> (std::istream& istr,Vec3& V) 00317 { 00318 istr 00319 >> V.data[0] 00320 >> V.data[1] 00321 >> V.data[2]; 00322 00323 return istr; 00324 } 00325 00326 #endif // __VEC3_HPP