based on https://github.com/raphlinus/font-rs/blob/master/src/raster.rs simplified and extended

raster.rs

use std::process::Command;
use std::cmp::min;
use std::io::prelude::*;
use std::fs::File;

#[derive(Copy, Clone)]
pub struct Point {
    pub x: f32,
    pub y: f32,
}

pub struct Raster {
    w: usize,
    h: usize,
    a: Vec<f32>,
}

impl Raster {
    pub fn new(w: usize, h: usize) -> Raster {
        Raster { w: w, h: h, a: vec![0.0; w * h + 4], }
    }

    pub fn draw_line(&mut self, p0: &Point, p1: &Point) {
        //println!("draw_line {} {}", p0, p1);
        if p0.y == p1.y {
            return;
        }
        let (dir, p0, p1) = if p0.y < p1.y {
            (1.0, p0, p1)
        } else {
            (-1.0, p1, p0)
        };
        let dxdy = (p1.x - p0.x) / (p1.y - p0.y);
        let mut x = p0.x;
        let y0 = p0.y as usize; // note: implicit max of 0 because usize (TODO: really true?)
        if p0.y < 0.0 {
            x -= p0.y * dxdy;
        }
        for y in y0..min(self.h, p1.y.ceil() as usize) {
            let linestart = y * self.w;
            let dy = ((y + 1) as f32).min(p1.y) - (y as f32).max(p0.y);
            let xnext = x + dxdy * dy;
            let d = dy * dir;
            let (x0, x1) = if x < xnext { (x, xnext) } else { (xnext, x) };
            let x0floor = x0.floor();
            let x0i = x0floor as i32;
            let x1ceil = x1.ceil();
            let x1i = x1ceil as i32;
            if x1i <= x0i + 1 {
                let xmf = 0.5 * (x + xnext) - x0floor;
                self.a[linestart + x0i as usize] += d - d * xmf;
                self.a[linestart + (x0i + 1) as usize] += d * xmf;
            } else {
                let s = 1.0 / (x1 - x0);
                let x0f = x0 - x0floor;
                let a0 = 0.5 * s * (1.0 - x0f) * (1.0 - x0f);
                let x1f = x1 - x1ceil + 1.0;
                let am = 0.5 * s * x1f * x1f;
                self.a[linestart + x0i as usize] += d * a0;
                if x1i == x0i + 2 {
                    self.a[linestart + (x0i + 1) as usize] += d * (1.0 - a0 - am);
                } else {
                    let a1 = s * (1.5 - x0f);
                    self.a[linestart + (x0i + 1) as usize] += d * (a1 - a0);
                    for xi in x0i + 2..x1i - 1 {
                        self.a[linestart + xi as usize] += d * s;
                    }
                    let a2 = a1 + (x1i - x0i - 3) as f32 * s;
                    self.a[linestart + (x1i - 1) as usize] += d * (1.0 - a2 - am);
                }
                self.a[linestart + x1i as usize] += d * am;
            }
            x = xnext;
        }
    }

    pub fn draw_quad(&mut self, p0: &Point, p1: &Point, p2: &Point) {
        //println!("draw_quad {} {} {}", p0, p1, p2);
        let devx = p0.x - 2.0 * p1.x + p2.x;
        let devy = p0.y - 2.0 * p1.y + p2.y;
        let devsq = devx * devx + devy * devy;
        if devsq < 0.333 {
            self.draw_line(p0, p2);
            return;
        }
        let tol = 3.0;
        let n = 1 + (tol * (devx * devx + devy * devy)).sqrt().sqrt().floor() as usize;
        //println!("n = {}", n);
        let mut p = *p0;
        let nrecip = 1.0 / (n as f32);
        let mut t = 0.0;
        for _i in 0..n - 1 {
            t += nrecip;

            // https://developer.roblox.com/en-us/articles/Bezier-curves
            let omt = 1.0 - t;
            let omt2 = omt * omt;
            let t2 = t * t;
            let pn = Point {
                x: omt2 * p0.x + 2.0 * omt * t * p1.x + t2 * p2.x,
                y: omt2 * p0.y + 2.0 * omt * t * p1.y + t2 * p2.y,
            };
            self.draw_line(&p, &pn);
            p = pn;
        }
        self.draw_line(&p, p2);
    }
    
    // Copied from draw_quad with little changes
    pub fn draw_cubic(&mut self, p0: &Point, p1: &Point, p2: &Point, p3: &Point) {
        //println!("draw_quad {} {} {}", p0, p1, p2);
        let devx = p0.x - 2.0 * (p1.x + p2.x) / 2.0 + p3.x; // just a made one, try to find a better way
        let devy = p0.y - 2.0 * (p1.x + p2.x) / 2.0 + p3.y;
        let devsq = devx * devx + devy * devy;
        if devsq < 0.333 {
            self.draw_line(p0, p3);
            return;
        }
        let tol = 3.0;
        let n = 1 + (tol * (devx * devx + devy * devy)).sqrt().sqrt().floor() as usize;
        //println!("n = {}", n);
        
        let mut p = *p0;
        let nrecip = 1.0 / (n as f32);
        let mut t = 0.0;
        for _i in 0..n - 1 {
            t += nrecip;
            
            // https://twitter.com/freyaholmer/status/1063633411246628864
            let omt = 1.0 - t;
            let omt2 = omt * omt;
            let t2 = t * t;
            let pn = Point {
                x: p0.x * (omt2 * omt) + p1.x * (3.0 * omt2 * t) + p2.x * (3.0 * omt * t2) + p3.x * (t2 * t),
                y: p0.y * (omt2 * omt) + p1.y * (3.0 * omt2 * t) + p2.y * (3.0 * omt * t2) + p3.y * (t2 * t),
            };
            self.draw_line(&p, &pn);
            p = pn;
        }
        self.draw_line(&p, p3);
    }

    pub fn get_bitmap(&self) -> Vec<u8> {
        let mut acc = 0.0;
        self.a[0..self.w * self.h].iter()
            .map(|c| {
                // This would translate really well to SIMD
                acc += c;
                let y = acc.abs();
                let y = if y < 1.0 { y } else { 1.0 };
                (255.0 * y) as u8
            })
            .collect()
    }
}

fn main() -> std::io::Result<()> {
    let mut buffer = File::create("out.pbm")?;
    
    let mut raster = Raster::new(100, 100);
    let x = Point { x: 5 as f32, y: 5 as f32, };
    let y = Point { x: 100 as f32, y: 5 as f32, };
    let z = Point { x: 30 as f32, y: 60 as f32, };
    let c1 = Point { x: 60 as f32, y: 10 as f32, };
    let c2 = Point { x: 55 as f32, y: 20 as f32, };
    raster.draw_line(&x, &y);
    // a.draw_quad(&y, &c2, &z);
    raster.draw_cubic(&y, &c1, &c2, &z);
    raster.draw_line(&z, &x);
    
    buffer.write_all(format!("P5 {} {} {}\n", 100, 100, 255).as_bytes())?;
    buffer.write_all(&raster.get_bitmap())?;
    Command::new("xdg-open").arg("out.pbm").output();
    
    Ok(())
}

main()