Correct voxel raycast, voxel types, reflectance
dkl9

dkl9 commited on 2025-175 11:55:45
Showing 10 changed files, with 192 additions and 74 deletions.

... ...
@@ -1,4 +1,5 @@
1 1
 *.swp
2 2
 *.log
3
-venv
3
+__pycache__/
4
+venv/
4 5
 tags
... ...
@@ -22,8 +22,8 @@ class Camera:
22 22
         return self.viewport[0] / self.viewport[1]
23 23
 
24 24
     def view_mat(self):
25
-        f = normalize(self.looking())
26
-        s = normalize(np.cross(f, self.up()))
25
+        f = normalise(self.looking())
26
+        s = normalise(np.cross(f, self.up()))
27 27
         u = np.cross(s, f)
28 28
         m = np.identity(4, dtype=FLOAT)
29 29
         m[:3, :3] = [s, u, -f]
... ...
@@ -50,10 +50,13 @@ class Camera:
50 50
         else:
51 51
             return self.persp_mat() @ self.view_mat() @ ent.translate() @ ent.orient
52 52
 
53
-    def draw(self, ent, loc_mvp):
53
+    def draw(self, ent):
54
+        if not ent.verts.shape[0]:
55
+            return
54 56
         if not ent.vao:
55 57
             ent.build_vao()
56
-        glUniformMatrix4fv(loc_mvp, 1, GL_TRUE, self.mvp(ent))
58
+        glUniformMatrix4fv(UNIFORMS["MVP"][2], 1, GL_TRUE, self.mvp(ent))
59
+        glUniform3fv(UNIFORMS["camPos"][2], 1, self.pos)
57 60
         glBindVertexArray(ent.vao)
58 61
         glDrawArrays(ent.draw_mode, 0, len(ent.verts))
59 62
 
... ...
@@ -1,51 +1,48 @@
1 1
 from libs import *
2 2
 
3
-# https://iquilezles.org/articles/intersectors/
4
-def intersect_tri(ro, rd, tv):
5
-    e1 = tv[1] - tv[0]
6
-    e2 = tv[2] - tv[0]
7
-    to = ro - tv[0]
8
-    n = np.cross(e1, e2)
9
-    q = np.cross(to, rd)
10
-    d = 1.0 / np.dot(rd, n)
11
-    u = d * np.dot(-q, e2)
12
-    v = d * np.dot(q, e1)
13
-    t = d * np.dot(-n, to)
14
-    if u > 0.0 and v > 0.0 and u + v < 1.0 and t > 1e-3:
15
-        return t
16
-
17 3
 class Entity:
18 4
     def __init__(self, verts, col, pos, **kwargs):
19 5
         self.verts = vec(verts)
20 6
         self.col = vec(col)
7
+        self.refls = vec(kwargs.get("refls", [0.5, 0.5, 0.0, 32.0]))
21 8
         self.pos = vec(pos)
22 9
         self.draw_mode = kwargs.get("mode", GL_TRIANGLE_STRIP)
23 10
         self.hud = kwargs.get("hud", False)
24 11
         self.orient = np.identity(4, dtype=FLOAT)
25 12
         self.onclick = lambda s, b, n: print(f"clicked {s}:{n} with button {b}")
26 13
         self.vao = None
27
-        self.vbos = [[None, 0], [None, 0]]
14
+        self.vbos = [[None, 0] for _ in VERTEX_ATTRS]
28 15
 
29 16
     def col_list(self):
30 17
         return np.repeat(vec([self.col]), len(self.verts), 0)
31 18
 
32
-    def update_attr(self, i, l=None):
33
-        if l is None:
34
-            l = self.verts if i == 0 else self.col_list()
35
-        if self.vao:
36
-            glBindBuffer(GL_ARRAY_BUFFER, self.vbos[i][0])
37
-            if l.shape[0] != self.vbos[i][1]:
38
-                glBufferData(GL_ARRAY_BUFFER, l.nbytes, l, GL_DYNAMIC_DRAW)
39
-                self.vbos[i][1] = l.shape[0]
19
+    def norm_list(self):
20
+        l = []
21
+        if self.draw_mode == GL_TRIANGLE_STRIP:
22
+            l.append(norm_tri(self.verts[:3]))
23
+            for i in range(1, len(self.verts) - 1):
24
+                l.append(norm_tri(self.verts[i - 1:i + 2]))
25
+            l.append(norm_tri(self.verts[-3:]))
26
+        elif self.draw_mode == GL_TRIANGLES:
27
+            for i in range(0, len(self.verts), 3):
28
+                for _ in range(3):
29
+                    l.append(norm_tri(self.verts[i:i + 3]))
30
+        elif self.draw_mode == GL_LINE_STRIP or self.draw_mode == GL_LINES:
31
+            for _ in range(len(self.verts)):
32
+                l.append([1, 0, 0])
40 33
         else:
41
-                glBufferSubData(GL_ARRAY_BUFFER, 0, l.nbytes, l)
34
+            raise NotImplementedError
35
+        return (self.orient[:3, :3] @ vec(l).T).T
36
+
37
+    def refls_list(self):
38
+        return np.repeat(vec([self.refls]), len(self.verts), 0)
42 39
 
43 40
     def build_vao(self):
44 41
         self.vao = glGenVertexArrays(1)
45 42
         glBindVertexArray(self.vao)
46
-        for (i, l) in enumerate((self.verts, self.col_list())):
43
+        for i in range(len(VERTEX_ATTRS)):
47 44
             self.vbos[i][0] = glGenBuffers(1)
48
-            self.update_attr(i, l)
45
+            l = update_attr(self, i)
49 46
             glEnableVertexAttribArray(i)
50 47
             glVertexAttribPointer(i, l.shape[1], GL_FLOAT, GL_FALSE, 0, None)
51 48
 
... ...
@@ -0,0 +1,51 @@
1
+from OpenGL.GL import *
2
+from OpenGL.GLU import *
3
+from OpenGL.GLUT import *
4
+import sys
5
+import time
6
+
7
+lf = time.time()
8
+fc = 0
9
+
10
+def reshape(width, height):
11
+    glViewport(0, 0, width, height)
12
+    glMatrixMode(GL_PROJECTION)
13
+    glLoadIdentity()
14
+    gluPerspective(45, width / float(height or 1), 1, 50)
15
+    glMatrixMode(GL_MODELVIEW)
16
+
17
+def draw_cube():
18
+    glutWireCube(1)
19
+
20
+def display():
21
+    global lf, fc
22
+    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
23
+    glLoadIdentity()
24
+    glTranslatef(0, 0, -5)
25
+    glRotatef((60 * time.time()) % 360, 2, 1, 0)
26
+    draw_cube()
27
+    glutSwapBuffers()
28
+    fc += 1
29
+    now = time.time()
30
+    if now - lf >= 1.0:
31
+        print(f"{fc} FPS")
32
+        fc = 0
33
+        lf = now
34
+
35
+def keyboard(key, x, y):
36
+    if key == b'\x1b':  # ESC
37
+        glutLeaveMainLoop()
38
+
39
+glutInit(sys.argv)
40
+glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGB | GLUT_DEPTH)
41
+glutInitWindowSize(500, 500)
42
+glutInitWindowPosition(100, 100)
43
+glutCreateWindow("PyOpenGL Cube")
44
+glEnable(GL_DEPTH_TEST)
45
+
46
+glutDisplayFunc(display)
47
+glutIdleFunc(display)
48
+glutReshapeFunc(reshape)
49
+glutKeyboardFunc(keyboard)
50
+
51
+glutMainLoop()
... ...
@@ -5,6 +5,8 @@ import pygame
5 5
 
6 6
 FLOAT = np.float32
7 7
 VOXEL_TYPE = np.uint8
8
+VERTEX_ATTRS = ["vec3 position", "vec4 col", "vec3 normal", "vec4 refls"]
9
+UNIFORMS = {"MVP": ["vert", "mat4", None], "lightPos": ["frag", "vec4", None], "camPos": ["frag", "vec3", None]}
8 10
 
9 11
 def vec(l):
10 12
     return np.array(l, dtype=FLOAT)
... ...
@@ -13,7 +15,7 @@ RED = vec([1.0, 0.0, 0.0, 1.0])
13 15
 YELLOW = vec([1.0, 1.0, 0.0, 1.0])
14 16
 GREEN = vec([0.0, 1.0, 0.0, 1.0])
15 17
 
16
-def normalize(v):
18
+def normalise(v):
17 19
     norm = np.linalg.norm(v)
18 20
     return v / norm if norm > 0 else v
19 21
 
... ...
@@ -28,3 +30,46 @@ def axis2offset(a, x=1):
28 30
     o = [0, 0, x * np.sign(a)]
29 31
     roll_vec3(o, abs(a))
30 32
     return o
33
+
34
+# https://iquilezles.org/articles/intersectors/
35
+def intersect_tri(ro, rd, tv):
36
+    e1 = tv[1] - tv[0]
37
+    e2 = tv[2] - tv[0]
38
+    to = ro - tv[0]
39
+    n = np.cross(e1, e2)
40
+    q = np.cross(to, rd)
41
+    d = 1.0 / np.dot(rd, n)
42
+    u = d * np.dot(-q, e2)
43
+    v = d * np.dot(q, e1)
44
+    t = d * np.dot(-n, to)
45
+    if u > 0.0 and v > 0.0 and u + v < 1.0 and t > 1e-3:
46
+        return t
47
+
48
+def norm_tri(tv):
49
+    assert tv.shape[0] == 3
50
+    return normalise(np.cross(tv[2] - tv[0], tv[1] - tv[0]))
51
+
52
+def update_attr(ent, n, l=None):
53
+    if isinstance(n, int):
54
+        i, name = n, VERTEX_ATTRS[n].split()[1]
55
+    elif isinstance(n, str):
56
+        # https://stackoverflow.com/a/30197797
57
+        i, name = next((i for (i, s) in enumerate(VERTEX_ATTRS) if n in s)), n
58
+    if l is None:
59
+        match name:
60
+            case "position":
61
+                l = ent.verts
62
+            case "col":
63
+                l = ent.col_list()
64
+            case "normal":
65
+                l = ent.norm_list()
66
+            case "refls":
67
+                l = ent.refls_list()
68
+    if ent.vao:
69
+        glBindBuffer(GL_ARRAY_BUFFER, ent.vbos[i][0])
70
+        if l.shape[0] != ent.vbos[i][1]:
71
+            glBufferData(GL_ARRAY_BUFFER, l.nbytes, l, GL_DYNAMIC_DRAW)
72
+            ent.vbos[i][1] = l.shape[0]
73
+        else:
74
+            glBufferSubData(GL_ARRAY_BUFFER, 0, l.nbytes, l)
75
+        return l
... ...
@@ -6,24 +6,29 @@ from mesh import SquareMesh
6 6
 from player import Player
7 7
 from region import WorldRegion
8 8
 
9
-VS = """
9
+VS = f"""
10 10
 #version 120
11
-attribute vec3 position;
12
-attribute vec4 colIn;
13
-varying vec4 colOut;
14
-uniform mat4 MVP;
15
-void main() {
11
+{"\n".join([f"attribute {a};" for a in VERTEX_ATTRS])}
12
+{"\n".join([f"varying {a}F;" for a in VERTEX_ATTRS])}
13
+{"\n".join([f"uniform {v[1]} {k};" for (k, v) in UNIFORMS.items() if v[0] == "vert"])}
14
+void main() {{
16 15
     gl_Position = MVP * vec4(position, 1.0);
17
-    colOut = colIn;
18
-}
16
+{"\n".join([f"    {a.split()[1]}F = {a.split()[1]};" for a in VERTEX_ATTRS])}
17
+}}
19 18
 """
20 19
 
21
-FS = """
20
+FS = f"""
22 21
 #version 120
23
-varying vec4 colOut;
24
-void main() {
25
-    gl_FragColor = colOut;
26
-}
22
+{"\n".join([f"varying {a}F;" for a in VERTEX_ATTRS])}
23
+{"\n".join([f"uniform {v[1]} {k};" for (k, v) in UNIFORMS.items() if v[0] == "frag"])}
24
+void main() {{
25
+    float ambient = reflsF.x;
26
+    vec3 lightDir = normalize(lightPos.xyz - positionF);
27
+    float lambert = dot(normalF, lightDir);
28
+    float diffuse = reflsF.y * lightPos.w * lambert / pow(distance(lightPos.xyz, positionF), 2);
29
+    float specular = reflsF.z * lightPos.w * pow(max(dot(normalize(camPos - positionF), reflect(-lightDir, normalF)), 0.0), reflsF.w);
30
+    gl_FragColor = vec4((ambient + diffuse + specular) * colF.xyz, colF.w);
31
+}}
27 32
 """
28 33
 
29 34
 FPS = 60
... ...
@@ -41,8 +46,8 @@ def click_mesh(self, button, coords):
41 46
             offset = axis2offset(coords[0])
42 47
             self.set_voxel(coords[1] + offset, 2)
43 48
     if button % 2 == 1:
44
-        self.update_attr(0)
45
-        self.update_attr(1)
49
+        for i in range(len(VERTEX_ATTRS)):
50
+            update_attr(self, i)
46 51
 
47 52
 def build_scene(mg):
48 53
     box = WorldRegion()
... ...
@@ -59,10 +64,12 @@ def main():
59 64
         shaders.compileShader(VS, GL_VERTEX_SHADER),
60 65
         shaders.compileShader(FS, GL_FRAGMENT_SHADER)
61 66
     )
62
-    loc_mvp = glGetUniformLocation(sp, "MVP")
67
+    for (k, v) in UNIFORMS.items():
68
+        v[2] = glGetUniformLocation(sp, k)
63 69
     glUseProgram(sp)
64 70
     glEnable(GL_DEPTH_TEST)
65 71
     glClearColor(0.1, 0.1, 0.1, 1.0)
72
+    glUniform4fv(UNIFORMS["lightPos"][2], 1, vec([5, 10, 10, 10]))
66 73
     mg = MapGen()
67 74
     scene = build_scene(mg)
68 75
     clock = pygame.time.Clock()
... ...
@@ -80,9 +87,9 @@ def main():
80 87
         cam.update_pos()
81 88
         glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
82 89
         for ent in scene:
83
-            cam.draw(ent, loc_mvp)
90
+            cam.draw(ent)
84 91
         perf = clock.get_fps() / FPS
85
-        cam.draw_hud(loc_mvp, scene[0], perf)
92
+        cam.draw_hud(scene[0], perf)
86 93
         pygame.display.flip()
87 94
         clock.tick(FPS)
88 95
 
... ...
@@ -20,8 +20,10 @@ class MapGen:
20 20
         dots = []
21 21
         for c in corners:
22 22
             ci = tuple(np.add(c, corner, dtype=int, casting="unsafe"))
23
-            dots.append(np.dot(fract - c, gt.setdefault(ci, np.random.rand(dim))))
24
-        return smerp(dots, fract)
23
+            cg = gt.setdefault(ci, normalise(np.random.normal(size=dim)))
24
+            assert abs(np.linalg.norm(cg) - 1) < 1e-3
25
+            dots.append(np.dot(fract - c, cg))
26
+        return 2 * smerp(dots, fract)
25 27
 
26 28
     def world_natural(self, wr):
27 29
         n = 0
... ...
@@ -29,8 +31,9 @@ class MapGen:
29 31
             for j in range(wr.r):
30 32
                 h = 8 + 4 * self.perlin("height", np.add(wr.pos[[0, 2]], (i, j)), np.array([10, 5]))
31 33
                 for k in range(int(h)):
32
-                    wr.grid[i, k, j] = 1
34
+                    wr.grid[i, k, j] = 3
33 35
                     n += 1
36
+                wr.grid[i, int(h), j] = 4
34 37
         return n
35 38
 
36 39
     def world_min(self, wr):
... ...
@@ -19,7 +19,8 @@ class SquareMesh(entity.Entity):
19 19
             [-w, -w, 0], [w, -w, 0], [-w, w, 0]
20 20
         ]), abs(axis), 1)
21 21
         if axis < 0:
22
-            sq = np.flip(sq, 0)
22
+            sq[1:3] = sq[2:0:-1]
23
+            sq[4:6] = sq[5:3:-1]
23 24
         self.verts = np.concatenate((self.verts, sq), dtype=FLOAT)
24 25
         for i in range(2):
25 26
             self.faces.append((axis, fd))
... ...
@@ -48,6 +49,10 @@ class SquareMesh(entity.Entity):
48 49
                         return top
49 50
                     case -2:
50 51
                         return bottom
52
+            case 3:
53
+                return vec([0.7, 0.7, 0.8, 1.0])
54
+            case 4:
55
+                return vec([0.5, 0.2, 0.0, 1.0])
51 56
 
52 57
     def col_list(self):
53 58
         cl = np.empty((3 * len(self.faces), 4), dtype=FLOAT)
... ...
@@ -79,7 +79,7 @@ class Player(camera.Camera):
79 79
             case pygame.WINDOWSIZECHANGED | pygame.WINDOWRESIZED:
80 80
                 self.viewport = reshape(ev.x, ev.y)
81 81
 
82
-    def draw_hud(self, loc_mvp, wr, perf):
82
+    def draw_hud(self, wr, perf):
83 83
         self.hud["skybox"].pos = self.pos
84 84
         c, ap = self.target([wr])
85 85
         if c:
... ...
@@ -88,6 +88,6 @@ class Player(camera.Camera):
88 88
         else:
89 89
             self.hud["shl"].pos = vec([0, 0, 0])
90 90
         self.hud["fps"].col = GREEN if 0.9 < perf < 1.1 else YELLOW if perf > 0.5 else RED
91
-        self.hud["fps"].update_attr(1)
91
+        update_attr(self.hud["fps"], "col")
92 92
         for ent in self.hud.values():
93
-            self.draw(ent, loc_mvp)
93
+            self.draw(ent)
... ...
@@ -15,13 +15,15 @@ def intersect_box(ro, rd, c, s):
15 15
 
16 16
 class WorldRegion(mesh.SquareMesh):
17 17
     def __init__(self, size=16):
18
-        super().__init__()
18
+        super().__init__(refls=[0.5, 0.5, 0.5, 32.0])
19 19
         self.grid = np.zeros((size, size, size), dtype=VOXEL_TYPE)
20 20
         self.r = size
21
-        self.dirty = 0
22 21
 
23 22
     def set_voxel(self, pos, vt=0):
24 23
         self.grid[*pos] = vt
24
+        if vt == 0:
25
+            self.remesh()
26
+        else:
25 27
             bp = bool(vt)
26 28
             for d in range(3):
27 29
                 for s in range(2):
... ...
@@ -31,12 +33,8 @@ class WorldRegion(mesh.SquareMesh):
31 33
                     o = np.array(axis2offset(a), dtype=np.int64)
32 34
                     bc = bool(self.grid[*(pos + o)])
33 35
                     if bp ^ bc:
34
-                    a = (bp - bc) * (d + 1)
36
+                        a *= bc - bp
35 37
                         self.add_sq(a, 1, pos + o / 2, pos + (0 if bp else o))
36
-        self.dirty += 1
37
-        if self.dirty > self.r:
38
-            self.dirty = 0
39
-            self.remesh()
40 38
 
41 39
     def remesh(self):
42 40
         offset = vec([0.5, 0.0, 0.0])
... ...
@@ -57,7 +55,10 @@ class WorldRegion(mesh.SquareMesh):
57 55
                                 [0, 1, 1], [0, -1, 1], [0, 1, -1],
58 56
                                 [0, -1, -1], [0, 1, -1], [0, -1, 1]
59 57
                             ])
60
-                            sq /= 2 * (bp - bc)
58
+                            sq /= 2
59
+                            if bp:
60
+                                sq[1:3] = sq[2:0:-1]
61
+                                sq[4:6] = sq[5:3:-1]
61 62
                             for v in sq:
62 63
                                 roll_vec3(v, d)
63 64
                             sq += p - offset / 2
... ...
@@ -68,19 +69,24 @@ class WorldRegion(mesh.SquareMesh):
68 69
                                 self.faces.append(((bp - bc) * (d + 1), p))
69 70
         self.verts = vec(vl)
70 71
 
72
+    def norm_list(self):
73
+        return np.astype(np.repeat([axis2offset(a) for (a, _) in self.faces], 3, 0), FLOAT)
74
+
71 75
     def raycast(self, ro, rd):
72 76
         vr = np.repeat(self.r, 3)
73 77
         tn, tf = intersect_box(ro, rd, self.pos + vr / 2 - 0.5, vr / 2)
74 78
         if tn is not None:
75
-            if tn < 0:
76
-                tn = 0
77
-            ids = [[(d, i) for d in ((np.arange(self.r) + 0.5) * np.sign(rd[i]) - np.mod(ro[i] + tn * rd[i], 1)) / rd[i]] for i in range(3)]
79
+            ids = []
80
+            for i in range(3):
81
+                k = np.arange(self.r + 1) if rd[i] >= 0 else np.arange(self.r + 1, 0, -1)
82
+                ids.append([(d, i) for d in (k - 0.5 - ro[i]) / rd[i]])
83
+            q, p = None, None
78 84
             for (te, i) in heapq.merge(*ids):
79
-                ts = tn + te
80
-                p = np.astype(np.round(ro + (ts + 1e-2) * rd - self.pos), np.int16)
81
-                if any(p < 0) or any(p >= self.r):
82
-                    break
83
-                elif self.grid[*p]:
85
+                if te < 0 or te < tn or te > tf - 1e-3:
86
+                    continue
87
+                q, p = p, np.astype(np.round(ro + (te + 1e-8) * rd - self.pos), np.int16)
88
+                assert q is None or abs(np.sum(p - q)) == 1
89
+                if self.grid[*p]:
84 90
                     sa = np.astype((i + 1) * -np.sign(rd[i]), np.int8)
85
-                    return ts, (sa, p)
91
+                    return te, (sa, p)
86 92
         return None, None
87 93