zig

fork of https://codeberg.org/ziglang/zig
Log | Files | Refs | README | LICENSE

windowsnumerics.impl.h (48101B) - Raw


      1 /**
      2  * This file has no copyright assigned and is placed in the Public Domain.
      3  * This file is part of the mingw-w64 runtime package.
      4  * No warranty is given; refer to the file DISCLAIMER.PD within this package.
      5  */
      6 
      7 /**
      8  * Normal users should include `windowsnumerics.h` instead of this header.
      9  * However, the cppwinrt headers set `_WINDOWS_NUMERICS_NAMESPACE_`,
     10  * `_WINDOWS_NUMERICS_BEGIN_NAMESPACE_` and `_WINDOWS_NUMERICS_END_NAMESPACE_`
     11  * to custom values and include `windowsnumerics.impl.h`. Therefore this shall
     12  * be considered a public header, and these macros are public API.
     13 */
     14 
     15 
     16 #ifdef min
     17 #  pragma push_macro("min")
     18 #  undef min
     19 #  define _WINDOWS_NUMERICS_IMPL_PUSHED_MIN_
     20 #endif
     21 
     22 #ifdef max
     23 #  pragma push_macro("max")
     24 #  undef max
     25 #  define _WINDOWS_NUMERICS_IMPL_PUSHED_MAX_
     26 #endif
     27 
     28 #include <algorithm>
     29 #include <cmath>
     30 
     31 #include "directxmath.h"
     32 
     33 
     34 // === Internal macros ===
     35 #define _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(_ty1, _op, _ty2) \
     36   inline _ty1 &operator _op ## =(_ty1 &val1, _ty2 val2) { \
     37     val1 = operator _op (val1, val2); \
     38     return val1; \
     39   }
     40 
     41 
     42 // === Internal functions ===
     43 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
     44   namespace _impl {
     45 
     46 #if 0 && defined(__cpp_lib_clamp)
     47     using std::clamp;
     48 #else
     49     constexpr const float &clamp(const float &val, const float &min, const float &max) {
     50       return val < min ? min : (val > max ? max : val);
     51     }
     52 #endif
     53 
     54 #if 0 && defined(__cpp_lib_interpolate)
     55     using std::lerp;
     56 #else
     57     constexpr float lerp(float val1, float val2, float amount) {
     58       // Don't do (val2 - val1) * amount + val1 as it has worse precision.
     59       return val2 * amount + val1 * (1.0f - amount);
     60     }
     61 #endif
     62 
     63   }
     64 } _WINDOWS_NUMERICS_END_NAMESPACE_
     65 
     66 
     67 // === Forward decls ===
     68 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
     69 
     70   struct float2;
     71   struct float3;
     72   struct float4;
     73   struct float3x2;
     74   struct float4x4;
     75   struct plane;
     76   struct quaternion;
     77 
     78 } _WINDOWS_NUMERICS_END_NAMESPACE_
     79 
     80 
     81 // === float2: Struct and function defs ===
     82 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
     83 
     84   struct float2 {
     85     float2() = default;
     86     constexpr float2(float x, float y)
     87       : x(x), y(y)
     88     {}
     89     constexpr explicit float2(float val)
     90       : x(val), y(val)
     91     {}
     92 
     93     static constexpr float2 zero() {
     94       return float2(0.0f);
     95     }
     96     static constexpr float2 one() {
     97       return float2(1.0f);
     98     }
     99     static constexpr float2 unit_x() {
    100       return { 1.0f, 0.0f };
    101     }
    102     static constexpr float2 unit_y() {
    103       return { 0.0f, 1.0f };
    104     }
    105 
    106     float x;
    107     float y;
    108   };
    109 
    110   // Forward decl functions
    111   inline float length(const float2 &val);
    112   inline float length_squared(const float2 &val);
    113   inline float distance(const float2 &val1, const float2 &val2);
    114   inline float distance_squared(const float2 &val1, const float2 &val2);
    115   inline float dot(const float2 &val1, const float2 &val2);
    116   inline float2 normalize(const float2 &val);
    117   inline float2 reflect(const float2 &vec, const float2 &norm);
    118   inline float2 min(const float2 &val1, const float2 &val2);
    119   inline float2 max(const float2 &val1, const float2 &val2);
    120   inline float2 clamp(const float2 &val, const float2 &min, const float2 &max);
    121   inline float2 lerp(const float2 &val1, const float2 &val2, float amount);
    122   inline float2 transform(const float2 &pos, const float3x2 &mat);
    123   inline float2 transform(const float2 &pos, const float4x4 &mat);
    124   inline float2 transform_normal(const float2 &norm, const float3x2 &mat);
    125   inline float2 transform_normal(const float2 &norm, const float4x4 &mat);
    126   inline float2 transform(const float2 &val, const quaternion &rot);
    127 
    128   // Define operators
    129 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    130   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    131     return { val1.x _op val2.x, val1.y _op val2.y }; \
    132   }
    133   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float2, +)
    134   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float2, -)
    135   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float2, *)
    136   inline float2 operator*(const float2 &val1, float val2) {
    137     return { val1.x * val2, val1.y * val2 };
    138   }
    139   inline float2 operator*(float val1, const float2 &val2) {
    140     return operator*(val2, val1);
    141   }
    142   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float2, /)
    143   inline float2 operator/(const float2 &val1, float val2) {
    144     return operator*(val1, 1.0f / val2);
    145   }
    146   inline float2 operator-(const float2 &val) {
    147     return { -val.x, -val.y };
    148   }
    149 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    150   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, +, const float2 &)
    151   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, -, const float2 &)
    152   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, *, const float2 &)
    153   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, *, float)
    154   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, /, const float2 &)
    155   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float2, /, float)
    156   inline bool operator==(const float2 &val1, const float2 &val2) {
    157     return val1.x == val2.x && val1.y == val2.y;
    158   }
    159   inline bool operator!=(const float2 &val1, const float2 &val2) {
    160     return !operator==(val1, val2);
    161   }
    162 
    163 } _WINDOWS_NUMERICS_END_NAMESPACE_
    164 
    165 
    166 // === float3: Struct and function defs ===
    167 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    168 
    169   struct float3 {
    170     float3() = default;
    171     constexpr float3(float x, float y, float z)
    172       : x(x), y(y), z(z)
    173     {}
    174     constexpr float3(float2 val, float z)
    175       : x(val.x), y(val.y), z(z)
    176     {}
    177     constexpr explicit float3(float val)
    178       : x(val), y(val), z(val)
    179     {}
    180 
    181     static constexpr float3 zero() {
    182       return float3(0.0);
    183     }
    184     static constexpr float3 one() {
    185       return float3(1.0);
    186     }
    187     static constexpr float3 unit_x() {
    188       return { 1.0f, 0.0f, 0.0f };
    189     }
    190     static constexpr float3 unit_y() {
    191       return { 0.0f, 1.0f, 0.0f };
    192     }
    193     static constexpr float3 unit_z() {
    194       return { 0.0f, 0.0f, 1.0f };
    195     }
    196 
    197     float x;
    198     float y;
    199     float z;
    200   };
    201 
    202   // Forward decl functions
    203   inline float length(const float3 &val);
    204   inline float length_squared(const float3 &val);
    205   inline float distance(const float3 &val1, const float3 &val2);
    206   inline float distance_squared(const float3 &val1, const float3 &val2);
    207   inline float dot(const float3 &val1, const float3 &val2);
    208   inline float3 normalize(const float3 &val);
    209   inline float3 cross(const float3 &val1, const float3 &val2);
    210   inline float3 reflect(const float3 &vec, const float3 &norm);
    211   inline float3 min(const float3 &val1, const float3 &val2);
    212   inline float3 max(const float3 &val1, const float3 &val2);
    213   inline float3 clamp(const float3 &val, const float3 &min, const float3 &max);
    214   inline float3 lerp(const float3 &val1, const float3 &val2, float amount);
    215   inline float3 transform(const float3 &pos, const float4x4 &mat);
    216   inline float3 transform_normal(const float3 &norm, const float4x4 &mat);
    217   inline float3 transform(const float3 &val, const quaternion &rot);
    218 
    219   // Define operators
    220 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    221   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    222     return { val1.x _op val2.x, val1.y _op val2.y, val1.z _op val2.z }; \
    223   }
    224   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3, +)
    225   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3, -)
    226   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3, *)
    227   inline float3 operator*(const float3 &val1, float val2) {
    228     return { val1.x * val2, val1.y * val2, val1.z * val2 };
    229   }
    230   inline float3 operator*(float val1, const float3 &val2) {
    231     return operator*(val2, val1);
    232   }
    233   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3, /)
    234   inline float3 operator/(const float3 &val1, float val2) {
    235     return operator*(val1, 1.0f / val2);
    236   }
    237   inline float3 operator-(const float3 &val) {
    238     return { -val.x, -val.y, -val.z };
    239   }
    240 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    241   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, +, const float3 &)
    242   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, -, const float3 &)
    243   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, *, const float3 &)
    244   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, *, float)
    245   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, /, const float3 &)
    246   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3, /, float)
    247   inline bool operator==(const float3 &val1, const float3 &val2) {
    248     return val1.x == val2.x && val1.y == val2.y && val1.z == val2.z;
    249   }
    250   inline bool operator!=(const float3 &val1, const float3 &val2) {
    251     return !operator==(val1, val2);
    252   }
    253 
    254 } _WINDOWS_NUMERICS_END_NAMESPACE_
    255 
    256 
    257 // === float4: Struct and function defs ===
    258 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    259 
    260   struct float4 {
    261     float4() = default;
    262     constexpr float4(float x, float y, float z, float w)
    263       : x(x), y(y), z(z), w(w)
    264     {}
    265     constexpr float4(float2 val, float z, float w)
    266       : x(val.x), y(val.y), z(z), w(w)
    267     {}
    268     constexpr float4(float3 val, float w)
    269       : x(val.x), y(val.y), z(val.z), w(w)
    270     {}
    271     constexpr explicit float4(float val)
    272       : x(val), y(val), z(val), w(val)
    273     {}
    274 
    275     static constexpr float4 zero() {
    276       return float4(0.0);
    277     }
    278     static constexpr float4 one() {
    279       return float4(1.0);
    280     }
    281     static constexpr float4 unit_x() {
    282       return { 1.0f, 0.0f, 0.0f, 0.0f };
    283     }
    284     static constexpr float4 unit_y() {
    285       return { 0.0f, 1.0f, 0.0f, 0.0f };
    286     }
    287     static constexpr float4 unit_z() {
    288       return { 0.0f, 0.0f, 1.0f, 0.0f };
    289     }
    290     static constexpr float4 unit_w() {
    291       return { 0.0f, 0.0f, 0.0f, 1.0f };
    292     }
    293 
    294     float x;
    295     float y;
    296     float z;
    297     float w;
    298   };
    299 
    300   // Forward decl functions
    301   inline float length(const float4 &val);
    302   inline float length_squared(const float4 &val);
    303   inline float distance(const float4 &val1, const float4 &val2);
    304   inline float distance_squared(const float4 &val1, const float4 &val2);
    305   inline float dot(const float4 &val1, const float4 &val2);
    306   inline float4 normalize(const float4 &val);
    307   inline float4 min(const float4 &val1, const float4 &val2);
    308   inline float4 max(const float4 &val1, const float4 &val2);
    309   inline float4 clamp(const float4 &val, const float4 &min, const float4 &max);
    310   inline float4 lerp(const float4 &val1, const float4 &val2, float amount);
    311   inline float4 transform(const float4 &pos, const float4x4 &mat);
    312   inline float4 transform4(const float3 &pos, const float4x4 &mat);
    313   inline float4 transform4(const float2 &pos, const float4x4 &mat);
    314   inline float4 transform(const float4 &val, const quaternion &rot);
    315   inline float4 transform4(const float3 &val, const quaternion &rot);
    316   inline float4 transform4(const float2 &val, const quaternion &rot);
    317 
    318   // Define operators
    319 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    320   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    321     return { val1.x _op val2.x, val1.y _op val2.y, val1.z _op val2.z, val1.w _op val2.w }; \
    322   }
    323   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4, +)
    324   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4, -)
    325   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4, *)
    326   inline float4 operator*(const float4 &val1, float val2) {
    327     return { val1.x * val2, val1.y * val2, val1.z * val2, val1.w * val2 };
    328   }
    329   inline float4 operator*(float val1, const float4 &val2) {
    330     return operator*(val2, val1);
    331   }
    332   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4, /)
    333   inline float4 operator/(const float4 &val1, float val2) {
    334     return operator*(val1, 1.0f / val2);
    335   }
    336   inline float4 operator-(const float4 &val) {
    337     return { -val.x, -val.y, -val.z, -val.w };
    338   }
    339 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    340   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, +, const float4 &)
    341   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, -, const float4 &)
    342   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, *, const float4 &)
    343   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, *, float)
    344   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, /, const float4 &)
    345   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4, /, float)
    346   inline bool operator==(const float4 &val1, const float4 &val2) {
    347     return val1.x == val2.x && val1.y == val2.y && val1.z == val2.z && val2.w == val2.w;
    348   }
    349   inline bool operator!=(const float4 &val1, const float4 &val2) {
    350     return !operator==(val1, val2);
    351   }
    352 
    353 } _WINDOWS_NUMERICS_END_NAMESPACE_
    354 
    355 
    356 // === float3x2: Struct and function defs ===
    357 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    358 
    359   struct float3x2 {
    360     float3x2() = default;
    361     constexpr float3x2(
    362       float m11, float m12,
    363       float m21, float m22,
    364       float m31, float m32
    365     )
    366       : m11(m11), m12(m12)
    367       , m21(m21), m22(m22)
    368       , m31(m31), m32(m32)
    369     {}
    370 
    371     static constexpr float3x2 identity() {
    372       return {
    373         1.0f, 0.0f,
    374         0.0f, 1.0f,
    375         0.0f, 0.0f
    376       };
    377     }
    378 
    379     float m11; float m12;
    380     float m21; float m22;
    381     float m31; float m32;
    382   };
    383 
    384   // Forward decl functions
    385   inline float3x2 make_float3x2_translation(const float2 &pos);
    386   inline float3x2 make_float3x2_translation(float xpos, float ypos);
    387   inline float3x2 make_float3x2_scale(float xscale, float yscale);
    388   inline float3x2 make_float3x2_scale(float xscale, float yscale, const float2 &center);
    389   inline float3x2 make_float3x2_scale(const float2 &xyscale);
    390   inline float3x2 make_float3x2_scale(const float2 &xyscale, const float2 &center);
    391   inline float3x2 make_float3x2_scale(float scale);
    392   inline float3x2 make_float3x2_scale(float scale, const float2 &center);
    393   inline float3x2 make_float3x2_skew(float xrad, float yrad);
    394   inline float3x2 make_float3x2_skew(float xrad, float yrad, const float2 &center);
    395   inline float3x2 make_float3x2_rotation(float rad);
    396   inline float3x2 make_float3x2_rotation(float rad, const float2 &center);
    397   inline bool is_identity(const float3x2 &val);
    398   inline float determinant(const float3x2 &val);
    399   inline float2 translation(const float3x2 &val);
    400   inline bool invert(const float3x2 &val, float3x2 *out);
    401   inline float3x2 lerp(const float3x2 &mat1, const float3x2 &mat2, float amount);
    402 
    403   // Define operators
    404 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    405   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    406     return { \
    407       val1.m11 _op val2.m11, val1.m12 _op val2.m12, \
    408       val1.m21 _op val2.m21, val1.m22 _op val2.m22, \
    409       val1.m31 _op val2.m31, val1.m32 _op val2.m32, \
    410     }; \
    411   }
    412   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3x2, +)
    413   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float3x2, -)
    414   inline float3x2 operator*(const float3x2 &val1, const float3x2 &val2) {
    415     // 2D transformation matrix has an implied 3rd column with (0, 0, 1)
    416     const float3 v1r1(val1.m11, val1.m12, 0.0f);
    417     const float3 v1r2(val1.m21, val1.m22, 0.0f);
    418     const float3 v1r3(val1.m31, val1.m32, 1.0f);
    419     const float3 v2c1(val2.m11, val2.m21, val2.m31);
    420     const float3 v2c2(val2.m12, val2.m22, val2.m32);
    421     // const float3 v2c3(0.0f, 0.0f, 1.0f);
    422     return {
    423       dot(v1r1, v2c1), dot(v1r1, v2c2),
    424       dot(v1r2, v2c1), dot(v1r2, v2c2),
    425       dot(v1r3, v2c1), dot(v1r3, v2c2)
    426     };
    427   }
    428   inline float3x2 operator*(const float3x2 &val1, float val2) {
    429     return {
    430       val1.m11 * val2, val1.m12 * val2,
    431       val1.m21 * val2, val1.m22 * val2,
    432       val1.m31 * val2, val1.m32 * val2
    433     };
    434   }
    435   inline float3x2 operator-(const float3x2 &val) {
    436     return {
    437       -val.m11, -val.m12,
    438       -val.m21, -val.m22,
    439       -val.m31, -val.m32
    440     };
    441   }
    442 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    443   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3x2, +, const float3x2 &)
    444   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3x2, -, const float3x2 &)
    445   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3x2, *, const float3x2 &)
    446   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float3x2, *, float)
    447   inline bool operator==(const float3x2 &val1, const float3x2 &val2) {
    448     return
    449       val1.m11 == val2.m11 && val1.m12 == val2.m12 &&
    450       val1.m21 == val2.m21 && val1.m22 == val2.m22 &&
    451       val1.m31 == val2.m31 && val1.m32 == val2.m32;
    452   }
    453   inline bool operator!=(const float3x2 &val1, const float3x2 &val2) {
    454     return !operator==(val1, val2);
    455   }
    456 
    457 } _WINDOWS_NUMERICS_END_NAMESPACE_
    458 
    459 
    460 // === float4x4: Struct and function defs ===
    461 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    462 
    463   struct float4x4 {
    464     float4x4() = default;
    465     constexpr float4x4(
    466       float m11, float m12, float m13, float m14,
    467       float m21, float m22, float m23, float m24,
    468       float m31, float m32, float m33, float m34,
    469       float m41, float m42, float m43, float m44
    470     )
    471       : m11(m11), m12(m12), m13(m13), m14(m14)
    472       , m21(m21), m22(m22), m23(m23), m24(m24)
    473       , m31(m31), m32(m32), m33(m33), m34(m34)
    474       , m41(m41), m42(m42), m43(m43), m44(m44)
    475     {}
    476 
    477     static constexpr float4x4 identity() {
    478       return {
    479         1.0f, 0.0f, 0.0f, 0.0f,
    480         0.0f, 1.0f, 0.0f, 0.0f,
    481         0.0f, 0.0f, 1.0f, 0.0f,
    482         0.0f, 0.0f, 0.0f, 1.0f
    483       };
    484     }
    485 
    486     float m11; float m12; float m13; float m14;
    487     float m21; float m22; float m23; float m24;
    488     float m31; float m32; float m33; float m34;
    489     float m41; float m42; float m43; float m44;
    490   };
    491 
    492   // Forward decl functions
    493   inline float4x4 make_float4x4_billboard(const float3 &objpos, const float3 &camerapos, const float3 &cameraup, const float3 &camerafwd);
    494   inline float4x4 make_float4x4_constrained_billboard(const float3 &objpos, const float3 &camerapos, const float3 &rotateaxis, const float3 &camerafwd, const float3 &objfwd);
    495   inline float4x4 make_float4x4_translation(const float3 &pos);
    496   inline float4x4 make_float4x4_translation(float xpos, float ypos, float zpos);
    497   inline float4x4 make_float4x4_scale(float xscale, float yscale, float zscale);
    498   inline float4x4 make_float4x4_scale(float xscale, float yscale, float zscale, const float3 &center);
    499   inline float4x4 make_float4x4_scale(const float3 &xyzscale);
    500   inline float4x4 make_float4x4_scale(const float3 &xyzscale, const float3 &center);
    501   inline float4x4 make_float4x4_scale(float scale);
    502   inline float4x4 make_float4x4_scale(float scale, const float3 &center);
    503   inline float4x4 make_float4x4_rotation_x(float rad);
    504   inline float4x4 make_float4x4_rotation_x(float rad, const float3 &center);
    505   inline float4x4 make_float4x4_rotation_y(float rad);
    506   inline float4x4 make_float4x4_rotation_y(float rad, const float3 &center);
    507   inline float4x4 make_float4x4_rotation_z(float rad);
    508   inline float4x4 make_float4x4_rotation_z(float rad, const float3 &center);
    509   inline float4x4 make_float4x4_from_axis_angle(const float3 &axis, float angle);
    510   inline float4x4 make_float4x4_perspective_field_of_view(float fov, float aspect, float nearplane, float farplane);
    511   inline float4x4 make_float4x4_perspective(float w, float h, float nearplane, float farplane);
    512   inline float4x4 make_float4x4_perspective_off_center(float left, float right, float bottom, float top, float nearplane, float farplane);
    513   inline float4x4 make_float4x4_orthographic(float w, float h, float znearplane, float zfarplane);
    514   inline float4x4 make_float4x4_orthographic_off_center(float left, float right, float bottom, float top, float znearplane, float zfarplane);
    515   inline float4x4 make_float4x4_look_at(const float3 &camerapos, const float3 &target, const float3 &cameraup);
    516   inline float4x4 make_float4x4_world(const float3 &pos, const float3 &fwd, const float3 &up);
    517   inline float4x4 make_float4x4_from_quaternion(const quaternion &quat);
    518   inline float4x4 make_float4x4_from_yaw_pitch_roll(float yaw, float pitch, float roll);
    519   inline float4x4 make_float4x4_shadow(const float3 &lightdir, const plane &plane);
    520   inline float4x4 make_float4x4_reflection(const plane &plane);
    521   inline bool is_identity(const float4x4 &val);
    522   inline float determinant(const float4x4 &val);
    523   inline float3 translation(const float4x4 &val);
    524   inline bool invert(const float4x4 &mat, float4x4 *out);
    525   inline bool decompose(const float4x4 &mat, float3 *out_scale, quaternion *out_rot, float3 *out_translate);
    526   inline float4x4 transform(const float4x4 &val, const quaternion &rot);
    527   inline float4x4 transpose(const float4x4 &val);
    528   inline float4x4 lerp(const float4x4 &val1, const float4x4 &val2, float amount);
    529 
    530   // Define operators
    531 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    532   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    533     return { \
    534       val1.m11 _op val2.m11, val1.m12 _op val2.m12, val1.m13 _op val2.m13, val1.m14 _op val2.m14, \
    535       val1.m21 _op val2.m21, val1.m22 _op val2.m22, val1.m23 _op val2.m23, val1.m24 _op val2.m24, \
    536       val1.m31 _op val2.m31, val1.m32 _op val2.m32, val1.m33 _op val2.m33, val1.m34 _op val2.m34, \
    537       val1.m41 _op val2.m41, val1.m42 _op val2.m42, val1.m43 _op val2.m43, val1.m44 _op val2.m44, \
    538     }; \
    539   }
    540   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4x4, +)
    541   _WINDOWS_NUMERICS_IMPL_BINARY_OP(float4x4, -)
    542   inline float4x4 operator*(const float4x4 &val1, const float4x4 &val2) {
    543     const float4 v1r1(val1.m11, val1.m12, val1.m13, val1.m14);
    544     const float4 v1r2(val1.m21, val1.m22, val1.m23, val1.m24);
    545     const float4 v1r3(val1.m31, val1.m32, val1.m33, val1.m34);
    546     const float4 v1r4(val1.m41, val1.m42, val1.m43, val1.m44);
    547     const float4 v2c1(val2.m11, val2.m21, val2.m31, val2.m41);
    548     const float4 v2c2(val2.m12, val2.m22, val2.m32, val2.m42);
    549     const float4 v2c3(val2.m13, val2.m23, val2.m33, val2.m43);
    550     const float4 v2c4(val2.m14, val2.m24, val2.m34, val2.m44);
    551     return {
    552       dot(v1r1, v2c1), dot(v1r1, v2c2), dot(v1r1, v2c3), dot(v1r1, v2c4),
    553       dot(v1r2, v2c1), dot(v1r2, v2c2), dot(v1r2, v2c3), dot(v1r2, v2c4),
    554       dot(v1r3, v2c1), dot(v1r3, v2c2), dot(v1r3, v2c3), dot(v1r3, v2c4),
    555       dot(v1r4, v2c1), dot(v1r4, v2c2), dot(v1r4, v2c3), dot(v1r4, v2c4)
    556     };
    557   }
    558   inline float4x4 operator*(const float4x4 &val1, float val2) {
    559     return {
    560       val1.m11 * val2, val1.m12 * val2, val1.m13 * val2, val1.m14 * val2,
    561       val1.m21 * val2, val1.m22 * val2, val1.m23 * val2, val1.m24 * val2,
    562       val1.m31 * val2, val1.m32 * val2, val1.m33 * val2, val1.m34 * val2,
    563       val1.m41 * val2, val1.m42 * val2, val1.m43 * val2, val1.m44 * val2
    564     };
    565   }
    566   inline float4x4 operator-(const float4x4 &val) {
    567     return {
    568       -val.m11, -val.m12, -val.m13, -val.m14,
    569       -val.m21, -val.m22, -val.m23, -val.m24,
    570       -val.m31, -val.m32, -val.m33, -val.m34,
    571       -val.m41, -val.m42, -val.m43, -val.m44
    572     };
    573   }
    574 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    575   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4x4, +, const float4x4 &)
    576   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4x4, -, const float4x4 &)
    577   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4x4, *, const float4x4 &)
    578   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(float4x4, *, float)
    579   inline bool operator==(const float4x4 &val1, const float4x4 &val2) {
    580     return
    581       val1.m11 == val2.m11 && val1.m12 == val2.m12 && val1.m13 == val2.m13 && val1.m14 == val2.m14 &&
    582       val1.m21 == val2.m21 && val1.m22 == val2.m22 && val1.m23 == val2.m23 && val1.m24 == val2.m24 &&
    583       val1.m31 == val2.m31 && val1.m32 == val2.m32 && val1.m33 == val2.m33 && val1.m34 == val2.m34 &&
    584       val1.m41 == val2.m41 && val1.m42 == val2.m42 && val1.m43 == val2.m43 && val1.m44 == val2.m44;
    585   }
    586   inline bool operator!=(const float4x4 &val1, const float4x4 &val2) {
    587     return !operator==(val1, val2);
    588   }
    589 
    590 } _WINDOWS_NUMERICS_END_NAMESPACE_
    591 
    592 
    593 // === plane: Struct and function defs ===
    594 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    595 
    596   struct plane {
    597     plane() = default;
    598     constexpr plane(float x, float y, float z, float d)
    599       : normal(float3(x, y, z)), d(d)
    600     {}
    601     constexpr plane(float3 normal, float d)
    602       : normal(normal), d(d)
    603     {}
    604     constexpr explicit plane(float4 val)
    605       : normal(float3(val.x, val.y, val.z)), d(val.w)
    606     {}
    607 
    608     float3 normal;
    609     float d;
    610   };
    611 
    612   // Forward decl functions
    613   inline plane make_plane_from_vertices(const float3 &pt1, const float3 &pt2, const float3 &pt3);
    614   inline plane normalize(const plane &val);
    615   inline plane transform(const plane &plane, const float4x4 &mat);
    616   inline plane transform(const plane &plane, const quaternion &rot);
    617   inline float dot(const plane &plane, const float4 &val);
    618   inline float dot_coordinate(const plane &plane, const float3 &val);
    619   inline float dot_normal(const plane &plane, const float3 &val);
    620 
    621   // Define operators
    622   inline bool operator==(const plane &val1, const plane &val2) {
    623     return val1.normal == val2.normal && val1.d == val2.d;
    624   }
    625   inline bool operator!=(const plane &val1, const plane &val2) {
    626     return !operator==(val1, val2);
    627   }
    628 
    629 } _WINDOWS_NUMERICS_END_NAMESPACE_
    630 
    631 
    632 // === quaternion: Struct and function defs ===
    633 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    634 
    635   struct quaternion {
    636     quaternion() = default;
    637     constexpr quaternion(float x, float y, float z, float w)
    638       : x(x), y(y), z(z), w(w)
    639     {}
    640     constexpr quaternion(float3 vecPart, float scalarPart)
    641       : x(vecPart.x), y(vecPart.y), z(vecPart.z), w(scalarPart)
    642     {}
    643 
    644     static constexpr quaternion identity() {
    645       return { 0.0f, 0.0f, 0.0f, 1.0f };
    646     }
    647 
    648     float x;
    649     float y;
    650     float z;
    651     float w;
    652   };
    653 
    654   // Forward decl functions
    655   inline quaternion make_quaternion_from_axis_angle(const float3 &axis, float angle);
    656   inline quaternion make_quaternion_from_yaw_pitch_roll(float yaw, float pitch, float roll);
    657   inline quaternion make_quaternion_from_rotation_matrix(const float4x4 &mat);
    658   inline bool is_identity(const quaternion &val);
    659   inline float length(const quaternion &val);
    660   inline float length_squared(const quaternion &val);
    661   inline float dot(const quaternion &val1, const quaternion &val2);
    662   inline quaternion normalize(const quaternion &val);
    663   inline quaternion conjugate(const quaternion &val);
    664   inline quaternion inverse(const quaternion &val);
    665   inline quaternion slerp(const quaternion &val1, const quaternion &val2, float amount);
    666   inline quaternion lerp(const quaternion &val1, const quaternion &val2, float amount);
    667   inline quaternion concatenate(const quaternion &val1, const quaternion &val2);
    668 
    669   // Define operators
    670 #define _WINDOWS_NUMERICS_IMPL_BINARY_OP(_ty, _op) \
    671   inline _ty operator _op(const _ty &val1, const _ty &val2) { \
    672     return { val1.x _op val2.x, val1.y _op val2.y, val1.z _op val2.z, val1.w _op val2.w }; \
    673   }
    674   _WINDOWS_NUMERICS_IMPL_BINARY_OP(quaternion, +)
    675   _WINDOWS_NUMERICS_IMPL_BINARY_OP(quaternion, -)
    676   inline quaternion operator*(const quaternion &val1, const quaternion &val2) {
    677     return {
    678       val1.w * val2.x + val1.x * val2.w + val1.y * val2.z - val1.z * val2.y,
    679       val1.w * val2.y - val1.x * val2.z + val1.y * val2.w + val1.z * val2.x,
    680       val1.w * val2.z + val1.x * val2.y - val1.y * val2.x + val1.z * val2.w,
    681       val1.w * val2.w - val1.x * val2.x - val1.y * val2.y - val1.z * val2.z
    682     }; }
    683   inline quaternion operator*(const quaternion &val1, float val2) {
    684     return { val1.x * val2, val1.y * val2, val1.z * val2, val1.w * val2 };
    685   }
    686   inline quaternion operator/(const quaternion &val1, const quaternion &val2) {
    687     return operator*(val1, inverse(val2));
    688   }
    689   inline quaternion operator-(const quaternion &val) {
    690     return { -val.x, -val.y, -val.z, -val.w };
    691   }
    692 #undef _WINDOWS_NUMERICS_IMPL_BINARY_OP
    693   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(quaternion, +, const quaternion &)
    694   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(quaternion, -, const quaternion &)
    695   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(quaternion, *, const quaternion &)
    696   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(quaternion, *, float)
    697   _WINDOWS_NUMERICS_IMPL_ASSIGN_OP(quaternion, /, const quaternion &)
    698   inline bool operator==(const quaternion &val1, const quaternion &val2) {
    699     return val1.x == val2.x && val1.y == val2.y && val1.z == val2.z && val2.w == val2.w;
    700   }
    701   inline bool operator!=(const quaternion &val1, const quaternion &val2) {
    702     return !operator==(val1, val2);
    703   }
    704 
    705 } _WINDOWS_NUMERICS_END_NAMESPACE_
    706 
    707 
    708 // === Function definitions ===
    709 
    710 // Define float2 functions
    711 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    712 
    713   inline float length(const float2 &val) {
    714     return ::std::sqrt(length_squared(val));
    715   }
    716   inline float length_squared(const float2 &val) {
    717     return val.x * val.x + val.y * val.y;
    718   }
    719   inline float distance(const float2 &val1, const float2 &val2) {
    720     return length(val2 - val1);
    721   }
    722   inline float distance_squared(const float2 &val1, const float2 &val2) {
    723     return length_squared(val2 - val1);
    724   }
    725   inline float dot(const float2 &val1, const float2 &val2) {
    726     return val1.x * val2.x + val1.y * val2.y;
    727   }
    728   inline float2 normalize(const float2 &val) {
    729     return val / length(val);
    730   }
    731   inline float2 reflect(const float2 &vec, const float2 &norm) {
    732     // norm is assumed to be normalized.
    733     return vec - 2.0f * dot(vec, norm) * norm;
    734   }
    735   inline float2 min(const float2 &val1, const float2 &val2) {
    736     return { ::std::min(val1.x, val2.x), ::std::min(val1.y, val2.y) };
    737   }
    738   inline float2 max(const float2 &val1, const float2 &val2) {
    739     return { ::std::max(val1.x, val2.x), ::std::max(val1.y, val2.y) };
    740   }
    741   inline float2 clamp(const float2 &val, const float2 &min, const float2 &max) {
    742     return { _impl::clamp(val.x, min.x, max.x), _impl::clamp(val.y, min.y, max.y) };
    743   }
    744   inline float2 lerp(const float2 &val1, const float2 &val2, float amount) {
    745     return { _impl::lerp(val1.x, val2.x, amount), _impl::lerp(val1.y, val2.y, amount) };
    746   }
    747   // TODO: impl
    748   inline float2 transform(const float2 &pos, const float3x2 &mat);
    749   inline float2 transform(const float2 &pos, const float4x4 &mat);
    750   inline float2 transform_normal(const float2 &norm, const float3x2 &mat);
    751   inline float2 transform_normal(const float2 &norm, const float4x4 &mat);
    752   inline float2 transform(const float2 &val, const quaternion &rot) {
    753     // See comments in the float3 transform function.
    754     quaternion result = rot * quaternion(float3(val.x, val.y, 0.0f), 0.0f) * inverse(rot);
    755     return { result.x, result.y };
    756   }
    757 
    758 } _WINDOWS_NUMERICS_END_NAMESPACE_
    759 
    760 
    761 // Define float3 functions
    762 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    763 
    764   inline float length(const float3 &val) {
    765     return ::std::sqrt(length_squared(val));
    766   }
    767   inline float length_squared(const float3 &val) {
    768     return val.x * val.x + val.y * val.y + val.z * val.z;
    769   }
    770   inline float distance(const float3 &val1, const float3 &val2) {
    771     return length(val2 - val1);
    772   }
    773   inline float distance_squared(const float3 &val1, const float3 &val2) {
    774     return length_squared(val2 - val1);
    775   }
    776   inline float dot(const float3 &val1, const float3 &val2) {
    777     return val1.x * val2.x + val1.y * val2.y + val1.z * val2.z;
    778   }
    779   inline float3 normalize(const float3 &val) {
    780     return val / length(val);
    781   }
    782   inline float3 cross(const float3 &val1, const float3 &val2) {
    783     return {
    784       val1.y * val2.z - val2.y * val1.z,
    785       val1.z * val2.x - val2.z * val1.x,
    786       val1.x * val2.y - val2.x * val1.y
    787     };
    788   }
    789   inline float3 reflect(const float3 &vec, const float3 &norm) {
    790     // norm is assumed to be normalized.
    791     return vec - 2.0f * dot(vec, norm) * norm;
    792   }
    793   inline float3 min(const float3 &val1, const float3 &val2) {
    794     return { ::std::min(val1.x, val2.x), ::std::min(val1.y, val2.y), ::std::min(val1.z, val2.z) };
    795   }
    796   inline float3 max(const float3 &val1, const float3 &val2) {
    797     return { ::std::max(val1.x, val2.x), ::std::max(val1.y, val2.y), ::std::max(val1.z, val2.z) };
    798   }
    799   inline float3 clamp(const float3 &val, const float3 &min, const float3 &max) {
    800     return { _impl::clamp(val.x, min.x, max.x), _impl::clamp(val.y, min.y, max.y), _impl::clamp(val.z, min.z, max.z) };
    801   }
    802   inline float3 lerp(const float3 &val1, const float3 &val2, float amount) {
    803     return { _impl::lerp(val1.x, val2.x, amount), _impl::lerp(val1.y, val2.y, amount), _impl::lerp(val1.z, val2.z, amount) };
    804   }
    805   // TODO: impl
    806   inline float3 transform(const float3 &pos, const float4x4 &mat);
    807   inline float3 transform_normal(const float3 &norm, const float4x4 &mat);
    808   inline float3 transform(const float3 &val, const quaternion &rot) {
    809     // https://en.wikipedia.org/wiki/Quaternions_and_spatial_rotation#Using_quaternions_as_rotations
    810     // If assuming rot is a unit quaternion, this could use
    811     // conjugate() instead of inverse() too.
    812     // This can be expressed as a matrix operation too, with
    813     // https://en.wikipedia.org/wiki/Quaternions_and_spatial_rotation#Quaternion-derived_rotation_matrix
    814     // (see make_float4x4_from_quaternion).
    815     quaternion result = rot * quaternion(val, 0.0f) * inverse(rot);
    816     return { result.x, result.y, result.z };
    817   }
    818 
    819 } _WINDOWS_NUMERICS_END_NAMESPACE_
    820 
    821 
    822 // Define float4 functions
    823 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    824 
    825   inline float length(const float4 &val) {
    826     return ::std::sqrt(length_squared(val));
    827   }
    828   inline float length_squared(const float4 &val) {
    829     return val.x * val.x + val.y * val.y + val.z * val.z + val.w * val.w;
    830   }
    831   inline float distance(const float4 &val1, const float4 &val2) {
    832     return length(val2 - val1);
    833   }
    834   inline float distance_squared(const float4 &val1, const float4 &val2) {
    835     return length_squared(val2 - val1);
    836   }
    837   inline float dot(const float4 &val1, const float4 &val2) {
    838     return val1.x * val2.x + val1.y * val2.y + val1.z * val2.z + val1.w * val2.w;
    839   }
    840   inline float4 normalize(const float4 &val) {
    841     return val / length(val);
    842   }
    843   inline float4 min(const float4 &val1, const float4 &val2) {
    844     return {
    845       ::std::min(val1.x, val2.x),
    846       ::std::min(val1.y, val2.y),
    847       ::std::min(val1.z, val2.z),
    848       ::std::min(val1.w, val2.w)
    849     };
    850   }
    851   inline float4 max(const float4 &val1, const float4 &val2) {
    852     return {
    853       ::std::max(val1.x, val2.x),
    854       ::std::max(val1.y, val2.y),
    855       ::std::max(val1.z, val2.z),
    856       ::std::max(val1.w, val2.w)
    857     };
    858   }
    859   inline float4 clamp(const float4 &val, const float4 &min, const float4 &max) {
    860     return {
    861       _impl::clamp(val.x, min.x, max.x),
    862       _impl::clamp(val.y, min.y, max.y),
    863       _impl::clamp(val.z, min.z, max.z),
    864       _impl::clamp(val.w, min.w, max.w)
    865     };
    866   }
    867   inline float4 lerp(const float4 &val1, const float4 &val2, float amount) {
    868     return {
    869       _impl::lerp(val1.x, val2.x, amount),
    870       _impl::lerp(val1.y, val2.y, amount),
    871       _impl::lerp(val1.z, val2.z, amount),
    872       _impl::lerp(val1.w, val2.w, amount)
    873     };
    874   }
    875   // TODO: impl
    876   inline float4 transform(const float4 &pos, const float4x4 &mat);
    877   inline float4 transform4(const float3 &pos, const float4x4 &mat);
    878   inline float4 transform4(const float2 &pos, const float4x4 &mat);
    879   inline float4 transform(const float4 &val, const quaternion &rot) {
    880     // See comments in the float3 transform function.
    881     quaternion result = rot * quaternion(float3(val.x, val.y, val.z), 0.0f) * inverse(rot);
    882     return { result.x, result.y, result.z, val.w };
    883   }
    884   inline float4 transform4(const float3 &val, const quaternion &rot) {
    885     quaternion result = rot * quaternion(val, 0.0f) * inverse(rot);
    886     return { result.x, result.y, result.z, 1.0f };
    887   }
    888   inline float4 transform4(const float2 &val, const quaternion &rot) {
    889     quaternion result = rot * quaternion(float3(val.x, val.y, 0.0f), 0.0f) * inverse(rot);
    890     return { result.x, result.y, result.z, 1.0f };
    891   }
    892 
    893 } _WINDOWS_NUMERICS_END_NAMESPACE_
    894 
    895 
    896 // Define float3x2 functions
    897 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    898 
    899   // TODO: impl
    900   inline float3x2 make_float3x2_translation(const float2 &pos);
    901   inline float3x2 make_float3x2_translation(float xpos, float ypos);
    902   inline float3x2 make_float3x2_scale(float xscale, float yscale);
    903   inline float3x2 make_float3x2_scale(float xscale, float yscale, const float2 &center);
    904   inline float3x2 make_float3x2_scale(const float2 &xyscale);
    905   inline float3x2 make_float3x2_scale(const float2 &xyscale, const float2 &center);
    906   inline float3x2 make_float3x2_scale(float scale);
    907   inline float3x2 make_float3x2_scale(float scale, const float2 &center);
    908   inline float3x2 make_float3x2_skew(float xrad, float yrad);
    909   inline float3x2 make_float3x2_skew(float xrad, float yrad, const float2 &center);
    910   inline float3x2 make_float3x2_rotation(float rad);
    911   inline float3x2 make_float3x2_rotation(float rad, const float2 &center);
    912   inline bool is_identity(const float3x2 &val) {
    913     return val == float3x2::identity();
    914   }
    915   inline float determinant(const float3x2 &val) {
    916     // 2D transformation matrix has an implied 3rd column with (0, 0, 1)
    917     // det = m11 * m22 * m33 + m12 * m23 * m31 + m13 * m21 * m32
    918     //     - m31 * m22 * m13 - m21 * m12 * m33 - m11 * m32 * m23;
    919     return val.m11 * val.m22 - val.m21 * val.m12;
    920   }
    921   inline float2 translation(const float3x2 &val) {
    922     return { val.m31, val.m32 };
    923   }
    924   inline bool invert(const float3x2 &val, float3x2 *out);
    925   inline float3x2 lerp(const float3x2 &mat1, const float3x2 &mat2, float amount);
    926 
    927 } _WINDOWS_NUMERICS_END_NAMESPACE_
    928 
    929 
    930 // Define float4x4 functions
    931 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
    932 
    933   // TODO: impl
    934   inline float4x4 make_float4x4_billboard(const float3 &objpos, const float3 &camerapos, const float3 &cameraup, const float3 &camerafwd);
    935   inline float4x4 make_float4x4_constrained_billboard(const float3 &objpos, const float3 &camerapos, const float3 &rotateaxis, const float3 &camerafwd, const float3 &objfwd);
    936   inline float4x4 make_float4x4_translation(const float3 &pos) {
    937     return {
    938       1.0f,  0.0f,  0.0f,  0.0f,
    939       0.0f,  1.0f,  0.0f,  0.0f,
    940       0.0f,  0.0f,  1.0f,  0.0f,
    941       pos.x, pos.y, pos.z, 1.0f
    942     };
    943   }
    944   inline float4x4 make_float4x4_translation(float xpos, float ypos, float zpos);
    945   inline float4x4 make_float4x4_scale(float xscale, float yscale, float zscale);
    946   inline float4x4 make_float4x4_scale(float xscale, float yscale, float zscale, const float3 &center);
    947   inline float4x4 make_float4x4_scale(const float3 &xyzscale) {
    948     return {
    949       xyzscale.x, 0.0f,       0.0f,       0.0f,
    950       0.0f,       xyzscale.y, 0.0f,       0.0f,
    951       0.0f,       0.0f,       xyzscale.z, 0.0f,
    952       0.0f,       0.0f,       0.0f,       1.0f
    953     };
    954   }
    955   inline float4x4 make_float4x4_scale(const float3 &xyzscale, const float3 &center);
    956   inline float4x4 make_float4x4_scale(float scale);
    957   inline float4x4 make_float4x4_scale(float scale, const float3 &center);
    958   inline float4x4 make_float4x4_rotation_x(float rad);
    959   inline float4x4 make_float4x4_rotation_x(float rad, const float3 &center);
    960   inline float4x4 make_float4x4_rotation_y(float rad);
    961   inline float4x4 make_float4x4_rotation_y(float rad, const float3 &center);
    962   inline float4x4 make_float4x4_rotation_z(float rad);
    963   inline float4x4 make_float4x4_rotation_z(float rad, const float3 &center);
    964   inline float4x4 make_float4x4_from_axis_angle(const float3 &axis, float angle);
    965   inline float4x4 make_float4x4_perspective_field_of_view(float fov, float aspect, float nearplane, float farplane);
    966   inline float4x4 make_float4x4_perspective(float w, float h, float nearplane, float farplane);
    967   inline float4x4 make_float4x4_perspective_off_center(float left, float right, float bottom, float top, float nearplane, float farplane);
    968   inline float4x4 make_float4x4_orthographic(float w, float h, float znearplane, float zfarplane);
    969   inline float4x4 make_float4x4_orthographic_off_center(float left, float right, float bottom, float top, float znearplane, float zfarplane);
    970   inline float4x4 make_float4x4_look_at(const float3 &camerapos, const float3 &target, const float3 &cameraup);
    971   inline float4x4 make_float4x4_world(const float3 &pos, const float3 &fwd, const float3 &up);
    972   inline float4x4 make_float4x4_from_quaternion(const quaternion &quat) {
    973     // https://en.wikipedia.org/wiki/Rotation_matrix#Quaternion
    974     float xx = quat.x * quat.x;
    975     float yy = quat.y * quat.y;
    976     float zz = quat.z * quat.z;
    977     float xy = quat.x * quat.y;
    978     float xz = quat.x * quat.z;
    979     float xw = quat.x * quat.w;
    980     float yz = quat.y * quat.z;
    981     float yw = quat.y * quat.w;
    982     float zw = quat.z * quat.w;
    983     return {
    984       1.0f - 2.0f*yy - 2.0f*zz, 2.0f*xy + 2.0f*zw,        2.0f*xz - 2.0f*yw,        0.0f,
    985       2.0f*xy - 2.0f*zw,        1.0f - 2.0f*xx - 2.0f*zz, 2.0f*yz + 2.0f*xw,        0.0f,
    986       2.0f*xz + 2.0f*yw,        2.0f*yz - 2.0f*xw,        1.0f - 2.0f*xx - 2.0f*yy, 0.0f,
    987       0.0f,                     0.0f,                     0.0f,                     1.0f
    988     };
    989   }
    990   inline float4x4 make_float4x4_from_yaw_pitch_roll(float yaw, float pitch, float roll);
    991   inline float4x4 make_float4x4_shadow(const float3 &lightdir, const plane &plane);
    992   inline float4x4 make_float4x4_reflection(const plane &plane);
    993   inline bool is_identity(const float4x4 &val) {
    994     return val == float4x4::identity();
    995   }
    996   inline float determinant(const float4x4 &val) {
    997     const float det_33_44 = (val.m33 * val.m44 - val.m34 * val.m43);
    998     const float det_32_44 = (val.m32 * val.m44 - val.m34 * val.m42);
    999     const float det_32_43 = (val.m32 * val.m43 - val.m33 * val.m42);
   1000     const float det_31_44 = (val.m31 * val.m44 - val.m34 * val.m41);
   1001     const float det_31_43 = (val.m31 * val.m43 - val.m33 * val.m41);
   1002     const float det_31_42 = (val.m31 * val.m42 - val.m32 * val.m41);
   1003     return val.m11 * (val.m22 * det_33_44 - val.m23 * det_32_44 + val.m24 * det_32_43)
   1004       - val.m12 * (val.m21 * det_33_44 - val.m23 * det_31_44 + val.m24 * det_31_43)
   1005       + val.m13 * (val.m21 * det_32_44 - val.m22 * det_31_44 + val.m24 * det_31_42)
   1006       - val.m14 * (val.m21 * det_32_43 - val.m22 * det_31_43 + val.m23 * det_31_42);
   1007   }
   1008   inline float3 translation(const float4x4 &val) {
   1009     return { val.m41, val.m42, val.m43 };
   1010   }
   1011   inline bool invert(const float4x4 &mat, float4x4 *out);
   1012   inline bool decompose(const float4x4 &mat, float3 *out_scale, quaternion *out_rot, float3 *out_translate);
   1013   inline float4x4 transform(const float4x4 &val, const quaternion &rot) {
   1014     return val * make_float4x4_from_quaternion(rot);
   1015   }
   1016   inline float4x4 transpose(const float4x4 &val) {
   1017     return  {
   1018       val.m11, val.m21, val.m31, val.m41,
   1019       val.m12, val.m22, val.m32, val.m42,
   1020       val.m13, val.m23, val.m33, val.m43,
   1021       val.m14, val.m24, val.m34, val.m44
   1022     };
   1023   }
   1024   inline float4x4 lerp(const float4x4 &val1, const float4x4 &val2, float amount);
   1025 
   1026 } _WINDOWS_NUMERICS_END_NAMESPACE_
   1027 
   1028 
   1029 // Define plane functions
   1030 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
   1031 
   1032   // TODO: impl
   1033   inline plane make_plane_from_vertices(const float3 &pt1, const float3 &pt2, const float3 &pt3);
   1034   inline plane normalize(const plane &val) {
   1035     const float invlen = 1.0f / length(val.normal);
   1036     return { val.normal * invlen, val.d * invlen };
   1037   }
   1038   inline plane transform(const plane &plane, const float4x4 &mat);
   1039   inline plane transform(const plane &plane, const quaternion &rot) {
   1040     quaternion result = rot * quaternion(plane.normal, 0.0f) * inverse(rot);
   1041     return { float3(result.x, result.y, result.z), plane.d };
   1042   }
   1043   inline float dot(const plane &plane, const float4 &val);
   1044   inline float dot_coordinate(const plane &plane, const float3 &val);
   1045   inline float dot_normal(const plane &plane, const float3 &val);
   1046 
   1047 } _WINDOWS_NUMERICS_END_NAMESPACE_
   1048 
   1049 
   1050 // Define quaternion functions
   1051 _WINDOWS_NUMERICS_BEGIN_NAMESPACE_ {
   1052 
   1053   inline quaternion make_quaternion_from_axis_angle(const float3 &axis, float angle) {
   1054     return quaternion(axis * ::std::sin(angle * 0.5f), ::std::cos(angle * 0.5f));
   1055   }
   1056   inline quaternion make_quaternion_from_yaw_pitch_roll(float yaw, float pitch, float roll) {
   1057     quaternion yq = make_quaternion_from_axis_angle(float3(0.0f, 1.0f, 0.0f), yaw);
   1058     quaternion pq = make_quaternion_from_axis_angle(float3(1.0f, 0.0f, 0.0f), pitch);
   1059     quaternion rq = make_quaternion_from_axis_angle(float3(0.0f, 0.0f, 1.0f), roll);
   1060     return concatenate(concatenate(rq, pq), yq);
   1061   }
   1062   inline quaternion make_quaternion_from_rotation_matrix(const float4x4 &mat) {
   1063     // https://en.wikipedia.org/wiki/Rotation_matrix#Quaternion
   1064     float t = mat.m11 + mat.m22 + mat.m33;
   1065     if (t > 0) {
   1066       float r = ::std::sqrt(1.0f + t);
   1067       float s = 1.0f / (2.0f * r);
   1068       return { (mat.m23 - mat.m32) * s, (mat.m31 - mat.m13) * s,
   1069                (mat.m12 - mat.m21) * s, r * 0.5f };
   1070     } else if (mat.m11 >= mat.m22 && mat.m11 >= mat.m33) {
   1071       float r = ::std::sqrt(1.0f + mat.m11 - mat.m22 - mat.m33);
   1072       float s = 1.0f / (2.0f * r);
   1073       return { r * 0.5f, (mat.m21 + mat.m12) * s,
   1074                (mat.m13 + mat.m31) * s, (mat.m23 - mat.m32) * s };
   1075     } else if (mat.m22 >= mat.m11 && mat.m22 >= mat.m33) {
   1076       float r = ::std::sqrt(1.0f + mat.m22 - mat.m11 - mat.m33);
   1077       float s = 1.0f / (2.0f * r);
   1078       return { (mat.m21 + mat.m12) * s, r * 0.5f,
   1079                (mat.m32 + mat.m23) * s, (mat.m31 - mat.m13) * s };
   1080     } else {
   1081       float r = ::std::sqrt(1.0f + mat.m33 - mat.m11 - mat.m22);
   1082       float s = 1.0f / (2.0f * r);
   1083       return { (mat.m13 + mat.m31) * s, (mat.m32 + mat.m23) * s,
   1084                r * 0.5f, (mat.m12 - mat.m21) * s };
   1085     }
   1086   }
   1087   inline bool is_identity(const quaternion &val) {
   1088     return val == quaternion::identity();
   1089   }
   1090   inline float length(const quaternion &val) {
   1091     return ::std::sqrt(length_squared(val));
   1092   }
   1093   inline float length_squared(const quaternion &val) {
   1094     return dot(val, val);
   1095   }
   1096   inline float dot(const quaternion &val1, const quaternion &val2) {
   1097     return val1.x * val2.x + val1.y * val2.y + val1.z * val2.z + val1.w * val2.w;
   1098   }
   1099   inline quaternion normalize(const quaternion &val) {
   1100     return operator*(val, 1.0f / length(val));
   1101   }
   1102   inline quaternion conjugate(const quaternion &val) {
   1103     return { -val.x, -val.y, -val.z, val.w};
   1104   }
   1105   inline quaternion inverse(const quaternion &val) {
   1106     return operator*(conjugate(val), 1.0f / length_squared(val));
   1107   }
   1108   inline quaternion slerp(const quaternion &val1, const quaternion &val2, float amount) {
   1109     // https://en.wikipedia.org/wiki/Slerp#Geometric_Slerp
   1110     float cosangle = dot(val1, val2);
   1111     quaternion end = val2;
   1112     if (cosangle < 0.0f) {
   1113       end = -val2;
   1114       cosangle = -cosangle;
   1115     }
   1116     float fact1, fact2;
   1117     const float epsilon = 1.0e-6f;
   1118     if (cosangle > 1.0f - epsilon) {
   1119       // Very small rotation range, or non-normalized input quaternions.
   1120       fact1 = (1.0f - amount);
   1121       fact2 = amount;
   1122     } else {
   1123       float angle = ::std::acos(cosangle);
   1124       float sinangle = ::std::sin(angle);
   1125       fact1 = ::std::sin((1.0f - amount) * angle) / sinangle;
   1126       fact2 = ::std::sin(amount * angle) / sinangle;
   1127     }
   1128     return val1 * fact1 + end * fact2;
   1129   }
   1130   inline quaternion lerp(const quaternion &val1, const quaternion &val2, float amount) {
   1131     quaternion end = val2;
   1132     if (dot(val1, val2) < 0.0f)
   1133       end = -val2;
   1134     return normalize(quaternion(
   1135       _impl::lerp(val1.x, end.x, amount), _impl::lerp(val1.y, end.y, amount),
   1136       _impl::lerp(val1.z, end.z, amount), _impl::lerp(val1.w, end.w, amount)
   1137     ));
   1138   }
   1139   inline quaternion concatenate(const quaternion &val1, const quaternion &val2) {
   1140     return val2 * val1;
   1141   }
   1142 
   1143 } _WINDOWS_NUMERICS_END_NAMESPACE_
   1144 
   1145 
   1146 /**
   1147  * FIXME: Implement interop functions with DirectXMath.
   1148  * This is where we are supposed to define the functions to convert between
   1149  * Windows::Foundation::Numerics types and XMVECTOR / XMMATRIX. But our
   1150  * directxmath.h does not contain the definitions for these types...
   1151  */
   1152 #if 0
   1153 // === DirectXMath Interop ===
   1154 namespace DirectX {
   1155 
   1156   // TODO: impl
   1157   XMVECTOR XMLoadFloat2(const _WINDOWS_NUMERICS_NAMESPACE_::float2 *src);
   1158   XMVECTOR XMLoadFloat3(const _WINDOWS_NUMERICS_NAMESPACE_::float3 *src);
   1159   XMVECTOR XMLoadFloat4(const _WINDOWS_NUMERICS_NAMESPACE_::float4 *src);
   1160   XMMATRIX XMLoadFloat3x2(const _WINDOWS_NUMERICS_NAMESPACE_::float3x2 *src);
   1161   XMMATRIX XMLoadFloat4x4(const _WINDOWS_NUMERICS_NAMESPACE_::float4x4 *src);
   1162   XMVECTOR XMLoadPlane(const _WINDOWS_NUMERICS_NAMESPACE_::plane *src);
   1163   XMVECTOR XMLoadQuaternion(const _WINDOWS_NUMERICS_NAMESPACE_::quaternion *src);
   1164   void XMStoreFloat2(_WINDOWS_NUMERICS_NAMESPACE_::float2 *out, FXMVECTOR in);
   1165   void XMStoreFloat3(_WINDOWS_NUMERICS_NAMESPACE_::float3 *out, FXMVECTOR in);
   1166   void XMStoreFloat4(_WINDOWS_NUMERICS_NAMESPACE_::float4 *out, FXMVECTOR in);
   1167   void XMStoreFloat3x2(_WINDOWS_NUMERICS_NAMESPACE_::float3x2 *out, FXMMATRIX in);
   1168   void XMStoreFloat4x4(_WINDOWS_NUMERICS_NAMESPACE_::float4x4 *out, FXMMATRIX in);
   1169   void XMStorePlane(_WINDOWS_NUMERICS_NAMESPACE_::plane *out, FXMVECTOR in);
   1170   void XMStoreQuaternion(_WINDOWS_NUMERICS_NAMESPACE_::quaternion *out, FXMVECTOR in);
   1171 
   1172 } /* namespace DirectX */
   1173 #endif
   1174 
   1175 
   1176 #undef _WINDOWS_NUMERICS_IMPL_ASSIGN_OP
   1177 
   1178 #ifdef _WINDOWS_NUMERICS_IMPL_PUSHED_MIN_
   1179 #  undef _WINDOWS_NUMERICS_IMPL_PUSHED_MIN_
   1180 #  pragma pop_macro("min")
   1181 #endif
   1182 
   1183 #ifdef _WINDOWS_NUMERICS_IMPL_PUSHED_MAX_
   1184 #  undef _WINDOWS_NUMERICS_IMPL_PUSHED_MAX_
   1185 #  pragma pop_macro("max")
   1186 #endif