Add image accumulation

All shader-related code has been removed; GPU entry points will be added once there's a better understanding of how the API should work.

.gitea/workflows/ci.yml

@@ -1,7 +1,6 @@
 name: "CI"
 on:
   push:
-  pull_request:
   workflow_dispatch:
 
 jobs:
@@ -22,7 +21,10 @@ jobs:
     steps:
       - uses: actions/checkout@v6
       - uses: actions-rust-lang/setup-rust-toolchain@v1
-      - run: cargo check
+        with:
+          components: clippy
+      - run: cargo check --all-targets
+      - run: cargo clippy --all-targets
       - run: cargo test
 
   test-gpu:
@@ -31,5 +33,5 @@ jobs:
     steps:
       - uses: actions/checkout@v6
       - uses: actions-rust-lang/setup-rust-toolchain@v1
-      - run: cargo install --git https://github.com/rust-gpu/rust-gpu cargo-gpu
-      - run: cargo gpu check -p enkou-shaders
+      - run: cargo install --git https://github.com/rust-gpu/rust-gpu cargo-gpu --rev 67f1ff2
+      - run: cargo gpu check --auto-install-rust-toolchain -p enkou-shaders

Cargo.toml

@@ -1,7 +1,6 @@
 [workspace]
 members = [
     "enkou-shaders",
-    "enkou-shaders-tests",
 ]
 resolver = "3"
 
@@ -13,6 +12,7 @@ license = "MIT"
 repository = ""
 
 [workspace.lints.rust]
+missing_docs = { level = "warn" }
 unexpected_cfgs = { level = "allow", check-cfg = ['cfg(target_arch, values("spirv"))'] }
 
 [workspace.dependencies]
@@ -21,6 +21,9 @@ spirv-std = { git = "https://github.com/Rust-GPU/rust-gpu.git", rev = "67f1ff2"
 
 anyhow = "1.0.102"
 bytemuck = { version = "1.25.0", features = ["derive"] }
-glam = { version = "0.33.1", default-features = false, features = ["libm"] }
-rspirv = "0.13.0"
-
+glam = { version = "0.33.1", default-features = false, features = ["bytemuck", "scalar-math"] }
+image = { version = "0.25.10", default-features = false, features = ["default-formats"]}
+libm = "0.2.16"
+rand = { version = "0.10.1", default-features = false }
+rand_xoshiro = "0.8.1"
+tempfile = "3.27.0"

enkou-shaders-tests/Cargo.toml

@@ -1,18 +0,0 @@
-[package]
-name = "enkou-shaders-tests"
-publish = false
-version.workspace = true
-authors.workspace = true
-edition.workspace = true
-license.workspace = true
-repository.workspace = true
-
-[lints]
-workspace = true
-
-[dependencies]
-rspirv.workspace = true
-
-[build-dependencies]
-anyhow.workspace = true
-cargo-gpu-install.workspace = true

enkou-shaders-tests/build.rs

@@ -1,26 +0,0 @@
-use cargo_gpu_install::install::Install;
-use cargo_gpu_install::spirv_builder::{ShaderPanicStrategy, SpirvMetadata};
-use std::path::PathBuf;
-
-pub fn main() -> anyhow::Result<()> {
-    let manifest_dir = env!("CARGO_MANIFEST_DIR");
-    let crate_path = [manifest_dir, "..", "enkou-shaders"]
-        .iter()
-        .copied()
-        .collect::<PathBuf>();
-
-    let mut install = Install::from_shader_crate(crate_path.clone());
-    install.build_script = true;
-    install.auto_install_rust_toolchain = true;
-    let install = install.run()?;
-
-    let mut builder = install.to_spirv_builder(crate_path, "spirv-unknown-vulkan1.3");
-    builder.build_script.defaults = true;
-    builder.shader_panic_strategy = ShaderPanicStrategy::SilentExit;
-    builder.spirv_metadata = SpirvMetadata::Full;
-
-    let compile_result = builder.build()?;
-    let spv_path = compile_result.module.unwrap_single();
-    println!("cargo::rustc-env=SHADER_SPV_PATH={}", spv_path.display());
-    Ok(())
-}

enkou-shaders-tests/src/lib.rs

@@ -1,67 +0,0 @@
-#[cfg(test)]
-mod test {
-    use rspirv::binary::parse_bytes;
-    use rspirv::dr::{Module, Operand};
-    use rspirv::spirv::ExecutionModel;
-    use std::sync::OnceLock;
-
-    static SHADER_MODULE: OnceLock<Module> = OnceLock::new();
-
-    fn shader() -> &'static Module {
-        SHADER_MODULE.get_or_init(|| {
-            let shader_bytes = include_bytes!(env!("SHADER_SPV_PATH"));
-
-            let mut loader = rspirv::dr::Loader::new();
-            parse_bytes(shader_bytes, &mut loader).expect("Unable to parse shader");
-
-            loader.module()
-        })
-    }
-
-    fn has_entry_point(execution_model: ExecutionModel, name: &str) -> bool {
-        for ref entry_point in shader().entry_points.iter() {
-            let operands: Vec<Operand> = entry_point
-                .operands
-                .iter()
-                .filter(|op| match op {
-                    Operand::ExecutionModel(_) | Operand::LiteralString(_) => true,
-                    _ => false,
-                })
-                .map(|op| op.clone())
-                .collect();
-
-            assert_eq!(operands.len(), 2);
-            match &operands[0] {
-                Operand::ExecutionModel(actual) => {
-                    if execution_model != *actual {
-                        continue;
-                    }
-                }
-                op => panic!("Unexpected operand; {}", op),
-            }
-
-            match &operands[1] {
-                Operand::LiteralString(actual) => {
-                    if name != actual {
-                        continue;
-                    }
-                }
-                op => panic!("Unexpected operand; {}", op),
-            }
-
-            return true;
-        }
-
-        false
-    }
-
-    #[test]
-    pub fn has_entry_main_fs() {
-        assert!(has_entry_point(ExecutionModel::Fragment, "main_fs"))
-    }
-
-    #[test]
-    pub fn has_entry_main_vs() {
-        assert!(has_entry_point(ExecutionModel::Vertex, "main_vs"))
-    }
-}

enkou-shaders/Cargo.toml

@@ -10,6 +10,14 @@ repository.workspace = true
 workspace = true
 
 [dependencies]
-spirv-std.workspace = true
-glam.workspace = true
 bytemuck.workspace = true
+glam.workspace = true
+libm.workspace = true
+rand.workspace = true
+rand_xoshiro.workspace = true
+spirv-std.workspace = true
+
+[dev-dependencies]
+anyhow.workspace = true
+image.workspace = true
+tempfile.workspace = true

enkou-shaders/examples/flame.rs

@@ -0,0 +1,106 @@
+//! # Fractal flame
+//!
+//! The `flam3` sources include example parameters (`test.flam3`) to demonstrate the fractal flame algorithm,
+//! which are copied here and used to generate an image that can be compared against the reference `flam3` renderer
+
+use anyhow::Result;
+use enkou_shaders::chaos_game::ChaosGame;
+use enkou_shaders::transform::Transform;
+use enkou_shaders::variation::{Variation, VariationKind};
+use glam::{UVec2, UVec3, UVec4, Vec2, Vec2Swizzles, Vec4};
+use rand::SeedableRng;
+use rand_xoshiro::Xoshiro256StarStar;
+use enkou_shaders::camera::Camera;
+use enkou_shaders::image::BlendMode;
+
+const ITERATIONS_DISCARD: usize = 100;
+const IMAGE_DIMENSIONS: UVec2 = UVec2::new(640, 480);
+const IMAGE_QUALITY: f32 = 10.0;
+
+// Palette data copied from `flam3-palettes.xml`.
+// Palettes are ASCII hexadecimal strings; two characters form a byte,
+// and bytes are grouped in triplets to get a single RGB color.
+const PALETTE: &'static str = "\
+00dadebc00eee6c500eef2ce00eef2cf00e6eee100eaeed800f2f1eb00f2f5d8
+00e6f2ce00deeac500d6dac600ced2bc00c2caa900becaa000ced6aa00dee2c5
+00eaedce00eaf2c500dee2c500c2caaa00aebeaa00a5b29600a2a98d0096a284
+008d8d7a0085897100858d710085856700797d6700797d670071795e00656d55
+004d5d4200344025003040250030381c002c3c1c002c341c00242c1200242400
+00242c0900283409003840120030401c0040502f0055694200657555006c7d5e
+00748d710074898400748d8400788d840079897a0079857100757d670071795e
+006c715e006d705e006c795e0068755e00697155006d7555006d755500697155
+0065715500696d550064715e006870670068706700686c67006c6c5e0071715e
+0079796700818571007d91710085927a0085927a007d9284007992840078928d
+00788d8d00748d84007492840075927a006c85670064795e0059694b00aa5700
+0038441c00303c1c002c3c1c003440250050614b005d6d5e0064715e0060715e
+0060755e0068755e006c795e006c795e0071796700707967006c7d6700687967
+006c7967006c75670071755e0071755e0075795e00757d5e00818d5e008d925e
+008d9267009a9a71009aa27a009aa27a009aa17a00929a71008992670081855e
+007d7d550069794b00616d4200444c250038441c0040512500454d2500716d42
+00797d4b00817d5500797955006d755500697d55006c795e006579540068795e
+00647967006479670068755e0064715e00646c5e00656d55004d584200344025
+002c381c0020281c001c1409001818000004140000081000000c1800001c2809
+00243012003c4425005d6555007579550085895e008991710096a271009aa27a
+009eaa7a009eaa7a00aaae7100a6aa7a00a2aa7a00a1a57a00969e7a0085967a
+0081927a0078927a0075927a00758d7a00708167007d7d670089896700929a71
+009eaa7a00aab68400b2b68d00b6ba9700c2ca9700b2be8d00b2b68d00aab28d
+00a2ae84009aa67a00929e7a00859a7a007d967a007d927a007d9284007d9284
+0081968400859684008596840081928400859a8400859a84008d9a8400929684
+009ea98400aeb28400aaba8400b2be8d00b6c2a000c6caa000c6ceaa00d6dab3
+00dae2c500d2d6bc00bec2a000aab68d009ea67a00929a710089897100817d67
+007d7d6700817867007d7d5e0079795e0079815e00817d6700817d6700818167
+008189710085917a0089927a00969d7a00969e7a0092968400969a8d00929284
+0089918400819284007d928d0078928d0074928d0078928d007896970081968d
+0081968d00819a8d00859a8d00899e8d00899e8d008da2970095a297008da297
+0096a68d009aa18d009ea984009ea67a00a2a571009ea671009aa67100959d71";
+
+/// Generate an image
+pub fn main() -> Result<()> {
+    let transforms: &[Transform] = &[];
+
+    let weights = [0.25; 4];
+
+    let variations = [Variation::new(
+        VariationKind::Spherical,
+        1.0,
+        [0.0; 4].into(),
+    )];
+
+    let camera = Camera::new(IMAGE_DIMENSIONS, Vec2::ZERO, 0.0, Vec2::ZERO, (IMAGE_DIMENSIONS / 2).as_vec2());
+
+    let palette_chars: String = PALETTE.chars().filter(|c| c.is_alphanumeric()).collect();
+    let palette_bytes: Vec<u8> = (0..palette_chars.len())
+        .step_by(2)
+        .map(|i| {
+            let s = &palette_chars[i..i + 2];
+            let b = u8::from_str_radix(s, 16).expect("illegal palette value");
+            b
+        })
+        .collect();
+    let palette_colors: Vec<Vec4> = palette_bytes
+        .chunks(3)
+        .map(|chunk| [chunk[0] as f32, chunk[1] as f32, chunk[2] as f32, 255.0].into())
+        .collect();
+
+    let mut rng = Xoshiro256StarStar::from_seed([4u8; 32]);
+    let chaos_game = ChaosGame::new(&mut rng, transforms, &weights, &variations);
+
+    let mut image_accum: Vec<Vec4> = Vec::new();
+    image_accum.resize(IMAGE_DIMENSIONS.element_product() as usize, Vec4::ZERO);
+
+    let iterations: usize = (IMAGE_DIMENSIONS.element_product() as f32 * IMAGE_QUALITY) as usize;
+    chaos_game.skip(ITERATIONS_DISCARD).take(iterations).for_each(|(ifs_point, ifs_color)| {
+        let pixel_point = camera.transform_point_to_image(ifs_point);
+        if pixel_point.is_none() {
+            return;
+        }
+
+        let pixel_point = pixel_point.unwrap();
+        let pixel_index = pixel_point.y * IMAGE_DIMENSIONS.x + pixel_point.x;
+
+        let pixel_color = BlendMode::Linear.ifs_to_rgb(ifs_color, palette_colors.as_ref());
+        image_accum[pixel_index as usize] += pixel_color;
+    });
+
+    Ok(())
+}

enkou-shaders/examples/gasket.rs

@@ -0,0 +1,86 @@
+//! # Sierpinski Gasket
+//!
+//! The Sierpinski Gasket is a simple iterated function system (IFS) that demonstrates
+//! the basic principles of using the "chaos game" to generate images. Using a set of three
+//! transforms, generate new points in the function system using the chaos game,
+//! then plot those points and display the resulting image.
+use anyhow::{Context, Result};
+use enkou_shaders::Coefficients2;
+use enkou_shaders::camera::Camera;
+use enkou_shaders::chaos_game::ChaosGame;
+use enkou_shaders::transform::Transform;
+use enkou_shaders::variation::Variation;
+use glam::{Affine2, UVec2, Vec2, uvec2};
+use image::{GrayImage, Luma};
+use rand::SeedableRng;
+use rand_xoshiro::Xoshiro256StarStar;
+use std::mem;
+use std::process::Command;
+use tempfile::NamedTempFile;
+
+const ITERATIONS_DISCARD: usize = 20;
+const IMAGE_DIMENSION: UVec2 = uvec2(600, 600);
+const IMAGE_QUALITY: f32 = 1.0;
+
+/// Build and display a simple fractal - the Sierpinski Gasket
+pub fn main() -> Result<()> {
+    let mut rng = Xoshiro256StarStar::from_seed([4u8; 32]);
+
+    let transforms = [
+        {
+            // F_0: (x / 2, y / 2)
+            let coefficients = Affine2::from_coefficients(0.5, 0.0, 0.0, 0.0, 0.5, 0.0);
+            Transform::new(coefficients, Affine2::IDENTITY, uvec2(0, 1), 0.0, 0.5)
+        },
+        {
+            // F_1: ((x + 1) / 2, y / 2)
+            let coefficients = Affine2::from_coefficients(0.5, 0.0, 0.5, 0.0, 0.5, 0.0);
+            Transform::new(coefficients, Affine2::IDENTITY, uvec2(0, 1), 0.5, 0.5)
+        },
+        {
+            // F_2: (x / 2, (y + 1) / 2)
+            let coefficients = Affine2::from_coefficients(0.5, 0.0, 0.0, 0.0, 0.5, 0.5);
+            Transform::new(coefficients, Affine2::IDENTITY, uvec2(0, 1), 1.0, 0.5)
+        },
+    ];
+    let weights = [1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0];
+    let variations = [Variation::IDENTITY];
+
+    // The gasket is defined on the range [0, 1] for both X and Y
+    let camera = Camera::new(
+        IMAGE_DIMENSION,
+        Vec2::ONE * 0.5,
+        0.0,
+        Vec2::ZERO,
+        IMAGE_DIMENSION.as_vec2(),
+    );
+
+    let mut image = GrayImage::new(IMAGE_DIMENSION.x, IMAGE_DIMENSION.y);
+    let chaos_game = ChaosGame::new(&mut rng, &transforms, &weights, &variations);
+
+    let iterations = (IMAGE_DIMENSION.element_product() as f32 * IMAGE_QUALITY) as usize;
+    chaos_game
+        .skip(ITERATIONS_DISCARD)
+        .take(iterations)
+        .for_each(|(ifs_point, _)| {
+            if let Some(pixel_point) = camera.transform_point_to_image(ifs_point) {
+                image.put_pixel(pixel_point.x, pixel_point.y, Luma([255u8]));
+            }
+        });
+
+    let temp = NamedTempFile::with_suffix(".png").context("Unable to create file for image")?;
+    image.save(temp.path()).context("Unable to save image")?;
+
+    let mut command = cfg_select! {
+        unix => Command::new("xdg-open").arg(temp.path()).spawn(),
+        windows => Command::new("PowerShell").arg("-Command").arg(format!("start {}", temp.path().display())).spawn(),
+        _ => Err(anyhow::anyhow!("No available program to open images"))?
+    }?;
+    command.wait()?;
+
+    // In case the image viewer forks and gives control back prior to reading the file,
+    // drop it and don't run the destructor
+    mem::forget(temp);
+
+    Ok(())
+}

enkou-shaders/src/camera.rs

@@ -0,0 +1,201 @@
+//! # Camera
+//!
+//! Map points from the IFS coordinate system to pixel coordinates. This is a lossy transformation.
+use bytemuck::{Pod, Zeroable};
+use glam::{Affine2, IVec2, UVec2, Vec2, vec2};
+use libm::powf;
+
+/// Settings used to map IFS coordinates to pixel coordinates.
+///
+/// The camera is itself an affine transformation, capable of zoom, rotation, and translation
+/// of the IFS coordinates before rendering to the final image.
+#[derive(Copy, Clone, Pod, Zeroable)]
+#[repr(C)]
+pub struct Camera {
+    dimensions: UVec2,
+    transform: Affine2,
+}
+
+impl Camera {
+    /// Construct a new camera for translating IFS coordinates to pixel coordinates.
+    ///
+    /// While the camera is implemented as a single affine transformation, it's helpful
+    /// to express the transform steps individually.
+    ///
+    /// # Arguments
+    ///
+    /// * `dimensions` - Width and height of the output image (in pixels).
+    /// * `center` - Location of the origin in IFS coordinates. Positive `x` shifts the image
+    ///   left, and positive `y` position shifts the image up.
+    /// * `rotate` - Rotation angle (in radians) of IFS coordinates. Rotation is applied after the
+    ///   `center` translation, so it is about the new origin.
+    /// * `zoom` - Zoom factor applied to IFS coordinates. IFS coordinates are scaled by
+    ///   `pow(2, zoom)`, so a zoom factor of 0 is the identity.
+    /// * `scale` - Pixels per unit of IFS coordinates. This parameter is usually chosen such
+    ///   that the largest dimension will cover the range `[-2, 2]`, but values higher or lower
+    ///   can be used as a secondary zoom.
+    pub fn new(dimensions: UVec2, center: Vec2, rotate: f32, zoom: Vec2, scale: Vec2) -> Camera {
+        let ifs_center_transform = Affine2::from_translation(-center);
+        let zoom_transform = Affine2::from_scale(vec2(powf(2.0, zoom.x), powf(2.0, zoom.y)));
+        let scale_transform = Affine2::from_scale(scale);
+        let rotate_transform = Affine2::from_angle(rotate);
+        let image_center_transform = Affine2::from_translation((dimensions / 2).as_vec2());
+
+        let transform = image_center_transform
+            * rotate_transform
+            * scale_transform
+            * zoom_transform
+            * ifs_center_transform;
+
+        Camera {
+            dimensions,
+            transform,
+        }
+    }
+
+    /// Map a point from IFS coordinates to pixel coordinates.
+    ///
+    /// ```
+    /// # use glam::{vec2, ivec2, uvec2, Vec2};
+    /// # use crate::enkou_shaders::camera::Camera;
+    /// // Output image is 600x600 pixels, centered at the origin, no rotation, no zoom,
+    /// // and scaled such that it covers the range [-2, 2].
+    /// // Use the origin as the IFS coordinate, so the pixel coordinate is the center of the image
+    /// let camera = Camera::new(
+    ///     uvec2(600, 600),
+    ///     Vec2::ZERO,
+    ///     0.0,
+    ///     Vec2::ZERO,
+    ///     vec2(150.0, 150.0)
+    /// );
+    /// assert_eq!(camera.transform_point(vec2(0.0, 0.0)), ivec2(300, 300));
+    /// ```
+    pub fn transform_point(&self, point: Vec2) -> IVec2 {
+        self.transform.transform_point2(point).as_ivec2()
+    }
+
+    /// Map a point from IFS coordinates to pixel coordinates (like [`transform_point`](Camera::transform_point)),
+    /// and check that the result is within the provided image dimensions.
+    pub fn transform_point_to_image(&self, point: Vec2) -> Option<UVec2> {
+        let pixel_coordinates = self.transform_point(point);
+        if pixel_coordinates.x < 0
+            || pixel_coordinates.y < 0
+            || (pixel_coordinates.x as u32) >= self.dimensions.x
+            || (pixel_coordinates.y as u32) >= self.dimensions.y
+        {
+            None
+        } else {
+            Some(pixel_coordinates.as_uvec2())
+        }
+    }
+}
+
+/// Shader entry point for running the camera transformation over a list of IFS coordinates
+pub mod entry {
+    use crate::camera::Camera;
+    use spirv_std::glam::{IVec2, Vec2};
+    use spirv_std::spirv;
+
+    /// Transform IFS coordinates to pixel coordinates
+    #[spirv(compute(entry_point_name = "main_camera", threads(1)))]
+    pub fn main_camera(
+        #[spirv(storage_buffer, descriptor_set = 0, binding = 0)] camera: &Camera,
+        #[spirv(storage_buffer, descriptor_set = 0, binding = 1)] coordinates_ifs: &[Vec2],
+        #[spirv(storage_buffer, descriptor_set = 1, binding = 0)] coordinates_pixel: &mut [IVec2],
+    ) {
+        for i in 0..coordinates_ifs.len() {
+            coordinates_pixel[i] = camera.transform_point(coordinates_ifs[i])
+        }
+    }
+}
+
+#[cfg(test)]
+mod test {
+    use crate::camera::Camera;
+    use glam::{Affine2, Vec2, ivec2, uvec2, vec2};
+    use libm::powf;
+
+    #[test]
+    pub fn manual_camera() {
+        let starting_point = vec2(1.0, 1.0);
+
+        // Move the origin; points move right and up by one unit, giving us (2.0, 2.0)
+        let center = vec2(-1.0, -1.0);
+        let point = starting_point - center;
+
+        // Rotate about the new origin; points move counter-clockwise, giving us (-2.0, 2.0)
+        let rotate = 90.0f32.to_radians();
+        let point = Affine2::from_angle(rotate).transform_point2(point);
+
+        // Zoom in by a factor of 1; points will be twice as far from the origin,
+        // giving us (-4.0, 4.0)
+        let zoom = vec2(1.0, 1.0);
+        let point = point * vec2(powf(2.0, zoom.x), powf(2.0, zoom.y));
+
+        // Apply scaling; scale 100 in a 1000 x 1000 image is an effective range
+        // of [-5, 5] in IFS coordinates.
+        // After scaling, the point is (-400.0, 400.0)
+        let scale = vec2(100.0, 100.0);
+        let point = point * scale;
+
+        // Move the origin from (0, 0) to image center,
+        // giving us (100.0, 900.0)
+        let dimensions = uvec2(1000, 1000);
+        let point = point.as_ivec2() + dimensions.as_ivec2() / 2;
+
+        // Check that the camera implementation ends up at the same point
+        let camera = Camera::new(dimensions, center, rotate, zoom, scale);
+
+        // The camera is implemented by composing affine transforms,
+        // which ends up with a slightly different result because of rounding.
+        let error = camera.transform_point(starting_point) - point;
+        assert!(error.x.abs() <= 1);
+        assert!(error.y.abs() <= 1);
+    }
+
+    #[test]
+    pub fn point_outside_camera() {
+        // Scale 250 for an image 1000 x 1000 gives an effective range of [-2, 2]
+        let camera = Camera::new(
+            uvec2(1000, 1000),
+            Vec2::ZERO,
+            0.0,
+            Vec2::ZERO,
+            vec2(250.0, 250.0),
+        );
+
+        // Converting a point outside the effective range is legal, but outside the image bounds
+        assert_eq!(camera.transform_point(vec2(3.0, 3.0)), ivec2(1250, 1250));
+    }
+
+    #[test]
+    pub fn point_outside_camera_negative() {
+        // Scale 250 for an image 1000 x 1000 gives an effective range of [-2, 2]
+        let camera = Camera::new(
+            uvec2(1000, 1000),
+            Vec2::ZERO,
+            0.0,
+            Vec2::ZERO,
+            vec2(250.0, 250.0),
+        );
+
+        // Converting a point outside the effective range is legal, but outside the image bounds
+        assert_eq!(camera.transform_point(vec2(-3.0, -3.0)), ivec2(-250, -250));
+    }
+
+    #[test]
+    pub fn aspect_ratio() {
+        // Scale 100 for an image 1600 x 900 gives an effective X range of [-8, 8],
+        // and effective Y range of [-4.5, 4.5]
+        let camera = Camera::new(
+            uvec2(1600, 900),
+            Vec2::ZERO,
+            0.0,
+            Vec2::ZERO,
+            vec2(100.0, 100.0),
+        );
+
+        // This point is inside the image width, but outside its height
+        assert_eq!(camera.transform_point(vec2(6.0, 6.0)), ivec2(1400, 1050));
+    }
+}

enkou-shaders/src/chaos_game.rs

@@ -0,0 +1,115 @@
+//! # Chaos Game
+//!
+//! Fractal flames are a class of
+//! [iterated function systems](https://en.wikipedia.org/wiki/Iterated_function_system)
+//! that generate images following a simple algorithm:
+//!
+//! - Pick a starting point `(x, y)`
+//! - Iterate:
+//!     - Pick a [`Transform`] from the set of available transforms
+//!     - Apply the current point to the chosen transform, generating a new point `(x, y)`
+//!     - Plot the new point `(x, y)`
+//!
+//! This algorithm is also known as the ["chaos game"](https://en.wikipedia.org/wiki/Chaos_game),
+//! and it forms the basic system for producing images.
+
+use crate::transform::Transform;
+use crate::variation::Variation;
+use rand::distr::{Distribution, StandardUniform};
+use rand::{Rng, RngExt};
+use spirv_std::glam::{Vec2, vec2};
+
+struct BiUnit;
+impl Distribution<f32> for BiUnit {
+    fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> f32 {
+        rng.sample::<f32, _>(StandardUniform) * 2.0 - 1.0
+    }
+}
+
+/// Iterate one step in the chaos game; choose the next transform, apply it,
+/// and return the resulting point. Also returns the transform index so that
+/// path-dependent weights (the "Xaos" table in Apophysis) can be chosen
+/// for the next iteration step.
+///
+/// # Arguments
+///
+/// * `weights` - Weights are assumed to be normalized; adding all elements together should return the value 1
+pub fn step_chaos_game<R: Rng>(
+    point: Vec2,
+    color: f32,
+    rng: &mut R,
+    transforms: &[Transform],
+    weights: &[f32],
+    variations: &[Variation],
+) -> (Vec2, f32) {
+    let mut choice_weight = rng.sample::<f32, _>(StandardUniform);
+    let mut transform_index: u32 = 0;
+
+    for i in 0..weights.len() {
+        choice_weight -= weights[i];
+        if choice_weight <= 0.0 {
+            break;
+        }
+
+        transform_index += 1;
+    }
+
+    let transform = &transforms[transform_index as usize];
+    (
+        transform.transform_point(rng, variations, point),
+        transform.transform_color(color),
+    )
+}
+
+/// Iterator for chaos game state. Holds the current point and references to all other data
+/// necessary to generate fractal flame images.
+///
+/// New points in the chaos game are produced by iterating on the chaos game.
+pub struct ChaosGame<'a, R: Rng> {
+    current_point: Vec2,
+    current_color: f32,
+    rng: &'a mut R,
+    transforms: &'a [Transform],
+    weights: &'a [f32],
+    variations: &'a [Variation],
+}
+
+impl<'a, R: Rng> ChaosGame<'a, R> {
+    /// Create a new chaos game iterator
+    pub fn new(
+        rng: &'a mut R,
+        transforms: &'a [Transform],
+        weights: &'a [f32],
+        variations: &'a [Variation],
+    ) -> Self {
+        let current_point = vec2(rng.sample(BiUnit), rng.sample(BiUnit));
+        let current_color = rng.sample(StandardUniform);
+        ChaosGame {
+            current_point,
+            current_color,
+            rng,
+            transforms,
+            weights,
+            variations,
+        }
+    }
+}
+
+impl<'a, R: Rng> Iterator for ChaosGame<'a, R> {
+    type Item = (Vec2, f32);
+
+    fn next(&mut self) -> Option<Self::Item> {
+        let (next_point, next_color) = step_chaos_game(
+            self.current_point,
+            self.current_color,
+            self.rng,
+            self.transforms,
+            self.weights,
+            self.variations,
+        );
+        self.current_point = next_point;
+        self.current_color = next_color;
+
+        Some((next_point, next_color))
+    }
+}

enkou-shaders/src/image.rs

@@ -0,0 +1,84 @@
+//! Image
+use glam::Vec4;
+use libm::floorf;
+
+/// Blending modes for mapping IFS color values (which are on a scale `[0, 1]`)
+/// to RGBA colors.
+#[derive(Copy, Clone, Default)]
+#[repr(u32)]
+pub enum BlendMode {
+    /// Map IFS color values to a linear blend of the nearest two palette colors
+    #[default]
+    Linear = 0,
+
+    /// Map IFS color values to the nearest single palette color
+    Step = 1,
+}
+
+impl BlendMode {
+    /// Map an IFS color value to RGBA color from the provided palette.
+    pub fn ifs_to_rgb(&self, color: f32, palette: &[Vec4]) -> Vec4 {
+        let colors_m_one = palette.len() - 1;
+        let period = 1.0 / colors_m_one as f32;
+        let index_lower = floorf(color / period) as usize;
+        let index_upper = (index_lower + 1).clamp(0, colors_m_one);
+        let rem = color % period / period;
+
+        match self {
+            BlendMode::Linear => palette[index_lower].lerp(palette[index_upper], rem),
+            BlendMode::Step => palette[index_lower],
+        }
+    }
+}
+
+// UNSAFE: Sound because enum has guaranteed layout (u32) and defined zero-value
+unsafe impl bytemuck::Zeroable for BlendMode {}
+// UNSAFE: Sound because enum has guaranteed layout (u32) and defined zero-value
+unsafe impl bytemuck::Pod for BlendMode {}
+
+#[cfg(test)]
+mod test {
+    use crate::image::BlendMode;
+    use glam::Vec4;
+
+    #[test]
+    fn blend_linear() {
+        let ifs_to_rgb = |color, palette| BlendMode::Linear.ifs_to_rgb(color, palette);
+
+        let palette = &[Vec4::splat(0.0), Vec4::splat(1.0)];
+        assert_eq!(ifs_to_rgb(0.0, palette), Vec4::splat(0.0));
+        assert_eq!(ifs_to_rgb(0.5, palette), Vec4::splat(0.5));
+        assert_eq!(ifs_to_rgb(1.0, palette), Vec4::splat(1.0));
+
+        let palette = &[Vec4::splat(1.0), Vec4::splat(2.0), Vec4::splat(3.0)];
+        assert_eq!(ifs_to_rgb(0.0, palette), Vec4::splat(1.0));
+        assert_eq!(ifs_to_rgb(0.5, palette), Vec4::splat(2.0));
+        assert_eq!(ifs_to_rgb(1.0, palette), Vec4::splat(3.0));
+
+        let palette = &[
+            Vec4::splat(1.0),
+            Vec4::splat(2.0),
+            Vec4::splat(3.0),
+            Vec4::splat(4.0),
+        ];
+        assert_eq!(ifs_to_rgb(0.0, palette), Vec4::splat(1.0));
+        assert_eq!(ifs_to_rgb(0.5, palette), Vec4::splat(2.5));
+        assert_eq!(ifs_to_rgb(1.0, palette), Vec4::splat(4.0));
+    }
+
+    #[test]
+    fn blend_step() {
+        let ifs_to_rgb = |color, palette| BlendMode::Step.ifs_to_rgb(color, palette);
+
+        let palette = &[Vec4::splat(0.0), Vec4::splat(1.0)];
+        assert_eq!(ifs_to_rgb(0.5, palette), Vec4::splat(0.0));
+
+        let palette = &[Vec4::splat(1.0), Vec4::splat(2.0), Vec4::splat(3.0)];
+        assert_eq!(ifs_to_rgb(0.0, palette), palette[0]);
+        assert_eq!(ifs_to_rgb(0.25, palette), palette[0]);
+        assert_eq!(ifs_to_rgb(0.4, palette), palette[0]);
+        assert_eq!(ifs_to_rgb(0.5, palette), palette[1]);
+        assert_eq!(ifs_to_rgb(0.7, palette), palette[1]);
+        assert_eq!(ifs_to_rgb(1.0, palette), palette[2]);
+    }
+}

enkou-shaders/src/lib.rs

@@ -1,36 +1,109 @@
+//! # Enkou
 #![no_std]
+#![allow(clippy::needless_range_loop)] // SPIR-V backend has issues with iteration over items
 
-use bytemuck::{Pod, Zeroable};
-use core::f32::consts::PI;
-use glam::{Vec3, Vec4, vec2, vec3};
-#[cfg(target_arch = "spirv")]
-use spirv_std::num_traits::Float;
-use spirv_std::spirv;
+pub mod camera;
+pub mod chaos_game;
+pub mod image;
+pub mod transform;
+pub mod variation;
 
-#[derive(Copy, Clone, Pod, Zeroable)]
-#[repr(C)]
-pub struct ShaderConstants {
-    pub width: u32,
-    pub height: u32,
-    pub time: f32,
+use glam::Affine2;
+
+/// Utility trait to convert between `flam3` notation and [`glam`].
+#[allow(missing_docs)]
+pub trait Coefficients2 {
+    /// Convert affine transformation coefficients to the [`glam`] representation.
+    /// Parameters use the following form:
+    ///
+    /// ```text
+    /// (a * x + b * y + c, d * x + e * y + f)
+    /// ```
+    ///
+    /// ```
+    /// # use glam::{Affine2, vec2};
+    /// # use crate::enkou_shaders::Coefficients2;
+    /// let coefs = Affine2::from_coefficients(1.0, 2.0, 3.0, 4.0, 5.0, 6.0);
+    /// let (x, y) = (7.0, 8.0);
+    /// assert_eq!(
+    ///     coefs.transform_point2(vec2(x, y)),
+    ///     vec2(
+    ///         coefs.a() * x + coefs.b() * y + coefs.c(),
+    ///         coefs.d() * x + coefs.e() * y + coefs.f()
+    ///     )
+    /// );
+    /// ```
+    fn from_coefficients(a: f32, b: f32, c: f32, d: f32, e: f32, f: f32) -> Affine2;
+
+    /// Convert affine transformation coefficients to the [`glam`] representation.
+    /// Parameters use the following form:
+    ///
+    /// ```text
+    /// (a * x + b * y + c, d * x + e * y + f)
+    /// ```
+    ///
+    /// ```
+    /// # use glam::{Affine2, vec2};
+    /// # use crate::enkou_shaders::Coefficients2;
+    /// let coefs = Affine2::from_coefficients_arr([1.0, 2.0, 3.0, 4.0, 5.0, 6.0]);
+    /// let (x, y) = (7.0, 8.0);
+    /// assert_eq!(
+    ///     coefs.transform_point2(vec2(x, y)),
+    ///     vec2(
+    ///         coefs.a() * x + coefs.b() * y + coefs.c(),
+    ///         coefs.d() * x + coefs.e() * y + coefs.f()
+    ///     )
+    /// );
+    /// ```
+    fn from_coefficients_arr(coefficients: [f32; 6]) -> Affine2;
+
+    fn a(&self) -> f32;
+    fn b(&self) -> f32;
+    fn c(&self) -> f32;
+    fn d(&self) -> f32;
+    fn e(&self) -> f32;
+    fn f(&self) -> f32;
 }
 
-#[spirv(fragment)]
-pub fn main_fs(vtx_color: Vec3, output: &mut Vec4) {
-    *output = Vec4::from((vtx_color, 1.));
+impl Coefficients2 for Affine2 {
+    #[inline]
+    fn from_coefficients(a: f32, b: f32, c: f32, d: f32, e: f32, f: f32) -> Affine2 {
+        Affine2::from_cols_array(&[a, d, b, e, c, f])
     }
 
-#[spirv(vertex)]
-pub fn main_vs(
-    #[spirv(vertex_index)] vert_id: i32,
-    #[spirv(descriptor_set = 0, binding = 0, storage_buffer)] constants: &ShaderConstants,
-    #[spirv(position)] vtx_pos: &mut Vec4,
-    vtx_color: &mut Vec3,
-) {
-    let speed = 0.4;
-    let time = constants.time * speed + vert_id as f32 * (2. * PI * 120. / 360.);
-    let position = vec2(f32::sin(time), f32::cos(time));
-    *vtx_pos = Vec4::from((position, 0.0, 1.0));
-
-    *vtx_color = [vec3(1., 0., 0.), vec3(0., 1., 0.), vec3(0., 0., 1.)][vert_id as usize % 3];
+    #[inline]
+    fn from_coefficients_arr(coefficients: [f32; 6]) -> Affine2 {
+        Affine2::from_coefficients(
+            coefficients[0],
+            coefficients[1],
+            coefficients[2],
+            coefficients[3],
+            coefficients[4],
+            coefficients[5],
+        )
+    }
+
+    fn a(&self) -> f32 {
+        self.matrix2.x_axis.x
+    }
+
+    fn b(&self) -> f32 {
+        self.matrix2.y_axis.x
+    }
+
+    fn c(&self) -> f32 {
+        self.translation.x
+    }
+
+    fn d(&self) -> f32 {
+        self.matrix2.x_axis.y
+    }
+
+    fn e(&self) -> f32 {
+        self.matrix2.y_axis.y
+    }
+
+    fn f(&self) -> f32 {
+        self.translation.y
+    }
 }

enkou-shaders/src/transform.rs

@@ -0,0 +1,138 @@
+//! # Transform
+//!
+//! Transforms are the "functions" in an iterated function system. They take in a point,
+//! and generate a new point. For fractal flames, transforms are always affine,
+//! but produce more interesting images once we add variations.
+use crate::variation::Variation;
+use bytemuck::{Pod, Zeroable};
+use glam::{Affine2, FloatExt, UVec2, Vec2};
+use rand::Rng;
+
+/// Affine transform for use in the [`chaos_game`](crate::chaos_game).
+#[derive(Copy, Clone, Pod, Zeroable)]
+#[repr(C)]
+pub struct Transform {
+    coefficients: Affine2,
+    coefficients_post: Affine2,
+    variation_range: UVec2,
+    color: f32,
+    color_speed: f32,
+}
+
+impl Transform {
+    /// Create a new transform from an affine transformation matrix
+    pub fn new(
+        coefficients: Affine2,
+        coefficients_post: Affine2,
+        variation_range: UVec2,
+        color: f32,
+        color_speed: f32,
+    ) -> Self {
+        Transform {
+            coefficients,
+            coefficients_post,
+            variation_range,
+            color,
+            color_speed,
+        }
+    }
+
+    /// Apply this transform to a point in IFS coordinates, producing a new point
+    pub fn transform_point<R: Rng>(
+        &self,
+        rng: &mut R,
+        variations: &[Variation],
+        point: Vec2,
+    ) -> Vec2 {
+        let point = self.coefficients.transform_point2(point);
+
+        let mut point_output = Vec2::ZERO;
+
+        let variation_start = self.variation_range.x;
+        let variation_end = self.variation_range.y;
+        for variation_index in variation_start..variation_end {
+            let variation = &variations[variation_index as usize];
+            point_output += variation.transform_point(point, rng, &self.coefficients)
+        }
+
+        self.coefficients_post.transform_point2(point)
+    }
+
+    /// Apply this transform to a color in IFS coordinates, producing a new color
+    pub fn transform_color(&self, color: f32) -> f32 {
+        self.color.lerp(color, self.color_speed)
+    }
+}
+
+#[cfg(test)]
+mod test {
+    use crate::transform::Transform;
+    use crate::variation::{Variation, VariationKind};
+    use glam::{Affine2, uvec2, vec2};
+    use rand::SeedableRng;
+    use rand_xoshiro::Xoshiro256StarStar;
+
+    #[test]
+    fn transform_scaling() {
+        let scale_coefficients = vec2(2.0, 0.5);
+        let transform = Transform::new(
+            Affine2::from_scale(scale_coefficients),
+            Affine2::IDENTITY,
+            uvec2(0, 1),
+            0.0,
+            0.0,
+        );
+
+        let mut rng = Xoshiro256StarStar::from_seed([0u8; 32]);
+        let variations = [Variation::IDENTITY];
+        let point = vec2(1.0, 1.0);
+
+        assert_eq!(
+            transform.transform_point(&mut rng, &variations, point),
+            scale_coefficients
+        );
+    }
+
+    #[test]
+    fn transform_scaling_post() {
+        let scale_coefficients = vec2(2.0, 0.5);
+        let transform_pdj =
+            Transform::new(Affine2::IDENTITY, Affine2::IDENTITY, uvec2(0, 1), 0.0, 0.0);
+        let transform_pdj_post = Transform::new(
+            Affine2::IDENTITY,
+            Affine2::from_scale(scale_coefficients),
+            uvec2(0, 1),
+            0.0,
+            0.0,
+        );
+
+        let mut rng = Xoshiro256StarStar::from_seed([0; 32]);
+        let variations = [Variation::new(VariationKind::Pdj, 1.0, [0.0f32; 4].into())];
+        let point = vec2(1.0, 1.0);
+
+        let point_pdj = transform_pdj.transform_point(&mut rng, &variations, point);
+        let point_pdj_post = transform_pdj_post.transform_point(&mut rng, &variations, point);
+
+        assert_eq!(point_pdj * scale_coefficients, point_pdj_post);
+    }
+
+    #[test]
+    fn transform_color() {
+        let starting_color = 0.0;
+        let transform_color = 1.0;
+        let transform_speed = 0.5;
+
+        let transform = Transform::new(
+            Affine2::IDENTITY,
+            Affine2::IDENTITY,
+            uvec2(0, 1),
+            starting_color,
+            transform_speed,
+        );
+
+        assert_eq!(
+            transform.transform_color(starting_color),
+            starting_color * (1.0 - transform_speed) + transform_color * transform_speed
+        );
+    }
+}

enkou-shaders/src/variation.rs

@@ -0,0 +1,139 @@
+//! # Variation
+//!
+//! Variations extend the fractal flame iterated function system
+//! with non-linear transforms (as opposed to [`Transform`]s,
+//! which are strictly affine transformations).
+use crate::Coefficients2;
+use bytemuck::{Pod, Zeroable};
+use core::f32::consts::PI;
+use glam::{Affine2, Vec2, vec2};
+use libm::{atan2f, cosf, powf, sinf, sqrtf, tanf};
+use rand::distr::StandardUniform;
+use rand::{Rng, RngExt};
+
+/// Generic variation parameters
+///
+/// Not all variations will use these parameters, but passing them
+/// as an array per variation allows shaders to use a consistent struct size
+/// no matter what the variation actually needs.
+#[derive(Copy, Clone, Pod, Zeroable)]
+#[repr(C)]
+pub struct VariationParams([f32; 4]);
+
+impl From<[f32; 4]> for VariationParams {
+    fn from(v: [f32; 4]) -> Self {
+        VariationParams(v)
+    }
+}
+
+/// Enum for all supported variation types
+///
+/// ID numbers are chosen to match the variation identifier also used by `flam3`
+#[derive(Copy, Clone)]
+#[repr(u32)]
+#[allow(missing_docs)]
+pub enum VariationKind {
+    /// Identity variation, returns the point as-is
+    Linear = 0,
+    Spherical = 2,
+    Julia = 13,
+    Popcorn = 17,
+    Pdj = 24,
+}
+
+// UNSAFE: Sound because enum has guaranteed layout (u32) and defined zero-value
+unsafe impl bytemuck::Zeroable for VariationKind {}
+// UNSAFE: Sound because enum has guaranteed layout (u32) and defined zero-value
+unsafe impl bytemuck::Pod for VariationKind {}
+
+/// Parameters required for shaders to run the variation function.
+///
+/// Not all variations use the [`VariationParams`], but using the struct
+/// makes it easy to provide parameters to the shader.
+#[derive(Copy, Clone, Pod, Zeroable)]
+#[repr(C)]
+pub struct Variation {
+    kind: VariationKind,
+    weight: f32,
+    params: VariationParams,
+}
+
+impl Variation {
+    /// Identity variation; calling [`transform_point`] will yield
+    /// the same point as the input.
+    pub const IDENTITY: Variation = Variation {
+        kind: VariationKind::Linear,
+        weight: 1.0,
+        params: VariationParams([0f32; 4]),
+    };
+
+    /// Create a new variation by providing the variation kind, weight, and parameters.
+    pub fn new(kind: VariationKind, weight: f32, params: VariationParams) -> Variation {
+        Variation {
+            kind,
+            weight,
+            params,
+        }
+    }
+
+    /// Transform a point by applying this variation.
+    ///
+    /// Output points are scaled by this variation's weight.
+    pub fn transform_point<R: Rng>(
+        &self,
+        point: Vec2,
+        rng: &mut R,
+        coefficients: &Affine2,
+    ) -> Vec2 {
+        (match self.kind {
+            VariationKind::Linear => transform_point_linear(point),
+            VariationKind::Spherical => transform_point_spherical(point),
+            VariationKind::Julia => transform_point_julia(point, rng),
+            VariationKind::Popcorn => transform_point_popcorn(point, coefficients),
+            VariationKind::Pdj => transform_point_pdj(point, &self.params),
+        }) * self.weight
+    }
+}
+
+fn transform_point_linear(point: Vec2) -> Vec2 {
+    point
+}
+
+fn transform_point_spherical(point: Vec2) -> Vec2 {
+    let r2 = (point * point).element_sum();
+    let r2_inv = 1.0 / r2;
+
+    point * Vec2::splat(r2_inv)
+}
+fn transform_point_julia<R: Rng>(point: Vec2, rng: &mut R) -> Vec2 {
+    let x2 = powf(point.x, 2.0);
+    let y2 = powf(point.y, 2.0);
+    let r = sqrtf(x2 + y2);
+
+    let theta = atan2f(point.x, point.y);
+    let omega = if rng.sample::<f32, _>(StandardUniform) > 0.5 {
+        PI
+    } else {
+        0.0
+    };
+
+    let sqrt_r = sqrtf(r);
+    let theta_val = theta / 2.0 + omega;
+
+    vec2(sqrt_r * cosf(theta_val), sqrt_r * sinf(theta_val))
+}
+
+fn transform_point_popcorn(point: Vec2, coefficients: &Affine2) -> Vec2 {
+    vec2(
+        point.x * coefficients.c() * sinf(tanf(3.0 * point.y)),
+        point.y + coefficients.f() * sinf(tanf(3.0 * point.x)),
+    )
+}
+
+fn transform_point_pdj(point: Vec2, params: &VariationParams) -> Vec2 {
+    let (pdj_a, pdj_b, pdj_c, pdj_d) = (params.0[0], params.0[1], params.0[2], params.0[3]);
+    vec2(
+        sinf(pdj_a * point.y) - cosf(pdj_b * point.x),
+        sinf(pdj_c * point.x) - cosf(pdj_d * point.y),
+    )
+}