/src/shared/vmath.cpp

  1#include "vmath.h"
  2
  3///////////////////////////////////////////////////////////////////////////
  4// SIMD
  5///////////////////////////////////////////////////////////////////////////
  6
  7const __m128 SIMD_EPSILON = _mm_set_ps1(EPSILON);
  8const __m128 SIMD_ZERO = _mm_setzero_ps();
  9const __m128 SIMD_ONE = _mm_set_ps1(1.0f);
 10const __m128 SIMD_1110 = _mm_set_ps(0.0f, 1.0f, 1.0, 1.0f);
 11const __m128 SIMD_1000 = _mm_set_ps(0.0f, 0.0f, 0.0, 1.0f);
 12const __m128 SIMD_0100 = _mm_set_ps(0.0f, 0.0f, 1.0, 0.0f);
 13const __m128 SIMD_0010 = _mm_set_ps(0.0f, 1.0f, 0.0, 0.0f);
 14const __m128 SIMD_0001 = _mm_set_ps(1.0f, 0.0f, 0.0, 0.0f);
 15
 16__declspec(align(16)) const int SIMD_INT_MASK_AND_XYZ0[4] = {0xFFFFFFFF, 0xFFFFFFFF, 0xFFFFFFFF, 0};
 17__declspec(align(16)) const int SIMD_INT_MASK_AND_ABS[4] =  {0x7FFFFFFF, 0x7FFFFFFF, 0x7FFFFFFF, 0x7FFFFFFF};
 18__declspec(align(16)) const int SIMD_INT_MASK_XOR_SIGN_XYZ0[4] = {0x80000000, 0x80000000, 0x80000000, 0x00000000};
 19__declspec(align(16)) const int SIMD_INT_MASK_XOR_SIGN_0Y00[4] = {0, 0x80000000, 0, 0};
 20__declspec(align(16)) const int SIMD_INT_MASK_XOR_SIGN_000W[4] = {0, 0, 0, 0x80000000};
 21
 22const __m128 SIMD_MASK_AND_XYZ0 = *reinterpret_cast<const __m128 *>(SIMD_INT_MASK_AND_XYZ0);
 23const __m128 SIMD_MASK_AND_ABS = *reinterpret_cast<const __m128 *>(SIMD_INT_MASK_AND_ABS);
 24const __m128 SIMD_MASK_XOR_SIGN_XYZ0 = *reinterpret_cast<const __m128 *>(SIMD_INT_MASK_XOR_SIGN_XYZ0);
 25const __m128 SIMD_MASK_XOR_SIGN_0Y00 = *reinterpret_cast<const __m128 *>(SIMD_INT_MASK_XOR_SIGN_0Y00);
 26const __m128 SIMD_MASK_XOR_SIGN_000W = *reinterpret_cast<const __m128 *>(SIMD_INT_MASK_XOR_SIGN_000W);
 27
 28const __m128 SIMD_MATRIX_INVERSE_SIGN_A = _mm_set_ps(1.0f, -1.0f, 1.0f, -1.0f);
 29const __m128 SIMD_MATRIX_INVERSE_SIGN_B = _mm_set_ps(-1.0f, 1.0f, -1.0f, 1.0f);
 30
 31///////////////////////////////////////////////////////////////////////////
 32// 4D vector
 33///////////////////////////////////////////////////////////////////////////
 34
 35const vec4 vec4::ZERO = vec4(SIMD_ZERO);
 36const vec4 vec4::ONE = vec4(SIMD_ONE);
 37const vec4 vec4::UNIT_X = vec4(SIMD_1000);
 38const vec4 vec4::UNIT_Y = vec4(SIMD_0100);
 39const vec4 vec4::UNIT_Z = vec4(SIMD_0010);
 40const vec4 vec4::UNIT_W = vec4(SIMD_0001);
 41
 42///////////////////////////////////////////////////////////////////////////
 43// Quaternion
 44///////////////////////////////////////////////////////////////////////////
 45
 46const vec4 QUAT_IDENTITY = vec4(SIMD_0001);
 47	
 48// multiply quaternions
 49vec4 quat_mult(const vec4 &lhs, const vec4 &rhs)
 50{
 51	__m128 a1 = _mm_shuffle_ps(lhs.xyzw, lhs.xyzw, _MM_SHUFFLE(1, 3, 2, 0));
 52	__m128 a2 = _mm_shuffle_ps(lhs.xyzw, lhs.xyzw, _MM_SHUFFLE(0, 2, 3, 1));
 53	__m128 a3 = _mm_shuffle_ps(lhs.xyzw, lhs.xyzw, _MM_SHUFFLE(3, 1, 0, 2));
 54	__m128 a4 = _mm_shuffle_ps(lhs.xyzw, lhs.xyzw, _MM_SHUFFLE(2, 0, 1, 3));
 55	__m128 b1 = _mm_shuffle_ps(rhs.xyzw, rhs.xyzw, _MM_SHUFFLE(1, 2, 0, 3));
 56	__m128 b2 = _mm_shuffle_ps(rhs.xyzw, rhs.xyzw, _MM_SHUFFLE(0, 3, 1, 2));
 57	__m128 b3 = _mm_shuffle_ps(rhs.xyzw, rhs.xyzw, _MM_SHUFFLE(3, 0, 2, 1));
 58	__m128 b4 = _mm_shuffle_ps(rhs.xyzw, rhs.xyzw, _MM_SHUFFLE(2, 1, 3, 0));
 59	return vec4(_mm_xor_ps(
 60		_mm_add_ps(
 61			_mm_add_ps(_mm_mul_ps(a2, b2), _mm_mul_ps(a1, b1)),
 62			_mm_sub_ps(_mm_mul_ps(a4, b4), _mm_mul_ps(a3, b3))),
 63		SIMD_MASK_XOR_SIGN_000W));
 64}
 65
 66// spherical linear interpolation
 67vec4 quat_slerp(const vec4 &v0, const vec4 &v1, float t)
 68{
 69	__m128 cosa = _mm_dot_ps(v0.xyzw, v1.xyzw);
 70	__m128 abscosa = _mm_and_ps(cosa, SIMD_MASK_AND_ABS);
 71	if (_mm_movemask_ps(_mm_cmpgt_ps(_mm_sub_ps(SIMD_ONE, abscosa), SIMD_EPSILON)) == 0xF)
 72	{
 73		float a;
 74		_mm_store_ss(&a, abscosa);
 75		a = acos(a);
 76		float sina = 1.0f / sin(a);
 77		float s0 = sin((1.0f - t) * a) * sina;
 78		float s1 = sin(t * a) * sina;
 79		__m128 scale0 = _mm_load_ps1(&s0);
 80		__m128 scale1 = _mm_load_ps1(&s1);
 81		if (_mm_movemask_ps(_mm_cmplt_ps(cosa, SIMD_ZERO)) == 0xF)
 82		{
 83			scale1 = _mm_sub_ps(SIMD_ZERO, scale1);
 84		}
 85		return vec4(_mm_add_ps(_mm_mul_ps(scale0, v0.xyzw), _mm_mul_ps(scale1, v1.xyzw)));
 86	}
 87	return v0;
 88}
 89
 90///////////////////////////////////////////////////////////////////////////
 91// 4x4 matrix
 92///////////////////////////////////////////////////////////////////////////
 93
 94const mat4 mat4::ZERO = mat4(SIMD_ZERO, SIMD_ZERO, SIMD_ZERO, SIMD_ZERO);
 95const mat4 mat4::IDENTITY = mat4(SIMD_1000, SIMD_0100, SIMD_0010, SIMD_0001);
 96
 97// mul
 98mat4 &operator*=(mat4 &lhs, const mat4 &rhs)
 99{
100	lhs = lhs * rhs;
101	return lhs;
102}
103
104const mat4 operator*(const mat4 &lhs, const mat4 &rhs)
105{
106	mat4 m;
107
108	m.row0 = _mm_add_ps(
109		_mm_add_ps(
110			_mm_mul_ps(rhs.row0, _mm_shuffle_ps(lhs.row0, lhs.row0, _MM_SHUFFLE(0, 0, 0, 0))),
111			_mm_mul_ps(rhs.row1, _mm_shuffle_ps(lhs.row0, lhs.row0, _MM_SHUFFLE(1, 1, 1, 1)))),
112		_mm_add_ps(
113			_mm_mul_ps(rhs.row2, _mm_shuffle_ps(lhs.row0, lhs.row0, _MM_SHUFFLE(2, 2, 2, 2))),
114			_mm_mul_ps(rhs.row3, _mm_shuffle_ps(lhs.row0, lhs.row0, _MM_SHUFFLE(3, 3, 3, 3)))));
115
116	m.row1 = _mm_add_ps(
117		_mm_add_ps(
118			_mm_mul_ps(rhs.row0, _mm_shuffle_ps(lhs.row1, lhs.row1, _MM_SHUFFLE(0, 0, 0, 0))),
119			_mm_mul_ps(rhs.row1, _mm_shuffle_ps(lhs.row1, lhs.row1, _MM_SHUFFLE(1, 1, 1, 1)))),
120		_mm_add_ps(
121			_mm_mul_ps(rhs.row2, _mm_shuffle_ps(lhs.row1, lhs.row1, _MM_SHUFFLE(2, 2, 2, 2))),
122			_mm_mul_ps(rhs.row3, _mm_shuffle_ps(lhs.row1, lhs.row1, _MM_SHUFFLE(3, 3, 3, 3)))));
123
124	m.row2 = _mm_add_ps(
125		_mm_add_ps(
126			_mm_mul_ps(rhs.row0, _mm_shuffle_ps(lhs.row2, lhs.row2, _MM_SHUFFLE(0, 0, 0, 0))),
127			_mm_mul_ps(rhs.row1, _mm_shuffle_ps(lhs.row2, lhs.row2, _MM_SHUFFLE(1, 1, 1, 1)))),
128		_mm_add_ps(
129			_mm_mul_ps(rhs.row2, _mm_shuffle_ps(lhs.row2, lhs.row2, _MM_SHUFFLE(2, 2, 2, 2))),
130			_mm_mul_ps(rhs.row3, _mm_shuffle_ps(lhs.row2, lhs.row2, _MM_SHUFFLE(3, 3, 3, 3)))));
131
132	m.row3 = _mm_add_ps(
133		_mm_add_ps(
134			_mm_mul_ps(rhs.row0, _mm_shuffle_ps(lhs.row3, lhs.row3, _MM_SHUFFLE(0, 0, 0, 0))),
135			_mm_mul_ps(rhs.row1, _mm_shuffle_ps(lhs.row3, lhs.row3, _MM_SHUFFLE(1, 1, 1, 1)))),
136		_mm_add_ps(
137			_mm_mul_ps(rhs.row2, _mm_shuffle_ps(lhs.row3, lhs.row3, _MM_SHUFFLE(2, 2, 2, 2))),
138			_mm_mul_ps(rhs.row3, _mm_shuffle_ps(lhs.row3, lhs.row3, _MM_SHUFFLE(3, 3, 3, 3)))));
139
140	return m;
141}
142
143// inverse matrix4x4
144mat4 matrix_inverse(const mat4 &v)
145{
146	// http://www.devmaster.net/forums/showthread.php?t=14569
147
148	__m128 Swp0a, Swp0b, Swp00, Swp01, Swp02, Swp03, Temp;
149
150	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(3, 3, 3, 3));
151	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(2, 2, 2, 2));
152	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(2, 2, 2, 2));
153	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
154	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
155	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(3, 3, 3, 3));
156	__m128 Fac0 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
157
158	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(3, 3, 3, 3));
159	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(1, 1, 1, 1));
160	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(1, 1, 1, 1));
161	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
162	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
163	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(3, 3, 3, 3));
164	__m128 Fac1 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
165
166	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(2, 2, 2, 2));
167	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(1, 1, 1, 1));
168	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(1, 1, 1, 1));
169	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
170	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
171	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(2, 2, 2, 2));
172	__m128 Fac2 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
173
174	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(3, 3, 3, 3));
175	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(0, 0, 0, 0));
176	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(0, 0, 0, 0));
177	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
178	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
179	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(3, 3, 3, 3));
180	__m128 Fac3 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
181
182	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(2, 2, 2, 2));
183	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(0, 0, 0, 0));
184	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(0, 0, 0, 0));
185	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
186	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
187	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(2, 2, 2, 2));
188	__m128 Fac4 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
189
190	Swp0a = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(1, 1, 1, 1));
191	Swp0b = _mm_shuffle_ps(v.row3, v.row2, _MM_SHUFFLE(0, 0, 0, 0));
192	Swp00 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(0, 0, 0, 0));
193	Swp01 = _mm_shuffle_ps(Swp0a, Swp0a, _MM_SHUFFLE(2, 0, 0, 0));
194	Swp02 = _mm_shuffle_ps(Swp0b, Swp0b, _MM_SHUFFLE(2, 0, 0, 0));
195	Swp03 = _mm_shuffle_ps(v.row2, v.row1, _MM_SHUFFLE(1, 1, 1, 1));
196	__m128 Fac5 = _mm_sub_ps(_mm_mul_ps(Swp00, Swp01), _mm_mul_ps(Swp02, Swp03));
197
198	Temp = _mm_shuffle_ps(v.row1, v.row0, _MM_SHUFFLE(0, 0, 0, 0));
199	__m128 Vec0 = _mm_shuffle_ps(Temp, Temp, _MM_SHUFFLE(2, 2, 2, 0));
200	Temp = _mm_shuffle_ps(v.row1, v.row0, _MM_SHUFFLE(1, 1, 1, 1));
201	__m128 Vec1 = _mm_shuffle_ps(Temp, Temp, _MM_SHUFFLE(2, 2, 2, 0));
202	Temp = _mm_shuffle_ps(v.row1, v.row0, _MM_SHUFFLE(2, 2, 2, 2));
203	__m128 Vec2 = _mm_shuffle_ps(Temp, Temp, _MM_SHUFFLE(2, 2, 2, 0));
204	Temp = _mm_shuffle_ps(v.row1, v.row0, _MM_SHUFFLE(3, 3, 3, 3));
205	__m128 Vec3 = _mm_shuffle_ps(Temp, Temp, _MM_SHUFFLE(2, 2, 2, 0));
206
207	__m128 Inv0 = _mm_mul_ps(
208		_mm_add_ps(
209			_mm_mul_ps(Vec3, Fac2),
210			_mm_sub_ps(
211				_mm_mul_ps(Vec1, Fac0),
212				_mm_mul_ps(Vec2, Fac1))),
213		SIMD_MATRIX_INVERSE_SIGN_B);
214
215	__m128 Inv1 = _mm_mul_ps(
216		_mm_add_ps(
217			_mm_mul_ps(Vec3, Fac4),
218			_mm_sub_ps(
219				_mm_mul_ps(Vec0, Fac0),
220				_mm_mul_ps(Vec2, Fac3))),
221		SIMD_MATRIX_INVERSE_SIGN_A);
222
223	__m128 Inv2 = _mm_mul_ps(
224		_mm_add_ps(
225			_mm_mul_ps(Vec3, Fac5),
226			_mm_sub_ps(
227				_mm_mul_ps(Vec0, Fac1),
228				_mm_mul_ps(Vec1, Fac3))),
229		SIMD_MATRIX_INVERSE_SIGN_B);
230
231	__m128 Inv3 = _mm_mul_ps(
232		_mm_add_ps(
233			_mm_mul_ps(Vec2, Fac5),
234			_mm_sub_ps(
235				_mm_mul_ps(Vec0, Fac2),
236				_mm_mul_ps(Vec1, Fac4))),
237		SIMD_MATRIX_INVERSE_SIGN_A);
238
239	__m128 Det = _mm_dot_ps(
240		_mm_shuffle_ps(
241			_mm_shuffle_ps(Inv0, Inv1, _MM_SHUFFLE(0, 0, 0, 0)),
242			_mm_shuffle_ps(Inv2, Inv3, _MM_SHUFFLE(0, 0, 0, 0)),
243			_MM_SHUFFLE(2, 0, 2, 0)),
244		v.row0);
245	Det = _mm_div_ps(SIMD_ONE, Det);
246	//Det = _mm_rcp_ps(Det);
247
248	//	Inverse /= Determinant;
249	mat4 m;
250	m.row0 = _mm_mul_ps(Inv0, Det);
251	m.row1 = _mm_mul_ps(Inv1, Det);
252	m.row2 = _mm_mul_ps(Inv2, Det);
253	m.row3 = _mm_mul_ps(Inv3, Det);
254	return m;
255}
256
257// scale matrix
258mat4 matrix_scaling(const vec4 &v)
259{
260	mat4 m;
261	m.row0 = _mm_mul_ps(v.xyzw, SIMD_1000);
262	m.row1 = _mm_mul_ps(v.xyzw, SIMD_0100);
263	m.row2 = _mm_mul_ps(v.xyzw, SIMD_0010);
264	m.row3 = SIMD_0001;
265	return m;
266}
267
268// translation matrix
269mat4 matrix_translation(const vec4 &v)
270{
271	mat4 m;
272	m.row0 = SIMD_1000;
273	m.row1 = SIMD_0100;
274	m.row2 = SIMD_0010;
275	m.row3 = _mm_merge_ps(v.xyzw, SIMD_1000);
276	return m;
277}
278
279// rotation matrix from quaternion
280mat4 matrix_rotation(const vec4 &v)
281{
282	__m128 v2 = _mm_and_ps(_mm_add_ps(v.xyzw, v.xyzw), SIMD_MASK_AND_XYZ0);
283		
284	__m128 wv = _mm_xor_ps(_mm_mul_ps(_mm_shuffle_ps(v.xyzw, v.xyzw, _MM_SHUFFLE(3, 3, 3, 3)), v2), SIMD_MASK_XOR_SIGN_0Y00);
285	__m128 xy_yz_xz = _mm_mul_ps(_mm_shuffle_ps(v.xyzw, v.xyzw, _MM_SHUFFLE(0, 0, 2, 1)), v2);
286	__m128 yz_xz_xy = _mm_shuffle_ps(xy_yz_xz, xy_yz_xz, _MM_SHUFFLE(3, 0, 2, 1));
287
288	__m128 a = _mm_add_ps(yz_xz_xy, wv);
289	__m128 b = _mm_sub_ps(yz_xz_xy, wv);
290
291	__m128 vv = _mm_mul_ps(v.xyzw, v2);
292	__m128 yy_xx_xx = _mm_shuffle_ps(vv, vv, _MM_SHUFFLE(3, 0, 0, 1));
293	__m128 zz_zz_yy = _mm_shuffle_ps(vv, vv, _MM_SHUFFLE(3, 1, 2, 2));
294
295	__m128 c = _mm_sub_ps(SIMD_1110, _mm_add_ps(yy_xx_xx, zz_zz_yy));
296
297	mat4 m;
298	m.row0 = _mm_shuffle_ps(c, a, _MM_SHUFFLE(1, 2, 3, 0));
299	m.row0 = _mm_shuffle_ps(m.row0, m.row0, _MM_SHUFFLE(1, 3, 2, 0));
300	m.row1 = _mm_shuffle_ps(_mm_shuffle_ps(b, c, _MM_SHUFFLE(1, 1, 2, 2)), a, _MM_SHUFFLE(3, 0, 2, 0));
301	m.row2 = _mm_shuffle_ps(b, c, _MM_SHUFFLE(3, 2, 0, 1));
302	m.row3 = SIMD_0001;
303	return m;
304}
305
306mat4 matrix_lookat(const vec4 &eye, const vec4 &at, const vec4 &up)
307{
308	vec4 zaxis = normalize(at - eye);
309	vec4 xaxis = normalize(cross(up, zaxis));
310	vec4 yaxis = cross(zaxis, xaxis);
311
312	__m128 tmp = _mm_hadd_ps(
313		_mm_hadd_ps(_mm_mul_ps(xaxis.xyzw, eye.xyzw), _mm_mul_ps(yaxis.xyzw, eye.xyzw)),
314		_mm_hadd_ps(_mm_mul_ps(zaxis.xyzw, eye.xyzw), _mm_setzero_ps()));
315	tmp = _mm_merge_ps(_mm_sub_ps(_mm_setzero_ps(), tmp), SIMD_1000);
316
317	mat4 m;
318	m.row0 = xaxis.xyzw;
319	m.row1 = yaxis.xyzw;
320	m.row2 = zaxis.xyzw;
321	m.row3 = _mm_setzero_ps();
322	m.transpose();
323	m.row3 = tmp;
324		
325	return m;
326}
327
328mat4 matrix_perspective(float fov, float aspect, float znear, float zfar)
329{
330	float y = 1.0f / tan(fov / 2.0f);
331	float x = y / aspect;
332	float a = zfar / (zfar - znear);
333	float b = -znear * zfar / (zfar - znear);
334	return mat4(
335		x, 0.0f, 0.0f, 0.0f,
336		0.0f, y, 0.0f, 0.0f,
337		0.0f, 0.0f, a, 1.0f,
338		0.0f, 0.0f, b, 0.0f);
339}
340
341mat4 matrix_ortho(float w, float h, float znear, float zfar)
342{
343	float x = w / 2.0f;
344	float y = h / 2.0f;
345	float a = 1.0f / (zfar - znear);
346	float b = -znear / (zfar - znear);
347	return mat4(
348		x, 0.0f, 0.0f, 0.0f,
349		0.0f, y, 0.0f, 0.0f,
350		0.0f, 0.0f, a, 0.0f,
351		0.0f, 0.0f, b, 1.0f);
352}
353
354///////////////////////////////////////////////////////////////////////////
355// Axis aligned bounding box
356///////////////////////////////////////////////////////////////////////////
357
358const bbox bbox::ZERO = bbox(vec4::ZERO, vec4::ZERO);
359const bbox bbox::DEGENERATE = bbox(vec4(INFINITY), vec4(-INFINITY));