## Saturday, February 14, 2015

### A C++ vec3 class

I'm teaching an intro graphics class as an adjunct at Westminster College (and really liking it-- it is a great environment and I would recommend it as a very good place to go as an undergrad).  The assignments are mostly ray tracing and I have been writing my solutions from scratch so I can understand what is needed for them.  My first decision was programming language.  The students are all either using Python or Java, and so I didn't want to use either of those.  I tried Swift which I am enamored with but I ran into a compiler bug early and decided some old battle-tested language would be better.  I almost did C but operator overloading is too much to give up.  So I went with C++.   In the past I have always went with the heuristic:

If it has the same representation (e.g., 3D Cartesian vectors and 3D rgb colors), that doesn't mean it is the same class.  What IS the thing?

This can be taken to an extreme where a length and a time are not floats, but are rather different types, and an operator length/time returns a velocity.  Jim Arvo and Brian Smits (and probably many others) experimented with this and found the C++ typing became a little too cumbersome to manage.  They both thought some real language support for SI units etc. might be a good idea, but rolling your own was probably too hard in practice.

On the opposite end of the spectrum are languages like GLSL which I have been using a lot lately where representation determines type and most graphics variables are vec3 or float.   If you add a color to a surface normal that is your problem.

I decided to go the GLSL route simply because I am liking it in practice, and verbose variable names make errors less likely.      So for the first time in at least a decade I wrote a new C++ vector class vec3.  First question was "float, double, or templated?"   I tried templated but it is too much typing.  For example, here is a cross product:

vec3 v = vec3::cross(u, w);   // ugly

So I just went with double (precision problems are a pain).  Some typedef REAL might be wise but I was too lazy.

I'll be curious to see if I consider this a mistake by the end of the semester!

Here's the resulting class:
```#ifndef VEC3_HPP
#define VEC3_HPP

#include

template  class Vec3
{
private:
// A Vec3 simply has three properties called x, y and z
T x, y, z;

public:
// ------------ Constructors ------------

// Default constructor
Vec3() { x = y = z = 0; };

// Three parameter constructor
Vec3(T xValue, T yValue, T zValue)
{
x = xValue;
y = yValue;
z = zValue;
}

// ------------ Getters and setters ------------

void set(const T &xValue, const T &yValue, const T &zValue)
{
x = xValue;
y = yValue;
z = zValue;
}

T getX() const { return x; }
T getY() const { return y; }
T getZ() const { return z; }

void setX(const T &xValue) { x = xValue; }
void setY(const T &yValue) { y = yValue; }
void setZ(const T &zValue) { z = zValue; }

// ------------ Helper methods ------------

// Method to reset a vector to zero
void zero()
{
x = y = z = 0;
}

// Method to normalise a vector
void normalise()
{
// Calculate the magnitude of our vector
T magnitude = sqrt((x * x) + (y * y) + (z * z));

// As long as the magnitude isn't zero, divide each element by the magnitude
// to get the normalised value between -1 and +1
if (magnitude != 0)
{
x /= magnitude;
y /= magnitude;
z /= magnitude;
}
}

// Static method to calculate and return the scalar dot product of two vectors
//
// Note: The dot product of two vectors tell us things about the angle between
// the vectors. That is, it tells us if they are pointing in the same direction
// (i.e. are they parallel? If so, the dot product will be 1), or if they're
// perpendicular (i.e. at 90 degrees to each other) the dot product will be 0,
// or if they're pointing in opposite directions then the dot product will be -1.
//
// Usage example: double foo = Vec3::dotProduct(vectorA, vectorB);
static T dotProduct(const Vec3 &vec1, const Vec3 &vec2)
{
return vec1.x * vec2.x + vec1.y * vec2.y + vec1.z * vec2.z;
}

// Non-static method to calculate and return the scalar dot product of this vector and another vector
//
// Usage example: double foo = vectorA.dotProduct(vectorB);
T dotProduct(const Vec3 &vec) const
{
return x * vec.x + y * vec.y + z * vec.z;
}

// Static method to calculate and return a vector which is the cross product of two vectors
//
// Note: The cross product is simply a vector which is perpendicular to the plane formed by
// the first two vectors. Think of a desk like the one your laptop or keyboard is sitting on.
// If you put one pencil pointing directly away from you, and then another pencil pointing to the
// right so they form a "L" shape, the vector perpendicular to the plane made by these two pencils
// points directly upwards.
//
// Whether the vector is perpendicularly pointing "up" or "down" depends on the "handedness" of the
// coordinate system that you're using.
//
// Further reading: http://en.wikipedia.org/wiki/Cross_product
//
// Usage example: Vec3 crossVect = Vec3::crossProduct(vectorA, vectorB);
static Vec3 crossProduct(const Vec3 &vec1, const Vec3 &vec2)
{
return Vec3(vec1.y * vec2.z - vec1.z * vec2.y, vec1.z * vec2.x - vec1.x * vec2.z, vec1.x * vec2.y - vec1.y * vec2.x);
}

void addX(T value) { x += value; }
void addY(T value) { y += value; }
void addZ(T value) { z += value; }

// Method to return the distance between two vectors in 3D space
//
// Note: This is accurate, but not especially fast - depending on your needs you might
// like to use the Manhattan Distance instead: http://en.wikipedia.org/wiki/Taxicab_geometry
// There's a good discussion of it here: http://stackoverflow.com/questions/3693514/very-fast-3d-distance-check
// The gist is, to find if we're within a given distance between two vectors you can use:
//
// bool within3DManhattanDistance(Vec3 c1, Vec3 c2, float distance)
// {
//      float dx = abs(c2.x - c1.x);
//      if (dx > distance) return false; // too far in x direction
//
//      float dy = abs(c2.y - c1.y);
//      if (dy > distance) return false; // too far in y direction
//
//      float dz = abs(c2.z - c1.z);
//      if (dz > distance) return false; // too far in z direction
//
//      return true; // we're within the cube
// }
//
// Or to just calculate the straight Manhattan distance you could use:
//
// float getManhattanDistance(Vec3 c1, Vec3 c2)
// {
//      float dx = abs(c2.x - c1.x);
//      float dy = abs(c2.y - c1.y);
//      float dz = abs(c2.z - c1.z);
//      return dx+dy+dz;
// }
//
static T getDistance(const Vec3 &v1, const Vec3 &v2)
{
T dx = v2.x - v1.x;
T dy = v2.y - v1.y;
T dz = v2.z - v1.z;

return sqrt(dx * dx + dy * dy + dz * dz);
}

// Method to display the vector so you can easily check the values
void display()
{
std::cout << "X: " << x << "\t Y: " << y << "\t Z: " << z << std::endl;
}

// ------------ Overloaded operators ------------

Vec3 operator+(const Vec3 &vector) const
{
return Vec3(x + vector.x, y + vector.y, z + vector.z);
}

// Overloaded add and asssign operator to add Vec3s together
void operator+=(const Vec3 &vector)
{
x += vector.x;
y += vector.y;
z += vector.z;
}

// Overloaded subtraction operator to subtract a Vec3 from another Vec3
Vec3 operator-(const Vec3 &vector) const
{
return Vec3(x - vector.x, y - vector.y, z - vector.z);
}

// Overloaded subtract and asssign operator to subtract a Vec3 from another Vec3
void operator-=(const Vec3 &vector)
{
x -= vector.x;
y -= vector.y;
z -= vector.z;
}

// Overloaded multiplication operator to multiply two Vec3s together
Vec3 operator*(const Vec3 &vector) const
{
return Vec3(x * vector.x, y * vector.y, z * vector.z);
}

// Overloaded multiply operator to multiply a vector by a scalar
Vec3 operator*(const T &value) const
{
return Vec3(x * value, y * value, z * value);
}

// Overloaded multiply and assign operator to multiply a vector by a scalar
void operator*=(const T &value)
{
x *= value;
y *= value;
z *= value;
}

// Overloaded multiply operator to multiply a vector by a scalar
Vec3 operator/(const T &value) const
{
return Vec3(x / value, y / value, z / value);
}

// Overloaded multiply and assign operator to multiply a vector by a scalar
void operator/=(const T &value)
{
x /= value;
y /= value;
z /= value;
}
};

#endif

```

Nick said...

Hi Pete, you might be interested in the header suite Apple ships with clang now. It enables GLSL/OpenCL style expressions at the C level.

Pseudonym said...

First off, you probably should be using floats, not doubles. :-)

Secondly, my classes always distinguish between vectors and points (and normals, if necessary). First off, they transform differently. But more importantly, it helps avoid simple type-related mistakes like trying to add two points together.

friedlinguini said...

If you're using GLSL as inspiration, you should definitely check out the GLM library, which aims to reproduce GLSL's types and functions. It even has some C++11 voodoo to implement swizzling without needing an explicit call and without increasing the size of the type.

Peter Shirley said...

Great leads I will definitely check those out! It appears I will stay with C family languages until I die :)

I agree on floats for perf but for a class I will always go doubles so I can do lazy programming and not have to think hard about tolerances.

Michael Taylor said...

Just a note, the code has the less than, and greater than symbols (< >) mangled, likely expecting it to be used for HTML tags.

The #include statements were the obvious case.

Peter Shirley said...

Thanks MT-- good eye! I just now put it in verbatim HTML mode. Thank you.

Serena Williams said...

Thank you, I was searching for this code.
Serena from Vector Art Service