diff --git a/.github/workflows/test-gglm.yml b/.github/workflows/test-gglm.yml
new file mode 100755
index 0000000..6e0bd9a
--- /dev/null
+++ b/.github/workflows/test-gglm.yml
@@ -0,0 +1,21 @@
+name: test-gglm
+
+on:
+  create:
+  workflow_dispatch:
+
+jobs:
+  test-gglm:
+    runs-on: [windows-latest, macos-latest, ubuntu-latest]
+    steps:
+      - name: Install golang
+        uses: actions/setup-go@v3
+        with:
+          go-version: ">=1.17"
+
+      - name: Clone gglm
+        run: git clone https://github.com/bloeys/gglm
+
+      - name: Test gglm
+        working-directory: nmage
+        run: cd gglm && go test ./... -v
diff --git a/gglm/quat.go b/gglm/quat.go
index 0eec549..55cea63 100755
--- a/gglm/quat.go
+++ b/gglm/quat.go
@@ -47,36 +47,15 @@ func (q *Quat) Axis() Vec3 {
 	}}
 }
 
-// Euler takes rotations in radians and produces a rotation that
+// NewQuatEulerVec takes rotations in radians and produces a rotation that
 // rotates around the z-axis, y-axis and lastly x-axis.
-func NewQuatEuler(v *Vec3) Quat {
-
-	//Some other common terminology: x=roll, y=pitch, z=yaw
-	sinX, cosX := Sincos32(v.Data[0] * 0.5)
-	sinY, cosY := Sincos32(v.Data[1] * 0.5)
-	sinZ, cosZ := Sincos32(v.Data[2] * 0.5)
-
-	//This produces a z->y->x multiply order, but its written as XYZ.
-	//This is due to XYZ meaning independent rotation matrices, so Z is applied
-	//first, then Y matrix and lastly X.
-	//See this for more info: https://github.com/godotengine/godot/issues/6816#issuecomment-254592170
-	//
-	//Note: On most conversion tools putting the multiply order (e.g. ZYX for us) is required.
-	return Quat{
-		Vec4: Vec4{
-			Data: [4]float32{
-				sinX*cosY*cosZ - cosX*sinY*sinZ,
-				cosX*sinY*cosZ + sinX*cosY*sinZ,
-				cosX*cosY*sinZ - sinX*sinY*cosZ,
-				cosX*cosY*cosZ + sinX*sinY*sinZ,
-			},
-		},
-	}
+func NewQuatEulerVec(v *Vec3) Quat {
+	return NewQuatEuler(v.X(), v.Y(), v.Z())
 }
 
-// Euler takes rotations in radians and produces a rotation that
+// NewQuatEuler takes rotations in radians and produces a rotation that
 // rotates around the z-axis, y-axis and lastly x-axis.
-func NewQuatEulerXYZ(x, y, z float32) Quat {
+func NewQuatEuler(x, y, z float32) Quat {
 
 	//Some other common terminology: x=roll, y=pitch, z=yaw
 	sinX, cosX := Sincos32(x * 0.5)
@@ -101,16 +80,21 @@ func NewQuatEulerXYZ(x, y, z float32) Quat {
 	}
 }
 
+// NewQuatAngleAxisVec produces a quaternion thats rotates rotRad radians around the *normalized* vector rotAxisNorm
+func NewQuatAngleAxisVec(rotRad float32, rotAxisNorm *Vec3) Quat {
+	return NewQuatAngleAxis(rotRad, rotAxisNorm.X(), rotAxisNorm.Y(), rotAxisNorm.Z())
+}
+
 // NewQuatAngleAxis produces a quaternion thats rotates rotRad radians around the *normalized* vector rotAxisNorm
-func NewQuatAngleAxis(rotRad float32, rotAxisNorm *Vec3) Quat {
+func NewQuatAngleAxis(rotRad float32, rotAxisNormX, rotAxisNormY, rotAxisNormZ float32) Quat {
 
 	s, c := Sincos32(rotRad * 0.5)
 	return Quat{
 		Vec4: Vec4{
 			Data: [4]float32{
-				rotAxisNorm.Data[0] * s,
-				rotAxisNorm.Data[1] * s,
-				rotAxisNorm.Data[2] * s,
+				rotAxisNormX * s,
+				rotAxisNormY * s,
+				rotAxisNormZ * s,
 				c,
 			},
 		},
@@ -124,3 +108,27 @@ func NewQuatId() Quat {
 		},
 	}
 }
+
+func NewQuat(x, y, z, w float32) Quat {
+	return Quat{
+		Vec4: Vec4{
+			Data: [4]float32{x, y, z, w},
+		},
+	}
+}
+
+func NewQuatArr(arr [4]float32) Quat {
+	return Quat{
+		Vec4: Vec4{
+			Data: arr,
+		},
+	}
+}
+
+func NewQuatVec(v *Vec4) Quat {
+	return Quat{
+		Vec4: Vec4{
+			Data: v.Data,
+		},
+	}
+}
diff --git a/gglm/quat_test.go b/gglm/quat_test.go
index 109da41..37a3ce7 100755
--- a/gglm/quat_test.go
+++ b/gglm/quat_test.go
@@ -10,14 +10,14 @@ func TestNewQuatEuler(t *testing.T) {
 
 	degs := gglm.NewVec3(180, 180, 180)
 	degs.Data = degs.AsRad().Data
-	q := gglm.NewQuatEuler(&degs)
+	q := gglm.NewQuatEulerVec(&degs)
 	ans := &gglm.Quat{Vec4: gglm.NewVec4(0, 0, 0, 1)}
 
 	if !gglm.EqF32(q.X(), ans.X()) || !gglm.EqF32(q.Y(), ans.Y()) || !gglm.EqF32(q.Z(), ans.Z()) || !gglm.EqF32(q.W(), ans.W()) {
 		t.Errorf("Got: %v; Expected: %v", q.String(), ans.String())
 	}
 
-	q = gglm.NewQuatEulerXYZ(180*gglm.Deg2Rad, 180*gglm.Deg2Rad, 180*gglm.Deg2Rad)
+	q = gglm.NewQuatEuler(180*gglm.Deg2Rad, 180*gglm.Deg2Rad, 180*gglm.Deg2Rad)
 
 	if !gglm.EqF32(q.X(), ans.X()) || !gglm.EqF32(q.Y(), ans.Y()) || !gglm.EqF32(q.Z(), ans.Z()) || !gglm.EqF32(q.W(), ans.W()) {
 		t.Errorf("Got: %v; Expected: %v", q.String(), ans.String())
@@ -27,7 +27,7 @@ func TestNewQuatEuler(t *testing.T) {
 func TestNewQuatAngleAxis(t *testing.T) {
 
 	rotAxis := gglm.NewVec3(0, 1, 0)
-	q := gglm.NewQuatAngleAxis(180*gglm.Deg2Rad, &rotAxis)
+	q := gglm.NewQuatAngleAxisVec(180*gglm.Deg2Rad, &rotAxis)
 	ans := &gglm.Quat{Vec4: gglm.NewVec4(0, 1, 0, 0)}
 
 	if !gglm.EqF32(q.X(), ans.X()) || !gglm.EqF32(q.Y(), ans.Y()) || !gglm.EqF32(q.Z(), ans.Z()) || !gglm.EqF32(q.W(), ans.W()) {
@@ -38,7 +38,7 @@ func TestNewQuatAngleAxis(t *testing.T) {
 func TestQuatAngle(t *testing.T) {
 
 	rotAxis := gglm.NewVec3(0, 1, 0)
-	quat := gglm.NewQuatAngleAxis(180*gglm.Deg2Rad, &rotAxis)
+	quat := gglm.NewQuatAngleAxisVec(180*gglm.Deg2Rad, &rotAxis)
 	a := quat.Angle()
 	var ans float32 = 180.0 * gglm.Deg2Rad
 
@@ -47,7 +47,7 @@ func TestQuatAngle(t *testing.T) {
 	}
 
 	rotAxis = gglm.NewVec3(1, 1, 0)
-	quat = gglm.NewQuatAngleAxis(90*gglm.Deg2Rad, rotAxis.Normalize())
+	quat = gglm.NewQuatAngleAxisVec(90*gglm.Deg2Rad, rotAxis.Normalize())
 	a = quat.Angle()
 	ans = 90 * gglm.Deg2Rad
 
@@ -56,7 +56,7 @@ func TestQuatAngle(t *testing.T) {
 	}
 
 	rotAxis = gglm.NewVec3(1, 1, 0)
-	quat = gglm.NewQuatAngleAxis(125*gglm.Deg2Rad, rotAxis.Normalize())
+	quat = gglm.NewQuatAngleAxisVec(125*gglm.Deg2Rad, rotAxis.Normalize())
 	a = quat.Angle()
 	ans = 125 * gglm.Deg2Rad
 
@@ -68,7 +68,7 @@ func TestQuatAngle(t *testing.T) {
 func TestQuatAxis(t *testing.T) {
 
 	rotAxis := gglm.NewVec3(0, 1, 0)
-	quat := gglm.NewQuatAngleAxis(1, &rotAxis)
+	quat := gglm.NewQuatAngleAxisVec(1, &rotAxis)
 	a := quat.Axis()
 	ans := gglm.NewVec3(0, 1, 0)
 
@@ -77,7 +77,7 @@ func TestQuatAxis(t *testing.T) {
 	}
 
 	rotAxis = gglm.NewVec3(1, 1, 0)
-	quat = gglm.NewQuatAngleAxis(1, rotAxis.Normalize())
+	quat = gglm.NewQuatAngleAxisVec(1, rotAxis.Normalize())
 	a = quat.Axis()
 	ans = gglm.NewVec3(1, 1, 0)
 	ans.Normalize()
@@ -87,7 +87,7 @@ func TestQuatAxis(t *testing.T) {
 	}
 
 	rotAxis = gglm.NewVec3(67, 46, 32)
-	quat = gglm.NewQuatAngleAxis(1, rotAxis.Normalize())
+	quat = gglm.NewQuatAngleAxisVec(1, rotAxis.Normalize())
 	a = quat.Axis()
 	ans = gglm.NewVec3(67, 46, 32)
 	ans.Normalize()
diff --git a/gglm/vec_test.go b/gglm/vec_test.go
index cf8fb1f..2339db0 100755
--- a/gglm/vec_test.go
+++ b/gglm/vec_test.go
@@ -209,7 +209,7 @@ func TestVecSwizzleSet(t *testing.T) {
 	// Test rot by quat
 	v32.SetXY(1, 0)
 	rotAxis := gglm.NewVec3(0, 1, 0)
-	quat := gglm.NewQuatAngleAxis(90*gglm.Deg2Rad, &rotAxis)
+	quat := gglm.NewQuatAngleAxisVec(90*gglm.Deg2Rad, &rotAxis)
 	v32.RotByQuat(&quat)
 	angleV3 = gglm.AngleVec3(&v3, &v32) * gglm.Rad2Deg
 	if angleV3 != 90 {
diff --git a/main.go b/main.go
index 32f0aee..c7935c5 100755
--- a/main.go
+++ b/main.go
@@ -143,11 +143,11 @@ func main() {
 
 	// Quaternion
 	vRot := &gglm.Vec3{Data: [3]float32{60, 30, 20}}
-	q := gglm.NewQuatEuler(vRot.AsRad())
+	q := gglm.NewQuatEulerVec(vRot.AsRad())
 	println("\n" + vRot.AsRad().String())
 	println(q.String(), "\n", q.Mag())
 
-	q = gglm.NewQuatAngleAxis(60*gglm.Deg2Rad, vRot.Normalize())
+	q = gglm.NewQuatAngleAxisVec(60*gglm.Deg2Rad, vRot.Normalize())
 	println("\n" + vRot.Normalize().String())
 	println(q.String())
 
@@ -155,7 +155,7 @@ func main() {
 	translationMat := gglm.NewTranslationMat(&gglm.Vec3{Data: [3]float32{1, 2, 3}})
 
 	rotDegs := gglm.NewVec3(60, 30, 20)
-	quat := gglm.NewQuatEuler(rotDegs.AsRad())
+	quat := gglm.NewQuatEulerVec(rotDegs.AsRad())
 	rotMat := gglm.NewRotMat(&quat)
 
 	scale := gglm.NewVec3(1, 1, 1)
@@ -186,10 +186,10 @@ func main() {
 
 	// Quat geo
 	q1Degs := gglm.NewVec3(180*gglm.Deg2Rad, 0, 0)
-	q1 := gglm.NewQuatEuler(&q1Degs)
+	q1 := gglm.NewQuatEulerVec(&q1Degs)
 
 	q2Degs := gglm.NewVec3(0, 180*gglm.Deg2Rad, 0)
-	q2 := gglm.NewQuatEuler(&q2Degs)
+	q2 := gglm.NewQuatEulerVec(&q2Degs)
 	println(gglm.AngleQuat(&q1, &q2) * gglm.Rad2Deg)
 
 	// LookAt