Initial start to disk renderer

examples/disk_viz.py

@@ -0,0 +1,22 @@
+import yt
+
+import yt_idv
+from yt_idv.scene_components.disk_slice import DiskSlice
+from yt_idv.scene_data.sliced_values import SlicedData
+
+ds = yt.load_sample("IsolatedGalaxy")
+dd = ds.all_data()
+
+rc = yt_idv.render_context(height=800, width=800, gui=True)
+sg = rc.add_scene(ds, None, no_ghost=True)
+
+sdata = SlicedData(data_source=dd)
+sren = DiskSlice(data=sdata)
+
+sg.data_objects.append(sdata)
+sdata.vertex_array
+sg.components.append(sren)
+
+print(sg.camera.near_plane, sg.camera.far_plane)
+# sg.camera.far_plane = 0.1
+rc.run()

yt_idv/scene_components/disk_slice.py

@@ -0,0 +1,65 @@
+import numpy as np
+import traitlets
+from OpenGL import GL
+
+from yt_idv.scene_components.base_component import SceneComponent
+from yt_idv.scene_data.sliced_values import SlicedData
+
+
+class DiskSlice(SceneComponent):
+    name = "disk_slice"
+    data = traitlets.Instance(SlicedData)
+    r_bounds = traitlets.Tuple((0.0, 1.0), trait=traitlets.CFloat())
+    theta_bounds = traitlets.Tuple((0.0, 2.0 * np.pi), trait=traitlets.CFloat())
+    phi_bounds = traitlets.Tuple((0.0, np.pi), trait=traitlets.CFloat())
+    center = traitlets.Tuple((0.5, 0.5, 0.5), trait=traitlets.CFloat())
+    disk_theta = traitlets.CFloat(0.0)
+    disk_phi = traitlets.CFloat(0.0)
+
+    def render_gui(self, imgui, renderer, scene):
+        changed = super(DiskSlice, self).render_gui(imgui, renderer, scene)
+        for val, (mi, ma) in (
+            ("r", (0.0, 1.0)),
+            ("theta", (0.0, 2.0 * np.pi)),
+            ("phi", (0.0, np.pi)),
+        ):
+            vals = getattr(self, f"{val}_bounds")
+            _, new_lower = imgui.slider_float2(val, vals[0], vals[1], mi, ma)
+            if _:
+                setattr(self, f"{val}_bounds", new_lower)
+                changed = True
+        _, val = imgui.slider_float("disk theta", self.disk_theta, 0, 2.0 * np.pi)
+        if _:
+            self.disk_theta = val
+            changed = True
+        _, val = imgui.slider_float("disk phi", self.disk_phi, 0.0, np.pi)
+        if _:
+            self.disk_phi = val
+            changed = True
+        return changed
+
+    def draw(self, scene, program):
+        each = self.data.vertex_array.each
+        GL.glEnable(GL.GL_CULL_FACE)
+        GL.glCullFace(GL.GL_BACK)
+        GL.glDrawArraysInstanced(GL.GL_TRIANGLE_STRIP, 0, each, self.data.size)
+
+    def _set_uniforms(self, scene, shader_program):
+        normal = np.array(
+            [
+                np.cos(self.disk_theta) * np.sin(self.disk_phi),
+                np.sin(self.disk_theta) * np.sin(self.disk_phi),
+                np.cos(self.disk_phi),
+            ]
+        )
+        cam = scene.camera
+        shader_program._set_uniform("r_bounds", np.array(self.r_bounds))
+        shader_program._set_uniform("theta_bounds", np.array(self.theta_bounds))
+        shader_program._set_uniform("phi_bounds", np.array(self.phi_bounds))
+        shader_program._set_uniform("disk_center", np.array(self.center))
+        shader_program._set_uniform("disk_normal", np.array(normal))
+        shader_program._set_uniform("projection", cam.projection_matrix)
+        shader_program._set_uniform("modelview", cam.view_matrix)
+        shader_program._set_uniform(
+            "inv_pmvm", np.linalg.inv(cam.projection_matrix @ cam.view_matrix)
+        )

yt_idv/scene_data/sliced_values.py

@@ -0,0 +1,41 @@
+import numpy as np
+import traitlets
+from yt.data_objects.data_containers import YTDataContainer
+
+from yt_idv.constants import aabb_triangle_strip
+from yt_idv.opengl_support import VertexArray, VertexAttribute
+from yt_idv.scene_data.base_data import SceneData
+
+
+class SlicedData(SceneData):
+    name = "sliced_data"
+    data_source = traitlets.Instance(YTDataContainer)
+    color_field = traitlets.Unicode("density", allow_none=True)
+    field_type = traitlets.Unicode("gas")
+    size = traitlets.CInt(-1)
+
+    @traitlets.default("vertex_array")
+    def _default_vertex_array(self):
+        va = VertexArray(name="data_positions", each=14)
+        va.attributes.append(
+            VertexAttribute(name="model_vertex", data=aabb_triangle_strip, divisor=0)
+        )
+        field = self.data_source[self.field_type, self.color_field].astype("f4").d
+        field.shape = (field.size, 1)
+        va.attributes.append(VertexAttribute(name="color_field", data=field, divisor=1))
+        self.size = field.size
+        positions = np.empty((self.size, 3), dtype="f4")
+        positions[:, 0] = self.data_source["index", "x"].in_units("unitary")
+        positions[:, 1] = self.data_source["index", "y"].in_units("unitary")
+        positions[:, 2] = self.data_source["index", "z"].in_units("unitary")
+        va.attributes.append(
+            VertexAttribute(name="position_field", data=positions, divisor=1)
+        )
+        widths = np.empty((self.size, 3), dtype="f4")
+        widths[:, 0] = self.data_source["index", "dx"].in_units("unitary")
+        widths[:, 1] = self.data_source["index", "dy"].in_units("unitary")
+        widths[:, 2] = self.data_source["index", "dz"].in_units("unitary")
+        va.attributes.append(
+            VertexAttribute(name="width_field", data=widths, divisor=1)
+        )
+        return va

yt_idv/shaders/disk.frag.glsl

@@ -0,0 +1,11 @@
+
+flat in vec4 v_model;
+flat in float field_value;
+out vec4 color;
+in vec2 UV;
+
+void main(){
+    //color = vec4(field_value * float(length(UV) < 1.0));
+    color = vec4(field_value);
+    color.a = 1.0;
+}

yt_idv/shaders/disk.vert.glsl

@@ -0,0 +1,55 @@
+in vec4 model_vertex;
+in vec3 position_field; // The location of the vertex in model space
+in vec3 width_field;
+in float color_field;
+
+flat out vec4 vradius;
+flat out mat4 inverse_proj;
+flat out mat4 inverse_mvm;
+flat out mat4 inverse_pmvm;
+flat out float field_value;
+flat out vec4 v_model;
+
+out vec2 UV;
+
+void main()
+{
+    vec3 v1, v2, v3;
+    v1 = disk_normal;
+    ortho_find(v1, v2, v3);
+    vec3 disk_position = position_field - disk_center;
+    float r = distance(disk_center, position_field);
+    float theta = atan(disk_position.y, disk_position.x) + PI;
+    float phi = acos(disk_position.z / r);
+
+    if ((r > r_bounds.y) || (r < r_bounds.x)
+        || (theta > theta_bounds.y) || (theta < theta_bounds.x)
+        || (phi > phi_bounds.y) || (phi < phi_bounds.x)) {
+        vradius = vec4(0.0);
+    } else {
+        vradius = vec4(width_field, 1.0);
+    }
+
+    v_model = vec4(position_field, 1.0) + vradius * (model_vertex - 0.5);
+    vec4 clip_pos = projection * modelview * v_model;
+    inverse_pmvm = inv_pmvm;
+    if (color_field == 0.0) {
+        field_value = 1.0;
+    } else {
+        field_value = color_field;
+    }
+
+    // since we have theta: [0, 2pi) and phi: [0, pi] we use
+    // r = sqrt(x^2 + y^2 + z^2)
+    // theta = arctan(y/x)
+    // phi = arccos(z/r)
+
+    //mat3 rotation = mat3(v1, v2, v3);
+    disk_position = (v_model.xyz - disk_center);
+    r = length(cross(disk_normal, disk_position));
+    theta = atan(dot(v3, disk_position), dot(v2, disk_position)) + PI;
+    phi = acos(dot(v1, disk_position) / r);
+
+    gl_Position = vec4(r * cos(theta), r * sin(theta), 1, 1);
+    UV = model_vertex.xy;
+}

yt_idv/shaders/known_uniforms.inc.glsl

@@ -46,3 +46,29 @@ uniform sampler3D ds_tex;
 
 // ray tracing control
 uniform float sample_factor;
+
+// disk slice renderer
+// in our convention, theta runs from 0 .. 2pi and phi from 0 .. pi
+uniform vec3 disk_center;
+uniform vec3 disk_normal;
+uniform vec2 theta_bounds; 
+uniform vec2 phi_bounds;
+uniform vec2 r_bounds;
+
+const float PI = 3.1415926535897932384626433832795;
+const float PI_2 = 1.57079632679489661923;
+const float PI_4 = 0.785398163397448309616;
+
+void ortho_find(inout vec3 v1, out vec3 v2, out vec3 v3) {
+    // Normalize
+    v1 = normalize(v1);
+    if (v1.z != 0) {
+        v2.xyz = vec3(1.0, 0.0, -(v1.x / v1.z));
+    } else if (v1.y != 0) {
+        v2.xyz = vec3(0.0, -(v1.z / v1.y), 1.0);
+    } else {
+        v2.xyz = vec3((-v1.y/v1.x), 1.0, 0.0);
+    }
+    v2 = normalize(v2);
+    v3.xyz = cross(v1, v2);
+}