// Emacs: -*- C++ -*-

static char rcsid[] = "$Id: Point.c.m4,v 1.13 1993/11/16 13:27:10 thoth Exp $";

/*
        Copyright 1993, Center for Computer Vision and Visualization,
        University of Florida.  All rights reserved.

// $Log: Point.c.m4,v $
// Revision 1.13  1993/11/16  13:27:10  thoth
// Pointcmp is now pointcmp
//
// Revision 1.12  1993/09/15  12:25:58  thoth
// Point is now a template with specializations only.
//
// Revision 1.11  93/05/28  15:43:31  thoth
// fix comment errors.
// 
// Revision 1.10  1993/05/26  16:59:57  rjj
// Added copyright notice
//
// Revision 1.9  93/04/06  22:28:17  rjj
// Changed name of file from "Point.C.m4" to "Point.c.m4"
// Fixed lattice bugs with != by defining functions to invoke ==
// Changed "const {int|double}" to "{int|double}"
// Fixed a leftover FloatPoint name ref in trunc()
// Changed bitwise unsigned parameters to int
// 
//Revision 1.8  93/04/05  23:18:00  rjj
//Many changes to support changes necessary to mutate to the
//  interface defined by newpoint.  Please see the list for
//  Point.h.m4 revison 2.10
//  Error in operator ~ () corrected
//  Lattice operator != is probably still wrong, need to check !
//  m4 defines *_*_RELATION changed to LATTICE_*_*
//  Postfix ++ and -- code added
//  file rearanged to suit my tastes, your milage my vary ;-)
//
//Revision 1.7  92/12/11  01:10:35  rjj
//Modified to support the IA_ prefix naming convention.
//Also put #ifdef IA_PRII around the #pragma stmt
//
//Revision 1.6  92/09/02  14:51:27  rjj
//Corrected code in FloatPoint( IntPoint ) to allocated array of type TYPE
//Changed enorm() cast from (float) to (double)
//
//Revision 1.5  92/09/01  21:08:17  rjj
//Changed max(POINT_TYPE,POINT_TYPE) to be supremum(POINT_TYPE,POINT_TYPE)
//Changed min(POINT_TYPE,POINT_TYPE) to be infimum(POINT_TYPE,POINT_TYPE)
//
//Revision 1.4  92/09/01  20:46:24  rjj
//Added support to ensure float op IntPoint and IntPoint op float yield FloatPoint.
//
//Revision 1.3  92/08/27  12:11:55  jnw
//Added C style comments around Log Messages to fixx gcc-cpp bug
//
//Revision 1.2  92/08/23  13:33:19  thoth
//Cfrontified
//didn't hint to inline for loop stuff any more (braindead compiler).
//hack-around for switch statements that CFront's flow analysis
//can't fully understand.
//throw changed to ia_throw.
//removed duplicate default values for formal parameters (ARM spec).
//
// Revision 1.1  92/08/20  22:24:54  rjj
// Initial revision
// 
// 
// 

*/

#include <math.h>
ifdef(`INTPOINT', `#include "IntPoint.h"')
ifdef(`FLOATPOINT', `#include "FloatPoint.h"')

#ifdef IA_PRII
#pragma implementation
#endif

////////////////////////////////////////////////////////////
// C O N S T R U C T O R S
////////////////////////////////////////////////////////////

POINT_TYPE extend_to_point (TYPE initial_value, unsigned length) {
    POINT_TYPE pt;

    pt.coordinates = new TYPE[pt.dimen = length];

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < length; i++)
	pt.coordinates[i] = initial_value;
#else
    switch (pt.dim()) {
    case 5: pt.coordinates[4] = initial_value;
    case 4: pt.coordinates[3] = initial_value;
    case 3: pt.coordinates[2] = initial_value;
    case 2: pt.coordinates[1] = initial_value;
    case 1: pt.coordinates[0] = initial_value;
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < length; i++)
		    pt.coordinates[i] = initial_value;
	    }
    }
#endif    
    return pt;
}


POINT_TYPE::POINT (const POINT_TYPE &ip)
{

// new TYPE[0] creates a pointer to one TYPE; we never use it.
// This is to keep from having to do "if (dim() == 0) ..." everywhere.
// The destructor will delete it, though.
    
    this->coordinates = new TYPE[this->dimen = ip.dim()];

#ifdef NO_LOOP_UNROLLING    
    for (unsigned i=0; i<this->dim(); i++)
	this->coordinates[i] = ip.coordinates[i];
#else
    switch (this->dim()) {
    case 5: this->coordinates[4] = ip.coordinates[4];
    case 4: this->coordinates[3] = ip.coordinates[3];
    case 3: this->coordinates[2] = ip.coordinates[2];
    case 2: this->coordinates[1] = ip.coordinates[1];
    case 1: this->coordinates[0] = ip.coordinates[0];
    case 0: break;
    default:
	    {
		for (unsigned i=0; i<this->dim(); i++)
		    this->coordinates[i] = ip.coordinates[i];
	    }
    }
#endif
}


POINT_TYPE::POINT (TYPE *vec, unsigned length, int giveaway)
{
    this->dimen = length;
    if (giveaway) {
	this->coordinates = vec;
	return;
    }
    
    this->coordinates = new TYPE[this->dim()];
    
#ifdef NO_LOOP_UNROLLING    
    for (unsigned length=0; length<this->dim(); length++)
	this->coordinates[length] = vec[length];
#else	
    switch (this->dim()) {
    case 5: this->coordinates[4] = vec[4];
    case 4: this->coordinates[3] = vec[3];
    case 3: this->coordinates[2] = vec[2];
    case 2: this->coordinates[1] = vec[1];
    case 1: this->coordinates[0] = vec[0];
    case 0: break;
    default:
	    {
		for (unsigned length=0; length<this->dim(); length++)
		    this->coordinates[length] = vec[length];
	    }
    }
#endif
}

ifdef(`FLOATPOINT',`
//  FloatPoint includes IntPoint.h so we can have a constructor
//  from IntPoint to FloatPoint: this gives us all FloatPoint/IntPoint
//  operations for free.

POINT_TYPE::POINT (const IPOINT &ipt) {

    this->coordinates = new TYPE[this->dimen = ipt.dim()];

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < ipt.dim(); i++)
	coordinates[i] = ipt[i];
#else
    switch (ipt.dim()) {
    case 5: coordinates[4] = ipt[4];
    case 4: coordinates[3] = ipt[3];
    case 3: coordinates[2] = ipt[2];
    case 2: coordinates[1] = ipt[1];
    case 1: coordinates[0] = ipt[0];
    case 0: break;     
    default:
	{
	    for (unsigned i = 0; i < ipt.dim(); i++)
		coordinates[i] = ipt[i];
	}
    }
#endif
}
')


// POINT_TYPE = concat (POINT_TYPE, POINT_TYPE)
//
// Given a combination of Points and numbers, create a new Point that is
// the concatenation of the arguments.  Lets the user construct larger
// points than the usual 5 arg constructor.

// TO DO -- These could be modified to not check for zero-length
// points (since we now know that "new foo[0]" returns a deletable
// pointer to one foo that we will never access), and also we could
// do some loop unrolling here as well.

POINT_TYPE concat(const POINT_TYPE& lhs, const POINT_TYPE& rhs)
{
    POINT_TYPE result_point;
    
    if (lhs.dim() == 0) {
	if (rhs.dim() == 0)
	    return result_point;
	else 
	    return rhs;
    } else if (rhs.dim() == 0)
	return lhs;

    unsigned int total_dimen = lhs.dim() + rhs.dim();

    result_point.coordinates = new TYPE[result_point.dimen = total_dimen];
    
    for (unsigned int i = 0; i < lhs.dim(); i++) 
	result_point.coordinates[i] = lhs.coordinates[i];

    for (unsigned int j = 0; j < rhs.dim(); j++) 
	result_point.coordinates[i+j] = rhs.coordinates[j];

    return result_point;
}


POINT_TYPE concat(const POINT_TYPE& lhs, TYPE rhs)
{
    POINT_TYPE result_point;

    if (lhs.dim() == 0) 
	return POINT_TYPE (rhs);
    
    unsigned int total_dimen = lhs.dim() + 1;

    result_point.coordinates = new TYPE[total_dimen];
    result_point.dimen = total_dimen;
    
    for (unsigned int i = 0; i < lhs.dim(); i++) 
	result_point.coordinates[i] = lhs.coordinates[i];
    
    result_point.coordinates[i] = rhs;
    
    return result_point;
}


POINT_TYPE concat(TYPE lhs, const POINT_TYPE& rhs)
{
    POINT_TYPE result_point;
    
    if (rhs.dim() == 0)
	return POINT_TYPE (lhs);

    unsigned int total_dimen = rhs.dim() + 1;

    result_point.coordinates = new TYPE[total_dimen];
    result_point.dimen = total_dimen;

    result_point.coordinates[0] = lhs;

    for (unsigned int i = 0; i < rhs.dim(); i++) 
	result_point.coordinates[i+1] = rhs.coordinates[i];

    return result_point;
}

//  cdr(POINT_TYPE, skip=1):
//
//  Given a point of dimension n, return a point made of the n-skip last
//  coordinates (or if skip is negative, the n+skip first coordinates),
//  where skip defaults to 1.


POINT_TYPE cdr(const POINT_TYPE& pt, int skip) {

    unsigned  abs_skip;

    abs_skip = (skip < 0) ? -skip : skip;
    
    if (pt.dim() < abs_skip) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    } else if (pt.dim() == skip)
	return POINT_TYPE ();

    if (skip < 0)
	return POINT_TYPE (pt.coordinates, pt.dim()+skip);
    else
	return POINT_TYPE (pt.coordinates+skip, pt.dim()-skip);
}


////////////////////////////////////////////////////////////
// M I S C    O P E R A T I O N S
////////////////////////////////////////////////////////////

POINT_TYPE &POINT_TYPE::operator = (const POINT_TYPE &point) {
    // return directly when self-assignment
    //  (not sure the if is worth it in general)

    if (this == &point)
	return *this;

    // We need to deallocate/reallocate if differing dimens.
    
    if (point.dim() != this->dim()) {
        delete[] this->coordinates;
        this->coordinates = new TYPE[this->dimen = point.dim()];
    }

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < this->dim(); i++)
	this->coordinates[i] = point.coordinates[i];
#else
    switch (this->dim()) {
    case 5: this->coordinates[4] = point.coordinates[4];
    case 4: this->coordinates[3] = point.coordinates[3];
    case 3: this->coordinates[2] = point.coordinates[2];
    case 2: this->coordinates[1] = point.coordinates[1];
    case 1: this->coordinates[0] = point.coordinates[0];
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < this->dim(); i++)
		    this->coordinates[i] = point.coordinates[i];
	    }
    }
#endif
    return *this;
}


FRIEND_SELECT_BINARY(inf, <)
FRIEND_SELECT_BINARY(sup, >)

FRIEND_SELECT_UNARY(min, <)
FRIEND_SELECT_UNARY(max, >)


ostream &operator << (ostream &stream, const POINT_TYPE& point)
{
    int last_index = point.dim() - 1;
    stream << "<";
    for  ( int i = 0; i < last_index; ++i ) {
	stream << point.coordinates[i] << " ";
    }
    if ( point.dim() > 0 ) {
	stream << point.coordinates[last_index];
    };
    stream << ">";
    return stream;
}


////////////////////////////////////////////////////////////
// I N D E X I N G    O P E R A T I O N S
////////////////////////////////////////////////////////////
// Index derangement: "POINT_TYPE (point_index, point_index..)"
//
// Construct new points from old points, possibly with different
// dimensions and indices rearanged.
// For instance, let  a=Point(10, 20, 30).
// Then b=a.derange(2,0) is equivalent to b=Point(30,10).


POINT_TYPE POINT_TYPE::derange (unsigned i0) const
{
    POINT_TYPE p(0);
    
    if (i0 >= this->dim()) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    }
    p.coordinates[0] = this->coordinates[i0];
    return p;
}


POINT_TYPE POINT_TYPE::derange (unsigned i0, unsigned i1) const
{
    POINT_TYPE p((TYPE)0, (TYPE)0);

    if (i0 >= this->dim() || i1 >= this->dim()) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    }
    p.coordinates[0] = this->coordinates[i0];
    p.coordinates[1] = this->coordinates[i1];
    return p;
}


POINT_TYPE POINT_TYPE::derange (unsigned i0, unsigned i1, unsigned i2) const
{
    POINT_TYPE p((TYPE)0, (TYPE)0, (TYPE)0);

    if (i0 >= this->dim() || i1 >= this->dim() || i2 >= this->dim()) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    }
    
    p.coordinates[0] = this->coordinates[i0];
    p.coordinates[1] = this->coordinates[i1];
    p.coordinates[2] = this->coordinates[i2];
    return p;
}

POINT_TYPE POINT_TYPE::derange (unsigned i0, unsigned i1, unsigned i2, unsigned i3) const
{
    POINT_TYPE p(0,0,0,0);

    if (i0 >= this->dim() || i1 >= this->dim() || i2 >= this->dim() || i3 >= this->dim()) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    }

    p.coordinates[0] = this->coordinates[i0];
    p.coordinates[1] = this->coordinates[i1];
    p.coordinates[2] = this->coordinates[i2];
    p.coordinates[3] = this->coordinates[i3];
    return p;
}

POINT_TYPE POINT_TYPE::derange (unsigned i0, unsigned i1, unsigned i2, unsigned i3, unsigned i4) const
{
    POINT_TYPE p(0,0,0,0,0);

    if (i0 >= this->dim() || i1 >= this->dim() || i2 >= this->dim() ||
	i3 >= this->dim() || i4 >= this->dim()) {
	ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__);
	return POINT_TYPE ();
    }

    p.coordinates[0] = this->coordinates[i0];
    p.coordinates[1] = this->coordinates[i1];
    p.coordinates[2] = this->coordinates[i2];
    p.coordinates[3] = this->coordinates[i3];
    p.coordinates[4] = this->coordinates[i4];
    return p;
}
    
POINT_TYPE POINT_TYPE::derange (unsigned *vec, const unsigned vecLen) const
{
    // We need to deallocate/reallocate if differing dimens.
    
    for (unsigned i=0; i<vecLen; i++) {
	if (vec[i] >= this->dim()) {
	    ia_throw(IA::POINT_INDEX_MISMATCH, __FILE__, __LINE__ );
	    return POINT_TYPE ();  // An empty point
	}
    }

    POINT_TYPE  point;  // Our new point
    point.coordinates = new TYPE[point.dimen = vecLen];

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < this->dim(); i++)
	this->coordinates[i] = point.coordinates[vec[i]];
#else
    switch (this->dim()) {
    case 5: this->coordinates[4] = point.coordinates[vec[4]];
    case 4: this->coordinates[3] = point.coordinates[vec[3]];
    case 3: this->coordinates[2] = point.coordinates[vec[2]];
    case 2: this->coordinates[1] = point.coordinates[vec[1]];
    case 1: this->coordinates[0] = point.coordinates[vec[0]];
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < this->dim(); i++)
		    this->coordinates[i] = point.coordinates[vec[i]];
	    }
    }
#endif
    return point;
}





////////////////////////////////////////////////////////////
// A R I T H M E T I C   O P E R A T I O N S
////////////////////////////////////////////////////////////

COMPOUND_ARITH_ASSIGN_POINT(+)
COMPOUND_ARITH_ASSIGN_POINT(-)
COMPOUND_ARITH_ASSIGN_POINT(*)
COMPOUND_ARITH_ASSIGN_POINT(/)


COMPOUND_ARITH_ASSIGN_NUMBER(+)
COMPOUND_ARITH_ASSIGN_NUMBER(-)
COMPOUND_ARITH_ASSIGN_NUMBER(*)
COMPOUND_ARITH_ASSIGN_NUMBER(/)
    

ARITH_POINT_POINT(+)
ARITH_POINT_POINT(-)
ARITH_POINT_POINT(*)
ARITH_POINT_POINT(/)
    
ARITH_NUMBER_POINT(+)
ARITH_NUMBER_POINT(-)
ARITH_NUMBER_POINT(*)
ARITH_NUMBER_POINT(/)
    
ARITH_POINT_NUMBER(+)
ARITH_POINT_NUMBER(-)
ARITH_POINT_NUMBER(*)
ARITH_POINT_NUMBER(/)
    
ifdef(`FLOATPOINT', `ARITH_FLOAT_INTPOINT(+)')
ifdef(`FLOATPOINT', `ARITH_FLOAT_INTPOINT(-)')
ifdef(`FLOATPOINT', `ARITH_FLOAT_INTPOINT(*)')
ifdef(`FLOATPOINT', `ARITH_FLOAT_INTPOINT(/)')

ifdef(`FLOATPOINT', `ARITH_INTPOINT_FLOAT(+)')
ifdef(`FLOATPOINT', `ARITH_INTPOINT_FLOAT(-)')
ifdef(`FLOATPOINT', `ARITH_INTPOINT_FLOAT(*)')
ifdef(`FLOATPOINT', `ARITH_INTPOINT_FLOAT(/)')


POINT_TYPE POINT_TYPE::operator - () const
{
    POINT_TYPE result_point;

    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return result_point;
    }

    result_point.coordinates = new TYPE[result_point.dimen = this->dim()];

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < result_point.dim(); i++)
	result_point.coordinates[i] =  -this->coordinates[i];
#else
    switch (result_point.dim()) {
    case 5: result_point.coordinates[4] = -this->coordinates[4];
    case 4: result_point.coordinates[3] = -this->coordinates[3];
    case 3: result_point.coordinates[2] = -this->coordinates[2];
    case 2: result_point.coordinates[1] = -this->coordinates[1];
    case 1: result_point.coordinates[0] = -this->coordinates[0];
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < result_point.dim(); i++)
		    result_point.coordinates[i] =  - this->coordinates[i];
	    }
    }
#endif
    return result_point;
}


// Prefix ++
POINT_TYPE& POINT_TYPE::operator ++ () {
    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return *this;
    }
    *this += 1;
    return  *this;
}


// Prefix --
POINT_TYPE& POINT_TYPE::operator -- () {
    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return *this;
    }
    *this -= 1;
    return  *this;
}


// Postfix ++
POINT_TYPE  POINT_TYPE::operator ++ (int) {
    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return *this;
    }
    POINT_TYPE  result = *this;
    *this += 1;
    return  result;
}


// Postfix --
POINT_TYPE  POINT_TYPE::operator -- (int) {
    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return *this;
    }
    POINT_TYPE  result = *this;
    *this -= 1;
    return  result;
}


// enorm: Euclidean distance from origin.

double POINT_TYPE::enorm () const
{
    double result_number = 0;
    
    for (unsigned i = 0; i < this->dim(); i++)
	result_number = result_number +
	    (double) (this->coordinates[i] * this->coordinates[i]);

    return sqrt(result_number);
}


ifdef(`FLOATPOINT', `
// Math rounding functions as friends.
//  
// All of the following routines perform some rounding function.  Except
// for irint, nint and trunc, the routines always return float points.
// irint, nint, and trunc alway pass back int points.  See rint (3)
// for details.')

ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_FLOAT(aint)')
ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_FLOAT(anint)')
ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_FLOAT(ceil)')
ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_FLOAT(floor)')
ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_FLOAT(rint)')

ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_INT(irint)')
ifdef(`FLOATPOINT', `ROUNDING_FLOAT_TO_INT(nint)')

ifdef(`FLOATPOINT', `
IPOINT trunc(const POINT_TYPE& p) 
{ 
    int	*iptr;
    iptr = new int[p.dim()];	    

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < p.dim(); i++) 
	iptr[i] = int(p.coordinates[i]);    
#else
    switch (p.dim()) {
    case 5: iptr[4] = int(p.coordinates[4]);
    case 4: iptr[3] = int(p.coordinates[3]);
    case 3: iptr[2] = int(p.coordinates[2]);
    case 2: iptr[1] = int(p.coordinates[1]);
    case 1: iptr[0] = int(p.coordinates[0]);
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < p.dim(); i++) 
		    iptr[i] = int(p.coordinates[i]);
	    }
    }
#endif    
    return IPOINT (iptr, p.dim(), 1);
}')


////////////////////////////////////////////////////////////
// L A T T I C E   O P E R A T I O N S
////////////////////////////////////////////////////////////

// Relational operators

// POINT_TYPE < POINT_TYPE    

LATTICE_POINT_POINT(<)
LATTICE_POINT_POINT(<=)
LATTICE_POINT_POINT(>)
LATTICE_POINT_POINT(>=)
LATTICE_POINT_POINT(==)

int POINT_TYPE::operator != (const POINT_TYPE& rhs) const
{
    return ! (*this == rhs);
}


// POINT_TYPE < TYPE

LATTICE_POINT_NUMBER(<)
LATTICE_POINT_NUMBER(<=)
LATTICE_POINT_NUMBER(>)
LATTICE_POINT_NUMBER(>=)
LATTICE_POINT_NUMBER(==)

int  POINT_TYPE::operator != (TYPE rhs) const
{
    return ! (*this == rhs);
}

// TYPE < POINT_TYPE 

LATTICE_NUMBER_POINT(<)
LATTICE_NUMBER_POINT(<=)
LATTICE_NUMBER_POINT(>)
LATTICE_NUMBER_POINT(>=)
LATTICE_NUMBER_POINT(==)

int operator != (TYPE lhs, const POINT_TYPE& rhs)  
{
    return !(rhs == lhs);
}


ifdef(`INTPOINT', `
// Intpoint op float
LATTICE_INTPOINT_FLOAT(<)
LATTICE_INTPOINT_FLOAT(<=)
LATTICE_INTPOINT_FLOAT(>)
LATTICE_INTPOINT_FLOAT(>=)
LATTICE_INTPOINT_FLOAT(==)

int POINT_TYPE::operator !=  (double rhs) const
{
    return ! (*this == rhs);
}


// float < INTPOINT
LATTICE_FLOAT_INTPOINT(<)
LATTICE_FLOAT_INTPOINT(<=)
LATTICE_FLOAT_INTPOINT(>)
LATTICE_FLOAT_INTPOINT(>=)
LATTICE_FLOAT_INTPOINT(==)

int operator != (double lhs, const POINT_TYPE& rhs)  
{
    return ! (rhs == lhs);
}
')

//  Point comparison: pointcmp(POINT_TYPE, POINT_TYPE)
//
// A sort of lexicographic comparison of Points, except that
// points are zero extended for comparison, e.g., <1 0> is equal
// to <1>;  <> is equal <0 0 0>, etc. Returns 0 for equality, -1
// if the first is less than the second, 1 if the first is greater
// than the second,

int pointcmp(const POINT_TYPE &lhs, const POINT_TYPE &rhs)
{
    if (lhs.dim() > rhs.dim()) {
	for (unsigned i=0; i < rhs.dim(); i++) {
	    if (lhs[i] < rhs[i])
		return -1;
	    else if (lhs[i] > rhs[i])
		return 1;
	}
	for ( ; i<lhs.dim(); i++) {
	    if (lhs[i] < 0)
		return -1;
	    else if (lhs[i] > 0)
		return 1;
	}
    } else {
	for (unsigned i=0; i < lhs.dim(); i++) {
	    if (lhs[i] < rhs[i])
		return -1;
	    else if (lhs[i] > rhs[i])
		return 1;
	}
	for ( ; i < rhs.dim(); i++) {
	    if (0 < rhs[i])
		return  -1;
	    else if (0 > rhs[i])
		return  1;
	}
    }
    return 0;

    // If this were to be used for a lexicographic sorter, then we would
    // need the following additional code:
    //
    // At this point the zero-extended Coordinates do not differ in value,
    // so check for number of dimensions; shorter Points are considered
    // less than longer Points.
    // 
    //     if (lhs.dim() == rhs.dim())
    // 		return 0;
    //     else if (lhs.dim() > rhs.dim())
    // 		return 1;
    //     else
    // 		return -1;

}


////////////////////////////////////////////////////////////
// B I T W I S E   O P E R A T I O N S
////////////////////////////////////////////////////////////

ifdef(`INTPOINT', `
POINT_TYPE POINT_TYPE::operator ~ () const
{
    POINT_TYPE result_point;
    
    if (this->dim() == 0) {
	ia_throw(IA::POINT_UNINITIALIZED, __FILE__, __LINE__);
	return *this;
    }

    result_point.coordinates = new int[result_point.dimen = this->dim()];

#ifdef NO_LOOP_UNROLLING
    for (unsigned i = 0; i < result_point.dim(); i++)
	result_point.coordinates[i] =  ~(this->coordinates[i]);
#else
    switch (result_point.dim()) {
    case 5: result_point.coordinates[4] = ~(this->coordinates[4]);
    case 4: result_point.coordinates[3] = ~(this->coordinates[3]);
    case 3: result_point.coordinates[2] = ~(this->coordinates[2]);
    case 2: result_point.coordinates[1] = ~(this->coordinates[1]);
    case 1: result_point.coordinates[0] = ~(this->coordinates[0]);
    case 0: break;
    default:
	    {
		for (unsigned i = 0; i < result_point.dim(); i++)
		    result_point.coordinates[i] = ~(this->coordinates[i]);
	    }
    }
#endif    
    return result_point;
}')


ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(%)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(<<)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(>>)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(&)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(^)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_POINT(|)')

ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(%)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(<<)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(>>)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(&)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(^)')
ifdef(`INTPOINT', `COMPOUND_BITWISE_ASSIGN_NUMBER(|)')

ifdef(`INTPOINT', `BITWISE_POINT_POINT(%)')
ifdef(`INTPOINT', `BITWISE_POINT_POINT(<<)')
ifdef(`INTPOINT', `BITWISE_POINT_POINT(>>)')
ifdef(`INTPOINT', `BITWISE_POINT_POINT(&)')
ifdef(`INTPOINT', `BITWISE_POINT_POINT(^)')
ifdef(`INTPOINT', `BITWISE_POINT_POINT(|)')

ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(%)')
ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(<<)')
ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(>>)')
ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(&)')
ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(^)')
ifdef(`INTPOINT', `BITWISE_NUMBER_POINT(|)')

ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(%)')
ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(<<)')
ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(>>)')
ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(&)')
ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(^)')
ifdef(`INTPOINT', `BITWISE_POINT_NUMBER(|)')

