From 88654d76f498424ed6cba9b156b407df4db206bb Mon Sep 17 00:00:00 2001
From: marco <marco@5a81b35b-ba03-0410-adc8-b2c5c5119f08>
Date: Mon, 12 Jul 2010 07:34:21 +0000
Subject: [PATCH] move small vector/matrix functions to header file

This allows for better inlining. In addition, I've changed the data
array of the vector classes to x,y,z,w member variables. This makes
accessing them easier. No indended functionality change.

First benchmarks suggest up to 4x speed improvement for vector/mat
heavy algorithms.

git-svn-id: https://dng.biozentrum.unibas.ch/svn/openstructure/trunk@2551 5a81b35b-ba03-0410-adc8-b2c5c5119f08
---
 modules/base/doc/generic.rst                |   2 +-
 modules/geom/src/composite2.cc              |  28 +--
 modules/geom/src/composite2_op.cc           |   8 +-
 modules/geom/src/mat2.cc                    |  14 --
 modules/geom/src/mat2.hh                    |  13 +-
 modules/geom/src/mat3.cc                    | 132 --------------
 modules/geom/src/mat3.hh                    | 182 ++++++++++++++++----
 modules/geom/src/quat.cc                    |  30 ++--
 modules/geom/src/vec2.cc                    | 155 -----------------
 modules/geom/src/vec2.hh                    | 140 +++++++++++----
 modules/geom/src/vec3.cc                    | 177 -------------------
 modules/geom/src/vec3.hh                    | 152 ++++++++++++----
 modules/geom/src/vec4.cc                    | 133 --------------
 modules/geom/src/vec4.hh                    | 157 +++++++++++++----
 modules/geom/src/vecmat2_op.cc              |  88 +---------
 modules/geom/src/vecmat2_op.hh              | 101 +++++++++--
 modules/geom/src/vecmat3_op.cc              | 117 -------------
 modules/geom/src/vecmat3_op.hh              | 132 +++++++++++---
 modules/geom/src/vecmat4_op.cc              | 121 +------------
 modules/geom/src/vecmat4_op.hh              | 113 ++++++++++--
 modules/gfx/src/impl/simple_renderer.cc     |   4 +-
 modules/gui/src/scene_selection.cc          |   6 +-
 modules/io/src/mol/dcd_io.cc                |   6 +-
 modules/io/src/mol/entity_io_crd_handler.cc |   6 +-
 24 files changed, 837 insertions(+), 1180 deletions(-)

diff --git a/modules/base/doc/generic.rst b/modules/base/doc/generic.rst
index e24ca06b9..919497f39 100644
--- a/modules/base/doc/generic.rst
+++ b/modules/base/doc/generic.rst
@@ -45,7 +45,7 @@ To store a float value with the key 'myfloatprop' in all atoms of an entity:
   
   import math
   for atom in entity.GetAtomList(): 
-    val=5*math.sin(0.4*atom.GetPos().GetX())
+    val=5*math.sin(0.4*atom.GetPos().x)
     atom.SetFloatProp("myfloatprop", val)
   
 If a GenericProp at a given level (i.e. atom, bond, residue, chain or entity) 
diff --git a/modules/geom/src/composite2.cc b/modules/geom/src/composite2.cc
index b8300fcfb..f6d7f6857 100644
--- a/modules/geom/src/composite2.cc
+++ b/modules/geom/src/composite2.cc
@@ -88,16 +88,16 @@ Vec2 Ellipse2::calc_(Real t) const
 Rectangle2 Ellipse2::GetBoundingBox() const
 {
   if(gamma_==0.0 || gamma_==M_PI){
-    return Rectangle2(Vec2(origin_.GetX()-a_,origin_.GetY()-b_),Vec2(origin_.GetX()+a_/2,origin_.GetY()+b_));
+    return Rectangle2(Vec2(origin_.x-a_,origin_.y-b_),Vec2(origin_.x+a_/2,origin_.y+b_));
   }
   if(gamma_==M_PI*0.5 || gamma_==3*M_PI*0.5){
-    return Rectangle2(Vec2(origin_.GetX()-b_,origin_.GetY()-a_),Vec2(origin_.GetX()+b_/2,origin_.GetY()+a_));
+    return Rectangle2(Vec2(origin_.x-b_,origin_.y-a_),Vec2(origin_.x+b_/2,origin_.y+a_));
   }
   Real t=atan(-b_/a_*tan(gamma_));
   Vec2 vx=calc_(t);
   t=atan(b_/a_/tan(gamma_));
   Vec2 vy=calc_(t);
-  return Rectangle2(Vec2(origin_.GetX()-vx.GetX(),origin_.GetY()-vy.GetY()),Vec2(origin_.GetX()+vx.GetX(),origin_.GetY()+vy.GetY()));
+  return Rectangle2(Vec2(origin_.x-vx.x,origin_.y-vy.y),Vec2(origin_.x+vx.x,origin_.y+vy.y));
 }
 
 Rectangle2::Rectangle2():
@@ -112,11 +112,11 @@ bottomright_(bottomright)
 }
 Real Rectangle2::GetWidth() const
 {
-    return bottomright_.GetX()-topleft_.GetX();
+    return bottomright_.x-topleft_.x;
 }
 Real Rectangle2::GetHeight() const
 {
-    return bottomright_.GetY()-topleft_.GetY();
+    return bottomright_.y-topleft_.y;
 }
 Real Rectangle2::GetArea() const
 {
@@ -207,7 +207,7 @@ Real Polygon2::GetArea() const
     unsigned int ip1=(i+1)%nodecount;
     Vec2 n1=GetNode(i);
     Vec2 n2=GetNode(ip1);
-    area+=n1.GetX()*n2.GetY()-n2.GetX()*n1.GetY();
+    area+=n1.x*n2.y-n2.x*n1.y;
   }
   return area/2.0;
 }
@@ -220,8 +220,8 @@ Vec2 Polygon2::GetCentroid() const
     unsigned int ip1=(i+1)%nodecount;
     Vec2 n1=GetNode(i);
     Vec2 n2=GetNode(ip1);
-    cx+=(n1.GetX()+n2.GetX())*(n1.GetX()*n2.GetY()-n2.GetX()*n1.GetY());
-    cy+=(n1.GetY()+n2.GetY())*(n1.GetX()*n2.GetY()-n2.GetX()*n1.GetY());
+    cx+=(n1.x+n2.x)*(n1.x*n2.y-n2.x*n1.y);
+    cy+=(n1.y+n2.y)*(n1.x*n2.y-n2.x*n1.y);
   }
   Real area=GetArea();
   return Vec2(cx/(6.0*area),cy/(6.0*area));
@@ -255,17 +255,17 @@ Rectangle2 Polygon2::GetBoundingBox() const
   {
     return Rectangle2(Vec2(0.0,0.0),Vec2(0.0,0.0));
   }
-  Real x_min=operator[](0).GetX();
+  Real x_min=operator[](0).x;
   Real x_max=x_min;
-  Real y_min=operator[](0).GetY();
+  Real y_min=operator[](0).y;
   Real y_max=y_min;
   unsigned int nodecount=GetNodeCount();
   for(unsigned int i=0;i<nodecount;++i){
     Vec2 n=operator[](i);
-    x_min=std::min<Real>(n.GetX(),x_min);
-    x_max=std::max<Real>(n.GetX(),x_max);
-    y_min=std::min<Real>(n.GetY(),y_min);
-    y_max=std::max<Real>(n.GetY(),y_max);
+    x_min=std::min<Real>(n.x,x_min);
+    x_max=std::max<Real>(n.x,x_max);
+    y_min=std::min<Real>(n.y,y_min);
+    y_max=std::max<Real>(n.y,y_max);
   }
   return Rectangle2(Vec2(x_min,y_min),Vec2(x_max,y_max));
 }
diff --git a/modules/geom/src/composite2_op.cc b/modules/geom/src/composite2_op.cc
index cc660f821..5044771a8 100644
--- a/modules/geom/src/composite2_op.cc
+++ b/modules/geom/src/composite2_op.cc
@@ -88,11 +88,11 @@ bool IsInPolygon(const Polygon2& p, const Vec2& v)
   }
   bool c = false;
   Vec2 old=p.GetNode(p.GetNodeCount()-1);
-  bool oldcomp= v.GetY()<old.GetY();
+  bool oldcomp= v.y<old.y;
   for (Polygon2::const_iterator it=p.begin();it!=p.end(); ++it) {
-    bool comp= v.GetY() < (*it).GetY();
+    bool comp= v.y < (*it).y;
     if( comp ^ oldcomp ) {
-      if( (v.GetX() < (old.GetX() - (*it).GetX()) * (v.GetY() - (*it).GetY()) / (old.GetY() - (*it).GetY()) + (*it).GetX()) ){
+      if( (v.x < (old.x - (*it).x) * (v.y - (*it).y) / (old.y - (*it).y) + (*it).x) ){
         c = !c;
       }
     }
@@ -123,7 +123,7 @@ DLLEXPORT bool IsInRectangle(const Rectangle2& r, const Vec2& v)
 {
   Vec2 start=r.GetStart();
   Vec2 end=r.GetEnd();
-  return v.GetX()>=start.GetX() && v.GetX()<=end.GetX() &&  v.GetY()>=start.GetY() && v.GetY()<=end.GetY();
+  return v.x>=start.x && v.x<=end.x &&  v.y>=start.y && v.y<=end.y;
 }
 DLLEXPORT bool IsInCircle(const Circle2& c, const Vec2& v)
 {
diff --git a/modules/geom/src/mat2.cc b/modules/geom/src/mat2.cc
index ac34069df..06e196d8c 100644
--- a/modules/geom/src/mat2.cc
+++ b/modules/geom/src/mat2.cc
@@ -65,20 +65,6 @@ bool Mat2::operator==(const Mat2& rhs) const
     data_[1][1] == rhs.data_[1][1];
 }
 
-
-Real& Mat2::operator()(std::size_t r, std::size_t c)
-{
-  if(r>1 || c>1) throw OutOfRangeException();
-  return data_[r][c];
-}
-
-const Real& Mat2::operator()(std::size_t r, std::size_t c) const
-{
-  if(r>1 || c>1) throw OutOfRangeException();
-  return data_[r][c];
-}
-
-
 Mat2& Mat2::operator+=(const Mat2& rhs)
 {
   data_[0][0]+=rhs(0,0);
diff --git a/modules/geom/src/mat2.hh b/modules/geom/src/mat2.hh
index 9231cb710..bc53416f2 100644
--- a/modules/geom/src/mat2.hh
+++ b/modules/geom/src/mat2.hh
@@ -21,6 +21,7 @@
 
 #include <cstddef> // for size_t
 #include <ostream>
+#include <cassert>
 
 #include <boost/operators.hpp>
 
@@ -63,9 +64,17 @@ public:
   bool operator==(const Mat2& rhs) const;
 
   //! element access
-  Real& operator()(std::size_t r, std::size_t c);
+  Real& operator()(std::size_t r, std::size_t c)
+  {
+    assert(r<=1 && c<=1);
+    return data_[r][c];
+  }
   //! const element access
-  const Real& operator()(std::size_t r, std::size_t c) const;
+  const Real& operator()(std::size_t r, std::size_t c) const
+  {
+    assert(r<=1 && c<=1);
+    return data_[r][c];
+  }
 
   //! addable op
   Mat2& operator+=(const Mat2& rhs);
diff --git a/modules/geom/src/mat3.cc b/modules/geom/src/mat3.cc
index f9492d3d4..a202f20a3 100644
--- a/modules/geom/src/mat3.cc
+++ b/modules/geom/src/mat3.cc
@@ -24,132 +24,6 @@
 
 namespace geom {
 
-Mat3::Mat3()
-{
-  this->set(1.0,0.0,0.0, 0.0,1.0,0.0, 0.0,0.0,1.0);
-}
-
-Mat3::Mat3(Real i00, Real i01, Real i02,Real i10, Real i11, Real i12,Real i20, Real i21, Real i22)
-{
-  this->set(i00,i01,i02,i10,i11,i12,i20,i21,i22);
-}
-
-Mat3::Mat3(const Mat3& m)
-{
-  this->set(m(0,0),m(0,1),m(0,2),m(1,0),m(1,1),m(1,2),m(2,0),m(2,1),m(2,2));
-}
-
-Mat3::Mat3(const Mat2& m)
-{
-  this->set(m(0,0),m(0,1),0.0, m(1,0),m(1,1),0.0, 0.0,0.0,1.0);
-}
-
-Mat3::Mat3(const Real arr[9])
-{
-  this->set(arr[0],arr[1],arr[2],arr[3],arr[4],arr[5],arr[6],arr[7],arr[8]);
-}
-
-Mat3& Mat3::operator=(const Mat3& m)
-{
-  if(&m!=this) {
-    this->set(m(0,0),m(0,1),m(0,2),m(1,0),m(1,1),m(1,2),m(2,0),m(2,1),m(2,2));
-  }
-  return *this;
-}
-
-Mat3 Mat3::Identity() 
-{
-  static Mat3 i(1.0,0.0,0.0,
-                0.0,1.0,0.0,
-                0.0,0.0,1.0);
-  return i;
-}
-
-bool Mat3::operator==(const Mat3& rhs) const
-{
-  return data_[0][0] ==  rhs.data_[0][0] &&
-    data_[1][0] == rhs.data_[1][0] &&
-    data_[2][0] == rhs.data_[2][0] &&
-    data_[0][1] == rhs.data_[0][1] &&
-    data_[1][1] ==  rhs.data_[1][1] &&
-    data_[2][1] ==  rhs.data_[2][1] &&
-    data_[0][2] ==  rhs.data_[0][2] &&
-    data_[1][2] ==  rhs.data_[1][2] &&
-    data_[2][2] ==  rhs.data_[2][2];
-}
-
-Real& Mat3::operator()(std::size_t r, std::size_t c)
-{
-  if(r>2 || c>2) throw OutOfRangeException();
-  return data_[r][c];
-}
-
-const Real& Mat3::operator()(std::size_t r, std::size_t c) const
-{
-  if(r>2 || c>2) throw OutOfRangeException();
-  return data_[r][c];
-}
-
-
-Mat3& Mat3::operator+=(const Mat3& rhs)
-{
-  data_[0][0]+=rhs(0,0);
-  data_[0][1]+=rhs(0,1);
-  data_[0][2]+=rhs(0,2);
-  data_[1][0]+=rhs(1,0);
-  data_[1][1]+=rhs(1,1);
-  data_[1][2]+=rhs(1,2);
-  data_[2][0]+=rhs(2,0);
-  data_[2][1]+=rhs(2,1);
-  data_[2][2]+=rhs(2,2);
-  return *this;
-}
-Mat3& Mat3::operator-=(const Mat3& rhs)
-{
-  data_[0][0]-=rhs(0,0);
-  data_[0][1]-=rhs(0,1);
-  data_[0][2]-=rhs(0,2);
-  data_[1][0]-=rhs(1,0);
-  data_[1][1]-=rhs(1,1);
-  data_[1][2]-=rhs(1,2);
-  data_[2][0]-=rhs(2,0);
-  data_[2][1]-=rhs(2,1);
-  data_[2][2]-=rhs(2,2);
-  return *this;
-}
-Mat3& Mat3::operator*=(const Real d)
-{
-  data_[0][0]*=d;
-  data_[0][1]*=d;
-  data_[0][2]*=d;
-  data_[1][0]*=d;
-  data_[1][1]*=d;
-  data_[1][2]*=d;
-  data_[2][0]*=d;
-  data_[2][1]*=d;
-  data_[2][2]*=d;
-  return *this;
-}
-Mat3& Mat3::operator/=(const Real d)
-{
-  data_[0][0]/=d;
-  data_[0][1]/=d;
-  data_[0][2]/=d;
-  data_[1][0]/=d;
-  data_[1][1]/=d;
-  data_[1][2]/=d;
-  data_[2][0]/=d;
-  data_[2][1]/=d;
-  data_[2][2]/=d;
-  return *this;
-}
-
-Mat3& Mat3::operator*=(const Mat3& m)
-{
-  (*this)=operator*(*this,m);
-  return *this;
-}
-
 geom::Vec3 Mat3::GetCol(int index) const
 {
   return geom::Vec3(data_[0][index], data_[1][index], data_[2][index]);  
@@ -163,12 +37,6 @@ geom::Vec3 Mat3::GetRow(int index) const
 
 ////////////////////////////////////////////
 // private member functions
-void Mat3::set(Real i00, Real i01, Real i02,Real i10, Real i11, Real i12,Real i20, Real i21, Real i22)
-{
-  data_[0][0]=i00; data_[0][1]=i01; data_[0][2]=i02;
-  data_[1][0]=i10; data_[1][1]=i11; data_[1][2]=i12;
-  data_[2][0]=i20; data_[2][1]=i21; data_[2][2]=i22;
-}
 
 std::ostream& operator<<(std::ostream& os, const Mat3& m)
 {
diff --git a/modules/geom/src/mat3.hh b/modules/geom/src/mat3.hh
index f7d90c3f3..158158633 100644
--- a/modules/geom/src/mat3.hh
+++ b/modules/geom/src/mat3.hh
@@ -21,10 +21,12 @@
 
 #include <cstddef> // for size_t
 #include <ostream>
+#include <cassert>
 
 #include <boost/operators.hpp>
 
 #include <ost/geom/module_config.hh>
+#include <ost/geom/mat2.hh>
 
 namespace geom {
 
@@ -36,12 +38,12 @@ class DLLEXPORT_OST_GEOM Mat3:
     private boost::additive1<Mat3>,
     private boost::multiplicative2<Mat3, Real>
 {
-public:
-  static Mat3 Identity();
-  
-  //! Default initialization, identity matrix
-  Mat3();
-
+public:  
+  //! Default initialization, identity matrix  
+  Mat3()
+  {
+    this->set(1.0,0.0,0.0, 0.0,1.0,0.0, 0.0,0.0,1.0);
+  }
   //! In with 9 values in row-major order
   /*!
     row-major order means that the matrix
@@ -52,37 +54,138 @@ public:
 
     is initialized with (a,b,c, d,e,f, g,h,i)
   */
-  Mat3(Real i00, Real i01, Real i02,Real i10, Real i11, Real i12,Real i20, Real i21, Real i22);
-
-  //! Copy ctor
-  Mat3(const Mat3& m);
-
-  // initialization with 2d matrix
-  explicit Mat3(const Mat2& m);
-
-  //! initialization from array
-  explicit Mat3(const Real[9]);
-
-  //! assignement op
-  Mat3& operator=(const Mat3& m);
-
-  //! comparable
-  bool operator==(const Mat3& rhs) const;
+  Mat3(Real i00, Real i01, Real i02,
+       Real i10, Real i11, Real i12,
+       Real i20, Real i21, Real i22)
+  {
+    this->set(i00,i01,i02,i10,i11,i12,i20,i21,i22);
+  }
+
+  Mat3(const Mat3& m)
+  {
+    this->set(m(0,0),m(0,1),m(0,2),m(1,0),m(1,1),m(1,2),m(2,0),m(2,1),m(2,2));
+  }
+
+  Mat3(const Mat2& m)
+  {
+    this->set(m(0, 0), m(0, 1), Real(0.0),
+              m(1, 0), m(1, 1), Real(0.0),
+              Real(0.0), Real(0.0), Real(1.0));
+  }
+  
+  explicit Mat3(const Real arr[9])
+  {
+    this->set(arr[0],arr[1],arr[2],arr[3],arr[4],arr[5],arr[6],arr[7],arr[8]);
+  }  
 
   //! element access
-  Real& operator()(std::size_t r, std::size_t c);
+  Real& operator()(std::size_t r, std::size_t c)
+  {
+    assert(r<=2 && c<=2);
+    return data_[r][c];
+  }
   //! const element access
-  const Real& operator()(std::size_t r, std::size_t c) const;
-
-  //! addable op
-  Mat3& operator+=(const Mat3& rhs);
-  //! subtractable op
-  Mat3& operator-=(const Mat3& rhs);
-
-  Mat3& operator*=(const Real d);
-  Mat3& operator/=(const Real d);
-
-  Mat3& operator*=(const Mat3& m);
+  const Real& operator()(std::size_t r, std::size_t c) const
+  {
+    assert(r<=2 && c<=2);
+    return data_[r][c];
+  }
+
+  Mat3& operator=(const Mat3& m)
+  {
+    if(&m!=this) {
+      this->set(m(0,0),m(0,1),m(0,2),m(1,0),m(1,1),m(1,2),m(2,0),m(2,1),m(2,2));
+    }
+    return *this;
+  }
+
+  static Mat3 Identity() 
+  {
+    static Mat3 i(1.0,0.0,0.0,
+                  0.0,1.0,0.0,
+                  0.0,0.0,1.0);
+    return i;
+  }
+
+  bool operator==(const Mat3& rhs) const
+  {
+    return data_[0][0] ==  rhs.data_[0][0] &&
+      data_[1][0] == rhs.data_[1][0] &&
+      data_[2][0] == rhs.data_[2][0] &&
+      data_[0][1] == rhs.data_[0][1] &&
+      data_[1][1] ==  rhs.data_[1][1] &&
+      data_[2][1] ==  rhs.data_[2][1] &&
+      data_[0][2] ==  rhs.data_[0][2] &&
+      data_[1][2] ==  rhs.data_[1][2] &&
+      data_[2][2] ==  rhs.data_[2][2];
+  }
+
+  Mat3& operator+=(const Mat3& rhs)
+  {
+    data_[0][0]+=rhs(0,0);
+    data_[0][1]+=rhs(0,1);
+    data_[0][2]+=rhs(0,2);
+    data_[1][0]+=rhs(1,0);
+    data_[1][1]+=rhs(1,1);
+    data_[1][2]+=rhs(1,2);
+    data_[2][0]+=rhs(2,0);
+    data_[2][1]+=rhs(2,1);
+    data_[2][2]+=rhs(2,2);
+    return *this;
+  }
+  Mat3& operator-=(const Mat3& rhs)
+  {
+    data_[0][0]-=rhs(0,0);
+    data_[0][1]-=rhs(0,1);
+    data_[0][2]-=rhs(0,2);
+    data_[1][0]-=rhs(1,0);
+    data_[1][1]-=rhs(1,1);
+    data_[1][2]-=rhs(1,2);
+    data_[2][0]-=rhs(2,0);
+    data_[2][1]-=rhs(2,1);
+    data_[2][2]-=rhs(2,2);
+    return *this;
+  }
+  Mat3& operator*=(const Real d)
+  {
+    data_[0][0]*=d;
+    data_[0][1]*=d;
+    data_[0][2]*=d;
+    data_[1][0]*=d;
+    data_[1][1]*=d;
+    data_[1][2]*=d;
+    data_[2][0]*=d;
+    data_[2][1]*=d;
+    data_[2][2]*=d;
+    return *this;
+  }
+  Mat3& operator/=(const Real d)
+  {
+    data_[0][0]/=d;
+    data_[0][1]/=d;
+    data_[0][2]/=d;
+    data_[1][0]/=d;
+    data_[1][1]/=d;
+    data_[1][2]/=d;
+    data_[2][0]/=d;
+    data_[2][1]/=d;
+    data_[2][2]/=d;
+    return *this;
+  }
+
+  Mat3& operator*=(const Mat3& m)
+  {
+    (*this)=Mat3((*this)(0,0)*m(0,0)+(*this)(0,1)*m(1,0)+(*this)(0,2)*m(2,0),
+                 (*this)(0,0)*m(0,1)+(*this)(0,1)*m(1,1)+(*this)(0,2)*m(2,1),
+                 (*this)(0,0)*m(0,2)+(*this)(0,1)*m(1,2)+(*this)(0,2)*m(2,2),
+                 (*this)(1,0)*m(0,0)+(*this)(1,1)*m(1,0)+(*this)(1,2)*m(2,0),
+                 (*this)(1,0)*m(0,1)+(*this)(1,1)*m(1,1)+(*this)(1,2)*m(2,1),
+                 (*this)(1,0)*m(0,2)+(*this)(1,1)*m(1,2)+(*this)(1,2)*m(2,2),                                   
+                 (*this)(2,0)*m(0,0)+(*this)(2,1)*m(1,0)+(*this)(2,2)*m(2,0),
+                 (*this)(2,0)*m(0,1)+(*this)(2,1)*m(1,1)+(*this)(2,2)*m(2,1),
+                 (*this)(2,0)*m(0,2)+(*this)(2,1)*m(1,2)+(*this)(2,2)*m(2,2));
+    return *this;
+  }
 
   Real* Data() {return &data_[0][0];}
   const Real* Data() const {return &data_[0][0];}
@@ -92,9 +195,16 @@ public:
 private:
   Real data_[3][3];
 
-  void set(Real i00, Real i01, Real i02,Real i10, Real i11, Real i12,Real i20, Real i21, Real i22);
+  void set(Real i00, Real i01, Real i02,
+           Real i10, Real i11, Real i12,
+           Real i20, Real i21, Real i22)
+  {
+    data_[0][0]=i00; data_[0][1]=i01; data_[0][2]=i02;
+    data_[1][0]=i10; data_[1][1]=i11; data_[1][2]=i12;
+    data_[2][0]=i20; data_[2][1]=i21; data_[2][2]=i22;
+  }
 };
-
+  
 DLLEXPORT_OST_GEOM std::ostream& operator<<(std::ostream& o, const Mat3& m);
 
 } // ns geom
diff --git a/modules/geom/src/quat.cc b/modules/geom/src/quat.cc
index 3d5f3f2f9..128d540f0 100644
--- a/modules/geom/src/quat.cc
+++ b/modules/geom/src/quat.cc
@@ -462,21 +462,21 @@ Vec3 Quat::Rotate(const Vec3& vec) const {
 Quat Grassmann(const Quat& lhs, const Quat& rhs)
 {
   return Quat(lhs.GetAngle()*rhs.GetAngle()-
-              lhs.GetAxis().GetX()*rhs.GetAxis().GetX()-
-              lhs.GetAxis().GetY()*rhs.GetAxis().GetY()-
-              lhs.GetAxis().GetZ()*rhs.GetAxis().GetZ(),
-                    lhs.GetAngle()*rhs.GetAxis().GetX()+
-                    lhs.GetAxis().GetX()*rhs.GetAngle()+
-                    lhs.GetAxis().GetY()*rhs.GetAxis().GetZ()-
-                    lhs.GetAxis().GetZ()*rhs.GetAxis().GetY(),
-               lhs.GetAngle()*rhs.GetAxis().GetY()-
-               lhs.GetAxis().GetX()*rhs.GetAxis().GetZ()+
-               lhs.GetAxis().GetY()*rhs.GetAngle()+
-               lhs.GetAxis().GetZ()*rhs.GetAxis().GetX(),
-                    lhs.GetAngle()*rhs.GetAxis().GetZ()+
-                    lhs.GetAxis().GetX()*rhs.GetAxis().GetY()-
-                    lhs.GetAxis().GetY()*rhs.GetAxis().GetX()+
-                    lhs.GetAxis().GetZ()*rhs.GetAngle());
+              lhs.GetAxis().x*rhs.GetAxis().x-
+              lhs.GetAxis().y*rhs.GetAxis().y-
+              lhs.GetAxis().z*rhs.GetAxis().z,
+                    lhs.GetAngle()*rhs.GetAxis().x+
+                    lhs.GetAxis().x*rhs.GetAngle()+
+                    lhs.GetAxis().y*rhs.GetAxis().z-
+                    lhs.GetAxis().z*rhs.GetAxis().y,
+               lhs.GetAngle()*rhs.GetAxis().y-
+               lhs.GetAxis().x*rhs.GetAxis().z+
+               lhs.GetAxis().y*rhs.GetAngle()+
+               lhs.GetAxis().z*rhs.GetAxis().x,
+                    lhs.GetAngle()*rhs.GetAxis().z+
+                    lhs.GetAxis().x*rhs.GetAxis().y-
+                    lhs.GetAxis().y*rhs.GetAxis().x+
+                    lhs.GetAxis().z*rhs.GetAngle());
 }
 
 std::ostream& operator<<(std::ostream& str, const Quat& q)
diff --git a/modules/geom/src/vec2.cc b/modules/geom/src/vec2.cc
index 129a2745c..4c30939b5 100644
--- a/modules/geom/src/vec2.cc
+++ b/modules/geom/src/vec2.cc
@@ -22,158 +22,3 @@
 
 #include "exc.hh"
 
-namespace geom {
-
-Vec2::Vec2()
-{
-  this->set(0.0,0.0);
-}
-
-Vec2::Vec2(Real x, Real y)
-{
-  this->set(x,y);
-}
-
-Vec2::Vec2(const Vec2& v)
-{
-  this->set(v.data_[0],v.data_[1]);
-}
-
-Vec2::Vec2(const Vec3& v)
-{
-  this->set(v[0],v[1]); // z is silenly swallowed
-}
-
-Vec2::Vec2(const Vec4& v) 
-{
-  if(v[3]==0.0) {
-    throw DivideByZeroException("zero w component in Vec4->Vec2 conversion");
-  }
-  this->set(v[0]/v[3],v[1]/v[3]); // z is silently swallowed
-}
-
-Vec2::Vec2(const Real d[2])
-{
-  this->set(d[0],d[1]);
-}
-
-/*Vec2::Vec2(const float d[2])
-{
-  this->set(d[0],d[1]);
-}
-*/
-Vec2& Vec2::operator=(const Vec2& v)
-{
-  if(&v!=this) {
-    this->set(v[0],v[1]);
-  }
-  return *this;
-}
-
-bool Vec2::operator==(const Vec2& v) const
-{
-  return data_[0]==v.data_[0] && 
-    data_[1]==v.data_[1];
-}
-
-
-Real& Vec2::operator[](std::size_t indx)
-{
-  if(indx>1) throw OutOfRangeException();
-  return data_[indx];
-}
-
-
-const Real& Vec2::operator[](std::size_t indx) const
-{
-  if(indx>1) throw OutOfRangeException();
-  return data_[indx];
-}
-
-
-Vec2 Vec2::operator-() const
-{
-  Vec2 nrvo(-data_[0],-data_[1]);
-  return nrvo;
-}
-
-Vec2& Vec2::operator+=(const Vec2& rhs)
-{
-  data_[0]+=rhs[0]; data_[1]+=rhs[1];
-  return *this;
-}
-
-Vec2& Vec2::operator-=(const Vec2& rhs)
-{
-  data_[0]-=rhs[0]; data_[1]-=rhs[1];
-  return *this;
-}
-
-Vec2& Vec2::operator*=(Real d)
-{
-  data_[0]*=d; data_[1]*=d;
-  return *this;
-}
-
-Vec2& Vec2::operator+=(Real d)
-{
-  data_[0]+=d; data_[1]+=d;
-  return *this;
-}
-
-Vec2& Vec2::operator-=(Real d)
-{
-  data_[0]-=d; data_[1]-=d;
-  return *this;
-}
-
-Vec2& Vec2::operator/=(Real d)
-{
-  data_[0]/=d; data_[1]/=d;
-  return *this;
-}
-
-////////////////////////////////////////////
-// private member functions
-
-void Vec2::set(Real x, Real y)
-{
-  data_[0]=x; data_[1]=y;
-}
-
-////////////////////////////////////////////
-// functions
-
-Vec2 operator/(Real d, const Vec2& v)
-{
-  Vec2 nrvo(d/v[0],d/v[1]);
-  return nrvo;
-}
-
-std::ostream& operator<<(std::ostream& os, const Vec2& v)
-{
-  os << "(" << v[0] << "," << v[1] << ")";
-  return os;
-}
-
-Real Vec2::GetX() const
-{
-  return data_[0];
-}
-
-Real Vec2::GetY() const
-{
-  return data_[1];
-}
-
-void Vec2::SetX(Real v)
-{
-  data_[0]=v;
-}
-void Vec2::SetY(Real v)
-{
-  data_[1]=v;
-}
-
-
-} // ns
diff --git a/modules/geom/src/vec2.hh b/modules/geom/src/vec2.hh
index 199b6a760..391e073ff 100644
--- a/modules/geom/src/vec2.hh
+++ b/modules/geom/src/vec2.hh
@@ -24,8 +24,9 @@
 #include <vector>
 #include <boost/operators.hpp>
 
-#include <ost/geom/module_config.hh>
+
 #include <ost/config.hh>
+#include <ost/geom/module_config.hh>
 
 namespace geom {
 
@@ -36,7 +37,7 @@ class Vec4;
 /*
   Two dimensional vector class, using Real precision.
 */
-class DLLEXPORT_OST_GEOM Vec2:
+class DLLEXPORT Vec2:
     private boost::equality_comparable<Vec2>,
     private boost::additive<Vec2>,
     private boost::additive<Vec2, Real>,
@@ -44,13 +45,13 @@ class DLLEXPORT_OST_GEOM Vec2:
 {
 public:
   //! Default initialization, all components are set to zero
-  Vec2();
+  Vec2(): x(0), y(0) { }
 
   //! Initialization with x, y and z component
-  Vec2(Real x, Real y);
+  Vec2(Real px, Real py): x(px), y(py) { }
 
   //! copy ctor
-  Vec2(const Vec2& v);
+  Vec2(const Vec2& v): x(v.x), y(v.y) { }
 
   //! explicit initialization with 3D vector
   explicit Vec2(const Vec3& v);
@@ -59,53 +60,120 @@ public:
   explicit Vec2(const Vec4& v);
 
   //! explicit initialization with an array of doubles
-  explicit Vec2(const Real[2]);
+  explicit Vec2(const float v[2]): x(v[0]), y(v[1]) { }
+
 
-#if OST_DOUBLE_PRECISION
   //! explicit initialization with an array of floats
-  explicit Vec2(const float[2]);
-#endif
-  //! assignement op
-  Vec2& operator=(const Vec2& v);
+  explicit Vec2(const double v[2]): x(v[0]), y(v[1]) { }
 
-  //! element access
-  Real& operator[](std::size_t indx);
-  //! const element access
-  const Real& operator[](std::size_t indx) const;
-  //! element access
-  Real GetX() const;
-  Real GetY() const;
-  void SetX(Real v);
-  void SetY(Real v);
+
+  Real GetX() const { return x; }
+  Real GetY() const { return y; }
+
+  void SetX(Real d) { x=d; }
+  void SetY(Real d) { y=d; }
 
   //! comparable
-  bool operator==(const Vec2& rhs) const;
+  bool operator==(const Vec2& rhs) const
+  {
+    return x==rhs.x && y==rhs.y;
+  }
 
+  //! element access
+  Real& operator[](std::size_t indx)
+  {
+    return (&x)[indx];
+  }
+  
+  //! const element access
+  const Real& operator[](std::size_t indx) const
+  {
+    
+    return (&x)[indx];
+  }
+  
   //! addable op
-  Vec2& operator+=(const Vec2& rhs);
-  Vec2& operator+=(Real d);
+  Vec2& operator+=(const Vec2& rhs)
+  {
+    x+=rhs.x;
+    y+=rhs.y;
+    return *this;
+  }
+  
+  Vec2& operator+=(Real d)
+  {
+    x+=d;
+    y+=d;
+    return *this;
+  }
+  
   //! subtractable op
-  Vec2& operator-=(const Vec2& rhs);
-  Vec2& operator-=(Real d);
+  Vec2& operator-=(const Vec2& rhs)
+  {
+    x-=rhs.x;
+    y-=rhs.y;
+    return *this;
+  }
+  
+  Vec2& operator-=(Real d)
+  {
+    x-=d;
+    y-=d;
+    return *this;
+  }
   //! negateable
-  Vec2 operator-() const;
+  Vec2 operator-() const
+  {
+    return Vec2(-x, -y);
+  }
+
   //! multipliable
-  Vec2& operator*=(Real d);
+  Vec2& operator*=(Real d)
+  {
+    x*=d;
+    y*=d;
+    return *this;
+  }
+
   //! dividable
-  Vec2& operator/=(Real d);
+  Vec2& operator/=(Real d)
+  {
+    Real one_over_d=Real(1.0)/d;
+    x*=one_over_d;
+    y*=one_over_d;
+    return *this;
+  }
+
+  Real* Data() {return &x;}
+  const Real* Data() const { return &x; }
+
+  Real x;
+  Real y;
+};
 
-  Real* Data() {return data_;}
-  const Real* Data() const {return data_;}
+inline Vec2 operator/(Real d, const Vec2& v)
+{
+  return Vec2(d/v.x, d/v.y);
+}
 
-private:
-  Real data_[2];
+inline std::ostream& operator<<(std::ostream& os, const Vec2& v)
+{
+  os << "(" << v[0] << "," << v[1] << ")";
+  return os;  
+}
+
+}
+
+#include <ost/geom/vec3.hh>
+#include <ost/geom/vec4.hh>
+
+namespace geom {
+
+inline Vec2::Vec2(const Vec3& v): x(v.x), y(v.y) { }
 
-  void set(Real x, Real y);
-};
 
-DLLEXPORT_OST_GEOM Vec2 operator/(Real d, const Vec2& v);
+inline Vec2::Vec2(const Vec4& v): x(v.x), y(v.y) { }
 
-DLLEXPORT_OST_GEOM std::ostream& operator<<(std::ostream&, const Vec2&);
 
 typedef std::vector<Vec2> Vec2List;
 
diff --git a/modules/geom/src/vec3.cc b/modules/geom/src/vec3.cc
index e84c2b64a..8df7803ce 100644
--- a/modules/geom/src/vec3.cc
+++ b/modules/geom/src/vec3.cc
@@ -22,181 +22,4 @@
 
 #include "exc.hh"
 #include <ost/config.hh>
-namespace geom {
-
-Vec3::Vec3()
-{
-  this->set(0.0,0.0,0.0);
-}
-
-Vec3::Vec3(Real x, Real y, Real z)
-{
-  this->set(x,y,z);
-}
-
-Vec3::Vec3(Real v)
-{
-  this->set(v, v, v);
-}
-
-
-//#if OST_DOUBLE_PRECISION
-//Vec3::Vec3(float x, float y, float z)
-//{
-//  this->set(x,y,z);
-//}
-//#endif
-
-Vec3::Vec3(const Vec3& v)
-{
-  this->set(v.data_[0],v.data_[1],v.data_[2]);
-}
-
-Vec3::Vec3(const Vec2& v)
-{
-  this->set(v[0],v[1],0.0);
-}
-
-Vec3::Vec3(const Vec4& v)
-{
-  if(v[3]==0.0) {
-    throw DivideByZeroException("zero w component in Vec4->Vec3 conversion");
-  }
-  this->set(v[0]/v[3],v[1]/v[3],v[2]/v[3]);
-}
-
-
-Vec3::Vec3(const Real d[3])
-{
-  this->set(d[0],d[1],d[2]);
-}
-#if OST_DOUBLE_PRECISION
-Vec3::Vec3(const float d[3])
-{
-  this->set(d[0],d[1],d[2]);
-}
-#endif
-
-Vec3& Vec3::operator=(const Vec3& v)
-{
-  if(&v!=this) {
-    this->set(v[0],v[1],v[2]);
-  }
-  return *this;
-}
-
-bool Vec3::operator==(const Vec3& v) const
-{
-  return data_[0]==v.data_[0] &&
-    data_[1]==v.data_[1] &&
-    data_[2]==v.data_[2];
-}
-
-Real& Vec3::operator[](std::size_t indx)
-{
-  if(indx>2) throw OutOfRangeException();
-  return data_[indx];
-}
-
-
-const Real& Vec3::operator[](std::size_t indx) const
-{
-  if(indx>2) throw OutOfRangeException();
-  return data_[indx];
-}
-
-Real Vec3::GetX() const
-{
-  return data_[0];
-}
-
-Real Vec3::GetY() const
-{
-  return data_[1];
-}
-
-Real Vec3::GetZ() const
-{
-  return data_[2];
-}
-
-void Vec3::SetX(Real v)
-{
-  data_[0]=v;
-}
-void Vec3::SetY(Real v)
-{
-  data_[1]=v;
-}
-void Vec3::SetZ(Real v)
-{
-  data_[2]=v;
-}
-
-Vec3 Vec3::operator-() const
-{
-  Vec3 nrvo(-data_[0],-data_[1],-data_[2]);
-  return nrvo;
-}
-
-Vec3& Vec3::operator+=(const Vec3& rhs)
-{
-  data_[0]+=rhs[0]; data_[1]+=rhs[1]; data_[2]+=rhs[2];
-  return *this;
-}
-
-Vec3& Vec3::operator-=(const Vec3& rhs)
-{
-  data_[0]-=rhs[0]; data_[1]-=rhs[1]; data_[2]-=rhs[2];
-  return *this;
-}
-
-Vec3& Vec3::operator*=(Real d)
-{
-  data_[0]*=d; data_[1]*=d; data_[2]*=d;
-  return *this;
-}
-
-Vec3& Vec3::operator+=(Real d)
-{
-  data_[0]+=d; data_[1]+=d; data_[2]+=d;
-  return *this;
-}
-
-Vec3& Vec3::operator-=(Real d)
-{
-  data_[0]-=d; data_[1]-=d; data_[2]-=d;
-  return *this;
-}
-
-Vec3& Vec3::operator/=(Real d)
-{
-  data_[0]/=d; data_[1]/=d; data_[2]/=d;
-  return *this;
-}
-
-////////////////////////////////////////////
-// private member functions
-
-void Vec3::set(Real x, Real y, Real z)
-{
-  data_[0]=x; data_[1]=y; data_[2]=z;
-}
-
-////////////////////////////////////////////
-// functions
-
-Vec3 operator/(Real d, const Vec3& v)
-{
-  Vec3 nrvo(d/v[0],d/v[1],d/v[2]);
-  return nrvo;
-}
-
-std::ostream& operator<<(std::ostream& os, const Vec3& v)
-{
-  os << "(" << v[0] << "," << v[1] << "," << v[2] << ")";
-  return os;
-}
-
-} // ns geom
 
diff --git a/modules/geom/src/vec3.hh b/modules/geom/src/vec3.hh
index ac64e76f7..627de4c1e 100644
--- a/modules/geom/src/vec3.hh
+++ b/modules/geom/src/vec3.hh
@@ -24,8 +24,9 @@
 #include <vector>
 #include <boost/operators.hpp>
 
-#include <ost/geom/module_config.hh>
+
 #include <ost/config.hh>
+#include <ost/geom/module_config.hh>
 
 namespace geom {
 
@@ -43,13 +44,13 @@ class DLLEXPORT_OST_GEOM Vec3:
 {
 public:
   //! Default initialization, all components are set to zero
-  Vec3();
+  Vec3(): x(0), y(0), z(0) {}
 
   //! Initialization with x, y and z component
-  Vec3(Real x, Real y, Real z);
+  Vec3(Real px, Real py, Real pz): x(px), y(py), z(pz) {}
 
   //! copy ctor
-  Vec3(const Vec3& v);
+  Vec3(const Vec3& v): x(v.x), y(v.y), z(v.z) { }
 
   //! (implicit) initialization with 2D vector
   Vec3(const Vec2& v);
@@ -60,62 +61,139 @@ public:
     to a 3D vector, resulting in (x/w,y/w,z/w)
   */
   explicit Vec3(const Vec4& v);
-
+  
+  explicit Vec3(Real v): x(v), y(v), z(v) { }
+  
   //! explicit initialization with an array of doubles
-  explicit Vec3(const Real[3]);
+  explicit Vec3(const double v[3]): x(v[0]), y(v[1]), z(v[2]) { }
 
-  explicit Vec3(Real v);
-  
-#if OST_DOUBLE_PRECISION
   //! explicit initialization with an array of floats
-  explicit Vec3(const float[3]);
-#endif
+  explicit Vec3(const float v[3]): x(v[0]), y(v[1]), z(v[2]) { }
+
   //! assignement op
-  Vec3& operator=(const Vec3& v);
+  Vec3& operator=(const Vec3& v)
+  {
+    x=v.x;
+    y=v.y;
+    z=v.z;
+    return *this;
+  }
 
   //! comparable
-  bool operator==(const Vec3& rhs) const;
+  bool operator==(const Vec3& rhs) const
+  {
+    return x==rhs.x && y==rhs.y && z==rhs.z;
+  }
 
   //! element access
-  Real& operator[](std::size_t indx);
+  Real& operator[](std::size_t indx)
+  {
+    return (&x)[indx];
+  }
+  
   //! const element access
-  const Real& operator[](std::size_t indx) const;
+  const Real& operator[](std::size_t indx) const
+  {
+    
+    return (&x)[indx];
+  }
   //! element access
-  Real GetX() const;
-  Real GetY() const;
-  Real GetZ() const;
-  void SetX(Real v);
-  void SetY(Real v);
-  void SetZ(Real v);
+  Real GetX() const { return x; }
+  Real GetY() const { return y; }
+  Real GetZ() const { return z; }
+  void SetX(Real v) { x=v; }
+  void SetY(Real v) { y=v; }
+  void SetZ(Real v) { z=v; }
 
   //! addable op
-  Vec3& operator+=(const Vec3& rhs);
-  Vec3& operator+=(Real d);
+  Vec3& operator+=(const Vec3& rhs)
+  {
+    x+=rhs.x;
+    y+=rhs.y;
+    z+=rhs.z;
+    return *this;
+  }
+  
+  Vec3& operator+=(Real d)
+  {
+    x+=d;
+    y+=d;
+    z+=d;
+    return *this;
+  }
+  
   //! subtractable op
-  Vec3& operator-=(const Vec3& rhs);
-  Vec3& operator-=(Real d);
+  Vec3& operator-=(const Vec3& rhs)
+  {
+    x-=rhs.x;
+    y-=rhs.y;
+    z-=rhs.z;
+    return *this;
+  }
+  
+  Vec3& operator-=(Real d)
+  {
+    x-=d;
+    y-=d;
+    z-=d;
+    return *this;
+  }
   //! negateable
-  Vec3 operator-() const;
+  Vec3 operator-() const
+  {
+    return Vec3(-x, -y, -z);
+  }
+  
   //! multipliable
-  Vec3& operator*=(Real d);
+  Vec3& operator*=(Real d)
+  {
+    x*=d;
+    y*=d;
+    z*=d;
+    return *this;
+  }
+  
   //! dividable
-  Vec3& operator/=(Real d);
-
-  Real* Data() {return data_;}
-  const Real* Data() const {return data_;}
+  Vec3& operator/=(Real d)
+  {
+    Real one_over_d=Real(1.0)/d;
+    x*=one_over_d;
+    y*=one_over_d;
+    z*=one_over_d;
+    return *this;
+  }
+
+  Real* Data() {return &x;}
+  const Real* Data() const {return &x;}
+
+  Real x;
+  Real y;
+  Real z;
+};
 
-private:
-  Real data_[3];
+inline Vec3 operator/(Real d, const Vec3& v)
+{
+  Vec3 nrvo(d/v[0],d/v[1],d/v[2]);
+  return nrvo;
+}
 
-  void set(Real x, Real y, Real z);
-};
+inline std::ostream& operator<<(std::ostream& os, const Vec3& v)
+{
+  os << "[" << v.x << ", " << v.y << "]";
+  return os;
+}
+}
 
-DLLEXPORT_OST_GEOM Vec3 operator/(Real, const Vec3& v);
+#include <ost/geom/vec2.hh>
+#include <ost/geom/vec4.hh>
 
-DLLEXPORT_OST_GEOM std::ostream& operator<<(std::ostream&, const Vec3&);
+namespace geom {
 
 typedef std::vector<Vec3> Vec3List;
 
+inline Vec3::Vec3(const Vec2& v): x(v.x), y(v.y), z(0.0) { }
+inline Vec3::Vec3(const Vec4& v): x(v.x/v.w), y(v.y/v.w), z(v.z/v.w) { }
+  
 } // namespace geom
 
 
diff --git a/modules/geom/src/vec4.cc b/modules/geom/src/vec4.cc
index 8f53869a9..4c30939b5 100644
--- a/modules/geom/src/vec4.cc
+++ b/modules/geom/src/vec4.cc
@@ -22,136 +22,3 @@
 
 #include "exc.hh"
 
-namespace geom {
-
-Vec4::Vec4()
-{
-  this->set(0.0,0.0,0.0,1.0);
-}
-
-Vec4::Vec4(Real x, Real y, Real z, Real w)
-{
-  this->set(x,y,z,w);
-}
-
-Vec4::Vec4(const Vec4& v)
-{
-  this->set(v.data_[0],v.data_[1],v.data_[2],v.data_[3]);
-}
-
-Vec4::Vec4(const Vec2& v)
-{
-  this->set(v[0],v[1],0.0,1.0);
-}
-
-Vec4::Vec4(const Vec3& v)
-{
-  this->set(v[0],v[1],v[2],1.0);
-}
-
-Vec4::Vec4(const Real d[4])
-{
-  this->set(d[0],d[1],d[2],d[3]);
-}
-
-#if OST_DOUBLE_PRECISION
-Vec4::Vec4(const float d[4])
-{
-  this->set(d[0],d[1],d[2],d[3]);
-}
-#endif
-
-Vec4& Vec4::operator=(const Vec4& v)
-{
-  if(&v!=this) {
-    this->set(v[0],v[1],v[2],v[3]);
-  }
-  return *this;
-}
-
-bool Vec4::operator==(const Vec4& v) const
-{
-  return data_[0]==v.data_[0] && 
-    data_[1]==v.data_[1] &&
-    data_[2]==v.data_[2] &&
-    data_[3]==v.data_[3];
-}
-
-Real& Vec4::operator[](std::size_t indx)
-{
-  if(indx>3) throw OutOfRangeException();
-  return data_[indx];
-}
-
-const Real& Vec4::operator[](std::size_t indx) const
-{
-  if(indx>3) throw OutOfRangeException();
-  return data_[indx];
-}
-
-Vec4 Vec4::operator-() const
-{
-  Vec4 nrvo(-data_[0],-data_[1],-data_[2],-data_[3]);
-  return nrvo;
-}
-
-Vec4& Vec4::operator+=(const Vec4& rhs)
-{
-  data_[0]+=rhs[0]; data_[1]+=rhs[1]; data_[2]+=rhs[2]; data_[3]+=rhs[3];
-  return *this;
-}
-
-Vec4& Vec4::operator-=(const Vec4& rhs)
-{
-  data_[0]-=rhs[0]; data_[1]-=rhs[1]; data_[2]-=rhs[2]; data_[3]-=rhs[3];
-  return *this;
-}
-
-Vec4& Vec4::operator*=(Real d)
-{
-  data_[0]*=d; data_[1]*=d; data_[2]*=d; data_[3]*=d;
-  return *this;
-}
-
-Vec4& Vec4::operator+=(Real d)
-{
-  data_[0]+=d; data_[1]+=d; data_[2]+=d; data_[3]+=d;
-  return *this;
-}
-
-Vec4& Vec4::operator-=(Real d)
-{
-  data_[0]-=d; data_[1]-=d; data_[2]-=d; data_[3]-=d;
-  return *this;
-}
-
-Vec4& Vec4::operator/=(Real d)
-{
-  data_[0]/=d; data_[1]/=d; data_[2]/=d; data_[3]/=d;
-  return *this;
-}
-
-////////////////////////////////////////////
-// private member functions
-
-void Vec4::set(Real x, Real y, Real z, Real w)
-{
-  data_[0]=x; data_[1]=y; data_[2]=z; data_[3]=w;
-}
-
-////////////////////////////////////////////
-// functions
-
-Vec4 operator/(Real d, const Vec4& v)
-{
-  Vec4 nrvo(d/v[0],d/v[1],d/v[2],d/v[3]);
-  return nrvo;
-}
-
-std::ostream& operator<<(std::ostream& os, const Vec4& v)
-{
-  os << "(" << v[0] << "," << v[1] << "," << v[2] << "," << v[3] << ")";
-  return os;
-}
-
-}
diff --git a/modules/geom/src/vec4.hh b/modules/geom/src/vec4.hh
index b114713dd..ad3367b97 100644
--- a/modules/geom/src/vec4.hh
+++ b/modules/geom/src/vec4.hh
@@ -32,12 +32,11 @@ namespace geom {
 // fw decl
 class Vec2;
 class Vec3;
-class Mat4;
 
 /*
   Four dimensional, homogeneous vector class, using Real precision.
 */
-class DLLEXPORT_OST_GEOM Vec4:
+class DLLEXPORT Vec4:
     private boost::equality_comparable<Vec4>,
     private boost::additive<Vec4>,
     private boost::additive<Vec4, Real>,
@@ -45,13 +44,13 @@ class DLLEXPORT_OST_GEOM Vec4:
 {
 public:
   //! Default initialization, all components are set to zero
-  Vec4();
+  Vec4(): x(0), y(0), z(0), w(0) { }
 
   //! Initialization with x, y and z component
-  Vec4(Real x, Real y, Real z, Real w);
+  Vec4(Real px, Real py, Real pz, Real pw): x(px), y(py), z(pz), w(pw) { }
 
   //! copy ctor
-  Vec4(const Vec4& v);
+  Vec4(const Vec4& v): x(v[0]), y(v[1]), z(v[2]), w(v[3]) { }
 
   //! (implicit) initialization with 2D vector
   Vec4(const Vec2& v);
@@ -59,50 +58,146 @@ public:
   //! (implicit) initialization with 3D vector
   Vec4(const Vec3& v);
 
-  //! explicit initialization with an array of doubles
-  explicit Vec4(const Real[4]);
-
-#if OST_DOUBLE_PRECISION
   //! explicit initialization with an array of floats
-  explicit Vec4(const float[4]);
-#endif
-
+  explicit Vec4(const float v[4]): x(v[0]), y(v[1]), z(v[2]), w(v[3]) { }
+  
+  //! explicit initialization with an array of doubles  
+  explicit Vec4(const double v[4]): x(v[0]), y(v[1]), z(v[2]), w(v[3]) { }
   //! assignement op
-  Vec4& operator=(const Vec4& v);
+  Vec4& operator=(const Vec4& v)
+  {
+    x=v.x;
+    y=v.y;
+    z=v.z;
+    w=v.w;
+    return *this;
+  }
 
   //! comparable
-  bool operator==(const Vec4& rhs) const;
+  bool operator==(const Vec4& rhs) const
+  {
+    return x==rhs.x && y==rhs.y && z==rhs.z && w==rhs.w;
+  }
 
   //! element access
-  Real& operator[](std::size_t indx);
+  Real& operator[](std::size_t indx)
+  {
+    return (&x)[indx];
+  }
+  
   //! const element access
-  const Real& operator[](std::size_t indx) const;
+  const Real& operator[](std::size_t indx) const
+  {
+    
+    return (&x)[indx];
+  }
+  //! element access
+  Real GetX() const { return x; }
+  Real GetY() const { return y; }
+  Real GetZ() const { return z; }
+  Real GetW() const { return w; }
+  
+  void SetX(Real v) { x=v; }
+  void SetY(Real v) { y=v; }
+  void SetZ(Real v) { z=v; }
+  void SetW(Real v) { w=v; }  
 
   //! addable op
-  Vec4& operator+=(const Vec4& rhs);
-  Vec4& operator+=(Real d);
+  Vec4& operator+=(const Vec4& rhs)
+  {
+    x+=rhs.x;
+    y+=rhs.y;
+    z+=rhs.z;
+    w+=rhs.w;
+    return *this;
+  }
+  
+  Vec4& operator+=(Real d)
+  {
+    x+=d;
+    y+=d;
+    z+=d;
+    w+=d;
+    return *this;
+  }
+  
   //! subtractable op
-  Vec4& operator-=(const Vec4& rhs);
-  Vec4& operator-=(Real d);
+  Vec4& operator-=(const Vec4& rhs)
+  {
+    x-=rhs.x;
+    y-=rhs.y;
+    z-=rhs.z;
+    w-=rhs.w;
+    return *this;
+  }
+  
+  Vec4& operator-=(Real d)
+  {
+    x-=d;
+    y-=d;
+    z-=d;
+    w-=d;
+    return *this;
+  }
   //! negateable
-  Vec4 operator-() const;
+  Vec4 operator-() const
+  {
+    return Vec4(-x, -y, -z, -w);
+  }
+  
   //! multipliable
-  Vec4& operator*=(Real d);
+  Vec4& operator*=(Real d)
+  {
+    x*=d;
+    y*=d;
+    z*=d;
+    w*=d;
+    return *this;
+  }
+  
   //! dividable
-  Vec4& operator/=(Real d);
+  Vec4& operator/=(Real d)
+  {
+    Real one_over_d=Real(1.0)/d;
+    x*=one_over_d;
+    y*=one_over_d;
+    z*=one_over_d;
+    w*=one_over_d;
+    return *this;
+  }
+
+  Real* Data() {return &x;}
+  const Real* Data() const {return &x;}
+
+  Real x;
+  Real y;
+  Real z;
+  Real w;
+};
 
-  Real* Data() {return data_;}
-  const Real* Data() const {return data_;}
+inline Vec4 operator/(Real d, const Vec4& v)
+{
+  Vec4 nrvo(d/v[0],d/v[1],d/v[2],d/v[3]);
+  return nrvo;
+}
 
-private:
-  Real data_[4];
+inline std::ostream& operator<<(std::ostream& os, const Vec4& v)
+{
+  os << "(" << v[0] << "," << v[1] << "," << v[2] << "," << v[3] << ")";
+  return os;
+}
+}
 
-  void set(Real x, Real y, Real z, Real w);
-};
+#include <ost/geom/vec2.hh>
+#include <ost/geom/vec3.hh>
 
-DLLEXPORT_OST_GEOM Vec4 operator/(Real, const Vec4& v);
+namespace geom {
+  
+//! (implicit) initialization with 2D vector
+inline Vec4::Vec4(const Vec2& v): x(v.x), y(v.y), z(0), w(1) { }
 
-DLLEXPORT_OST_GEOM std::ostream& operator<<(std::ostream&, const Vec4&);
+//! (implicit) initialization with 3D vector
+inline Vec4::Vec4(const Vec3& v): x(v.x), y(v.y), z(v.z), w(1.0) { }
 
 } // namespace geom
 
diff --git a/modules/geom/src/vecmat2_op.cc b/modules/geom/src/vecmat2_op.cc
index 6bda1a12a..c23866293 100644
--- a/modules/geom/src/vecmat2_op.cc
+++ b/modules/geom/src/vecmat2_op.cc
@@ -26,70 +26,6 @@
 #include "mat2.hh"
 
 namespace geom {
-
-Real Length(const Vec2& v)
-{
-  return std::sqrt(Length2(v));
-}
-
-Real Length2(const Vec2& v)
-{
-  return v[0]*v[0]+v[1]*v[1];
-}
-
-bool Equal(const Vec2& v1, const Vec2& v2, Real ephilon)
-{
-  return std::fabs(v1[0]-v2[0])<ephilon &&
-    std::fabs(v1[1]-v2[1])<ephilon;
-}
-bool Equal(const Mat2& m1, const Mat2& m2, Real ephilon)
-{
-  return std::fabs(m1(0,0)-m2(0,0))<ephilon &&
-    std::fabs(m1(0,1)-m2(0,1))<ephilon &&
-    std::fabs(m1(1,0)-m2(1,0))<ephilon &&
-    std::fabs(m1(1,1)-m2(1,1))<ephilon;
-}
-
-Real Dot(const Vec2& v1, const Vec2& v2)
-{
-  return v1[0]*v2[0]+v1[1]*v2[1];
-}
-
-Vec2 Normalize(const Vec2& v)
-{
-  Real l=Length(v);
-  if(l==0.0) {
-    return v;
-  }
-  return v/l;
-}
-
-
-Vec2 CompMultiply(const Vec2& v1, const Vec2& v2)
-{
-  Vec2 nrvo(v1[0]*v2[0],v1[1]*v2[1]);
-  return nrvo;
-}
-
-Vec2 CompDivide(const Vec2& v1, const Vec2& v2)
-{
-  Vec2 nrvo(v1[0]/v2[0],v1[1]/v2[1]);
-  return nrvo;
-}
-
-Vec2 operator*(const Vec2& v,const Mat2& m)
-{
-  Vec2 nrvo(v[0]*m(0,0)+v[1]*m(1,0),
-      v[0]*m(0,1)+v[1]*m(1,1));
-  return nrvo;
-}
-
-Vec2 operator*(const Mat2& m, const Vec2& v)
-{
-  Vec2 nrvo(v[0]*m(0,0)+v[1]*m(0,1),
-      v[0]*m(1,0)+v[1]*m(1,1));
-  return nrvo;
-}
 // determinant
 Real Det(const Mat2& m)
 {
@@ -116,14 +52,6 @@ Real Angle(const Vec2& v1, const Vec2& v2)
 {
   return std::acos(Dot(v1,v2)/Length(v1)/Length(v2));
 }
-Mat2 operator*(const Mat2& m1, const Mat2& m2)
-{
-  Mat2 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0),
-      m1(0,0)*m2(0,1)+m1(0,1)*m2(1,1),
-      m1(1,0)*m2(0,0)+m1(1,1)*m2(1,0),
-      m1(1,0)*m2(0,1)+m1(1,1)*m2(1,1));
-  return nrvo;
-}
 
 Real SignedAngle(const Vec2& v1, const Vec2& v2)
 {
@@ -131,21 +59,7 @@ Real SignedAngle(const Vec2& v1, const Vec2& v2)
   vc=Cross(Vec3(v1),Vec3(v2));
   if(Length(v1)==0.0 || Length(v2)==0.0 || Length(vc)==0.0)
     return 0.0;
-  return acos(Dot(v1,v2)/Length(v1)/Length(v2))*vc.GetZ()/std::fabs(vc.GetZ());
-}
-
-Vec2 Min(const Vec2& v1, const Vec2& v2)
-{
-  Vec2 nrvo(std::min(v1[0],v2[0]),
-            std::min(v1[1],v2[1]));
-  return nrvo;
-}
-
-Vec2 Max(const Vec2& v1, const Vec2& v2)
-{
-  Vec2 nrvo(std::max(v1[0],v2[0]),
-            std::max(v1[1],v2[1]));
-  return nrvo;
+  return acos(Dot(v1,v2)/Length(v1)/Length(v2))*vc.z/std::fabs(vc.z);
 }
 
 Vec2 Rotate(const Vec2& v,Real ang)
diff --git a/modules/geom/src/vecmat2_op.hh b/modules/geom/src/vecmat2_op.hh
index cf2051367..31ae6430a 100644
--- a/modules/geom/src/vecmat2_op.hh
+++ b/modules/geom/src/vecmat2_op.hh
@@ -23,33 +23,71 @@
 
 #include <ost/geom/module_config.hh>
 
+#include <ost/geom/vec2.hh>
+#include <ost/geom/mat2.hh>
+
 namespace geom {
 
 class Vec2;
 class Mat2;
 
-//! returns length of vector
-Real DLLEXPORT_OST_GEOM Length(const Vec2& v);
-
 //! returns squared length of vector
-Real DLLEXPORT_OST_GEOM Length2(const Vec2& v);
+inline Real Length2(const Vec2& v)
+{
+  return v[0]*v[0]+v[1]*v[1];
+}
+
+//! returns length of vector
+inline Real Length(const Vec2& v)
+{
+  return std::sqrt(Length2(v));
+}
 
 //! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Vec2& v1, const Vec2& v2, Real ephilon=EPSILON);
+inline bool Equal(const Vec2& v1, const Vec2& v2, Real ephilon=EPSILON)
+{
+  return std::fabs(v1[0]-v2[0])<ephilon &&
+    std::fabs(v1[1]-v2[1])<ephilon;
+}
+
 //! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Mat2& m1, const Mat2& m2, Real ephilon=EPSILON);
+inline bool Equal(const Mat2& m1, const Mat2& m2, Real ephilon=EPSILON)
+{
+  return std::fabs(m1(0,0)-m2(0,0))<ephilon &&
+    std::fabs(m1(0,1)-m2(0,1))<ephilon &&
+    std::fabs(m1(1,0)-m2(1,0))<ephilon &&
+    std::fabs(m1(1,1)-m2(1,1))<ephilon;
+}
 
 //! vector dot product
-Real DLLEXPORT_OST_GEOM Dot(const Vec2& v1, const Vec2& v2);
+inline Real Dot(const Vec2& v1, const Vec2& v2)
+{
+    return v1[0]*v2[0]+v1[1]*v2[1];
+}
 
 //! Normalize vector
-Vec2 DLLEXPORT_OST_GEOM Normalize(const Vec2& v);
+inline Vec2 Normalize(const Vec2& v)
+{
+  Real l=Length(v);
+  if(l==0.0) {
+    return v;
+  }
+  return v/l;
+}
 
 //! multiply each component of v1 with that of v2
-Vec2 DLLEXPORT_OST_GEOM CompMultiply(const Vec2& v1, const Vec2& v2);
+inline Vec2 CompMultiply(const Vec2& v1, const Vec2& v2)
+{
+  Vec2 nrvo(v1[0]*v2[0],v1[1]*v2[1]);
+  return nrvo;
+}
 
 //! divide each component of v1 by that of v2
-Vec2 DLLEXPORT_OST_GEOM CompDivide(const Vec2& v1, const Vec2& v2);
+inline Vec2 CompDivide(const Vec2& v1, const Vec2& v2)
+{
+  Vec2 nrvo(v1[0]/v2[0],v1[1]/v2[1]);
+  return nrvo;
+}
 
 //! vector matrix multiplication
 /*!
@@ -58,7 +96,12 @@ Vec2 DLLEXPORT_OST_GEOM CompDivide(const Vec2& v1, const Vec2& v2);
           |a b|
   (x,y) * |c d| = (ax+cy, bx+dy)
 */
-Vec2 DLLEXPORT_OST_GEOM operator*(const Vec2& v,const Mat2& m);
+inline Vec2 operator*(const Vec2& v,const Mat2& m)
+{
+  Vec2 nrvo(v[0]*m(0,0)+v[1]*m(1,0),
+      v[0]*m(0,1)+v[1]*m(1,1));
+  return nrvo;
+}
 
 //! vector matrix multiplication
 /*!
@@ -67,7 +110,12 @@ Vec2 DLLEXPORT_OST_GEOM operator*(const Vec2& v,const Mat2& m);
    |a b|   (x)   (ax+by)
    |c d| * (y)  =(cx+dy)
 */
-Vec2 DLLEXPORT_OST_GEOM operator*(const Mat2& m, const Vec2& v);
+inline Vec2 operator*(const Mat2& m, const Vec2& v)
+{
+  Vec2 nrvo(v[0]*m(0,0)+v[1]*m(0,1),
+      v[0]*m(1,0)+v[1]*m(1,1));
+  return nrvo;
+}
 
 //! determinant
 Real DLLEXPORT_OST_GEOM Det(const Mat2& m);
@@ -83,13 +131,28 @@ Real DLLEXPORT_OST_GEOM Angle(const Vec2& v1, const Vec2& v2);
 Real  DLLEXPORT_OST_GEOM SignedAngle(const Vec2& v1, const Vec2& v2);
 
 //! matrix matrix multiplication
-Mat2 DLLEXPORT_OST_GEOM operator*(const Mat2& m1, const Mat2& m2);
-
-//! returns std::min of each component
-Vec2 DLLEXPORT_OST_GEOM Min(const Vec2& v1, const Vec2& v2);
-
-//! returns std::max of each component
-Vec2 DLLEXPORT_OST_GEOM Max(const Vec2& v1, const Vec2& v2);
+inline Mat2 operator*(const Mat2& m1, const Mat2& m2)
+{
+  Mat2 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0),
+      m1(0,0)*m2(0,1)+m1(0,1)*m2(1,1),
+      m1(1,0)*m2(0,0)+m1(1,1)*m2(1,0),
+      m1(1,0)*m2(0,1)+m1(1,1)*m2(1,1));
+  return nrvo;
+}
+
+inline Vec2 Min(const Vec2& v1, const Vec2& v2)
+{
+  Vec2 nrvo(std::min(v1[0],v2[0]),
+            std::min(v1[1],v2[1]));
+  return nrvo;
+}
+
+inline Vec2 Max(const Vec2& v1, const Vec2& v2)
+{
+  Vec2 nrvo(std::max(v1[0],v2[0]),
+            std::max(v1[1],v2[1]));
+  return nrvo;
+}
 
 //! Rotation
 DLLEXPORT Vec2 Rotate(const Vec2& v,Real ang);
diff --git a/modules/geom/src/vecmat3_op.cc b/modules/geom/src/vecmat3_op.cc
index 5487c5d57..8ef740e37 100644
--- a/modules/geom/src/vecmat3_op.cc
+++ b/modules/geom/src/vecmat3_op.cc
@@ -29,102 +29,6 @@
 
 namespace geom {
 
-Real Length(const Vec3& v)
-{
-  return std::sqrt(Length2(v));
-}
-
-Real Length2(const Vec3& v)
-{
-  return v[0]*v[0]+v[1]*v[1]+v[2]*v[2];
-}
-
-bool Equal(const Vec3& v1, const Vec3& v2, Real ephilon)
-{
-  return std::fabs(v1[0]-v2[0])<ephilon &&
-    std::fabs(v1[1]-v2[1])<ephilon &&
-    std::fabs(v1[2]-v2[2])<ephilon;
-}
-
-bool Equal(const Mat3& m1, const Mat3& m2, Real ephilon)
-{
-  return std::fabs(m1(0,0)-m2(0,0))<ephilon &&
-    std::fabs(m1(0,1)-m2(0,1))<ephilon &&
-    std::fabs(m1(0,2)-m2(0,2))<ephilon &&
-    std::fabs(m1(1,0)-m2(1,0))<ephilon &&
-    std::fabs(m1(1,1)-m2(1,1))<ephilon &&
-    std::fabs(m1(1,2)-m2(1,2))<ephilon &&
-    std::fabs(m1(2,0)-m2(2,0))<ephilon &&
-    std::fabs(m1(2,1)-m2(2,1))<ephilon &&
-    std::fabs(m1(2,2)-m2(2,2))<ephilon;
-}
-
-Real Dot(const Vec3& v1, const Vec3& v2)
-{
-  return v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2];
-}
-
-Vec3 Cross(const Vec3& v1, const Vec3& v2)
-{
-  Vec3 nrvo(v1[1]*v2[2]-v2[1]*v1[2],
-      v1[2]*v2[0]-v2[2]*v1[0],
-      v1[0]*v2[1]-v2[0]*v1[1]);
-  return nrvo;
-}
-
-Vec3 Normalize(const Vec3& v)
-{
-  Real l=Length(v);
-  if(l==0.0) {
-    return v;
-  }
-  return v/l;
-}
-
-Vec3 CompMultiply(const Vec3& v1, const Vec3& v2)
-{
-  Vec3 nrvo(v1[0]*v2[0],v1[1]*v2[1],v1[2]*v2[2]);
-  return nrvo;
-}
-
-Vec3 CompDivide(const Vec3& v1, const Vec3& v2)
-{
-  Vec3 nrvo(v1[0]/v2[0],v1[1]/v2[1],v1[2]/v2[2]);
-  return nrvo;
-}
-
-Vec3 operator*(const Vec3& v,const Mat3& m)
-{
-  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(1,0)+v[2]*m(2,0),
-      v[0]*m(0,1)+v[1]*m(1,1)+v[2]*m(2,1),
-      v[0]*m(0,2)+v[1]*m(1,2)+v[2]*m(2,2));
-  return nrvo;
-}
-
-Vec3 operator*(const Mat3& m, const Vec3& v)
-{
-  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(0,1)+v[2]*m(0,2),
-      v[0]*m(1,0)+v[1]*m(1,1)+v[2]*m(1,2),
-      v[0]*m(2,0)+v[1]*m(2,1)+v[2]*m(2,2));
-  return nrvo;
-}
-
-Mat3 operator*(const Mat3& m1, const Mat3& m2)
-{
-  Mat3 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0)+m1(0,2)*m2(2,0),
-      m1(0,0)*m2(0,1)+m1(0,1)*m2(1,1)+m1(0,2)*m2(2,1),
-      m1(0,0)*m2(0,2)+m1(0,1)*m2(1,2)+m1(0,2)*m2(2,2),
-
-      m1(1,0)*m2(0,0)+m1(1,1)*m2(1,0)+m1(1,2)*m2(2,0),
-      m1(1,0)*m2(0,1)+m1(1,1)*m2(1,1)+m1(1,2)*m2(2,1),
-      m1(1,0)*m2(0,2)+m1(1,1)*m2(1,2)+m1(1,2)*m2(2,2),
-
-      m1(2,0)*m2(0,0)+m1(2,1)*m2(1,0)+m1(2,2)*m2(2,0),
-      m1(2,0)*m2(0,1)+m1(2,1)*m2(1,1)+m1(2,2)*m2(2,1),
-      m1(2,0)*m2(0,2)+m1(2,1)*m2(1,2)+m1(2,2)*m2(2,2));
-  return nrvo;
-}
-
 Real Comp(const Mat3& m,unsigned int i_in, unsigned int j_in)
 {
   return Minor(m,i_in,j_in)* static_cast<Real>( ((i_in+j_in)&0x1) ? -1 : 1);
@@ -264,25 +168,4 @@ Mat3 AxisRotation(const Vec3& axis, Real angle)
           xz-ca*xz-sa*y, yz-ca*yz+sa*x,zz+ca-ca*zz);
 }
 
-Vec3 Min(const Vec3& v1, const Vec3& v2)
-{
-  Vec3 nrvo(std::min(v1[0],v2[0]),
-            std::min(v1[1],v2[1]),
-            std::min(v1[2],v2[2]));
-  return nrvo;
-}
-
-Vec3 Max(const Vec3& v1, const Vec3& v2)
-{
-  Vec3 nrvo(std::max(v1[0],v2[0]),
-            std::max(v1[1],v2[1]),
-            std::max(v1[2],v2[2]));
-  return nrvo;
-}
-
-Real  Distance(const Vec3& p1, const Vec3& p2)
-{
-  return Length(p1-p2);
-}
-
 } // ns
diff --git a/modules/geom/src/vecmat3_op.hh b/modules/geom/src/vecmat3_op.hh
index 44c835768..e4ca5a65b 100644
--- a/modules/geom/src/vecmat3_op.hh
+++ b/modules/geom/src/vecmat3_op.hh
@@ -23,48 +23,121 @@
 #include "constants.hh"
 
 #include <ost/geom/module_config.hh>
+#include <ost/geom/vec3.hh>
+#include <ost/geom/mat3.hh>
 
 namespace geom {
 
-class Vec3;
-class Mat3;
-
-//! returns length of vector
-Real DLLEXPORT_OST_GEOM Length(const Vec3& v);
 
 //! returns squared length of vector
-Real DLLEXPORT_OST_GEOM Length2(const Vec3& v);
+inline Real Length2(const Vec3& v)
+{
+  return v[0]*v[0]+v[1]*v[1]+v[2]*v[2];
+}
 
-//! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Vec3& v1, const Vec3& v2, 
-                              Real ephilon=EPSILON);
+//! returns length of vector
+inline Real Length(const Vec3& v)
+{
+  return std::sqrt(Length2(v));
+}
+
+//! return true if components differ not more than epsilon
+inline bool Equal(const Vec3& v1, const Vec3& v2, 
+                  Real ephilon=EPSILON)
+{
+  return std::fabs(v1[0]-v2[0])<ephilon &&
+    std::fabs(v1[1]-v2[1])<ephilon &&
+    std::fabs(v1[2]-v2[2])<ephilon;
+}
 
 //! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Mat3& m1, const Mat3& m2, 
-                              Real ephilon=EPSILON);
+inline bool Equal(const Mat3& m1, const Mat3& m2, 
+                  Real ephilon=EPSILON)
+{
+  return std::fabs(m1(0,0)-m2(0,0))<ephilon &&
+    std::fabs(m1(0,1)-m2(0,1))<ephilon &&
+    std::fabs(m1(0,2)-m2(0,2))<ephilon &&
+    std::fabs(m1(1,0)-m2(1,0))<ephilon &&
+    std::fabs(m1(1,1)-m2(1,1))<ephilon &&
+    std::fabs(m1(1,2)-m2(1,2))<ephilon &&
+    std::fabs(m1(2,0)-m2(2,0))<ephilon &&
+    std::fabs(m1(2,1)-m2(2,1))<ephilon &&
+    std::fabs(m1(2,2)-m2(2,2))<ephilon;
+}
 
 //! vector dot product
-Real DLLEXPORT_OST_GEOM Dot(const Vec3& v1, const Vec3& v2);
+inline Real Dot(const Vec3& v1, const Vec3& v2)
+{
+  return v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2];
+}
 
 //! Normalize
-Vec3 DLLEXPORT_OST_GEOM Normalize(const Vec3& v);
+inline Vec3 Normalize(const Vec3& v)
+{
+  Real l=Length(v);
+  if(l==0.0) {
+    return v;
+  }
+  return v/l;
+}
 
 //! vector cross product
-Vec3 DLLEXPORT_OST_GEOM Cross(const Vec3& v1, const Vec3& v2);
+inline Vec3 Cross(const Vec3& v1, const Vec3& v2)
+{
+  Vec3 nrvo(v1[1]*v2[2]-v2[1]*v1[2],
+      v1[2]*v2[0]-v2[2]*v1[0],
+      v1[0]*v2[1]-v2[0]*v1[1]);
+  return nrvo;
+}
 
 //! multiply each component of v1 with that of v2
-Vec3 DLLEXPORT_OST_GEOM CompMultiply(const Vec3& v1, const Vec3& v2);
+inline Vec3 CompMultiply(const Vec3& v1, const Vec3& v2)
+{
+  Vec3 nrvo(v1[0]*v2[0],v1[1]*v2[1],v1[2]*v2[2]);
+  return nrvo;
+}
 
 //! divide each component of v1 by that of v2
-Vec3 DLLEXPORT_OST_GEOM CompDivide(const Vec3& v1, const Vec3& v2);
+inline Vec3 CompDivide(const Vec3& v1, const Vec3& v2)
+{
+  Vec3 nrvo(v1[0]/v2[0],v1[1]/v2[1],v1[2]/v2[2]);
+  return nrvo;
+}
 
 //! vector matrix multiplication
-Vec3 DLLEXPORT_OST_GEOM operator*(const Vec3& v,const Mat3& m);
+inline Vec3 operator*(const Vec3& v,const Mat3& m)
+{
+  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(1,0)+v[2]*m(2,0),
+      v[0]*m(0,1)+v[1]*m(1,1)+v[2]*m(2,1),
+      v[0]*m(0,2)+v[1]*m(1,2)+v[2]*m(2,2));
+  return nrvo;
+}
+
 //! vector matrix multiplication
-Vec3 DLLEXPORT_OST_GEOM operator*(const Mat3& m, const Vec3& v);
+inline Vec3 operator*(const Mat3& m, const Vec3& v)
+{
+  Vec3 nrvo(v[0]*m(0,0)+v[1]*m(0,1)+v[2]*m(0,2),
+      v[0]*m(1,0)+v[1]*m(1,1)+v[2]*m(1,2),
+      v[0]*m(2,0)+v[1]*m(2,1)+v[2]*m(2,2));
+  return nrvo;
+}
 
 //! matrix matrix multiplication
-Mat3 DLLEXPORT_OST_GEOM operator*(const Mat3& m1, const Mat3& m2);
+inline Mat3 operator*(const Mat3& m1, const Mat3& m2)
+{
+  Mat3 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0)+m1(0,2)*m2(2,0),
+      m1(0,0)*m2(0,1)+m1(0,1)*m2(1,1)+m1(0,2)*m2(2,1),
+      m1(0,0)*m2(0,2)+m1(0,1)*m2(1,2)+m1(0,2)*m2(2,2),
+
+      m1(1,0)*m2(0,0)+m1(1,1)*m2(1,0)+m1(1,2)*m2(2,0),
+      m1(1,0)*m2(0,1)+m1(1,1)*m2(1,1)+m1(1,2)*m2(2,1),
+      m1(1,0)*m2(0,2)+m1(1,1)*m2(1,2)+m1(1,2)*m2(2,2),
+
+      m1(2,0)*m2(0,0)+m1(2,1)*m2(1,0)+m1(2,2)*m2(2,0),
+      m1(2,0)*m2(0,1)+m1(2,1)*m2(1,1)+m1(2,2)*m2(2,1),
+      m1(2,0)*m2(0,2)+m1(2,1)*m2(1,2)+m1(2,2)*m2(2,2));
+  return nrvo;
+}
 
 Mat3 DLLEXPORT_OST_GEOM Invert(const Mat3& m);
 Mat3 DLLEXPORT_OST_GEOM Transpose(const Mat3& m);
@@ -91,13 +164,28 @@ Mat3 DLLEXPORT_OST_GEOM AxisRotation(const Vec3& axis, Real angle);
 Vec3 DLLEXPORT_OST_GEOM OrthogonalVector(const Vec3& axis);
 
 //! returns std::min of each component
-Vec3 DLLEXPORT_OST_GEOM Min(const Vec3& v1, const Vec3& v2);
+inline Vec3 Min(const Vec3& v1, const Vec3& v2)
+{
+  Vec3 nrvo(std::min(v1[0],v2[0]),
+            std::min(v1[1],v2[1]),
+            std::min(v1[2],v2[2]));
+  return nrvo;
+}
 
 //! returns std::max of each component
-Vec3 DLLEXPORT_OST_GEOM Max(const Vec3& v1, const Vec3& v2);
+inline Vec3 Max(const Vec3& v1, const Vec3& v2)
+{
+  Vec3 nrvo(std::max(v1[0],v2[0]),
+            std::max(v1[1],v2[1]),
+            std::max(v1[2],v2[2]));
+  return nrvo;
+}
 
 //! returns the distance between two points
-Real DLLEXPORT_OST_GEOM Distance(const Vec3& p1, const Vec3& p2);
+inline Real Distance(const Vec3& p1, const Vec3& p2)
+{
+    return Length(p1-p2);
+}
 
 } // ns
 
diff --git a/modules/geom/src/vecmat4_op.cc b/modules/geom/src/vecmat4_op.cc
index a677c7882..456757541 100644
--- a/modules/geom/src/vecmat4_op.cc
+++ b/modules/geom/src/vecmat4_op.cc
@@ -23,93 +23,9 @@
 #include "vecmat3_op.hh"
 
 #include "exc.hh"
-#include "vec4.hh"
-#include "mat4.hh"
-#include "mat3.hh"
 
-namespace geom {
-
-Real Length(const Vec4& v)
-{
-  return std::sqrt(Length2(v));
-}
-
-Real Length2(const Vec4& v)
-{
-  return v[0]*v[0]+v[1]*v[1]+v[2]*v[2]+v[3]*v[3];
-}
-
-bool Equal(const Vec4& v1, const Vec4& v2, Real epsilon)
-{
-  return std::fabs(v1[0]-v2[0])<epsilon &&
-    std::fabs(v1[1]-v2[1])<epsilon &&
-    std::fabs(v1[2]-v2[2])<epsilon &&
-    std::fabs(v1[3]-v2[3])<epsilon;
-}
-
-bool Equal(const Mat4& m1, const Mat4& m2, Real epsilon)
-{
-  return std::fabs(m1(0,0)-m2(0,0))<epsilon &&
-    std::fabs(m1(0,1)-m2(0,1))<epsilon &&
-    std::fabs(m1(0,2)-m2(0,2))<epsilon &&
-    std::fabs(m1(0,3)-m2(0,3))<epsilon &&
-    std::fabs(m1(1,0)-m2(1,0))<epsilon &&
-    std::fabs(m1(1,1)-m2(1,1))<epsilon &&
-    std::fabs(m1(1,2)-m2(1,2))<epsilon &&
-    std::fabs(m1(1,3)-m2(1,3))<epsilon &&
-    std::fabs(m1(2,0)-m2(2,0))<epsilon &&
-    std::fabs(m1(2,1)-m2(2,1))<epsilon &&
-    std::fabs(m1(2,2)-m2(2,2))<epsilon &&
-    std::fabs(m1(2,3)-m2(2,3))<epsilon &&
-    std::fabs(m1(3,0)-m2(3,0))<epsilon &&
-    std::fabs(m1(3,1)-m2(3,1))<epsilon &&
-    std::fabs(m1(3,2)-m2(3,2))<epsilon &&
-    std::fabs(m1(3,3)-m2(3,3))<epsilon;
-}
-
-Real Dot(const Vec4& v1, const Vec4& v2)
-{
-  return v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2]+v1[3]*v2[3];
-}
-
-Vec4 Normalize(const Vec4& v)
-{
-  Real l=Length(v);
-  if(l==0.0) {
-    return v;
-  }
-  return v/l;
-}
-
-Vec4 CompMultiply(const Vec4& v1, const Vec4& v2)
-{
-  Vec4 nrvo(v1[0]*v2[0],v1[1]*v2[1],v1[2]*v2[2],v1[3]*v2[3]);
-  return nrvo;
-}
-
-Vec4 CompDivide(const Vec4& v1, const Vec4& v2)
-{
-  Vec4 nrvo(v1[0]/v2[0],v1[1]/v2[1],v1[2]/v2[2],v1[3]/v2[3]);
-  return nrvo;
-}
 
-Vec4 operator*(const Vec4& v,const Mat4& m)
-{
-  Vec4 nrvo(v[0]*m(0,0)+v[1]*m(1,0)+v[2]*m(2,0)+v[3]*m(3,0),
-      v[0]*m(0,1)+v[1]*m(1,1)+v[2]*m(2,1)+v[3]*m(3,1),
-      v[0]*m(0,2)+v[1]*m(1,2)+v[2]*m(2,2)+v[3]*m(3,2),
-      v[0]*m(0,3)+v[1]*m(1,3)+v[2]*m(2,3)+v[3]*m(3,3));
-  return nrvo;
-}
-
-Vec4 operator*(const Mat4& m, const Vec4& v)
-{
-  Vec4 nrvo(v[0]*m(0,0)+v[1]*m(0,1)+v[2]*m(0,2)+v[3]*m(0,3),
-      v[0]*m(1,0)+v[1]*m(1,1)+v[2]*m(1,2)+v[3]*m(1,3),
-      v[0]*m(2,0)+v[1]*m(2,1)+v[2]*m(2,2)+v[3]*m(2,3),
-      v[0]*m(3,0)+v[1]*m(3,1)+v[2]*m(3,2)+v[3]*m(3,3));
-  return nrvo;
-}
+namespace geom {
 
 Real Comp(const Mat4& m, unsigned int i_in, unsigned int j_in)
 {
@@ -158,7 +74,6 @@ Mat4 Transpose(const Mat4& m)
   return r;
 }
 
-#if 1
 Mat4 Invert(const Mat4& in)
 {
   Mat4 m(in);
@@ -221,22 +136,6 @@ Mat4 Invert(const Mat4& in)
   return m;
 }
 
-#else
-
-Mat4 Invert(const Mat4& m)
-{
-  Mat4 r;
-  Real d=Det(m);
-  for (int i=0;i<4;i++){
-    for (int j=0;j<4;j++){
-      r(j,i)=Comp(m,i,j)/d;
-    }
-  }
-  return r;
-}
-
-#endif
-
 Mat4 operator*(const Mat4& m1, const Mat4& m2)
 {
   Mat4 nrvo(m1(0,0)*m2(0,0)+m1(0,1)*m2(1,0)+m1(0,2)*m2(2,0)+m1(0,3)*m2(3,0),
@@ -267,22 +166,4 @@ Real Angle(const Vec4& v1, const Vec4& v2)
   return std::acos(t);
 }
 
-Vec4 Min(const Vec4& v1, const Vec4& v2)
-{
-  Vec4 nrvo(std::min(v1[0],v2[0]),
-            std::min(v1[1],v2[1]),
-            std::min(v1[2],v2[2]),
-            std::min(v1[3],v2[3]));
-  return nrvo;
-}
-
-Vec4 Max(const Vec4& v1, const Vec4& v2)
-{
-  Vec4 nrvo(std::max(v1[0],v2[0]),
-            std::max(v1[1],v2[1]),
-            std::max(v1[2],v2[2]),
-            std::max(v1[3],v2[3]));
-  return nrvo;
-}
-
 } // ns
diff --git a/modules/geom/src/vecmat4_op.hh b/modules/geom/src/vecmat4_op.hh
index 501a09ebb..7f7d6a2a9 100644
--- a/modules/geom/src/vecmat4_op.hh
+++ b/modules/geom/src/vecmat4_op.hh
@@ -22,39 +22,104 @@
 #include "constants.hh"
 
 #include <ost/geom/module_config.hh>
+#include <ost/geom/vec4.hh>
+#include <ost/geom/vec3.hh>
+#include <ost/geom/mat4.hh>
+#include <ost/geom/mat3.hh>
 
 namespace geom {
 
-class Vec4;
-class Mat4;
-
-//! returns length of vector
-Real DLLEXPORT_OST_GEOM Length(const Vec4& v);
 
 //! returns squared length of vector
-Real DLLEXPORT_OST_GEOM Length2(const Vec4& v);
+inline Real Length2(const Vec4& v)
+{
+    return v[0]*v[0]+v[1]*v[1]+v[2]*v[2]+v[3]*v[3];
+}
+
+//! returns length of vector
+inline Real Length(const Vec4& v)
+{
+  return std::sqrt(Length2(v));
+}
 
 //! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Vec4& v1, const Vec4& v2, Real ephilon=EPSILON);
+inline bool Equal(const Vec4& v1, const Vec4& v2, Real epsilon=EPSILON)
+{
+  return std::fabs(v1[0]-v2[0])<epsilon &&
+    std::fabs(v1[1]-v2[1])<epsilon &&
+    std::fabs(v1[2]-v2[2])<epsilon &&
+    std::fabs(v1[3]-v2[3])<epsilon;
+}
 
 //! return true if components differ not more than ephilon
-bool DLLEXPORT_OST_GEOM Equal(const Mat4& m1, const Mat4& v2, Real ephilon=EPSILON);
+inline bool Equal(const Mat4& m1, const Mat4& m2, Real epsilon=EPSILON)
+{
+  return std::fabs(m1(0,0)-m2(0,0))<epsilon &&
+    std::fabs(m1(0,1)-m2(0,1))<epsilon &&
+    std::fabs(m1(0,2)-m2(0,2))<epsilon &&
+    std::fabs(m1(0,3)-m2(0,3))<epsilon &&
+    std::fabs(m1(1,0)-m2(1,0))<epsilon &&
+    std::fabs(m1(1,1)-m2(1,1))<epsilon &&
+    std::fabs(m1(1,2)-m2(1,2))<epsilon &&
+    std::fabs(m1(1,3)-m2(1,3))<epsilon &&
+    std::fabs(m1(2,0)-m2(2,0))<epsilon &&
+    std::fabs(m1(2,1)-m2(2,1))<epsilon &&
+    std::fabs(m1(2,2)-m2(2,2))<epsilon &&
+    std::fabs(m1(2,3)-m2(2,3))<epsilon &&
+    std::fabs(m1(3,0)-m2(3,0))<epsilon &&
+    std::fabs(m1(3,1)-m2(3,1))<epsilon &&
+    std::fabs(m1(3,2)-m2(3,2))<epsilon &&
+    std::fabs(m1(3,3)-m2(3,3))<epsilon;
+}
 
 //! vector dot product
-Real DLLEXPORT_OST_GEOM Dot(const Vec4& v1, const Vec4& v2);
-
-Vec4 DLLEXPORT_OST_GEOM Normalize(const Vec4& v);
+inline Real Dot(const Vec4& v1, const Vec4& v2)
+{
+    return v1[0]*v2[0]+v1[1]*v2[1]+v1[2]*v2[2]+v1[3]*v2[3];
+}
+
+inline Vec4 Normalize(const Vec4& v)
+{
+  Real l=Length(v);
+  if(l==0.0) {
+    return v;
+  }
+  return v/l;
+}
 
 //! multiply each component of v1 with that of v2
-Vec4 DLLEXPORT_OST_GEOM CompMultiply(const Vec4& v1, const Vec4& v2);
+inline Vec4 CompMultiply(const Vec4& v1, const Vec4& v2)
+{
+  Vec4 nrvo(v1[0]*v2[0],v1[1]*v2[1],v1[2]*v2[2],v1[3]*v2[3]);
+  return nrvo;
+}
 
 //! divide each component of v1 by that of v2
-Vec4 DLLEXPORT_OST_GEOM CompDivide(const Vec4& v1, const Vec4& v2);
+inline Vec4 CompDivide(const Vec4& v1, const Vec4& v2)
+{
+  Vec4 nrvo(v1[0]/v2[0],v1[1]/v2[1],v1[2]/v2[2],v1[3]/v2[3]);
+  return nrvo;
+}
 
 //! vector matrix multiplication
-Vec4 DLLEXPORT_OST_GEOM operator*(const Vec4& v,const Mat4& m);
+inline Vec4 operator*(const Vec4& v,const Mat4& m)
+{
+  Vec4 nrvo(v[0]*m(0,0)+v[1]*m(1,0)+v[2]*m(2,0)+v[3]*m(3,0),
+      v[0]*m(0,1)+v[1]*m(1,1)+v[2]*m(2,1)+v[3]*m(3,1),
+      v[0]*m(0,2)+v[1]*m(1,2)+v[2]*m(2,2)+v[3]*m(3,2),
+      v[0]*m(0,3)+v[1]*m(1,3)+v[2]*m(2,3)+v[3]*m(3,3));
+  return nrvo;
+}
+
 //! vector matrix multiplication
-Vec4 DLLEXPORT_OST_GEOM operator*(const Mat4& m, const Vec4& v);
+inline Vec4 operator*(const Mat4& m, const Vec4& v)
+{
+  Vec4 nrvo(v[0]*m(0,0)+v[1]*m(0,1)+v[2]*m(0,2)+v[3]*m(0,3),
+      v[0]*m(1,0)+v[1]*m(1,1)+v[2]*m(1,2)+v[3]*m(1,3),
+      v[0]*m(2,0)+v[1]*m(2,1)+v[2]*m(2,2)+v[3]*m(2,3),
+      v[0]*m(3,0)+v[1]*m(3,1)+v[2]*m(3,2)+v[3]*m(3,3));
+  return nrvo;
+}
 
 // algebraic complement
 Real DLLEXPORT_OST_GEOM Comp(const Mat4& m, unsigned int i, unsigned int j);
@@ -69,10 +134,24 @@ Mat4 DLLEXPORT_OST_GEOM operator*(const Mat4& m1, const Mat4& m2);
 Real DLLEXPORT_OST_GEOM Angle(const Vec4& v1, const Vec4& v2);
 
 //! returns std::min of each component
-Vec4 DLLEXPORT_OST_GEOM Min(const Vec4& v1, const Vec4& v2);
+inline Vec4 Min(const Vec4& v1, const Vec4& v2)
+{
+  Vec4 nrvo(std::min(v1[0],v2[0]),
+            std::min(v1[1],v2[1]),
+            std::min(v1[2],v2[2]),
+            std::min(v1[3],v2[3]));
+  return nrvo;
+}
 
 //! returns std::max of each component
-Vec4 DLLEXPORT_OST_GEOM Max(const Vec4& v1, const Vec4& v2);
+inline Vec4  Max(const Vec4& v1, const Vec4& v2)
+{
+  Vec4 nrvo(std::max(v1[0],v2[0]),
+            std::max(v1[1],v2[1]),
+            std::max(v1[2],v2[2]),
+            std::max(v1[3],v2[3]));
+  return nrvo;
+}
 
 } // ns
 
diff --git a/modules/gfx/src/impl/simple_renderer.cc b/modules/gfx/src/impl/simple_renderer.cc
index 4e228f0dc..bfb30d313 100644
--- a/modules/gfx/src/impl/simple_renderer.cc
+++ b/modules/gfx/src/impl/simple_renderer.cc
@@ -99,7 +99,7 @@ geom::Vec3 SimpleRenderer::GetDoubleBondPlane(mol::BondHandle b)
   GetBondPartnerNormal(vec, n, bond_vec, partner2, a1, a2);
   if(n==0) {
     //set vec into x direction if no partners are present
-    vec.SetX(1);
+    vec.x=1;
   }
   vec=geom::Normalize(geom::Cross(vec,bond_vec));
   return vec;
@@ -159,7 +159,7 @@ void SimpleRenderer::PrepareRendering(GfxView& view, IndexedVertexArray& va)
         if(it->bond.GetBondOrder()==2) {
           vec=GetDoubleBondPlane(it->bond);
         } else {
-          vec.SetX(1);
+          vec.x=1;
         }
         Real dbd=options_->GetBondOrderDistance();
         VertexID id00 = va.Add(p0+dbd*vec,geom::Vec3(),it->atom1->color);
diff --git a/modules/gui/src/scene_selection.cc b/modules/gui/src/scene_selection.cc
index 72e8f4f55..b408907f4 100644
--- a/modules/gui/src/scene_selection.cc
+++ b/modules/gui/src/scene_selection.cc
@@ -110,9 +110,9 @@ void SceneSelection::CenterOnObjects() {
     }
   }
   if(changed){
-    gfx::Scene::Instance().SetCenter(geom::Vec3(((max.GetX()+min.GetX())/2),
-                                                ((max.GetY()+min.GetY())/2),
-                                                ((max.GetZ()+min.GetZ())/2)));
+    gfx::Scene::Instance().SetCenter(geom::Vec3(((max.x+min.x)/2),
+                                                ((max.y+min.y)/2),
+                                                ((max.z+min.z)/2)));
   }
 }
 
diff --git a/modules/io/src/mol/dcd_io.cc b/modules/io/src/mol/dcd_io.cc
index 0ac5a7456..29cbb7bd2 100644
--- a/modules/io/src/mol/dcd_io.cc
+++ b/modules/io/src/mol/dcd_io.cc
@@ -163,7 +163,7 @@ mol::CoordGroupHandle load_dcd(const mol::AtomHandleList& alist2,
     if(swap_flag) swap_float(&xlist[0],xlist.size());
     if(gap_flag) ff.read(dummy,sizeof(dummy));
     for(uint j=0;j<clist.size();++j) {
-      clist[j].SetX(xlist[j]);
+      clist[j].x=xlist[j];
     }
 
     // y coord
@@ -172,7 +172,7 @@ mol::CoordGroupHandle load_dcd(const mol::AtomHandleList& alist2,
     if(swap_flag) swap_float(&xlist[0],xlist.size());
     if(gap_flag) ff.read(dummy,sizeof(dummy));
     for(uint j=0;j<clist.size();++j) {
-      clist[j].SetY(xlist[j]);
+      clist[j].y=xlist[j];
     }
 
     // z coord
@@ -181,7 +181,7 @@ mol::CoordGroupHandle load_dcd(const mol::AtomHandleList& alist2,
     if(swap_flag) swap_float(&xlist[0],xlist.size());
     if(gap_flag) ff.read(dummy,sizeof(dummy));
     for(uint j=0;j<clist.size();++j) {
-      clist[j].SetZ(xlist[j]);
+      clist[j].z=xlist[j];
     }
 
     cg.AddFrame(clist);
diff --git a/modules/io/src/mol/entity_io_crd_handler.cc b/modules/io/src/mol/entity_io_crd_handler.cc
index 2fd39293f..1b093be94 100644
--- a/modules/io/src/mol/entity_io_crd_handler.cc
+++ b/modules/io/src/mol/entity_io_crd_handler.cc
@@ -186,9 +186,9 @@ bool CRDWriter::VisitAtom(const mol::AtomHandle& atom)
               << format("%5i") % res.GetNumber() << " "
               << format("%4s") % res.GetKey() << " "
               << format("%-4s") % atom.GetName()
-              << format("%10.5f") % atom.GetPos().GetX()
-              << format("%10.5f") % atom.GetPos().GetY()
-              << format("%10.5f") % atom.GetPos().GetZ() << " "
+              << format("%10.5f") % atom.GetPos().x
+              << format("%10.5f") % atom.GetPos().y
+              << format("%10.5f") % atom.GetPos().z << " "
               << format("%-4s") % e_name << " "
               << format("%-5i") % res.GetNumber() << " "
               << format("%8.5f") % atom.GetBFactor()
-- 
GitLab