power_play/src/base/base_math.c

////////////////////////////////////////////////////////////
//~ Min / max

//- Min
u8 MinU8(u8 a, u8 b) { return a <= b ? a : b; }
i8 MinI8(i8 a, i8 b) { return a <= b ? a : b; }
u32 MinU32(u32 a, u32 b) { return a <= b ? a : b; }
i32 MinI32(i32 a, i32 b) { return a <= b ? a : b; }
f32 MinF32(f32 a, f32 b) { return a <= b ? a : b; }
u64 MinU64(u64 a, u64 b) { return a <= b ? a : b; }
i64 MinI64(i64 a, i64 b) { return a <= b ? a : b; }
f64 MinF64(f64 a, f64 b) { return a <= b ? a : b; }

//- Max
u8 MaxU8(u8 a, u8 b) { return a >= b ? a : b; }
i8 MaxI8(i8 a, i8 b) { return a >= b ? a : b; }
u32 MaxU32(u32 a, u32 b) { return a >= b ? a : b; }
i32 MaxI32(i32 a, i32 b) { return a >= b ? a : b; }
f32 MaxF32(f32 a, f32 b) { return a >= b ? a : b; }
u64 MaxU64(u64 a, u64 b) { return a >= b ? a : b; }
i64 MaxI64(i64 a, i64 b) { return a >= b ? a : b; }
f64 MaxF64(f64 a, f64 b) { return a >= b ? a : b; }

//- Clamp
u32 ClampU32(u32 v, u32 min, u32 max) { return v < min ? min : v > max ? max : v; }
i32 ClampI32(i32 v, i32 min, i32 max) { return v < min ? min : v > max ? max : v; }
f32 ClampF32(f32 v, f32 min, f32 max) { return v < min ? min : v > max ? max : v; }
u64 ClampU64(u64 v, u64 min, u64 max) { return v < min ? min : v > max ? max : v; }
i64 ClampI64(i64 v, i64 min, i64 max) { return v < min ? min : v > max ? max : v; }
f64 ClampF64(f64 v, f64 min, f64 max) { return v < min ? min : v > max ? max : v; }

//- Saturate
u32 SaturateU32(u32 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }
i32 SaturateI32(i32 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }
f32 SaturateF32(f32 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }
u64 SaturateU64(u64 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }
i64 SaturateI64(i64 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }
f64 SaturateF64(f64 v) { return v < 0 ? 0 : v > 1 ? 1 : v; }

////////////////////////////////////////////////////////////
//~ Numeric ops

//- Sign

i32 SignF32(f32 v)
{
  return (v >= 0) - (v < 0);
}

i32 SignF64(f64 v)
{
  return (v >= 0) - (v < 0);
}

//- Normalize

f32 Norm8(u32 v)
{
  return (v & 0xFF) / (f32)0x100;
}

f32 Norm16(u32 v)
{
  return (v & 0xFFFF) / (f32)0x10000;
}

f32 Norm24(u32 v)
{
  return (v & 0xFFFFFF) / (f32)0x1000000;
}

f64 Norm53(u64 v)
{
  return (v & 0x1FFFFFFFFFFFFFull) / (f64)0x20000000000000ull;
}

////////////////////////////////////////////////////////////
//~ Exponential ops

//- Pow u64

// Taken from https://gist.github.com/orlp/3551590
u64 PowU64(u64 v, u8 exp)
{
  PERSIST const u8 highest_bit_set[] = {
    0, 1, 2, 2, 3, 3, 3, 3,
    4, 4, 4, 4, 4, 4, 4, 4,
    5, 5, 5, 5, 5, 5, 5, 5,
    5, 5, 5, 5, 5, 5, 5, 5,
    6, 6, 6, 6, 6, 6, 6, 6,
    6, 6, 6, 6, 6, 6, 6, 6,
    6, 6, 6, 6, 6, 6, 6, 6,
    6, 6, 6, 6, 6, 6, 6, 255, // Anything past 63 is a guaranteed overflow with v > 1
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
    255, 255, 255, 255, 255, 255, 255, 255,
  };

  u64 result = 1;

  switch (highest_bit_set[exp])
  {
    case 255:
    {
      // 255 = overflow, return 0
      if (v == 1)
      {
        return 1;
      }
      // if (v == -1)
      // {
      //     return 1 - 2 * (exp & 1);
      // }
      return 0;
    } break;
    case 6:
    {
      if (exp & 1) result *= v;
      exp >>= 1;
      v *= v;
    } FALLTHROUGH;
    case 5:
    {
      if (exp & 1) result *= v;
      exp >>= 1;
      v *= v;
    } FALLTHROUGH;
    case 4:
    {
      if (exp & 1) result *= v;
      exp >>= 1;
      v *= v;
    } FALLTHROUGH;
    case 3:
    {
      if (exp & 1) result *= v;
      exp >>= 1;
      v *= v;
    } FALLTHROUGH;
    case 2:
    {
      if (exp & 1) result *= v;
      exp >>= 1;
      v *= v;
    } FALLTHROUGH;
    case 1:
    {
      if (exp & 1) result *= v;
    } FALLTHROUGH;
    default: return result;
  }
}

////////////////////////////////////////////////////////////
//~ Align

u64 AlignU64(u64 x, u64 align)
{
  align = MaxU64(align, 1);
  u64 result = (x + (align - 1));
  result -= result % align;
  return result;
}

u64 NextPow2U64(u64 x)
{
  u64 result = 0;
  if (x > 0)
  {
    result = x - 1;
    result |= result >> 1;
    result |= result >> 2;
    result |= result >> 4;
    result |= result >> 8;
    result |= result >> 16;
    result |= result >> 32;
    ++result;
  }
  return result;
}

////////////////////////////////////////////////////////////
//~ Angle unwind

// Returns angle in range [-Pi, Pi]
f32 UnwindAngleF32(f32 a)
{
  f32 d = ModF32(a, Tau);
  return ModF32(2.0f * d, Tau) - d;
}

////////////////////////////////////////////////////////////
//~ Float lerp

f32 LerpF32(f32 a, f32 b, f32 t)
{
  return a + ((b - a) * t);
}

f64 LerpF64(f64 a, f64 b, f64 t)
{
  return a + ((b - a) * t);
}

f32 LerpAngleF32(f32 a, f32 b, f32 t)
{
  f32 diff = UnwindAngleF32(b - a);
  return a + diff * t;
}

////////////////////////////////////////////////////////////
//~ Int lerp

i32 LerpI32(i32 a, i32 b, f32 t)
{
  return a + RoundF32((f32)(b - a) * t);
}

i64 LerpI64(i64 a, i64 b, f64 t)
{
  return a + RoundF64((f64)(b - a) * t);
}

i32 LerpU32(u32 a, u32 b, f32 t)
{
  return a + RoundF32((f32)(b - a) * t);
}

i64 LerpU64(u64 a, u64 b, f64 t)
{
  return a + RoundF64((f64)(b - a) * t);
}

////////////////////////////////////////////////////////////
//~ Smoothstep

f32 SmoothstepF32(f32 a, f32 b, f32 t)
{
  f32 result = 0;
  if (a == b)
  {
    result = (t < a) ? 0 : 1;
  }
  else
  {
    t = ClampF32((t - a) / (b - a), 0, 1);
    result = t * t * (3.0f - 2.0f * t);
  }
  return result;
}

f64 SmoothstepF64(f64 a, f64 b, f64 t)
{
  f64 result = 0;
  if (a == b)
  {
    result = (t < a) ? 0 : 1;
  }
  else
  {
    t = ClampF64((t - a) / (b - a), 0, 1);
    result = t * t * (3.0f - 2.0f * t);
  }
  return result;
}

////////////////////////////////////////////////////////////
//~ Color

u8 MulNormalizedU8(u8 a, u8 b)
{
  u32 t = ((u32)a * (u32)b) + 0x80;
  return ((t >> 8) + t) >> 8;
}

f32 SrgbFromLinearF32(f32 lin)
{
  f32 result = 0;
  if (lin <= 0.0031308f)
  {
    result = lin * 12.92f;
  }
  else
  {
    result = 1.055f * PowF32(lin, 1.0f/2.4f) - 0.055f;
  }
  return result;
}

f32 LinearFromSrgbF32(f32 srgb)
{
  f32 result = 0;
  if (srgb <= 0.04045f)
  {
    result = srgb / 12.92f;
  }
  else
  {
    result = PowF32((srgb + 0.055f) / 1.055f, 2.4f);
  }
  return result;
}

Vec4 LinearFromSrgb(Vec4 srgb)
{
  Vec4 result = Zi;
  result.x = LinearFromSrgbF32(srgb.x);
  result.y = LinearFromSrgbF32(srgb.y);
  result.z = LinearFromSrgbF32(srgb.z);
  result.w = srgb.w;
  return result;
}

Vec4 SrgbFromLinear(Vec4 lin)
{
  Vec4 result = Zi;
  result.x = SrgbFromLinearF32(lin.x);
  result.y = SrgbFromLinearF32(lin.y);
  result.z = SrgbFromLinearF32(lin.z);
  result.w = lin.w;
  return result;
}

Vec4 SrgbFromHsv(f32 h, f32 s, f32 v)
{
  h = ModF32(h, 360.0f);
  if (h < 0)
  {
    h += 360.0f;
  }
  f32 c = v * s;
  f32 x = c * (1.0f - AbsF32(ModF32(h / 60.0f, 2.0f) - 1.0f));
  f32 m = v - c;
  Vec4 result = VEC4(0, 0, 0, 1);
  if      (h < 60.0f)  { result.r = c; result.g = x; result.b = 0; }
  else if (h < 120.0f) { result.r = x; result.g = c; result.b = 0; }
  else if (h < 180.0f) { result.r = 0; result.g = c; result.b = x; }
  else if (h < 240.0f) { result.r = 0; result.g = x; result.b = c; }
  else if (h < 300.0f) { result.r = x; result.g = 0; result.b = c; }
  else                 { result.r = c; result.g = 0; result.b = x; }
  return result;
}

Vec4 PremulFromLinear(Vec4 lin)
{
  Vec4 result = Zi;
  result.x = lin.x * lin.w;
  result.y = lin.y * lin.w;
  result.z = lin.z * lin.w;
  result.w = lin.w;
  return result;
}

Vec4 PremulFromSrgb(Vec4 srgb)
{
  Vec4 lin = LinearFromSrgb(srgb);
  Vec4 premul = PremulFromLinear(lin);
  return premul;
}

Vec4 LerpSrgb(Vec4 a, Vec4 b, f32 t)
{
  Vec4 result = Zi;
  Vec4 a_lin = LinearFromSrgb(a);
  Vec4 b_lin = LinearFromSrgb(b);
  Vec4 lerp_lin = LerpVec4(a_lin, b_lin, t);
  result = SrgbFromLinear(lerp_lin);
  return result;
}

////////////////////////////////////////////////////////////
//~ Vec2

b32 IsVec2Zero(Vec2 a)
{
  return a.x == 0 && a.y == 0;
}

b32 MatchVec2(Vec2 a, Vec2 b)
{
  return a.x == b.x && a.y == b.y;
}

//- Mul

Vec2 MulVec2(Vec2 a, f32 s)
{
  return VEC2(a.x * s, a.y * s);
}

Vec2 MulVec2Vec2(Vec2 a, Vec2 b)
{
  return VEC2(a.x * b.x, a.y * b.y);
}

Vec2 NegVec2(Vec2 a)
{
  return VEC2(-a.x, -a.y);
}

//- Div

Vec2 DivVec2(Vec2 a, f32 s)
{
  f32 d = 1 / s;
  return VEC2(a.x * d, a.y * d);
}

Vec2 DivVec2Vec2(Vec2 a, Vec2 b)
{
  return VEC2(a.x * (1 / b.x), a.y * (1 / b.y));
}

Vec2 RecipVec2(Vec2 a)
{
  return VEC2(1.0 / a.x, 1.0 / a.y);
};

//- Add

Vec2 AddVec2(Vec2 a, Vec2 b)
{
  return VEC2(a.x + b.x, a.y + b.y);
}

Vec2 SubVec2(Vec2 a, Vec2 b)
{
  return VEC2(a.x - b.x, a.y - b.y);
}

//- Len

f32 Vec2Len(Vec2 a)
{
  return SqrtF32(a.x * a.x + a.y * a.y);
}

f32 Vec2LenSq(Vec2 a)
{
  return a.x * a.x + a.y * a.y;
}

Vec2 Vec2WithLen(Vec2 a, f32 len)
{
  f32 l_sq = a.x * a.x + a.y * a.y;
  if (l_sq > 0)
  {
    f32 denom = len / SqrtF32(l_sq);
    a.x *= denom;
    a.y *= denom;
  }
  else
  {
    a = VEC2(len, 0);
  }
  return a;
}

Vec2 ClampVec2Len(Vec2 a, f32 min, f32 max)
{
  f32 l_sq = a.x * a.x + a.y * a.y;
  if (l_sq > 0)
  {
    if (l_sq > max * max)
    {
      f32 denom = max / SqrtF32(l_sq);
      a.x *= denom;
      a.y *= denom;
    }
    else if (l_sq < min * min)
    {
      f32 denom = min / SqrtF32(l_sq);
      a.x *= denom;
      a.y *= denom;
    }
  }
  else
  {
    a = VEC2(min, 0);
  }
  return a;
}

f32 Vec2Distance(Vec2 a, Vec2 b)
{
  f32 dx = b.x - a.x;
  f32 dy = b.y - a.y;
  return SqrtF32(dx * dx + dy * dy);
}

Vec2 NormVec2(Vec2 a)
{
  return Vec2WithLen(a, 1.f);
}

//- Clamp

Vec2 ClampVec2(Vec2 v, Rng2 r)
{
  Vec2 result = Zi;
  result.x = ClampF32(v.x, r.p0.x, r.p1.x);
  result.y = ClampF32(v.y, r.p0.y, r.p1.y);
  return result;
}

Vec2 SaturateVec2(Vec2 v)
{
  Vec2 result = Zi;
  result.x = ClampF32(v.x, 0, 1);
  result.y = ClampF32(v.y, 0, 1);
  return result;
}

//- Dot

f32 DotVec2(Vec2 a, Vec2 b)
{
  return a.x * b.x + a.y * b.y;
}

// Returns signed area between vectors (positive in clockwise direction)
// AKA perpendicular dot product
// AKA 2d cross-product
f32 WedgeVec2(Vec2 a, Vec2 b)
{
  return a.x * b.y - a.y * b.x;
}

// Clockwise (right) perpendicular vector
Vec2 PerpVec2(Vec2 a)
{
  return VEC2(-a.y, a.x);
}

// Clockwise (right) perpendicular vector scaled by s
Vec2 MulPerpVec2(Vec2 a, f32 s)
{
  return VEC2(s * -a.y, s * a.x);
}

Vec2 PerpVec2TowardsDir(Vec2 v, Vec2 dir)
{
  f32 wedge = WedgeVec2(v, dir);
  return MulPerpVec2(v, (wedge >= 0) - (wedge < 0));
}

//- Round / floor / ceil

Vec2 RoundVec2(Vec2 a)
{
  return VEC2(RoundF32(a.x), RoundF32(a.y));
}

Vec2 FloorVec2(Vec2 a)
{
  return VEC2(FloorF32(a.x), FloorF32(a.y));
}

Vec2 CeilVec2(Vec2 a)
{
  return VEC2(CeilF32(a.x), CeilF32(a.y));
}

//- Rotation

// Returns 1 if winding between vectors a & b is clockwise or straight, -1 if counter-clockwise
i32 WindingFromVec2(Vec2 a, Vec2 b)
{
  f32 w = WedgeVec2(a, b);
  return (w >= 0) - (w < 0);
}

Vec2 RotateVec2Angle(Vec2 v, f32 a)
{
  f32 c = CosF32(a);
  f32 s = SinF32(a);
  return VEC2(v.x * c - v.y * s, v.x * s + v.y * c);
}

Vec2 RotateVec2(Vec2 a, Vec2 b)
{
  return VEC2(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);
}

Vec2 InvertRot(Vec2 r)
{
  r.y = -r.y;
  return r;
}

Vec2 Vec2FromAngle(f32 a)
{
  return VEC2(CosF32(a), SinF32(a));
}

f32 AngleFromVec2(Vec2 v)
{
  return ArcTan2F32(v.y, v.x);
}

//- Area

f32 AreaFromVec2(Vec2 v)
{
  return AbsF32(v.x * v.y);
}

//- Closest point

Vec2 ClosestPointFromRay(Vec2 ray_pos, Vec2 ray_dir_norm, Vec2 p)
{
  Vec2 ray_p_dir = SubVec2(p, ray_pos);
  f32 dot = DotVec2(ray_dir_norm, ray_p_dir);
  Vec2 ray_dir_closest = MulVec2(ray_dir_norm, dot);
  return AddVec2(ray_pos, ray_dir_closest);
}

//- Lerp

// Interpolate position vectors
Vec2 LerpVec2(Vec2 a, Vec2 b, f32 t)
{
  return VEC2(LerpF32(a.x, b.x, t), LerpF32(a.y, b.y, t));
}

Vec2 LerpVec2Vec2(Vec2 a, Vec2 b, Vec2 t)
{
  return VEC2(LerpF32(a.x, b.x, t.x), LerpF32(a.y, b.y, t.y));
}

// Interpolate direction vectors (spherical lerp)
Vec2 SlerpVec2(Vec2 a, Vec2 b, f32 t)
{
  f32 rot = LerpAngleF32(AngleFromVec2(a), AngleFromVec2(b), t);
  f32 len = LerpF32(Vec2Len(a), Vec2Len(b), t);
  return MulVec2(Vec2FromAngle(rot), len);
}

////////////////////////////////////////////////////////////
//~ Vec2I32

b32 MatchVec2I32(Vec2I32 a, Vec2I32 b)
{
  return a.x == b.x && a.y == b.y;
}

Vec2I32 NegVec2I32(Vec2I32 a)
{
  return VEC2I32(-a.x, -a.y);
}

Vec2I32 AddVec2I32(Vec2I32 a, Vec2I32 b)
{
  return VEC2I32(a.x + b.x, a.y + b.y);
}

Vec2I32 SubVec2I32(Vec2I32 a, Vec2I32 b)
{
  return VEC2I32(a.x - b.x, a.y - b.y);
}

Vec2I32 MulVec2I32Vec2I32(Vec2I32 a, Vec2I32 b)
{
  return VEC2I32(a.x * b.x, a.y * b.y);
}

Vec2I32 DivVec2I32Vec2I32(Vec2I32 a, Vec2I32 b)
{
  return VEC2I32(a.x / b.x, a.y / b.y);
}

////////////////////////////////////////////////////////////
//~ Vec4

//- Mul

Vec4 MulVec4(Vec4 v, f32 s)
{
  return VEC4(v.x * s, v.y * s, v.z * s, v.w * s);
}

Vec4 MulVec4Vec4(Vec4 a, Vec4 b)
{
  return VEC4(a.x * b.x, a.y * b.y, a.z * b.z, a.w * b.w);
}

//- Lerp

Vec4 LerpVec4(Vec4 a, Vec4 b, f32 t)
{
  Vec4 result = Zi;
  result.x = LerpF32(a.x, b.x, t);
  result.y = LerpF32(a.y, b.y, t);
  result.z = LerpF32(a.z, b.z, t);
  result.w = LerpF32(a.w, b.w, t);
  return result;
}

//- Conversion

Vec4 Vec4FromU32(u32 v)
{
  Vec4 result = Zi;
  result.x = ((v >> 0) & 0xFF) / 255.0;
  result.y = ((v >> 8) & 0xFF) / 255.0;
  result.z = ((v >> 16) & 0xFF) / 255.0;
  result.w = ((v >> 24) & 0xFF) / 255.0;
  return result;
}

u32 U32FromVec4(Vec4 v)
{
  u32 result = 0;
  result |= (((u32)(v.x * 255.0)) & 0xFF) << 0;
  result |= (((u32)(v.y * 255.0)) & 0xFF) << 8;
  result |= (((u32)(v.z * 255.0)) & 0xFF) << 16;
  result |= (((u32)(v.w * 255.0)) & 0xFF) << 24;
  return result;
}

////////////////////////////////////////////////////////////
//~ Range

//- Rng1

f32 NormRng(Rng r, f32 v)
{
  return (r.max - r.min) / v;
}

//- Rng2

b32 IsRng2Empty(Rng2 r)
{
  return r.p0.x >= r.p1.x || r.p0.y >= r.p1.y;
}

Vec2 DimsFromRng2(Rng2 r)
{
  Vec2 result = Zi;
  result.x = r.p1.x - r.p0.x;
  result.y = r.p1.y - r.p0.y;
  return result;
}

Vec2 CenterFromRng2(Rng2 r)
{
  Vec2 result = Zi;
  Vec2 dims = DimsFromRng2(r);
  result.x = r.p0.x + dims.x * 0.5;
  result.y = r.p0.y + dims.y * 0.5;
  return result;
}

Vec2 NormRng2(Rng2 r, Vec2 v)
{
  Vec2 result = Zi;
  Vec2 dims = SubVec2(r.p1, r.p0);
  result.x = v.x / dims.x;
  result.y = v.y / dims.y;
  return result;
}

Rng2 UnionRng2(Rng2 a, Rng2 b)
{
  Rng2 result = Zi;
  result.p0.x = MinF32(a.p0.x, b.p0.x);
  result.p0.y = MinF32(a.p0.y, b.p0.y);
  result.p1.x = MaxF32(a.p1.x, b.p1.x);
  result.p1.y = MaxF32(a.p1.y, b.p1.y);
  return result;
}

Rng2 IntersectRng2(Rng2 a, Rng2 b)
{
  Rng2 result = Zi;
  result.p0.x = MaxF32(a.p0.x, b.p0.x);
  result.p0.y = MaxF32(a.p0.y, b.p0.y);
  result.p1.x = MinF32(a.p1.x, b.p1.x);
  result.p1.y = MinF32(a.p1.y, b.p1.y);
  return result;
}

b32 IsIntersectingRng2(Rng2 a, Rng2 b)
{
  return !IsRng2Empty(IntersectRng2(a, b));
}

Rng2 AddRng2Vec2(Rng2 r, Vec2 v)
{
  Rng2 result = Zi;
  result.p0 = AddVec2(r.p0, v);
  result.p1 = AddVec2(r.p1, v);
  return result;
}

Rng2 MulRng2Vec2(Rng2 r, Vec2 v)
{
  Rng2 result = Zi;
  result.p0 = MulVec2Vec2(r.p0, v);
  result.p1 = MulVec2Vec2(r.p1, v);
  return result;
}

Rng2 DivRng2Vec2(Rng2 r, Vec2 v)
{
  Rng2 result = Zi;
  result.p0 = DivVec2Vec2(r.p0, v);
  result.p1 = DivVec2Vec2(r.p1, v);
  return result;
}

//- Rng2I32

b32 IsRng2I32Empty(Rng2I32 r)
{
  return r.p0.x >= r.p1.x || r.p0.y >= r.p1.y;
}

Vec2I32 DimsFromRng2I32(Rng2I32 r)
{
  Vec2I32 result = Zi;
  result.x = r.p1.x - r.p0.x;
  result.y = r.p1.y - r.p0.y;
  return result;
}

Vec2I32 CenterFromRng2I32(Rng2I32 r)
{
  Vec2I32 result = Zi;
  Vec2I32 dims = DimsFromRng2I32(r);
  result.x = r.p0.x + dims.x * 0.5;
  result.y = r.p0.y + dims.y * 0.5;
  return result;
}

Rng2I32 UnionRng2I32(Rng2I32 a, Rng2I32 b)
{
  Rng2I32 result = Zi;
  result.p0.x = MinI32(a.p0.x, b.p0.x);
  result.p0.y = MinI32(a.p0.y, b.p0.y);
  result.p1.x = MaxI32(a.p1.x, b.p1.x);
  result.p1.y = MaxI32(a.p1.y, b.p1.y);
  return result;
}

Rng2I32 IntersectRng2I32(Rng2I32 a, Rng2I32 b)
{
  Rng2I32 result = Zi;
  result.p0.x = MaxI32(a.p0.x, b.p0.x);
  result.p0.y = MaxI32(a.p0.y, b.p0.y);
  result.p1.x = MinI32(a.p1.x, b.p1.x);
  result.p1.y = MinI32(a.p1.y, b.p1.y);
  return result;
}

b32 IsIntersectingRng2I32(Rng2I32 a, Rng2I32 b)
{
  return !IsRng2I32Empty(IntersectRng2I32(a, b));
}

Rng2I32 AddRng2I32Vec2I32(Rng2I32 r, Vec2I32 v)
{
  Rng2I32 result = Zi;
  result.p0 = AddVec2I32(r.p0, v);
  result.p1 = AddVec2I32(r.p1, v);
  return result;
}

Rng2I32 MulRng2I32Vec2I32(Rng2I32 r, Vec2I32 v)
{
  Rng2I32 result = Zi;
  result.p0 = MulVec2I32Vec2I32(r.p0, v);
  result.p1 = MulVec2I32Vec2I32(r.p1, v);
  return result;
}

Rng2I32 DivRng2I32Vec2I32(Rng2I32 r, Vec2I32 v)
{
  Rng2I32 result = Zi;
  result.p0 = DivVec2I32Vec2I32(r.p0, v);
  result.p1 = DivVec2I32Vec2I32(r.p1, v);
  return result;
}

////////////////////////////////////////////////////////////
//~ Affine

b32 MatchAffine(Affine af0, Affine af1)
{
  return MatchVec2(af0.og, af1.og) && MatchVec2(af0.bx, af1.bx) && MatchVec2(af0.by, af1.by);
}

//- Initialization

Affine AffineFromTranslation(Vec2 v)
{
  Affine af;
  af.bx = VEC2(1, 0);
  af.by = VEC2(0, 1);
  af.og = v;
  return af;
}

Affine AffineFromRot(Vec2 r)
{
  Affine result;
  result.bx = VEC2(r.x, r.y);
  result.by = VEC2(-r.y, r.x);
  result.og = VEC2(0, 0);
  return result;
}

Affine AffineFromAngle(f32 r)
{
  Affine result;
  f32 c = CosF32(r);
  f32 s = SinF32(r);
  result.bx = VEC2(c, s);
  result.by = VEC2(-s, c);
  result.og = VEC2(0, 0);
  return result;
}

Affine AffineFromScale(Vec2 scale)
{
  Affine result;
  result.bx = VEC2(scale.x, 0);
  result.by = VEC2(0, scale.y);
  result.og = VEC2(0, 0);
  return result;
}

//- Translation

Affine TranslateAffine(Affine af, Vec2 v)
{
  af.og = VEC2(af.bx.x * v.x + af.by.x * v.y + af.og.x, af.bx.y * v.x + af.by.y * v.y + af.og.y);
  return af;
}

Affine PreTranslateAffine(Affine af, Vec2 v)
{
  af.og = AddVec2(af.og, v);
  return af;
}

//- Rotation

Affine RotateAffine(Affine af, Vec2 r)
{
  return MulAffine(af, AffineFromRot(r));
}

Affine PreRotateAffine(Affine af, Vec2 r)
{
  return MulAffine(AffineFromRot(r), af);
}

//- Scale

Affine ScaleAffine(Affine af, Vec2 scale)
{
  af.bx = MulVec2(af.bx, scale.x);
  af.by = MulVec2(af.by, scale.y);
  return af;
}

Affine PreScaleAffine(Affine af, Vec2 scale)
{
  Affine result;
  result.bx = MulVec2Vec2(af.bx, scale);
  result.by = MulVec2Vec2(af.by, scale);
  result.og = MulVec2Vec2(af.og, scale);
  return result;
}

//- Lerp

Affine LerpAffine(Affine a, Affine b, f32 t)
{
  Affine result;
  result.bx = SlerpVec2(a.bx, b.bx, t);
  result.by = SlerpVec2(a.by, b.by, t);
  result.og = LerpVec2(a.og, b.og, t);
  return result;
}

//- Invert

Affine InvertAffine(Affine af)
{
  f32 det = DeterminantFromAffine(af);
  f32 inv_det = 1.0f / det;

  f32 old_bx_x = af.bx.x;
  af.bx.x = af.by.y;
  af.by.y = old_bx_x;

  af.bx = MulVec2Vec2(af.bx, VEC2(inv_det, -inv_det));
  af.by = MulVec2Vec2(af.by, VEC2(-inv_det, inv_det));

  af.og = MulAffineBasisVec2(af, NegVec2(af.og));

  return af;
}

//- Mul

Affine MulAffine(Affine a, Affine b)
{
  Affine result;
  result.bx.x = a.bx.x * b.bx.x + a.by.x * b.bx.y;
  result.bx.y = a.bx.y * b.bx.x + a.by.y * b.bx.y;
  result.by.x = a.bx.x * b.by.x + a.by.x * b.by.y;
  result.by.y = a.bx.y * b.by.x + a.by.y * b.by.y;
  result.og = MulAffineVec2(a, b.og);
  return result;
}

Vec2 MulAffineBasisVec2(Affine af, Vec2 v)
{
  return VEC2(
    af.bx.x * v.x + af.by.x * v.y,
    af.bx.y * v.x + af.by.y * v.y
  );
}

Vec2 MulAffineVec2(Affine af, Vec2 v)
{
  Vec2 result = MulAffineBasisVec2(af, v);
  result = AddVec2(result, af.og);
  return result;
}

Affine MulAffineXform(Affine af, Xform xf)
{
  Affine result;
  result.bx.x = af.bx.x * xf.r.x + af.by.x * xf.r.y;
  result.bx.y = af.bx.y * xf.r.x + af.by.y * xf.r.y;
  result.by.x = af.bx.x * -xf.r.y + af.by.x * xf.r.x;
  result.by.y = af.bx.y * -xf.r.y + af.by.y * xf.r.x;
  result.og = MulAffineVec2(af, xf.t);
  return result;
}

//- Helpers

Affine BasisFromAffine(Affine af)
{
  Affine result = Zi;
  result.bx = af.bx;
  result.by = af.by;
  return result;
}

f32 DeterminantFromAffine(Affine af)
{
  return WedgeVec2(af.bx, af.by);
}

Vec2 RightFromAffine(Affine af)
{
  return af.bx;
}

Vec2 LeftFromAffine(Affine af)
{
  return NegVec2(af.bx);
}

Vec2 UpFromAffine(Affine af)
{
  return NegVec2(af.by);
}

Vec2 DownFromAffine(Affine af)
{
  return af.by;
}

f32 AngleFromAffine(Affine af)
{
  return AngleFromVec2(af.bx);
}

Vec2 ScaleFromAffine(Affine af)
{
  f32 det = DeterminantFromAffine(af);
  i32 det_sign = (det >= 0) - (det < 0);
  return VEC2(Vec2Len(af.bx), det_sign * Vec2Len(af.by));
}

////////////////////////////////////////////////////////////
//~ Xform

Xform XformTR(Vec2 t, Vec2 r)
{
  Xform xf;
  xf.t = t;
  xf.r = r;
  return xf;
}

Xform XformRT(Vec2 r, Vec2 t)
{
  Xform xf;
  xf.t = VEC2(0, 0);
  xf.r = r;
  xf = TranslateXform(xf, t);
  return xf;
}

Xform MulXform(Xform a, Xform b)
{
  Xform result;
  result.r.x = a.r.x * b.r.x + -a.r.y * b.r.y;
  result.r.y = a.r.y * b.r.x + a.r.x * b.r.y;
  result.t.x = a.r.x * b.t.x + -a.r.y * b.t.y;
  result.t.y = a.r.y * b.t.x + a.r.x * b.t.y;
  result.t = AddVec2(result.t, a.t);
  return result;
}

Xform TranslateXform(Xform xf, Vec2 t)
{
  xf.t = VEC2(xf.r.x * t.x + -xf.r.y * t.y + xf.t.x, xf.r.y * t.x + xf.r.x * t.y + xf.t.y);
  return xf;
}

Xform RotateXform(Xform xf, Vec2 r)
{
  xf.r = RotateVec2(xf.r, r);
  return xf;
}

Xform InvertXform(Xform xf)
{
  xf.t = NegVec2(xf.t);
  xf.r.y = -xf.r.y;
  return xf;
}

Vec2 MulXformVec2(Xform xf, Vec2 v)
{
  Vec2 result = Zi;
  result = AddVec2(RotateVec2(v, xf.r), xf.t);
  return result;
}

Xform NormXform(Xform xf)
{
  xf.r = NormVec2(xf.r);
  return xf;
}

////////////////////////////////////////////////////////////
//~ Line

Vec2 IntersectLines(Vec2 a0, Vec2 b0, Vec2 a1, Vec2 b1)
{
  Vec2 vab0 = SubVec2(b0, a0);
  Vec2 vab1 = SubVec2(b1, a1);
  Vec2 va0a1 = SubVec2(a1, a0);
  f32 denom = WedgeVec2(vab0, vab1);
  f32 t0 = WedgeVec2(va0a1, vab1) / denom;
  Vec2 result = AddVec2(a0, MulVec2(vab0, t0));
  return result;
}

////////////////////////////////////////////////////////////
//~ Spring

// https://box2d.org/files/ErinCatto_SoftConstraints_GDC2011.pdf
SoftSpring MakeSpring(f32 hertz, f32 damping_ratio, f32 dt)
{
  SoftSpring result;
  if (hertz == 0.0f)
  {
    result.bias_rate = 0;
    result.mass_scale = 1;
    result.impulse_scale = 0;
  }
  else
  {
    f32 angular_frequency = Tau * hertz;
    f32 a = 2 * damping_ratio + angular_frequency * dt;
    f32 b = angular_frequency * a * dt;
    f32 c = 1 / (b + 1);
    result.bias_rate = angular_frequency / a;
    result.mass_scale = b * c;
    result.impulse_scale = c;
  }
  return result;
}

////////////////////////////////////////////////////////////
//~ Mat4x4

Mat4x4 Mat4x4FromAffine(Affine af)
{
  return (Mat4x4)
  {
    .e = {
      {af.bx.x, af.bx.y, 0, 0},
      {af.by.x, af.by.y, 0, 0},
      {0,       0,       1, 0},
      {af.og.x, af.og.y, 0, 1},
    }
  };
}

Mat4x4 Mat4x4FromOrtho(f32 left, f32 right, f32 bottom, f32 top, f32 near_z, f32 far_z)
{
  Mat4x4 m = Zi;

  f32 rl = 1.0f / (right - left);
  f32 tb = 1.0f / (top - bottom);
  f32 fn = -1.0f / (far_z - near_z);

  m.e[0][0] = 2.0f * rl;
  m.e[1][1] = 2.0f * tb;
  m.e[2][2] = 2.0f * fn;
  m.e[3][0] = -(right + left) * rl;
  m.e[3][1] = -(top + bottom) * tb;
  m.e[3][2] = (far_z + near_z) * fn;
  m.e[3][3] = 1.0f;

  return m;
}

Mat4x4 MulMat4x4(Mat4x4 m0, Mat4x4 m1)
{
  f32 a00 = m0.e[0][0], a01 = m0.e[0][1], a02 = m0.e[0][2], a03 = m0.e[0][3],
  a10 = m0.e[1][0], a11 = m0.e[1][1], a12 = m0.e[1][2], a13 = m0.e[1][3],
  a20 = m0.e[2][0], a21 = m0.e[2][1], a22 = m0.e[2][2], a23 = m0.e[2][3],
  a30 = m0.e[3][0], a31 = m0.e[3][1], a32 = m0.e[3][2], a33 = m0.e[3][3],

  b00 = m1.e[0][0], b01 = m1.e[0][1], b02 = m1.e[0][2], b03 = m1.e[0][3],
  b10 = m1.e[1][0], b11 = m1.e[1][1], b12 = m1.e[1][2], b13 = m1.e[1][3],
  b20 = m1.e[2][0], b21 = m1.e[2][1], b22 = m1.e[2][2], b23 = m1.e[2][3],
  b30 = m1.e[3][0], b31 = m1.e[3][1], b32 = m1.e[3][2], b33 = m1.e[3][3];

  Mat4x4 result;
  result.e[0][0] = a00 * b00 + a10 * b01 + a20 * b02 + a30 * b03;
  result.e[0][1] = a01 * b00 + a11 * b01 + a21 * b02 + a31 * b03;
  result.e[0][2] = a02 * b00 + a12 * b01 + a22 * b02 + a32 * b03;
  result.e[0][3] = a03 * b00 + a13 * b01 + a23 * b02 + a33 * b03;
  result.e[1][0] = a00 * b10 + a10 * b11 + a20 * b12 + a30 * b13;
  result.e[1][1] = a01 * b10 + a11 * b11 + a21 * b12 + a31 * b13;
  result.e[1][2] = a02 * b10 + a12 * b11 + a22 * b12 + a32 * b13;
  result.e[1][3] = a03 * b10 + a13 * b11 + a23 * b12 + a33 * b13;
  result.e[2][0] = a00 * b20 + a10 * b21 + a20 * b22 + a30 * b23;
  result.e[2][1] = a01 * b20 + a11 * b21 + a21 * b22 + a31 * b23;
  result.e[2][2] = a02 * b20 + a12 * b21 + a22 * b22 + a32 * b23;
  result.e[2][3] = a03 * b20 + a13 * b21 + a23 * b22 + a33 * b23;
  result.e[3][0] = a00 * b30 + a10 * b31 + a20 * b32 + a30 * b33;
  result.e[3][1] = a01 * b30 + a11 * b31 + a21 * b32 + a31 * b33;
  result.e[3][2] = a02 * b30 + a12 * b31 + a22 * b32 + a32 * b33;
  result.e[3][3] = a03 * b30 + a13 * b31 + a23 * b32 + a33 * b33;

  return result;
}

Mat4x4 ProjectMat4x4View(Affine view, f32 viewport_width, f32 viewport_height)
{
  Mat4x4 projection = Mat4x4FromOrtho(0.0, viewport_width, viewport_height, 0.0, -1.0, 1.0);
  Mat4x4 view4x4 = Mat4x4FromAffine(view);
  return MulMat4x4(projection, view4x4);
}