/*
  Copyright (c) 2025 Alex Diener
  
  This software is provided 'as-is', without any express or implied
  warranty. In no event will the authors be held liable for any damages
  arising from the use of this software.
  
  Permission is granted to anyone to use this software for any purpose,
  including commercial applications, and to alter it and redistribute it
  freely, subject to the following restrictions:
  
  1. The origin of this software must not be misrepresented; you must not
     claim that you wrote the original software. If you use this software
     in a product, an acknowledgment in the product documentation would be
     appreciated but is not required.
  2. Altered source versions must be plainly marked as such, and must not be
     misrepresented as being the original software.
  3. This notice may not be removed or altered from any source distribution.
  
  Alex Diener alex@ludobloom.com
*/

#include "calculator/calculator_shared.h"
#include "gamemath/FixedPoint.h"
#include "gamemath/Matrix3x3f.h"
#include "gamemath/Matrix4x4f.h"
#include "gamemath/Vector2f.h"
#include "gamemath/Vector3f.h"
#include "gamemath/Vector4f.h"
#include "gamemath/Quaternionf.h"
#include "utilities/IOUtilities.h"
#include <stdarg.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <unistd.h>

static unsigned int namedMatrix3x3Count;
static struct {
	Matrix3x3f matrix;
	const char * name;
} * namedMatrices3x3;

static unsigned int namedMatrix4x4Count;
static struct {
	Matrix4x4f matrix;
	const char * name;
} * namedMatrices4x4;

static unsigned int namedQuaternionCount;
static struct {
	Quaternionf quaternion;
	const char * name;
} * namedQuaternions;

static unsigned int namedVector2Count;
static struct {
	Vector2f vector;
	const char * name;
} * namedVector2s;

static unsigned int namedVector3Count;
static struct {
	Vector3f vector;
	const char * name;
} * namedVector3s;

static unsigned int namedVector4Count;
static struct {
	Vector4f vector;
	const char * name;
} * namedVector4s;

static Matrix3x3f readMatrix3x3fArgv(int * ioArgIndex, int argc, const char ** argv) {
	Matrix3x3f matrix;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a matrix\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedMatrix3x3Count; nameIndex++) {
			if (!strcmp(namedMatrices3x3[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedMatrices3x3[nameIndex].matrix;
			}
		}
		fprintf(stderr, "Matrix with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "identity")) {
		matrix = MATRIX3x3f_IDENTITY;
		*ioArgIndex += 1;
		
	} else if (!strcmp(argv[argIndex], "-")) {
		for (int numberIndex = 0; numberIndex < 9; numberIndex++) {
			matrix.m[numberIndex] = readFloatFromFile(stdin);
		}
		*ioArgIndex += 1;
		
	} else {
		for (int numberIndex = 0; numberIndex < 9; numberIndex++) {
			matrix.m[numberIndex] = readFloatArgv(ioArgIndex, argc, argv);
		}
	}
	Matrix3x3f_transpose(&matrix);
	return matrix;
}

static void printMatrix3x3f(Matrix3x3f matrix) {
	Matrix3x3f_transpose(&matrix);
	calculatorPrintf("{%f, %f, %f,\n %f, %f, %f,\n %f, %f, %f}\n", matrix.m[0], matrix.m[1], matrix.m[2], matrix.m[3], matrix.m[4], matrix.m[5], matrix.m[6], matrix.m[7], matrix.m[8]);
}

static Matrix4x4f readMatrix4x4fArgv(int * ioArgIndex, int argc, const char ** argv) {
	Matrix4x4f matrix;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a matrix\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedMatrix4x4Count; nameIndex++) {
			if (!strcmp(namedMatrices4x4[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedMatrices4x4[nameIndex].matrix;
			}
		}
		fprintf(stderr, "Matrix with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "identity")) {
		matrix = MATRIX4x4f_IDENTITY;
		*ioArgIndex += 1;
		
	} else if (!strcmp(argv[argIndex], "-")) {
		for (int numberIndex = 0; numberIndex < 16; numberIndex++) {
			matrix.m[numberIndex] = readFloatFromFile(stdin);
		}
		*ioArgIndex += 1;
		
	} else {
		for (int numberIndex = 0; numberIndex < 16; numberIndex++) {
			matrix.m[numberIndex] = readFloatArgv(ioArgIndex, argc, argv);
		}
	}
	Matrix4x4f_transpose(&matrix);
	return matrix;
}

static void printMatrix4x4f(Matrix4x4f matrix) {
	Matrix4x4f_transpose(&matrix);
	calculatorPrintf("{%f, %f, %f, %f,\n %f, %f, %f, %f,\n %f, %f, %f, %f,\n %f, %f, %f, %f}\n", matrix.m[0], matrix.m[1], matrix.m[2], matrix.m[3], matrix.m[4], matrix.m[5], matrix.m[6], matrix.m[7], matrix.m[8], matrix.m[9], matrix.m[10], matrix.m[11], matrix.m[12], matrix.m[13], matrix.m[14], matrix.m[15]);
}

static Vector3f readVector3fArgv(int * ioArgIndex, int argc, const char ** argv) {
	Vector3f vector;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a vector3\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedVector3Count; nameIndex++) {
			if (!strcmp(namedVector3s[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedVector3s[nameIndex].vector;
			}
		}
		fprintf(stderr, "Vector3 with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "barycenter")) {
		*ioArgIndex += 1;
		Vector3f barycenter = readVector3fArgv(ioArgIndex, argc, argv);
		Vector3f vertex0 = readVector3fArgv(ioArgIndex, argc, argv);
		Vector3f vertex1 = readVector3fArgv(ioArgIndex, argc, argv);
		Vector3f vertex2 = readVector3fArgv(ioArgIndex, argc, argv);
		vector = Vector3f_fromBarycenter(barycenter, vertex0, vertex1, vertex2);
		
	} else if (!strcmp(argv[argIndex], "-")) {
		vector.x = readFloatFromFile(stdin);
		vector.y = readFloatFromFile(stdin);
		vector.z = readFloatFromFile(stdin);
		*ioArgIndex += 1;
		
	} else {
		vector.x = readFloatArgv(ioArgIndex, argc, argv);
		vector.y = readFloatArgv(ioArgIndex, argc, argv);
		vector.z = readFloatArgv(ioArgIndex, argc, argv);
	}
	return vector;
}

static Quaternionf readQuaternionfArgv(int * ioArgIndex, int argc, const char ** argv) {
	Quaternionf quaternion;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a quaternion\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedQuaternionCount; nameIndex++) {
			if (!strcmp(namedQuaternions[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedQuaternions[nameIndex].quaternion;
			}
		}
		fprintf(stderr, "Quaternion with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "identity")) {
		quaternion = QUATERNIONf_IDENTITY;
		*ioArgIndex += 1;
		
	} else if (!strcmp(argv[argIndex], "axis_angle")) {
		++*ioArgIndex;
		Vector3f axis = readVector3fArgv(ioArgIndex, argc, argv);
		float angle = readFloatAngleArgv(ioArgIndex, argc, argv);
		quaternion = Quaternionf_fromAxisAngle(axis, angle);
		
	} else if (!strcmp(argv[argIndex], "forward")) {
		++*ioArgIndex;
		Vector3f forward = readVector3fArgv(ioArgIndex, argc, argv);
		quaternion = Quaternionf_forwardDirect(forward);
		
	} else if (!strcmp(argv[argIndex], "forward_level")) {
		++*ioArgIndex;
		Vector3f forward = readVector3fArgv(ioArgIndex, argc, argv);
		quaternion = Quaternionf_forwardLevel(forward);
		
	} else if (!strcmp(argv[argIndex], "-")) {
		quaternion.x = readFloatFromFile(stdin);
		quaternion.y = readFloatFromFile(stdin);
		quaternion.z = readFloatFromFile(stdin);
		quaternion.w = readFloatFromFile(stdin);
		*ioArgIndex += 1;
		
	} else {
		quaternion.x = readFloatArgv(ioArgIndex, argc, argv);
		quaternion.y = readFloatArgv(ioArgIndex, argc, argv);
		quaternion.z = readFloatArgv(ioArgIndex, argc, argv);
		quaternion.w = readFloatArgv(ioArgIndex, argc, argv);
	}
	return quaternion;
}

static void printQuaternionf(Quaternionf quaternion) {
	calculatorPrintf("{%f, %f, %f, %f}\n", quaternion.x, quaternion.y, quaternion.z, quaternion.w);
}

static void printVector3f(Vector3f vector) {
	calculatorPrintf("{%f, %f, %f}\n", vector.x, vector.y, vector.z);
}

static Vector4f readVector4fArgv(int * ioArgIndex, int argc, const char ** argv) {
	Vector4f vector;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a vector4\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedVector4Count; nameIndex++) {
			if (!strcmp(namedVector4s[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedVector4s[nameIndex].vector;
			}
		}
		fprintf(stderr, "Vector4 with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "barycenter")) {
		*ioArgIndex += 1;
		Vector3f barycenter = readVector3fArgv(ioArgIndex, argc, argv);
		Vector4f vertex0 = readVector4fArgv(ioArgIndex, argc, argv);
		Vector4f vertex1 = readVector4fArgv(ioArgIndex, argc, argv);
		Vector4f vertex2 = readVector4fArgv(ioArgIndex, argc, argv);
		vector = Vector4f_fromBarycenter(barycenter, vertex0, vertex1, vertex2);
		
	} else if (!strcmp(argv[argIndex], "-")) {
		vector.x = readFloatFromFile(stdin);
		vector.y = readFloatFromFile(stdin);
		vector.z = readFloatFromFile(stdin);
		vector.w = readFloatFromFile(stdin);
		*ioArgIndex += 1;
		
	} else {
		vector.x = readFloatArgv(ioArgIndex, argc, argv);
		vector.y = readFloatArgv(ioArgIndex, argc, argv);
		vector.z = readFloatArgv(ioArgIndex, argc, argv);
		vector.w = readFloatArgv(ioArgIndex, argc, argv);
	}
	return vector;
}

static unsigned int swizzleVector(float * vector, unsigned int dimensions, float * outComponents, const char * letters) {
	unsigned int length = strlen(letters);
	if (length < 2 || length > 4) {
		fprintf(stderr, "Invalid swizzle: \"%s\"\n", letters);
		exit(EXIT_FAILURE);
	}
	for (unsigned int letterIndex = 0; letterIndex < length; letterIndex++) {
		switch (letters[letterIndex]) {
			case 'x':
			case 'X':
				outComponents[letterIndex] = vector[0];
				break;
			case 'y':
			case 'Y':
				outComponents[letterIndex] = vector[1];
				break;
			case 'z':
			case 'Z':
				if (dimensions < 3) {
					fprintf(stderr, "Swizzle \"%s\" attempts to reference z component in a vector with only %u dimensions\n", letters, dimensions);
					exit(EXIT_FAILURE);
				}
				outComponents[letterIndex] = vector[2];
				break;
			case 'w':
			case 'W':
				if (dimensions < 4) {
					fprintf(stderr, "Swizzle \"%s\" attempts to reference w component in a vector with only %u dimensions\n", letters, dimensions);
					exit(EXIT_FAILURE);
				}
				outComponents[letterIndex] = vector[3];
				break;
			default:
				fprintf(stderr, "Swizzle \"%s\" contains invalid characters (only x, y, z, and w are allowed)\n", letters);
				exit(EXIT_FAILURE);
		}
	}
	return length;
}

static void printVector4f(Vector4f vector) {
	calculatorPrintf("{%f, %f, %f, %f}\n", vector.x, vector.y, vector.z, vector.w);
}

static void printAxisAngle(Vector3f axis, float angle) {
	calculatorPrintf("{%f, %f, %f}, %f\n", axis.x, axis.y, axis.z, angle);
}

float readFloatArgv(int * ioArgIndex, int argc, const char ** argv) {
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting more numbers\n");
		exit(EXIT_FAILURE);
	}
	unsigned int hexValue;
	if (sscanf(argv[argIndex], "0x%X", &hexValue) == 1) {
		*ioArgIndex += 1;
		return *(float *) &hexValue;
	}
	float floatValue;
	if (sscanf(argv[argIndex], "%f", &floatValue) == 1) {
		*ioArgIndex += 1;
		return floatValue;
	}
	fprintf(stderr, "Couldn't parse \"%s\" as a number\n", argv[argIndex]);
	exit(EXIT_FAILURE);
}

float readFloatFromFile(FILE * file) {
	unsigned int hexValue;
	if (fscanf(file, "0x%X", &hexValue) == 1 || fscanf(file, " 0x%X", &hexValue) == 1) {
		return *(float *) &hexValue;
	}
	float floatValue;
	if (fscanf(file, "%f", &floatValue) == 1 || fscanf(file, " %f", &floatValue) == 1) {
		return floatValue;
	}
	fprintf(stderr, "Couldn't parse next number from input file\n");
	exit(EXIT_FAILURE);
}

float readFloatAngleArgv(int * ioArgIndex, int argc, const char ** argv) {
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting more numbers\n");
		exit(EXIT_FAILURE);
	}
	float floatValue;
	if (sscanf(argv[argIndex], "deg:%f", &floatValue) == 1) {
		*ioArgIndex += 1;
		floatValue *= M_PI / 180.0f;
		return floatValue;
	}
	if (sscanf(argv[argIndex], "pi:%f", &floatValue) == 1) {
		*ioArgIndex += 1;
		floatValue *= M_PI;
		return floatValue;
	}
	return readFloatArgv(ioArgIndex, argc, argv);
}

Vector2f readVector2fArgv(int * ioArgIndex, int argc, const char ** argv) {
	Vector2f vector;
	int argIndex = *ioArgIndex + 1;
	if (argIndex >= argc) {
		fprintf(stderr, "Reached end of argv when expecting a vector2\n");
		exit(EXIT_FAILURE);
	}
	if (argv[argIndex][0] == '@') {
		for (unsigned int nameIndex = 0; nameIndex < namedVector2Count; nameIndex++) {
			if (!strcmp(namedVector2s[nameIndex].name, argv[argIndex] + 1)) {
				*ioArgIndex += 1;
				return namedVector2s[nameIndex].vector;
			}
		}
		fprintf(stderr, "Vector2 with name %s was not found\n", argv[argIndex]);
		exit(EXIT_FAILURE);
	}
	if (!strcmp(argv[argIndex], "barycenter")) {
		*ioArgIndex += 1;
		Vector3f barycenter = readVector3fArgv(ioArgIndex, argc, argv);
		Vector2f vertex0 = readVector2fArgv(ioArgIndex, argc, argv);
		Vector2f vertex1 = readVector2fArgv(ioArgIndex, argc, argv);
		Vector2f vertex2 = readVector2fArgv(ioArgIndex, argc, argv);
		vector = Vector2f_fromBarycenter(barycenter, vertex0, vertex1, vertex2);
		
	} else if (!strcmp(argv[argIndex], "angle")) {
		*ioArgIndex += 1;
		float angle = readFloatAngleArgv(ioArgIndex, argc, argv);
		vector = VECTOR2f(cosf(angle), sinf(angle));
		
	} else if (!strcmp(argv[argIndex], "-")) {
		vector.x = readFloatFromFile(stdin);
		vector.y = readFloatFromFile(stdin);
		*ioArgIndex += 1;
		
	} else {
		vector.x = readFloatArgv(ioArgIndex, argc, argv);
		vector.y = readFloatArgv(ioArgIndex, argc, argv);
	}
	return vector;
}

void calculatorPrintf(const char * format, ...) {
	unsigned int length = strlen(format);
	va_list args;
	va_start(args, format);
	bool firstValue = true;
	for (unsigned int charIndex = 0; charIndex < length; charIndex++) {
		if (format[charIndex] == '%' && charIndex < length - 1) {
			if (outputCompact && !firstValue) {
				putchar(' ');
			}
			firstValue = false;
			if (format[charIndex + 1] == 'f') {
				double value = va_arg(args, double);
				if (outputHex) {
					float valuef = value;
					printf("0x%X", *(unsigned int *) &valuef);
				} else {
					char decimalString[DECIMAL_STRING_MAX];
					formatDecimalString(value, decimalString, sizeof(decimalString), outputFloatPrecisionMin, outputFloatPrecisionMax);
					printf("%s", decimalString);
				}
				
			} else if (format[charIndex + 1] == 'x') {
				fixed16_16 value = va_arg(args, fixed16_16);
				if (outputFloat) {
					float valuef = xtof(value);
					if (outputHex) {
						printf("0x%X", *(unsigned int *) &valuef);
					} else {
						char decimalString[DECIMAL_STRING_MAX];
						formatDecimalString(valuef, decimalString, sizeof(decimalString), outputFloatPrecisionMin, outputFloatPrecisionMax);
						printf("%s", decimalString);
					}
				} else {
					printf("0x%05X", value);
				}
			}
			charIndex++;
			
		} else if (!outputCompact) {
			putchar(format[charIndex]);
		}
	}
	va_end(args);
	if (outputCompact) {
		putchar('\n');
	}
}

void printFloat(float value) {
	calculatorPrintf("%f\n", value);
}

void printVector2f(Vector2f vector) {
	calculatorPrintf("{%f, %f}\n", vector.x, vector.y);
}

static void matrix3x3Loop(int argc, const char ** argv, int * ioArgIndex, Matrix3x3f matrix);
static void matrix4x4Loop(int argc, const char ** argv, int * ioArgIndex, Matrix4x4f matrix);
static void quaternionLoop(int argc, const char ** argv, int * ioArgIndex, Quaternionf quaternion);
static void vector2Loop(int argc, const char ** argv, int * ioArgIndex, Vector2f vector);
static void vector3Loop(int argc, const char ** argv, int * ioArgIndex, Vector3f vector);
static void vector4Loop(int argc, const char ** argv, int * ioArgIndex, Vector4f vector);

static void restartLoopWithVectorSwizzle(int argc, const char ** argv, int * ioArgIndex, float * vectorComponents, unsigned int dimensionCount) {
	float components[4];
	int argIndex = *ioArgIndex + 1;
	unsigned int componentCount = swizzleVector(vectorComponents, dimensionCount, components, argv[argIndex]);
	switch (componentCount) {
		case 2: {
			Vector2f vector2 = {components[0], components[1]};
			vector2Loop(argc, argv, &argIndex, vector2);
			break;
		}
		case 3: {
			Vector3f vector3 = {components[0], components[1], components[2]};
			vector3Loop(argc, argv, &argIndex, vector3);
			break;
		}
		case 4: {
			Vector4f vector4 = {components[0], components[1], components[2], components[3]};
			vector4Loop(argc, argv, &argIndex, vector4);
			break;
		}
	}
	*ioArgIndex = argIndex;
}

static void matrix3x3Loop(int argc, const char ** argv, int * ioArgIndex, Matrix3x3f matrix) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "multiply")) {
			Matrix3x3f matrix2 = readMatrix3x3fArgv(&argIndex, argc, argv);
			Matrix3x3f_multiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "left_multiply")) {
			Matrix3x3f matrix2 = readMatrix3x3fArgv(&argIndex, argc, argv);
			Matrix3x3f_leftMultiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "multiply_quaternion")) {
			Quaternionf quaternion = readQuaternionfArgv(&argIndex, argc, argv);
			Matrix3x3f matrix2 = Matrix3x3f_fromQuaternionf(quaternion);
			Matrix3x3f_multiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "scale")) {
			Vector3f scale = readVector3fArgv(&argIndex, argc, argv);
			Matrix3x3f_scale(&matrix, scale);
			
		} else if (!strcmp(argv[argIndex], "rotate")) {
			Vector3f axis = readVector3fArgv(&argIndex, argc, argv);
			float angle = readFloatAngleArgv(&argIndex, argc, argv);
			Matrix3x3f_rotate(&matrix, axis, angle);
			
		} else if (!strcmp(argv[argIndex], "shear_x")) {
			float y = readFloatArgv(&argIndex, argc, argv);
			float z = readFloatArgv(&argIndex, argc, argv);
			Matrix3x3f_shearX(&matrix, y, z);
			
		} else if (!strcmp(argv[argIndex], "shear_y")) {
			float x = readFloatArgv(&argIndex, argc, argv);
			float z = readFloatArgv(&argIndex, argc, argv);
			Matrix3x3f_shearY(&matrix, x, z);
			
		} else if (!strcmp(argv[argIndex], "shear_z")) {
			float x = readFloatArgv(&argIndex, argc, argv);
			float y = readFloatArgv(&argIndex, argc, argv);
			Matrix3x3f_shearZ(&matrix, x, y);
			
		} else if (!strcmp(argv[argIndex], "transpose")) {
			Matrix3x3f_transpose(&matrix);
			
		} else if (!strcmp(argv[argIndex], "interpolate")) {
			Matrix3x3f matrix2 = readMatrix3x3fArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			Matrix3x3f_interpolate(&matrix, matrix2, value);
			
		} else if (!strcmp(argv[argIndex], "normalize")) {
			Matrix3x3f_normalize(&matrix);
			
		} else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
			Matrix3x3f_invert(&matrix);
			
		} else if (!strcmp(argv[argIndex], "determinant")) {
			float determinant = Matrix3x3f_determinant(matrix);
			printFloat(determinant);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "to_quaternion")) {
			Quaternionf quaternion = Quaternionf_fromMatrix3x3f(matrix);
			*ioArgIndex = argIndex;
			quaternionLoop(argc, argv, ioArgIndex, quaternion);
			return;
			
		} else if (!strcmp(argv[argIndex], "to_matrix4x4")) {
			Matrix4x4f matrix4x4 = Matrix4x4f_fromMatrix3x3f(matrix);
			*ioArgIndex = argIndex;
			matrix4x4Loop(argc, argv, ioArgIndex, matrix4x4);
			return;
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printMatrix3x3f(matrix);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedMatrixIndex = 0; namedMatrixIndex < namedMatrix3x3Count; namedMatrixIndex++) {
				if (!strcmp(namedMatrices3x3[namedMatrixIndex].name, argv[argIndex + 1])) {
					namedMatrices3x3[namedMatrixIndex].matrix = matrix;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedMatrices3x3 = realloc(namedMatrices3x3, (namedMatrix3x3Count + 1) * sizeof(*namedMatrices3x3));
			namedMatrices3x3[namedMatrix3x3Count].name = argv[++argIndex];
			namedMatrices3x3[namedMatrix3x3Count++].matrix = matrix;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized matrix operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printMatrix3x3f(matrix);
	exit(EXIT_SUCCESS);
}

static void matrix4x4Loop(int argc, const char ** argv, int * ioArgIndex, Matrix4x4f matrix) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "multiply")) {
			Matrix4x4f matrix2 = readMatrix4x4fArgv(&argIndex, argc, argv);
			Matrix4x4f_multiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "left_multiply")) {
			Matrix4x4f matrix2 = readMatrix4x4fArgv(&argIndex, argc, argv);
			Matrix4x4f_leftMultiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "multiply_quaternion")) {
			Quaternionf quaternion = readQuaternionfArgv(&argIndex, argc, argv);
			Matrix4x4f matrix2 = Matrix4x4f_fromQuaternionf(quaternion);
			Matrix4x4f_multiply(&matrix, matrix2);
			
		} else if (!strcmp(argv[argIndex], "translate")) {
			Vector3f translation = readVector3fArgv(&argIndex, argc, argv);
			Matrix4x4f_translate(&matrix, translation);
			
		} else if (!strcmp(argv[argIndex], "scale")) {
			Vector3f scale = readVector3fArgv(&argIndex, argc, argv);
			Matrix4x4f_scale(&matrix, scale);
			
		} else if (!strcmp(argv[argIndex], "rotate")) {
			Vector3f axis = readVector3fArgv(&argIndex, argc, argv);
			float angle = readFloatAngleArgv(&argIndex, argc, argv);
			Matrix4x4f_rotate(&matrix, axis, angle);
			
		} else if (!strcmp(argv[argIndex], "shear_x")) {
			float y = readFloatArgv(&argIndex, argc, argv);
			float z = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_shearX(&matrix, y, z);
			
		} else if (!strcmp(argv[argIndex], "shear_y")) {
			float x = readFloatArgv(&argIndex, argc, argv);
			float z = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_shearY(&matrix, x, z);
			
		} else if (!strcmp(argv[argIndex], "shear_z")) {
			float x = readFloatArgv(&argIndex, argc, argv);
			float y = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_shearZ(&matrix, x, y);
			
		} else if (!strcmp(argv[argIndex], "perspective")) {
			float fovYDegrees = readFloatArgv(&argIndex, argc, argv);
			float aspect = readFloatArgv(&argIndex, argc, argv);
			float zNear = readFloatArgv(&argIndex, argc, argv);
			float zFar = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_applyPerspective(&matrix, fovYDegrees, aspect, zNear, zFar);
			
		} else if (!strcmp(argv[argIndex], "ortho")) {
			float left = readFloatArgv(&argIndex, argc, argv);
			float right = readFloatArgv(&argIndex, argc, argv);
			float bottom = readFloatArgv(&argIndex, argc, argv);
			float top = readFloatArgv(&argIndex, argc, argv);
			float zNear = readFloatArgv(&argIndex, argc, argv);
			float zFar = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_applyOrtho(&matrix, left, right, bottom, top, zNear, zFar);
			
		} else if (!strcmp(argv[argIndex], "transpose")) {
			Matrix4x4f_transpose(&matrix);
			
		} else if (!strcmp(argv[argIndex], "interpolate")) {
			Matrix4x4f matrix2 = readMatrix4x4fArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			Matrix4x4f_interpolate(&matrix, matrix2, value);
			
		} else if (!strcmp(argv[argIndex], "normalize")) {
			Matrix4x4f_normalize(&matrix);
			
		} else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
			Matrix4x4f_invert(&matrix);
			
		} else if (!strcmp(argv[argIndex], "determinant")) {
			float determinant = Matrix4x4f_determinant(matrix);
			printFloat(determinant);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "to_matrix3x3")) {
			Matrix3x3f matrix3x3 = Matrix3x3f_fromMatrix4x4f(matrix);
			*ioArgIndex = argIndex;
			matrix3x3Loop(argc, argv, ioArgIndex, matrix3x3);
			return;
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printMatrix4x4f(matrix);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedMatrixIndex = 0; namedMatrixIndex < namedMatrix4x4Count; namedMatrixIndex++) {
				if (!strcmp(namedMatrices4x4[namedMatrixIndex].name, argv[argIndex + 1])) {
					namedMatrices4x4[namedMatrixIndex].matrix = matrix;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedMatrices4x4 = realloc(namedMatrices4x4, (namedMatrix4x4Count + 1) * sizeof(*namedMatrices4x4));
			namedMatrices4x4[namedMatrix4x4Count].name = argv[++argIndex];
			namedMatrices4x4[namedMatrix4x4Count++].matrix = matrix;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized matrix operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printMatrix4x4f(matrix);
	exit(EXIT_SUCCESS);
}

static void quaternionLoop(int argc, const char ** argv, int * ioArgIndex, Quaternionf quaternion) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "multiply")) {
			Quaternionf quaternion2 = readQuaternionfArgv(&argIndex, argc, argv);
			Quaternionf_multiply(&quaternion, quaternion2);
			
    } else if (!strcmp(argv[argIndex], "left_multiply")) {
			Quaternionf quaternion2 = readQuaternionfArgv(&argIndex, argc, argv);
			Quaternionf_leftMultiply(&quaternion, quaternion2);
			
    } else if (!strcmp(argv[argIndex], "rotate")) {
			Vector3f axis = readVector3fArgv(&argIndex, argc, argv);
			float angle = readFloatAngleArgv(&argIndex, argc, argv);
			Quaternionf_rotate(&quaternion, axis, angle);
			
    } else if (!strcmp(argv[argIndex], "magnitude")) {
    	float magnitude = Quaternionf_magnitude(quaternion);
			printFloat(magnitude);
			exit(EXIT_SUCCESS);
			
    } else if (!strcmp(argv[argIndex], "normalize")) {
			Quaternionf_normalize(&quaternion);
			
    } else if (!strcmp(argv[argIndex], "slerp")) {
			Quaternionf quaternion2 = readQuaternionfArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			quaternion = Quaternionf_slerp(quaternion, quaternion2, value);
			
    } else if (!strcmp(argv[argIndex], "squad")) {
			Quaternionf quaternion2 = readQuaternionfArgv(&argIndex, argc, argv);
			Quaternionf quaternion3 = readQuaternionfArgv(&argIndex, argc, argv);
			Quaternionf quaternion4 = readQuaternionfArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			quaternion = Quaternionf_squad(quaternion, quaternion2, quaternion3, quaternion4, value);
			
    } else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
    	Quaternionf_invert(&quaternion);
			
    } else if (!strcmp(argv[argIndex], "level_horizon")) {
    	quaternion = Quaternionf_levelHorizon(quaternion);
			
    } else if (!strcmp(argv[argIndex], "to_matrix3x3")) {
			Matrix3x3f matrix = Matrix3x3f_fromQuaternionf(quaternion);
			*ioArgIndex = argIndex;
			matrix3x3Loop(argc, argv, ioArgIndex, matrix);
			return;
			
    } else if (!strcmp(argv[argIndex], "to_matrix4x4") || !strcmp(argv[argIndex], "to_matrix")) {
			Matrix4x4f matrix = Matrix4x4f_fromQuaternionf(quaternion);
			*ioArgIndex = argIndex;
			matrix4x4Loop(argc, argv, ioArgIndex, matrix);
			return;
			
    } else if (!strcmp(argv[argIndex], "to_axis_angle")) {
			Vector3f axis;
			float angle;
			Quaternionf_toAxisAngle(quaternion, &axis, &angle);
			printAxisAngle(axis, angle);
			exit(EXIT_SUCCESS);
			
    } else if (!strcmp(argv[argIndex], "get_angle")) {
			Quaternionf quaternion2 = readQuaternionfArgv(&argIndex, argc, argv);
			float angle = Quaternionf_getAngle(quaternion, quaternion2);
			printFloat(angle);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printQuaternionf(quaternion);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedQuaternionIndex = 0; namedQuaternionIndex < namedQuaternionCount; namedQuaternionIndex++) {
				if (!strcmp(namedQuaternions[namedQuaternionIndex].name, argv[argIndex + 1])) {
					namedQuaternions[namedQuaternionIndex].quaternion = quaternion;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedQuaternions = realloc(namedQuaternions, (namedQuaternionCount + 1) * sizeof(*namedQuaternions));
			namedQuaternions[namedQuaternionCount].name = argv[++argIndex];
			namedQuaternions[namedQuaternionCount++].quaternion = quaternion;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized quaternion operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printQuaternionf(quaternion);
	exit(EXIT_SUCCESS);
}

static void vector2Loop(int argc, const char ** argv, int * ioArgIndex, Vector2f vector) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "normalize")) {
			Vector2f_normalize(&vector);
			
		} else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
			Vector2f_invert(&vector);
			
		} else if (!strcmp(argv[argIndex], "add")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "subtract")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_subtract(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "scale_to")) {
			float value = readFloatArgv(&argIndex, argc, argv);
			Vector2f_scaleTo(&vector, value);
			
		} else if (!strcmp(argv[argIndex], "multiply_scalar")) {
			float multiplier = readFloatArgv(&argIndex, argc, argv);
			vector = Vector2f_multiplyScalar(vector, multiplier);
			
		} else if (!strcmp(argv[argIndex], "multiply_vector")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_multiplyComponents(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "divide_scalar")) {
			float divisor = readFloatArgv(&argIndex, argc, argv);
			vector = Vector2f_divideScalar(vector, divisor);
			
		} else if (!strcmp(argv[argIndex], "divide_vector")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "interpolate")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			vector = Vector2f_interpolate(vector, vector2, value);
			
		} else if (!strcmp(argv[argIndex], "round")) {
			vector = Vector2f_round(vector);
			
		} else if (!strcmp(argv[argIndex], "reflect")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_reflect(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "project")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_project(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "line_position")) {
			Vector2f origin = readVector2fArgv(&argIndex, argc, argv);
			Vector2f normal = readVector2fArgv(&argIndex, argc, argv);
			vector = Vector2f_positionOnLine(vector, origin, normal);
			
		} else if (!strcmp(argv[argIndex], "transpose")) {
			vector = Vector2f_transpose(vector);
			
		} else if (!strcmp(argv[argIndex], "dot")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			float dot = Vector2f_dot(vector, vector2);
			printFloat(dot);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "cross")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			float cross = Vector2f_cross(vector, vector2);
			printFloat(cross);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "intersect")) {
			Vector2f direction = readVector2fArgv(&argIndex, argc, argv);
			Vector2f normal = readVector2fArgv(&argIndex, argc, argv);
			Vector2f intersection;
			bool success = Vector2f_intersectLine(vector, direction, normal, &intersection, NULL);
			if (!success) {
				fprintf(stderr, "Couldn't find intersection point on parallel lines ({%f, %f} -> {%f, %f}, {%f, %f})\n", vector.x, vector.y, direction.x, direction.y, normal.x, normal.y);
				exit(EXIT_FAILURE);
			}
			vector = intersection;
			
		} else if (!strcmp(argv[argIndex], "intersect_segments")) {
			Vector2f point1 = readVector2fArgv(&argIndex, argc, argv);
			Vector2f point2 = readVector2fArgv(&argIndex, argc, argv);
			Vector2f point3 = readVector2fArgv(&argIndex, argc, argv);
			Vector2f intersection;
			bool success = Vector2f_intersectLineSegments(vector, point1, point2, point3, &intersection);
			if (!success) {
				fprintf(stderr, "Couldn't find intersection point on line segments {%f, %f} -> {%f, %f} and {%f, %f} -> {%f, %f}\n", vector.x, vector.y, point1.x, point1.y, point2.x, point2.y, point3.x, point3.y);
				exit(EXIT_FAILURE);
			}
			vector = intersection;
			
		} else if (!strcmp(argv[argIndex], "rotate")) {
			float radians = readFloatAngleArgv(&argIndex, argc, argv);
			vector = Vector2f_rotate(vector, radians);
			
		} else if (!strcmp(argv[argIndex], "magnitude")) {
			float magnitude = Vector2f_magnitude(vector);
			printFloat(magnitude);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "distance")) {
			Vector2f vector2 = readVector2fArgv(&argIndex, argc, argv);
			float distance = Vector2f_distance(vector, vector2);
			printFloat(distance);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "to_barycenter")) {
			Vector2f vertex0 = readVector2fArgv(&argIndex, argc, argv);
			Vector2f vertex1 = readVector2fArgv(&argIndex, argc, argv);
			Vector2f vertex2 = readVector2fArgv(&argIndex, argc, argv);
			Vector3f barycenter = Vector2f_toBarycenter(vector, vertex0, vertex1, vertex2);
			printVector3f(barycenter);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix3x3")) {
			Matrix3x3f matrix = readMatrix3x3fArgv(&argIndex, argc, argv);
			vector = Matrix3x3f_multiplyVector2f(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix4x4") || !strcmp(argv[argIndex], "multiply_matrix")) {
			Matrix4x4f matrix = readMatrix4x4fArgv(&argIndex, argc, argv);
			vector = Matrix4x4f_multiplyVector2f(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix4x4_rotation_only") || !strcmp(argv[argIndex], "multiply_matrix_rotation_only")) {
			Matrix4x4f matrix = readMatrix4x4fArgv(&argIndex, argc, argv);
			vector = Matrix4x4f_multiplyVector2f_rotationOnly(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "swizzle")) {
			restartLoopWithVectorSwizzle(argc, argv, &argIndex, (float *) &vector, 2);
			*ioArgIndex = argIndex;
			return;
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printVector2f(vector);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedVector2Index = 0; namedVector2Index < namedVector2Count; namedVector2Index++) {
				if (!strcmp(namedVector2s[namedVector2Index].name, argv[argIndex + 1])) {
					namedVector2s[namedVector2Index].vector = vector;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedVector2s = realloc(namedVector2s, (namedVector2Count + 1) * sizeof(*namedVector2s));
			namedVector2s[namedVector2Count].name = argv[++argIndex];
			namedVector2s[namedVector2Count++].vector = vector;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized vector2 operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printVector2f(vector);
	exit(EXIT_SUCCESS);
}

static void vector3Loop(int argc, const char ** argv, int * ioArgIndex, Vector3f vector) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "normalize")) {
			Vector3f_normalize(&vector);
			
		} else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
			Vector3f_invert(&vector);
			
		} else if (!strcmp(argv[argIndex], "add")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "subtract")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_subtract(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "scale_to")) {
			float value = readFloatArgv(&argIndex, argc, argv);
			Vector3f_scaleTo(&vector, value);
			
		} else if (!strcmp(argv[argIndex], "multiply_scalar")) {
			float multiplier = readFloatArgv(&argIndex, argc, argv);
			vector = Vector3f_multiplyScalar(vector, multiplier);
			
		} else if (!strcmp(argv[argIndex], "multiply_vector")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_multiplyComponents(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "divide_scalar")) {
			float divisor = readFloatArgv(&argIndex, argc, argv);
			vector = Vector3f_divideScalar(vector, divisor);
			
		} else if (!strcmp(argv[argIndex], "divide_vector")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "rotate")) {
			if (argc <= argIndex + 2) {
				fprintf(stderr, "Not enough arguments to rotate command (requires <x/y/z>, <cw/ccw/180>)\n");
				exit(EXIT_FAILURE);
			}
			argIndex++;
			unsigned int axisNumber = 0;
			if (!strcmp(argv[argIndex], "x")) {
				axisNumber = 0;
			} else if (!strcmp(argv[argIndex], "y")) {
				axisNumber = 1;
			} else if (!strcmp(argv[argIndex], "z")) {
				axisNumber = 2;
			} else {
				fprintf(stderr, "Unrecognized axis argument \"%s\" to rotate command; must be \"x\", \"y\", or \"z\")\n", argv[argIndex]);
				exit(EXIT_FAILURE);
			}
			argIndex++;
			if (!strcmp(argv[argIndex], "cw")) {
				switch (axisNumber) {
					case 0:
						vector = Vector3f_rotate90XCW(vector);
						break;
					case 1:
						vector = Vector3f_rotate90YCW(vector);
						break;
					case 2:
						vector = Vector3f_rotate90ZCW(vector);
						break;
				}
			} else if (!strcmp(argv[argIndex], "ccw")) {
				switch (axisNumber) {
					case 0:
						vector = Vector3f_rotate90XCCW(vector);
						break;
					case 1:
						vector = Vector3f_rotate90YCCW(vector);
						break;
					case 2:
						vector = Vector3f_rotate90ZCCW(vector);
						break;
				}
			} else if (!strcmp(argv[argIndex], "180")) {
				switch (axisNumber) {
					case 0:
						vector = Vector3f_rotate180X(vector);
						break;
					case 1:
						vector = Vector3f_rotate180Y(vector);
						break;
					case 2:
						vector = Vector3f_rotate180Z(vector);
						break;
				}
			} else {
				fprintf(stderr, "Unrecognized angle argument \"%s\" to rotate command; must be \"cw\", \"ccw\", or \"180\")\n", argv[argIndex]);
				exit(EXIT_FAILURE);
			}
			
		} else if (!strcmp(argv[argIndex], "interpolate")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			vector = Vector3f_interpolate(vector, vector2, value);
			
		} else if (!strcmp(argv[argIndex], "round")) {
			vector = Vector3f_round(vector);
			
		} else if (!strcmp(argv[argIndex], "reflect")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_reflect(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "project")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_project(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "line_position")) {
			Vector3f origin = readVector3fArgv(&argIndex, argc, argv);
			Vector3f normal = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_positionOnLine(vector, origin, normal);
			
		} else if (!strcmp(argv[argIndex], "plane_position")) {
			Vector3f origin = readVector3fArgv(&argIndex, argc, argv);
			Vector3f normal = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_positionOnPlane(vector, origin, normal);
			
		} else if (!strcmp(argv[argIndex], "closest_line_to_line")) {
			Vector3f normal0 = readVector3fArgv(&argIndex, argc, argv);
			Vector3f origin1 = readVector3fArgv(&argIndex, argc, argv);
			Vector3f normal1 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_closestLinePointToLine(vector, normal0, origin1, normal1);
			
		} else if (!strcmp(argv[argIndex], "perpendicular")) {
			vector = Vector3f_perpendicular(vector);
			
		} else if (!strcmp(argv[argIndex], "dot")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			float dot = Vector3f_dot(vector, vector2);
			printFloat(dot);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "cross")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			vector = Vector3f_cross(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "intersect")) {
			Vector3f direction = readVector3fArgv(&argIndex, argc, argv);
			Vector3f normal = readVector3fArgv(&argIndex, argc, argv);
			Vector3f intersection;
			bool success = Vector3f_intersectPlane(vector, direction, normal, &intersection, NULL);
			if (!success) {
				fprintf(stderr, "Couldn't find intersection point on line parallel to plane ({%f, %f, %f} -> {%f, %f, %f}, {%f, %f, %f})\n", vector.x, vector.y, vector.z, direction.x, direction.y, direction.z, normal.x, normal.y, normal.z);
				exit(EXIT_FAILURE);
			}
			vector = intersection;
			
		} else if (!strcmp(argv[argIndex], "magnitude")) {
			float magnitude = Vector3f_magnitude(vector);
			printFloat(magnitude);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "distance")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			float distance = Vector3f_distance(vector, vector2);
			printFloat(distance);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "rotation_to")) {
			Vector3f vector2 = readVector3fArgv(&argIndex, argc, argv);
			if (argc <= argIndex + 1) {
				fprintf(stderr, "Not enough arguments to rotation_to command (requires \"print\" or \"quaternion\" after vector argument)\n");
				exit(EXIT_FAILURE);
			}
			Vector3f axis;
			float angle;
			Vector3f_getAxisAngleForRotation(vector, vector2, &axis, &angle);
			argIndex++;
			if (!strcmp(argv[argIndex], "quaternion")) {
				*ioArgIndex = argIndex;
				quaternionLoop(argc, argv, ioArgIndex, Quaternionf_fromAxisAngle(axis, angle));
				return;
			}
			if (!strcmp(argv[argIndex], "print")) {
				printAxisAngle(axis, angle);
				exit(EXIT_SUCCESS);
			}
			fprintf(stderr, "Unrecognized argument \"%s\" to rotation_to command; must be \"print\" or \"quaternion\")\n", argv[argIndex]);
			exit(EXIT_FAILURE);
			
		} else if (!strcmp(argv[argIndex], "to_barycenter")) {
			Vector3f vertex0 = readVector3fArgv(&argIndex, argc, argv);
			Vector3f vertex1 = readVector3fArgv(&argIndex, argc, argv);
			Vector3f vertex2 = readVector3fArgv(&argIndex, argc, argv);
			Vector3f barycenter = Vector3f_toBarycenter(vector, vertex0, vertex1, vertex2);
			printVector3f(barycenter);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix3x3")) {
			Matrix3x3f matrix = readMatrix3x3fArgv(&argIndex, argc, argv);
			vector = Matrix3x3f_multiplyVector3f(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix4x4") || !strcmp(argv[argIndex], "multiply_matrix")) {
			Matrix4x4f matrix = readMatrix4x4fArgv(&argIndex, argc, argv);
			vector = Matrix4x4f_multiplyVector3f(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix4x4_rotation_only") || !strcmp(argv[argIndex], "multiply_matrix_rotation_only")) {
			Matrix4x4f matrix = readMatrix4x4fArgv(&argIndex, argc, argv);
			vector = Matrix4x4f_multiplyVector3f_rotationOnly(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_quaternion")) {
			Quaternionf quaternion = readQuaternionfArgv(&argIndex, argc, argv);
			vector = Quaternionf_multiplyVector3f(quaternion, vector);
			
		} else if (!strcmp(argv[argIndex], "swizzle")) {
			restartLoopWithVectorSwizzle(argc, argv, &argIndex, (float *) &vector, 3);
			*ioArgIndex = argIndex;
			return;
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printVector3f(vector);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedVector3Index = 0; namedVector3Index < namedVector3Count; namedVector3Index++) {
				if (!strcmp(namedVector3s[namedVector3Index].name, argv[argIndex + 1])) {
					namedVector3s[namedVector3Index].vector = vector;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedVector3s = realloc(namedVector3s, (namedVector3Count + 1) * sizeof(*namedVector3s));
			namedVector3s[namedVector3Count].name = argv[++argIndex];
			namedVector3s[namedVector3Count++].vector = vector;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized vector3 operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printVector3f(vector);
	exit(EXIT_SUCCESS);
}

static void vector4Loop(int argc, const char ** argv, int * ioArgIndex, Vector4f vector) {
	for (int argIndex = *ioArgIndex + 1; argIndex < argc; argIndex++) {
		if (!strcmp(argv[argIndex], "normalize")) {
			Vector4f_normalize(&vector);
			
		} else if (!strcmp(argv[argIndex], "normalize_w")) {
			Vector4f_normalizeW(&vector);
			
		} else if (!strcmp(argv[argIndex], "inverse") || !strcmp(argv[argIndex], "invert")) {
			Vector4f_invert(&vector);
			
		} else if (!strcmp(argv[argIndex], "add")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "subtract")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_subtract(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "scale_to")) {
			float value = readFloatArgv(&argIndex, argc, argv);
			Vector4f_scaleTo(&vector, value);
			
		} else if (!strcmp(argv[argIndex], "multiply_scalar")) {
			float multiplier = readFloatArgv(&argIndex, argc, argv);
			vector = Vector4f_multiplyScalar(vector, multiplier);
			
		} else if (!strcmp(argv[argIndex], "multiply_vector")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_multiplyComponents(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "divide_scalar")) {
			float divisor = readFloatArgv(&argIndex, argc, argv);
			vector = Vector4f_divideScalar(vector, divisor);
			
		} else if (!strcmp(argv[argIndex], "divide_vector")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_add(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "interpolate")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			float value = readFloatArgv(&argIndex, argc, argv);
			vector = Vector4f_interpolate(vector, vector2, value);
			
		} else if (!strcmp(argv[argIndex], "round")) {
			vector = Vector4f_round(vector);
			
		} else if (!strcmp(argv[argIndex], "reflect")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_reflect(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "project")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_project(vector, vector2);
			
		} else if (!strcmp(argv[argIndex], "line_position")) {
			Vector4f origin = readVector4fArgv(&argIndex, argc, argv);
			Vector4f normal = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_positionOnLine(vector, origin, normal);
			
		} else if (!strcmp(argv[argIndex], "plane_position")) {
			Vector4f origin = readVector4fArgv(&argIndex, argc, argv);
			Vector4f normal = readVector4fArgv(&argIndex, argc, argv);
			vector = Vector4f_positionOnPlane(vector, origin, normal);
			
		} else if (!strcmp(argv[argIndex], "dot")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			float dot = Vector4f_dot(vector, vector2);
			printFloat(dot);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "magnitude")) {
			float magnitude = Vector4f_magnitude(vector);
			printFloat(magnitude);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "distance")) {
			Vector4f vector2 = readVector4fArgv(&argIndex, argc, argv);
			float distance = Vector4f_distance(vector, vector2);
			printFloat(distance);
			exit(EXIT_SUCCESS);
			
		} else if (!strcmp(argv[argIndex], "multiply_matrix4x4") || !strcmp(argv[argIndex], "multiply_matrix")) {
			Matrix4x4f matrix = readMatrix4x4fArgv(&argIndex, argc, argv);
			vector = Matrix4x4f_multiplyVector4f(matrix, vector);
			
		} else if (!strcmp(argv[argIndex], "multiply_quaternion")) {
			Quaternionf quaternion = readQuaternionfArgv(&argIndex, argc, argv);
			vector = Quaternionf_multiplyVector4f(quaternion, vector);
			
		} else if (!strcmp(argv[argIndex], "swizzle")) {
			restartLoopWithVectorSwizzle(argc, argv, &argIndex, (float *) &vector, 4);
			*ioArgIndex = argIndex;
			return;
			
		} else if (!strcmp(argv[argIndex], "print")) {
			printVector4f(vector);
			
		} else if (!strcmp(argv[argIndex], "assign")) {
			if (argc <= argIndex + 1) {
				fprintf(stderr, "No name specified after assign command\n");
				printUsage();
				exit(EXIT_FAILURE);
			}
			for (unsigned int namedVector4Index = 0; namedVector4Index < namedVector4Count; namedVector4Index++) {
				if (!strcmp(namedVector4s[namedVector4Index].name, argv[argIndex + 1])) {
					namedVector4s[namedVector4Index].vector = vector;
					*ioArgIndex = argIndex + 1;
					return;
				}
			}
			namedVector4s = realloc(namedVector4s, (namedVector4Count + 1) * sizeof(*namedVector4s));
			namedVector4s[namedVector4Count].name = argv[++argIndex];
			namedVector4s[namedVector4Count++].vector = vector;
			*ioArgIndex = argIndex;
			return;
			
		} else {
			fprintf(stderr, "Unrecognized vector4 operation: %s\n", argv[argIndex]);
			break;
		}
	}
	printVector4f(vector);
	exit(EXIT_SUCCESS);
}

void startMatrix3x3fLoop(int argc, const char ** argv, int * ioArgIndex) {
	Matrix3x3f matrix = readMatrix3x3fArgv(ioArgIndex, argc, argv);
	matrix3x3Loop(argc, argv, ioArgIndex, matrix);
}

void startMatrix4x4fLoop(int argc, const char ** argv, int * ioArgIndex) {
	Matrix4x4f matrix = readMatrix4x4fArgv(ioArgIndex, argc, argv);
	matrix4x4Loop(argc, argv, ioArgIndex, matrix);
}

void startQuaternionfLoop(int argc, const char ** argv, int * ioArgIndex) {
	Quaternionf quaternion = readQuaternionfArgv(ioArgIndex, argc, argv);
	quaternionLoop(argc, argv, ioArgIndex, quaternion);
}

void startVector2fLoop(int argc, const char ** argv, int * ioArgIndex) {
	Vector2f vector = readVector2fArgv(ioArgIndex, argc, argv);
	vector2Loop(argc, argv, ioArgIndex, vector);
}

void startVector3fLoop(int argc, const char ** argv, int * ioArgIndex) {
	Vector3f vector = readVector3fArgv(ioArgIndex, argc, argv);
	vector3Loop(argc, argv, ioArgIndex, vector);
}

void startVector4fLoop(int argc, const char ** argv, int * ioArgIndex) {
	Vector4f vector = readVector4fArgv(ioArgIndex, argc, argv);
	vector4Loop(argc, argv, ioArgIndex, vector);
}