diff --git a/include/boost/histogram/arithmetic_operators.hpp b/include/boost/histogram/arithmetic_operators.hpp
index 7372cce7..023f0483 100644
--- a/include/boost/histogram/arithmetic_operators.hpp
+++ b/include/boost/histogram/arithmetic_operators.hpp
@@ -13,78 +13,64 @@ namespace boost {
 namespace histogram {
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator+(H<A, S> &&a,
-          const H<A, S> &b) {
+H<A, S> &&operator+(H<A, S> &&a, const H<A, S> &b) {
   a += b;
   return std::move(a);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator+(H<A, S> &&a,
-          H<A, S> &&b) {
+H<A, S> &&operator+(H<A, S> &&a, H<A, S> &&b) {
   a += b;
   return std::move(a);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator+(const H<A, S> &a,
-          H<A, S> &&b) {
+H<A, S> &&operator+(const H<A, S> &a, H<A, S> &&b) {
   b += a;
   return std::move(b);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S>
-operator+(const H<A, S> &a,
-          const H<A, S> &b) {
+H<A, S> operator+(const H<A, S> &a, const H<A, S> &b) {
   H<A, S> r(a);
   r += b;
   return r;
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator*(H<A, S> &&a, const double x) {
+H<A, S> &&operator*(H<A, S> &&a, const double x) {
   a *= x;
   return std::move(a);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator*(const double x, H<A, S> &&b) {
+H<A, S> &&operator*(const double x, H<A, S> &&b) {
   b *= x;
   return std::move(b);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S>
-operator*(const H<A, S> &a, const double x) {
+H<A, S> operator*(const H<A, S> &a, const double x) {
   H<A, S> r(a);
   r *= x;
   return r;
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S>
-operator*(const double x, const H<A, S> &b) {
+H<A, S> operator*(const double x, const H<A, S> &b) {
   H<A, S> r(b);
   r *= x;
   return r;
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S> &&
-operator/(H<A, S> &&a, const double x) {
+H<A, S> &&operator/(H<A, S> &&a, const double x) {
   a /= x;
   return std::move(a);
 }
 
 template <template <class, class> class H, typename A, typename S>
-H<A, S>
-operator/(const H<A, S> &a, const double x) {
+H<A, S> operator/(const H<A, S> &a, const double x) {
   H<A, S> r(a);
   r /= x;
   return r;
diff --git a/include/boost/histogram/axis/any.hpp b/include/boost/histogram/axis/any.hpp
index 573e5406..a2029f16 100644
--- a/include/boost/histogram/axis/any.hpp
+++ b/include/boost/histogram/axis/any.hpp
@@ -8,10 +8,10 @@
 #define _BOOST_HISTOGRAM_AXIS_ANY_HPP_
 
 #include <boost/histogram/axis/iterator.hpp>
-#include <boost/histogram/interval.hpp>
-#include <boost/utility/string_view.hpp>
 #include <boost/histogram/detail/axis_visitor.hpp>
+#include <boost/histogram/interval.hpp>
 #include <boost/mpl/contains.hpp>
+#include <boost/utility/string_view.hpp>
 #include <boost/variant.hpp>
 #include <stdexcept>
 #include <type_traits>
@@ -42,13 +42,15 @@ struct uoflow : public static_visitor<bool> {
 };
 
 struct get_label : public static_visitor<string_view> {
-  template <typename A> ::boost::string_view operator()(const A& a) const { return a.label(); }
+  template <typename A>::boost::string_view operator()(const A &a) const {
+    return a.label();
+  }
 };
 
 struct set_label : public static_visitor<void> {
   const ::boost::string_view label;
   set_label(const ::boost::string_view x) : label(x) {}
-  template <typename A> void operator()(A& a) const { a.label(label); }
+  template <typename A> void operator()(A &a) const { a.label(label); }
 };
 
 template <typename T> struct index : public static_visitor<int> {
@@ -57,11 +59,12 @@ template <typename T> struct index : public static_visitor<int> {
   template <typename Axis> int operator()(const Axis &a) const {
     return impl(std::is_convertible<T, typename Axis::value_type>(), a);
   }
-  template <typename Axis> int impl(std::true_type, const Axis& a) const {
+  template <typename Axis> int impl(std::true_type, const Axis &a) const {
     return a.index(t);
   }
-  template <typename Axis> int impl(std::false_type, const Axis&) const {
-    throw std::runtime_error("index argument not convertible to axis value type");
+  template <typename Axis> int impl(std::false_type, const Axis &) const {
+    throw std::runtime_error(
+        "index argument not convertible to axis value type");
   }
 };
 
@@ -83,9 +86,9 @@ struct bin : public static_visitor<axis::interval<double>> {
 } // namespace detail
 
 /// Polymorphic axis type
-template <typename Axes>
-class any : public make_variant_over<Axes>::type {
+template <typename Axes> class any : public make_variant_over<Axes>::type {
   using base_type = typename make_variant_over<Axes>::type;
+
 public:
   using types = typename base_type::types;
   using value_type = double;
@@ -93,43 +96,34 @@ public:
   using const_iterator = axis_iterator<any>;
 
   any() = default;
-  any(const any& t) = default;
-  any(any&& t) = default;
-  any& operator=(const any& t) = default;
-  any& operator=(any&& t) = default;
+  any(const any &t) = default;
+  any(any &&t) = default;
+  any &operator=(const any &t) = default;
+  any &operator=(any &&t) = default;
 
   template <typename T, typename = typename std::enable_if<
-    mpl::contains<Axes, T>::value
-  >::type>
-  any(const T& t) : base_type(t) {}
+                            mpl::contains<Axes, T>::value>::type>
+  any(const T &t) : base_type(t) {}
 
   template <typename T, typename = typename std::enable_if<
-    mpl::contains<Axes, T>::value
-  >::type>
-  any& operator=(const T& t) {
+                            mpl::contains<Axes, T>::value>::type>
+  any &operator=(const T &t) {
     // ugly workaround for compiler bug
-    return reinterpret_cast<any&>(base_type::operator=(t));
+    return reinterpret_cast<any &>(base_type::operator=(t));
   }
 
   template <typename T, typename = typename std::enable_if<
-    mpl::contains<Axes, T>::value
-  >::type>
-  any& operator=(T&& t) {
+                            mpl::contains<Axes, T>::value>::type>
+  any &operator=(T &&t) {
     // ugly workaround for compiler bug
-    return reinterpret_cast<any&>(base_type::operator=(std::move(t)));
+    return reinterpret_cast<any &>(base_type::operator=(std::move(t)));
   }
 
-  int size() const {
-    return apply_visitor(detail::size(), *this);
-  }
+  int size() const { return apply_visitor(detail::size(), *this); }
 
-  int shape() const {
-    return apply_visitor(detail::shape(), *this);
-  }
+  int shape() const { return apply_visitor(detail::shape(), *this); }
 
-  bool uoflow() const {
-    return apply_visitor(detail::uoflow(), *this);
-  }
+  bool uoflow() const { return apply_visitor(detail::uoflow(), *this); }
 
   // note: this only works for axes with compatible value type
   int index(const value_type x) const {
@@ -150,8 +144,8 @@ public:
     return apply_visitor(detail::bin(i), *this);
   }
 
-  bool operator==(const any& rhs) const {
-    return base_type::operator==(static_cast<const base_type&>(rhs));
+  bool operator==(const any &rhs) const {
+    return base_type::operator==(static_cast<const base_type &>(rhs));
   }
 
   const_iterator begin() const { return const_iterator(*this, 0); }
@@ -160,34 +154,49 @@ public:
 
 private:
   friend class ::boost::serialization::access;
-  template <typename Archive> void serialize(Archive&, unsigned);
+  template <typename Archive> void serialize(Archive &, unsigned);
 };
 
 // dynamic casts
 template <typename T, typename Axes>
-typename std::add_lvalue_reference<T>::type cast(any<Axes>& any) { return get<T>(any); }
+typename std::add_lvalue_reference<T>::type cast(any<Axes> &any) {
+  return get<T>(any);
+}
 
 template <typename T, typename Axes>
-const typename std::add_lvalue_reference<T>::type cast(const any<Axes>& any) { return get<T>(any); }
+const typename std::add_lvalue_reference<T>::type cast(const any<Axes> &any) {
+  return get<T>(any);
+}
 
 template <typename T, typename Axes>
-typename std::add_pointer<T>::type cast(any<Axes>* any) { return get<T>(&any); }
+typename std::add_pointer<T>::type cast(any<Axes> *any) {
+  return get<T>(&any);
+}
 
 template <typename T, typename Axes>
-const typename std::add_pointer<T>::type cast(const any<Axes>* any) { return get<T>(&any); }
+const typename std::add_pointer<T>::type cast(const any<Axes> *any) {
+  return get<T>(&any);
+}
 
-// pass-through versions for generic programming, i.e. when you switch to static histogram
-template <typename T>
-typename std::add_lvalue_reference<T>::type cast(T& t) { return t; }
+// pass-through versions for generic programming, i.e. when you switch to static
+// histogram
+template <typename T> typename std::add_lvalue_reference<T>::type cast(T &t) {
+  return t;
+}
 
 template <typename T>
-const typename std::add_lvalue_reference<T>::type cast(const T& t) { return t; }
+const typename std::add_lvalue_reference<T>::type cast(const T &t) {
+  return t;
+}
+
+template <typename T> typename std::add_pointer<T>::type cast(T *t) {
+  return t;
+}
 
 template <typename T>
-typename std::add_pointer<T>::type cast(T* t) { return t; }
-
-template <typename T>
-const typename std::add_pointer<T>::type cast(const T* t) { return t; }
+const typename std::add_pointer<T>::type cast(const T *t) {
+  return t;
+}
 
 } // namespace axis
 } // namespace histogram
diff --git a/include/boost/histogram/axis/axis.hpp b/include/boost/histogram/axis/axis.hpp
index c011c329..15d8ad20 100644
--- a/include/boost/histogram/axis/axis.hpp
+++ b/include/boost/histogram/axis/axis.hpp
@@ -163,6 +163,7 @@ struct pow {
   bool operator==(const pow &other) const noexcept {
     return value == other.value;
   }
+
 private:
   friend ::boost::serialization::access;
   template <class Archive> void serialize(Archive &, unsigned);
diff --git a/include/boost/histogram/axis/iterator.hpp b/include/boost/histogram/axis/iterator.hpp
index 7e4908e9..8010322f 100644
--- a/include/boost/histogram/axis/iterator.hpp
+++ b/include/boost/histogram/axis/iterator.hpp
@@ -10,7 +10,9 @@
 #include <boost/iterator/iterator_facade.hpp>
 #include <utility>
 
-namespace boost { namespace histogram { namespace axis {
+namespace boost {
+namespace histogram {
+namespace axis {
 
 template <typename Axis>
 class axis_iterator
@@ -21,8 +23,8 @@ class axis_iterator
 public:
   explicit axis_iterator(const Axis &axis, int idx) : axis_(axis), idx_(idx) {}
 
-  axis_iterator(const axis_iterator &o) = default;
-  axis_iterator &operator=(const axis_iterator &o) = default;
+  axis_iterator(const axis_iterator &) = default;
+  axis_iterator &operator=(const axis_iterator &) = default;
 
 private:
   void increment() noexcept { ++idx_; }
@@ -37,11 +39,13 @@ private:
   std::pair<int, typename Axis::bin_type> dereference() const {
     return std::make_pair(idx_, axis_[idx_]);
   }
-  const Axis& axis_;
+  const Axis &axis_;
   int idx_;
-  friend class boost::iterator_core_access;
+  friend class ::boost::iterator_core_access;
 };
 
-}}}
+} // namespace axis
+} // namespace histogram
+} // namespace boost
 
 #endif
diff --git a/include/boost/histogram/detail/axis_visitor.hpp b/include/boost/histogram/detail/axis_visitor.hpp
index bfafc681..a736e3ea 100644
--- a/include/boost/histogram/detail/axis_visitor.hpp
+++ b/include/boost/histogram/detail/axis_visitor.hpp
@@ -14,7 +14,6 @@
 #include <boost/fusion/include/size.hpp>
 #include <boost/fusion/support/is_sequence.hpp>
 #include <boost/histogram/interval.hpp>
-#include <boost/utility/string_view.hpp>
 #include <boost/mpl/bool.hpp>
 #include <boost/mpl/contains.hpp>
 #include <boost/variant/get.hpp>
@@ -83,7 +82,7 @@ template <typename Iterator> struct fusion_assign_axis2 {
   template <typename T> void operator()(const T &t) const { *(iter++) = t; }
 };
 
-struct field_count : public static_visitor<void> {
+struct field_count_visitor : public static_visitor<void> {
   mutable std::size_t value = 1;
   template <typename T> void operator()(const T &t) const {
     value *= t.shape();
diff --git a/include/boost/histogram/detail/meta.hpp b/include/boost/histogram/detail/meta.hpp
index 0fec7dd0..b5a40694 100644
--- a/include/boost/histogram/detail/meta.hpp
+++ b/include/boost/histogram/detail/meta.hpp
@@ -11,10 +11,10 @@
 #include <boost/fusion/include/fold.hpp>
 #include <boost/mpl/at.hpp>
 #include <boost/mpl/back_inserter.hpp>
-#include <boost/mpl/int.hpp>
 #include <boost/mpl/contains.hpp>
 #include <boost/mpl/copy_if.hpp>
 #include <boost/mpl/for_each.hpp>
+#include <boost/mpl/int.hpp>
 #include <boost/mpl/not.hpp>
 #include <boost/mpl/sort.hpp>
 #include <boost/mpl/transform.hpp>
@@ -33,24 +33,23 @@ namespace detail {
 template <typename T, typename = decltype(std::declval<T &>().size(),
                                           std::declval<T &>().increase(0),
                                           std::declval<T &>().value(0))>
-struct is_storage {};
+struct requires_storage {};
 
-template <typename T>
-struct has_variance_support
-{
-  template<typename U, typename T::value_type (U::*)(std::size_t) const> struct SFINAE {};
-  template<typename U> static std::true_type Test(SFINAE<U, &U::variance>*);
-  template<typename U> static std::false_type Test(...);
+template <typename T> struct has_variance_support {
+  template <typename U, typename T::value_type (U::*)(std::size_t) const>
+  struct SFINAE {};
+  template <typename U> static std::true_type Test(SFINAE<U, &U::variance> *);
+  template <typename U> static std::false_type Test(...);
   using type = decltype(Test<T>(nullptr));
 };
 
 template <typename T>
 using has_variance_support_t = typename has_variance_support<T>::type;
 
-template <typename S>
-using requires_variance_support = typename std::enable_if<
-    has_variance_support_t<S>::value, typename S::value_type
-  >::type;
+template <typename S, typename ReturnValue = typename S::value_type>
+using requires_variance_support =
+    typename std::enable_if<has_variance_support_t<S>::value,
+                            ReturnValue>::type;
 
 template <typename T,
           typename = decltype(*std::declval<T &>(), ++std::declval<T &>())>
@@ -99,12 +98,12 @@ void axes_assign_subset(Axes1 &axes1, const Axes &axes) {
 
 template <typename Ns>
 using unique_sorted =
-  typename mpl::unique<typename mpl::sort<Ns>::type,
-                       std::is_same<mpl::_1, mpl::_2>>::type;
+    typename mpl::unique<typename mpl::sort<Ns>::type,
+                         std::is_same<mpl::_1, mpl::_2>>::type;
 
 template <typename Axes, typename Numbers>
 using axes_select =
-  typename mpl::transform<Numbers, mpl::at<Axes, mpl::_>>::type;
+    typename mpl::transform<Numbers, mpl::at<Axes, mpl::_>>::type;
 
 } // namespace detail
 } // namespace histogram
diff --git a/include/boost/histogram/detail/utility.hpp b/include/boost/histogram/detail/utility.hpp
index 86900ea2..a52df239 100644
--- a/include/boost/histogram/detail/utility.hpp
+++ b/include/boost/histogram/detail/utility.hpp
@@ -8,8 +8,8 @@
 #define _BOOST_HISTOGRAM_DETAIL_UTILITY_HPP_
 
 #include <algorithm>
-#include <boost/utility/string_view.hpp>
 #include <boost/histogram/detail/meta.hpp>
+#include <boost/utility/string_view.hpp>
 #include <ostream>
 #include <vector>
 
@@ -18,20 +18,20 @@ namespace histogram {
 namespace detail {
 
 template <typename Storage>
-void storage_add_impl(std::true_type, Storage& result, const Storage& input,
+void storage_add_impl(std::true_type, Storage &result, const Storage &input,
                       std::size_t dst, std::size_t src) {
   result.add(dst, input.value(src), input.variance(src));
 }
 
 template <typename Storage>
-void storage_add_impl(std::false_type, Storage& result, const Storage& input,
+void storage_add_impl(std::false_type, Storage &result, const Storage &input,
                       std::size_t dst, std::size_t src) {
   result.add(dst, input.value(src));
 }
 
 template <typename Storage>
-void storage_add(Storage& result, const Storage& input,
-                 std::size_t dst, std::size_t src) {
+void storage_add(Storage &result, const Storage &input, std::size_t dst,
+                 std::size_t src) {
   storage_add_impl(has_variance_support_t<Storage>(), result, input, dst, src);
 }
 
diff --git a/include/boost/histogram/detail/weight_counter.hpp b/include/boost/histogram/detail/weight_counter.hpp
deleted file mode 100644
index 50592c9e..00000000
--- a/include/boost/histogram/detail/weight_counter.hpp
+++ /dev/null
@@ -1,87 +0,0 @@
-// Copyright 2015-2016 Hans Dembinski
-//
-// Distributed under the Boost Software License, version 1.0.
-// (See accompanying file LICENSE_1_0.txt
-// or copy at http://www.boost.org/LICENSE_1_0.txt)
-
-#ifndef _BOOST_HISTOGRAM_DETAIL_WEIGHT_HPP_
-#define _BOOST_HISTOGRAM_DETAIL_WEIGHT_HPP_
-
-namespace boost {
-namespace histogram {
-namespace detail {
-
-/// Double counter which holds a sum of weights and a sum of squared weights
-struct weight_counter {
-  double w, w2;
-  weight_counter() = default;
-  weight_counter(const weight_counter &) = default;
-  weight_counter(weight_counter &&) = default;
-  weight_counter &operator=(const weight_counter &) = default;
-  weight_counter &operator=(weight_counter &&) = default;
-
-  weight_counter(double value, double variance) : w(value), w2(variance) {}
-
-  weight_counter &operator++() {
-    ++w;
-    ++w2;
-    return *this;
-  }
-
-  weight_counter &operator+=(const weight_counter &rhs) {
-    w += rhs.w;
-    w2 += rhs.w2;
-    return *this;
-  }
-
-  weight_counter &operator*=(const double x) {
-    w *= x;
-    w2 *= x*x;
-    return *this;
-  }
-
-  bool operator==(const weight_counter &rhs) const {
-    return w == rhs.w && w2 == rhs.w2;
-  }
-  bool operator!=(const weight_counter &rhs) const { return !operator==(rhs); }
-  template <typename T> bool operator==(const T &rhs) const {
-    return w == static_cast<double>(rhs) && w2 == static_cast<double>(rhs);
-  }
-  template <typename T> bool operator!=(const T &rhs) const {
-    return !operator==(rhs);
-  }
-
-  weight_counter &operator+=(double t) {
-    w += t;
-    w2 += t * t;
-    return *this;
-  }
-
-  template <typename T>
-  explicit weight_counter(const T &t)
-      : w(static_cast<double>(t)), w2(static_cast<double>(t)) {}
-
-  template <typename T> weight_counter &operator=(const T &t) {
-    w = static_cast<double>(t);
-    w2 = static_cast<double>(t);
-    return *this;
-  }
-  template <typename T> weight_counter &operator+=(const T &t) {
-    w += static_cast<double>(t);
-    w2 += static_cast<double>(t);
-    return *this;
-  }
-};
-
-template <typename T> bool operator==(const T &t, const weight_counter &w) {
-  return w == t;
-}
-
-template <typename T> bool operator!=(const T &t, const weight_counter &w) {
-  return !(w == t);
-}
-} // namespace detail
-} // namespace histogram
-} // namespace boost
-
-#endif
diff --git a/include/boost/histogram/dynamic_histogram.hpp b/include/boost/histogram/dynamic_histogram.hpp
index 4d5db349..5f6b1eee 100644
--- a/include/boost/histogram/dynamic_histogram.hpp
+++ b/include/boost/histogram/dynamic_histogram.hpp
@@ -9,13 +9,14 @@
 
 #include <algorithm>
 #include <boost/config.hpp>
-#include <boost/histogram/axis/axis.hpp>
 #include <boost/histogram/axis/any.hpp>
+#include <boost/histogram/axis/axis.hpp>
 #include <boost/histogram/detail/axis_visitor.hpp>
 #include <boost/histogram/detail/meta.hpp>
 #include <boost/histogram/detail/utility.hpp>
 #include <boost/histogram/histogram_fwd.hpp>
 #include <boost/histogram/storage/operators.hpp>
+#include <boost/histogram/value_iterator.hpp>
 #include <boost/mpl/count.hpp>
 #include <boost/mpl/empty.hpp>
 #include <boost/mpl/int.hpp>
@@ -43,12 +44,13 @@ class access;
 namespace boost {
 namespace histogram {
 
-template <typename Axes, typename Storage>
-class dynamic_histogram {
+template <typename Axes, typename Storage> class dynamic_histogram {
   static_assert(!mpl::empty<Axes>::value, "at least one axis required");
+
 public:
   using any_axis_type = axis::any<Axes>;
   using value_type = typename Storage::value_type;
+  using value_iterator = value_iterator<Storage>;
 
 private:
   using axes_type = std::vector<any_axis_type>;
@@ -61,7 +63,8 @@ public:
   dynamic_histogram &operator=(dynamic_histogram &&) = default;
 
   template <typename... Axes1>
-  explicit dynamic_histogram(const Axes1 &... axes) : axes_({any_axis_type(axes)...}) {
+  explicit dynamic_histogram(const Axes1 &... axes)
+      : axes_({any_axis_type(axes)...}) {
     storage_ = Storage(bincount_from_axes());
   }
 
@@ -73,12 +76,14 @@ public:
   }
 
   template <typename A, typename S>
-  explicit dynamic_histogram(const static_histogram<A, S> &rhs) : storage_(rhs.storage_) {
+  explicit dynamic_histogram(const static_histogram<A, S> &rhs)
+      : storage_(rhs.storage_) {
     detail::axes_assign(axes_, rhs.axes_);
   }
 
   template <typename A, typename S>
-  explicit dynamic_histogram(const dynamic_histogram<A, S> &rhs) : storage_(rhs.storage_) {
+  explicit dynamic_histogram(const dynamic_histogram<A, S> &rhs)
+      : storage_(rhs.storage_) {
     detail::axes_assign(axes_, rhs.axes_);
   }
 
@@ -166,15 +171,17 @@ public:
         (n_count::value + n_weight::value) <= 1,
         "arguments may contain at most one instance of type count or weight");
     if (dim() != sizeof...(args) - n_count::value - n_weight::value)
-      throw std::invalid_argument("fill arguments does not match histogram dimension");
+      throw std::invalid_argument(
+          "fill arguments does not match histogram dimension");
     fill_impl(mpl::int_<(n_count::value + 2 * n_weight::value)>(),
               std::forward<Args>(args)...);
   }
 
   template <typename Iterator, typename = detail::is_iterator<Iterator>>
   void fill(Iterator begin, Iterator end) {
-    if(dim() != std::distance(begin, end))
-      throw std::invalid_argument("fill iterator range does not match histogram dimension");
+    if (dim() != std::distance(begin, end))
+      throw std::invalid_argument(
+          "fill iterator range does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
     xlin_iter(idx, stride, begin);
     if (stride) {
@@ -184,8 +191,9 @@ public:
 
   template <typename Iterator, typename = detail::is_iterator<Iterator>>
   void fill(Iterator begin, Iterator end, const count n) {
-    if(dim() != std::distance(begin, end))
-      throw std::invalid_argument("fill iterator range does not match histogram dimension");
+    if (dim() != std::distance(begin, end))
+      throw std::invalid_argument(
+          "fill iterator range does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
     xlin_iter(idx, stride, begin);
     if (stride) {
@@ -195,20 +203,22 @@ public:
 
   template <typename Iterator, typename = detail::is_iterator<Iterator>>
   void fill(Iterator begin, Iterator end, const weight w) {
-    if(dim() != std::distance(begin, end))
-      throw std::invalid_argument("fill iterator range does not match histogram dimension");
+    if (dim() != std::distance(begin, end))
+      throw std::invalid_argument(
+          "fill iterator range does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
     xlin_iter(idx, stride, begin);
     if (stride) {
-      storage_.weighted_increase(idx, w.value);
+      storage_.increase_by_weight(idx, w.value);
     }
   }
 
   template <typename... Indices> value_type value(Indices &&... indices) const {
-    if(dim() != sizeof...(indices))
-      throw std::invalid_argument("value arguments does not match histogram dimension");
+    if (dim() != sizeof...(indices))
+      throw std::invalid_argument(
+          "value arguments does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
-    lin<0>(idx, stride, std::forward<Indices>(indices)...);
+    lin<0>(idx, stride, indices...);
     if (stride == 0)
       throw std::out_of_range("invalid index");
     return storage_.value(idx);
@@ -216,8 +226,9 @@ public:
 
   template <typename Iterator, typename = detail::is_iterator<Iterator>>
   value_type value(Iterator begin, Iterator end) const {
-    if(dim() != std::distance(begin, end))
-      throw std::invalid_argument("value iterator range does not match histogram dimension");
+    if (dim() != std::distance(begin, end))
+      throw std::invalid_argument(
+          "value iterator range does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
     lin_iter(idx, stride, begin);
     if (stride == 0)
@@ -226,23 +237,25 @@ public:
   }
 
   template <typename S = Storage, typename... Indices>
-  detail::requires_variance_support<S>
-  variance(Indices &&... indices) const {
-    if(dim() != sizeof...(indices))
-      throw std::invalid_argument("variance arguments does not match histogram dimension");
+  detail::requires_variance_support<S> variance(Indices &&... indices) const {
+    if (dim() != sizeof...(indices))
+      throw std::invalid_argument(
+          "variance arguments does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
-    lin<0>(idx, stride, std::forward<Indices>(indices)...);
+    lin<0>(idx, stride, indices...);
     if (stride == 0) {
       throw std::out_of_range("invalid index");
     }
     return storage_.variance(idx);
   }
 
-  template <typename S = Storage, typename Iterator, typename = detail::is_iterator<Iterator>>
-  detail::requires_variance_support<S>
-  variance(Iterator begin, Iterator end) const {
-    if(dim() != std::distance(begin, end))
-      throw std::invalid_argument("variance iterator range does not match histogram dimension");
+  template <typename S = Storage, typename Iterator,
+            typename = detail::is_iterator<Iterator>>
+  detail::requires_variance_support<S> variance(Iterator begin,
+                                                Iterator end) const {
+    if (dim() != std::distance(begin, end))
+      throw std::invalid_argument(
+          "variance iterator range does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
     lin_iter(idx, stride, begin);
     if (stride == 0) {
@@ -272,33 +285,36 @@ public:
 
   /// Return axis \a i
   any_axis_type &axis(unsigned i = 0) {
-    if(i >= dim())
+    if (i >= dim())
       throw std::out_of_range("axis index out of range");
     return axes_[i];
   }
 
   /// Return axis \a i (const version)
   const any_axis_type &axis(unsigned i = 0) const {
-    if(i >= dim())
+    if (i >= dim())
       throw std::out_of_range("axis index out of range");
     return axes_[i];
   }
 
   /// Apply unary functor/function to each axis
-  template <typename Unary> void for_each_axis(Unary &unary) const {
+  template <typename Unary> void for_each_axis(Unary &&unary) const {
     for (const auto &a : axes_) {
       apply_visitor(detail::unary_visitor<Unary>(unary), a);
     }
   }
 
   /// Return a lower dimensional histogram
-  template <int N, typename... Rest> dynamic_histogram reduce_to(mpl::int_<N>, Rest...) const {
-    const auto b = detail::bool_mask<mpl::vector<mpl::int_<N>, Rest...>>(dim(), true);
+  template <int N, typename... Rest>
+  dynamic_histogram reduce_to(mpl::int_<N>, Rest...) const {
+    const auto b =
+        detail::bool_mask<mpl::vector<mpl::int_<N>, Rest...>>(dim(), true);
     return reduce_impl(b);
   }
 
   /// Return a lower dimensional histogram
-  template <typename... Rest> dynamic_histogram reduce_to(int n, Rest... rest) const {
+  template <typename... Rest>
+  dynamic_histogram reduce_to(int n, Rest... rest) const {
     std::vector<bool> b(dim(), false);
     for (const auto &i : {n, rest...})
       b[i] = true;
@@ -314,14 +330,22 @@ public:
     return reduce_impl(b);
   }
 
+  value_iterator begin() const noexcept {
+    return value_iterator(*this, storage_);
+  }
+
+  value_iterator end() const noexcept {
+    return value_iterator(storage_);
+  }
+
 private:
   axes_type axes_;
   Storage storage_;
 
   std::size_t bincount_from_axes() const noexcept {
-    detail::field_count fc;
-    for_each_axis(fc);
-    return fc.value;
+    detail::field_count_visitor v;
+    for_each_axis(v);
+    return v.value;
   }
 
   template <typename... Args>
@@ -349,7 +373,7 @@ private:
     double w = 0.0;
     xlin_w<0>(idx, stride, w, std::forward<Args>(args)...);
     if (stride) {
-      storage_.weighted_increase(idx, w);
+      storage_.increase_by_weight(idx, w);
     }
   }
 
@@ -374,14 +398,13 @@ private:
       impl(std::is_convertible<Value, typename Axis::value_type>(), a);
     }
 
-    template <typename Axis>
-    void impl(std::true_type, const Axis& a) const {
+    template <typename Axis> void impl(std::true_type, const Axis &a) const {
       detail::xlin(idx, stride, a, val);
     }
 
-    template <typename Axis>
-    void impl(std::false_type, const Axis&) const {
-      throw std::runtime_error("fill argument not convertible to axis value type");
+    template <typename Axis> void impl(std::false_type, const Axis &) const {
+      throw std::runtime_error(
+          "fill argument not convertible to axis value type");
     }
   };
 
@@ -436,9 +459,9 @@ private:
 
   template <unsigned D, typename First, typename... Rest>
   inline typename disable_if<is_same<First, count>>::type
-  xlin_n(std::size_t &idx, std::size_t &stride, unsigned &x, First &&f,
+  xlin_n(std::size_t &idx, std::size_t &stride, unsigned &x, First &&first,
          Rest &&... rest) const {
-    apply_visitor(xlin_visitor<First>{idx, stride, std::forward<First>(f)},
+    apply_visitor(xlin_visitor<First>{idx, stride, std::forward<First>(first)},
                   axes_[D]);
     return xlin_n<D + 1>(idx, stride, x, std::forward<Rest>(rest)...);
   }
@@ -491,31 +514,37 @@ dynamic_histogram<detail::combine_t<axis::builtins, mpl::vector<Axes...>>>
 make_dynamic_histogram(Axes &&... axes) {
 
   return dynamic_histogram<
-                   detail::combine_t<axis::builtins, mpl::vector<Axes...>>>(
+      detail::combine_t<axis::builtins, mpl::vector<Axes...>>>(
       std::forward<Axes>(axes)...);
 }
 
 template <typename Storage, typename... Axes>
 dynamic_histogram<detail::combine_t<axis::builtins, mpl::vector<Axes...>>,
-                 Storage>
+                  Storage>
 make_dynamic_histogram_with(Axes &&... axes) {
   return dynamic_histogram<
-    detail::combine_t<axis::builtins, mpl::vector<Axes...>>, Storage>(
+      detail::combine_t<axis::builtins, mpl::vector<Axes...>>, Storage>(
       std::forward<Axes>(axes)...);
 }
 
 template <typename Iterator, typename = detail::is_iterator<Iterator>>
-dynamic_histogram<detail::combine_t<axis::builtins, typename Iterator::value_type::types>>
+dynamic_histogram<
+    detail::combine_t<axis::builtins, typename Iterator::value_type::types>>
 make_dynamic_histogram(Iterator begin, Iterator end) {
-  return dynamic_histogram<detail::combine_t<axis::builtins, typename Iterator::value_type::types>>(
+  return dynamic_histogram<
+      detail::combine_t<axis::builtins, typename Iterator::value_type::types>>(
       begin, end);
 }
 
-template <typename Storage, typename Iterator, typename = detail::is_iterator<Iterator>>
-dynamic_histogram<detail::combine_t<axis::builtins, typename Iterator::value_type::types>, Storage>
+template <typename Storage, typename Iterator,
+          typename = detail::is_iterator<Iterator>>
+dynamic_histogram<
+    detail::combine_t<axis::builtins, typename Iterator::value_type::types>,
+    Storage>
 make_dynamic_histogram_with(Iterator begin, Iterator end) {
-  return dynamic_histogram<detail::combine_t<axis::builtins, typename Iterator::value_type::types>, Storage>(
-      begin, end);
+  return dynamic_histogram<
+      detail::combine_t<axis::builtins, typename Iterator::value_type::types>,
+      Storage>(begin, end);
 }
 
 } // namespace histogram
diff --git a/include/boost/histogram/histogram.hpp b/include/boost/histogram/histogram.hpp
index 4745278e..87bba5d8 100644
--- a/include/boost/histogram/histogram.hpp
+++ b/include/boost/histogram/histogram.hpp
@@ -7,8 +7,8 @@
 #ifndef _BOOST_HISTOGRAM_HISTOGRAM_HPP_
 #define _BOOST_HISTOGRAM_HISTOGRAM_HPP_
 
+#include <boost/histogram/arithmetic_operators.hpp>
 #include <boost/histogram/dynamic_histogram.hpp>
 #include <boost/histogram/static_histogram.hpp>
-#include <boost/histogram/arithmetic_operators.hpp>
 
 #endif
diff --git a/include/boost/histogram/histogram_fwd.hpp b/include/boost/histogram/histogram_fwd.hpp
index 2a4d76e2..ce765269 100644
--- a/include/boost/histogram/histogram_fwd.hpp
+++ b/include/boost/histogram/histogram_fwd.hpp
@@ -18,14 +18,15 @@ class adaptive_storage;
 namespace axis {
 
 namespace transform {
-  struct identity;
-  struct log;
-  struct sqrt;
-  struct cos;
-  struct pow;
-}
+struct identity;
+struct log;
+struct sqrt;
+struct cos;
+struct pow;
+} // namespace transform
 
-template <typename RealType = double, typename Transform = transform::identity> class regular;
+template <typename RealType = double, typename Transform = transform::identity>
+class regular;
 template <typename RealType = double> class circular;
 template <typename RealType = double> class variable;
 template <typename IntType = int> class integer;
@@ -36,19 +37,16 @@ using builtins =
                 axis::regular<double, axis::transform::sqrt>,
                 axis::regular<double, axis::transform::cos>,
                 axis::regular<double, axis::transform::pow>, axis::circular<>,
-                axis::variable<>, axis::integer<>,
-                axis::category<>, axis::category<std::string>>;
+                axis::variable<>, axis::integer<>, axis::category<>,
+                axis::category<std::string>>;
 
-template <typename Axes = builtins>
-class any;
+template <typename Axes = builtins> class any;
 
 } // namespace axis
 
-template <class Axes, class Storage = adaptive_storage>
-class static_histogram;
+template <class Axes, class Storage = adaptive_storage> class static_histogram;
 
-template <class Axes, class Storage = adaptive_storage>
-class dynamic_histogram;
+template <class Axes, class Storage = adaptive_storage> class dynamic_histogram;
 
 struct weight {
   weight(double w) : value(w) {}
diff --git a/include/boost/histogram/ostream_operators.hpp b/include/boost/histogram/ostream_operators.hpp
index a8fc8f81..ec4f5c20 100644
--- a/include/boost/histogram/ostream_operators.hpp
+++ b/include/boost/histogram/ostream_operators.hpp
@@ -26,8 +26,7 @@ struct axis_ostream_visitor {
 template <template <class, class> class H, typename A, typename S>
 inline std::ostream &operator<<(std::ostream &os, const H<A, S> &h) {
   os << "histogram(";
-  detail::axis_ostream_visitor sh(os);
-  h.for_each_axis(sh);
+  h.for_each_axis(detail::axis_ostream_visitor(os));
   os << (h.dim() ? "\n)" : ")");
   return os;
 }
diff --git a/include/boost/histogram/serialization.hpp b/include/boost/histogram/serialization.hpp
index a266e30e..c67c7264 100644
--- a/include/boost/histogram/serialization.hpp
+++ b/include/boost/histogram/serialization.hpp
@@ -10,10 +10,10 @@
 #include <boost/fusion/algorithm/iteration/for_each.hpp>
 #include <boost/fusion/include/for_each.hpp>
 #include <boost/histogram/detail/utility.hpp>
-#include <boost/histogram/detail/weight_counter.hpp>
 #include <boost/histogram/histogram.hpp>
 #include <boost/histogram/storage/adaptive_storage.hpp>
 #include <boost/histogram/storage/array_storage.hpp>
+#include <boost/histogram/storage/weight_counter.hpp>
 #include <boost/serialization/array.hpp>
 #include <boost/serialization/unique_ptr.hpp>
 #include <boost/serialization/variant.hpp>
@@ -28,21 +28,20 @@ namespace boost {
 namespace histogram {
 
 namespace detail {
-
-template <class Archive>
-void serialize(Archive &ar, weight_counter &wt, unsigned /* version */) {
-  ar &wt.w;
-  ar &wt.w2;
-}
-
 template <typename Archive> struct serialize_helper {
   Archive &ar_;
   explicit serialize_helper(Archive &ar) : ar_(ar) {}
   template <typename T> void operator()(T &t) const { ar_ &t; }
 };
-
 } // namespace detail
 
+template <typename RealType>
+template <class Archive>
+void weight_counter<RealType>::serialize(Archive &ar, unsigned /* version */) {
+  ar &w;
+  ar &w2;
+}
+
 template <class Archive, typename Container>
 void serialize(Archive &ar, array_storage<Container> &store,
                unsigned /* version */) {
@@ -105,13 +104,11 @@ void adaptive_storage::serialize(Archive &ar, unsigned /* version */) {
       tid = 4u;
       ar &tid;
       ar &serialization::make_array(a->begin(), size);
-    } else if (auto *a =
-                   get<array<detail::mp_int>>(&buffer_)) {
+    } else if (auto *a = get<array<detail::mp_int>>(&buffer_)) {
       tid = 5u;
       ar &tid;
       ar &serialization::make_array(a->begin(), size);
-    } else if (auto *a =
-                   get<array<detail::weight_counter>>(&buffer_)) {
+    } else if (auto *a = get<array<detail::weight_counter>>(&buffer_)) {
       tid = 6u;
       ar &tid;
       ar &serialization::make_array(a->begin(), size);
diff --git a/include/boost/histogram/static_histogram.hpp b/include/boost/histogram/static_histogram.hpp
index 8d6b935d..a273b81f 100644
--- a/include/boost/histogram/static_histogram.hpp
+++ b/include/boost/histogram/static_histogram.hpp
@@ -25,6 +25,7 @@
 #include <boost/histogram/detail/utility.hpp>
 #include <boost/histogram/histogram_fwd.hpp>
 #include <boost/histogram/storage/operators.hpp>
+#include <boost/histogram/value_iterator.hpp>
 #include <boost/mpl/count.hpp>
 #include <boost/mpl/empty.hpp>
 #include <boost/mpl/int.hpp>
@@ -41,14 +42,14 @@ class access;
 namespace boost {
 namespace histogram {
 
-template <typename Axes, typename Storage>
-class static_histogram {
+template <typename Axes, typename Storage> class static_histogram {
   static_assert(!mpl::empty<Axes>::value, "at least one axis required");
   using axes_size = typename fusion::result_of::size<Axes>::type;
 
 public:
   using axes_type = typename fusion::result_of::as_vector<Axes>::type;
   using value_type = typename Storage::value_type;
+  using value_iterator = value_iterator<Storage>;
 
   static_histogram() = default;
   static_histogram(const static_histogram &rhs) = default;
@@ -66,12 +67,14 @@ public:
   }
 
   template <typename S>
-  explicit static_histogram(const static_histogram<Axes, S> &rhs) : storage_(rhs.storage_) {
+  explicit static_histogram(const static_histogram<Axes, S> &rhs)
+      : storage_(rhs.storage_) {
     detail::axes_assign(axes_, rhs.axes_);
   }
 
   template <typename A, typename S>
-  explicit static_histogram(const dynamic_histogram<A, S> &rhs) : storage_(rhs.storage_) {
+  explicit static_histogram(const dynamic_histogram<A, S> &rhs)
+      : storage_(rhs.storage_) {
     detail::axes_assign(axes_, rhs.axes_);
   }
 
@@ -168,9 +171,8 @@ public:
     return storage_.value(idx);
   }
 
-  template <typename S=Storage, typename... Indices>
-  detail::requires_variance_support<S>
-  variance(Indices &&... indices) const {
+  template <typename S = Storage, typename... Indices>
+  detail::requires_variance_support<S> variance(Indices &&... indices) const {
     static_assert(sizeof...(indices) == axes_size::value,
                   "number of arguments does not match histogram dimension");
     std::size_t idx = 0, stride = 1;
@@ -223,35 +225,44 @@ public:
   }
 
   // Get first axis (convenience for 1-d histograms)
-  typename fusion::result_of::value_at_c<axes_type, 0>::type &
-  axis() {
+  typename fusion::result_of::value_at_c<axes_type, 0>::type &axis() {
     return fusion::at_c<0>(axes_);
   }
 
   /// Apply unary functor/function to each axis
-  template <typename Unary> void for_each_axis(Unary &unary) const {
+  template <typename Unary> void for_each_axis(Unary &&unary) const {
     fusion::for_each(axes_, unary);
   }
 
   /// Returns a lower-dimensional histogram
   template <int N, typename... Rest>
-  auto reduce_to(mpl::int_<N>, Rest...) const
-      -> static_histogram<detail::axes_select<Axes, mpl::vector<mpl::int_<N>, Rest...>>, Storage> {
-    using HR = static_histogram<detail::axes_select<Axes, mpl::vector<mpl::int_<N>, Rest...>>, Storage>;
+  auto reduce_to(mpl::int_<N>, Rest...) const -> static_histogram<
+      detail::axes_select<Axes, mpl::vector<mpl::int_<N>, Rest...>>, Storage> {
+    using HR = static_histogram<
+        detail::axes_select<Axes, mpl::vector<mpl::int_<N>, Rest...>>, Storage>;
     typename HR::axes_type axes;
     detail::axes_assign_subset<mpl::vector<mpl::int_<N>, Rest...>>(axes, axes_);
     auto hr = HR(std::move(axes));
-    const auto b = detail::bool_mask<mpl::vector<mpl::int_<N>, Rest...>>(dim(), true);
+    const auto b =
+        detail::bool_mask<mpl::vector<mpl::int_<N>, Rest...>>(dim(), true);
     reduce_impl(hr, b);
     return hr;
   }
 
+  value_iterator begin() const noexcept {
+    return value_iterator(*this, storage_);
+  }
+
+  value_iterator end() const noexcept {
+    return value_iterator(storage_);
+  }
+
 private:
   axes_type axes_;
   Storage storage_;
 
   std::size_t bincount_from_axes() const noexcept {
-    detail::field_count fc;
+    detail::field_count_visitor fc;
     fusion::for_each(axes_, fc);
     return fc.value;
   }
@@ -279,7 +290,7 @@ private:
     double w = 0.0;
     xlin_w<0>(idx, stride, w, std::forward<Args>(args)...);
     if (stride)
-      storage_.weighted_increase(idx, w);
+      storage_.increase_by_weight(idx, w);
   }
 
   template <unsigned D> inline void lin(std::size_t &, std::size_t &) const {}
@@ -341,7 +352,7 @@ private:
     return xlin_n<D>(idx, stride, x, std::forward<Rest>(rest)...);
   }
 
-  struct shape_assign_helper {
+  struct shape_assign_visitor {
     mutable std::vector<unsigned>::iterator ni;
     template <typename Axis> void operator()(const Axis &a) const {
       *ni = a.shape();
@@ -352,8 +363,7 @@ private:
   template <typename H>
   void reduce_impl(H &h, const std::vector<bool> &b) const {
     std::vector<unsigned> n(dim());
-    auto helper = shape_assign_helper{n.begin()};
-    for_each_axis(helper);
+    for_each_axis(shape_assign_visitor{n.begin()});
     detail::index_mapper m(n, b);
     do {
       detail::storage_add(h.storage_, storage_, m.second, m.first);
diff --git a/include/boost/histogram/storage/adaptive_storage.hpp b/include/boost/histogram/storage/adaptive_storage.hpp
index 99c496bc..b08a1030 100644
--- a/include/boost/histogram/storage/adaptive_storage.hpp
+++ b/include/boost/histogram/storage/adaptive_storage.hpp
@@ -10,7 +10,7 @@
 #include <algorithm>
 #include <boost/cstdint.hpp>
 #include <boost/histogram/detail/meta.hpp>
-#include <boost/histogram/detail/weight_counter.hpp>
+#include <boost/histogram/storage/weight_counter.hpp>
 #include <boost/mpl/int.hpp>
 #include <boost/multiprecision/cpp_int.hpp>
 #include <boost/variant.hpp>
@@ -18,11 +18,10 @@
 #include <memory>
 #include <type_traits>
 #ifdef BOOST_HISTOGRAM_TRACE_ALLOCS
-#include <iostream>
 #include <boost/core/typeinfo.hpp>
+#include <iostream>
 #endif
 
-
 // forward declaration for serialization
 namespace boost {
 namespace serialization {
@@ -42,13 +41,14 @@ namespace histogram {
 
 namespace detail {
 
-using mp_int = multiprecision::cpp_int;
+using mp_int = ::boost::multiprecision::cpp_int;
+using weight_counter = ::boost::histogram::weight_counter<double>;
 
-template <typename T>
-inline T* alloc(std::size_t s) {
+template <typename T> inline T *alloc(std::size_t s) {
 #ifdef BOOST_HISTOGRAM_TRACE_ALLOCS
-  boost::core::typeinfo const & ti = BOOST_CORE_TYPEID(T);
-  std::cerr << "alloc " << boost::core::demangled_name( ti ) << "[" << s << "]" << std::endl;
+  boost::core::typeinfo const &ti = BOOST_CORE_TYPEID(T);
+  std::cerr << "alloc " << boost::core::demangled_name(ti) << "[" << s << "]"
+            << std::endl;
 #endif
   return new T[s];
 }
@@ -230,17 +230,19 @@ struct wincrease_visitor : public static_visitor<void> {
 
   template <typename T> void operator()(array<T> &lhs) const {
     array<weight_counter> a(lhs);
-    a[idx] += rhs;
+    a[idx].increase_by_weight(rhs);
     lhs_any = std::move(a);
   }
 
   void operator()(array<void> &lhs) const {
     array<weight_counter> a(lhs.size);
-    a[idx] += rhs;
+    a[idx].increase_by_weight(rhs);
     lhs_any = std::move(a);
   }
 
-  void operator()(array<weight_counter> &lhs) const { lhs[idx] += rhs; }
+  void operator()(array<weight_counter> &lhs) const {
+    lhs[idx].increase_by_weight(rhs);
+  }
 };
 
 struct value_visitor : public static_visitor<double> {
@@ -253,7 +255,9 @@ struct value_visitor : public static_visitor<double> {
 
   double operator()(const array<void> & /*b*/) const { return 0; }
 
-  double operator()(const array<weight_counter> &b) const { return b[idx].w; }
+  double operator()(const array<weight_counter> &b) const {
+    return b[idx].value();
+  }
 };
 
 struct variance_visitor : public static_visitor<double> {
@@ -266,7 +270,9 @@ struct variance_visitor : public static_visitor<double> {
 
   double operator()(const array<void> & /*b*/) const { return 0; }
 
-  double operator()(const array<weight_counter> &b) const { return b[idx].w2; }
+  double operator()(const array<weight_counter> &b) const {
+    return b[idx].variance();
+  }
 };
 
 template <typename RHS> struct radd_visitor : public static_visitor<void> {
@@ -294,7 +300,8 @@ template <typename RHS> struct radd_visitor : public static_visitor<void> {
     lhs[idx] += static_cast<mp_int>(rhs);
   }
 
-  void operator()(array<weight_counter> &lhs) const { lhs[idx] += rhs; }
+  void operator()(array<weight_counter> &lhs) const {
+    lhs[idx].increase_by_count(static_cast<double>(rhs)); }
 };
 
 template <> struct radd_visitor<weight_counter> : public static_visitor<void> {
@@ -384,7 +391,7 @@ public:
   adaptive_storage(adaptive_storage &&) = default;
   adaptive_storage &operator=(adaptive_storage &&) = default;
 
-  template <typename RHS, typename = detail::is_storage<RHS>>
+  template <typename RHS>
   explicit adaptive_storage(const RHS &rhs)
       : buffer_(detail::array<void>(rhs.size())) {
     using T = typename RHS::value_type;
@@ -432,8 +439,9 @@ public:
     }
   }
 
-  void weighted_increase(std::size_t i, const double weight_counter) {
-    apply_visitor(detail::wincrease_visitor(buffer_, i, weight_counter), buffer_);
+  void increase_by_weight(std::size_t i, const value_type weight_counter) {
+    apply_visitor(detail::wincrease_visitor(buffer_, i, weight_counter),
+                  buffer_);
   }
 
   value_type value(std::size_t i) const {
@@ -460,7 +468,7 @@ public:
   }
 
   // precondition: storages have same size
-  template <typename RHS, typename = detail::is_storage<RHS>>
+  template <typename RHS>
   adaptive_storage &operator+=(const RHS &rhs) {
     for (auto i = 0ul, n = size(); i < n; ++i)
       apply_visitor(detail::radd_visitor<typename RHS::value_type>(
diff --git a/include/boost/histogram/storage/array_storage.hpp b/include/boost/histogram/storage/array_storage.hpp
index 63444213..4620194a 100644
--- a/include/boost/histogram/storage/array_storage.hpp
+++ b/include/boost/histogram/storage/array_storage.hpp
@@ -23,12 +23,19 @@ class access;
 namespace boost {
 namespace histogram {
 
-namespace detail {} // namespace detail
+namespace detail {
+template <typename T>
+struct counter_traits {
+  using value_type = T;
+  static value_type value(const T& t) noexcept { return t; }
+  static void increase_by_count(T& lhs, const T& n) noexcept { lhs += n; }
+  static void increase_by_weight(T& lhs, const T& w) noexcept { lhs += w; }
+};
+} // namespace detail
 
 template <typename T> class array_storage {
 public:
-  static_assert(is_arithmetic<T>::value, "T must support arithmetic operators");
-  using value_type = T;
+  using value_type = typename detail::counter_traits<T>::value_type;
 
   explicit array_storage(std::size_t s) { init(s); }
 
@@ -56,7 +63,7 @@ public:
     return *this;
   }
 
-  template <typename S, typename = detail::is_storage<S>>
+  template <typename S>
   explicit array_storage(const S &other) {
     reset(other.size());
     for (decltype(size_) i = 0u; i < size_; ++i) {
@@ -73,10 +80,30 @@ public:
   }
 
   std::size_t size() const noexcept { return size_; }
-  void increase(std::size_t i) noexcept { ++array_[i]; }
-  void add(std::size_t i, const value_type &n) noexcept { array_[i] += n; }
 
-  value_type value(std::size_t i) const noexcept { return array_[i]; }
+  void increase(std::size_t i) noexcept { ++array_[i]; }
+
+  void increase_by_weight(std::size_t i, const value_type w) noexcept {
+    detail::counter_traits<T>::increase_by_weight(array_[i], w);
+  }
+
+  void add(std::size_t i, const value_type &n) noexcept {
+    detail::counter_traits<T>::increase_by_count(array_[i], n);
+  }
+
+  template <typename U=T, typename = decltype(U(0, 0))>
+  void add(std::size_t i, const value_type& value, const value_type& variance) noexcept {
+    array_[i] += T(value, variance);
+  }
+
+  value_type value(std::size_t i) const noexcept {
+    return detail::counter_traits<T>::value(array_[i]);
+  }
+
+  template <typename U=T, typename = decltype(std::declval<U &>().variance())>
+  value_type variance(std::size_t i) const noexcept {
+    return array_[i].variance();
+  }
 
   template <typename U>
   array_storage &operator+=(const array_storage<U> &rhs) noexcept {
diff --git a/include/boost/histogram/storage/operators.hpp b/include/boost/histogram/storage/operators.hpp
index a58f278d..5554cb20 100644
--- a/include/boost/histogram/storage/operators.hpp
+++ b/include/boost/histogram/storage/operators.hpp
@@ -12,53 +12,47 @@
 namespace boost {
 namespace histogram {
 namespace detail {
-  template <typename S1, typename S2>
-  bool equal_impl(std::true_type, std::true_type, const S1& s1, const S2& s2)
-  {
-    for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
-      if (s1.value(i) != s2.value(i) || s1.variance(i) != s2.variance(i))
-        return false;
-    return true;
-  }
-
-  template <typename S1, typename S2>
-  bool equal_impl(std::true_type, std::false_type, const S1& s1, const S2& s2)
-  {
-    for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
-      if (s1.value(i) != s2.value(i) || s1.value(i) != s1.variance(i))
-        return false;
-    return true;
-  }
-
-  template <typename S1, typename S2>
-  bool equal_impl(std::false_type, std::true_type, const S1& s1, const S2& s2)
-  {
-    return equal_impl(std::true_type(), std::false_type(), s2, s1);
-  }
-
-  template <typename S1, typename S2>
-  bool equal_impl(std::false_type, std::false_type, const S1& s1, const S2& s2)
-  {
-    for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
-      if (s1.value(i) != s2.value(i))
-        return false;
-    return true;
-  }
+template <typename S1, typename S2>
+bool equal_impl(std::true_type, std::true_type, const S1 &s1, const S2 &s2) {
+  for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
+    if (s1.value(i) != s2.value(i) || s1.variance(i) != s2.variance(i))
+      return false;
+  return true;
 }
 
-template <typename S1, typename S2, typename = detail::is_storage<S1>,
-          typename = detail::is_storage<S2>>
+template <typename S1, typename S2>
+bool equal_impl(std::true_type, std::false_type, const S1 &s1, const S2 &s2) {
+  for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
+    if (s1.value(i) != s2.value(i) || s1.value(i) != s1.variance(i))
+      return false;
+  return true;
+}
+
+template <typename S1, typename S2>
+bool equal_impl(std::false_type, std::true_type, const S1 &s1, const S2 &s2) {
+  return equal_impl(std::true_type(), std::false_type(), s2, s1);
+}
+
+template <typename S1, typename S2>
+bool equal_impl(std::false_type, std::false_type, const S1 &s1, const S2 &s2) {
+  for (decltype(s1.size()) i = 0, n = s1.size(); i < n; ++i)
+    if (s1.value(i) != s2.value(i))
+      return false;
+  return true;
+}
+} // namespace detail
+
+template <typename S1, typename S2, typename = detail::requires_storage<S1>,
+          typename = detail::requires_storage<S2>>
 bool operator==(const S1 &s1, const S2 &s2) noexcept {
   if (s1.size() != s2.size())
     return false;
-  return detail::equal_impl(
-    detail::has_variance_support_t<S1>(),
-    detail::has_variance_support_t<S2>(),
-    s1, s2);
+  return detail::equal_impl(detail::has_variance_support_t<S1>(),
+                            detail::has_variance_support_t<S2>(), s1, s2);
 }
 
-template <typename S1, typename S2, typename = detail::is_storage<S1>,
-          typename = detail::is_storage<S2>>
+template <typename S1, typename S2, typename = detail::requires_storage<S1>,
+          typename = detail::requires_storage<S2>>
 bool operator!=(const S1 &s1, const S2 &s2) noexcept {
   return !operator==(s1, s2);
 }
diff --git a/include/boost/histogram/storage/weight_counter.hpp b/include/boost/histogram/storage/weight_counter.hpp
new file mode 100644
index 00000000..0f79a09c
--- /dev/null
+++ b/include/boost/histogram/storage/weight_counter.hpp
@@ -0,0 +1,120 @@
+// Copyright 2015-2016 Hans Dembinski
+//
+// Distributed under the Boost Software License, version 1.0.
+// (See accompanying file LICENSE_1_0.txt
+// or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+#ifndef _BOOST_HISTOGRAM_STORAGE_WEIGHT_COUNTER_HPP_
+#define _BOOST_HISTOGRAM_STORAGE_WEIGHT_COUNTER_HPP_
+
+namespace boost {
+
+namespace serialization { class access; } // forward declare
+
+namespace histogram {
+
+/// Double counter which holds a sum of weights and a sum of squared weights
+template <typename RealType> class weight_counter {
+public:
+  /// Beware: For performance reasons counters are not initialized
+  weight_counter() = default;
+  weight_counter(const weight_counter &) = default;
+  weight_counter(weight_counter &&) = default;
+  weight_counter &operator=(const weight_counter &) = default;
+  weight_counter &operator=(weight_counter &&) = default;
+
+  weight_counter(RealType value, RealType variance) : w(value), w2(variance) {}
+
+  template <typename T>
+  explicit weight_counter(const T &t)
+      : w(static_cast<RealType>(t)), w2(static_cast<RealType>(t)) {}
+
+  template <typename T> weight_counter &operator=(const T &t) {
+    w = static_cast<RealType>(t);
+    w2 = static_cast<RealType>(t);
+    return *this;
+  }
+
+  weight_counter &operator++() {
+    ++w;
+    ++w2;
+    return *this;
+  }
+
+  template <typename U>
+  weight_counter &operator+=(const weight_counter<U> &rhs) {
+    w += rhs.w;
+    w2 += rhs.w2;
+    return *this;
+  }
+
+  weight_counter &operator*=(const RealType x) {
+    w *= x;
+    w2 *= x * x;
+    return *this;
+  }
+
+  weight_counter &increase_by_count(const RealType n) {
+    w += n;
+    w2 += n;
+    return *this;
+  }
+
+  weight_counter &increase_by_weight(const RealType x) {
+    w += x;
+    w2 += x * x;
+    return *this;
+  }
+
+  template <typename U>
+  bool operator==(const weight_counter<U> &rhs) const {
+    return w == rhs.w && w2 == rhs.w2;
+  }
+
+  template <typename U>
+  bool operator!=(const weight_counter<U> &rhs) const { return !operator==(rhs); }
+
+  template <typename T> bool operator==(const T &rhs) const {
+    return w == w2 && w == static_cast<RealType>(rhs);
+  }
+
+  template <typename T> bool operator!=(const T &rhs) const {
+    return !operator==(rhs);
+  }
+
+  RealType value() const noexcept { return w; }
+  RealType variance() const noexcept { return w2; }
+
+private:
+  friend class ::boost::serialization::access;
+
+  template <class Archive> void serialize(Archive &ar, unsigned /* version */);
+
+  RealType w, w2;
+};
+
+template <typename T, typename U> bool operator==(const T &t, const weight_counter<U> &w) {
+  return w == t;
+}
+
+template <typename T, typename U> bool operator!=(const T &t, const weight_counter<U> &w) {
+  return !(w == t);
+}
+
+namespace detail {
+  template <typename T> struct counter_traits; // fwd declaration
+
+  // specialization
+  template <typename T>
+  struct counter_traits<weight_counter<T>> {
+    using value_type = T;
+    static value_type value(const weight_counter<T>& w) noexcept { return w.value(); }
+    static void increase_by_count(weight_counter<T>& lhs, const T& n) noexcept { lhs.increase_by_count(n); }
+    static void increase_by_weight(weight_counter<T>& lhs, const T& w) noexcept { lhs.increase_by_weight(w); }
+  };
+}
+
+} // namespace histogram
+} // namespace boost
+
+#endif
diff --git a/include/boost/histogram/value_iterator.hpp b/include/boost/histogram/value_iterator.hpp
new file mode 100644
index 00000000..f6f045cc
--- /dev/null
+++ b/include/boost/histogram/value_iterator.hpp
@@ -0,0 +1,82 @@
+// Copyright 2015-2016 Hans Dembinski
+//
+// Distributed under the Boost Software License, Version 1.0.
+// (See accompanying file LICENSE_1_0.txt
+// or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+#ifndef _BOOST_HISTOGRAM_VALUE_ITERATOR_HPP_
+#define _BOOST_HISTOGRAM_VALUE_ITERATOR_HPP_
+
+#include <boost/iterator/iterator_facade.hpp>
+#include <limits>
+#include <vector>
+
+namespace boost {
+namespace histogram {
+
+template <typename Storage>
+class value_iterator
+    : public iterator_facade<
+          value_iterator<Storage>, typename Storage::value_type,
+          forward_traversal_tag, typename Storage::value_type> {
+  struct dim_t {
+    int idx, size;
+    std::size_t stride;
+  };
+
+  struct dim_visitor {
+    mutable std::size_t stride;
+    std::vector<dim_t> &dims;
+    template <typename Axis> void operator()(const Axis &a) const {
+      dims.emplace_back(dim_t{0, a.size(), stride});
+      stride *= a.shape();
+    }
+  };
+
+public:
+  /// begin iterator
+  template <typename Histogram>
+  value_iterator(const Histogram &h, const Storage &s) : s_(s), idx_(0) {
+    dims_.reserve(h.dim());
+    h.for_each_axis(dim_visitor{1, dims_});
+  }
+
+  /// end iterator
+  explicit value_iterator(const Storage &s)
+      : s_(s), idx_(std::numeric_limits<std::size_t>::max()) {}
+
+  value_iterator(const value_iterator &) = default;
+  value_iterator &operator=(const value_iterator &) = default;
+
+  int idx(int dim = 0) const noexcept { return dims_[dim].idx; }
+
+private:
+  void increment() noexcept {
+    auto iter = dims_.begin();
+    while (iter->idx == (iter->size - 1)) {
+      iter->idx = 0;
+      ++iter;
+    }
+    if (iter == dims_.end())
+      idx_ = std::numeric_limits<std::size_t>::max();
+    else {
+      ++(iter->idx);
+      idx_ = 0;
+      for (; iter != dims_.end(); ++iter)
+        idx_ += iter->idx * iter->stride;
+    }
+  }
+  bool equal(const value_iterator &other) const noexcept {
+    return &s_ == &(other.s_) && idx_ == other.idx_;
+  }
+  typename Storage::value_type dereference() const { return s_.value(idx_); }
+
+  const Storage &s_;
+  std::size_t idx_;
+  std::vector<dim_t> dims_;
+  friend class ::boost::iterator_core_access;
+};
+} // namespace histogram
+} // namespace boost
+
+#endif
\ No newline at end of file
diff --git a/test/adaptive_storage_test.cpp b/test/adaptive_storage_test.cpp
index 498172d0..b10e9b6a 100644
--- a/test/adaptive_storage_test.cpp
+++ b/test/adaptive_storage_test.cpp
@@ -17,7 +17,7 @@ namespace boost {
 
 namespace histogram {
 
-template <typename T> adaptive_storage prepare(unsigned n = 1) {
+template <typename T> adaptive_storage prepare(std::size_t n = 1) {
   adaptive_storage s(n);
   s.increase(0);
   const auto tmax = std::numeric_limits<T>::max();
@@ -27,18 +27,18 @@ template <typename T> adaptive_storage prepare(unsigned n = 1) {
   return s;
 }
 
-template <> adaptive_storage prepare<void>(unsigned n) {
+template <> adaptive_storage prepare<void>(std::size_t n) {
   adaptive_storage s(n);
   return s;
 }
 
-template <> adaptive_storage prepare<detail::weight_counter>(unsigned n) {
+template <> adaptive_storage prepare<detail::weight_counter>(std::size_t n) {
   adaptive_storage s(n);
-  s.weighted_increase(0, 1.0);
+  s.increase_by_weight(0, 1.0);
   return s;
 }
 
-template <> adaptive_storage prepare<detail::mp_int>(unsigned n) {
+template <> adaptive_storage prepare<detail::mp_int>(std::size_t n) {
   adaptive_storage s(n);
   s.increase(0);
   auto tmax = static_cast<double>(std::numeric_limits<uint64_t>::max());
@@ -56,7 +56,7 @@ namespace python { // cheating to get access
 class access {
 public:
   template <typename T>
-  static histogram::adaptive_storage set_value(unsigned n, T x) {
+  static histogram::adaptive_storage set_value(std::size_t n, T x) {
     histogram::adaptive_storage s = histogram::prepare<T>(n);
     get<histogram::detail::array<T>>(s.buffer_)[0] = x;
     return s;
@@ -104,7 +104,7 @@ template <typename T> void serialization_impl() {
 }
 
 template <> void serialization_impl<void>() {
-  adaptive_storage a(1);
+  adaptive_storage a(std::size_t(1));
   std::ostringstream os;
   std::string buf;
   {
@@ -124,7 +124,7 @@ template <> void serialization_impl<void>() {
 }
 
 template <typename T> void equal_impl() {
-  adaptive_storage a(1);
+  adaptive_storage a(std::size_t(1));
   auto b = python::access::set_value(1, T(0));
   BOOST_TEST_EQ(a.value(0), 0.0);
   BOOST_TEST_EQ(a.variance(0), 0.0);
@@ -132,7 +132,7 @@ template <typename T> void equal_impl() {
   b.increase(0);
   BOOST_TEST(!(a == b));
 
-  array_storage<unsigned> c(1);
+  array_storage<unsigned> c(std::size_t(1));
   auto d = python::access::set_value(1, T(0));
   BOOST_TEST(c == d);
   c.increase(0);
@@ -140,10 +140,10 @@ template <typename T> void equal_impl() {
 }
 
 template <> void equal_impl<void>() {
-  adaptive_storage a(1);
+  adaptive_storage a(std::size_t(1));
   auto b = python::access::set_value(1, uint8_t(0));
   auto c = python::access::set_value(2, uint8_t(0));
-  auto d = array_storage<unsigned>(1);
+  auto d = array_storage<unsigned>(std::size_t(1));
   BOOST_TEST_EQ(a.value(0), 0.0);
   BOOST_TEST_EQ(a.variance(0), 0.0);
   BOOST_TEST(a == b);
@@ -169,7 +169,7 @@ template <typename T> void increase_and_grow_impl() {
   n.increase(0);
   n.increase(0);
 
-  adaptive_storage x(2);
+  adaptive_storage x(std::size_t(2));
   x.increase(0);
   n2.add(0, x.value(0));
   n2.add(0, x.value(0));
@@ -183,7 +183,7 @@ template <typename T> void increase_and_grow_impl() {
 }
 
 template <> void increase_and_grow_impl<void>() {
-  adaptive_storage s(2);
+  adaptive_storage s(std::size_t(2));
   s.increase(0);
   BOOST_TEST_EQ(s.value(0), 1.0);
   BOOST_TEST_EQ(s.value(1), 0.0);
@@ -191,7 +191,7 @@ template <> void increase_and_grow_impl<void>() {
 
 template <typename T> void convert_array_storage_impl() {
   const auto aref = python::access::set_value(1, T(0));
-  array_storage<uint8_t> s(1);
+  array_storage<uint8_t> s(std::size_t(1));
   s.increase(0);
 
   auto a = aref;
@@ -213,7 +213,7 @@ template <typename T> void convert_array_storage_impl() {
   BOOST_TEST(c == s);
   BOOST_TEST(s == c);
 
-  array_storage<float> t(1);
+  array_storage<float> t(std::size_t(1));
   t.increase(0);
   while (t.value(0) < 1e20)
     t.add(0, t.value(0));
@@ -240,7 +240,7 @@ template <typename T> void convert_array_storage_impl() {
   BOOST_TEST(g == s);
   BOOST_TEST(s == g);
 
-  array_storage<uint8_t> u(2);
+  array_storage<uint8_t> u(std::size_t(2));
   u.increase(0);
   auto h = aref;
   BOOST_TEST(!(h == u));
@@ -249,9 +249,9 @@ template <typename T> void convert_array_storage_impl() {
 }
 
 template <> void convert_array_storage_impl<void>() {
-  const auto aref = adaptive_storage(1);
+  const auto aref = adaptive_storage(std::size_t(1));
   BOOST_TEST_EQ(aref.value(0), 0.0);
-  array_storage<uint8_t> s(1);
+  array_storage<uint8_t> s(std::size_t(1));
   s.increase(0);
 
   auto a = aref;
@@ -267,7 +267,7 @@ template <> void convert_array_storage_impl<void>() {
   BOOST_TEST(c == s);
   BOOST_TEST(s == c);
 
-  array_storage<uint8_t> t(2);
+  array_storage<uint8_t> t(std::size_t(2));
   t.increase(0);
   auto d = aref;
   BOOST_TEST(!(d == t));
@@ -348,27 +348,27 @@ int main() {
 
   // add_and_grow
   {
-    adaptive_storage a(1);
+    adaptive_storage a(std::size_t(1));
     a.increase(0);
     double x = 1.0;
-    adaptive_storage y(1);
+    adaptive_storage y(std::size_t(1));
     BOOST_TEST_EQ(y.value(0), 0.0);
     a.add(0, y.value(0));
     BOOST_TEST_EQ(a.value(0), x);
     for (unsigned i = 0; i < 80; ++i) {
       a.add(0, a.value(0));
       x += x;
-      adaptive_storage b(1);
+      adaptive_storage b(std::size_t(1));
       b.add(0, a.value(0));
       BOOST_TEST_EQ(a.value(0), x);
       BOOST_TEST_EQ(a.variance(0), x);
       BOOST_TEST_EQ(b.value(0), x);
       BOOST_TEST_EQ(b.variance(0), x);
-      b.weighted_increase(0, 0.0);
+      b.increase_by_weight(0, 0.0);
       BOOST_TEST_EQ(b.value(0), x);
       BOOST_TEST_EQ(b.variance(0), x);
-      adaptive_storage c(1);
-      c.weighted_increase(0, a.value(0));
+      adaptive_storage c(std::size_t(1));
+      c.increase_by_weight(0, a.value(0));
       BOOST_TEST_EQ(c.value(0), x);
       BOOST_TEST_EQ(c.variance(0), x * x);
     }
@@ -376,7 +376,7 @@ int main() {
 
   // multiply
   {
-    adaptive_storage a(2);
+    adaptive_storage a(std::size_t(2));
     a.increase(0);
     a *= 3;
     BOOST_TEST_EQ(a.value(0), 3);
diff --git a/test/array_storage_test.cpp b/test/array_storage_test.cpp
index 908eb3c4..c9fbe72a 100644
--- a/test/array_storage_test.cpp
+++ b/test/array_storage_test.cpp
@@ -12,24 +12,23 @@
 #include <boost/histogram/storage/adaptive_storage.hpp>
 #include <boost/histogram/storage/array_storage.hpp>
 #include <boost/histogram/storage/operators.hpp>
-#include <deque>
+#include <boost/histogram/storage/weight_counter.hpp>
 #include <limits>
-#include <vector>
 
 int main() {
   using namespace boost::histogram;
 
   // ctor
   {
-    array_storage<unsigned> a(1);
+    array_storage<unsigned> a(std::size_t(1));
     BOOST_TEST_EQ(a.size(), 1u);
     BOOST_TEST_EQ(a.value(0), 0u);
   }
 
   // increase
   {
-    array_storage<unsigned> a(1), b(1);
-    array_storage<unsigned char> c(1), d(2);
+    array_storage<unsigned> a(std::size_t(1)), b(std::size_t(1));
+    array_storage<unsigned char> c(std::size_t(1)), d(std::size_t(2));
     a.increase(0);
     b.increase(0);
     c.increase(0);
@@ -49,12 +48,12 @@ int main() {
 
   // multiply
   {
-    array_storage<unsigned> a(2);
+    array_storage<unsigned> a(std::size_t(2));
     a.increase(0);
     a *= 3;
     BOOST_TEST_EQ(a.value(0), 3.0);
     BOOST_TEST_EQ(a.value(1), 0.0);
-    a.add(1, 2.0);
+    a.add(1, 2);
     BOOST_TEST_EQ(a.value(0), 3.0);
     BOOST_TEST_EQ(a.value(1), 2.0);
     a *= 3;
@@ -64,9 +63,9 @@ int main() {
 
   // copy
   {
-    array_storage<unsigned> a(1);
+    array_storage<unsigned> a(std::size_t(1));
     a.increase(0);
-    decltype(a) b(2);
+    decltype(a) b(std::size_t(2));
     BOOST_TEST(!(a == b));
     b = a;
     BOOST_TEST(a == b);
@@ -78,7 +77,7 @@ int main() {
     BOOST_TEST_EQ(c.size(), 1u);
     BOOST_TEST_EQ(c.value(0), 1u);
 
-    array_storage<unsigned char> d(1);
+    array_storage<unsigned char> d(std::size_t(1));
     BOOST_TEST(!(a == d));
     d = a;
     BOOST_TEST(a == d);
@@ -88,7 +87,7 @@ int main() {
 
   // move
   {
-    array_storage<unsigned> a(1);
+    array_storage<unsigned> a(std::size_t(1));
     a.increase(0);
     decltype(a) b;
     BOOST_TEST(!(a == b));
@@ -102,5 +101,18 @@ int main() {
     BOOST_TEST_EQ(b.size(), 0u);
   }
 
+  // with weight_counter
+  {
+    array_storage<weight_counter<double>> a(std::size_t(1));
+    a.increase(0);
+    a.add(0, 1);
+    a.add(0, 1, 0);
+    BOOST_TEST_EQ(a.value(0), 3);
+    BOOST_TEST_EQ(a.variance(0), 2);
+    a.increase_by_weight(0, 2);
+    BOOST_TEST_EQ(a.value(0), 5);
+    BOOST_TEST_EQ(a.variance(0), 6);
+  }
+
   return boost::report_errors();
 }
diff --git a/test/axis_test.cpp b/test/axis_test.cpp
index 6e0076fd..40959d89 100644
--- a/test/axis_test.cpp
+++ b/test/axis_test.cpp
@@ -8,12 +8,12 @@
 #include <boost/fusion/container/generation/make_vector.hpp>
 #include <boost/fusion/include/make_vector.hpp>
 #include <boost/fusion/include/vector.hpp>
-#include <boost/histogram/axis/axis.hpp>
 #include <boost/histogram/axis/any.hpp>
-#include <boost/histogram/histogram_fwd.hpp>
+#include <boost/histogram/axis/axis.hpp>
 #include <boost/histogram/axis/ostream_operators.hpp>
 #include <boost/histogram/detail/axis_visitor.hpp>
 #include <boost/histogram/detail/utility.hpp>
+#include <boost/histogram/histogram_fwd.hpp>
 #include <boost/math/constants/constants.hpp>
 #include <boost/variant.hpp>
 #include <limits>
diff --git a/test/detail_test.cpp b/test/detail_test.cpp
index e4c9dbd0..600b96f0 100644
--- a/test/detail_test.cpp
+++ b/test/detail_test.cpp
@@ -7,43 +7,10 @@
 #include <boost/core/lightweight_test.hpp>
 #include <boost/histogram/detail/axis_visitor.hpp>
 #include <boost/histogram/detail/utility.hpp>
-#include <boost/histogram/detail/weight_counter.hpp>
 #include <boost/variant.hpp>
-#include <cstring>
-#include <sstream>
 using namespace boost::histogram::detail;
 
-namespace boost { namespace histogram { namespace detail {
-std::ostream& operator<<(std::ostream& os, const weight_counter& w) {
-  os << "[ " << w.w << ", " << w.w2 << "]";
-  return os;
-}
-}}}
-
 int main() {
-  // weight_counter
-  {
-    BOOST_TEST_EQ(weight_counter(0), weight_counter());
-    weight_counter w(1);
-    BOOST_TEST_EQ(w, weight_counter(1));
-    BOOST_TEST_NE(w, weight_counter());
-    BOOST_TEST_EQ(1, w);
-    BOOST_TEST_EQ(w, 1);
-    BOOST_TEST_NE(2, w);
-    BOOST_TEST_NE(w, 2);
-    w += 2.0;
-    BOOST_TEST_EQ(w, weight_counter(3, 5));
-    // consistency: a weighted counter increased by weight 1 multiplied
-    // by 2 must be the same as a weighted counter increased by weight 2
-    weight_counter u(0);
-    u += 1.0;
-    u *= 2;
-    BOOST_TEST_EQ(u, weight_counter(2, 4));
-    weight_counter v(0);
-    v += 2.0;
-    BOOST_TEST_EQ(u, v);
-  }
-
   // escape0
   {
     std::ostringstream os;
diff --git a/test/histogram_test.cpp b/test/histogram_test.cpp
index a2e6994e..0b084b1d 100644
--- a/test/histogram_test.cpp
+++ b/test/histogram_test.cpp
@@ -8,16 +8,17 @@
 #include <boost/archive/text_oarchive.hpp>
 #include <boost/core/lightweight_test.hpp>
 #include <boost/histogram/histogram.hpp>
-#include <boost/histogram/ostream_operators.hpp>
 #include <boost/histogram/literals.hpp>
+#include <boost/histogram/ostream_operators.hpp>
 #include <boost/histogram/serialization.hpp>
 #include <boost/histogram/storage/adaptive_storage.hpp>
 #include <boost/histogram/storage/array_storage.hpp>
+#include <boost/histogram/storage/weight_counter.hpp>
 #include <boost/mpl/int.hpp>
+#include <cstdlib>
 #include <limits>
 #include <sstream>
 #include <vector>
-#include <cstdlib>
 
 using namespace boost::histogram;
 using namespace boost::histogram::literals; // to get _c suffix
@@ -52,12 +53,11 @@ template <typename Type> void run_tests() {
 
   // init_0
   {
-    auto h =
-        static_histogram<mpl::vector<axis::integer<>>, adaptive_storage>();
+    auto h = static_histogram<mpl::vector<axis::integer<>>, adaptive_storage>();
     BOOST_TEST_EQ(h.dim(), 1);
     BOOST_TEST_EQ(h.bincount(), 0);
     auto h2 = static_histogram<mpl::vector<axis::integer<>>,
-                        array_storage<unsigned>>();
+                               array_storage<unsigned>>();
     BOOST_TEST_EQ(h2, h);
     auto h3 =
         static_histogram<mpl::vector<axis::regular<>>, adaptive_storage>();
@@ -140,7 +140,7 @@ template <typename Type> void run_tests() {
     auto h2 = decltype(h)(h);
     BOOST_TEST(h2 == h);
     auto h3 = static_histogram<mpl::vector<axis::integer<>, axis::integer<>>,
-                        array_storage<unsigned>>(h);
+                               array_storage<unsigned>>(h);
     BOOST_TEST_EQ(h3, h);
   }
 
@@ -157,7 +157,7 @@ template <typename Type> void run_tests() {
     h2 = h2;
     BOOST_TEST_EQ(h, h2);
     auto h3 = static_histogram<mpl::vector<axis::integer<>, axis::integer<>>,
-                        array_storage<unsigned>>();
+                               array_storage<unsigned>>();
     h3 = h;
     BOOST_TEST_EQ(h, h3);
   }
@@ -185,8 +185,9 @@ template <typename Type> void run_tests() {
 
   // axis methods
   {
-    enum { A=3, B=5 };
-    auto a = make_histogram<adaptive_storage>(Type(), axis::regular<>(1, 1, 2, "foo"));
+    enum { A = 3, B = 5 };
+    auto a = make_histogram<adaptive_storage>(Type(),
+                                              axis::regular<>(1, 1, 2, "foo"));
     BOOST_TEST_EQ(a.axis().size(), 1);
     BOOST_TEST_EQ(a.axis().shape(), 3);
     BOOST_TEST_EQ(a.axis().index(1.0), 0);
@@ -326,22 +327,22 @@ template <typename Type> void run_tests() {
     auto h =
         make_histogram<adaptive_storage>(Type(), axis::regular<>(2, -1, 1));
     h.fill(0);
-    h.fill(weight(2.0), -1.0);
-    h.fill(-1.0);
-    h.fill(-2.0);
-    h.fill(weight(5.0), 10.0);
+    h.fill(weight(0.5), -1);
+    h.fill(-1);
+    h.fill(-2);
+    h.fill(weight(5), 10);
 
-    BOOST_TEST_EQ(h.sum(), 10);
+    BOOST_TEST_EQ(h.sum(), 8.5);
 
-    BOOST_TEST_EQ(h.value(-1), 1.0);
-    BOOST_TEST_EQ(h.value(0), 3.0);
-    BOOST_TEST_EQ(h.value(1), 1.0);
-    BOOST_TEST_EQ(h.value(2), 5.0);
+    BOOST_TEST_EQ(h.value(-1), 1);
+    BOOST_TEST_EQ(h.value(0), 1.5);
+    BOOST_TEST_EQ(h.value(1), 1);
+    BOOST_TEST_EQ(h.value(2), 5);
 
-    BOOST_TEST_EQ(h.variance(-1), 1.0);
-    BOOST_TEST_EQ(h.variance(0), 5.0);
-    BOOST_TEST_EQ(h.variance(1), 1.0);
-    BOOST_TEST_EQ(h.variance(2), 25.0);
+    BOOST_TEST_EQ(h.variance(-1), 1);
+    BOOST_TEST_EQ(h.variance(0), 1.25);
+    BOOST_TEST_EQ(h.variance(1), 1);
+    BOOST_TEST_EQ(h.variance(2), 25);
   }
 
   // d1w2
@@ -349,17 +350,39 @@ template <typename Type> void run_tests() {
     auto h =
         make_histogram<array_storage<float>>(Type(), axis::regular<>(2, -1, 1));
     h.fill(0);
-    h.fill(count(2.0), -1.0);
+    h.fill(weight(0.5), -1);
+    h.fill(-1);
+    h.fill(-2);
+    h.fill(count(5), 10);
+
+    BOOST_TEST_EQ(h.sum(), 8.5);
+
+    BOOST_TEST_EQ(h.value(-1), 1);
+    BOOST_TEST_EQ(h.value(0), 1.5);
+    BOOST_TEST_EQ(h.value(1), 1);
+    BOOST_TEST_EQ(h.value(2), 5);
+  }
+
+  // d1w3
+  {
+    auto h =
+        make_histogram<array_storage<weight_counter<double>>>(Type(), axis::regular<>(2, -1, 1));
+    h.fill(0);
+    h.fill(weight(0.5), -1.0);
     h.fill(-1.0);
     h.fill(-2.0);
-    h.fill(count(5.0), 10.0);
+    h.fill(count(5), 10.0);
 
-    BOOST_TEST_EQ(h.sum(), 10);
+    BOOST_TEST_EQ(h.sum(), 8.5);
 
-    BOOST_TEST_EQ(h.value(-1), 1.0);
-    BOOST_TEST_EQ(h.value(0), 3.0);
-    BOOST_TEST_EQ(h.value(1), 1.0);
-    BOOST_TEST_EQ(h.value(2), 5.0);
+    BOOST_TEST_EQ(h.value(-1), 1);
+    BOOST_TEST_EQ(h.variance(-1), 1);
+    BOOST_TEST_EQ(h.value(0), 1.5);
+    BOOST_TEST_EQ(h.variance(0), 1.25);
+    BOOST_TEST_EQ(h.value(1), 1);
+    BOOST_TEST_EQ(h.variance(1), 1);
+    BOOST_TEST_EQ(h.value(2), 5);
+    BOOST_TEST_EQ(h.variance(2), 5);
   }
 
   // d2
@@ -752,16 +775,14 @@ template <typename Type> void run_tests() {
   // custom axis
   {
     struct custom_axis : public axis::integer<> {
-      using value_type = const char*; // type that is fed to the axis
+      using value_type = const char *; // type that is fed to the axis
 
       using integer::integer; // inherit ctors of base
 
       // the customization point
       // - accept const char* and convert to int
       // - then call index method of base class
-      int index(value_type s) const {
-        return integer::index(std::atoi(s));
-      }
+      int index(value_type s) const { return integer::index(std::atoi(s)); }
     };
 
     auto h = make_histogram<adaptive_storage>(Type(), custom_axis(0, 3));
@@ -776,6 +797,47 @@ template <typename Type> void run_tests() {
     BOOST_TEST_EQ(h.value(1), 1);
     BOOST_TEST_EQ(h.value(2), 0);
   }
+
+  // value iterator
+  {
+    auto h = make_histogram<adaptive_storage>(Type(), axis::integer<>(0, 3),
+                                              axis::integer<>(2, 4));
+    h.fill(0, 2);
+    h.fill(1, 2);
+    h.fill(1, 3);
+    auto it = h.begin();
+    BOOST_TEST_EQ(it.idx(0), 0);
+    BOOST_TEST_EQ(it.idx(1), 0);
+    BOOST_TEST_EQ(*it, 1);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST_EQ(it.idx(0), 1);
+    BOOST_TEST_EQ(it.idx(1), 0);
+    BOOST_TEST_EQ(*it, 1);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST_EQ(it.idx(0), 2);
+    BOOST_TEST_EQ(it.idx(1), 0);
+    BOOST_TEST_EQ(*it, 0);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST_EQ(it.idx(0), 0);
+    BOOST_TEST_EQ(it.idx(1), 1);
+    BOOST_TEST_EQ(*it, 0);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST_EQ(it.idx(0), 1);
+    BOOST_TEST_EQ(it.idx(1), 1);
+    BOOST_TEST_EQ(*it, 1);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST_EQ(it.idx(0), 2);
+    BOOST_TEST_EQ(it.idx(1), 1);
+    BOOST_TEST_EQ(*it, 0);
+    BOOST_TEST_EQ(h.value(it.idx(0), it.idx(1)), *it);
+    ++it;
+    BOOST_TEST(it == h.end());
+  }
 }
 
 template <typename T1, typename T2> void run_mixed_tests() {
@@ -829,11 +891,11 @@ int main() {
   // using iterator ranges
   {
     auto h = make_dynamic_histogram(axis::regular<>(2, -1, 1),
-                                    axis::regular<>(2,  2, 4));
+                                    axis::regular<>(2, 2, 4));
     auto v = std::vector<double>(2);
     v = {-0.5, 2.5};
     h.fill(v.begin(), v.end());
-    v = { 0.5, 3.5};
+    v = {0.5, 3.5};
     h.fill(v.begin(), v.end());
     auto i = std::vector<int>(2);
     i = {0, 0};
diff --git a/test/meta_test.cpp b/test/meta_test.cpp
index f74882f9..54051d30 100644
--- a/test/meta_test.cpp
+++ b/test/meta_test.cpp
@@ -25,7 +25,8 @@ int main() {
     value_type variance(std::size_t) const;
   };
 
-  BOOST_TEST_EQ(typename has_variance_support<no_variance_method>::type(), false);
+  BOOST_TEST_EQ(typename has_variance_support<no_variance_method>::type(),
+                false);
   BOOST_TEST_EQ(typename has_variance_support<variance_method>::type(), true);
 
   return boost::report_errors();
diff --git a/test/weight_counter_test.cpp b/test/weight_counter_test.cpp
new file mode 100644
index 00000000..0305fa4a
--- /dev/null
+++ b/test/weight_counter_test.cpp
@@ -0,0 +1,52 @@
+// Copyright 2015-2016 Hans Dembinski
+//
+// Distributed under the Boost Software License, Version 1.0.
+// (See accompanying file LICENSE_1_0.txt
+// or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+#include <boost/core/lightweight_test.hpp>
+#include <boost/histogram/storage/weight_counter.hpp>
+#include <ostream>
+
+namespace boost {
+namespace histogram {
+template <typename T>
+std::ostream &operator<<(std::ostream &os, const weight_counter<T> &w) {
+  os << "[ " << w.value() << ", " << w.variance() << "]";
+  return os;
+}
+}
+}
+
+int main() {
+  using weight_counter = boost::histogram::weight_counter<double>;
+
+  weight_counter w(1);
+  BOOST_TEST_EQ(w, weight_counter(1));
+  BOOST_TEST_NE(w, weight_counter(0));
+  BOOST_TEST_EQ(1, w);
+  BOOST_TEST_EQ(w, 1);
+  BOOST_TEST_NE(2, w);
+  BOOST_TEST_NE(w, 2);
+
+  w.increase_by_weight(2);
+  BOOST_TEST_EQ(w.value(), 3);
+  BOOST_TEST_EQ(w.variance(), 5);
+
+  // consistency: a weighted counter increased by weight 1 multiplied
+  // by 2 must be the same as a weighted counter increased by weight 2
+  weight_counter u(0);
+  u.increase_by_weight(1);
+  u *= 2;
+  BOOST_TEST_EQ(u, weight_counter(2, 4));
+
+  weight_counter v(0);
+  v.increase_by_weight(2);
+  BOOST_TEST_EQ(u, v);
+
+  weight_counter x(0);
+  x.increase_by_count(2);
+  BOOST_TEST_EQ(x, weight_counter(2, 2));
+
+  return boost::report_errors();
+}