diff --git a/Readme.md b/Readme.md
index 0e40590bf99013515b2cfbd96ff24adfd8f45607..9a747c70fee7591ca92c1b485b0ca5b392f76131 100644
--- a/Readme.md
+++ b/Readme.md
@@ -8,10 +8,14 @@ Best used with glm and glow.
 
 * Properties
 * Algorithms
-* Loader/Writer
 * Tests
 * std::less and std::hash for _index (and maybe _handle)
 * attribute transformations (also between different types)
 * lambda to attribute (from attribute to attribute or from make_attribute to attribute)
 * Debug: store compactify generation in handles to check for invalidation
-* Debug: insert is_removed assertions into handle access
\ No newline at end of file
+* Debug: insert is_removed assertions into handle access
+* Switch primitives and valid_primitives, check if compact flag is inlined
+* Test self-adjacent faces
+* smart ranges: average, min, max, any, all, first, last ...
+* mid-level topo API: edge-rotate-next/prev, edge-split, edge-collapse, halfedge-split, halfedge-collapse, vertex-collapse
+* annotate property preservation for mid-level topo API
diff --git a/src/polymesh/Mesh.hh b/src/polymesh/Mesh.hh
index 6f8ce8ad830b2c2960423d2e37df4a6799d02946..2505ff1519250300511b4935ddf98a00d464cf56 100644
--- a/src/polymesh/Mesh.hh
+++ b/src/polymesh/Mesh.hh
@@ -8,6 +8,9 @@
 #include "cursors.hh"
 #include "ranges.hh"
 
+// often used helper
+#include "attribute_collection.hh"
+
 namespace polymesh
 {
 using SharedMesh = std::shared_ptr<Mesh>;
@@ -77,6 +80,9 @@ public:
     vertex_handle operator[](vertex_index idx) const { return handle_of(idx); }
     halfedge_handle operator[](halfedge_index idx) const { return handle_of(idx); }
 
+    /// deletes all faces, vertices, edges, and halfedges
+    void clear();
+
     // helper
 public:
     /// Returns true if the mesh is guaranteed compact, otherwise call compactify() to be sure
@@ -104,10 +110,10 @@ public:
     // internal helper
 private:
     // reserves a certain number of primitives
-    void reserve_faces(size_t capacity);
-    void reserve_vertices(size_t capacity);
-    void reserve_edges(size_t capacity);
-    void reserve_halfedges(size_t capacity);
+    void reserve_faces(int capacity);
+    void reserve_vertices(int capacity);
+    void reserve_edges(int capacity);
+    void reserve_halfedges(int capacity);
 
     int size_faces() const { return (int)mFaces.size(); }
     int size_vertices() const { return (int)mVertices.size(); }
@@ -159,10 +165,10 @@ private:
     /// Adds a face consisting of N vertices
     /// The vertices must already be sorted in CCW order
     /// (note: trying to add already existing halfedges triggers assertions)
-    face_index add_face(vertex_handle const *v_handles, size_t vcnt);
-    face_index add_face(vertex_index const *v_indices, size_t vcnt);
-    face_index add_face(halfedge_handle const *half_loop, size_t vcnt);
-    face_index add_face(halfedge_index const *half_loop, size_t vcnt);
+    face_index add_face(vertex_handle const *v_handles, int vcnt);
+    face_index add_face(vertex_index const *v_indices, int vcnt);
+    face_index add_face(halfedge_handle const *half_loop, int vcnt);
+    face_index add_face(halfedge_index const *half_loop, int vcnt);
 
     /// Adds an edge between two existing, distinct vertices
     /// if edge already exists, returns it
@@ -389,38 +395,38 @@ inline vertex_index Mesh::add_vertex()
     return idx;
 }
 
-inline face_index Mesh::add_face(const vertex_handle *v_handles, size_t vcnt)
+inline face_index Mesh::add_face(const vertex_handle *v_handles, int vcnt)
 {
     mFaceInsertCache.resize(vcnt);
-    for (auto i = 0u; i < vcnt; ++i)
+    for (auto i = 0; i < vcnt; ++i)
         mFaceInsertCache[i] = add_or_get_halfedge(v_handles[i].idx, v_handles[(i + 1) % vcnt].idx);
     return add_face(mFaceInsertCache.data(), vcnt);
 }
 
-inline face_index Mesh::add_face(const vertex_index *v_indices, size_t vcnt)
+inline face_index Mesh::add_face(const vertex_index *v_indices, int vcnt)
 {
     mFaceInsertCache.resize(vcnt);
-    for (auto i = 0u; i < vcnt; ++i)
+    for (auto i = 0; i < vcnt; ++i)
         mFaceInsertCache[i] = add_or_get_halfedge(v_indices[i], v_indices[(i + 1) % vcnt]);
     return add_face(mFaceInsertCache.data(), vcnt);
 }
 
-inline face_index Mesh::add_face(const halfedge_handle *half_loop, size_t vcnt)
+inline face_index Mesh::add_face(const halfedge_handle *half_loop, int vcnt)
 {
     mFaceInsertCache.resize(vcnt);
-    for (auto i = 0u; i < vcnt; ++i)
+    for (auto i = 0; i < vcnt; ++i)
         mFaceInsertCache[i] = half_loop[i].idx;
     return add_face(mFaceInsertCache.data(), vcnt);
 }
 
-inline face_index Mesh::add_face(const halfedge_index *half_loop, size_t vcnt)
+inline face_index Mesh::add_face(const halfedge_index *half_loop, int vcnt)
 {
     assert(vcnt >= 3 && "no support for less-than-triangular faces");
 
     auto fidx = face_index((int)mFaces.size());
 
     // ensure that half-edges are adjacent at each vertex
-    for (auto i = 0u; i < vcnt; ++i)
+    for (auto i = 0; i < vcnt; ++i)
     {
         auto h0 = half_loop[i];
         auto h1 = half_loop[(i + 1) % vcnt];
@@ -438,7 +444,7 @@ inline face_index Mesh::add_face(const halfedge_index *half_loop, size_t vcnt)
     }
 
     // fix boundary states
-    for (auto i = 0u; i < vcnt; ++i)
+    for (auto i = 0; i < vcnt; ++i)
     {
         auto h = half_loop[i];
         auto v = halfedge(h).to_vertex;
@@ -1015,12 +1021,12 @@ inline void face_collection::reserve(int capacity) const
     mesh->reserve_faces(capacity);
 }
 
-inline face_handle face_collection::add(const vertex_handle *v_handles, size_t vcnt) const
+inline face_handle face_collection::add(const vertex_handle *v_handles, int vcnt) const
 {
     return mesh->handle_of(mesh->add_face(v_handles, vcnt));
 }
 
-inline face_handle face_collection::add(const halfedge_handle *half_loop, size_t vcnt) const
+inline face_handle face_collection::add(const halfedge_handle *half_loop, int vcnt) const
 {
     return mesh->handle_of(mesh->add_face(half_loop, vcnt));
 }
@@ -1072,7 +1078,7 @@ template <size_t N>
 inline face_handle face_collection::add(const vertex_handle (&v_handles)[N]) const
 {
     halfedge_index hs[N];
-    for (auto i = 0u; i < N; ++i)
+    for (auto i = 0; i < N; ++i)
         hs[i] = mesh->find_halfedge(v_handles[i].idx, v_handles[(i + 1) % N].idx);
     return mesh->handle_of(mesh->add_face(hs, N));
 }
@@ -1081,7 +1087,7 @@ template <size_t N>
 inline face_handle face_collection::add(const halfedge_handle (&half_loop)[N]) const
 {
     halfedge_index hs[N];
-    for (auto i = 0u; i < N; ++i)
+    for (auto i = 0; i < N; ++i)
         hs[i] = half_loop[i].idx;
     return mesh->handle_of(mesh->add_face(hs, N));
 }
@@ -1375,7 +1381,22 @@ inline void Mesh::compactify()
     mCompact = true;
 }
 
-inline void Mesh::reserve_faces(size_t capacity)
+/// ======== OTHER IMPLEMENTATION ========
+
+inline void Mesh::clear()
+{
+    for (auto &v : mVertices)
+        v.set_removed();
+    for (auto &h : mHalfedges)
+        h.set_removed();
+    for (auto &f : mFaces)
+        f.set_removed();
+
+    mCompact = false;
+    compactify();
+}
+
+inline void Mesh::reserve_faces(int capacity)
 {
     for (auto a = mFaceAttrs; a; a = a->mNextAttribute)
         a->resize(capacity, false);
@@ -1383,7 +1404,7 @@ inline void Mesh::reserve_faces(size_t capacity)
     mFaces.reserve(capacity);
 }
 
-inline void Mesh::reserve_vertices(size_t capacity)
+inline void Mesh::reserve_vertices(int capacity)
 {
     for (auto a = mVertexAttrs; a; a = a->mNextAttribute)
         a->resize(capacity, false);
@@ -1391,7 +1412,7 @@ inline void Mesh::reserve_vertices(size_t capacity)
     mVertices.reserve(capacity);
 }
 
-inline void Mesh::reserve_edges(size_t capacity)
+inline void Mesh::reserve_edges(int capacity)
 {
     for (auto a = mEdgeAttrs; a; a = a->mNextAttribute)
         a->resize(capacity, false);
@@ -1401,7 +1422,7 @@ inline void Mesh::reserve_edges(size_t capacity)
     mHalfedges.reserve(capacity * 2);
 }
 
-inline void Mesh::reserve_halfedges(size_t capacity)
+inline void Mesh::reserve_halfedges(int capacity)
 {
     for (auto a = mHalfedgeAttrs; a; a = a->mNextAttribute)
         a->resize(capacity, false);
@@ -1624,49 +1645,49 @@ inline vertex_vertex_ring vertex_handle::adjacent_vertices() const
 
 /// ======== attributes IMPLEMENTATION ========
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT> vertex_collection::make_attribute(const AttrT &def_value)
 {
     return vertex_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT> const_vertex_collection::make_attribute(const AttrT &def_value)
 {
     return vertex_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT> face_collection::make_attribute(const AttrT &def_value)
 {
     return face_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT> const_face_collection::make_attribute(const AttrT &def_value)
 {
     return face_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT> edge_collection::make_attribute(const AttrT &def_value)
 {
     return edge_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT> const_edge_collection::make_attribute(const AttrT &def_value)
 {
     return edge_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT> halfedge_collection::make_attribute(const AttrT &def_value)
 {
     return halfedge_attribute<AttrT>(mesh, def_value);
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT> const_halfedge_collection::make_attribute(const AttrT &def_value)
 {
     return halfedge_attribute<AttrT>(mesh, def_value);
@@ -1860,103 +1881,114 @@ inline void halfedge_attribute_base::register_attr()
     mMesh->register_attr(this);
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT>::vertex_attribute(const Mesh *mesh, const AttrT &def_value)
   : vertex_attribute_base(mesh), mDefaultValue(def_value)
 {
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT>::face_attribute(const Mesh *mesh, const AttrT &def_value)
   : face_attribute_base(mesh), mDefaultValue(def_value)
 {
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT>::edge_attribute(const Mesh *mesh, const AttrT &def_value)
   : edge_attribute_base(mesh), mDefaultValue(def_value)
 {
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT>::halfedge_attribute(const Mesh *mesh, const AttrT &def_value)
   : halfedge_attribute_base(mesh), mDefaultValue(def_value)
 {
     register_attr();
 }
 
-template <typename AttrT>
-size_t vertex_attribute<AttrT>::size() const
+template <class AttrT>
+int vertex_attribute<AttrT>::size() const
 {
     return mMesh->vertices().size();
 }
 
-template <typename AttrT>
+template <class AttrT>
 void vertex_attribute<AttrT>::clear(AttrT const &value)
 {
     mData.clear();
     mData.resize(mMesh->vertices().size(), value);
 }
-template <typename AttrT>
+template <class AttrT>
 void vertex_attribute<AttrT>::clear()
 {
     clear(mDefaultValue);
 }
 
-template <typename AttrT>
-size_t face_attribute<AttrT>::size() const
+template <class AttrT>
+int face_attribute<AttrT>::size() const
 {
     return mMesh->vertices().size();
 }
 
-template <typename AttrT>
+template <class AttrT>
 void face_attribute<AttrT>::clear(AttrT const &value)
 {
     mData.clear();
     mData.resize(mMesh->vertices().size(), value);
 }
-template <typename AttrT>
+template <class AttrT>
 void face_attribute<AttrT>::clear()
 {
     clear(mDefaultValue);
 }
 
-template <typename AttrT>
-size_t edge_attribute<AttrT>::size() const
+template <class AttrT>
+int edge_attribute<AttrT>::size() const
 {
     return mMesh->vertices().size();
 }
 
-template <typename AttrT>
+template <class AttrT>
 void edge_attribute<AttrT>::clear(AttrT const &value)
 {
     mData.clear();
     mData.resize(mMesh->vertices().size(), value);
 }
-template <typename AttrT>
+template <class AttrT>
 void edge_attribute<AttrT>::clear()
 {
     clear(mDefaultValue);
 }
 
-template <typename AttrT>
-size_t halfedge_attribute<AttrT>::size() const
+template <class AttrT>
+int halfedge_attribute<AttrT>::size() const
 {
     return mMesh->vertices().size();
 }
 
-template <typename AttrT>
+template <class AttrT>
 void halfedge_attribute<AttrT>::clear(AttrT const &value)
 {
     mData.clear();
     mData.resize(mMesh->vertices().size(), value);
 }
-template <typename AttrT>
+template <class AttrT>
 void halfedge_attribute<AttrT>::clear()
 {
     clear(mDefaultValue);
 }
+
+template <class AttrT>
+template <class FuncT>
+auto vertex_attribute<AttrT>::map(FuncT f) const -> tmp::result_type_of<FuncT, AttrT>
+{
+    auto attr = mMesh->vertices().make_attribute<tmp::result_type_of<FuncT, AttrT>>();
+    auto s = size();
+    for (auto i = 0; i < s; ++i)
+        attr.mData[i] = f(mData[i]);
+    return attr; // copy elison
+}
 }
diff --git a/src/polymesh/algorithms.hh b/src/polymesh/algorithms.hh
index 78e921081e9fc4334d1b29c183dfab66a8cd4881..3fa1f9bea2d301842589d8224c1f81eb6acd99d5 100644
--- a/src/polymesh/algorithms.hh
+++ b/src/polymesh/algorithms.hh
@@ -16,9 +16,7 @@
 
 // Basic mesh operations, including:
 // - elementary subdivision
-// - edge splits
 // - intersections
-// - collapses
 #include "algorithms/operations.hh"
 
 // TODO:
diff --git a/src/polymesh/algorithms/operations.hh b/src/polymesh/algorithms/operations.hh
index 8a83d60bf874f6fabff8305a377805f4db6f1a43..16e96f4b7c2b55245d1ccf3a942224395af10156 100644
--- a/src/polymesh/algorithms/operations.hh
+++ b/src/polymesh/algorithms/operations.hh
@@ -6,9 +6,7 @@
 
 // Basic mesh operations, including:
 // - elementary subdivision
-// - edge splits
 // - intersections
-// - collapses
 
 namespace polymesh
 {
diff --git a/src/polymesh/algorithms/properties.hh b/src/polymesh/algorithms/properties.hh
index 7a8694edf79970d49c00fa6e233a3f040fe85df6..4a7fe637daed428e97ff17d2dd884c7b20162c81 100644
--- a/src/polymesh/algorithms/properties.hh
+++ b/src/polymesh/algorithms/properties.hh
@@ -3,6 +3,7 @@
 #include <glm/glm.hpp>
 
 #include "../Mesh.hh"
+#include "../fields.hh"
 
 // Derived mesh properties, including:
 // - valences
@@ -20,16 +21,20 @@ namespace polymesh
 int valence_of(vertex_handle v);
 
 /// returns the area of the (flat) polygonal face
-float face_area(face_handle f, vertex_attribute<glm::vec3> const& position);
+template <class Vec3>
+typename field_3d<Vec3>::Scalar face_area(face_handle f, vertex_attribute<Vec3> const& position);
 
 /// returns the center of gravity for a given (flat) polygonal face
-glm::vec3 face_centroid(face_handle f, vertex_attribute<glm::vec3> const& position);
+template <class Vec3>
+Vec3 face_centroid(face_handle f, vertex_attribute<Vec3> const& position);
 
 /// returns the area of a given triangle
-float triangle_area(face_handle f, vertex_attribute<glm::vec3> const& position);
+template <class Vec3>
+typename field_3d<Vec3>::Scalar triangle_area(face_handle f, vertex_attribute<Vec3> const& position);
 
 /// returns the center of gravity for a given triangle
-glm::vec3 triangle_centroid(face_handle f, vertex_attribute<glm::vec3> const& position);
+template <class Vec3>
+Vec3 triangle_centroid(face_handle f, vertex_attribute<Vec3> const& position);
 
 /// ======== IMPLEMENTATION ========
 
@@ -38,29 +43,32 @@ inline int valence_of(vertex_handle v)
     return v.adjacent_vertices().size();
 }
 
-inline float triangle_area(face_handle f, vertex_attribute<glm::vec3> const& position)
+template <class Vec3>
+typename field_3d<Vec3>::Scalar triangle_area(face_handle f, vertex_attribute<Vec3> const& position)
 {
     auto h = f.any_halfedge();
     auto p0 = position[h.vertex_from()];
     auto p1 = position[h.vertex_to()];
     auto p2 = position[h.next().vertex_to()];
 
-    return 0.5f * length(cross(p0 - p1, p0 - p2));
+    return field_3d<Vec3>::length(field_3d<Vec3>::cross(p0 - p1, p0 - p2)) * field_3d<Vec3>::scalar(0.5f);
 }
 
-inline glm::vec3 triangle_centroid(face_handle f, vertex_attribute<glm::vec3> const& position)
+template <class Vec3>
+Vec3 triangle_centroid(face_handle f, vertex_attribute<Vec3> const& position)
 {
     auto h = f.any_halfedge();
     auto p0 = position[h.vertex_from()];
     auto p1 = position[h.vertex_to()];
     auto p2 = position[h.next().vertex_to()];
 
-    return (p0 + p1 + p2) / 3.0f;
+    return (p0 + p1 + p2) / field_3d<Vec3>::scalar(3);
 }
 
-inline float face_area(face_handle f, vertex_attribute<glm::vec3> const& position)
+template <class Vec3>
+typename field_3d<Vec3>::Scalar face_area(face_handle f, vertex_attribute<Vec3> const& position)
 {
-    glm::vec3 varea;
+    auto varea = field_3d<Vec3>::zero();
 
     auto h = f.any_halfedge();
 
@@ -74,22 +82,23 @@ inline float face_area(face_handle f, vertex_attribute<glm::vec3> const& positio
     {
         auto p_curr = h.vertex_to()[position];
 
-        varea += cross(p_prev - p0, p_curr - p0);
+        varea += field_3d<Vec3>::cross(p_prev - p0, p_curr - p0);
 
         // circulate
         h = h.next();
         p_prev = p_curr;
     } while (h.vertex_to() != v0);
 
-    return length(varea) * 0.5f;
+    return field_3d<Vec3>::length(varea) * 0.5f;
 }
 
-inline glm::vec3 face_centroid(face_handle f, vertex_attribute<glm::vec3> const& position)
+template <class Vec3>
+Vec3 face_centroid(face_handle f, vertex_attribute<Vec3> const& position)
 {
     // TODO: make correct for non-convex polygons!
 
-    float area = 0.0f;
-    glm::vec3 centroid;
+    auto area = field_3d<Vec3>::scalar(0);
+    auto centroid = field_3d<Vec3>::zero();
 
     auto h = f.any_halfedge();
 
@@ -103,7 +112,7 @@ inline glm::vec3 face_centroid(face_handle f, vertex_attribute<glm::vec3> const&
     {
         auto p_curr = h.vertex_to()[position];
 
-        auto a = length(cross(p_prev - p0, p_curr - p0));
+        auto a = field_3d<Vec3>::length(field_3d<Vec3>::cross(p_prev - p0, p_curr - p0));
         area += a;
         centroid += (p_prev + p_curr + p0) * a;
 
diff --git a/src/polymesh/attribute_base.hh b/src/polymesh/attribute_base.hh
index 00f62bb5dbe6fd417078181721b826a3cfad3438..d3052897cabeb19008ba5571780bfe56e6d56e7f 100644
--- a/src/polymesh/attribute_base.hh
+++ b/src/polymesh/attribute_base.hh
@@ -9,10 +9,10 @@ namespace polymesh
 {
 class Mesh;
 
-template <typename DataT>
+template <class DataT>
 struct attribute_data
 {
-    size_t size = 0;
+    int size = 0;
     DataT* data = nullptr;
 
     DataT& operator[](int i) { return data[i]; }
@@ -24,7 +24,7 @@ struct attribute_data
         size = rhs.size;
         data = new DataT[size];
 
-        for (size_t i = 0; i < size; ++i)
+        for (int i = 0; i < size; ++i)
             data[i] = rhs.data[i];
     }
     attribute_data(attribute_data<DataT>&& rhs) // move
@@ -42,7 +42,7 @@ struct attribute_data
         size = rhs.size;
         data = new DataT[size];
 
-        for (size_t i = 0; i < size; ++i)
+        for (int i = 0; i < size; ++i)
             data[i] = rhs.data[i];
 
         return *this;
@@ -61,21 +61,21 @@ struct attribute_data
     }
     ~attribute_data() { delete[] data; }
 
-    void resize(size_t new_size, DataT const& default_value)
+    void resize(int new_size, DataT const& default_value)
     {
         auto new_data = new DataT[new_size];
 
         if (new_size < size)
         {
-            for (size_t i = 0; i < new_size; ++i)
+            for (int i = 0; i < new_size; ++i)
                 new_data[i] = data[i];
         }
         else
         {
-            for (size_t i = 0; i < size; ++i)
+            for (int i = 0; i < size; ++i)
                 new_data[i] = data[i];
 
-            for (size_t i = size; i < new_size; ++i)
+            for (int i = size; i < new_size; ++i)
                 new_data[i] = default_value;
         }
 
@@ -91,7 +91,7 @@ private:
     vertex_attribute_base* mNextAttribute = nullptr;
     vertex_attribute_base* mPrevAttribute = nullptr;
 
-    void resize(size_t newSize, bool force)
+    void resize(int newSize, bool force)
     {
         if (force)
         {
@@ -108,11 +108,11 @@ private:
     }
 
 protected:
-    size_t mDataSize = 0;
+    int mDataSize = 0;
     Mesh const* mMesh;
     vertex_attribute_base(Mesh const* mesh);
     virtual ~vertex_attribute_base() { deregister_attr(); }
-    virtual void on_resize(size_t newSize) = 0;
+    virtual void on_resize(int newSize) = 0;
     virtual void apply_remapping(std::vector<int> const& map) = 0;
     void register_attr();
     void deregister_attr();
@@ -125,7 +125,7 @@ private:
     face_attribute_base* mNextAttribute = nullptr;
     face_attribute_base* mPrevAttribute = nullptr;
 
-    void resize(size_t newSize, bool force)
+    void resize(int newSize, bool force)
     {
         if (force)
         {
@@ -142,11 +142,11 @@ private:
     }
 
 protected:
-    size_t mDataSize = 0;
+    int mDataSize = 0;
     Mesh const* mMesh;
     face_attribute_base(Mesh const* mesh);
     virtual ~face_attribute_base() { deregister_attr(); }
-    virtual void on_resize(size_t newSize) = 0;
+    virtual void on_resize(int newSize) = 0;
     virtual void apply_remapping(std::vector<int> const& map) = 0;
     void register_attr();
     void deregister_attr();
@@ -159,7 +159,7 @@ private:
     edge_attribute_base* mNextAttribute = nullptr;
     edge_attribute_base* mPrevAttribute = nullptr;
 
-    void resize(size_t newSize, bool force)
+    void resize(int newSize, bool force)
     {
         if (force)
         {
@@ -176,11 +176,11 @@ private:
     }
 
 protected:
-    size_t mDataSize = 0;
+    int mDataSize = 0;
     Mesh const* mMesh;
     edge_attribute_base(Mesh const* mesh);
     virtual ~edge_attribute_base() { deregister_attr(); }
-    virtual void on_resize(size_t newSize) = 0;
+    virtual void on_resize(int newSize) = 0;
     virtual void apply_remapping(std::vector<int> const& map) = 0;
     void register_attr();
     void deregister_attr();
@@ -193,7 +193,7 @@ private:
     halfedge_attribute_base* mNextAttribute = nullptr;
     halfedge_attribute_base* mPrevAttribute = nullptr;
 
-    void resize(size_t newSize, bool force)
+    void resize(int newSize, bool force)
     {
         if (force)
         {
@@ -210,11 +210,11 @@ private:
     }
 
 protected:
-    size_t mDataSize = 0;
+    int mDataSize = 0;
     Mesh const* mMesh;
     halfedge_attribute_base(Mesh const* mesh);
     virtual ~halfedge_attribute_base() { deregister_attr(); }
-    virtual void on_resize(size_t newSize) = 0;
+    virtual void on_resize(int newSize) = 0;
     virtual void apply_remapping(std::vector<int> const& map) = 0;
     void register_attr();
     void deregister_attr();
diff --git a/src/polymesh/attribute_collection.hh b/src/polymesh/attribute_collection.hh
new file mode 100644
index 0000000000000000000000000000000000000000..6d21a7377bb1dde1d57da5446eca3fa0c7f99995
--- /dev/null
+++ b/src/polymesh/attribute_collection.hh
@@ -0,0 +1,16 @@
+#pragma once
+
+#include <map>
+#include <string>
+
+#include "attributes.hh"
+
+namespace polymesh
+{
+struct attribute_collection
+{
+private:
+    std::map<std::string, vertex_attribute_base> mVertexAttrs;
+    // TODO
+};
+}
diff --git a/src/polymesh/attributes.hh b/src/polymesh/attributes.hh
index 61d332b7f73d4fe33d6e95fdbaa6060fb78e6e2b..36ef96c43ccfd8563f3d28e2122a4c2670d9fac8 100644
--- a/src/polymesh/attributes.hh
+++ b/src/polymesh/attributes.hh
@@ -5,6 +5,7 @@
 
 #include "attribute_base.hh"
 #include "cursors.hh"
+#include "tmp.hh"
 
 /** Attributes
  *
@@ -22,7 +23,7 @@
 
 namespace polymesh
 {
-template <typename AttrT>
+template <class AttrT>
 struct vertex_attribute : vertex_attribute_base
 {
     // data access
@@ -34,19 +35,26 @@ public:
 
     AttrT* data() { return mData; }
     AttrT const* data() const { return mData; }
-    size_t size() const;
+    int size() const;
 
     // methods
 public:
     void clear(AttrT const& value);
     void clear();
 
+    /// returns a new attribute where the given function was applied to each entry
+    template <class FuncT>
+    auto map(FuncT f) const -> tmp::result_type_of<FuncT, AttrT>;
+    /// applies to given function to each attribute entry
+    template <class FuncT>
+    void apply(FuncT f);
+
     // data
 private:
     attribute_data<AttrT> mData;
     AttrT mDefaultValue;
 
-    void on_resize(size_t newSize) override { mData.resize(newSize, mDefaultValue); }
+    void on_resize(int newSize) override { mData.resize(newSize, mDefaultValue); }
     void apply_remapping(std::vector<int> const& map) override;
 
     // ctor
@@ -62,7 +70,7 @@ public:
     vertex_attribute& operator=(vertex_attribute&&);
 };
 
-template <typename AttrT>
+template <class AttrT>
 struct face_attribute : face_attribute_base
 {
     // data access
@@ -74,7 +82,7 @@ public:
 
     AttrT* data() { return mData.data(); }
     AttrT const* data() const { return mData.data(); }
-    size_t size() const;
+    int size() const;
 
     // methods
 public:
@@ -86,7 +94,7 @@ private:
     attribute_data<AttrT> mData;
     AttrT mDefaultValue;
 
-    void on_resize(size_t newSize) override { mData.resize(newSize, mDefaultValue); }
+    void on_resize(int newSize) override { mData.resize(newSize, mDefaultValue); }
     void apply_remapping(std::vector<int> const& map) override;
 
     // ctor
@@ -102,7 +110,7 @@ public:
     face_attribute& operator=(face_attribute&&);
 };
 
-template <typename AttrT>
+template <class AttrT>
 struct edge_attribute : edge_attribute_base
 {
     // data access
@@ -114,7 +122,7 @@ public:
 
     AttrT* data() { return mData.data(); }
     AttrT const* data() const { return mData.data(); }
-    size_t size() const;
+    int size() const;
 
     // methods
 public:
@@ -126,7 +134,7 @@ private:
     attribute_data<AttrT> mData;
     AttrT mDefaultValue;
 
-    void on_resize(size_t newSize) override { mData.resize(newSize, mDefaultValue); }
+    void on_resize(int newSize) override { mData.resize(newSize, mDefaultValue); }
     void apply_remapping(std::vector<int> const& map) override;
 
     // ctor
@@ -142,7 +150,7 @@ public:
     edge_attribute& operator=(edge_attribute&&);
 };
 
-template <typename AttrT>
+template <class AttrT>
 struct halfedge_attribute : halfedge_attribute_base
 {
     // data access
@@ -154,14 +162,14 @@ public:
 
     AttrT* data() { return mData.data(); }
     AttrT const* data() const { return mData.data(); }
-    size_t size() const;
+    int size() const;
 
     // methods
 public:
     void clear(AttrT const& value);
     void clear();
 
-    void on_resize(size_t newSize) override { mData.resize(newSize, mDefaultValue); }
+    void on_resize(int newSize) override { mData.resize(newSize, mDefaultValue); }
     void apply_remapping(std::vector<int> const& map) override;
 
     // data
@@ -184,32 +192,32 @@ public:
 
 /// ======== IMPLEMENTATION ========
 
-template <typename AttrT>
+template <class AttrT>
 void vertex_attribute<AttrT>::apply_remapping(const std::vector<int>& map)
 {
     for (auto i = 0u; i < map.size(); ++i)
         mData[i] = mData[map[i]];
 }
-template <typename AttrT>
+template <class AttrT>
 void face_attribute<AttrT>::apply_remapping(const std::vector<int>& map)
 {
     for (auto i = 0u; i < map.size(); ++i)
         mData[i] = mData[map[i]];
 }
-template <typename AttrT>
+template <class AttrT>
 void edge_attribute<AttrT>::apply_remapping(const std::vector<int>& map)
 {
     for (auto i = 0u; i < map.size(); ++i)
         mData[i] = mData[map[i]];
 }
-template <typename AttrT>
+template <class AttrT>
 void halfedge_attribute<AttrT>::apply_remapping(const std::vector<int>& map)
 {
     for (auto i = 0u; i < map.size(); ++i)
         mData[i] = mData[map[i]];
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT>::vertex_attribute(vertex_attribute const& rhs) : vertex_attribute_base(rhs.mMesh) // copy
 {
     mDefaultValue = rhs.mDefaultValue;
@@ -219,7 +227,7 @@ vertex_attribute<AttrT>::vertex_attribute(vertex_attribute const& rhs) : vertex_
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT>::vertex_attribute(vertex_attribute&& rhs) : vertex_attribute_base(rhs.mMesh) // move
 {
     mDefaultValue = std::move(rhs.mDefaultValue);
@@ -230,7 +238,7 @@ vertex_attribute<AttrT>::vertex_attribute(vertex_attribute&& rhs) : vertex_attri
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT>& vertex_attribute<AttrT>::operator=(vertex_attribute const& rhs) // copy
 {
     deregister_attr();
@@ -243,7 +251,7 @@ vertex_attribute<AttrT>& vertex_attribute<AttrT>::operator=(vertex_attribute con
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 vertex_attribute<AttrT>& vertex_attribute<AttrT>::operator=(vertex_attribute&& rhs) // move
 {
     deregister_attr();
@@ -257,7 +265,7 @@ vertex_attribute<AttrT>& vertex_attribute<AttrT>::operator=(vertex_attribute&& r
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT>::face_attribute(face_attribute const& rhs) : face_attribute_base(rhs.mMesh) // copy
 {
     mDefaultValue = rhs.mDefaultValue;
@@ -267,7 +275,7 @@ face_attribute<AttrT>::face_attribute(face_attribute const& rhs) : face_attribut
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT>::face_attribute(face_attribute&& rhs) : face_attribute_base(rhs.mMesh) // move
 {
     mDefaultValue = std::move(rhs.mDefaultValue);
@@ -278,7 +286,7 @@ face_attribute<AttrT>::face_attribute(face_attribute&& rhs) : face_attribute_bas
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT>& face_attribute<AttrT>::operator=(face_attribute const& rhs) // copy
 {
     deregister_attr();
@@ -291,7 +299,7 @@ face_attribute<AttrT>& face_attribute<AttrT>::operator=(face_attribute const& rh
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 face_attribute<AttrT>& face_attribute<AttrT>::operator=(face_attribute&& rhs) // move
 {
     deregister_attr();
@@ -305,7 +313,7 @@ face_attribute<AttrT>& face_attribute<AttrT>::operator=(face_attribute&& rhs) //
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT>::edge_attribute(edge_attribute const& rhs) : edge_attribute_base(rhs.mMesh) // copy
 {
     mDefaultValue = rhs.mDefaultValue;
@@ -315,7 +323,7 @@ edge_attribute<AttrT>::edge_attribute(edge_attribute const& rhs) : edge_attribut
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT>::edge_attribute(edge_attribute&& rhs) : edge_attribute_base(rhs.mMesh) // move
 {
     mDefaultValue = std::move(rhs.mDefaultValue);
@@ -326,7 +334,7 @@ edge_attribute<AttrT>::edge_attribute(edge_attribute&& rhs) : edge_attribute_bas
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT>& edge_attribute<AttrT>::operator=(edge_attribute const& rhs) // copy
 {
     deregister_attr();
@@ -339,7 +347,7 @@ edge_attribute<AttrT>& edge_attribute<AttrT>::operator=(edge_attribute const& rh
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 edge_attribute<AttrT>& edge_attribute<AttrT>::operator=(edge_attribute&& rhs) // move
 {
     deregister_attr();
@@ -353,7 +361,7 @@ edge_attribute<AttrT>& edge_attribute<AttrT>::operator=(edge_attribute&& rhs) //
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT>::halfedge_attribute(halfedge_attribute const& rhs)
   : halfedge_attribute_base(rhs.mMesh) // copy
 {
@@ -364,7 +372,7 @@ halfedge_attribute<AttrT>::halfedge_attribute(halfedge_attribute const& rhs)
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT>::halfedge_attribute(halfedge_attribute&& rhs) : halfedge_attribute_base(rhs.mMesh) // move
 {
     mDefaultValue = std::move(rhs.mDefaultValue);
@@ -375,7 +383,7 @@ halfedge_attribute<AttrT>::halfedge_attribute(halfedge_attribute&& rhs) : halfed
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT>& halfedge_attribute<AttrT>::operator=(halfedge_attribute const& rhs) // copy
 {
     deregister_attr();
@@ -388,7 +396,7 @@ halfedge_attribute<AttrT>& halfedge_attribute<AttrT>::operator=(halfedge_attribu
     register_attr();
 }
 
-template <typename AttrT>
+template <class AttrT>
 halfedge_attribute<AttrT>& halfedge_attribute<AttrT>::operator=(halfedge_attribute&& rhs) // move
 {
     deregister_attr();
@@ -404,87 +412,172 @@ halfedge_attribute<AttrT>& halfedge_attribute<AttrT>::operator=(halfedge_attribu
 
 /// ======== CURSOR IMPLEMENTATION ========
 
-template <typename AttrT>
+template <class AttrT>
 AttrT& face_index::operator[](face_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& face_index::operator[](face_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT& face_handle::operator[](face_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& face_handle::operator[](face_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
 
-template <typename AttrT>
+template <class AttrT>
 AttrT& vertex_index::operator[](vertex_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& vertex_index::operator[](vertex_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT& vertex_handle::operator[](vertex_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& vertex_handle::operator[](vertex_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
 
-template <typename AttrT>
+template <class AttrT>
 AttrT& edge_index::operator[](edge_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& edge_index::operator[](edge_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT& edge_handle::operator[](edge_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& edge_handle::operator[](edge_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
 
-template <typename AttrT>
+template <class AttrT>
 AttrT& halfedge_index::operator[](halfedge_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& halfedge_index::operator[](halfedge_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT& halfedge_handle::operator[](halfedge_attribute<AttrT>& attr) const
 {
     return attr[*this];
 }
-template <typename AttrT>
+template <class AttrT>
 AttrT const& halfedge_handle::operator[](halfedge_attribute<AttrT> const& attr) const
 {
     return attr[*this];
 }
+
+template <class AttrT>
+AttrT& face_index::operator[](face_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& face_index::operator[](face_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT& face_handle::operator[](face_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& face_handle::operator[](face_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+
+template <class AttrT>
+AttrT& vertex_index::operator[](vertex_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& vertex_index::operator[](vertex_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT& vertex_handle::operator[](vertex_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& vertex_handle::operator[](vertex_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+
+template <class AttrT>
+AttrT& edge_index::operator[](edge_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& edge_index::operator[](edge_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT& edge_handle::operator[](edge_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& edge_handle::operator[](edge_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+
+template <class AttrT>
+AttrT& halfedge_index::operator[](halfedge_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& halfedge_index::operator[](halfedge_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT& halfedge_handle::operator[](halfedge_attribute<AttrT>* attr) const
+{
+    return (*attr)[*this];
+}
+template <class AttrT>
+AttrT const& halfedge_handle::operator[](halfedge_attribute<AttrT> const* attr) const
+{
+    return (*attr)[*this];
+}
+
 }
diff --git a/src/polymesh/cursors.hh b/src/polymesh/cursors.hh
index bbd2a2157440eff8f4bea86baec128df98b73871..508109dafa59b54dd9622ce2d4958eb7c19fe40c 100644
--- a/src/polymesh/cursors.hh
+++ b/src/polymesh/cursors.hh
@@ -6,13 +6,13 @@ namespace polymesh
 {
 class Mesh;
 
-template <typename PropT>
+template <class PropT>
 struct vertex_attribute;
-template <typename PropT>
+template <class PropT>
 struct face_attribute;
-template <typename PropT>
+template <class PropT>
 struct edge_attribute;
-template <typename PropT>
+template <class PropT>
 struct halfedge_attribute;
 
 struct vertex_handle;
@@ -47,10 +47,14 @@ struct face_index
     bool operator==(face_index const& rhs) const { return value == rhs.value; }
     bool operator!=(face_index const& rhs) const { return value != rhs.value; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](face_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](face_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](face_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](face_attribute<PropT> const* prop) const;
 };
 
 struct vertex_index
@@ -67,10 +71,14 @@ struct vertex_index
     bool operator==(vertex_index const& rhs) const { return value == rhs.value; }
     bool operator!=(vertex_index const& rhs) const { return value != rhs.value; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](vertex_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](vertex_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](vertex_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](vertex_attribute<PropT> const* prop) const;
 };
 
 struct edge_index
@@ -87,10 +95,14 @@ struct edge_index
     bool operator==(edge_index const& rhs) const { return value == rhs.value; }
     bool operator!=(edge_index const& rhs) const { return value != rhs.value; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](edge_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](edge_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](edge_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](edge_attribute<PropT> const* prop) const;
 };
 
 struct halfedge_index
@@ -107,10 +119,14 @@ struct halfedge_index
     bool operator==(halfedge_index const& rhs) const { return value == rhs.value; }
     bool operator!=(halfedge_index const& rhs) const { return value != rhs.value; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](halfedge_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](halfedge_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](halfedge_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](halfedge_attribute<PropT> const* prop) const;
 };
 
 // ======================== HANDLES ========================
@@ -128,10 +144,14 @@ struct face_handle
     bool operator==(face_handle const& rhs) const { return mesh == rhs.mesh && idx == rhs.idx; }
     bool operator!=(face_handle const& rhs) const { return mesh != rhs.mesh || idx != rhs.idx; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](face_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](face_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](face_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](face_attribute<PropT> const* prop) const;
 
     bool is_valid() const { return idx.is_valid(); }    ///< valid idx (but could be deleted in some iterators)
     bool is_invalid() const { return !idx.is_valid(); } ///< invalid idx
@@ -161,10 +181,14 @@ struct vertex_handle
     bool operator==(vertex_handle const& rhs) const { return mesh == rhs.mesh && idx == rhs.idx; }
     bool operator!=(vertex_handle const& rhs) const { return mesh != rhs.mesh || idx != rhs.idx; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](vertex_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](vertex_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](vertex_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](vertex_attribute<PropT> const* prop) const;
 
     bool is_valid() const { return idx.is_valid(); }    ///< valid idx (but could be deleted in some iterators)
     bool is_invalid() const { return !idx.is_valid(); } ///< invalid idx
@@ -199,10 +223,14 @@ struct edge_handle
     bool operator==(edge_handle const& rhs) const { return mesh == rhs.mesh && idx == rhs.idx; }
     bool operator!=(edge_handle const& rhs) const { return mesh != rhs.mesh || idx != rhs.idx; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](edge_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](edge_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](edge_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](edge_attribute<PropT> const* prop) const;
 
     bool is_valid() const { return idx.is_valid(); }    ///< valid idx (but could be deleted in some iterators)
     bool is_invalid() const { return !idx.is_valid(); } ///< invalid idx
@@ -232,10 +260,14 @@ struct halfedge_handle
     bool operator==(halfedge_handle const& rhs) const { return mesh == rhs.mesh && idx == rhs.idx; }
     bool operator!=(halfedge_handle const& rhs) const { return mesh != rhs.mesh || idx != rhs.idx; }
 
-    template <typename PropT>
+    template <class PropT>
     PropT& operator[](halfedge_attribute<PropT>& prop) const;
-    template <typename PropT>
+    template <class PropT>
     PropT const& operator[](halfedge_attribute<PropT> const& prop) const;
+    template <class PropT>
+    PropT& operator[](halfedge_attribute<PropT>* prop) const;
+    template <class PropT>
+    PropT const& operator[](halfedge_attribute<PropT> const* prop) const;
 
     bool is_valid() const { return idx.is_valid(); }    ///< valid idx (but could be deleted in some iterators)
     bool is_invalid() const { return !idx.is_valid(); } ///< invalid idx
diff --git a/src/polymesh/fields.hh b/src/polymesh/fields.hh
new file mode 100644
index 0000000000000000000000000000000000000000..231845c4fe15c207b5e9a13db61197098aaad28b
--- /dev/null
+++ b/src/polymesh/fields.hh
@@ -0,0 +1,41 @@
+#pragma once
+
+#include <cmath>
+#include <cstddef>
+#include <utility>
+
+namespace polymesh
+{
+/// Type trait for 3D vector types
+template <class Vec3>
+struct field_3d
+{
+    using Point = Vec3;
+    using Scalar = typename std::decay<decltype(std::declval<Vec3>()[0])>::type;
+
+    constexpr static Point make(Scalar x, Scalar y, Scalar z) { return Point(x, y, z); }
+    constexpr static Point zero() { return make(0, 0, 0); }
+
+    constexpr static Scalar dot(Vec3 const& a, Vec3 const& b)
+    {
+        return a[0] * b[0] + //
+               a[1] * b[1] + //
+               a[2] * b[2];
+    }
+
+    constexpr static Point cross(Vec3 const& a, Vec3 const& b)
+    {
+        return Point(a[1] * b[2] - a[2] * b[1], //
+                     a[2] * b[0] - a[0] * b[2], //
+                     a[0] * b[1] - a[1] * b[0]);
+    }
+
+    constexpr static Scalar length(Vec3 const& a) { return std::sqrt(dot(a, a)); }
+
+    template <class T>
+    constexpr static Scalar scalar(T const& t)
+    {
+        return static_cast<Scalar>(t);
+    }
+};
+}
diff --git a/src/polymesh/formats/obj.cc b/src/polymesh/formats/obj.cc
new file mode 100644
index 0000000000000000000000000000000000000000..d5f470877f325e9f36212c117bb9b5b292c925eb
--- /dev/null
+++ b/src/polymesh/formats/obj.cc
@@ -0,0 +1,296 @@
+#include "obj.hh"
+
+#include <fstream>
+#include <sstream>
+
+using namespace polymesh;
+
+obj_writer::obj_writer(const std::string &filename)
+{
+    tmp_out = new std::ofstream(filename);
+    out = tmp_out;
+}
+
+obj_writer::obj_writer(std::ostream &out)
+{
+    this->out = &out;
+}
+
+obj_writer::~obj_writer()
+{
+    delete tmp_out;
+}
+
+void obj_writer::write_object_name(std::string object_name)
+{
+    *out << "o " << object_name << "\n";
+}
+
+void obj_writer::write_mesh(const Mesh &mesh,
+                            vertex_attribute<glm::vec3> const &position,
+                            vertex_attribute<glm::vec2> const *tex_coord,
+                            vertex_attribute<glm::vec3> const *normal)
+{
+    auto base_v = vertex_idx;
+    auto base_t = texture_idx;
+    auto base_n = normal_idx;
+
+    for (auto v : mesh.vertices())
+    {
+        auto pos = v[position];
+        *out << "v " << pos.x << " " << pos.y << " " << pos.z << "\n";
+        ++vertex_idx;
+    }
+
+    if (tex_coord)
+        for (auto v : mesh.vertices())
+        {
+            auto t = v[*tex_coord];
+            *out << "vt " << t.x << " " << t.y << "\n";
+            ++texture_idx;
+        }
+
+    if (normal)
+        for (auto v : mesh.vertices())
+        {
+            auto n = v[*normal];
+            *out << "vn " << n.x << " " << n.y << " " << n.z << "\n";
+            ++normal_idx;
+        }
+
+    for (auto f : mesh.valid_faces())
+    {
+        *out << "f";
+        for (auto v : f.vertices())
+        {
+            auto i = v.idx.value;
+            *out << " ";
+            *out << base_v + i;
+            if (tex_coord || normal)
+                *out << "/";
+            if (tex_coord)
+                *out << base_t + i;
+            if (normal)
+            {
+                *out << base_n + i;
+                *out << "/";
+            }
+        }
+        *out << "\n";
+    }
+}
+
+void write_obj(const std::string &filename,
+               const Mesh &mesh,
+               const vertex_attribute<glm::vec3> &position,
+               const vertex_attribute<glm::vec2> *tex_coord,
+               const vertex_attribute<glm::vec3> *normal)
+{
+    obj_writer obj(filename);
+    obj.write_mesh(mesh, position, tex_coord, normal);
+}
+
+obj_reader::obj_reader(const std::string &filename, Mesh &mesh)
+  : positions(mesh.vertices().make_attribute<glm::vec4>()),
+    tex_coords(mesh.halfedges().make_attribute<glm::vec3>()),
+    normals(mesh.halfedges().make_attribute<glm::vec3>())
+{
+    std::ifstream file(filename);
+    if (!file.good())
+        std::cerr << "Cannot read from file `" << filename << "'" << std::endl;
+    else
+        parse(file, mesh);
+}
+
+obj_reader::obj_reader(std::istream &in, Mesh &mesh)
+  : positions(mesh.vertices().make_attribute<glm::vec4>()),
+    tex_coords(mesh.halfedges().make_attribute<glm::vec3>()),
+    normals(mesh.halfedges().make_attribute<glm::vec3>())
+{
+    parse(in, mesh);
+}
+
+void obj_reader::parse(std::istream &in, Mesh &mesh)
+{
+    mesh.clear();
+
+    std::vector<glm::vec3> raw_tex_coords;
+    std::vector<glm::vec3> raw_normals;
+
+    struct face
+    {
+        int v = 0;
+        int t = 0;
+        int n = 0;
+        vertex_handle vh;
+    };
+
+    std::vector<face> poly;
+    std::vector<halfedge_handle> poly_hs;
+    std::string fs;
+
+    std::string line_s;
+    auto line_nr = 0;
+    while (std::getline(in, line_s))
+    {
+        ++line_nr;
+        std::istringstream line(line_s);
+        std::string type;
+
+        line >> type;
+
+        // empty lines
+        if (type.empty())
+            continue;
+
+        // comments
+        else if (type[0] == '#')
+            continue;
+
+        // vertices
+        else if (type == "v")
+        {
+            auto v = mesh.vertices().add();
+
+            glm::vec4 p;
+            p.w = 1.0f;
+
+            line >> p.x;
+            line >> p.y;
+            line >> p.z;
+            float w;
+            line >> w;
+            if (line.good())
+                p.w = w;
+
+            positions[v] = p;
+        }
+
+        // textures
+        else if (type == "vt")
+        {
+            glm::vec3 t;
+            t.z = 1.0f;
+
+            line >> t.x;
+            line >> t.y;
+            float z;
+            line >> z;
+            if (line.good())
+                t.z = z;
+
+            raw_tex_coords.push_back(t);
+        }
+
+        // normals
+        else if (type == "vn")
+        {
+            glm::vec3 n;
+            line >> n.x;
+            line >> n.y;
+            line >> n.z;
+            raw_normals.push_back(n);
+        }
+
+        // faces
+        else if (type == "f")
+        {
+            poly_hs.clear();
+            poly.clear();
+
+            while (line.good())
+            {
+                fs.clear();
+                line >> fs;
+                int sc = 0;
+                auto first_s = fs.find_first_of('/');
+                auto last_s = fs.find_last_of('/');
+                for (auto &c : fs)
+                    if (c == '/')
+                    {
+                        c = ' ';
+                        ++sc;
+                    }
+
+                std::istringstream ss(fs);
+                face f;
+                switch (sc)
+                {
+                case 0:
+                    ss >> f.v;
+                    break;
+
+                case 1:
+                    ss >> f.v;
+                    ss >> f.t;
+                    break;
+
+                case 2:
+                    ss >> f.v;
+                    if (first_s + 1 != last_s) // "1//2"
+                        ss >> f.t;
+                    ss >> f.n;
+                    break;
+                }
+
+                f.vh = mesh.handle_of(vertex_index(f.v - 1));
+                poly.push_back(f);
+            }
+
+            if (poly.size() < 3)
+            {
+                std::cerr << "faces with less than 3 vertices are not supported. Use lines instead." << std::endl;
+                continue;
+            }
+
+            poly_hs.resize(poly.size());
+            for (auto i = 0u; i < poly.size(); ++i)
+            {
+                auto const &f0 = poly[i];
+                auto const &f1 = poly[(i + 1) % poly.size()];
+                auto v0 = f0.vh;
+                auto v1 = f1.vh;
+                auto hh = mesh.halfedges().add_or_get(v0, v1);
+                poly_hs[i] = hh;
+
+                if (f1.t > 0)
+                    tex_coords[hh] = raw_tex_coords[f1.t - 1];
+                if (f1.n > 0)
+                    normals[hh] = raw_normals[f1.n - 1];
+            }
+
+            mesh.faces().add(poly_hs);
+        }
+
+        // lines
+        else if (type == "l")
+        {
+            int i0, i1;
+            line >> i0;
+            line >> i1;
+            while (line.good())
+            {
+                mesh.edges().add_or_get(mesh[vertex_index(i0 - 1)], mesh[vertex_index(i1 - 1)]);
+                i0 = i1;
+                line >> i1;
+            }
+        }
+
+        // not implemented
+        else if (type == "s")
+            continue;
+        else if (type == "o")
+            continue;
+        else if (type == "g")
+            continue;
+        else if (type == "usemtl")
+            continue;
+        else if (type == "mtllib")
+            continue;
+
+        else
+        {
+            std::cerr << "Unable to parse line " << line_nr << ": " << line_s << std::endl;
+        }
+    }
+}
diff --git a/src/polymesh/formats/obj.hh b/src/polymesh/formats/obj.hh
new file mode 100644
index 0000000000000000000000000000000000000000..0e81ff20c87925c5bd39ed20ad225c696035366d
--- /dev/null
+++ b/src/polymesh/formats/obj.hh
@@ -0,0 +1,63 @@
+#pragma once
+
+#include <glm/glm.hpp>
+
+#include <iostream>
+#include <string>
+
+#include "../Mesh.hh"
+
+namespace polymesh
+{
+void write_obj(std::string const& filename,
+               Mesh const& mesh,
+               vertex_attribute<glm::vec3> const& position,
+               vertex_attribute<glm::vec2> const* tex_coord = nullptr,
+               vertex_attribute<glm::vec3> const* normal = nullptr);
+
+struct obj_writer
+{
+    obj_writer(std::string const& filename);
+    obj_writer(std::ostream& out);
+    ~obj_writer();
+
+    void write_object_name(std::string object_name);
+    void write_mesh(Mesh const& mesh,
+                    vertex_attribute<glm::vec3> const& position,
+                    vertex_attribute<glm::vec2> const* tex_coord = nullptr,
+                    vertex_attribute<glm::vec3> const* normal = nullptr);
+
+    // TODO: tex coords and normals as half-edge attributes
+
+private:
+    std::ostream* tmp_out = nullptr;
+    std::ostream* out = nullptr;
+
+    int vertex_idx = 1;
+    int texture_idx = 1;
+    int normal_idx = 1;
+};
+
+// clears the given mesh before adding data
+// obj must be manifold
+// no negative indices
+struct obj_reader
+{
+    obj_reader(std::string const& filename, Mesh& mesh);
+    obj_reader(std::istream& in, Mesh& mesh);
+
+    // get properties of the obj
+    // NOTE: these always return fresh copies of the attribute!
+public:
+    vertex_attribute<glm::vec4> positions_vec4() const { return positions; }
+    halfedge_attribute<glm::vec3> tex_coords_vec3() const { return tex_coords; }
+    halfedge_attribute<glm::vec3> normals_vec3() const { return tex_coords; }
+
+private:
+    void parse(std::istream& in, Mesh& mesh);
+
+    vertex_attribute<glm::vec4> positions;
+    halfedge_attribute<glm::vec3> tex_coords;
+    halfedge_attribute<glm::vec3> normals;
+};
+}
diff --git a/src/polymesh/obj_writer.cc b/src/polymesh/obj_writer.cc
deleted file mode 100644
index 9620abc60f1e10eeaff47d540eda918689e0e7c5..0000000000000000000000000000000000000000
--- a/src/polymesh/obj_writer.cc
+++ /dev/null
@@ -1,82 +0,0 @@
-#include "obj_writer.hh"
-
-#include <fstream>
-
-#include "Mesh.hh"
-
-using namespace polymesh;
-
-obj_writer::obj_writer(const std::string &filename)
-{
-    tmp_out = new std::ofstream(filename);
-    out = tmp_out;
-}
-
-obj_writer::obj_writer(std::ostream &out)
-{
-    this->out = &out;
-}
-
-obj_writer::~obj_writer()
-{
-    delete tmp_out;
-}
-
-void obj_writer::write_object_name(std::string object_name)
-{
-    *out << "o " << object_name << "\n";
-}
-
-void obj_writer::write_mesh(const Mesh &mesh,
-                            vertex_attribute<glm::vec3> const &position,
-                            vertex_attribute<glm::vec2> const *tex_coord,
-                            vertex_attribute<glm::vec3> const *normal)
-{
-    auto base_v = vertex_idx;
-    auto base_t = texture_idx;
-    auto base_n = normal_idx;
-
-    for (auto v : mesh.vertices())
-    {
-        auto pos = v[position];
-        *out << "v " << pos.x << " " << pos.y << " " << pos.z << "\n";
-        ++vertex_idx;
-    }
-
-    if (tex_coord)
-        for (auto v : mesh.vertices())
-        {
-            auto t = v[*tex_coord];
-            *out << "vt " << t.x << " " << t.y << "\n";
-            ++texture_idx;
-        }
-
-    if (normal)
-        for (auto v : mesh.vertices())
-        {
-            auto n = v[*normal];
-            *out << "vn " << n.x << " " << n.y << " " << n.z << "\n";
-            ++normal_idx;
-        }
-
-    for (auto f : mesh.valid_faces())
-    {
-        *out << "f";
-        for (auto v : f.vertices())
-        {
-            auto i = v.idx.value;
-            *out << " ";
-            *out << base_v + i;
-            if (tex_coord || normal)
-                *out << "/";
-            if (tex_coord)
-                *out << base_t + i;
-            if (normal)
-            {
-                *out << base_n + i;
-                *out << "/";
-            }
-        }
-        *out << "\n";
-    }
-}
diff --git a/src/polymesh/obj_writer.hh b/src/polymesh/obj_writer.hh
deleted file mode 100644
index 5df42d0a8f2ec9c96222909992555a39b3e6b2f5..0000000000000000000000000000000000000000
--- a/src/polymesh/obj_writer.hh
+++ /dev/null
@@ -1,32 +0,0 @@
-#pragma once
-
-#include <glm/glm.hpp>
-
-#include <iostream>
-#include <string>
-
-#include "Mesh.hh"
-
-namespace polymesh
-{
-struct obj_writer
-{
-    obj_writer(std::string const& filename);
-    obj_writer(std::ostream& out);
-    ~obj_writer();
-
-    void write_object_name(std::string object_name);
-    void write_mesh(Mesh const& mesh,
-                    vertex_attribute<glm::vec3> const& position,
-                    vertex_attribute<glm::vec2> const* tex_coord = nullptr,
-                    vertex_attribute<glm::vec3> const* normal = nullptr);
-
-private:
-    std::ostream* tmp_out = nullptr;
-    std::ostream* out = nullptr;
-
-    int vertex_idx = 1;
-    int texture_idx = 1;
-    int normal_idx = 1;
-};
-}
diff --git a/src/polymesh/ranges.hh b/src/polymesh/ranges.hh
index e08e5f7e85024c2cbd661394112f77418f3e7b5b..ce442f40416eab9433bc8913791b15c0e25c9723 100644
--- a/src/polymesh/ranges.hh
+++ b/src/polymesh/ranges.hh
@@ -28,7 +28,7 @@ struct vertex_collection
     void remove(vertex_handle v) const;
 
     /// Creates a new vertex attribute
-    template <typename PropT>
+    template <class PropT>
     vertex_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -46,7 +46,7 @@ struct const_vertex_collection
     int size() const;
 
     /// Creates a new vertex attribute
-    template <typename PropT>
+    template <class PropT>
     vertex_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -90,20 +90,20 @@ struct face_collection
     face_handle add(vertex_handle v0, vertex_handle v1, vertex_handle v2) const;
     face_handle add(vertex_handle v0, vertex_handle v1, vertex_handle v2, vertex_handle v3) const;
     face_handle add(std::vector<vertex_handle> const& v_handles) const;
-    face_handle add(vertex_handle const* v_handles, size_t vcnt) const;
+    face_handle add(vertex_handle const* v_handles, int vcnt) const;
     template <size_t N>
     face_handle add(const halfedge_handle (&half_loop)[N]) const;
     face_handle add(halfedge_handle h0, halfedge_handle h1, halfedge_handle h2) const;
     face_handle add(halfedge_handle h0, halfedge_handle h1, halfedge_handle h2, halfedge_handle h3) const;
     face_handle add(std::vector<halfedge_handle> const& half_loop) const;
-    face_handle add(halfedge_handle const* half_loop, size_t vcnt) const;
+    face_handle add(halfedge_handle const* half_loop, int vcnt) const;
 
     /// Removes a face (adjacent edges and vertices are NOT removed)
     /// (marks it as removed, compactify mesh to actually remove it)
     void remove(face_handle f) const;
 
     /// Creates a new face attribute
-    template <typename PropT>
+    template <class PropT>
     face_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -121,7 +121,7 @@ struct const_face_collection
     int size() const;
 
     /// Creates a new face attribute
-    template <typename PropT>
+    template <class PropT>
     face_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -166,7 +166,7 @@ struct edge_collection
     void remove(edge_handle e) const;
 
     /// Creates a new edge attribute
-    template <typename PropT>
+    template <class PropT>
     edge_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -184,7 +184,7 @@ struct const_edge_collection
     int size() const;
 
     /// Creates a new edge attribute
-    template <typename PropT>
+    template <class PropT>
     edge_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -230,7 +230,7 @@ struct halfedge_collection
     void remove_edge(halfedge_handle h) const;
 
     /// Creates a new half-edge attribute
-    template <typename PropT>
+    template <class PropT>
     halfedge_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
@@ -248,7 +248,7 @@ struct const_halfedge_collection
     int size() const;
 
     /// Creates a new half-edge attribute
-    template <typename PropT>
+    template <class PropT>
     halfedge_attribute<PropT> make_attribute(PropT const& def_value = PropT());
 
     // Iteration:
diff --git a/src/polymesh/tmp.hh b/src/polymesh/tmp.hh
new file mode 100644
index 0000000000000000000000000000000000000000..d30a4776c52c815ca67967bef9ffaa80f33b8d16
--- /dev/null
+++ b/src/polymesh/tmp.hh
@@ -0,0 +1,19 @@
+#pragma once
+
+#include <utility>
+
+namespace polymesh
+{
+// small template metaprogramming
+namespace tmp
+{
+template <class FuncT, class ArgT>
+struct result_of
+{
+    using type = typename std::decay<decltype((std::declval<FuncT>())(std::declval<ArgT>()))>::type;
+};
+
+template <class FuncT, class ArgT>
+using result_type_of = typename result_of<FuncT, ArgT>::type;
+}
+}