hittekaart 0.2.0

//! Actual rendering functions for heatmaps.
//!
//! We begin the rendering by having `Renderer::prepare()` to turn a list of GPX tracks into a
//! [`HeatCounter`], which is basically a grayscale heatmap, where each pixel represents the number
//! of tracks that goes through this pixel.
//!
//! We then render the colored heatmap tiles in `Renderer::colorize()`, which provides us with
//! colorful PNG data.
use colorgrad::Gradient;
use crossbeam_channel::Sender;
use image::{ImageBuffer, Luma, Pixel, Rgba, RgbaImage};
use nalgebra::{vector, Vector2};
use rayon::iter::ParallelIterator;
use std::iter;

use super::{
    super::{
        error::{Error, Result},
        gpx::Coordinates,
        layer::{self, TileLayer},
    },
    RenderedTile,
};

/// Type for the intermediate heat counters.
pub type HeatCounter = TileLayer<Luma<u8>>;

fn render_circle<P: Pixel>(layer: &mut TileLayer<P>, center: (u64, u64), radius: u64, pixel: P) {
    let topleft = (center.0 - radius, center.1 - radius);
    let rad_32: u32 = radius.try_into().unwrap();
    let mut circle = ImageBuffer::<P, Vec<P::Subpixel>>::new(rad_32 * 2 + 1, rad_32 * 2 + 1);
    imageproc::drawing::draw_filled_circle_mut(
        &mut circle,
        (
            i32::try_from(radius).unwrap(),
            i32::try_from(radius).unwrap(),
        ),
        radius.try_into().unwrap(),
        pixel,
    );
    layer.blit_nonzero(topleft.0, topleft.1, &circle);
}

fn direction_vector(a: (u64, u64), b: (u64, u64)) -> Vector2<f64> {
    let dx = if b.0 > a.0 {
        (b.0 - a.0) as f64
    } else {
        -((a.0 - b.0) as f64)
    };
    let dy = if b.1 > a.1 {
        (b.1 - a.1) as f64
    } else {
        -((a.1 - b.1) as f64)
    };
    vector![dx, dy]
}

fn render_line<P: Pixel>(
    layer: &mut TileLayer<P>,
    start: (u64, u64),
    end: (u64, u64),
    thickness: u64,
    pixel: P,
) {
    use imageproc::point::Point;

    if start == end {
        return;
    }

    fn unsigned_add(a: Vector2<u64>, b: Vector2<i32>) -> Vector2<u64> {
        let x = if b[0] < 0 {
            a[0] - u64::from(b[0].unsigned_abs())
        } else {
            a[0] + u64::try_from(b[0]).unwrap()
        };
        let y = if b[1] < 0 {
            a[1] - u64::from(b[1].unsigned_abs())
        } else {
            a[1] + u64::try_from(b[1]).unwrap()
        };
        vector![x, y]
    }

    let r = direction_vector(start, end);
    let normal = vector![r[1], -r[0]].normalize();

    let start = vector![start.0, start.1];
    let end = vector![end.0, end.1];

    let displacement = normal * thickness as f64;
    let displacement = displacement.map(|x| x as i32);
    if displacement == vector![0, 0] {
        return;
    }
    let polygon = [
        unsigned_add(start, displacement),
        unsigned_add(end, displacement),
        unsigned_add(end, -displacement),
        unsigned_add(start, -displacement),
    ];
    let min_x = polygon.iter().map(|p| p[0]).min().unwrap();
    let min_y = polygon.iter().map(|p| p[1]).min().unwrap();
    let max_x = polygon.iter().map(|p| p[0]).max().unwrap();
    let max_y = polygon.iter().map(|p| p[1]).max().unwrap();

    let mut overlay = ImageBuffer::<P, Vec<P::Subpixel>>::new(
        (max_x - min_x).try_into().unwrap(),
        (max_y - min_y).try_into().unwrap(),
    );
    let adjusted_poly = polygon
        .into_iter()
        .map(|p| Point::new((p[0] - min_x) as i32, (p[1] - min_y) as i32))
        .collect::<Vec<_>>();
    imageproc::drawing::draw_polygon_mut(&mut overlay, &adjusted_poly, pixel);

    layer.blit_nonzero(min_x, min_y, &overlay);
}

fn merge_heat_counter(base: &mut HeatCounter, overlay: &HeatCounter) {
    for (tx, ty, source) in overlay.enumerate_tiles() {
        let target = base.tile_mut(tx, ty);
        for (t, s) in target.iter_mut().zip(source.iter()) {
            *t += *s;
        }
    }
}

fn colorize_tile(tile: &ImageBuffer<Luma<u8>, Vec<u8>>, lut: &[Rgba<u8>]) -> RgbaImage {
    let mut result = ImageBuffer::from_pixel(tile.width(), tile.height(), [0, 0, 0, 0].into());
    for (pixel, target) in tile.pixels().zip(result.pixels_mut()) {
        if pixel[0] > 0 {
            *target = lut[usize::from(pixel[0])];
        }
    }
    result
}

fn prepare_lut(max: u8) -> Vec<Rgba<u8>> {
    let gradient = colorgrad::GradientBuilder::new()
        .css("rgba(0, 0, 0, 1) 0%, rgba(164, 0, 0, 1) 25%, rgba(224, 5, 0, 1) 50%, rgba(255, 187, 15, 1) 75%, rgba(255, 255, 255, 1) 100%")
        .build::<colorgrad::LinearGradient>()
        .unwrap();
    iter::once([0, 0, 0, 0].into())
        .chain((1..=max).map(|count| {
            let alpha = count as f32 / max as f32;
            // If we simply use alpha here, we get a linear mapping, and a jump from 1 to 2
            // repetitions is worth as much as a jump from 9 to 10. By transforming alpha via the
            // cubic function given below, we give more weight to the "early" repetitions by having
            // a steeper slope, and less weight to later repetitions.
            // There is no science behind the power 7, others work as well (if you adjust the
            // sign). It seemed like a decent trade-off though.
            let alpha = (alpha - 1.0).powi(7) + 1.0;
            let color = gradient.at(alpha);
            color.to_rgba8().into()
        }))
        .collect()
}

#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
pub struct Renderer;

impl super::Renderer for Renderer {
    type Prepared = HeatCounter;

    /// Renders the heat counter for the given zoom level and track points.
    ///
    /// The given callback will be called when a track has been rendered and merged into the
    /// accumulator, to allow for UI feedback. The passed parameter is the number of tracks that have
    /// been rendered since the last call.
    fn prepare(
        &self,
        zoom: u32,
        tracks: &[Vec<Coordinates>],
        tick: Sender<()>,
    ) -> Result<HeatCounter> {
        let mut heatcounter = TileLayer::from_pixel([0].into());

        for track in tracks {
            let mut layer = TileLayer::from_pixel([0].into());

            let points = track
                .iter()
                .map(|coords| coords.web_mercator(zoom))
                .collect::<Vec<_>>();

            for point in points.iter() {
                render_circle(&mut layer, *point, (zoom as u64 / 4).max(2) - 1, [1].into());
            }

            for (a, b) in points.iter().zip(points.iter().skip(1)) {
                render_line(&mut layer, *a, *b, (zoom as u64 / 4).max(1), [1].into());
            }

            merge_heat_counter(&mut heatcounter, &layer);
            tick.send(()).unwrap();
        }
        Ok(heatcounter)
    }

    /// Lazily colorizes a [`HeatCounter`] by colorizing it tile-by-tile and saving a tile before
    /// rendering the next one.
    ///
    /// This function calls the given callback with each rendered tile, and the function is responsible
    /// for saving it. If the callback returns an `Err(...)`, the error is passed through.
    ///
    /// Note that this function internally uses `rayon` for parallization. If you want to limit the
    /// number of threads used, set up the global [`rayon::ThreadPool`] first.
    fn colorize(&self, layer: HeatCounter, tx: Sender<RenderedTile>) -> Result<()> {
        let max = layer.pixels().map(|l| l.0[0]).max().unwrap_or_default();
        if max == 0 {
            return Ok(());
        }

        // max is a u8, so at most we have 256 colors to compute. We can handle that.
        let color_lut = prepare_lut(max);

        layer
            .into_parallel_tiles()
            .try_for_each_with(tx, |tx, (tile_x, tile_y, tile)| {
                let colorized = colorize_tile(&tile, &color_lut);
                let data = layer::compress_png_as_bytes(&colorized)?;
                tx.send(RenderedTile {
                    x: tile_x,
                    y: tile_y,
                    data,
                })?;
                Ok::<(), Error>(())
            })
    }

    fn tile_count(&self, layer: &HeatCounter) -> Result<u64> {
        Ok(layer.tile_count().try_into().unwrap())
    }
}