/*
  Copyright (C) 2018 David Robillard <d@drobilla.net>
  Copyright (C) 2006-2007 Chris Hamilton <chamilton@cs.dal.ca>

  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, see <https://www.gnu.org/licenses/>.
*/

#ifndef CHILBERT_GRAYCODERANK_HPP
#define CHILBERT_GRAYCODERANK_HPP

#include "chilbert/BigBitVec.hpp"
#include "chilbert/FixBitVec.hpp"

#include <cassert>

#define MODSPLIT(r, b, k)  \
	{                      \
		b = (k);           \
		r = b / FBV_BITS;  \
		b -= r * FBV_BITS; \
	}

namespace chilbert {

// This is the bulk of the cost in calculating
// a Compact Hilbert Index.  It compresses a previously
// calculated index when provided the rotation
// at each level of precision.
template <class H, class HC>
inline void
compactIndex(const size_t* const ms,
             const size_t* const ds,
             const size_t        n,
             const size_t        m,
             H&                  h,
             HC&                 hc)
{
	hc = 0;

	size_t hi  = 0;
	size_t hci = 0;

	// Run through the levels of precision
	for (size_t i = 0; i < m; i++) {
		// Run through the dimensions
		size_t j = ds[i];
		do {
			// This dimension contributes a bit?
			if (ms[j] > i) {
				if (testBit(h, hi)) {
					setBit(hc, hci);
				}
				++hci;
			}

			if (++j == n) {
				j = 0;
			}
			++hi;
		} while (j != ds[i]);
	}
}

template <class I>
inline void
grayCodeRank(const I& mask, const I& gi, const size_t n, I& r)
{
	r.reset();

	int      jr = 0;
	FBV_UINT jm = 1;
	int      ir = 0;
	FBV_UINT im = 1;
	for (size_t i = 0; i < n; ++i) {
		if (mask.racks()[ir] & im) {
			if (gi.racks()[ir] & im) {
				r.racks()[jr] |= jm;
			}
			jm <<= 1;
			if (jm == 0) {
				jm = 1;
				++jr;
			}
		}

		im <<= 1;
		if (im == 0) {
			im = 1;
			++ir;
		}
	}
}

template <class I>
inline void
grayCodeRankInv(const I&     mask,
                const I&     ptrn,
                const I&     r,
                const size_t n,
                const size_t b,
                I&           g,
                I&           gi)
{
	g.reset();
	gi.reset();

	size_t   ir, jr;
	FBV_UINT im, jm;

	intptr_t i = static_cast<intptr_t>(n - 1);
	MODSPLIT(ir, im, (n - 1));
	im = (FBV1 << im);

	size_t j = b - 1;
	MODSPLIT(jr, jm, j);
	jm = (FBV1 << jm);

	FBV_UINT gi0, gi1, g0;
	gi1 = gi0 = g0 = 0;

	for (; i >= 0; --i) {
		// Unconstrained bit?
		if (mask.racks()[ir] & im) {
			gi1 = gi0;
			gi0 = (r.racks()[jr] & jm) > 0;
			g0  = gi0 ^ gi1;
			if (gi0) {
				gi.racks()[ir] |= im;
			}
			if (g0) {
				g.racks()[ir] |= im;
			}
			jm >>= 1;
			if (jm == 0) {
				jm = (FBV_UINT{1}) << (FBV_BITS - 1);
				--jr;
			}
		} else {
			g0  = (ptrn.racks()[ir] & im) > 0;
			gi1 = gi0;
			gi0 = g0 ^ gi1;
			if (gi0) {
				gi.racks()[ir] |= im;
			}
			if (g0) {
				g.racks()[ir] |= im;
			}
		}

		im >>= 1;
		if (im == 0) {
			im = (FBV_UINT{1}) << (FBV_BITS - 1);
			--ir;
		}
	}
}

template <class I>
inline void
extractMask(const size_t* const ms,
            const size_t        n,
            const size_t        d,
            const size_t        i,
            I&                  mask,
            size_t&             b)
{
	assert(0 <= d && d < n);

	mask.reset();
	b = 0;

	FBV_UINT jm = 1;
	size_t   jr = 0;
	size_t   j  = d; // #D j = (d==n-1) ? 0 : d+1;
	do {
		if (ms[j] > i) {
			mask.racks()[jr] |= jm;
			++b;
		}

		jm <<= 1;
		if (jm == 0) {
			jm = 1;
			++jr;
		}
		if (++j == n) {
			j = 0;
		}
	} while (j != d);
}

} // namespace chilbert

#endif