/*
 * 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, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 */

#ifndef _GRAYCODERANK_HPP_
#define _GRAYCODERANK_HPP_


#include <Hilbert/BigBitVec.hpp>


namespace Hilbert
{
	// 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 int *ms,
		const int *ds,
		int n,
		int m,
		H &h,
		HC &hc
		)
	{
		int i, j, hr, hcr;
		FBV_UINT hm, hcm;

		hc.reset();

		hr = hcr = 0;
		hm = hcm = 1;

		// Run through the levels of precision
		for ( i = 0; i < m; i++ )
		{
			// Run through the dimensions
			j = ds[i];
			do
			{
				// This dimension contributes a bit?
				if ( ms[j] > i )
				{
					if ( h.racks()[hr] & hm )
						hc.racks()[hcr] |= hcm;
					hcm<<=1; if (hcm==0) {hcm=1;++hcr;}
				}

				if ( ++j == n ) j = 0;
				hm<<=1; if (hm==0) {hm=1;++hr;}
			}
			while ( j != ds[i] );
		}

		return;
	}

	template<class I>
	inline
	void
	grayCodeRank(
		const I &mask,
		const I &gi,
		int n,
		I &r
		)
	{
		r.reset();
		int i, ir, jr;
		FBV_UINT im, jm;

		jr = 0; jm = 1;
		ir = 0; im = 1;
		for ( 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; }
		}

		return;
	}


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

		int i, j, ir, jr;
		FBV_UINT im, jm;

		i = n-1;
		BBV_MODSPLIT(ir,im,i);
		im = (FBV1<<im);

		j = b-1;
		BBV_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; }
		}

		return;
	}

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

		mask.reset();
		b = 0;
		
		jm = 1; jr = 0;
		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 );

		return;
	}
};


#endif