/* -----------------------------------------------------------------------------
 *
 *  Copyright (C) 1997-2010 Krzysztof M. Gorski, Eric Hivon,
 *                          Benjamin D. Wandelt, Anthony J. Banday,
 *                          Matthias Bartelmann,
 *                          Reza Ansari & Kenneth M. Ganga
 *
 *
 *  This file is part of HEALPix.
 *
 *  HEALPix 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.
 *
 *  HEALPix 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 HEALPix; if not, write to the Free Software
 *  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
 *
 *  For more information about HEALPix see http://healpix.jpl.nasa.gov
 *
 *----------------------------------------------------------------------------- */
/* mk_pix2xy.c
 *
 */

/* Standard Includes */
#include <math.h>
#include "chealpix.h"
void
mk_pix2xy(int *pix2x, int *pix2y)
{

	/*
	 * =======================================================================
	 * subroutine mk_pix2xy
	 * =======================================================================
	 * constructs the array giving x and y in the face from pixel number for
	 * the nested (quad-cube like) ordering of pixels
	 *
	 * the bits corresponding to x and y are interleaved in the pixel number
	 * one breaks up the pixel number by even and odd bits
	 * =======================================================================
	 */

	int			i,
				kpix,
				jpix,
				IX,
				IY,
				IP,
				ID;

	for (i = 0; i < 1023; i++)
		pix2x[i] = 0;

	for (kpix = 0; kpix < 1024; kpix++)
	{
		jpix = kpix;
		IX = 0;
		IY = 0;
		IP = 1;
		//!bit position(in x and y)
			while (jpix != 0)
		{
			//!go through all the bits
				ID = (int) fmod(jpix, 2);
			//!bit value(in kpix), goes in ix
				jpix = jpix / 2;
			IX = ID * IP + IX;

			ID = (int) fmod(jpix, 2);
			//!bit value(in kpix), goes in iy
				jpix = jpix / 2;
			IY = ID * IP + IY;

			IP = 2 * IP;
			//!next bit(in x and y)
		}

		pix2x[kpix] = IX;
		//!in 0, 31
			pix2y[kpix] = IY;
		//!in 0, 31
	}

	/* Later */
	return;
}