Skip to content

Density-Wave Galaxy

This example renders an animated density-wave spiral galaxy.

Preview

Run And Adapt

Commands below assume a Datoviz source checkout and start at the repository root. Use your configured build environment; Python routes additionally require local bindings.

Route Availability Command or action
C Canonical native source just example-c showcases/galaxy (build and run), or rerun ./build/examples/c/showcases/galaxy
Python Available python3 -m examples.python.gallery.showcases.galaxy
Browser Live WebGPU route Open live example

Prepared data required

This example intentionally fails when its prepared input is absent; it does not substitute synthetic data. Prepare it from the repository root with deterministic runtime simulation; no external data payload.

Use this example as capability or integration evidence, not as a minimal copy-paste template. Start from the nearest supported, copy-safe example and add this feature after verifying the linked API reference.

What To Look For

Small temperature-colored stars form a warm central bulge while large soft dust sprites reveal nested blue spiral arms. Additive blending lets overlapping particles accumulate light without turning the dense core into an opaque disc. Drag to orbit through the thin stellar and dust layers, and use the wheel to inspect the core or the full disk. A sparse world-space star shell provides depth without competing with the simulated particles.

The density-wave equations and rendering composition are adapted from Nicolas P. Rougier's Glumpy galaxy example and Ingo Berg's galaxy simulation. See galaxy_model.c for the BSD notice and model details.

Source

#!/usr/bin/env python3
"""Animated NumPy density-wave spiral galaxy.

The density-wave model is adapted from Nicolas P. Rougier's Glumpy example.

Copyright (c) 2014 Nicolas P. Rougier.

Redistribution and use in source and binary forms, with or without modification, are permitted
provided that the following conditions are met:

* Redistributions of source code must retain the above copyright notice, this list of conditions
  and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright notice, this list of
  conditions and the following disclaimer in the documentation and/or other materials provided
  with the distribution.
* Neither Nicolas P. Rougier's name nor the names of contributors may be used to endorse or promote
  products derived from this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER
IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
"""

from __future__ import annotations

import ctypes
from dataclasses import dataclass

import numpy as np

import datoviz as dvz
from examples.python.gallery import common as ex


STAR_COUNT = 28_000
DUST_COUNT = 18_000
BACKGROUND_COUNT = 1_200
INITIAL_ANGLES = (ctypes.c_float * 3)(-0.302710, +0.044938, -0.017917)
BLACK_BODY = np.array(
    [[1.00, 0.23, 0.01], [1.00, 0.55, 0.25], [1.00, 0.78, 0.62], [1.00, 0.94, 0.99], [0.76, 0.80, 1.00], [0.60, 0.69, 1.00]],
    dtype=np.float32,
)


@dataclass
class Galaxy:
    radius: np.ndarray
    ratio: np.ndarray
    orientation: np.ndarray
    phase: np.ndarray
    height: np.ndarray
    colors: np.ndarray
    sizes: np.ndarray

    def positions(self, elapsed: float) -> np.ndarray:
        theta = self.phase + 0.095 * elapsed
        x = self.radius * np.cos(theta)
        y = self.radius * self.ratio * np.sin(theta)
        c, s = np.cos(self.orientation), np.sin(self.orientation)
        return np.column_stack((c * x - s * y, self.height, s * x + c * y)).astype(np.float32)


def _model() -> Galaxy:
    rng = np.random.default_rng(0x6A09E667)
    star_radius = np.abs(0.47 * rng.normal(size=STAR_COUNT))
    dust_radius = np.sqrt(rng.uniform(0.03, 1.0, DUST_COUNT))
    radius = np.concatenate((star_radius, dust_radius)).astype(np.float32)
    ratio = np.clip(1.0 - 0.10 * np.minimum(radius / 0.30, 1.0), 0.90, 1.0)
    orientation = np.where(np.arange(len(radius)) < STAR_COUNT, np.pi / 2 - 5.2 * radius, 5.2 * radius)
    phase = rng.uniform(0.0, 2.0 * np.pi, len(radius)).astype(np.float32)
    scale_height = np.where(np.arange(len(radius)) < STAR_COUNT, 0.012 + 0.055 * np.exp(-radius / 0.20), 0.007)
    height = (scale_height * rng.normal(size=len(radius))).astype(np.float32)

    temperature = np.concatenate((rng.uniform(3000, 9000, STAR_COUNT), 6000 + 3000 * dust_radius))
    coordinate = np.clip((temperature - 3000) / 6000 * (len(BLACK_BODY) - 1), 0, len(BLACK_BODY) - 1)
    lo = coordinate.astype(np.int32)
    hi = np.minimum(lo + 1, len(BLACK_BODY) - 1)
    rgb = (1.0 - (coordinate - lo)[:, None]) * BLACK_BODY[lo] + (coordinate - lo)[:, None] * BLACK_BODY[hi]
    alpha = np.concatenate((rng.uniform(45, 125, STAR_COUNT), rng.uniform(3, 18, DUST_COUNT)))
    colors = np.column_stack((np.clip(255 * rgb, 0, 255).astype(np.uint8), alpha.astype(np.uint8)))
    sizes = np.concatenate((rng.uniform(1.2, 3.2, STAR_COUNT), rng.uniform(5.0, 13.0, DUST_COUNT))).astype(np.float32)
    return Galaxy(radius, ratio, orientation.astype(np.float32), phase, height, colors, sizes)


def _background():
    rng = np.random.default_rng(42)
    directions = rng.normal(size=(BACKGROUND_COUNT, 3))
    directions /= np.linalg.norm(directions, axis=1, keepdims=True)
    positions = (8.0 * directions).astype(np.float32)
    brightness = rng.uniform(0.35, 1.0, BACKGROUND_COUNT)
    colors = np.column_stack((150 * brightness, 175 * brightness, 220 * brightness, np.full(BACKGROUND_COUNT, 255))).astype(np.uint8)
    sizes = rng.uniform(0.7, 2.0, BACKGROUND_COUNT).astype(np.float32)
    return positions, colors, sizes


def _build_scene():
    model = _model()
    scene, figure, panel = ex.scene_panel()
    dvz.dvz_panel_set_background_color(panel, dvz.DvzColor(0, 0, 8, 255))
    camera = dvz.dvz_camera_desc()
    camera.view.eye[:] = (0.0, -3.0, +4.8)
    camera.view.target[:] = (0.0, 0.0, 0.0)
    camera.view.up[:] = (0.0, 1.0, 0.0)
    camera.projection.fov_y = 0.55
    camera.projection.near_clip = 0.05
    camera.projection.far_clip = 100.0
    dvz.dvz_panel_set_camera_desc(panel, ctypes.byref(camera))

    background = dvz.dvz_point(scene, 0)
    positions, colors, sizes = _background()
    dvz.dvz_visual_set_data_many(background, {"position": positions, "color": colors, "diameter_px": sizes})
    ex.set_filled_point_style(background)
    dvz.dvz_visual_set_depth_test(background, False)
    ex.add_visual(panel, background)

    particles = dvz.dvz_point(scene, 0)
    dvz.dvz_visual_set_data_many(particles, {"position": model.positions(0.0), "color": model.colors, "diameter_px": model.sizes})
    ex.set_filled_point_style(particles)
    dvz.dvz_visual_set_alpha_mode(particles, dvz.DVZ_ALPHA_BLENDED)
    dvz.dvz_visual_set_blend_mode(particles, dvz.DVZ_BLEND_ADDITIVE)
    dvz.dvz_visual_set_depth_test(particles, False)
    ex.add_visual(panel, particles)
    return scene, figure, panel, particles, model


def main() -> None:
    scene, figure, panel, particles, model = _build_scene()

    def configure(view) -> None:
        arcball = dvz.dvz_view_arcball(view, panel, None)
        if not arcball or dvz.dvz_arcball_set(arcball, INITIAL_ANGLES) != 0:
            raise RuntimeError("arcball setup failed")
        dvz.dvz_arcball_zoom(arcball, 7.389051)

    def on_frame(_view, _frame_index: int, elapsed: float) -> None:
        if dvz.dvz_visual_set_data(particles, "position", model.positions(elapsed)) != 0:
            raise RuntimeError("galaxy position update failed")

    ex.run_with_frame_callback(scene, figure, "Density-Wave Galaxy", on_frame, configure)


if __name__ == "__main__":
    main()
/*
 * Copyright (c) 2021 Cyrille Rossant and contributors. All rights reserved.
 * Licensed under the MIT license. See LICENSE file in the project root for details.
 * SPDX-License-Identifier: MIT
 */

/* galaxy - This example renders an animated density-wave spiral galaxy.
 *
 * What to look for: small temperature-colored stars form a warm central bulge while large soft
 * dust sprites reveal nested blue spiral arms. Additive blending lets overlapping particles
 * accumulate light without turning the dense core into an opaque disc. Drag to orbit through the
 * thin stellar and dust layers, and use the wheel to inspect the core or the full disk. A sparse
 * world-space star shell provides depth without competing with the simulated particles.
 *
 * The density-wave equations and rendering composition are adapted from Nicolas P. Rougier's
 * Glumpy galaxy example and Ingo Berg's galaxy simulation. See galaxy_model.c for the retained BSD
 * notice and model details.
 *
 * Scenario: showcases_galaxy
 * Style: showcase, astronomical, interactive 3D
 *
 * Build:  just example-c showcases/galaxy
 * Run:    ./build/examples/c/showcases/galaxy --live
 * Smoke:  ./build/examples/c/showcases/galaxy --png
 */



/*************************************************************************************************/
/*  Includes                                                                                     */
/*************************************************************************************************/

#include <math.h>
#include <stdbool.h>
#include <stdint.h>

#include "_alloc.h"
#include "_assertions.h"
#include "datoviz/scene.h"
#include "example_common.h"
#include "example_tuner.h"
#include "galaxy_model.h"
#include "runner/scenario_runner.h"



/*************************************************************************************************/
/*  Constants                                                                                    */
/*************************************************************************************************/

#define WIDTH  EXAMPLE_WINDOW_WIDTH
#define HEIGHT EXAMPLE_WINDOW_HEIGHT

#define SPRITE_SIZE            64u
#define BACKGROUND_STAR_COUNT  1200u
#define BACKGROUND_STAR_RADIUS 48.0f

#define GALAXY_INITIAL_YEARS                    100000.0
#define GALAXY_DEFAULT_MILLION_YEARS_PER_SECOND 1.0f



/*************************************************************************************************/
/*  Forward declarations                                                                         */
/*************************************************************************************************/

DvzScenarioSpec dvz_showcase_galaxy_scenario(void);



/*************************************************************************************************/
/*  Structs                                                                                      */
/*************************************************************************************************/

typedef struct GalaxyState
{
    ExampleTuner tuner;
    GalaxyModel model;
    DvzVisual* visual;
    DvzVisual* background_stars;
    DvzArcball* arcball;
    float* angles;
    DvzSymbolId* symbols;
    double elapsed_years;
    float million_years_per_second;
    bool animate;
    bool show_background_stars;
} GalaxyState;



/*************************************************************************************************/
/*  Helpers                                                                                      */
/*************************************************************************************************/

/**
 * Clamp a normalized scalar.
 *
 * @param value input value
 * @return value clamped to [0, 1]
 */
static float _saturate(float value) { return fminf(fmaxf(value, 0.0f), 1.0f); }



/**
 * Return one deterministic pseudo-random integer for the background star shell.
 *
 * @param state generator state
 * @return next pseudo-random integer
 */
static uint32_t _background_random_u32(uint32_t* state)
{
    ANN(state);
    *state = 1664525u * (*state) + 1013904223u;
    return *state;
}



/**
 * Return one deterministic pseudo-random scalar in [0, 1).
 *
 * @param state generator state
 * @return next pseudo-random scalar
 */
static float _background_random_f32(uint32_t* state)
{
    return (float)(_background_random_u32(state) & 0x00FFFFFFu) / (float)0x01000000u;
}



/**
 * Create a sparse, deterministic world-space background star shell.
 *
 * This follows the point-shell technique used by the textured-planet showcase, but deliberately
 * omits its prominent sun and uses dimmer colors so the galaxy remains the focal layer.
 *
 * @param scene scene
 * @return background point visual, or NULL on failure
 */
static DvzVisual* _create_background_stars(DvzScene* scene)
{
    ANN(scene);
    vec3* positions = dvz_calloc(BACKGROUND_STAR_COUNT, sizeof(*positions));
    DvzColor* colors = dvz_calloc(BACKGROUND_STAR_COUNT, sizeof(*colors));
    float* sizes = dvz_calloc(BACKGROUND_STAR_COUNT, sizeof(*sizes));
    if (positions == NULL || colors == NULL || sizes == NULL)
        goto fail;

    uint32_t random = 0xA54FF53Au;
    for (uint32_t i = 0; i < BACKGROUND_STAR_COUNT; i++)
    {
        const float z = 2.0f * _background_random_f32(&random) - 1.0f;
        const float phi = 2.0f * (float)DVZ_PI * _background_random_f32(&random);
        const float radius_xy = sqrtf(fmaxf(0.0f, 1.0f - z * z));
        positions[i][0] = BACKGROUND_STAR_RADIUS * radius_xy * cosf(phi);
        positions[i][1] = BACKGROUND_STAR_RADIUS * radius_xy * sinf(phi);
        positions[i][2] = BACKGROUND_STAR_RADIUS * z;

        const float brightness = 0.35f + 0.65f * _background_random_f32(&random);
        const float warmth = _background_random_f32(&random);
        colors[i] = (DvzColor){
            (uint8_t)(brightness * (150.0f + 45.0f * warmth)),
            (uint8_t)(brightness * (165.0f + 28.0f * warmth)),
            (uint8_t)(brightness * (205.0f - 18.0f * warmth)), 255};
        const float size_random = _background_random_f32(&random);
        sizes[i] = 0.8f + 1.8f * size_random * size_random;
    }

    DvzVisual* stars = dvz_point(scene, 0);
    if (stars == NULL)
        goto fail;
    DvzVisualDataUpdate updates[] = {
        {.attr_name = "position", .data = positions, .item_count = BACKGROUND_STAR_COUNT},
        {.attr_name = "color", .data = colors, .item_count = BACKGROUND_STAR_COUNT},
        {.attr_name = "size", .data = sizes, .item_count = BACKGROUND_STAR_COUNT},
    };
    if (dvz_visual_set_data_many(stars, updates, 3) != DVZ_OK ||
        dvz_visual_set_depth_test(stars, false) != DVZ_OK)
    {
        goto fail;
    }

    dvz_free(sizes);
    dvz_free(colors);
    dvz_free(positions);
    return stars;

fail:
    dvz_free(sizes);
    dvz_free(colors);
    dvz_free(positions);
    return NULL;
}



/**
 * Fill a soft white particle texture matching the original 64-pixel Glumpy sprite profile.
 *
 * RGB carries the radial intensity while alpha stays opaque. The particle's retained color alpha
 * therefore controls additive brightness without squaring the sampled intensity.
 *
 * @param rgba output RGBA8 sprite
 */
static void _fill_particle_sprite(DvzColor rgba[SPRITE_SIZE * SPRITE_SIZE])
{
    ANN(rgba);
    for (uint32_t y = 0; y < SPRITE_SIZE; y++)
    {
        for (uint32_t x = 0; x < SPRITE_SIZE; x++)
        {
            const float px = 2.0f * ((float)x + 0.5f) / (float)SPRITE_SIZE - 1.0f;
            const float py = 2.0f * ((float)y + 0.5f) / (float)SPRITE_SIZE - 1.0f;
            const float radius2 = px * px + py * py;
            const float core = 0.70f * expf(-radius2 / 0.080f);
            const float glow = 0.17f * expf(-radius2 / 0.450f);
            const float halo = 0.05f * expf(-radius2 / 0.800f);
            const uint8_t intensity = (uint8_t)(255.0f * _saturate(core + glow + halo) + 0.5f);
            rgba[y * SPRITE_SIZE + x] = (DvzColor){intensity, intensity, intensity, 255};
        }
    }
}



/**
 * Free scenario-owned model and attribute arrays.
 *
 * @param state galaxy scenario state
 */
static void _free_state(GalaxyState* state)
{
    if (state == NULL)
        return;
    example_tuner_detach(&state->tuner);
    dvz_free(state->symbols);
    dvz_free(state->angles);
    galaxy_model_destroy(&state->model);
    dvz_free(state);
}



/**
 * Restore the live galaxy controls to their showcase defaults.
 *
 * @param user galaxy scenario state
 */
static void _galaxy_controls_reset(void* user)
{
    GalaxyState* state = (GalaxyState*)user;
    if (state == NULL)
        return;
    state->animate = true;
    state->show_background_stars = true;
    state->million_years_per_second = GALAXY_DEFAULT_MILLION_YEARS_PER_SECOND;
    if (state->background_stars != NULL)
        (void)dvz_visual_set_visible(state->background_stars, true);
}



/**
 * Draw the live simulation and background controls.
 *
 * @param gui GUI
 * @param user galaxy scenario state
 * @return whether a control changed
 */
static bool _galaxy_controls_gui(DvzGui* gui, void* user)
{
    GalaxyState* state = (GalaxyState*)user;
    if (gui == NULL || state == NULL)
        return false;

    bool changed = false;
    changed |= dvz_gui_checkbox(gui, "Animate galaxy", &state->animate);
    changed |= dvz_gui_slider_float_format(
        gui, "Simulation speed", &state->million_years_per_second, 0.0f, 6.0f,
        "%.2f million years/s");
    changed |= dvz_gui_checkbox(gui, "Background stars", &state->show_background_stars);
    return changed;
}



/**
 * Apply live galaxy controls that affect retained visual state.
 *
 * @param user galaxy scenario state
 */
static void _galaxy_controls_apply(void* user)
{
    GalaxyState* state = (GalaxyState*)user;
    if (state == NULL || state->background_stars == NULL)
        return;
    (void)dvz_visual_set_visible(state->background_stars, state->show_background_stars);
}



/*************************************************************************************************/
/*  Scenario callbacks                                                                           */
/*************************************************************************************************/

/**
 * Initialize the retained density-wave galaxy scenario.
 *
 * @param ctx scenario context
 * @param out_user output scenario state
 * @return whether initialization succeeded
 */
static bool _scenario_init(DvzScenarioContext* ctx, void** out_user)
{
    if (ctx == NULL)
        return false;
    if (out_user != NULL)
        *out_user = NULL;

    bool ok = false;
    GalaxyState* state = dvz_calloc(1, sizeof(*state));
    if (state == NULL)
        return false;
    if (out_user != NULL)
        *out_user = state;
    state->tuner = example_tuner("Galaxy settings");
    _galaxy_controls_reset(state);

    ok = galaxy_model_init(&state->model, 0x6A09E667u);
    EXAMPLE_CHECK(ok, "galaxy model initialization failed");
    state->angles = dvz_calloc(state->model.particle_count, sizeof(*state->angles));
    state->symbols = dvz_calloc(state->model.particle_count, sizeof(*state->symbols));
    EXAMPLE_CHECK(
        state->angles != NULL && state->symbols != NULL,
        "galaxy marker attribute allocation failed");

    ctx->figure = dvz_figure(ctx->scene, ctx->width, ctx->height, 0);
    EXAMPLE_CHECK(ctx->figure != NULL, "dvz_figure() failed");
    example_tuner_figure(&state->tuner, ctx->figure);
    DvzPanel* panel = dvz_panel_full(ctx->figure);
    EXAMPLE_CHECK(panel != NULL, "dvz_panel_full() failed");
    dvz_panel_set_background_color(panel, dvz_color_rgba(0, 0, 8, 255));

    DvzCameraDesc camera = dvz_camera_desc();
    camera.view.eye[0] = 0.0f;
    camera.view.eye[1] = -3.0f;
    camera.view.eye[2] = +4.8f;
    camera.view.target[0] = 0.0f;
    camera.view.target[1] = 0.0f;
    camera.view.target[2] = 0.0f;
    camera.view.up[0] = 0.0f;
    camera.view.up[1] = 1.0f;
    camera.view.up[2] = 0.0f;
    camera.projection.fov_y = 0.55f;
    camera.projection.near_clip = 0.05f;
    camera.projection.far_clip = 100.0f;
    DvzResult rc = dvz_panel_set_camera_desc(panel, &camera);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy camera setup failed");

    state->background_stars = _create_background_stars(ctx->scene);
    EXAMPLE_CHECK(state->background_stars != NULL, "galaxy background star creation failed");
    rc = dvz_panel_add_visual(panel, state->background_stars, NULL);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy background star attachment failed");

    DvzSymbolSet* symbol_set = dvz_symbol_set(ctx->scene, 0);
    EXAMPLE_CHECK(symbol_set != NULL, "dvz_symbol_set() failed");
    DvzColor sprite[SPRITE_SIZE * SPRITE_SIZE] = {{0}};
    _fill_particle_sprite(sprite);
    DvzSymbolId sprite_id = dvz_symbol_bitmap(
        symbol_set, "galaxy_soft_particle", sprite, SPRITE_SIZE, SPRITE_SIZE, NULL);
    EXAMPLE_CHECK(sprite_id != DVZ_SYMBOL_ID_INVALID, "galaxy particle symbol creation failed");
    for (uint32_t i = 0; i < state->model.particle_count; i++)
        state->symbols[i] = sprite_id;

    state->visual = dvz_marker(ctx->scene, 0);
    EXAMPLE_CHECK(state->visual != NULL, "dvz_marker() failed");
    rc = dvz_marker_set_symbols(state->visual, symbol_set);
    EXAMPLE_CHECK(rc == DVZ_OK, "dvz_marker_set_symbols() failed");
    DvzMarkerStyle style = dvz_marker_style();
    style.aspect = DVZ_SHAPE_ASPECT_FILLED;
    style.stroke_width_px = 0.0f;
    rc = dvz_marker_set_style(state->visual, &style);
    EXAMPLE_CHECK(rc == DVZ_OK, "dvz_marker_set_style() failed");

    DvzVisualDataUpdate updates[] = {
        {.attr_name = "position",
         .data = state->model.positions,
         .item_count = state->model.particle_count},
        {.attr_name = "color",
         .data = state->model.colors,
         .item_count = state->model.particle_count},
        {.attr_name = "diameter_px",
         .data = state->model.sizes,
         .item_count = state->model.particle_count},
        {.attr_name = "angle", .data = state->angles, .item_count = state->model.particle_count},
        {.attr_name = "symbol", .data = state->symbols, .item_count = state->model.particle_count},
    };
    rc = dvz_visual_set_data_many(state->visual, updates, 5);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy marker upload failed");
    rc = dvz_visual_set_alpha_mode(state->visual, DVZ_ALPHA_BLENDED);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy alpha mode setup failed");
    rc = dvz_visual_set_blend_mode(state->visual, DVZ_BLEND_ADDITIVE);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy additive blend setup failed");
    rc = dvz_visual_set_depth_test(state->visual, false);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy depth setup failed");
    rc = dvz_panel_add_visual(panel, state->visual, NULL);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy visual attachment failed");

    DvzController* controller = dvz_arcball(ctx->scene, NULL);
    EXAMPLE_CHECK(controller != NULL, "dvz_arcball() failed");
    state->arcball = dvz_controller_arcball(controller);
    EXAMPLE_CHECK(state->arcball != NULL, "dvz_controller_arcball() failed");
    rc = dvz_scenario_bind_controller(ctx, panel, controller, DVZ_DIM_MASK_XYZ);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy arcball binding failed");
    vec3 initial_angles = {-0.302710f, +0.044938f, -0.017917f};
    vec2 initial_pan = {0.0f, 0.0f};
    rc = dvz_arcball_initial(state->arcball, initial_angles);
    EXAMPLE_CHECK(rc == DVZ_OK, "galaxy arcball initialization failed");
    example_tuner_arcball(
        &state->tuner, "Arcball", state->arcball, initial_angles, 7.389051f, initial_pan);
    example_tuner_camera(&state->tuner, "Camera", panel, &camera);
    (void)example_tuner_add_component(
        &state->tuner, "Simulation", state, NULL, _galaxy_controls_gui, _galaxy_controls_apply,
        _galaxy_controls_reset, NULL);

    state->elapsed_years = GALAXY_INITIAL_YEARS;
    ok = true;
cleanup:
    if (!ok && out_user == NULL)
        _free_state(state);
    return ok;
}



/**
 * Advance the density-wave orbital phases and update dynamic marker attributes.
 *
 * @param ctx scenario context
 * @param user_data galaxy scenario state
 */
static void _scenario_frame(DvzScenarioContext* ctx, void* user_data)
{
    GalaxyState* state = (GalaxyState*)user_data;
    if (ctx == NULL || state == NULL || state->visual == NULL)
        return;

    if (ctx->preview_mode)
    {
        state->elapsed_years = GALAXY_INITIAL_YEARS + 1000000.0 *
                                                          GALAXY_DEFAULT_MILLION_YEARS_PER_SECOND *
                                                          dvz_scenario_preview_time(ctx);
    }
    else if (state->animate)
    {
        const double dt = fmin(fmax(ctx->dt, 0.0), 0.1);
        state->elapsed_years += 1000000.0 * (double)state->million_years_per_second * dt;
    }
    galaxy_model_update(&state->model, state->elapsed_years);
    DvzVisualDataUpdate updates[] = {
        {.attr_name = "position",
         .data = state->model.positions,
         .item_count = state->model.particle_count},
        {.attr_name = "diameter_px",
         .data = state->model.sizes,
         .item_count = state->model.particle_count},
    };
    (void)dvz_visual_set_data_many(state->visual, updates, 2);
}



/**
 * Destroy scenario-owned CPU state.
 *
 * @param ctx scenario context
 * @param user_data galaxy scenario state
 */
static void _scenario_destroy(DvzScenarioContext* ctx, void* user_data)
{
    (void)ctx;
    _free_state((GalaxyState*)user_data);
}



/**
 * Attach the live-only galaxy tuner to the native view.
 *
 * @param ctx scenario context
 * @param app app
 * @param view native view
 * @param user_data galaxy scenario state
 * @return whether the tuner was attached or intentionally skipped
 */
static bool
_scenario_native_view(DvzScenarioContext* ctx, DvzApp* app, DvzView* view, void* user_data)
{
    (void)app;
    GalaxyState* state = (GalaxyState*)user_data;
    if (ctx == NULL || ctx->presentation != DVZ_RUNNER_PRESENT_GLFW || state == NULL ||
        view == NULL)
    {
        return true;
    }
    return example_tuner_attach(&state->tuner, view);
}



/*************************************************************************************************/
/*  Functions                                                                                    */
/*************************************************************************************************/

/**
 * Return the density-wave galaxy scenario specification.
 *
 * @return scenario specification
 */
DvzScenarioSpec dvz_showcase_galaxy_scenario(void)
{
    return (DvzScenarioSpec){
        .id = "showcases_galaxy",
        .title = "Density-Wave Galaxy",
        .width = WIDTH,
        .height = HEIGHT,
        .fps = 60.0,
        .requirements = DVZ_SCENARIO_REQ_POINT_VISUAL | DVZ_SCENARIO_REQ_MARKER_VISUAL |
                        DVZ_SCENARIO_REQ_FRAME_CALLBACKS | DVZ_SCENARIO_REQ_CONTROLLER |
                        DVZ_SCENARIO_REQ_ARCBALL | DVZ_SCENARIO_REQ_CONTINUOUS_FRAMES,
        .init = _scenario_init,
        .frame = _scenario_frame,
        .destroy = _scenario_destroy,
    };
}



#ifndef DVZ_EXAMPLE_NO_MAIN
int main(int argc, char** argv)
{
    DvzScenarioSpec spec = dvz_showcase_galaxy_scenario();
    if (example_cli_wants_live_gui(argc, argv))
        spec.native_view = _scenario_native_view;
    return dvz_scenario_run_native_cli(&spec, argc, argv) == DVZ_OK ? 0 : 1;
}
#endif
/*
 * Copyright (c) 2021 Cyrille Rossant and contributors. All rights reserved.
 * Licensed under the MIT license. See LICENSE file in the project root for details.
 * SPDX-License-Identifier: MIT
 *
 * The density-wave model is adapted from Glumpy's galaxy_simulation.py:
 * Copyright (c) 2014 Nicolas P. Rougier.
 *
 * Redistribution and use in source and binary forms, with or without modification, are permitted
 * provided that the following conditions are met:
 *
 * - Redistributions of source code must retain the above copyright notice, this list of conditions
 *   and the following disclaimer.
 * - Redistributions in binary form must reproduce the above copyright notice, this list of
 *   conditions and the following disclaimer in the documentation and/or other materials provided
 *   with the distribution.
 * - Neither the name of Nicolas P. Rougier nor the names of its contributors may be used to
 * endorse or promote products derived from this software without specific prior written
 * permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
 * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
 * FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
 * OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 * POSSIBILITY OF SUCH DAMAGE.
 *
 * Stellar colors reproduce the public-domain specrend black-body reference values used by Glumpy.
 */



/*************************************************************************************************/
/*  Includes                                                                                     */
/*************************************************************************************************/

#include "galaxy_model.h"

#include <math.h>
#include <stdint.h>

#include "_alloc.h"
#include "_assertions.h"
#include "_compat.h"



/*************************************************************************************************/
/*  Constants                                                                                    */
/*************************************************************************************************/

#define GALAXY_RADIUS         13000.0f
#define GALAXY_CORE_RADIUS    4000.0f
#define GALAXY_NORMALIZATION  27000.0f
#define GALAXY_ANGULAR_OFFSET 0.0004f
#define GALAXY_ECCENTRICITY   0.90f

#define GALAXY_STAR_DISK_HEIGHT  180.0f
#define GALAXY_STAR_BULGE_HEIGHT 950.0f
#define GALAXY_BULGE_RADIUS      2800.0f
#define GALAXY_DUST_HEIGHT       90.0f
#define GALAXY_HII_HEIGHT        70.0f

static const float PI = 3.14159265358979323846f;

static const float BLACK_BODY_RGB[19][3] = {
    {1.000f, 0.007f, 0.000f}, {1.000f, 0.126f, 0.000f}, {1.000f, 0.234f, 0.010f},
    {1.000f, 0.349f, 0.067f}, {1.000f, 0.454f, 0.151f}, {1.000f, 0.549f, 0.254f},
    {1.000f, 0.635f, 0.370f}, {1.000f, 0.710f, 0.493f}, {1.000f, 0.778f, 0.620f},
    {1.000f, 0.837f, 0.746f}, {1.000f, 0.890f, 0.869f}, {1.000f, 0.937f, 0.988f},
    {0.907f, 0.888f, 1.000f}, {0.827f, 0.839f, 1.000f}, {0.762f, 0.800f, 1.000f},
    {0.711f, 0.766f, 1.000f}, {0.668f, 0.738f, 1.000f}, {0.632f, 0.714f, 1.000f},
    {0.602f, 0.693f, 1.000f},
};



/*************************************************************************************************/
/*  Structs                                                                                      */
/*************************************************************************************************/

typedef struct GalaxyRandom
{
    uint32_t state;
    bool has_spare;
    float spare;
} GalaxyRandom;



/*************************************************************************************************/
/*  Helpers                                                                                      */
/*************************************************************************************************/

/**
 * Return a deterministic pseudorandom unsigned integer.
 *
 * @param random random generator state
 * @return pseudorandom integer
 */
static uint32_t _random_u32(GalaxyRandom* random)
{
    ANN(random);
    uint32_t x = random->state;
    x ^= x << 13;
    x ^= x >> 17;
    x ^= x << 5;
    random->state = x;
    return x;
}



/**
 * Return a deterministic uniform sample in [0, 1).
 *
 * @param random random generator state
 * @return uniform sample
 */
static float _random_uniform(GalaxyRandom* random)
{
    return (float)(_random_u32(random) >> 8) * (1.0f / 16777216.0f);
}



/**
 * Return a deterministic standard-normal sample.
 *
 * @param random random generator state
 * @return standard-normal sample
 */
static float _random_normal(GalaxyRandom* random)
{
    ANN(random);
    if (random->has_spare)
    {
        random->has_spare = false;
        return random->spare;
    }
    const float u = fmaxf(_random_uniform(random), 1e-7f);
    const float v = _random_uniform(random);
    const float radius = sqrtf(-2.0f * logf(u));
    const float angle = 2.0f * PI * v;
    random->spare = radius * sinf(angle);
    random->has_spare = true;
    return radius * cosf(angle);
}



/**
 * Return the density-wave ellipse eccentricity at one radius.
 *
 * @param radius galactic radius
 * @return ellipse minor-to-major ratio
 */
static float _eccentricity(float radius)
{
    if (radius < GALAXY_CORE_RADIUS)
    {
        return 1.0f + (radius / GALAXY_CORE_RADIUS) * (GALAXY_ECCENTRICITY - 1.0f);
    }
    if (radius <= GALAXY_RADIUS)
    {
        return GALAXY_ECCENTRICITY;
    }
    if (radius < 2.0f * GALAXY_RADIUS)
    {
        const float t = (radius - GALAXY_RADIUS) / GALAXY_RADIUS;
        return GALAXY_ECCENTRICITY + t * (1.0f - GALAXY_ECCENTRICITY);
    }
    return 1.0f;
}



/**
 * Interpolate the black-body color table used by the original Glumpy example.
 *
 * @param temperature temperature in kelvin
 * @param rgb output normalized RGB
 */
static void _black_body_rgb(float temperature, float rgb[3])
{
    ANN(rgb);
    const float x = fminf(fmaxf((temperature - 1000.0f) / 500.0f, 0.0f), 18.0f);
    const uint32_t i0 = (uint32_t)floorf(x);
    const uint32_t i1 = i0 < 18u ? i0 + 1u : i0;
    const float t = x - (float)i0;
    for (uint32_t channel = 0; channel < 3; channel++)
        rgb[channel] = (1.0f - t) * BLACK_BODY_RGB[i0][channel] + t * BLACK_BODY_RGB[i1][channel];
}



/**
 * Convert normalized color and particle brightness to the retained RGBA style.
 *
 * @param particle source particle
 * @return retained color
 */
static DvzColor _particle_color(const GalaxyParticle* particle)
{
    ANN(particle);
    float rgb[3] = {0};
    _black_body_rgb(particle->temperature, rgb);
    float brightness = 1.8f * particle->brightness;
    if (particle->type == GALAXY_PARTICLE_HII_GLOW)
    {
        rgb[0] = fminf(2.0f * rgb[0], 1.0f);
        brightness *= 2.0f;
    }
    else if (particle->type == GALAXY_PARTICLE_HII_CORE)
    {
        rgb[0] = rgb[1] = rgb[2] = 0.90f;
        brightness *= 1.6f;
    }
    return dvz_color_rgba(
        (uint8_t)(255.0f * rgb[0] + 0.5f), (uint8_t)(255.0f * rgb[1] + 0.5f),
        (uint8_t)(255.0f * rgb[2] + 0.5f), (uint8_t)(255.0f * fminf(brightness, 1.0f) + 0.5f));
}



/**
 * Initialize one particle's shared density-wave parameters.
 *
 * @param particle particle to initialize
 * @param radius ellipse major radius
 * @param orientation ellipse orientation in radians
 * @param phase_degrees orbital phase in degrees
 */
static void
_particle_orbit(GalaxyParticle* particle, float radius, float orientation, float phase_degrees)
{
    ANN(particle);
    particle->major_radius = radius;
    particle->minor_radius = radius * _eccentricity(radius);
    particle->orientation = orientation;
    particle->phase_degrees = phase_degrees;
    particle->angular_velocity = 0.000005f;
}



/*************************************************************************************************/
/*  Functions                                                                                    */
/*************************************************************************************************/

/**
 * Initialize a deterministic density-wave galaxy.
 *
 * @param model output model
 * @param seed deterministic nonzero seed
 * @return whether initialization succeeded
 */
bool galaxy_model_init(GalaxyModel* model, uint32_t seed)
{
    ANN(model);
    dvz_memset(model, sizeof(*model), 0, sizeof(*model));
    model->particle_count = GALAXY_PARTICLE_COUNT;
    model->particles = dvz_calloc(model->particle_count, sizeof(*model->particles));
    model->positions = dvz_calloc(model->particle_count, sizeof(*model->positions));
    model->colors = dvz_calloc(model->particle_count, sizeof(*model->colors));
    model->sizes = dvz_calloc(model->particle_count, sizeof(*model->sizes));
    if (model->particles == NULL || model->positions == NULL || model->colors == NULL ||
        model->sizes == NULL)
    {
        galaxy_model_destroy(model);
        return false;
    }

    GalaxyRandom random = {.state = seed != 0 ? seed : 0x91E10DA5u};
    for (uint32_t i = 0; i < GALAXY_STAR_COUNT; i++)
    {
        GalaxyParticle* particle = &model->particles[i];
        const float radius = 0.5f * GALAXY_RADIUS * _random_normal(&random);
        _particle_orbit(
            particle, radius, 90.0f - radius * GALAXY_ANGULAR_OFFSET,
            360.0f * _random_uniform(&random));
        particle->type = GALAXY_PARTICLE_STAR;
        particle->temperature = 3000.0f + 6000.0f * _random_uniform(&random);
        particle->brightness = 0.05f + 0.20f * _random_uniform(&random);
        particle->base_size_px = 3.0f;
        const float scale_height =
            GALAXY_STAR_DISK_HEIGHT +
            GALAXY_STAR_BULGE_HEIGHT * expf(-fabsf(radius) / GALAXY_BULGE_RADIUS);
        particle->height = scale_height * _random_normal(&random);
    }

    const uint32_t dust_first = GALAXY_STAR_COUNT;
    for (uint32_t i = 0; i < GALAXY_DUST_COUNT; i++)
    {
        GalaxyParticle* particle = &model->particles[dust_first + i];
        const float x = 2.0f * GALAXY_RADIUS * _random_uniform(&random);
        const float y = GALAXY_RADIUS * (2.0f * _random_uniform(&random) - 1.0f);
        const float radius = sqrtf(x * x + y * y);
        _particle_orbit(
            particle, radius, radius * GALAXY_ANGULAR_OFFSET, 360.0f * _random_uniform(&random));
        particle->type = GALAXY_PARTICLE_DUST;
        particle->temperature = 6000.0f + 0.25f * radius;
        particle->brightness = 0.01f + 0.01f * _random_uniform(&random);
        particle->base_size_px = 64.0f;
        particle->height = GALAXY_DUST_HEIGHT * _random_normal(&random);
    }

    const uint32_t hii_first = dust_first + GALAXY_DUST_COUNT;
    for (uint32_t i = 0; i < GALAXY_HII_COUNT; i++)
    {
        const float x = GALAXY_RADIUS * (2.0f * _random_uniform(&random) - 1.0f);
        const float y = GALAXY_RADIUS * (2.0f * _random_uniform(&random) - 1.0f);
        const float radius = sqrtf(x * x + y * y);
        const float phase = 360.0f * _random_uniform(&random);
        const float temperature = 3000.0f + 6000.0f * _random_uniform(&random);
        const float brightness = 0.005f + 0.005f * _random_uniform(&random);
        GalaxyParticle* glow = &model->particles[hii_first + i];
        GalaxyParticle* core = &model->particles[hii_first + GALAXY_HII_COUNT + i];
        _particle_orbit(glow, radius, radius * GALAXY_ANGULAR_OFFSET, phase);
        _particle_orbit(core, radius + 1000.0f, radius * GALAXY_ANGULAR_OFFSET, phase);
        glow->type = GALAXY_PARTICLE_HII_GLOW;
        core->type = GALAXY_PARTICLE_HII_CORE;
        glow->temperature = core->temperature = temperature;
        glow->brightness = core->brightness = brightness;
        glow->height = core->height = GALAXY_HII_HEIGHT * _random_normal(&random);
    }

    for (uint32_t i = 0; i < model->particle_count; i++)
        model->colors[i] = _particle_color(&model->particles[i]);
    galaxy_model_update(model, 100000.0);
    return true;
}



/**
 * Update all particle positions at an absolute elapsed simulation time.
 *
 * @param model galaxy model
 * @param elapsed_years elapsed simulation time in years
 */
void galaxy_model_update(GalaxyModel* model, double elapsed_years)
{
    ANN(model);
    ANN(model->particles);
    ANN(model->positions);
    ANN(model->sizes);

    for (uint32_t i = 0; i < model->particle_count; i++)
    {
        const GalaxyParticle* particle = &model->particles[i];
        const float theta =
            (particle->phase_degrees + particle->angular_velocity * (float)elapsed_years) * PI /
            180.0f;
        const float beta = -particle->orientation;
        const float cos_theta = cosf(theta);
        const float sin_theta = sinf(theta);
        const float cos_beta = cosf(beta);
        const float sin_beta = sinf(beta);
        const float x = particle->major_radius * cos_theta * cos_beta -
                        particle->minor_radius * sin_theta * sin_beta;
        const float y = particle->major_radius * cos_theta * sin_beta +
                        particle->minor_radius * sin_theta * cos_beta;
        model->positions[i][0] = x / GALAXY_NORMALIZATION;
        model->positions[i][1] = y / GALAXY_NORMALIZATION;
        model->positions[i][2] = particle->height / GALAXY_NORMALIZATION;
        model->sizes[i] = particle->base_size_px;
    }

    const uint32_t hii_first = GALAXY_STAR_COUNT + GALAXY_DUST_COUNT;
    for (uint32_t i = 0; i < GALAXY_HII_COUNT; i++)
    {
        const uint32_t glow_index = hii_first + i;
        const uint32_t core_index = hii_first + GALAXY_HII_COUNT + i;
        const float dx = (model->positions[glow_index][0] - model->positions[core_index][0]) *
                         GALAXY_NORMALIZATION;
        const float dy = (model->positions[glow_index][1] - model->positions[core_index][1]) *
                         GALAXY_NORMALIZATION;
        const float distance = sqrtf(dx * dx + dy * dy);
        const float scale = fmaxf(1.0f, (1000.0f - distance) / 10.0f - 50.0f);
        const float glow_size = 2.0f * scale;
        const float core_size = scale / 6.0f;
        model->sizes[glow_index] = glow_size > 2.0f ? fminf(glow_size, 64.0f) : 0.0f;
        model->sizes[core_index] = core_size > 2.0f ? fminf(core_size, 64.0f) : 0.0f;
    }
}



/**
 * Destroy all model-owned particle arrays.
 *
 * @param model galaxy model
 */
void galaxy_model_destroy(GalaxyModel* model)
{
    if (model == NULL)
        return;
    dvz_free(model->sizes);
    dvz_free(model->colors);
    dvz_free(model->positions);
    dvz_free(model->particles);
    dvz_memset(model, sizeof(*model), 0, sizeof(*model));
}
/*
 * Copyright (c) 2021 Cyrille Rossant and contributors. All rights reserved.
 * Licensed under the MIT license. See LICENSE file in the project root for details.
 * SPDX-License-Identifier: MIT
 */

/* Deterministic density-wave galaxy model used by the galaxy showcase. */

#pragma once



/*************************************************************************************************/
/*  Includes                                                                                     */
/*************************************************************************************************/

#include <stdbool.h>
#include <stdint.h>

#include "datoviz/scene.h"



/*************************************************************************************************/
/*  Constants                                                                                    */
/*************************************************************************************************/

#define GALAXY_STAR_COUNT     35000u
#define GALAXY_DUST_COUNT     26250u
#define GALAXY_HII_COUNT      200u
#define GALAXY_PARTICLE_COUNT (GALAXY_STAR_COUNT + GALAXY_DUST_COUNT + 2u * GALAXY_HII_COUNT)



/*************************************************************************************************/
/*  Enums                                                                                        */
/*************************************************************************************************/

typedef enum GalaxyParticleType
{
    GALAXY_PARTICLE_STAR = 0,
    GALAXY_PARTICLE_DUST,
    GALAXY_PARTICLE_HII_GLOW,
    GALAXY_PARTICLE_HII_CORE,
} GalaxyParticleType;



/*************************************************************************************************/
/*  Structs                                                                                      */
/*************************************************************************************************/

typedef struct GalaxyParticle
{
    float phase_degrees;
    float angular_velocity;
    float orientation;
    float major_radius;
    float minor_radius;
    float temperature;
    float brightness;
    float base_size_px;
    float height;
    GalaxyParticleType type;
} GalaxyParticle;


typedef struct GalaxyModel
{
    uint32_t particle_count;
    GalaxyParticle* particles;
    vec3* positions;
    DvzColor* colors;
    float* sizes;
} GalaxyModel;



/*************************************************************************************************/
/*  Functions                                                                                    */
/*************************************************************************************************/

bool galaxy_model_init(GalaxyModel* model, uint32_t seed);

void galaxy_model_update(GalaxyModel* model, double elapsed_years);

void galaxy_model_destroy(GalaxyModel* model);
Example details

Tags

scientific, simulation, astronomy, particles, point, marker, bitmap-symbol, additive-blending, 3d, interactive, arcball, gui, starfield, animation, capture

Data

Field Value
kind simulated

Dataset

Field Value
name Density-wave spiral galaxy simulation
source https://github.com/glumpy/glumpy/tree/master/examples
source_url https://glumpy.github.io/gallery.html
license BSD-3-Clause for the adapted Glumpy simulation; specrend color reference is public domain
citation Nicolas P. Rougier, Glumpy galaxy example; Ingo Berg, Simulating a Galaxy with the density wave theory
preprocessing deterministic runtime simulation; no external data payload
provenance The orbital model is adapted from Glumpy's galaxy_simulation.py and preserves its BSD notice in galaxy_model.c. Stellar colors use the public-domain specrend reference values.