The reason I wrote a class on quaternions was simply for the reason of having a deeper understanding of how quaternions worked at the time. I was really interested on the fact that there is no gimbal lock and how matrix multiplications can insanely optimize your code from the generic Rodriguez Rotation Formula.

Quaternion.h

/***************************************
Author: Amir Azmi
Date: 8/30/2019
File Name: Quaternion.h
Description: Header for Quaternion Class
***************************************/
#pragma once
#include <cmath>
#include <algorithm>
#include <SFML/Graphics/Glsl.hpp>

/*
  Quaternions Rules

  associative but not communitive
  (a *b)^-1 = b^-1 *a^-1

  q * p * q^-1
*/

namespace math
{
  float ConvertToRadians(float degrees); //helper function

  struct Vec3 //didnt want use glm so made my own vec3 struct
  {
    public:
    float x, y, z;
    Vec3();
    Vec3(float x, float y, float z);
  };

  //think of quaternions as ordered pairs
  class Quaternion
  {
  public:
    //constructors
    Quaternion(); //identity
    Quaternion(float w, float x, float y, float z); //set quaternion
    Quaternion(const Vec3 &axis, float angle); // axis angle qauternion

    //operators
    const Quaternion operator+(const Quaternion& b) const;
    const Quaternion operator-(const Quaternion& b) const;
    const Quaternion operator*(const Quaternion& b); //quaternion cross product
    const Quaternion operator*(float scalar);
    Quaternion operator/(float scalar);

    //convert to matrix
    const sf::Transform QuaternionToMatrix();

    //operations
    Quaternion& normalized();
    Quaternion& conjugateQuaternion();
    Quaternion& inverse();
    Quaternion& projection(Quaternion& b);
    float magnitudeQuaternionSquared();
    float dotProduct(Quaternion& b);
    void negateQuaternion();
    void setToRotatAboutX(float theta);
    void setToRotatAboutY(float theta);
    void setToRotatAboutZ(float theta);
    void setToRotatAboutAxis(Vec3& axis, float theta);

    float x, y, z, w;

  private:
    //s is the w component
    //q = [s,v] where s is an element of all real numbers and v is an element of r3
    //quat = [s, xi + yj + zk] where i,j,k are elements of real numbers
  };
  
  template<class T>
  const T& clamp(const T& x, const T& upper, const T& lower);

  template<class T>
  const T& clamp(const T& x, const T& upper, const T& lower)
  {
    return std::min(upper, std::max(x, lower));
  }

  //more quat functions but differnt ways of doing the same thing
  float dotProduct(Quaternion a, Quaternion b);
  Quaternion Slerp(Quaternion& a, Quaternion& b, float time);
  const Quaternion operator*(float scalar, const Quaternion& a);
}

Quaternion.cpp

/***************************************
Author: Amir Azmi
Date: 8/30/2019
File Name: Quaternion.cpp
Description: Quaternion implementation
             becasue I wanted to learn
             what quaternions were.
***************************************/

#include "Quaternion.h" //header file
#define _USE_MATH_DEFINES //allowed to use pi
#include <math.h> //pi define, sinf, sqrtf, cosf
using namespace math;

//helper function to convert to radians
float math::ConvertToRadians(float degrees)
{
  return degrees * (M_PI / 180.0f);
}

//quaternion identity constructor
Quaternion::Quaternion() :w(1.0f), x(0.0f), y(0.0f), z(0.0f)
{
}

//set quaternion constructor
Quaternion::Quaternion(float w, float x, float y, float z) : w(w), x(x), y(y), z(z)
{
}

//axis angle qauternion constructor
math::Quaternion::Quaternion(const Vec3& axis, float angle)
{
  angle = ConvertToRadians(angle);
  w = cosf(0.5f * angle);

  float thetaOver2 = sinf(0.5f * angle);
  x = axis.x * thetaOver2;
  y = axis.y * thetaOver2;
  z = axis.z * thetaOver2;
}

//adding two quaternions
const Quaternion Quaternion::operator+(const Quaternion& b) const
{
  Quaternion c;
  c.w = this->w + b.w;
  c.x = this->x + b.x;
  c.y = this->y + b.y;
  c.z = this->z + b.z;

  return c;
}

//subtractiong two qauternions
const Quaternion Quaternion::operator-(const Quaternion& b) const
{
  Quaternion c;
  c.w = this->w - b.w;
  c.x = this->x - b.x;
  c.y = this->y - b.y;
  c.z = this->z - b.z;

  return c;
}

//cross product of quaternions
const Quaternion Quaternion::operator*(const Quaternion& b)
{
  Quaternion c;

  c.w = this->w * b.w - (this->x * b.x) - (this->y * b.y) - (this->z * b.z);
  c.x = this->w * b.x + (this->x * b.w) + (this->z * b.y) - (this->y * b.z);
  c.y = this->w * b.y + (this->y * b.w) + (this->x * b.z) - (this->z * b.x);
  c.z = this->w * b.z + (this->z * b.w) + (this->y * b.x) - (this->x * b.y);

  return c;
}

//scaling a quaternion
const Quaternion Quaternion::operator*(float scalar)
{
  Quaternion b;

  b.w *= scalar;
  b.x *= scalar;
  b.y *= scalar;
  b.z *= scalar;
  return b;
}

//dividing the componenets of a quaternion
Quaternion Quaternion::operator/(float scalar)
{
  return Quaternion(this->w / scalar, this->x / scalar, this->y / scalar, this->z / scalar);
}

//converting a quaternion into a matrix
//column major order for sfml
//returning a transform here is similar to returning a mat4 in SFML
const sf::Transform math::Quaternion::QuaternionToMatrix()
{
  return sf::Transform(
    1 - (2 * std::pow(this->y, 2)) - (2 * std::pow(this->z, 2)), (2 * this->x * this->y) - (2 * this->w * this->z), (2 * this->x * this->z) - (2 * this->w * this->y),
    (2 * this->x * this->y) + (2 * this->w * this->z), 1 - (2 * std::pow(this->x, 2)) - (2 * std::pow(this->z, 2)), (2 * this->y * this->z) + (2 * this->w * this->x),
    (2 * this->x * this->z) + (2 * this->w * this->y), (2 * this->y * this->z) - (2 * this->w * this->x), 1 - (2 * std::pow(this->x, 2)) - (2 * std::pow(this->y, 2)));
}

//normalizing a qauternion
Quaternion& math::Quaternion::normalized()
{
  float len_inv = 1.0f / sqrtf(this->w * this->w + this->x * this->x + this->y * this->y + this->z * this->z); //sqrtf is a heavy operation but cant be avoided

  this->w *= len_inv;
  this->x *= len_inv;
  this->y *= len_inv;
  this->z *= len_inv;

  return *this;
}

Quaternion& Quaternion::conjugateQuaternion() //conjugate quaternion
{
  this->x *= -1.0f;
  this->y *= -1.0f;
  this->z *= -1.0f;
  //notice how we don't negate w component

  return *this;
}

Quaternion& Quaternion::inverse()
{
  Quaternion b = *this;

  *this = b.conjugateQuaternion() / magnitudeQuaternionSquared();

  return *this;

}

Quaternion& math::Quaternion::projection(Quaternion& b)
{
  const float dot = dotProduct(b);
  this->w = dot * b.w;
  this->x = dot * b.x;
  this->y = dot * b.y;
  this->z = dot * b.z;

  return *this;
}

//squaring a quaternion
float Quaternion::magnitudeQuaternionSquared()
{
  return  std::pow(w, 2) + std::pow(x, 2) + std::pow(y, 2) + std::pow(z, 2);
}

//component wise multiplication to get a float aka dotproduct
float math::Quaternion::dotProduct(Quaternion& b)
{
  return this->w * b.w + this->x * b.x + this->y * b.y + this->z * b.z;
}

//quaternion negation
void Quaternion::negateQuaternion()
{
  this->w *= -1.0f;
  this->x *= -1.0f;
  this->y *= -1.0f;
  this->z *= -1.0f;
}

//rotate around X axis
void Quaternion::setToRotatAboutX(float theta)
{
  float thetaOver2 = sinf(theta / 2.0f);

  this->w = cosf(theta * 0.5f);
  this->x *= thetaOver2;
  this->y = 0.0f;
  this->z = 0.0f;

}

//rotate around Y axis
void Quaternion::setToRotatAboutY(float theta)
{
  float thetaOver2 = sinf(theta / 2.0f);

  this->w = cosf(theta * 0.5f);
  this->x = 0.0f;
  this->y *= thetaOver2;
  this->z = 0.0f;
}

//rotate around Z axis
void Quaternion::setToRotatAboutZ(float theta)
{
  float thetaOver2 = sinf(theta / 2.0f);

  this->w = cosf(theta * 0.5f);
  this->y = 0.0f;
  this->x = 0.0f;
  this->z *= thetaOver2;
}

//rotate around any axis
void Quaternion::setToRotatAboutAxis(Vec3& axis, float theta)
{
  float thetaOver2 = theta * 0.5f;
  float sinThetaOver2 = sinf(thetaOver2);

  this->w = cosf(thetaOver2);
  this->x = axis.x * sinThetaOver2;
  this->y = axis.y * sinThetaOver2;
  this->z = axis.z * sinThetaOver2;
}

// useful duplicate functions to have
//dot product of quaternion
float math::dotProduct(Quaternion a, Quaternion b)
{
  return a.w * b.w + a.x * b.x + a.y * b.y + a.z * b.z;
}

//slerp for quaternions ie go from one rotation to the other in a set amount of time (Spherical Linear Interpolation)
Quaternion math::Slerp(Quaternion& a, Quaternion& b, float time)
{
  const float omega = acosf(math::clamp(a.dotProduct(b), -1.0f, 1.0f));
  const float sin_inv = 1.0f / sinf(omega);

  return sinf((1.0f - time) * omega) * sin_inv * a + sin(time * omega) * sin_inv * b;
}

const Quaternion math::operator*(float scalar, const Quaternion& a)
{
  return scalar * a;
}

//constructor for vec3 identity
math::Vec3::Vec3() :x(1.0f), y(0.0f), z(0.0f)
{
}

//constructor for setting vec3
math::Vec3::Vec3(float x, float y, float z) : x(x), y(y), z(z)
{
}

Main.cpp

/***************************************
Author: Amir Azmi
Date: 8/30/2019
File Name: main.cpp
Description: testing my quaterion class
***************************************/
#include <SFML/Graphics.hpp> //rect, and window functions
#include "Quaternion.h" //header file

int main()
{
  math::Quaternion a(math::Vec3(0.0f,0.0f,1.0f), 45.f); //axis angle quaternion

  sf::RenderWindow window(sf::VideoMode(640, 480 ), "SFML Window" ); //create the window in SFML

  while (window.isOpen())
  {
    sf::RectangleShape rect; //create empty rectangle shape
    
    //compare SFML built in functions with my Quaternion rotation function
    rect.setSize(sf::Vector2f(100,50));
    //rect.setScale(2, 1);
    //rect.setRotation(10);
    //rect.setPosition(100,100);

    /**************************************************************************/
    //quaternion to matrix conversion
    
    sf::Transform t1 = a.QuaternionToMatrix(); //rotation matrix

    sf::Transform t5 = { //scale matrix
       2,0,0,
       0,2,0,
       0,0,1 };

    sf::Transform t6 = { //translate matrix
       1,0,100,
       0,1,100,
       0,0,1 };
    /**************************************************************************/

    sf::Event event;

    while (window.pollEvent(event))
    {
      // "close requested" event: we close the window
      if (event.type == sf::Event::Closed)
        window.close();
    }

    window.clear(); //clear the window -> probalby GlClearColor()

    window.draw(rect, t6 * t1 * t5); //draw rectangle with applied transformations
    window.display(); //display -> maybe swapping back buffer here, I dont know
  }
}