Misc. fixes & additions

[architektonas] / src / vector.cpp

//
// vector.cpp: Various structures used for 3 dimensional imaging
//
// by James Hammons
// (C) 2006 Underground Software
//
// JLH = James L. Hammons <jlhamm@acm.org>
//
// WHO  WHEN        WHAT
// ---  ----------  ------------------------------------------------------------
// JLH  09/19/2006  Created this file
// JLH  03/22/2011  Moved implementation of constructor from header to here
// JLH  04/02/2011  Fixed divide-by-zero bug in Unit(), added Angle() function
// JLH  08/04/2013  Added Parameter() function
//

#include "vector.h"

#include <math.h>                                                               // For sqrt()
#include "mathconstants.h"

// Vector implementation

Vector::Vector(double xx/*= 0*/, double yy/*= 0*/, double zz/*= 0*/): x(xx), y(yy), z(zz)
{
}


Vector::Vector(Vector tail, Vector head): x(head.x - tail.x), y(head.y - tail.y), z(head.z - tail.z)
{
}


Vector Vector::operator=(Vector const v)
{
        x = v.x, y = v.y, z = v.z;

        return *this;
}


Vector Vector::operator+(Vector const v)
{
        return Vector(x + v.x, y + v.y, z + v.z);
}


Vector Vector::operator-(Vector const v)
{
        return Vector(x - v.x, y - v.y, z - v.z);
}


// Unary negation

Vector Vector::operator-(void)
{
        return Vector(-x, -y, -z);
}


// Vector x constant

Vector Vector::operator*(double const v)
{
        return Vector(x * v, y * v, z * v);
}


// Vector x constant

Vector Vector::operator*(float const v)
{
        return Vector(x * v, y * v, z * v);
}


// Vector / constant

Vector Vector::operator/(double const v)
{
        return Vector(x / v, y / v, z / v);
}


// Vector / constant

Vector Vector::operator/(float const v)
{
        return Vector(x / v, y / v, z / v);
}


// Vector (cross) product

Vector Vector::operator*(Vector const v)
{
        // a x b = [a2b3 - a3b2, a3b1 - a1b3, a1b2 - a2b1]
        return Vector((y * v.z) - (z * v.y), (z * v.x) - (x * v.z), (x * v.y) - (y * v.x));
}


// Dot product

double Vector::Dot(Vector const v)
{
        return (x * v.x) + (y * v.y) + (z * v.z);
}


// Vector x constant, self assigned

Vector& Vector::operator*=(double const v)
{
        x *= v, y *= v, z *= v;

        return *this;
}


// Vector / constant, self assigned

Vector& Vector::operator/=(double const v)
{
        x /= v, y /= v, z /= v;

        return *this;
}

// Vector + vector, self assigned

Vector& Vector::operator+=(Vector const v)
{
        x += v.x, y += v.y, z += v.z;

        return *this;
}


// Vector + constant, self assigned

Vector& Vector::operator+=(double const v)
{
        x += v, y += v, z += v;

        return *this;
}


// Vector - vector, self assigned

Vector& Vector::operator-=(Vector const v)
{
        x -= v.x, y -= v.y, z -= v.z;

        return *this;
}


// Vector - constant, self assigned

Vector& Vector::operator-=(double const v)
{
        x -= v, y -= v, z -= v;

        return *this;
}


// Check for equality
bool Vector::operator==(Vector const v)
{
        return (x == v.x && y == v.y && z == v.z ? true : false);
}


// Check for inequality
bool Vector::operator!=(Vector const v)
{
        return (x != v.x || y != v.y || z != v.z ? true : false);
}


Vector Vector::Unit(void)
{
        double mag = Magnitude();

        // If the magnitude of the vector is zero, then the Unit vector is undefined...
        if (mag == 0)
                return Vector(0, 0, 0);

        return Vector(x / mag, y / mag, z / mag);
}


double Vector::Magnitude(void)
{
        return sqrt(x * x + y * y + z * z);
}


double Vector::Angle(void)
{
        // acos returns a value between zero and PI, which means we don't know which
        // quadrant the angle is in... Though, if the y-coordinate of the vector is
        // negative, that means that the angle is in quadrants III - IV.
        double rawAngle = acos(Unit().x);
        double correctedAngle = (y < 0 ? (2.0 * PI) - rawAngle : rawAngle);

        return correctedAngle;
}


bool Vector::isZero(double epsilon/*= 1e-6*/)
{
        return (fabs(x) < epsilon && fabs(y) < epsilon && fabs(z) < epsilon ? true : false);
}


// Class methods

double Vector::Dot(Vector v1, Vector v2)
{
        return (v1.x * v2.x) + (v1.y * v2.y) + (v1.z * v2.z);
}


double Vector::Magnitude(Vector v1, Vector v2)
{
        double xx = v1.x - v2.x;
        double yy = v1.y - v2.y;
        double zz = v1.z - v2.z;
        return sqrt((xx * xx) + (yy * yy) + (zz * zz));
}


// Returns the parameter of a point in space to this vector. If the parameter
// is between 0 and 1, the normal of the vector to the point is on the vector.
// Note: v1 is the tail, v2 is the head of the line (vector).
double Vector::Parameter(Vector tail, Vector head, Vector p)
{
        // Geometric interpretation:
        // The parameterized point on the vector lineSegment is where the normal of
        // the lineSegment to the point intersects lineSegment. If the pp < 0, then
        // the perpendicular lies beyond the 1st endpoint. If pp > 1, then the
        // perpendicular lies beyond the 2nd endpoint.

        Vector lineSegment = head - tail;
        double magnitude = lineSegment.Magnitude();
        Vector pointSegment = p - tail;
        double t = lineSegment.Dot(pointSegment) / (magnitude * magnitude);
        return t;
}


// Return the normal to the linesegment formed by the passed in points.
// (Not sure which is head or tail, or which hand the normal lies)
// [v1 should be the tail, v2 should be the head, in which case the normal should
//  rotate anti-clockwise.]
///*static*/ Vector Vector::Normal(Vector v1, Vector v2)
/*static*/ Vector Vector::Normal(Vector tail, Vector head)
{
//      Vector v = (v1 - v2).Unit();
        Vector v = (head - tail).Unit();
        return Vector(-v.y, v.x);
}