libcruft-util/point.hpp

/*
 * This Source Code Form is subject to the terms of the Mozilla Public
 * License, v. 2.0. If a copy of the MPL was not distributed with this
 * file, You can obtain one at http://mozilla.org/MPL/2.0/.
 *
 * Copyright 2011-2018 Danny Robson <danny@nerdcruft.net>
 */

#ifndef CRUFT_UTIL_POINT_HPP
#define CRUFT_UTIL_POINT_HPP

#include "vector.hpp"
#include "coord.hpp"
#include "maths.hpp"
#include "view.hpp"

#include <algorithm>
#include <type_traits>

namespace cruft {
    /// An n-dimensional position in space.
    ///
    /// \tparam S number of dimensions
    /// \tparam T the underlying per-dimension datatype
    template <size_t S, typename T>
    struct point : public coord::base<S,T,point<S,T>>
    {
        using coord::base<S,T,point<S,T>>::base;

        // use a forwarding assignment operator so that we can let the base
        // take care of the many different types of parameters. otherwise we
        // have to deal with scalar, vector, initializer_list, ad nauseum.
        template <typename Arg>
        point&
        operator= (Arg&&arg)
        {
            coord::base<S,T,point<S,T>>::operator=(std::forward<Arg> (arg));
            return *this;
        }

        vector<S,T> to (point dst) const { return dst - *this; }
        vector<S,T> from (point src) const { return *this - src; }

        /// expand point to use homogenous coordinates of a higher dimension.
        /// ie, fill with (0,..,0,1)
        point<S+1,T>
        homog (void) const
        {
            return this->template redim<S+1> (1);
        }

        ///////////////////////////////////////////////////////////////////////
        static constexpr
        auto min (void)
        {
            return point { std::numeric_limits<T>::lowest () };
        }

        //-------------------------------------------------------------------
        static constexpr
        auto max (void)
        {
            return point { std::numeric_limits<T>::max () };
        }


        //-------------------------------------------------------------------
        static constexpr
        point<S,T> origin (void)
        {
            return point<S,T> {0};
        }


        ///////////////////////////////////////////////////////////////////////
        void sanity (void) const;
    };


    // Convenience typedefs
    template <typename T> using point1 = point<1,T>;
    template <typename T> using point2 = point<2,T>;
    template <typename T> using point3 = point<3,T>;
    template <typename T> using point4 = point<4,T>;

    template <size_t S> using pointi = point<S,int>;
    template <size_t S> using pointf = point<S,float>;

    typedef point1<float> point1f;
    typedef point2<float> point2f;
    typedef point3<float> point3f;
    typedef point4<float> point4f;

    typedef point2<double> point2d;
    typedef point3<double> point3d;
    typedef point4<double> point4d;

    typedef point1<unsigned> point1u;
    typedef point2<unsigned> point2u;
    typedef point3<unsigned> point3u;
    typedef point4<unsigned> point4u;

    typedef point2<int> point2i;
    typedef point3<int> point3i;
    typedef point4<int> point4i;


    ///////////////////////////////////////////////////////////////////////////
    // distance operators

    /// computes the exact euclidean distance between two points.
    template <size_t S, typename T, typename U>
    typename std::common_type<T,U>::type
    distance (point<S,T> a, point<S,U> b)
    {
        using type_t = typename std::common_type<T,U>::type;
        static_assert (std::is_floating_point<type_t>::value,
                       "sqrt likely requires fractional types");

        return std::sqrt (distance2 (a, b));
    }


    /// computes the squared euclidean distance between two points.
    ///
    /// useful if you just need to compare distances because it avoids a sqrt
    /// operation.
    template <size_t S, typename T, typename U>
    constexpr typename std::common_type<T,U>::type
    distance2 (point<S,T> a, point<S,U> b)
    {
        return sum (pow (a - b, 2));
    }

    /// computes the octile distance between two points. that is, the shortest
    /// distance between `a' and `b' where travel is only allowed beween the 8
    /// grid neighbours and cost for diagonals is proportionally larger than
    /// cardinal movement. see also: chebyshev.
    template <typename T, typename U>
    typename std::common_type<T,U>::type
    octile (point<2,T> a, point<2,U> b)
    {
        using type_t = typename std::common_type<T,U>::type;
        static_assert (!std::is_integral<type_t>::value,
                       "octile requires more than integer precision");

        const type_t D1 = 1;
        const type_t D2 = std::sqrt (type_t {2});

        auto diff = cruft::abs (a - b);

        // distance for axis-aligned walks
        auto axis = D1 * (diff.x + diff.y);

        // the savings from diagonal walks
        auto diag = (D2 - 2 * D1) * cruft::min (diff);

        return axis + diag;
    }


    /// computes the manhattan distance between two points. that is, the
    /// distance where travel is only allowed along cardinal directions.
    template <size_t S, typename T, typename U>
    constexpr typename std::common_type<T,U>::type
    manhattan (point<S,T> a, point<S,U> b)
    {
        return sum (abs (a - b));
    }


    /// computes the cheyvshev distance between two points. that is, the
    /// shortest distance between `a' and `b' where travel is only allowed
    /// between the 8 grid neighbours and cost for diagonals is the same as
    /// cardinal movement. see also: octile.
    template <size_t S, typename T, typename U>
    constexpr typename std::common_type<T,U>::type
    chebyshev (point<S,T> a, point<S,U> b)
    {
        return cruft::max (abs (a - b));
    }


    ///////////////////////////////////////////////////////////////////////////
    // returns the most distant pair of points in a set
    //
    // performance has no guarantees. in fact it's probably spectacularly slow.
    //
    // especially given we have nothing to accelerate lookups with. if you
    // want it to be fast it may be an idea to construct a bounding volume and
    // pass those vertices instead.
    template <size_t S, typename T>
    std::pair<
        cruft::point<S,T>,
        cruft::point<S,T>
    >
    furthest (cruft::view<const cruft::point<S,T>*>);


    ///////////////////////////////////////////////////////////////////////////
    /// computes the mean point across a view of points
    template <typename InputT>
    auto
    center (cruft::view<InputT> points)
    {
        CHECK_NEZ (points.size ());

        using point_type = typename std::iterator_traits<InputT>::value_type;
        using value_type = typename point_type::value_type;

        cruft::vector<point_type::elements,value_type> accum = 0;

        for (auto const &i: points)
            accum += i.template as<cruft::vector> ();

        return (accum / points.size ()).template as<cruft::point> ();
    }
}

#endif // __UTIL_POINT_HPP