diff options
Diffstat (limited to 'src')
39 files changed, 3164 insertions, 0 deletions
diff --git a/src/core/color.c b/src/core/color.c new file mode 100644 index 0000000..45ca555 --- /dev/null +++ b/src/core/color.c @@ -0,0 +1,20 @@ +#include "vert.h" +#include "device.h" + +Color color32_tocolor(Color32* c) { + return ssr_color(0xff * c->r, 0xff * c->g, 0xff * c->b, 0xff * c->a); +} + +void color_tocolor32(Color c, Color32* out) { + out->r = COLOR_R(c) / (float)0xff; + out->g = COLOR_G(c) / (float)0xff; + out->b = COLOR_B(c) / (float)0xff; + out->a = COLOR_A(c) / (float)0xff; +} + +void color32_saturate(Color32* c) { + c->r = clamp(c->r, 0, 1); + c->g = clamp(c->g, 0, 1); + c->b = clamp(c->b, 0, 1); + c->a = clamp(c->a, 0, 1); +} diff --git a/src/core/color.h b/src/core/color.h new file mode 100644 index 0000000..1c9f45c --- /dev/null +++ b/src/core/color.h @@ -0,0 +1,7 @@ +#ifndef _SOFTSHADEROOM_COLOR_H_ +#define _SOFTSHADEROOM_COLOR_H_ + +#define COLOR_WHITE 0xffffffff +#define COLOR_BLACK 0xff000000 + +#endif diff --git a/src/core/device.c b/src/core/device.c new file mode 100644 index 0000000..f5ed4eb --- /dev/null +++ b/src/core/device.c @@ -0,0 +1,339 @@ +#include "device.h" +#include "rasterizer.h" +#include "../util/assert.h" + +ssr_Config config; + +typedef enum ssr_VertexAttrMask { + VERTEXATTR_POS = 1, + VERTEXATTR_COLOR = 1 << 1, + VERTEXATTR_UV = 1 << 2, + VERTEXATTR_NORMAL = 1 << 3, +} ssr_VertexAttrMask; + +typedef struct VertexAttr { + int start; + int offset; +} VertexAttr; + +/* +** 状态结构 +*/ +struct { + /* MVP矩阵栈 */ + Mat4 matrices[3][MATRIXDEPTH]; + uint matrixtop[3]; + ssr_MatrixMode matrixmode; + + /* 屏幕缓冲区 */ + Color* target; + Color* framebuffer; + uint buffersize; /*framebuffer size in bytess*/ + + /* zbuffer 深度会被映射到0~1的非线性区间*/ + uint* zbuffer; + + /* 顶点数据 */ + Vert** verts; uint nverts; + uint* indices; uint nprims; + + /* shader program */ + Program* program; + UniformCollection uniforms; + + uint enable; + +} state; + +/* +** MVP矩阵栈 +*/ +#define MATRIXTOP state.matrixtop[state.matrixmode] +#define MATRIX state.matrices[state.matrixmode][MATRIXTOP] +#define MATRIXSTACK state.matrices[state.matrixmode] + +#define GETMATRIX(MODE) state.matrices[MODE][state.matrixtop[MODE]] + +#define BUFFER (state.framebuffer) + +#define contains(x, y, l, r, t, b) (x >= l && x <= r && y <= b && y >= t) + +void ssr_init(ssr_Config* conf) { + config = *conf; + ssrM_zero(state.matrixtop, sizeof(state.matrixtop)); + state.target = conf->target; + if (config.dbuffer) { + state.framebuffer = ssrM_newvector(Color, config.width * config.height); + } + else { + state.framebuffer = conf->target; + } + state.buffersize = sizeof(Color) * config.width * config.height; + state.zbuffer = ssrM_newvector(uint, config.width * config.height); + memset(state.zbuffer, 0xff, sizeof(uint)*config.width*config.height); +} + +int ssr_getframebufferw() { + return config.width; +} + +int ssr_getframebufferh() { + return config.height; +} + +void ssr_matrixmode(ssr_MatrixMode mode) { + state.matrixmode = mode; +} + +void ssr_loadidentity() { + mat4_setidentity(&MATRIX); +} + +/* push进去之后会将最新的拷贝一份放在栈顶 */ +void ssr_pushmatrix() { + ssr_assert(MATRIXTOP < MATRIXDEPTH - 1); + ++MATRIXTOP; + MATRIX = MATRIXSTACK[MATRIXTOP - 1]; +} + +void ssr_popmatrix() { + MATRIXTOP = clamp(--MATRIXTOP, 0, MATRIXDEPTH - 1); +} + +/*用来做world->view的变换矩阵*/ +/*http://warmcat.org/chai/blog/?p=559*/ +void ssr_lookat(Vec3* pos, Vec3* target, Vec3* up) { + ssr_assert(pos && target && up); + Vec3 z; + vec3_minus(pos, target, &z); vec3_normalize(&z, &z); /*去除缩放影响*/ + Vec3 x; + vec3_cross(up, &z, &x); vec3_normalize(&x, &x); + Vec3 y; + vec3_cross(&z, &x, &y); vec3_normalize(&y, &y); + Mat4 m = { /*注意这里是列主序*/ + x.x, y.x, z.x, 0, + x.y, y.y, z.y, 0, + x.z, y.z, z.z, 0, + -vec3_dot(&x, pos), -vec3_dot(&y, pos), -vec3_dot(&z, pos), 1 + }; + mat4_multiply(&MATRIX, &m, &MATRIX); +} + +/* scalematrix * matrix */ +void ssr_scale(float sx, float sy, float sz) { + Vec3 scale = { sx, sy, sz }; + mat4_scale(&MATRIX, &scale, &MATRIX); +} + +void ssr_rotate(float angle, float x, float y, float z) { + Vec3 axis = { x, y, z }; + vec3_normalize(&axis, &axis); + mat4_rotate(&MATRIX, angle, &axis, &MATRIX); +} + +void ssr_translate(float x, float y, float z) { + Vec3 trans = { x, y, z }; + mat4_translate(&MATRIX, &trans, &MATRIX); +} + +void ssr_loadmatrix(Mat4* m) { + ssr_assert(m); + MATRIX = *m; +} + +void ssr_multmatrix(Mat4* m) { + ssr_assert(m); + mat4_multiply(&MATRIX, m, &MATRIX); +} + +void ssr_getmvp(Mat4* out) { + ssr_assert(out); + mat4_multiply(&GETMATRIX(MATRIX_VIEW), &GETMATRIX(MATRIX_MODEL), out); + mat4_multiply(&GETMATRIX(MATRIX_PROJECTION), out, out); +} + +void ssr_getmv(Mat4* out) { + ssr_assert(out); + mat4_multiply(&GETMATRIX(MATRIX_VIEW), &GETMATRIX(MATRIX_MODEL), out); +} + +void ssr_enable(uint mask) { + state.enable |= mask; +} + +void ssr_disable(uint mask) { + state.enable &= (~mask); +} + +bool ssr_isenable(uint mask) { + return state.enable & mask; +} + +void ssr_viewport(float l, float r, float b, float t) { +} + +void ssr_frustum(float l, float r, float b, float t, float n, float f) { + Mat4 m; + mat4_setfrustum(l, r, b, t, n, f, &m); + mat4_multiply(&MATRIX, &m, &MATRIX); +} + +void ssr_perspective(float fov, float aspect, float n, float f) { + Mat4 m; + mat4_setperspective(fov, aspect, n, f, &m); + mat4_multiply(&MATRIX, &m, &MATRIX); +} + +void ssr_present() { + if (config.dbuffer) { + memcpy(state.target, state.framebuffer, state.buffersize); + } +} + +void ssr_clearcolor(Color color) { + if (color == 0x00) { + memset(state.framebuffer, 0, state.buffersize); + } + else { + for (int i = 0; i < config.width*config.height; ++i) + state.framebuffer[i] = color; + } +} + +void ssr_cleardepth() { + memset(state.zbuffer, 0xff, sizeof(uint)*config.width*config.height); +} + +void ssr_putpoint(uint screenx, uint screeny, Color color) { + if (!contains(screenx, screeny, 0, config.width - 1, 0, config.height - 1)) + return; + BUFFER[screeny * config.width + screenx] = color; +} + +Color ssr_color(unsigned char r, unsigned char g, unsigned char b, unsigned char a) { + unsigned int c = (a << 24) | (r << 16) | (g << 8) | b; + return c; +} + +bool ssr_testdepth(uint x, uint y, uint depth) { + if (!contains(x, y, 0, config.width - 1, 0, config.height - 1)) + return 0; + uint off = x + y * config.width; + if (state.zbuffer[off] < depth) + return 0; + state.zbuffer[off] = depth; + return 1; +} + +bool ssr_testdepthf(uint x, uint y, float depth) { + uint d = UINT_MAX * (double)depth; + return ssr_testdepth(x, y, d); +} + +void ssrU_viewport(Vec2* p, Vec2* out) { + ssr_assert(p && out); + int halfw = config.width >> 1, halfh = config.height >> 1; + out->x = p->x * halfw + halfw + 0.5f; + out->y = halfh - p->y * halfh + 0.5f; +} + +void ssr_bindvertices(Vert** verts, int nverts, uint* indices, int nprims) { + ssr_assert(verts && indices); + state.verts = verts; + state.nverts = nverts; + state.indices = indices; + state.nprims = nprims; +} + +void ssr_useprogram(Program* program) { + ssr_assert(program); + state.program = program; +} + +void ssr_unuseprogram() { + state.program = NULL; +} + +static Vec4* homos = 0; +static uint tri[3]; /*三角形三个顶点的索引*/ +static Vec3 ab, ac; +static Mat4 mvp, mv; +void ssr_draw(ssr_PrimitiveType primitive) { + ssr_assert(state.verts && state.indices); + + ssr_getmvp(&mvp); + ssr_getmv(&mv); + state.uniforms.model = &GETMATRIX(MATRIX_MODEL); + state.uniforms.view = &GETMATRIX(MATRIX_VIEW); + state.uniforms.projection = &GETMATRIX(MATRIX_PROJECTION); + state.uniforms.mvp = &mvp; + state.uniforms.mv = &mv; + + VertexShaderIn vertIn = { NULL }; + if (!homos) + homos = ssrM_newvector(Vec4, state.nverts); + for (int i = 0; i < state.nverts; ++i) { + Vert* vert = state.verts[i]; + vertIn.vertex = vert; + state.program->vertexshader(&state.uniforms, &vertIn, &homos[i]); + } + + if (primitive == PRIMITIVE_TRIANGLE) { + for (int i = 0; i < state.nprims; ++i) { + tri[0] = state.indices[i * 3]; + tri[1] = state.indices[i * 3 + 1]; + tri[2] = state.indices[i * 3 + 2]; + + /* back face culling */ + if (ssr_isenable(ENABLEMASK_BACKFACECULL)) { + vec3_minus(&homos[tri[1]], &homos[tri[0]], &ab); + vec3_minus(&homos[tri[2]], &homos[tri[0]], &ac); + if (ab.x * ac.y - ab.y * ac.x <= 0) { + continue; + } + } + + ssrR_triangle( + &homos[tri[0]], + &homos[tri[1]], + &homos[tri[2]], + state.verts[tri[0]], + state.verts[tri[1]], + state.verts[tri[2]], + state.program, + &state.uniforms + ); + } + } + else if (primitive == PRIMITIVE_POINT) { + + } + else if (primitive == PRIMITIVE_LINE) { + + } + //ssrM_free(homos); +} + +void ssr_setuniformmat4(uint idx, Mat4* src) { + ssr_assert(src); + if (idx < 0 || idx > 7) return; + state.uniforms.var_mat4[idx] = *src; +} + +void ssr_setuniformvec4(uint idx, Vec4* src) { + ssr_assert(src); + if (idx < 0 || idx > 7) return; + state.uniforms.var_vec4[idx] = *src; +} + +void ssr_setuniformvec3(uint idx, Vec3* src) { + ssr_assert(src); + if (idx < 0 || idx > 7) return; + state.uniforms.var_vec3[idx] = *src; +} + +void ssr_setuniformvec2(uint idx, Vec2* src) { + ssr_assert(src); + if (idx < 0 || idx > 7) return; + state.uniforms.var_vec2[idx] = *src; +} diff --git a/src/core/device.h b/src/core/device.h new file mode 100644 index 0000000..aa6e73f --- /dev/null +++ b/src/core/device.h @@ -0,0 +1,98 @@ +#ifndef _SOFTSHADEROOM_DEVICE_H_ +#define _SOFTSHADEROOM_DEVICE_H_ + +#include "../math/math.h" +#include "../util/type.h" +#include "mem.h" +#include "limits.h" +#include "shader.h" +#include "rasterizer.h" +#include "vert.h" + +typedef enum ssr_MatrixMode { + MATRIX_MODEL = 0, + MATRIX_VIEW = 1, + MATRIX_PROJECTION = 2, +} ssr_MatrixMode; + +typedef struct ssr_Config { + int width, height; + bool dbuffer;/* double buffer? */ + Color* target; /* screen target buffer */ +} ssr_Config; + +typedef enum ssr_PrimitiveType { + PRIMITIVE_POINT, + PRIMITIVE_LINE, + PRIMITIVE_TRIANGLE, +} ssr_PrimitiveType; + +typedef enum ssr_EnableMask { + ENABLEMASK_BACKFACECULL = 1 , + ENABLEMASK_DEPTHTEST = 1 << 1, +} ssr_EnableMask; + +void ssr_init(ssr_Config* config); + +int ssr_getframebufferw(); +int ssr_getframebufferh(); + +void ssr_matrixmode(ssr_MatrixMode mode); +void ssr_loadidentity(); +void ssr_pushmatrix(); +void ssr_popmatrix(); +void ssr_lookat(Vec3* pos, Vec3* target, Vec3* up); + +/*后乘，意味着这些操作会被率先应用到物体上，应该按照反向顺序调用这些函数*/ +void ssr_translate(float x, float y, float z); +void ssr_rotate(float angle, float x, float y, float z); +void ssr_scale(float sx, float sy, float sz); + +void ssr_multmatrix(Mat4* m);/*右乘矩阵*/ +void ssr_loadmatrix(Mat4* m); +void ssr_frustum(float l, float r, float b, float t, float n, float f); +void ssr_perspective(float fov, float aspect, float n, float f); + +void ssr_viewport(float l, float r, float b, float t); + +void ssr_getmv(Mat4* out); +void ssr_getmvp(Mat4* out); + +void ssr_enable(uint mask); +void ssr_disable(uint mask); +bool ssr_isenable(uint mask); + +/* 绑定顶点数据 */ +void ssr_bindvertices(Vert** verts, int nverts, uint* indices, int nprims); +void ssr_unbindvertices(); + +void ssr_useprogram(Program* program); +void ssr_unuseprogram(); + +void ssr_setuniformmat4(uint idx, Mat4* src); +void ssr_setuniformvec4(uint idx, Vec4* src); +void ssr_setuniformvec3(uint idx, Vec3* src); +void ssr_setuniformvec2(uint idx, Vec2* src); + +void ssr_draw(ssr_PrimitiveType primitive); +void ssr_clearcolor(Color color); +void ssr_cleardepth(); + +/*直接绘制到缓冲区*/ +void ssr_putpoint(uint screenx, uint screeny, Color color); /*直接修改某个位置的颜色*/ + +void ssr_present(); + +Color ssr_color(unsigned char r, unsigned char g, unsigned char b, unsigned char a); + +//bool ssr_cutline(int* x0, int* y0, int* x1, int* y1); + +bool ssr_testdepth(uint x, uint y, uint depth); /*尝试写入深度，如果可绘制返回true，否则discard*/ +bool ssr_testdepthf(uint x, uint y, float depth01); + +/* +** Utils +*/ +void ssrU_viewport(Vec2* p, Vec2* out); + +#endif
\ No newline at end of file diff --git a/src/core/limits.h b/src/core/limits.h new file mode 100644 index 0000000..4ef0b75 --- /dev/null +++ b/src/core/limits.h @@ -0,0 +1,6 @@ +#ifndef _SOFTSHADEROOM_LIMITS_H_ +#define _SOFTSHADEROOM_LIMITS_H_ + +#define MATRIXDEPTH 32 + +#endif
\ No newline at end of file diff --git a/src/core/mem.c b/src/core/mem.c new file mode 100644 index 0000000..e69de29 --- /dev/null +++ b/src/core/mem.c diff --git a/src/core/mem.h b/src/core/mem.h new file mode 100644 index 0000000..1a763f9 --- /dev/null +++ b/src/core/mem.h @@ -0,0 +1,20 @@ +#ifndef _SOFTSHADEROOM_MEM_H_ +#define _SOFTSHADEROOM_MEM_H_ + +#include <memory.h> + +/* +** 内存分配相关 +*/ + +void ssrM_freemem(); +void ssrM_realloc(); + +#define ssrM_zero(obj, size) memset(obj, 0, size) +#define ssrM_new(T) (T*)malloc(sizeof(T)) +#define ssrM_newvector(T, c) ((T*)calloc(c, sizeof(T))) +#define ssrM_free(obj) free(obj) + +#define cast(v, T) (*(T*)(&v)) + +#endif
\ No newline at end of file diff --git a/src/core/rasterizer.c b/src/core/rasterizer.c new file mode 100644 index 0000000..0c8c6eb --- /dev/null +++ b/src/core/rasterizer.c @@ -0,0 +1,180 @@ +#include "rasterizer.h" +#include "device.h" +#include "shader.h" +#include "../math/math.h" +#include "../util/assert.h" + +void ssrR_putline(int x0, int y0, int x1, int y1, Color color) { + int steep = 0; + if (abs(x0 - x1) < abs(y0 - y1)) { + swapi(x0, y0); + swapi(x1, y1); + steep = 1; + } + if (x0 > x1) { + swapi(x0, x1); + swapi(y0, y1); + } + int dx = x1 - x0; + int dy = y1 - y0; + int derror2 = abs(dy) << 1; + int error2 = 0; + int y = y0; + for (int x = x0; x <= x1; ++x) { + if (steep) { + ssr_putpoint(y, x, color); + } + else { + ssr_putpoint(x, y, color); + } + error2 += derror2; + if (error2 > dx) { + y += (y1 > y0 ? 1 : -1); + error2 -= (dx << 1); + } + } +} + +/*使用辅助梯形计算三角形面积，《3D数学基础》P245*/ +float ssrR_area(Vec2* v1, Vec2* v2, Vec2* v3) { + ssr_assert(v1 && v2 && v3); + float area = 0.5f * ((v1->y - v2->y) * (v3->x - v2->x) + (v3->y - v2->y) * (v2->x - v1->x)); + return area; +} + +/*from https://github.com/ssloy/tinyrenderer*/ +static Vec3 s[2], u; +bool ssrR_barycentric(Vec2* A, Vec2* B, Vec2* C, Vec2* p, Vec3* out) { + ssr_assert(A && B && C && p && out); + s[0].x = C->x - A->x; s[0].y = B->x - A->x; s[0].z = A->x - p->x; + s[1].x = C->y - A->y; s[1].y = B->y - A->y; s[1].z = A->y - p->y; + vec3_cross(&s[0], &s[1], &u); + if (compare(u.z, 0)) { + return 0; + } else { + float uz = 1.f / u.z; + out->x = 1 - (u.x + u.y) * uz; + out->y = u.y * uz; + out->z = u.x * uz; + return 1; + } +} + +bool ssrR_ispointintriangle(Vec2* A, Vec2* B, Vec2* C, Vec2* p) { + ssr_assert(A && B && C && p); + Vec3 c; ssrR_barycentric(A, B, C, p, &c); + return c.x >= 0 && c.y >= 0 && c.z >= 0; +} + +void ssrR_center(Vec2* A, Vec2* B, Vec2* C, Vec2* out) { + ssr_assert(A && B && C && out); + float k = 1 / 3.f; + out->x = k * (A->x + B->x + C->x); + out->y = k * (A->y + B->y + C->y); +} + +/* +** clipc 按列排放的三角形三个顶点的齐次裁剪坐标 +** program shader程序 +*/ +static Vec3 SA, SB, SC; // screen coords +static Vec3 bcp; // 重心坐标 +void _ssrR_triangle(Vec4* CA, Vec4* CB, Vec4* CC, Vert* A, Vert* B, Vert* C, Program* program, UniformCollection* uniforms) { + ssr_assert(CA && CB && CC && program); + vec4_dividew(CA, &SA); ssrU_viewport(&SA, &SA); + vec4_dividew(CB, &SB); ssrU_viewport(&SB, &SB); + vec4_dividew(CC, &SC); ssrU_viewport(&SC, &SC); + /*根据屏幕空间的y值排序, A的y最小，B居中，C最大*/ + Vec3 *sa = &SA, *sb = &SB, *sc = &SC, *tmp; Vec4* v4tmp; Vert* vtemp; + if (sb->y < sa->y) { sa = &SB; sb = &SA; v4tmp = CA; CA = CB; CB = v4tmp; vtemp = A; A = B; B = vtemp; } + if (sc->y < sb->y) { tmp = sc; sc = sb; sb = tmp; v4tmp = CC; CC = CB; CB = v4tmp; vtemp = B; B = C; C = vtemp; } + if (sc->y < sa->y) { tmp = sc; sc = sa; sa = tmp; v4tmp = CC; CC = CA; CA = v4tmp; vtemp = C; C = A; A = vtemp; } + + float invkAC = (SC.x - SA.x) / (SC.y - SA.y + EPSILON); + float invkAB = (SB.x - SA.x) / (SB.y - SA.y + EPSILON); + float invkBC = (SC.x - SB.x) / (SC.y - SB.y + EPSILON); + + int order = 1; if (((sb->y - sa->y) / (sb->x - sa->x) - invkAC) * invkAC > 0) { order = -1; } + + FragmentShaderIn in = {A, B, C, 0}; + Color color; + Vec2 l, r, p; + float depth; Vec3 cdepth = {SA.z, SB.z, SC.z}; + float CAw = 1.f / CA->w, CBw = 1.f / CB->w, CCw = 1.f / CC->w; + for (p.y = sa->y; p.y <= sc->y; ++p.y) { + l.y = p.y; + r.y = p.y; + if (compare(SC.y, SA.y)) r.x = sc->x; + else r.x = sa->x + invkAC * (p.y - sa->y); /*AC*/ + if (p.y <= sb->y) { + if (compare(SB.y, SA.y)) l.x = sb->x; + else l.x = sa->x + invkAB * (p.y - sa->y); + } else { + if (compare(SC.y, SB.y)) l.x = sc->x; + else l.x = sb->x + invkBC * (p.y - sb->y); + } + for (p.x = l.x; order * (r.x - p.x) >= 0; p.x += order) { + /*计算线性的（齐次裁剪空间）的重心坐标，分三步*/ + if (!ssrR_barycentric(sa, sb, sc, &p, &bcp)) { + continue; /*discard fragment*/ + } + bcp.x *= CAw; bcp.y *= CBw; bcp.z *= CCw; + vec3_scale(&bcp, 1.f / (bcp.x + bcp.y + bcp.z), &bcp); + /*depth test*/ + depth = vec3_dot(&bcp, &cdepth) * 0.5f + 0.5f; + if (!ssr_testdepthf(p.x, p.y, depth)) { + continue; /*discard fragment*/ + } + /*enter fragment shader*/ + in.bc = &bcp; + program->fragmentshader(uniforms, &in, &color); + ssr_putpoint(p.x, p.y, color); + } + } +} + +static Vec2 bboxmin, bboxmax; +void ssrR_triangle(Vec4* CA, Vec4* CB, Vec4* CC, Vert* A, Vert* B, Vert* C, Program* program, UniformCollection* uniforms) { + ssr_assert(CA && CB && CC && program); + vec4_dividew(CA, &SA); ssrU_viewport(&SA, &SA); + vec4_dividew(CB, &SB); ssrU_viewport(&SB, &SB); + vec4_dividew(CC, &SC); ssrU_viewport(&SC, &SC); + bboxmin.x = min(SA.x, min(SB.x, SC.x)); + bboxmax.x = max(SA.x, max(SB.x, SC.x)); + bboxmin.y = min(SA.y, min(SB.y, SC.y)); + bboxmax.y = max(SA.y, max(SB.y, SC.y)); + float CAw = 1.f / CA->w, CBw = 1.f / CB->w, CCw = 1.f / CC->w; + float depth; Vec3 cdepth = { SA.z, SB.z, SC.z }; + FragmentShaderIn in = { A, B, C, 0 }; + Vec2 p; Color color; + for (p.y = bboxmin.y; p.y <= bboxmax.y; ++p.y) { + for (p.x = bboxmin.x; p.x <= bboxmax.x; ++p.x) { + /*计算线性的（齐次裁剪空间）的重心坐标，分三步*/ + if (!ssrR_barycentric(&SA, &SB, &SC, &p, &bcp)) { + continue; /*discard fragment*/ + } + if (bcp.x < 0 || bcp.y < 0 || bcp.z < 0) { + continue; /*discard fragment*/ + } + bcp.x *= CAw; bcp.y *= CBw; bcp.z *= CCw; + vec3_scale(&bcp, 1.f / (bcp.x + bcp.y + bcp.z), &bcp); + + /*深度测试*/ + if (ssr_isenable(ENABLEMASK_DEPTHTEST)) { + depth = vec3_dot(&bcp, &cdepth) * 0.5f + 0.5f; + if (!ssr_testdepthf(p.x, p.y, depth)) { + continue; /*discard fragment*/ + } + } + + /*enter fragment shader*/ + in.bc = &bcp; + bool discad = program->fragmentshader(uniforms, &in, &color); + if (!discad) + ssr_putpoint(p.x, p.y, color); + } + } + //ssrR_putline(SA.x, SA.y, SB.x, SB.y, 0xffff0000); + //ssrR_putline(SA.x, SA.y, SC.x, SC.y, 0xffff0000); + //ssrR_putline(SC.x, SC.y, SB.x, SB.y, 0xffff0000); +}
\ No newline at end of file diff --git a/src/core/rasterizer.h b/src/core/rasterizer.h new file mode 100644 index 0000000..e35b24a --- /dev/null +++ b/src/core/rasterizer.h @@ -0,0 +1,17 @@ +#ifndef _SOFTSHADEROOM_GEOMETRY_H_ +#define _SOFTSHADEROOM_GEOMETRY_H_ + +#include "../math/math.h" +#include "shader.h" + +void ssrR_putline(int x0, int y0, int x1, int y1, Color color); + +float ssrR_area(Vec2* v1, Vec2* v2, Vec2* v3); /*计算三角形面积，注意v1v2v3按逆时针顺序定义*/ + +bool ssrR_barycentric(Vec2* A, Vec2* B, Vec2* C, Vec2* p, Vec3* out); /*计算重心坐标*/ +void ssrR_center(Vec2* A, Vec2* B, Vec2* C, Vec2* out); /*获得重心*/ +bool ssrR_ispointintriangle(Vec2* A, Vec2* B, Vec2* C, Vec2* p); + +void ssrR_triangle(Vec4* CA, Vec4* CB, Vec4* CC, Vert* A, Vert* B, Vert* C, Program* program, UniformCollection* uniforms); /*绘制三角形图元*/ + +#endif
\ No newline at end of file diff --git a/src/core/shader.c b/src/core/shader.c new file mode 100644 index 0000000..fe25408 --- /dev/null +++ b/src/core/shader.c @@ -0,0 +1,34 @@ +#include "shader.h" +#include "vert.h" +#include "device.h" + +void ssrS_bcpcolor(Vec3* bc, Color A, Color B, Color C, Color* out) { + ssr_assert(bc && out); + *out = ssr_color( + bc->A * COLOR_R(A) + bc->B * COLOR_R(B) + bc->C * COLOR_R(C), + bc->A * COLOR_G(A) + bc->B * COLOR_G(B) + bc->C * COLOR_G(C), + bc->A * COLOR_B(A) + bc->B * COLOR_B(B) + bc->C * COLOR_B(C), + bc->A * COLOR_A(A) + bc->B * COLOR_A(B) + bc->C * COLOR_A(C) + ); +} + +void ssrS_bcpvec2(Vec3* bc, Vec2* A, Vec2* B, Vec2* C, Vec2* out) { + ssr_assert(bc && A && B && C && out); + out->x = bc->A * A->x + bc->B * B->x + bc->C * C->x; + out->y = bc->A * A->y + bc->B * B->y + bc->C * C->y; +} + +void ssrS_bcpvec3(Vec3* bc, Vec3* A, Vec3* B, Vec3* C, Vec3* out) { + ssr_assert(bc && A && B && C && out); + out->x = bc->A * A->x + bc->B * B->x + bc->C * C->x; + out->y = bc->A * A->y + bc->B * B->y + bc->C * C->y; + out->z = bc->A * A->z + bc->B * B->z + bc->C * C->z; +} + +void ssrS_bcpvec4(Vec3* bc, Vec4* A, Vec4* B, Vec4* C, Vec4* out) { + ssr_assert(bc && A && B && C && out); + out->x = bc->A * A->x + bc->B * B->x + bc->C * C->x; + out->y = bc->A * A->y + bc->B * B->y + bc->C * C->y; + out->z = bc->A * A->z + bc->B * B->z + bc->C * C->z; + out->w = bc->A * A->w + bc->B * B->w + bc->C * C->w; +} diff --git a/src/core/shader.h b/src/core/shader.h new file mode 100644 index 0000000..8da8ae3 --- /dev/null +++ b/src/core/shader.h @@ -0,0 +1,89 @@ +#ifndef _SOFTSHADEROOM_SHADER_H_ +#define _SOFTSHADEROOM_SHADER_H_ + +#include "../math/math.h" +#include "vert.h" + +typedef struct UniformCollection { + /*built in varaibles*/ + Mat4* model; + Mat4* view; + Mat4* projection; + Mat4* mvp; + Mat4* mv; + /*userside variables*/ + Mat4 var_mat4[8]; + Vec4 var_vec4[8]; + Vec3 var_vec3[8]; + Vec2 var_vec2[8]; +} UniformCollection; + +#define ssrum4(i) (uniforms->var_mat4[i]) +#define ssruv2(i) (uniforms->var_vec2[i]) +#define ssruv3(i) (uniforms->var_vec3[i]) +#define ssruv4(i) (uniforms->var_vec4[i]) + +typedef struct VertexShaderIn { + Vert* vertex; +} VertexShaderIn; + +typedef void(*VertexShader)(UniformCollection* uniforms, VertexShaderIn* in, Vec4* homocoord); + +typedef struct FragmentShaderIn { + Vert *A, *B, *C; + Vec3* bc; +} FragmentShaderIn; + +typedef bool(*FragmentShader)(UniformCollection* uniforms, FragmentShaderIn* in, Color* color); + +typedef struct Program { + VertexShader vertexshader; + FragmentShader fragmentshader; +} Program; + +void ssrS_bcpcolor(Vec3* bc, Color A, Color B, Color C, Color* out); +void ssrS_bcpvec2(Vec3* bc, Vec2* A, Vec2* B, Vec2* C, Vec2* out); +void ssrS_bcpvec3(Vec3* bc, Vec3* A, Vec3* B, Vec3* C, Vec3* out); +void ssrS_bcpvec4(Vec3* bc, Vec4* A, Vec4* B, Vec4* C, Vec4* out); + +/* +** 顶点数据外，还提供额外的寄存器 +*/ +enum RegisterType { + REGTYPE_VEC2 = 1, + REGTYPE_VEC4 = 2, + REGTYPE_VEC3 = 3, +}; + +typedef struct Register { + uint length; + uint type; + union { + Vec4* v4; + Vec2* v2; + Vec3* v3; + }; +} Register; + +void ssrS_setregvec4_01i(uint idx, Vec4* value); +void ssrS_setregvec4_02i(uint idx, Vec4* value); +void ssrS_setregvec2_01i(uint idx, Vec2* value); +void ssrS_setregvec2_02i(uint idx, Vec2* value); +void ssrS_setregvec3_01i(uint idx, Vec3* value); +void ssrS_setregvec3_02i(uint idx, Vec3* value); + +Vec4* ssrS_getregvec4_01i(uint idx); +Vec4* ssrS_getregvec4_02i(uint idx); +Vec2* ssrS_getregvec2_01i(uint idx); +Vec2* ssrS_getregvec2_02i(uint idx); +Vec3* ssrS_getregvec3_01i(uint idx); +Vec3* ssrS_getregvec3_02i(uint idx); + +Vec4* ssrS_getregvec4_01(Vert* vert); +Vec4* ssrS_getregvec4_02(Vert* vert); +Vec2* ssrS_getregvec2_01(Vert* vert); +Vec2* ssrS_getregvec2_02(Vert* vert); +Vec3* ssrS_getregvec3_01(Vert* vert); +Vec3* ssrS_getregvec3_02(Vert* vert); + +#endif
\ No newline at end of file diff --git a/src/core/tris.c b/src/core/tris.c new file mode 100644 index 0000000..e69de29 --- /dev/null +++ b/src/core/tris.c diff --git a/src/core/tris.h b/src/core/tris.h new file mode 100644 index 0000000..1f4208f --- /dev/null +++ b/src/core/tris.h @@ -0,0 +1,17 @@ +#ifndef _SOFTSHADEROOM_TRIS_H_ +#define _SOFTSHADEROOM_TRIS_H_ + +#include "../util/type.h" + +/* 用索引表示的三角形，在一系列顶点中索引 */ +typedef union Triangle +{ + uint v[3]; + struct { + uint v1, v2, v3; + }; +} Triangle; + + + +#endif
\ No newline at end of file diff --git a/src/core/vert.c b/src/core/vert.c new file mode 100644 index 0000000..0e62087 --- /dev/null +++ b/src/core/vert.c @@ -0,0 +1,42 @@ +#include "vert.h" + +Vert* vert_new(uint comp) { + Vert* vert = ssrM_new(Vert); + vert->comp = comp; + vert->position = ssrM_new(Vec3); + vert->normal = ssrM_new(Vec3); + vert->tangent = ssrM_new(Vec3); + vert->uv = ssrM_new(Vec2); + return vert; +} + +void vert_init(Vert* v, uint comp) { + ssr_assert(v); + v->comp = comp; + v->index = 0; + if ((comp & VERTMASK_POSITION) && !v->position) + v->position = ssrM_new(Vec3); + else if (!(comp & VERTMASK_POSITION) && v->position) + ssrM_free(v->position); + if ((comp & VERTMASK_NORMAL) && !v->normal) + v->normal = ssrM_new(Vec3); + else if (!(comp & VERTMASK_NORMAL) && v->normal) + ssrM_free(v->normal); + if ((comp & VERTMASK_TANGENT) && !v->tangent) + v->tangent = ssrM_new(Vec3); + else if (!(comp & VERTMASK_TANGENT) && v->tangent) + ssrM_free(v->tangent); + if ((comp & VERTMASK_UV) && !v->uv) + v->uv = ssrM_new(Vec2); + else if (!(comp & VERTMASK_UV) && v->uv) + ssrM_free(v->uv); +} + +void vert_free(Vert* v) { + ssr_assert(v); + if (v->position) ssrM_free(v->position); + if (v->normal) ssrM_free(v->normal); + if (v->tangent) ssrM_free(v->tangent); + if (v->uv) ssrM_free(v->uv); + ssrM_free(v); +}
\ No newline at end of file diff --git a/src/core/vert.h b/src/core/vert.h new file mode 100644 index 0000000..a2fdb74 --- /dev/null +++ b/src/core/vert.h @@ -0,0 +1,48 @@ +#ifndef _SOFTSHADEROOM_VERT_H_ +#define _SOFTSHADEROOM_VERT_H_ + +#include "mem.h" +#include "color.h" +#include "../math/math.h" +#include "../util/type.h" +#include "../util/assert.h" + +//typedef uint Color; // ARGB by default + +typedef uint Color; // ARGB by default + +#define COLOR_A(c) ((c >> 24) & 0xff) +#define COLOR_R(c) ((c >> 16) & 0xff) +#define COLOR_G(c) ((c >> 8) & 0xff) +#define COLOR_B(c) (c & 0xff) + +typedef Vec4 Color32; // 4个float的颜色 + +Color color32_tocolor(Color32* c); +void color_tocolor32(Color c, Color32* out); +void color32_saturate(Color32* c); + +enum VertMask { + VERTMASK_POSITION = 1 , + VERTMASK_NORMAL = 1 << 1, + VERTMASK_COLOR = 1 << 2, + VERTMASK_TANGENT = 1 << 3, + VERTMASK_UV = 1 << 4, +}; + +/*标准的顶点结构，禁止修改*/ +typedef struct Vert { + uint index; + uint comp; + Vec3* position; + Vec3* normal; + Vec3* tangent; + Vec2* uv; + Color color; +} Vert; + +Vert* vert_new(uint comp); /*新建vert*/ +void vert_init(Vert* v, uint comp); /*如果某些component没创建，创建它*/ +void vert_free(Vert* v); /*释放vert*/ + +#endif
\ No newline at end of file diff --git a/src/example/example.h b/src/example/example.h new file mode 100644 index 0000000..7879257 --- /dev/null +++ b/src/example/example.h @@ -0,0 +1,21 @@ +#ifndef _SSR_EXAMPLE_H_ +#define _SSR_EXAMPLE_H_ + +#include "../ssr.h" +#include "../util/type.h" +#include "SDL2/SDL.h" + +#define EXAMPLE(i)\ +extern void onload##i(void*);\ +extern void onevent##i(void*);\ +extern void onupdate##i(void*);\ +extern void ondraw##i(void*); + +#define EXAMPLECUR texture + +EXAMPLE(cube); +EXAMPLE(line); +EXAMPLE(dot); +EXAMPLE(texture); + +#endif
\ No newline at end of file diff --git a/src/example/example_cube.c b/src/example/example_cube.c new file mode 100644 index 0000000..caedd64 --- /dev/null +++ b/src/example/example_cube.c @@ -0,0 +1,73 @@ +#include "example.h" + +/*正方体*/ +Vec3 verts[] = { + // front face + {1, 1, 1}, {-1, 1, 1}, {-1, -1, 1}, {1, -1, 1}, + // back face + {1, 1, -1}, {-1, 1, -1}, {-1, -1, -1}, {1, -1, -1}, +}; + +Color colors[] = { + 0xffff0000, 0xff00ff00, 0xffff00ff, 0xff00ffff, + 0xff0000ff, 0xff000000, 0xffffff00, 0xffffffff, +}; + +int cube[] = { + 0, 2, 1, 0, 3, 2, + 1, 2, 5, 2, 6, 5, + 4, 5, 6, 4, 6, 7, + 0, 4, 7, 0, 7, 3, + 0, 1, 4, 1, 5, 4, + 2, 3, 6, 3, 7, 6 +}; +Mat4 m; +void onloadcube(void* data) { +} + +void oneventcube(void* data) { + SDL_Event* e = (SDL_Event*)data; +} + +float _t = 0; + +void onupdatecube(void*data) { + uint dt = *(uint*)data; + ssr_matrixmode(MATRIX_MODEL); + ssr_loadidentity(); + ssr_translate(0, 0, -3); + ssr_rotate(360 * sin(_t += 0.001f), 1, 1, 1); + ssr_matrixmode(MATRIX_PROJECTION); + ssr_loadidentity(); + ssr_perspective(100 + 20 * sin(_t * 10), 1.25f, -0.1f, -100); + ssr_matrixmode(MATRIX_VIEW); + ssr_loadidentity(); + Vec3 pos = { 0,0,0 }, target = { 0,0,-1 }, up = { 0,1,0 }; + ssr_lookat(&pos, &target, &up); + ssr_getmvp(&m); +} + +void ondrawcube(void*data) { + ssr_clearcolor(0); + Vec2 proj[8]; + + for (int i = 0; i < 8; ++i) { + Vec4 v = { verts[i].x, verts[i].y ,verts[i].z ,1 }, temp; + mat4_applytovec4(&m, &v, &temp); + temp.x /= temp.w; + temp.y /= temp.w; + temp.z /= temp.w; + //vec4_print(&temp); + proj[i].x = temp.x; + proj[i].y = temp.y; + } + for (int j = 1; j < sizeof(cube) / sizeof(int); ++j) { + int fromx = proj[cube[j]].x * 250.f + 250, fromy = 200 - proj[cube[j]].y * 200.f; + int tox = proj[cube[j - 1]].x * 250.f + 250, toy = 200 - proj[cube[j - 1]].y * 200.f; + ssrR_putline(fromx, fromy , tox, toy, 0xffff0000); + } + + Vec2 v1 = { 0, 0 }, v2 = { 3, 1 }, v3 = {1, 5}; + float area = ssrR_area(&v1, &v2, &v3); + printf("%f\n", area); +} diff --git a/src/example/example_dot.c b/src/example/example_dot.c new file mode 100644 index 0000000..0a88084 --- /dev/null +++ b/src/example/example_dot.c @@ -0,0 +1,36 @@ +#include "example.h" +#include "../core/rasterizer.h" + +Vec3 pos = { 0, 0, 2 }, target = { 0,0,-1 }, up = { 0, 1, 0 }; + +float dot[] = { 0, 0, -2 }; + +void onloaddot(void* data) { + +} + +void oneventdot(void* data) { + SDL_Event* e = (SDL_Event*)data; +} + +void onupdatedot(void*data) { + uint dt = *(uint*)data; + ssr_matrixmode(MATRIX_VIEW); + ssr_lookat(&pos, &target, &up); + ssr_matrixmode(MATRIX_PROJECTION); +} + +float j = 0; + +void ondrawdot(void*data) { + ssr_clearcolor(0xffffffff); + + //for (int i = -50; i < 200; ++i) { + // ssr_putpoint(i * sin(j += 0.1f), 300 * cos(j), ssr_color(0xff, 0, 0xff, 0)); + //} + + ssrR_putline(-100, 0, 100, 200, 0xffff0000); + +// printf("%u\n", (UINT_MAX * (double)0.5F)); + +} diff --git a/src/example/example_texture.c b/src/example/example_texture.c new file mode 100644 index 0000000..5411a12 --- /dev/null +++ b/src/example/example_texture.c @@ -0,0 +1,101 @@ +#include "example.h" + +Vert** quad; +static Vec3 verts[] = { + // front face + {1, 1, 1}, {-1, 1, 1}, {-1, -1, 1}, {1, -1, 1}, + // back face + {1, 1, -1}, {-1, 1, -1}, {-1, -1, -1}, {1, -1, -1}, +}; +static Color colors[] = { + 0xffff0000, 0xff00ff00, 0xff0000ff, 0xffff00ff, + 0xffaa28aa, 0xffFFC58E, 0xffA100FF, 0xffFAFF00, +}; + +int face[] = { + 0, 1, 2, 0, 2, 3, + 1, 5, 2, 2, 5, 6, + 4, 6, 5, 4, 7, 6, + 0, 7, 4, 0, 3, 7, + 0, 4, 1, 1, 4, 5, + 2, 6, 3, 3, 6, 7 +}; + +static Mat4 m; + +void vert(UniformCollection* uniforms, VertexShaderIn* in, Vec4* homocoord) { + static Vec4 p; p.xyz = *in->vertex->position; p.w = 1; + mat4_applytovec4(uniforms->mvp, &p, homocoord); +} +bool frag(UniformCollection* uniforms, FragmentShaderIn* in, Color* color) { + ssrS_bcpcolor(in->bc, in->A->color, in->B->color, in->C->color, color); + //ssrum4(1); + //Vec2 uv; + //ssrS_bcpvec2(in->bc, in->A->uv, in->B->uv, in->C->uv, &uv); + //int x = uv.x * 10; + //int y = uv.y * 10; + //if (x % 2 && y % 2) { + // *color = 0xffffffff; + //} else { + // *color = 0xff111111; + //} + return 0; +} +Program program = { vert, frag }; + +void onloadtexture(void* data) { + ssr_matrixmode(MATRIX_PROJECTION); + ssr_loadidentity(); + ssr_perspective(100, 1.25f, 0.1f, 10); + ssr_matrixmode(MATRIX_VIEW); + ssr_loadidentity(); + //Vec3 pos = { 0,0,0 }, target = { 0,0,-1 }, up = { 0,1,0 }; + //ssr_lookat(&pos, &target, &up); + + /*设置顶点数据*/ + quad = ssrM_newvector(Vert*, 8); + quad[0] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[1] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[2] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[3] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[4] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[5] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[6] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + quad[7] = vert_new(VERTMASK_POSITION | VERTMASK_COLOR); + *quad[0]->position = verts[0]; quad[0]->color = colors[0]; quad[0]->index = 0; + *quad[1]->position = verts[1]; quad[1]->color = colors[1]; quad[1]->index = 1; + *quad[2]->position = verts[2]; quad[2]->color = colors[2]; quad[2]->index = 2; + *quad[3]->position = verts[3]; quad[3]->color = colors[3]; quad[3]->index = 3; + *quad[4]->position = verts[4]; quad[4]->color = colors[4]; quad[4]->index = 4; + *quad[5]->position = verts[5]; quad[5]->color = colors[5]; quad[5]->index = 5; + *quad[6]->position = verts[6]; quad[6]->color = colors[6]; quad[6]->index = 6; + *quad[7]->position = verts[7]; quad[7]->color = colors[7]; quad[7]->index = 7; + + ssr_bindvertices(quad, 8, face, 12); + + ssr_useprogram(&program); + ssr_enable(ENABLEMASK_BACKFACECULL); + ssr_enable(ENABLEMASK_DEPTHTEST); +} + +void oneventtexture(void* data) { + SDL_Event* e = (SDL_Event*)data; +} + +static float _t = 0; + +void onupdatetexture(void*data) { + uint dt = *(uint*)data; + ssr_matrixmode(MATRIX_MODEL); + ssr_loadidentity(); + ssr_translate(0, 0, -3); +// ssr_rotate(-30, 1, 1, 0); + ssr_rotate(360 * sin(_t += dt/5000.f), 1, 1, 0); +} + +void ondrawtexture(void*data) { + ssr_clearcolor(0x00); + ssr_cleardepth(); + + ssr_draw(PRIMITIVE_TRIANGLE); +} diff --git a/src/extend/camera.c b/src/extend/camera.c new file mode 100644 index 0000000..e69de29 --- /dev/null +++ b/src/extend/camera.c diff --git a/src/extend/camera.h b/src/extend/camera.h new file mode 100644 index 0000000..4f4ba47 --- /dev/null +++ b/src/extend/camera.h @@ -0,0 +1,10 @@ +#ifndef _SOFTSHADEROOM_CAMERA_H_ +#define _SOFTSHADEROOM_CAMERA_H_ + +#include "../math/math.h" + +typedef struct Camera { + Vec3 front; +} Camera; + +#endif
\ No newline at end of file diff --git a/src/extend/transform.c b/src/extend/transform.c new file mode 100644 index 0000000..e69de29 --- /dev/null +++ b/src/extend/transform.c diff --git a/src/extend/transform.h b/src/extend/transform.h new file mode 100644 index 0000000..b366305 --- /dev/null +++ b/src/extend/transform.h @@ -0,0 +1,14 @@ +#include "../math/math.h" + +/* +** Transform要注意按 scale -> rotation -> position 顺序计算 +*/ +typedef struct Transform { + Vec3 position; + Vec3 scale; + Quat rotation; +} Transform; + +void transformtranslate(Transform* trans, Vec4* v, Vec4* out); + +Vec4 transform_translate(Transform* trans, Vec4* v); diff --git a/src/main.c b/src/main.c new file mode 100644 index 0000000..d06803a --- /dev/null +++ b/src/main.c @@ -0,0 +1,74 @@ +#include <SDL2/SDL.h> +#include <stdio.h> +#include <windows.h> +#include "math/math.h" +#include "util/assert.h" +#include "ssr.h" +#include "example/example.h" + +SDL_Surface* suf; + +#define SCREEN_WIDTH 500 +#define SCREEN_HEIGHT 400 + +typedef void(*F)(void*); +F onload; +F onupdate; +F onevent; +F ondraw; + +/*macro这里需要绕一下弯*/ +/*https://stackoverflow.com/questions/1489932/how-to-concatenate-twice-with-the-c-preprocessor-and-expand-a-macro-as-in-arg*/ +#define SETEXAMPLEF(f, e) \ +f = f##e; + +#define SETEXAMPLE(i) \ +SETEXAMPLEF(onload, i)\ +SETEXAMPLEF(ondraw, i)\ +SETEXAMPLEF(onevent, i)\ +SETEXAMPLEF(onupdate, i) + +int main(int argc, char* argv[]) +{ + if (SDL_Init(SDL_INIT_VIDEO) < 0) + return 1; + SDL_Window* wnd = SDL_CreateWindow("Test", SDL_WINDOWPOS_UNDEFINED, SDL_WINDOWPOS_UNDEFINED, SCREEN_WIDTH, SCREEN_HEIGHT, SDL_WINDOW_SHOWN); + SDL_Event e; + suf = SDL_GetWindowSurface(wnd); + /* ARGB format */ + ssr_assert(suf->format->BitsPerPixel == 32); + ssr_assert(suf->format->Rshift == 16); + ssr_assert(suf->format->Gshift == 8); + ssr_assert(suf->format->Bshift == 0); + /* init ssr */ + ssr_Config config = { + SCREEN_WIDTH, SCREEN_HEIGHT,/* screen size */ + 0, /* double buffer */ + suf->pixels /* screen buffer */ + }; + ssr_init(&config); + SETEXAMPLE(EXAMPLECUR); + onload(0); + /* main loop */ + uint previous = SDL_GetTicks(); + uint dt = 0; + while (1) { + while (SDL_PollEvent(&e)) { + if (e.type == SDL_QUIT) { + goto quit; + } else { + onevent(&e); + } + } + dt = SDL_GetTicks() - previous; + previous = dt + previous; + onupdate(&dt); + ondraw(0); + ssr_present(); + SDL_UpdateWindowSurface(wnd); + Sleep(10); /*100fps limit*/ + } + quit: + SDL_Quit(); + return 0; +} diff --git a/src/math/math.c b/src/math/math.c new file mode 100644 index 0000000..4cca29e --- /dev/null +++ b/src/math/math.c @@ -0,0 +1,21 @@ +#include "math.h" + +char printbuffer[2048] = { 0 }; + +float rsqrt(float number) { + long i; + float x2, y; + const float threehalfs = 1.5F; + + x2 = number * 0.5F; + y = number; + i = *(long *)&y; + i = 0x5f3759df - (i >> 1); + y = *(float *)&i; + y = y * (threehalfs - (x2 * y * y)); + return y; +} + +float lerp(float a, float b, float t) { + return a * (1 - t) + b * t; +}
\ No newline at end of file diff --git a/src/math/math.h b/src/math/math.h new file mode 100644 index 0000000..4498920 --- /dev/null +++ b/src/math/math.h @@ -0,0 +1,306 @@ +#ifndef _SOFTSHADEROOM_MATH_H_ +#define _SOFTSHADEROOM_MATH_H_ + +#include <stdio.h> +#include <string.h> +#include <math.h> +#include <limits.h> + +#include "../util/type.h" + +#define PI 3.141592653f +#define RAD2DEG 57.295779523f /*180.f/PI*/ +#define DEG2RAG 0.0174532925f /*PI/180.f*/ +#define EPSILON 0.000001f + +/* 用来打印的公共buffer */ +extern char printbuffer[2048]; + +/* +** 数学函数 +*/ +#define min(a, b) ((a) < (b) ? (a) : (b)) +#define max(a, b) ((a) > (b) ? (a) : (b)) +#define clamp(v, l, h) ((v) > (l) ? ((v) < (h) ? (v) : (h)) : (l)) +#define absf(v) ((v )> 0 ? (v ): -(v)) +#define radians(angle) (angle * DEG2RAG) +#define degree(rad) (rad * RAD2DEG) +#define compare(v1, v2) (absf((v1) - (v2)) < EPSILON) +#define swapi(a, b) {int temp = a; a = b; b = temp;} +float rsqrt(float n); +float lerp(float from, float to, float t); + +/* +** 二维向量，用来做屏幕上的一些计算 +*/ +typedef struct Vec2 { + float x, y; +} Vec2; + +/* +** 三维向量，用来做三维空间的计算 +*/ +typedef union Vec3 { + struct { + float x, y, z; + }; + struct { + float A, B, C; /*重心坐标*/ + }; + Vec2 xy; +} Vec3; + +/* +** 齐次坐标，列主项，平移变换和透视投影需要 +*/ +typedef union Vec4 { + struct { + float x, y, z, w; + }; + struct { + float r, g, b, a; + }; + Vec3 xyz; +} Vec4; + +/* +** 用来可视化四元数，欧拉角默认使用角度存储，用euler_deg2rad()转弧度 +*/ +typedef union Euler { + struct { + float x, y, z; + }; + struct { + float pitch, yaw, roll; + }; +} Euler; + +/* +** 四元数，用来做旋转变换。在进行变换复合以及插值的时候用，但最终还是需要通过quat_mat4转换成矩阵和其他变换矩阵 +** 一起对向量进行变换 +*/ +typedef struct Quat { + float x, y, z, w; +} Quat; + +/* +** 4x4矩阵，列主项，用来做平移和缩放变换。之所以用列主序存储，是为了快速读取矩阵的基向量 +*/ +typedef union Mat4 { + float l[16]; + float m[4][4]; + struct { + float + e00, e10, e20, e30, /*colum 0*/ + e01, e11, e21, e31, + e02, e12, e22, e32, + e03, e13, e23, e33; + }; + struct { + Vec4 x;/*colum 0*/ + Vec4 y; + Vec4 z; + Vec4 w; + } axis; /*轴*/ + struct { + Vec4 x;/*colum 0*/ + Vec4 y; + Vec4 z; + Vec4 w; + } basis; /*基向量*/ + struct { + Vec4 axisx; + Vec4 axisy; + Vec4 axisz; + Vec4 pos; + }; + Vec4 colums[4]; +} Mat4; + +typedef union Mat3 { + struct { + float + e00, e10, e20, /*colum 0*/ + e01, e11, e21, + e02, e12, e22; + }; +} Mat3; + +typedef union Mat23 { + struct { + float + e00, e10, /*colum 0*/ + e01, e11, + e02, e12; + }; +} Mat23; + +typedef union Mat43 { + struct { + float + e00, e10, e20, e30, /*colum 0*/ + e01, e11, e21, e31, + e02, e12, e22, e32; + }; + struct { /*三个齐次裁剪坐标*/ + Vec4 p1; + Vec4 p2; + Vec4 p3; + }; +} Mat43; + +//#define MAT(m, r, c) (m->l[r + (c<<2)]) +#define MAT(M, r, c) (M->m[c][r]) +/************************************************************************/ +/* Vec */ +/************************************************************************/ + +void vec2_scale(Vec2* v, float k, Vec2* out); +void vec2_plus(Vec2* v1, Vec2* v2, Vec2* out); +void vec2_offset(Vec2* v, float offset, Vec2* out); +void vec2_rotate(Vec2* v, float angle, Vec2* out); + +float vec2_dot(Vec2* v1, Vec2* v2); + +void vec2_tostring(Vec2* v, char buf[]); +void vec2_print(Vec2* v); + +#define vec3_xy(v) (v->xy) + +extern Vec3 vec3forward; /*(0,0,1)*/ +extern Vec3 vec3up; /*(0,1,0)*/ +extern Vec3 vec3left;/*(1,0,0)*/ + +void vec3_tostring(Vec3* v, char buf[]); +void vec3_print(Vec3* v); + +float vec3_intersection(Vec3* v1, Vec3* v2); /*夹角*/ +void vec3_projection(Vec3* v1, Vec3* v2, Vec3* out);/*v1在v2上的投影*/ +void vec3_scale(Vec3* v, float k, Vec3* out); +void vec3_plus(Vec3* v1, Vec3* v2, Vec3* out); +void vec3_offset(Vec3* v, float offset, Vec3* out); +void vec3_normalize(Vec3* v, Vec3* out); +void vec3_vec4(float w, Vec4* out); + +void vec3_minus(Vec3* v1, Vec3* v2, Vec3* out); +float vec3_dot(Vec3* v1, Vec3* v2); +void vec3_cross(Vec3* v1, Vec3* v2, Vec3* out); +void vec3_multiply(Vec3* v1, Vec3* v2, Quat* out);// 向量的乘法，st=sxt-s*t，结果是一个四元数 + +float vec3_magnitude(Vec3* v1); +float vec3_magnitude2(Vec3* v1); + +void vec3_lerp(Vec3* v1, Vec3* v2, float t, Vec3* out); +void vec3_slerp(Vec3* v1, Vec3* v2, float t, Vec3* out); + +void vec4_dividew(Vec4* v, Vec3* out); + +void vec4_tostring(Vec4* v, char buf[]); +void vec4_print(Vec4* v); + +/************************************************************************/ +/* Matrix */ +/************************************************************************/ + +extern Mat4 mat4identity; + +void mat4_tostring(Mat4* m, char str[]); +void mat4_print(Mat4* m); + +void mat4_zero(Mat4* out); +void mat4_setidentity(Mat4* out); +void mat4_setfrustum(float l, float r, float b, float t, float n, float f, Mat4* out); +void mat4_setperspective(float fov, float aspect, float near, float far, Mat4* out); +void mat4_setscale(float kx, float ky, float kz, Mat4* out); +void mat4_setposition(float x, float y, float z, Mat4* out); +void mat4_setrotatez(float angle, Mat4* out); +void mat4_setrotatex(float angle, Mat4* out); +void mat4_setrotatey(float angle, Mat4* out); +void mat4_setrotate(float angleX, float angleY, float angleZ, Mat4* out);/*RyRxRz*/ +void mat4_setaxisangle(Vec3* axis, float angle, Mat4* out); + +bool mat4_setlookrotation(Vec3* view, Vec3* up, Mat4* out); +void mat4_setorthonormalbias(Vec3* x, Vec3* y, Vec3* z, Mat4* out); /*正交的三个轴*/ + +void mat4_orthogonalize(Mat4* in, Mat4* out); /*解决矩阵蠕变，对左上角3x3矩阵进行正交化，结果是右手系的正交矩阵*/ +bool mat4_isorthogonal(Mat4* m); /*判断是不是正交矩阵*/ +bool mat4_isidentity(Mat4* m); + +void mat4_settr(Vec3* pos, Quat* rot, Mat4* out); /*用旋转和平移初始化mat4*/ +void mat4_settrs(Vec3* pos, Quat* rot, Vec3* scale, Mat4* out); +void mat4_settrinverse(Vec3* pos, Quat* rot, Mat4* out); + +void mat4_multiply(Mat4* m1, Mat4* m2, Mat4* out); /* m1的行乘m2的列，意义是用m1变换m2 */ + +void mat4_transpose(Mat4* m, Mat4* out); + +void mat4_scale(Mat4* m, Vec3* scale, Mat4* out);/* 后乘post-multiply scale */ +void mat4_translate(Mat4* m, Vec3* pos, Mat4* out); /* 后乘post-multiply translate */ +void mat4_rotate(Mat4*m, float angle, Vec3* rot, Mat4* out);/*后乘绕任意轴向量旋转矩阵*/ + +bool mat4_invertfull(Mat4* in, Mat4* out); /* 并不是所有矩阵都能求逆 */ +bool mat4_invertgeneral3d(Mat4* in, Mat4* out); /* 对scale rotate translate求逆 */ +void mat4_invertscale(Mat4* scale, Mat4* out); /* 对缩放矩阵求逆 */ +void mat4_invertrot(Mat4* rot, Mat4* out); /* 对旋转矩阵求逆 */ +void mat4_invertpos(Mat4* pos, Mat4* out); /* 对平移矩阵求逆 */ + +void mat4_decomposetrs(Mat4* src, Vec3* pos, Quat* quat, Vec3* scale); /*分解trs矩阵*/ + +void mat4_applytovec4(Mat4* m, Vec4* v, Vec4* out); + +bool mat4_toeuler(Mat4* in, Euler* out); /* 计算YXZ旋转矩阵的欧拉角 */ +void mat4_toquat(Mat4* in, Quat* out); /*in是正交矩阵*/ + +#define ROWMAT(A, ...)\ +Mat4 A={__VA_ARGS__};mat4_transpose(&A, &A); + +void mat3_applytovec3(Mat3* m, Vec3* v, Vec3* out); +void mat23_applytovec3(Mat23* m, Vec3* v, Vec2* out); +void mat43_applytovec3(Mat43* m, Vec3* v, Vec4* out); + +/************************************************************************/ +/* Quat */ +/************************************************************************/ + +void quat_tostring(Quat* q, char str[]); +void quat_print(Quat* q); + +void euler_toquat(Euler* e, Quat* out); +void euler_deg2rad(Euler* in, Euler* out); +void euler_rad2deg(Euler* in, Euler* out); + +void euler_tostring(Euler* v, char buf[]); +void euler_print(Euler* v); + +void quat_fromaxisangle(Vec3* axis, float angle, Quat* out); /*轴角转四元数*/ +void quat_fromeuler(Euler* euler, Quat* out); /*按照zxy顺序*/ + +void quat_tomat4(Quat* q, Mat4* out); +void quat_toeuler(Quat*q, Euler* out); + +void quat_normalize(Quat* q, Quat* out); /*解决蠕变，保持四元数合法*/ + +void quat_scale(Quat* q, float scale, Quat* out); +void quat_rotate(); + +void quat_minus(Quat* q1, Quat* q2, Quat* out); +void quat_slerp(Quat* start, Quat* end, float t, Quat* out); +void quat_lerp(Quat* start, Quat* end, float t, Quat* out); +void quat_translate(Quat* q, Vec4* v, Vec4* out); +void quat_invert(Quat* q, Quat* out); +float quat_dot(Quat* q1, Quat* q2); +void quat_multiply(Quat* q1, Quat* q2, Quat* out); +void quat_devide(Quat* q, float k, Quat* out); +void quat_negtive(Quat* in, Quat* out); +bool quat_isidentity(Quat* q); + +void quat_applytovec3(Quat* q, Vec3* v, Vec3* out); /*用四元数直接旋转向量*/ + +void quat_conjugate(Quat* in, Quat* out); + +bool quat_setlookrotation(Vec3* view, Vec3* up, Quat* out); + +float quat_magnitude(Quat* q); +float quat_magnitude2(Quat* q); + +#endif
\ No newline at end of file diff --git a/src/math/matrix.c b/src/math/matrix.c new file mode 100644 index 0000000..458f432 --- /dev/null +++ b/src/math/matrix.c @@ -0,0 +1,741 @@ +#include <math.h> +#include <stdio.h> +#include <string.h> + +#include "math.h" +#include "../util/assert.h" +#include "../core/mem.h" + + +static Mat4 sharedMat; +static Mat4 sharedMat2; +static Vec4 sharedVec4; + +Mat4 mat4identity = { + 1,0,0,0, + 0,1,0,0, + 0,0,1,0, + 0,0,0,1 +}; + +#define shrmat(p) \ +do{\ +sharedMat = *p;\ +p = &sharedMat;\ +}while(0) + +#define shrmat2(p) \ +do{\ +sharedMat2 = *p;\ +p = &sharedMat2;\ +}while(0) + +void mat4_tostring(Mat4* m, char str[]) { + ssrM_zero(str, sizeof(str)); + for (int r = 0; r < 4; ++r) { + for (int c = 0; c < 4; ++c) { + sprintf(str, "%8.3f ", MAT(m, r, c) == -0 ? +0 : MAT(m, r, c)); + str += strlen(str); + } + if(r != 3) sprintf(str, "\n"); + str += strlen(str); + } +} + +void mat4_print(Mat4* m) { + mat4_tostring(m, printbuffer); + printf("\n%s\n", printbuffer); +} + +void mat4_zero(Mat4* out) { + ssr_assert(out); + ssrM_zero(out, sizeof(Mat4)); +} + +void mat4_setidentity(Mat4* out) { + ssr_assert(out); + mat4_zero(out); + out->e00 = 1; + out->e11 = 1; + out->e22 = 1; + out->e33 = 1; +} + +void mat4_setfrustum(float l, float r, float b, float t, float n, float f, Mat4* out) { + ssr_assert(out); + mat4_zero(out); + out->e00 = (2.f * n) / (r - l); + out->e02 = (r + l) / (r - l); + out->e11 = 2.f * n / (t - b); + out->e12 = (t + b) / (t - b); + out->e22 = -(f + n) / (f - n); + out->e23 = -2.f * f * n / (f - n); + out->e32 = -1; +} + +void mat4_setperspective(float _fov, float aspect, float near, float far, Mat4* out) { + float fov = _fov * PI / 180.f; + float tanf = tan(fov * 0.5); + mat4_setfrustum( + -near*tanf*aspect, + near*tanf*aspect, + -near*tanf, + near*tanf, + near, + far, + out + ); +} + +static float _mul(float* r, float* c) { + return c[0] * r[0] + c[1] * r[4] + c[2] * r[8] + c[3] * r[12]; +} + +#define mul(r, c) _mul(&MAT(m1,r,0), &MAT(m2,0,c)) + +void mat4_multiply(Mat4* m1, Mat4* m2, Mat4* out) { + ssr_assert(m1 && m2 && out); + if (mat4_isidentity(m1)) { if(m2 != out) *out = *m2; return; } + if (mat4_isidentity(m2)) { if(m1 != out) *out = *m1; return; } + if (m1 == out) shrmat(m1); + if (m2 == out) shrmat2(m2); + + out->e00 = mul(0, 0); out->e01 = mul(0, 1); out->e02 = mul(0, 2); out->e03 = mul(0, 3); + out->e10 = mul(1, 0); out->e11 = mul(1, 1); out->e12 = mul(1, 2); out->e13 = mul(1, 3); + out->e20 = mul(2, 0); out->e21 = mul(2, 1); out->e22 = mul(2, 2); out->e23 = mul(2, 3); + out->e30 = mul(3, 0); out->e31 = mul(3, 1); out->e32 = mul(3, 2); out->e33 = mul(3, 3); +} + +void mat4_setscale(float kx, float ky, float kz, Mat4* out) { + ssr_assert(out); + mat4_zero(out); + out->e00 = kx; + out->e11 = ky; + out->e22 = kz; + out->e33 = 1; +} + +void mat4_setposition(float x, float y, float z, Mat4* out) { + ssr_assert(out); + mat4_setidentity(out); + out->e03 = x; + out->e13 = y; + out->e23 = z; +} + +void mat4_setrotatez(float angle, Mat4* out) { + ssr_assert(out); + mat4_setidentity(out); + angle = radians(angle); + float s = sin(angle), c = cos(angle); + out->e00 = c; out->e01 = -s; + out->e10 = s; out->e11 = c; +} + +void mat4_setrotatex(float angle, Mat4* out) { + ssr_assert(out); + mat4_setidentity(out); + angle = radians(angle); + float s = sin(angle), c = cos(angle); + out->e11 = c; out->e12 = -s; + out->e21 = s; out->e22 = c; +} + +void mat4_setrotatey(float angle, Mat4* out) { + ssr_assert(out); + mat4_setidentity(out); + angle = radians(angle); + float s = sin(angle), c = cos(angle); + out->e00 = c; out->e02 = s; + out->e20 = -s; out->e22 = c; +} + +/*https://www.geometrictools.com/Documentation/EulerAngles.pdf*/ +void mat4_setrotate(float angleX, float angleY, float angleZ, Mat4* out) { + ssr_assert(out); + mat4_setidentity(out); + angleX = radians(angleX); angleY = radians(angleY); angleZ = radians(angleZ); + float sx = sin(angleX), cx = cos(angleX); + float sy = sin(angleY), cy = cos(angleY); + float sz = sin(angleZ), cz = cos(angleZ); + out->e00 = cy * cz + sx * sy * sz; out->e01 = cz * sx*sy - cy * sz; out->e02 = cx * sy; + out->e10 = cx * sz; out->e11 = cx * cz; out->e12 = -sx; + out->e20 = -cz * sy + cy * sx * sz; out->e21 = cy * cz*sx + sy * sz; out->e22 = cx * cy; +} + +void mat4_setaxisangle(Vec3* ax, float angle, Mat4* out) { + ssr_assert(ax && out); + + float a = radians(angle); + float c = cos(a); + float s = sin(a); + + Vec3 axis = *ax; + Vec3 temp; + vec3_normalize(&axis, &axis); + vec3_scale(&axis, 1 - c, &temp); + + /* + rotation matrix 推导过程 https://zhuanlan.zhihu.com/p/56587491 + X^2(1-c)+c, XY(1-c)-Zs, XZ(1-c)+Ys, 0 + XY(1-c)+Zs, Y^2(1-c)+c, YZ(1-c)-Xs, 0 + XZ(1-c)-Ys, YZ(1-c)+Xs, Z^2(1-c)+c, 0 + 0, 0, 0, 1 + */ + + mat4_setidentity(out); + out->m[0][0] = c + temp.x * axis.x; + out->m[0][1] = 0 + temp.x * axis.y + s * axis.z; + out->m[0][2] = 0 + temp.x * axis.z - s * axis.y; + + out->m[1][0] = 0 + temp.y * axis.x - s * axis.z; + out->m[1][1] = c + temp.y * axis.y; + out->m[1][2] = 0 + temp.y * axis.z + s * axis.x; + + out->m[2][0] = 0 + temp.z * axis.x + s * axis.y; + out->m[2][1] = 0 + temp.z * axis.y - s * axis.x; + out->m[2][2] = c + temp.z * axis.z; + +} + +void mat4_setorthonormalbias(Vec3* x, Vec3* y, Vec3* z, Mat4* out) { + ssr_assert(x && y && z); + mat4_setidentity(out); + Vec4 asix = { x->x, x->y, x->z, 0 }; + Vec4 asiy = { y->x, y->y, y->z, 0 }; + Vec4 asiz = { z->x, z->y, z->z, 0 }; + out->colums[0] = asix; + out->colums[1] = asiy; + out->colums[2] = asiz; +} + +bool mat4_isidentity(Mat4* m) { + ssr_assert(m); + //return memcmp(m, &mat4identity, sizeof(Mat4)) == 0; + return + compare(m->axisx.x, 1) && compare(m->axisx.y, 0) && compare(m->axisx.z,0) && compare(m->axisx.w, 0) && + compare(m->axisy.x, 0) && compare(m->axisy.y, 1) && compare(m->axisy.z,0) && compare(m->axisy.w, 0) && + compare(m->axisz.x, 0) && compare(m->axisz.y, 0) && compare(m->axisz.z,1) && compare(m->axisz.w, 0) && + compare(m->pos.x, 0 ) && compare(m->pos.y, 0 ) && compare(m->pos.z, 0 ) &&compare( m->pos.w, 1); +} + +bool mat4_isorthogonal(Mat4* m) { + ssr_assert(m); + Mat4 trans = {0}, res = { 0 }; + mat4_transpose(m, &trans); + mat4_multiply(m, &trans, &res); + return mat4_isidentity(&res); +} + +/* +** 以z轴为准进行正交化，分为施密特正交化和叉乘正交化，施密特过程更加普遍，叉乘适用于三维空间，两种方法实际上等价 +** 如果用叉乘的方法，只需要关注yz，x通过叉乘得到 +*/ +void mat4_orthogonalize(Mat4* in, Mat4* out) { + ssr_assert(in && out); + if (in == out) { + shrmat(in); + } + + mat4_setidentity(out); + Vec4 z = in->basis.z; + vec3_normalize(&z, &z); + Vec4 y = in->basis.y; + Vec4 x = {0}; + vec3_cross(&y, &z, &x); + vec3_normalize(&x, &x); + vec3_cross(&z, &x, &y); + out->basis.x = x; + out->basis.y = y; + out->basis.z = z; + + /* + mat4_setidentity(out); + + Vec4 x = in->basis.x; + Vec4 y = in->basis.y; + Vec4 z = in->basis.z; + Vec3 temp, temp2; + + vec3_normalize(&z, &z); + out->basis.z = z; + + float dot = vec3_dot(&y, &z); + vec3_scale(&z, dot, &temp); + vec3_minus(&y, &temp, &y); + vec3_normalize(&y, &y); + out->basis.y = y; + + vec3_cross(&y, &z, &out->basis.x); + */ + /*针对右手系调整basis.x的方向*/ + /*https://math.stackexchange.com/questions/1847465/why-to-use-gram-schmidt-process-to-orthonormalise-a-basis-instead-of-cross-produ*/ + /*由于需要针对右手系，这里不这样计算，因为可能要对结果进行翻转 + dot = vec3_dot(&x, &z); + vec3_scale(&z, dot, &temp); + vec3_minus(&x, &temp, &temp2); + dot = vec3_dot(&x, &y); + vec3_scale(&y, dot, &temp); + vec3_minus(&temp2, &temp, &x); + vec3_normalize(&x, &x); + out->basis.x = x; + */ +} + +bool mat4_setlookrotation(Vec3* view, Vec3* up, Mat4* out) { + ssr_assert(view && up && out); + + /*正交化*/ + float mag = vec3_magnitude(view); + if (mag < EPSILON) return 0; + Vec3 z; + vec3_scale(view, 1.f / mag, &z); + + Vec3 x; + vec3_cross(up, &z, &x); + mag = vec3_magnitude(&x); + if (mag < EPSILON) return 0; + vec3_scale(&x, 1.f / mag, &x); + + Vec3 y; + vec3_cross(&z, &x, &y); + mag = vec3_magnitude(&y); + if (!compare(mag, 1)) return 0; + + mat4_setorthonormalbias(&x, &y, &z, out); /*xyz正交*/ + + return 1; +} + +void mat4_applytovec4(Mat4* mat, Vec4* v, Vec4* out) { + ssr_assert(mat && v && out); + if (v == out) { + sharedVec4 = *v; + v = &sharedVec4; + } + out->x = mat->e00 * v->x + mat->e01 * v->y + mat->e02 * v->z + mat->e03 * v->w; + out->y = mat->e10 * v->x + mat->e11 * v->y + mat->e12 * v->z + mat->e13 * v->w; + out->z = mat->e20 * v->x + mat->e21 * v->y + mat->e22 * v->z + mat->e23 * v->w; + out->w = mat->e30 * v->x + mat->e31 * v->y + mat->e32 * v->z + mat->e33 * v->w; +} + +#define trans(r, c) out->e##r##c = m->e##c##r + +void mat4_transpose(Mat4* m, Mat4* out) { + ssr_assert(m && out); + if (m == out) shrmat(m); + + trans(0, 0); trans(0, 1); trans(0, 2); trans(0, 3); + trans(1, 0); trans(1, 1); trans(1, 2); trans(1, 3); + trans(2, 0); trans(2, 1); trans(2, 2); trans(2, 3); + trans(3, 0); trans(3, 1); trans(3, 2); trans(3, 3); +} + +/* +** 使用高斯消元法计算任意矩阵的逆矩阵。针对不含投影的3D变换矩阵，应该使用 +** mat4_invertgeneral3d() +** 更快一些 +*/ +bool mat4_invertfull(Mat4* m, Mat4* out) { + ssr_assert(m && out); + +#define _m(r, c) MAT(m, r, c) + float wtmp[4][8] = { + { /*M*/ _m(0,0), _m(0, 1), _m(0, 2), _m(0, 3), /*I*/ 1, 0, 0, 0 }, + { /*M*/ _m(1,0), _m(1, 1), _m(1, 2), _m(1, 3), /*I*/ 0, 1, 0, 0 }, + { /*M*/ _m(2,0), _m(2, 1), _m(2, 2), _m(2, 3), /*I*/ 0, 0, 1, 0 }, + { /*M*/ _m(3,0), _m(3, 1), _m(3, 2), _m(3, 3), /*I*/ 0, 0, 0, 1 }, + }; +#undef _m + float m0, m1, m2, m3, s; + float *r0, *r1, *r2, *r3; + r0 = wtmp[0], r1 = wtmp[1], r2 = wtmp[2], r3 = wtmp[3]; +#define SWAP_ROWS(a, b) { float *_tmp = a; (a)=(b); (b)=_tmp; } + + //#define optimize(block) if(s!=0.f){block} +#define optimize(block) block + + /* choose pivot - or die */ + if (absf(r3[0]) > absf(r2[0])) SWAP_ROWS(r3, r2); + if (absf(r2[0]) > absf(r1[0])) SWAP_ROWS(r2, r1); + if (absf(r1[0]) > absf(r0[0])) SWAP_ROWS(r1, r0); + if (0.0f == r0[0]) return 0; + + /* eliminate first variable */ + m1 = r1[0] / r0[0]; m2 = r2[0] / r0[0]; m3 = r3[0] / r0[0]; + s = r0[1]; optimize(r1[1] -= m1 * s; r2[1] -= m2 * s; r3[1] -= m3 * s; ) + s = r0[2]; optimize(r1[2] -= m1 * s; r2[2] -= m2 * s; r3[2] -= m3 * s; ) + s = r0[3]; optimize(r1[3] -= m1 * s; r2[3] -= m2 * s; r3[3] -= m3 * s; ) + s = r0[4]; optimize(r1[4] -= m1 * s; r2[4] -= m2 * s; r3[4] -= m3 * s; ) + s = r0[5]; optimize(r1[5] -= m1 * s; r2[5] -= m2 * s; r3[5] -= m3 * s; ) + s = r0[6]; optimize(r1[6] -= m1 * s; r2[6] -= m2 * s; r3[6] -= m3 * s; ) + s = r0[7]; optimize(r1[7] -= m1 * s; r2[7] -= m2 * s; r3[7] -= m3 * s; ) + + /* choose pivot - or die */ + if (absf(r3[1]) > absf(r2[1])) SWAP_ROWS(r3, r2); + if (absf(r2[1]) > absf(r1[1])) SWAP_ROWS(r2, r1); + if (0.0F == r1[1]) return 0; + + /* eliminate second variable */ + m2 = r2[1] / r1[1]; m3 = r3[1] / r1[1]; + s = r1[2]; optimize(r2[2] -= m2 * s; r3[2] -= m3 * s; ) + s = r1[3]; optimize(r2[3] -= m2 * s; r3[3] -= m3 * s; ) + s = r1[4]; optimize(r2[4] -= m2 * s; r3[4] -= m3 * s; ) + s = r1[5]; optimize(r2[5] -= m2 * s; r3[5] -= m3 * s; ) + s = r1[6]; optimize(r2[6] -= m2 * s; r3[6] -= m3 * s; ) + s = r1[7]; optimize(r2[7] -= m2 * s; r3[7] -= m3 * s; ) + + /* choose pivot - or die */ + if (absf(r3[2])>absf(r2[2])) SWAP_ROWS(r3, r2); + if (0.0F == r2[2]) return 0; + + /* eliminate third variable */ + m3 = r3[2] / r2[2]; + s = r2[3]; optimize(r3[3] -= m3 * s; ) + s = r2[4]; optimize(r3[4] -= m3 * s; ) + s = r2[5]; optimize(r3[5] -= m3 * s; ) + s = r2[6]; optimize(r3[6] -= m3 * s; ) + s = r2[7]; optimize(r3[7] -= m3 * s; ) + +#undef optimize + + /* last check */ + if (0.0F == r3[3]) return 0; + + s = 1.0F / r3[3]; /* now back substitute row 3 */ + r3[4] *= s; r3[5] *= s; r3[6] *= s; r3[7] *= s; + + m2 = r2[3]; /* now back substitute row 2 */ + s = 1.0F / r2[2]; + r2[4] = s * (r2[4] - r3[4] * m2), r2[5] = s * (r2[5] - r3[5] * m2), + r2[6] = s * (r2[6] - r3[6] * m2), r2[7] = s * (r2[7] - r3[7] * m2); + m1 = r1[3]; + r1[4] -= r3[4] * m1, r1[5] -= r3[5] * m1, + r1[6] -= r3[6] * m1, r1[7] -= r3[7] * m1; + m0 = r0[3]; + r0[4] -= r3[4] * m0, r0[5] -= r3[5] * m0, + r0[6] -= r3[6] * m0, r0[7] -= r3[7] * m0; + + m1 = r1[2]; /* now back substitute row 1 */ + s = 1.0F / r1[1]; + r1[4] = s * (r1[4] - r2[4] * m1), r1[5] = s * (r1[5] - r2[5] * m1), + r1[6] = s * (r1[6] - r2[6] * m1), r1[7] = s * (r1[7] - r2[7] * m1); + m0 = r0[2]; + r0[4] -= r2[4] * m0, r0[5] -= r2[5] * m0, + r0[6] -= r2[6] * m0, r0[7] -= r2[7] * m0; + + m0 = r0[1]; /* now back substitute row 0 */ + s = 1.0F / r0[0]; + r0[4] = s * (r0[4] - r1[4] * m0), r0[5] = s * (r0[5] - r1[5] * m0), + r0[6] = s * (r0[6] - r1[6] * m0), r0[7] = s * (r0[7] - r1[7] * m0); + + out->e00 = r0[4]; out->e01 = r0[5]; out->e02 = r0[6]; out->e03 = r0[7]; + out->e10 = r1[4]; out->e11 = r1[5]; out->e12 = r1[6]; out->e13 = r1[7]; + out->e20 = r2[4]; out->e21 = r2[5]; out->e22 = r2[6]; out->e23 = r2[7]; + out->e30 = r3[4]; out->e31 = r3[5]; out->e32 = r3[6]; out->e33 = r3[7]; + +#undef SWAP_ROWS + + return 1; +} + +/* +** 对只包含基本3D变换的矩阵进行变换，先计算左上角3x3的RS矩阵的逆（通过伴随矩阵），然后 +** 乘上平移矩阵的逆矩阵，即 +** M^-1 = (T(RS))^-1 = (RS)^-1 * T^-1 +*/ +bool mat4_invertgeneral3d(Mat4* in, Mat4* out) { + ssr_assert(in && out); + if (in == out) shrmat(in); + + mat4_setidentity(out); + + /*计算左上角3x3矩阵的行列式*/ + float pos = 0, neg = 0, t; + float det; + + t = in->e00 * in->e11 * in->e22; + if (t >= 0) pos += t; else neg += t; + t = in->e10 * in->e21 * in->e02; + if (t >= 0) pos += t; else neg += t; + t = in->e20 * in->e01 * in->e12; + if (t >= 0) pos += t; else neg += t; + + t = -in->e20 * in->e11 * in->e02; + if (t >= 0) pos += t; else neg += t; + t = -in->e10 * in->e01 * in->e22; + if (t >= 0) pos += t; else neg += t; + t = -in->e00 * in->e21 * in->e12; + if (t >= 0) pos += t; else neg += t; + + det = pos + neg; + + if (det * det < 1e-25) + return 0; /*行列式为0*/ + + det = 1.f / det; + MAT(out, 0, 0) = ((MAT(in, 1, 1)*MAT(in, 2, 2) - MAT(in, 2, 1)*MAT(in, 1, 2))*det); + MAT(out, 0, 1) = (-(MAT(in, 0, 1)*MAT(in, 2, 2) - MAT(in, 2, 1)*MAT(in, 0, 2))*det); + MAT(out, 0, 2) = ((MAT(in, 0, 1)*MAT(in, 1, 2) - MAT(in, 1, 1)*MAT(in, 0, 2))*det); + MAT(out, 1, 0) = (-(MAT(in, 1, 0)*MAT(in, 2, 2) - MAT(in, 2, 0)*MAT(in, 1, 2))*det); + MAT(out, 1, 1) = ((MAT(in, 0, 0)*MAT(in, 2, 2) - MAT(in, 2, 0)*MAT(in, 0, 2))*det); + MAT(out, 1, 2) = (-(MAT(in, 0, 0)*MAT(in, 1, 2) - MAT(in, 1, 0)*MAT(in, 0, 2))*det); + MAT(out, 2, 0) = ((MAT(in, 1, 0)*MAT(in, 2, 1) - MAT(in, 2, 0)*MAT(in, 1, 1))*det); + MAT(out, 2, 1) = (-(MAT(in, 0, 0)*MAT(in, 2, 1) - MAT(in, 2, 0)*MAT(in, 0, 1))*det); + MAT(out, 2, 2) = ((MAT(in, 0, 0)*MAT(in, 1, 1) - MAT(in, 1, 0)*MAT(in, 0, 1))*det); + + // 乘T^-1 + MAT(out, 0, 3) = -(MAT(in, 0, 3) * MAT(out, 0, 0) + + MAT(in, 1, 3) * MAT(out, 0, 1) + + MAT(in, 2, 3) * MAT(out, 0, 2)); + MAT(out, 1, 3) = -(MAT(in, 0, 3) * MAT(out, 1, 0) + + MAT(in, 1, 3) * MAT(out, 1, 1) + + MAT(in, 2, 3) * MAT(out, 1, 2)); + MAT(out, 2, 3) = -(MAT(in, 0, 3) * MAT(out, 2, 0) + + MAT(in, 1, 3) * MAT(out, 2, 1) + + MAT(in, 2, 3) * MAT(out, 2, 2)); + + return 1; +} + +void mat4_invertpos(Mat4* in, Mat4* out) { + +} + +void mat4_invertscale(Mat4* in, Mat4* out) { + +} + +void mat4_invertrot(Mat4* in, Mat4* out) { + ssr_assert(in && out); + mat4_transpose(in, out); +} + +void mat4_settr(Vec3* pos, Quat* rot, Mat4* out) { + ssr_assert(pos && rot && out); + mat4_zero(out); + quat_tomat4(rot, out); + out->e03 = pos->x; + out->e13 = pos->y; + out->e23 = pos->z; +} + +void mat4_settrs(Vec3* pos, Quat* rot, Vec3* scale, Mat4* out) { + ssr_assert(pos && rot && scale && out); + mat4_zero(out); + quat_tomat4(rot, out); /*pos*rot*scale的顺序*/ + out->e00 *= scale->x; out->e01 *= scale->y; out->e02 *= scale->z; + out->e10 *= scale->x; out->e11 *= scale->y; out->e12 *= scale->z; + out->e20 *= scale->x; out->e21 *= scale->y; out->e22 *= scale->z; + out->e03 = pos->x; + out->e13 = pos->y; + out->e23 = pos->z; +} + +void mat4_settrinverse(Vec3* pos, Quat* rot, Mat4* out) { + ssr_assert(pos && rot && out); + mat4_zero(out); + quat_invert(rot, rot); + quat_tomat4(rot, out); + Vec3 reverse = { -pos->x, -pos->y, -pos->z}; + mat4_translate(out, &reverse, out); /* (TR)^-1 = R^-1*T^-1所以这里是右乘*/ +} + +void mat4_scale(Mat4* m, Vec3* scale, Mat4* out) { + ssr_assert(m && scale && out); + if (out != m) { + *out = *m; + } + /* + scale matrix + x, 0, 0, 0, + 0, y, 0, 0, + 0, 0, z, 0, + 0, 0, 0, 1 + */ + out->e00 *= scale->x; + out->e10 *= scale->x; + out->e20 *= scale->x; + out->e30 *= scale->x; + + out->e01 *= scale->y; + out->e11 *= scale->y; + out->e21 *= scale->y; + out->e31 *= scale->y; + + out->e02 *= scale->z; + out->e12 *= scale->z; + out->e22 *= scale->z; + out->e32 *= scale->z; +} + +void mat4_translate(Mat4* m, Vec3* pos, Mat4* out) { + ssr_assert(m && pos && out); + if (out != m) { + *out = *m; + } + /* + translate matrix + 1, 0, 0, x, + 0, 1, 0, y, + 0, 0, 1, z, + 0, 0, 0, 1, + */ + out->e03 = out->e00 * pos->x + out->e01 * pos->y + out->e02 * pos->z + out->e03; + out->e13 = out->e10 * pos->x + out->e11 * pos->y + out->e12 * pos->z + out->e13; + out->e23 = out->e20 * pos->x + out->e21 * pos->y + out->e22 * pos->z + out->e23; + out->e33 = out->e30 * pos->x + out->e31 * pos->y + out->e32 * pos->z + out->e33; +} + +void mat4_rotate(Mat4* m, float angle, Vec3* ax, Mat4* out) { + ssr_assert(m && ax && out); + Mat4 rot; + mat4_setaxisangle(ax, angle, &rot); + mat4_multiply(m, &rot, out); +} + +void mat4_decomposetrs(Mat4* src, Vec3* pos, Quat* quat, Vec3* scale) { + ssr_assert(src && pos && quat && scale); + + Vec3* x = &src->colums[0]; + Vec3* y = &src->colums[1]; + Vec3* z = &src->colums[2]; + Vec3* w = &src->colums[3]; + + *pos = *w; + + quat_setlookrotation(z, y, quat); + + scale->x = vec3_magnitude(x); + scale->y = vec3_magnitude(y); + scale->z = vec3_magnitude(z); +} + +static void MakePositive(Euler* euler) {/*弧度制欧拉角*/ + const float negativeFlip = -0.0001F; + const float positiveFlip = (PI * 2.0F) - 0.0001F; + + if (euler->x < negativeFlip) + euler->x += 2.0 * PI; + else if (euler->x > positiveFlip) + euler->x -= 2.0 * PI; + + if (euler->y < negativeFlip) + euler->y += 2.0 * PI; + else if (euler->y > positiveFlip) + euler->y -= 2.0 * PI; + + if (euler->z < negativeFlip) + euler->z += 2.0 * PI; + else if (euler->z > positiveFlip) + euler->z -= 2.0 * PI; +} + +static void SanitizeEuler(Euler* e) {/*弧度制欧拉角*/ + MakePositive(e); +} + +/*from unity src*/ +bool mat4_toeuler(Mat4* in, Euler* out) { + ssr_assert(in && out); + // from http://www.geometrictools.com/Documentation/EulerAngles.pdf + // YXZ order + if (MAT(in, 1, 2) < 0.999F) // some fudge for imprecision + { + if (MAT(in, 1, 2) > -0.999F) // some fudge for imprecision + { + out->x = asin(-MAT(in, 1, 2)); + out->y = atan2(MAT(in, 0, 2), MAT(in, 2, 2)); + out->z = atan2(MAT(in, 1, 0), MAT(in, 1, 1)); + //euler_rad2deg(out, out); + SanitizeEuler(out); + euler_rad2deg(out, out); + return 1; + } + else + { + // WARNING. Not unique. YA - ZA = atan2(r01,r00) + out->x = PI * 0.5F; + out->y = atan2(MAT(in, 0, 1), MAT(in, 0, 0)); + out->z = 0.0F; + //euler_rad2deg(out, out); + SanitizeEuler(out); + euler_rad2deg(out, out); + return 0; + } + } + else + { + // WARNING. Not unique. YA + ZA = atan2(-r01,r00) + out->x = -PI * 0.5F; + out->y = atan2(-MAT(in, 0, 1), MAT(in, 0, 0)); + out->z = 0.0F; + //euler_rad2deg(out, out); + SanitizeEuler(out); + euler_rad2deg(out, out); + return 0; + } +} + +/*from unity src*/ +/*https://www.euclideanspace.com/maths/geometry/rotations/conversions/matrixToQuaternion/*/ +void mat4_toquat(Mat4* in, Quat* out) { + // Algorithm in Ken Shoemake's article in 1987 SIGGRAPH course notes + // article "Quaternionf Calculus and Fast Animation". + float fTrace = MAT(in, 0, 0) + MAT(in, 1, 1) + MAT(in, 2, 2); + float fRoot; + + if (fTrace > 0.0f) + { + // |w| > 1/2, may as well choose w > 1/2 + fRoot = sqrt(fTrace + 1.0f); // 2w + out->w = 0.5f*fRoot; + fRoot = 0.5f / fRoot; // 1/(4w) + out->x = (MAT(in, 2, 1) - MAT(in, 1, 2))*fRoot; + out->y = (MAT(in, 0, 2) - MAT(in, 2, 0))*fRoot; + out->z = (MAT(in, 1, 0) - MAT(in, 0, 1))*fRoot; + } + else + { + // |w| <= 1/2 + int s_iNext[3] = { 1, 2, 0 }; + int i = 0; + if (MAT(in, 1, 1) > MAT(in, 0, 0)) + i = 1; + if (MAT(in, 2, 2) > MAT(in, i, i)) + i = 2; + int j = s_iNext[i]; + int k = s_iNext[j]; + + fRoot = sqrt(MAT(in, i, i) - MAT(in, j, j) - MAT(in, k, k) + 1.0f); + float* apkQuat[3] = { &out->x, &out->y, &out->z }; + ssr_assert(fRoot >= EPSILON); + *apkQuat[i] = 0.5f*fRoot; + fRoot = 0.5f / fRoot; + out->w = (MAT(in, k, j) - MAT(in, j, k)) * fRoot; + *apkQuat[j] = (MAT(in, j, i) + MAT(in, i, j))*fRoot; + *apkQuat[k] = (MAT(in, k, i) + MAT(in, i, k))*fRoot; + } + quat_normalize(out, out); +} + +void mat3_applytovec3(Mat3* m, Vec3* v, Vec3* out) { + ssr_assert(m && v && out); + out->x = m->e00 * v->x + m->e01 * v->y + m->e02 * v->z; + out->y = m->e10 * v->x + m->e11 * v->y + m->e12 * v->z; + out->z = m->e20 * v->x + m->e21 * v->y + m->e22 * v->z; +} + +void mat23_applytovec3(Mat23* m, Vec3* v, Vec2* out) { + ssr_assert(m && v && out); + out->x = m->e00 * v->x + m->e01 * v->y + m->e02 * v->z; + out->y = m->e10 * v->x + m->e11 * v->y + m->e12 * v->z; +} + +void mat43_applytovec3(Mat43* m, Vec3* v, Vec4* out) { + ssr_assert(m && v && out); + out->x = m->e00 * v->x + m->e01 * v->y + m->e02 * v->z; + out->y = m->e10 * v->x + m->e11 * v->y + m->e12 * v->z; + out->z = m->e20 * v->x + m->e21 * v->y + m->e22 * v->z; + out->w = m->e30 * v->x + m->e31 * v->y + m->e32 * v->z; +}
\ No newline at end of file diff --git a/src/math/quat.c b/src/math/quat.c new file mode 100644 index 0000000..e68e254 --- /dev/null +++ b/src/math/quat.c @@ -0,0 +1,314 @@ +#include <math.h> +#include <stdio.h> +#include <string.h> + +#include "math.h" +#include "../util/assert.h" +#include "../core/mem.h" + +static Vec3 sharedVec3; +static Quat sharedQuat; +static Quat sharedQuat2; +static Euler sharedEuler; + +#define shrquat(q) \ +do{\ +sharedQuat = *q;\ +q = &sharedQuat;\ +}while(0) + +#define shrquat2(q) \ +do{\ +sharedQuat2 = *q;\ +q = &sharedQuat2;\ +}while(0) + +#define shrvec3(v) \ +do{\ +sharedVec3 = *v;\ +v = &sharedVec3;\ +}while(0) + +#define shreuler(v) \ +do{\ +sharedEuler = *v;\ +v = &sharedEuler;\ +}while(0) + +void quat_tostring(Quat* v, char buf[]) { + ssr_assert(v); + sprintf(buf, "%8.3f %8.3f %8.3f", v->x, v->y, v->z); +} + +void quat_print(Quat* q) { + ssr_assert(q); + quat_tostring(q, printbuffer); + printf("\n%s\n", printbuffer); +} + +void euler_tostring(Euler* v, char buf[]) { + sprintf(buf, "%8.3f %8.3f %8.3f", v->x, v->y, v->z); +} + +void euler_print(Euler* v) { + euler_tostring(v, printbuffer); + printf("\n%s\n", printbuffer); +} + +void euler_deg2rad(Euler* in, Euler* out) { + ssr_assert(in && out); + out->x = radians(in->x); + out->y = radians(in->y); + out->z = radians(in->z); +} + +void euler_rad2deg(Euler* in, Euler* out) { + ssr_assert(in && out); + out->x = degree(in->x); + out->y = degree(in->y); + out->z = degree(in->z); +} + +void euler_toquat(Euler* euler, Quat* out) { + ssr_assert(euler && out); + quat_fromeuler(euler, out); +} + +bool quat_isidentity(Quat* q) { + return compare(quat_magnitude(q), 1.f); +} + +void quat_fromaxisangle(Vec3* axis, float angle, Quat* out) { + ssr_assert(compare(vec3_magnitude2(axis), 1.f)); + angle = radians(angle); + float half = angle * 0.5; + float s = sin(half); + out->w = cos(half); + out->x = s * axis->x; + out->y = s * axis->y; + out->z = s * axis->z; +} + +void quat_fromeuler(Euler* el, Quat* out) { + ssr_assert(el && out); + Euler euler; + euler_deg2rad(el, &euler); + float cx = cos(euler.x * 0.5f); + float sx = sin(euler.x * 0.5f); + float cy = cos(euler.y * 0.5f); + float sy = sin(euler.y * 0.5f); + float cz = cos(euler.z * 0.5f); + float sz = sin(euler.z * 0.5f); + Quat qx = {sx, 0, 0,cx}; + Quat qy = {0, sy, 0,cy}; + Quat qz = {0, 0,sz,cz}; + /*按ZXY顺序*/ + quat_multiply(&qx, &qz, &qx); + quat_multiply(&qy, &qx, out); + ssr_assert(quat_isidentity(out)); +} + +/* +** 四元数直接对向量进行旋转和先转成矩阵形式再旋转计算量一样 +*/ +void quat_applytovec3(Quat* q, Vec3* v, Vec3* out) { + ssr_assert(q && v && out); + if (v == out) { + shrvec3(v); + } + /* q[v,0]q^-1 */ + /* https://gamedev.stackexchange.com/questions/28395/rotating-vector3-by-a-quaternion */ + float x = q->x * 2.0F; + float y = q->y * 2.0F; + float z = q->z * 2.0F; + float xx = q->x * x; + float yy = q->y * y; + float zz = q->z * z; + float xy = q->x * y; + float xz = q->x * z; + float yz = q->y * z; + float wx = q->w * x; + float wy = q->w * y; + float wz = q->w * z; + + /*从这里能看到矩阵形式*/ + out->x = (1.0f - (yy + zz)) * v->x + (xy - wz) * v->y + (xz + wy) * v->z; + out->y = (xy + wz) * v->x + (1.0f - (xx + zz)) * v->y + (yz - wx) * v->z; + out->z = (xz - wy) * v->x + (yz + wx) * v->y + (1.0f - (xx + yy)) * v->z; +} + +void quat_tomat4(Quat* q, Mat4* out) { + ssr_assert(q && out); + ssr_assert(quat_isidentity(q)); + mat4_setidentity(out); + + float x = q->x * 2.0F; /*从quat_applytovec3能得到矩阵形式*/ + float y = q->y * 2.0F; + float z = q->z * 2.0F; + float xx = q->x * x; + float yy = q->y * y; + float zz = q->z * z; + float xy = q->x * y; + float xz = q->x * z; + float yz = q->y * z; + float wx = q->w * x; + float wy = q->w * y; + float wz = q->w * z; + + /*和mat4_setaxisangle实际上结果一样*/ + out->l[0] = 1.0f - (yy + zz); + out->l[1] = xy + wz; + out->l[2] = xz - wy; + + out->l[4] = xy - wz; + out->l[5] = 1.0f - (xx + zz); + out->l[6] = yz + wx; + + out->l[8] = xz + wy; + out->l[9] = yz - wx; + out->l[10] = 1.0f - (xx + yy); +} + +void quat_toeuler(Quat* q, Euler* out) { + ssr_assert(q && out); + Mat4 mat; + quat_tomat4(q, &mat); /*计算变换矩阵*/ + mat4_toeuler(&mat, out); /*mat是按照RyRxRz顺序(z->y->x)的旋转矩阵*/ +} + +void quat_normalize(Quat* q, Quat* out) { + ssr_assert(q && out); + float mag = quat_magnitude(q); + quat_scale(q, 1.f/mag, out); +} + +void quat_scale(Quat* q, float scale, Quat* out) { + ssr_assert(q && out); + out->x = q->x * scale; + out->y = q->y * scale; + out->z = q->z * scale; + out->w = q->w * scale; +} + +void quat_negtive(Quat* in, Quat* out) { + ssr_assert(in && out); + out->w = -in->w; + out->x = -in->x; + out->y = -in->y; + out->z = -in->z; +} + +void quat_devide(Quat* q, float k, Quat* out) { + ssr_assert(q && out); + k = 1 / k; + out->w = q->w * k; + out->x = q->x * k; + out->y = q->y * k; + out->z = q->z * k; +} + +void quat_invert(Quat* q, Quat* out) { + ssr_assert(q && out); + float mag2 = quat_magnitude2(q); + quat_conjugate(q, out); + quat_devide(out, mag2, out); + /*实际上如果是单位四元数，共轭就是逆*/ +} + +bool quat_setlookrotation(Vec3* view, Vec3* up, Quat* out) { + ssr_assert(view && up && out); + Mat4 m; + mat4_setlookrotation(view, up, &m); /*先以view为准构建正交矩阵*/ + mat4_toquat(&m, out); + return 1; +} + +/*q1 * q2^-1*/ +void quat_minus(Quat* q1, Quat* q2, Quat* out) { + ssr_assert(q1 && q2 && out); + Quat q2i; quat_invert(q2, &q2i); + quat_multiply(q1, &q2i, out); +} + +void quat_multiply(Quat* q1, Quat* q2, Quat* out) { + ssr_assert(q1 && q2 && out); + + float w1 = q1->w, x1 = q1->x, y1 = q1->y, z1 = q1->z, + w2 = q2->w, x2 = q2->x, y2 = q2->y, z2 = q2->z; + out->x = x1 * w2 + x2 * w1 + y1 * z2 - z1 * y2; + out->y = w1 * y2 + w2 * y1 + z1 * x2 - z2 * x1; + out->z = w1 * z2 + w2 * z1 + x1 * y2 - x2 * y1; + out->w = w1 * w2 - x1 * x2 - y1 * y2 - z1 * z2; +} + +/* +** 共轭的几何意义是颠倒旋转轴方向，代表了和q相反的角位移 +*/ +void quat_conjugate(Quat* in, Quat* out) { + ssr_assert(in && out); + out->w = in->w; + out->x = -in->x; + out->y = -in->y; + out->z = -in->z; +} + +float quat_dot(Quat* q1, Quat* q2) { + ssr_assert(q1 && q2); + return (q1->x*q2->x + q1->y*q2->y + q1->z*q2->z + q1->w*q2->w); +} + +float quat_magnitude(Quat* q) { + ssr_assert(q); + return sqrt(quat_dot(q, q)); +} + +float quat_magnitude2(Quat* q) { + ssr_assert(q); + return quat_dot(q, q); +} + +void quat_slerp(Quat* q1, Quat* q2, float t, Quat* out) { + ssr_assert(q1 && q2 && out); + ssr_assert(quat_isidentity(q1) && quat_isidentity(q2)); /*适用于单位四元数*/ + float dot = quat_dot(q1, q2); + Quat temp; + if (dot < 0.0f) { + dot = -dot; + temp.x = -q2->x; + temp.y = -q2->y; + temp.z = -q2->z; + temp.w = -q2->w; + } else { + temp = *q2; + } + if (dot < 0.95f) { + float angle = acos(dot); + float sinadiv, sinat, sinaomt; + sinadiv = 1.0f / sin(angle); + sinat = sin(angle*t); + sinaomt = sin(angle*(1.0f - t)); + out->x = (q1->x*sinaomt + temp.x*sinat)*sinadiv; + out->y = (q1->y*sinaomt + temp.y*sinat)*sinadiv; + out->z = (q1->z*sinaomt + temp.z*sinat)*sinadiv; + out->w = (q1->w*sinaomt + temp.w*sinat)*sinadiv; + } else { /*小范围时使用lerp*/ + quat_lerp(q1, &temp, t, out); + } +} + +void quat_lerp(Quat* q1, Quat* q2, float t, Quat* out) { + ssr_assert(q1 && q2 && out); + float t2 = 1 - t; + if (quat_dot(q1, q2) < 0.f) { /*点乘不能是负的，翻转一个四元数的符号，使得落回360°*/ + out->x = q1->x * t2 - t * q2->x; + out->y = q1->y * t2 - t * q2->y; + out->z = q1->z * t2 - t * q2->z; + out->w = q1->w * t2 - t * q2->w; + } else { + out->x = q1->x * t2 + t * q2->x; + out->y = q1->y * t2 + t * q2->y; + out->z = q1->z * t2 + t * q2->z; + out->w = q1->w * t2 + t * q2->w; + } + quat_normalize(out, out); +} diff --git a/src/math/vec2.c b/src/math/vec2.c new file mode 100644 index 0000000..86192dc --- /dev/null +++ b/src/math/vec2.c @@ -0,0 +1,36 @@ +#include "math.h" +#include "../util/assert.h" +#include "../core/mem.h" + +void vec2_scale(Vec2* v, float k, Vec2* out) { + ssr_assert(v && out); + out->x = v->x * k; + out->y = v->y * k; +} + +void vec2_plus(Vec2* v1, Vec2* v2, Vec2* out) { + ssr_assert(v1 && v2 && out); + out->x = v1->x + v2->x; + out->y = v1->y + v2->y; +} + +void vec2_offset(Vec2* v, float offset, Vec2* out) { + ssr_assert(v && out); + out->x = v->x + offset; + out->y = v->y + offset; +} + +float vec2_dot(Vec2* v1, Vec2* v2) { + ssr_assert(v1 && v2); + float d = v1->x * v2->x + v1->y * v2->y; + return d; +} + +void vec2_tostring(Vec2* v, char buf[]) { + sprintf(buf, "%8.3f %8.3f", v->x, v->y); +} + +void vec2_print(Vec2* v) { + vec2_tostring(v, printbuffer); + printf("%s", printbuffer); +}
\ No newline at end of file diff --git a/src/math/vec3.c b/src/math/vec3.c new file mode 100644 index 0000000..4d77b96 --- /dev/null +++ b/src/math/vec3.c @@ -0,0 +1,120 @@ +#include <math.h> + +#include "math.h" +#include "../util/assert.h" +#include "../core/mem.h" + +Vec3 vec3forward = {0,0,1}; +Vec3 vec3up = {0, 1, 0}; +Vec3 vec3left = {1, 0, 0}; + +void vec3_cross(Vec3* v1, Vec3* v2, Vec3* out) { + ssr_assert(v1 && v2 && out); + out->x = v1->y*v2->z - v1->z*v2->y; + out->y = v1->z*v2->x - v1->x*v2->z; + out->z = v1->x*v2->y - v1->y*v2->x; +} + +void vec3_scale(Vec3* v, float k, Vec3* out) { + ssr_assert(v && out); + out->x = v->x * k; + out->y = v->y * k; + out->z = v->z * k; +} + +void vec3_plus(Vec3* v1, Vec3* v2, Vec3* out) { + ssr_assert(v1 && v2 && out); + out->x = v1->x + v2->x; + out->y = v1->y + v2->y; + out->z = v1->z + v2->z; +} + +void vec3_offset(Vec3* v, float offset, Vec3* out) { + ssr_assert(v && out); + out->x = v->x + offset; + out->y = v->y + offset; + out->z = v->z + offset; +} + +float vec3_magnitude(Vec3* v) { + ssr_assert(v); + float l2 = vec3_magnitude2(v); + return sqrt(l2); +} + +float vec3_magnitude2(Vec3* v) { + ssr_assert(v); + return v->x*v->x + v->y*v->y + v->z*v->z; +} + +void vec3_lerp(Vec3* v1, Vec3* v2, float t, Vec3* out) { + ssr_assert(v1 && v2 && out); + float t2 = 1 - t; + out->x = v1->x * t2 + v2->x * t; + out->y = v1->y * t2 + v2->y * t; + out->z = v1->z * t2 + v2->z * t; +} + +void vec3_slerp(Vec3* v1, Vec3* v2, float t, Vec3* out) { + ssr_assert(v1 && v2 && out); + float mag1 = vec3_magnitude(v1); + float mag2 = vec3_magnitude(v2); + if (mag1 < EPSILON || mag2 < EPSILON) { + vec3_lerp(v1, v2, t, out); + return; + } + float lerplen = lerp(mag1, mag2, t); + /*分为长度上的lerp和角度上的slerp*/ + float dot = vec3_dot(v1, v2) / (mag1 * mag2); + if (dot > 1.f - EPSILON) { /*几乎是同向*/ + vec3_lerp(v1, v2, t, out); + return; + } else if (dot < -1 + EPSILON) { /*几乎是反向*/ + + } else { + Vec3 axis; vec3_cross(v1, v2, &axis); vec3_normalize(&axis, &axis); + Vec3 v1n; vec3_normalize(v1, &v1n); + float angle = acos(dot) * t; + //Mat4 m; mat4_setaxisangle(&axis, angle, &m); + //mat4_applytovec4(&m, &v1n, &v1n); + Quat q; quat_fromaxisangle(&axis, angle, &q); + quat_applytovec3(&q, &v1n, &v1n); + vec3_scale(&v1n, lerplen, out); + } +} + +float vec3_dot(Vec3* v1, Vec3* v2) { + ssr_assert(v1 && v2); + return v1->x * v2->x + v1->y * v2->y + v1->z * v2->z; +} + +void vec3_multiply(Vec3* v1, Vec3* v2, Quat* out) { + ssr_assert(v1 && v2 && out); + vec3_cross(v1, v2, out); + out->w = -vec3_dot(v1, v2); +} + +void vec3_normalize(Vec3* v, Vec3* out) { + ssr_assert(v && out); + //float mag = rsqrt(v->x * v->x + v->y * v->y + v->z * v->z); + float mag = 1.f / vec3_magnitude(v); + out->x = v->x * mag; + out->y = v->y * mag; + out->z = v->z * mag; +} + +void vec3_tostring(Vec3* v, char buf[]) { + sprintf(buf, "%8.3f %8.3f %8.3f", v->x, v->y, v->z); +} + +void vec3_print(Vec3* v) { + vec3_tostring(v, printbuffer); + printf("\n%s\n", printbuffer); +} + +void vec3_minus(Vec3* v1, Vec3* v2, Vec3* out) { + ssr_assert(v1 && v2 && out); + out->x = v1->x - v2->x; + out->y = v1->y - v2->y; + out->z = v1->z - v2->z; +}
\ No newline at end of file diff --git a/src/math/vec4.c b/src/math/vec4.c new file mode 100644 index 0000000..4216a08 --- /dev/null +++ b/src/math/vec4.c @@ -0,0 +1,24 @@ +#include <math.h> + +#include "math.h" +#include "../util/assert.h" +#include "../core/mem.h" + + + +void vec4_dividew(Vec4* v, Vec3* out) { + ssr_assert(out && v); + float w = 1.f / v->w; + out->x = v->x * w ; + out->y = v->y * w; + out->z = v->z * w; +} + +void vec4_tostring(Vec4* v, char buf[]) { + sprintf(buf, "%8.3f %8.3f %8.3f %8.3f", v->x, v->y, v->z, v->w); +} + +void vec4_print(Vec4* v) { + vec4_tostring(v, printbuffer); + printf("\n%s\n", printbuffer); +}
\ No newline at end of file diff --git a/src/ssr.h b/src/ssr.h new file mode 100644 index 0000000..5f80e92 --- /dev/null +++ b/src/ssr.h @@ -0,0 +1,6 @@ +#ifndef _SSR_H_ +#define _SSR_H_ + +#include "core/device.h" + +#endif
\ No newline at end of file diff --git a/src/test/test.h b/src/test/test.h new file mode 100644 index 0000000..4e7d77e --- /dev/null +++ b/src/test/test.h @@ -0,0 +1,15 @@ +#ifndef _SOFTSHADEROOM_TEST_H_ +#define _SOFTSHADEROOM_TEST_H_ + +#include <stdlib.h> +#include <stdio.h> +#include <string.h> + +#define print(str) printf("%s", str) +#define TEST(f) void f(){print(#f);print("\n"); +#define END } + +void test_mat4(); +void test_quat(); + +#endif
\ No newline at end of file diff --git a/src/test/test_mat4.c b/src/test/test_mat4.c new file mode 100644 index 0000000..e7a65fc --- /dev/null +++ b/src/test/test_mat4.c @@ -0,0 +1,208 @@ +#include "test.h" +#include "../math/math.h" +#include "../util/time.h" + +Mat4 m1 = { + 1, 5, 9,13, + 2, 6,10,14, + 3, 7,11,15, + 4, 8,12,16 +}; + +Mat4 m2 = { + 16,12, 8, 4, + 15,11, 7, 3, + 14,10, 6, 2, + 13, 9, 5, 1 +}; + +static void test_multiply() { + print("test_multiply\n"); + Mat4 res; + mat4_multiply(&m1, &m2, &res); + mat4_print(&res); +} + +static void test_invert() { + print("test_invert\n"); + Mat4 m = { + 0, 3, 0, 0, + 2, 7, 0, 1, + 8, 1, 1, 0, + 6, 0, 2, 1, + }; + Mat4 res, res2; + mat4_invertfull(&m, &res); + mat4_print(&res); + mat4_multiply(&m1, &m, &res2); + mat4_multiply(&res2, &res, &res2); + mat4_print(&res2); + /* + -0.250 0.333 1.667 -1.833 + 0.083 0.000 -0.667 0.833 + 0.167 -0.000 -0.333 -0.333 + -0.083 0.000 0.667 0.167 + */ +} + +static void test_transpose() { + print("test_transpose\n"); + Mat4 res; + mat4_transpose(&m1, &res); + mat4_print(&res); +} + +static void test_scale() { + print("test_scale\n"); + Mat4 scale; + mat4_setscale(1, 2, 3, &scale); + Mat4 res; + mat4_multiply(&scale, &m1, &res); + mat4_print(&scale); +} + +static void test_position() { + print("test_position\n"); + Mat4 pos; + mat4_setposition(1, 2, 3, &pos); + mat4_print(&pos); +} + +static void test_rotation() { + print("test_rotation\n"); + Mat4 i; + mat4_setidentity(&i); + Vec3 axis = { 0, 0, 1 }; + mat4_rotate(&i, 90, &axis, &i); + mat4_print(&i); +} + +static void test_invertgeneral3d() { + print("test_invertgeneral3d\n"); + Mat4 trans; + mat4_setidentity(&trans); + print("original matrix"); mat4_print(&m1); + Vec3 v; + v.x = 1; v.y = 3.3; v.z = 0; mat4_translate(&trans, &v, &trans); + v.x = 1; v.y = 0; v.z = 0; mat4_rotate(&trans, 23, &v, &trans); + v.x = 1; v.y = 1.2; v.z = 1; mat4_scale(&trans, &v, &trans); + v.x = 1; v.y = 3; v.z = 0; mat4_rotate(&trans, 23, &v, &trans); + v.x = 2; v.y = 3.3; v.z = 0; mat4_translate(&trans, &v, &trans); + print("transform matrix"); mat4_print(&trans); + Mat4 res; + mat4_multiply(&trans, &m1, &res); + mat4_invertgeneral3d(&trans, &trans); + print("inverse transform matrix"); mat4_print(&trans); + mat4_multiply(&trans, &res, &res); + print("inverted result"); mat4_print(&res); +} + +void test_invertrot() { + print("test_invertrot\n"); + print("original matrix"); mat4_print(&m1); + Mat4 rot; mat4_setidentity(&rot); + Vec3 axis = { 1,2,3 }; mat4_rotate(&rot, 90, &axis, &rot); + print("rotation matrix"); mat4_print(&rot); + Mat4 res; + mat4_multiply(&rot, &m1, &res); + print("result matrix"); mat4_print(&res); + mat4_invertgeneral3d(&rot, &rot); + print("inverse rotation matrix"); mat4_print(&rot); + mat4_multiply(&rot, &res, &res); + print("result matrix"); mat4_print(&res); +} + +TEST(test_orthogonalize) + ROWMAT(A0, + 1, 2, 4, 0, + 0, 0, 5, 0, + 0, 3, 6, 0, + 0, 0, 0, 0 + ); + Mat4 A; + + TIME_STAMP("mat4_orthogonalize time used") + for (int i = 0; i < 1000000; ++i) { + mat4_orthogonalize(&A0, &A); + } + TIME_STAMP_END + + TIME_STAMP("mat4_isorthogonal time used") + for (int i = 0; i < 1000000; ++i) { + mat4_isorthogonal(&A0); + } + TIME_STAMP_END + + mat4_print(&A); + float res = vec3_dot(&A.basis.y, &A.basis.x); + printf("%f\n", res); + printf("%d\n", mat4_isorthogonal(&A)); + Mat4 rev; + mat4_transpose(&A, &rev); + mat4_multiply(&A, &rev, &rev); + mat4_print(&rev); +END + +TEST(test_toeuler) + Euler euler = { 90, 0, 33 }; + Quat rot; + quat_fromeuler(&euler, &rot); + Mat4 m; + quat_tomat4(&rot, &m); + Euler res; + quat_toeuler(&rot, &res); + euler_print(&res); +END + +TEST(test_rotatematrix) + Euler euler = { 10, 20, 30 }; + Quat rot; + quat_fromeuler(&euler, &rot); + Mat4 m; + quat_tomat4(&rot, &m); + mat4_print(&m); + Mat4 mz, mx, my; + mat4_setrotatez(30, &mz); + mat4_setrotatex(10, &mx); + mat4_setrotatey(20, &my); + mat4_multiply(&mx, &mz, &mx); + mat4_multiply(&my, &mx, &m); + mat4_print(&m); + mat4_setrotate(10, 20, 30, &m); + mat4_print(&m); +END + +TEST(test_isidentity) + ROWMAT(mat, + 1, 0, 0, 0, + 0, 1, 0, 0, + 0, 0, 1, 0, + 0, 0, 0, 1 + ) + bool res; + TIME_STAMP_1000000("test is isidentity") + res = mat4_isidentity(&mat); + TIME_STAMP_END + printf("%d\n", res); +END + +void test_mat4() { + print("================================================\n"); + print("test_mat4\n"); + print("================================================\n"); + + //test_transpose(); + //test_multiply(); + //test_invert(); + //test_scale(); + //test_position(); + //test_rotation(); + //test_invertgeneral3d(); + + //test_invertrot(); + //test_orthogonalize(); + //test_toeuler(); + test_rotatematrix(); + + //test_isidentity(); +} diff --git a/src/test/test_quat.c b/src/test/test_quat.c new file mode 100644 index 0000000..c932082 --- /dev/null +++ b/src/test/test_quat.c @@ -0,0 +1,82 @@ +#include "test.h" +#include "../math/math.h" + +TEST(test_euler) + +Euler e = {90, 90, 0}; +Vec3 p = {2,3,4}, res; +Quat qut; +quat_fromeuler(&e, &qut); +quat_applytovec3(&qut, &p, &res); +vec3_print(&res); + +/*用轴角定义*/ +Vec3 axis = { 1, 0,0 }; +Quat qut2; +quat_fromaxisangle(&axis, 90, &qut2); +quat_applytovec3(&qut, &p, &res); +vec3_print(&res); + +Mat4 mat; +quat_tomat4(&qut, &mat); +mat4_print(&mat); + +Vec4 p2 = { 2,3,4,1 }; +Vec4 res2; +mat4_applytovec4(&mat, &p2, &res2); +vec4_print(&res2); + +END + +TEST(test_quat2euler) +Euler e = { 10, 20, 30}; + +Quat q; +euler_toquat(&e, &q); +euler_deg2rad(&e, &e); +quat_print(&q); +euler_print(&e); + +Euler res; +quat_toeuler(&q, &res); + +euler_print(&res); + +quat_fromeuler(&res, &q); +quat_print(&q); + +END + +TEST(test_conjugate) +Quat rot; +Vec3 axis = {2,3,4}; +vec3_normalize(&axis, &axis); +quat_fromaxisangle(&axis, 46, &rot); +Vec3 p = {12, 3, 55}; +vec3_print(&p); +Vec3 res; +quat_applytovec3(&rot, &p, &res); +vec3_print(&res); +Quat rev; +quat_invert(&rot, &rev); +quat_applytovec3(&rev, &res, &res); +vec3_print(&res); +END + +TEST(test_lookrotation) +Vec3 view = {1, 0, 0}; +Vec3 up = {0, 1, 0}; +Quat rot; +quat_setlookrotation(&view, &up, &rot); +quat_print(&rot); +Euler euler; +quat_toeuler(&rot, &euler); +euler_print(&euler); +END + +void test_quat() { + //test_euler(); + //test_quat2euler(); + //test_conjugate(); + //test_lookrotation(); +}
\ No newline at end of file diff --git a/src/util/assert.h b/src/util/assert.h new file mode 100644 index 0000000..563e991 --- /dev/null +++ b/src/util/assert.h @@ -0,0 +1,9 @@ +#ifndef _SOFTSHADEROOM_ASSERT_H_ +#define _SOFTSHADEROOM_ASSERT_H_ + +#include <assert.h> + +#define ssr_assert assert + + +#endif
\ No newline at end of file diff --git a/src/util/time.c b/src/util/time.c new file mode 100644 index 0000000..5bfc936 --- /dev/null +++ b/src/util/time.c @@ -0,0 +1,5 @@ +#include "time.h" + +clock_t sharedTick; + +const char* sharedTickHint; diff --git a/src/util/time.h b/src/util/time.h new file mode 100644 index 0000000..6b850b4 --- /dev/null +++ b/src/util/time.h @@ -0,0 +1,24 @@ +#ifndef _SOFTSHADEROOM_TIME_H_ +#define _SOFTSHADEROOM_TIME_H_ + +#include <time.h> +#include "assert.h" + +extern clock_t sharedTick; +extern const char* sharedTickHint; + +#define TIME_STAMP(hint) \ +ssr_assert(!sharedTickHint);\ +sharedTickHint = hint;\ +sharedTick = clock();{ +#define TIME_STAMP_END \ +}\ +sharedTick= clock()-sharedTick;\ +printf("%s : %.10fs %ldtick\n", sharedTickHint, sharedTick / (double)CLOCKS_PER_SEC, sharedTick);\ +sharedTickHint = 0; +#define TIME_STAMP_i(hint, i)\ +TIME_STAMP(hint);}\ +for(int _i = 0; _i < i; ++_i){ +#define TIME_STAMP_1000000(hint) TIME_STAMP_i(hint, 1000000) + +#endif
\ No newline at end of file diff --git a/src/util/type.h b/src/util/type.h new file mode 100644 index 0000000..709cb32 --- /dev/null +++ b/src/util/type.h @@ -0,0 +1,7 @@ +#ifndef _SOFTSHADEROOM_TYPE_H_ +#define _SOFTSHADEROOM_TYPE_H_ + +typedef unsigned int uint; +typedef int bool; + +#endif
\ No newline at end of file |