Taper.h

/*
Hyperspace Explorer - visualizing higher-dimensional geometry
Copyright (C) 2009-2010  Lene Preuss <lene.preuss@gmail.com>

This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.

*/

#ifndef TAPER_H
#define TAPER_H

#include "Rotope.h"


template <unsigned D>
    class taper_base: public VertexData<D> {

        public:
            taper_base(const VertexData<D> &v):
                VertexData<D>(v), _previous_dim(-1), _pre_previous_dim (-1) { }


            void taper(unsigned d);

        private:
            bool alreadyTapered() { return _previous_dim >= 0; }

            int _previous_dim;
            int _pre_previous_dim;
    };


template <unsigned D, unsigned Dmin, unsigned Dmax>
    class Taper: public Taper<D, Dmin, Dmax-1> {
        public:

            Taper(const VertexData<D> &v): Taper<D, Dmin, Dmax-1>(v) {
                taper_base<D>::taper(Dmax);
            }
    };


template <unsigned D, unsigned Dmin>
    class Taper<D, Dmin, Dmin>: public taper_base<D> {
        public:

            Taper(const VertexData<D> &v): taper_base<D>(v) {
                taper_base<D>::taper(Dmin);
            }
    };

template <unsigned D> void taper_base<D>::taper(unsigned d) {

    if (d >= D) {
        throw std::logic_error(
            "taper_base::taper() called on a higher dimension than the "
            "vector space allows");
    }
    if (VertexData<D>::raw_vertices().size() < 2) {
        throw std::logic_error(
            "taper_base::taper() can only operate on at least two vertices");
    }

    Vector<D> xnew;     //  connect all vertices to this new vertex

    if (!alreadyTapered()) {
        double displacement = sqrt(1.-0.25);
        for(unsigned i = 0; i < VertexData<D>::raw_vertices().size(); ++i) {
            xnew += VertexData<D>::raw_vertices()[i];
            VertexData<D>::raw_vertices()[i][d] -= displacement;
        }
        xnew *= 1./VertexData<D>::raw_vertices().size();

        xnew[d] = displacement;
    } else {
        xnew = VertexData<D>::raw_vertices().back();
        double displacement;

        if (_pre_previous_dim < 0) {
            displacement = sqrt(xnew[_previous_dim]*xnew[_previous_dim]-.25);
            xnew[d] = displacement;

            xnew[_previous_dim] /= 2.;
        } else {
            displacement = sqrt(fabs(xnew[_previous_dim]*xnew[_previous_dim]-
                        xnew[_pre_previous_dim]*xnew[_pre_previous_dim]));
            xnew[d] = displacement;

            xnew[_previous_dim] /= 2.;
        }
        for(unsigned i = 0; i < VertexData<D>::raw_vertices().size(); ++i) {
            VertexData<D>::raw_vertices()[i][d] -= displacement;
        }

    }

    VertexData<D>::raw_vertices().push_back(xnew);

    VertexData<D>::realm() = VertexData<D>::realm().tapered(VertexData<D>::raw_vertices().size()-1);

    _pre_previous_dim = _previous_dim;
    _previous_dim = d;

    VertexData<D>::dimension()++;   //  object is now one dimension higher
}

#endif