From 7e61ccca09f5dcc072d0f7af389de71c48d7ee21 Mon Sep 17 00:00:00 2001 From: Kenzo Lespagnol Date: Thu, 30 Nov 2023 10:41:07 +0100 Subject: [PATCH] Initial commit for poisson sampling --- .gitignore | 1 + .vscode/launch.json | 64 ++ Cargo.lock | 75 +++ Cargo.toml | 10 + notebook/test.ipynb | 1547 +++++++++++++++++++++++++++++++++++++++++++ rustfmt.toml | 4 + src/lib.rs | 243 +++++++ src/main.rs | 9 + 8 files changed, 1953 insertions(+) create mode 100644 .gitignore create mode 100644 .vscode/launch.json create mode 100644 Cargo.lock create mode 100644 Cargo.toml create mode 100644 notebook/test.ipynb create mode 100644 rustfmt.toml create mode 100644 src/lib.rs create mode 100644 src/main.rs diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..ea8c4bf --- /dev/null +++ b/.gitignore @@ -0,0 +1 @@ +/target diff --git a/.vscode/launch.json b/.vscode/launch.json new file mode 100644 index 0000000..3f90ed4 --- /dev/null +++ b/.vscode/launch.json @@ -0,0 +1,64 @@ +{ + // Use IntelliSense to learn about possible attributes. + // Hover to view descriptions of existing attributes. + // For more information, visit: https://go.microsoft.com/fwlink/?linkid=830387 + "version": "0.2.0", + "configurations": [ + { + "type": "lldb", + "request": "launch", + "name": "Debug unit tests in library 'poisson_sampling'", + "cargo": { + "args": [ + "test", + "--no-run", + "--lib", + "--package=poisson_sampling" + ], + "filter": { + "name": "poisson_sampling", + "kind": "lib" + } + }, + "args": [], + "cwd": "${workspaceFolder}" + }, + { + "type": "lldb", + "request": "launch", + "name": "Debug executable 'poisson_sampling'", + "cargo": { + "args": [ + "build", + "--bin=poisson_sampling", + "--package=poisson_sampling" + ], + "filter": { + "name": "poisson_sampling", + "kind": "bin" + } + }, + "args": [], + "cwd": "${workspaceFolder}" + }, + { + "type": "lldb", + "request": "launch", + "name": "Debug unit tests in executable 'poisson_sampling'", + "cargo": { + "args": [ + "test", + "--no-run", + "--bin=poisson_sampling", + "--package=poisson_sampling" + ], + "filter": { + "name": "poisson_sampling", + "kind": "bin" + } + }, + "args": [], + "cwd": "${workspaceFolder}" + } + ] +} \ No newline at end of file diff --git a/Cargo.lock b/Cargo.lock new file mode 100644 index 0000000..c8c305c --- /dev/null +++ b/Cargo.lock @@ -0,0 +1,75 @@ +# This file is automatically @generated by Cargo. +# It is not intended for manual editing. +version = 3 + +[[package]] +name = "cfg-if" +version = "1.0.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "baf1de4339761588bc0619e3cbc0120ee582ebb74b53b4efbf79117bd2da40fd" + +[[package]] +name = "getrandom" +version = "0.2.11" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "fe9006bed769170c11f845cf00c7c1e9092aeb3f268e007c3e760ac68008070f" +dependencies = [ + "cfg-if", + "libc", + "wasi", +] + +[[package]] +name = "libc" +version = "0.2.150" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "89d92a4743f9a61002fae18374ed11e7973f530cb3a3255fb354818118b2203c" + +[[package]] +name = "poisson_sampling" +version = "0.1.0" +dependencies = [ + "rand", +] + +[[package]] +name = "ppv-lite86" +version = "0.2.17" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5b40af805b3121feab8a3c29f04d8ad262fa8e0561883e7653e024ae4479e6de" + +[[package]] +name = "rand" +version = "0.8.5" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "34af8d1a0e25924bc5b7c43c079c942339d8f0a8b57c39049bef581b46327404" +dependencies = [ + "libc", + "rand_chacha", + "rand_core", +] + +[[package]] +name = "rand_chacha" +version = "0.3.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e6c10a63a0fa32252be49d21e7709d4d4baf8d231c2dbce1eaa8141b9b127d88" +dependencies = [ + "ppv-lite86", + "rand_core", +] + +[[package]] +name = "rand_core" +version = "0.6.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ec0be4795e2f6a28069bec0b5ff3e2ac9bafc99e6a9a7dc3547996c5c816922c" +dependencies = [ + "getrandom", +] + +[[package]] +name = "wasi" +version = "0.11.0+wasi-snapshot-preview1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "9c8d87e72b64a3b4db28d11ce29237c246188f4f51057d65a7eab63b7987e423" diff --git a/Cargo.toml b/Cargo.toml new file mode 100644 index 0000000..08ce89f --- /dev/null +++ b/Cargo.toml @@ -0,0 +1,10 @@ +[package] +name = "poisson_sampling" +version = "0.1.0" +edition = "2021" +publish = ["artifactory"] + +# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html + +[dependencies] +rand = { version = "0.8.5", registry = "crates-io" } diff --git a/notebook/test.ipynb b/notebook/test.ipynb new file mode 100644 index 0000000..bba1231 --- /dev/null +++ b/notebook/test.ipynb @@ -0,0 +1,1547 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Poisson sampling 2D\n", + "Visualize poisson distribution on a 2D disk based on Bridson's algorithm.\n", + "\n", + "[Original paper from Bridson](https://www.cs.ubc.ca/~rbridson/docs/bridson-siggraph07-poissondisk.pdf)\n", + "\n", + "[Detailed implementation](https://sighack.com/post/poisson-disk-sampling-bridsons-algorithm)" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "vscode": { + "languageId": "rust" + } + }, + "outputs": [], + "source": [ + ":dep plotters = { version = \"^0.3.1\", default_features = false, features = [\"evcxr\", \"all_series\"] }\n", + ":dep poisson_sampling = { path = \"D:\\\\Square\\\\lordOznek\\\\poisson_sampling\\\\\" }" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": { + "vscode": { + "languageId": "rust" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Number of samples 275\n" + ] + }, + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + "Poisson sampling\n", + "\n", + "\n", + "\n", + "-1.0\n", + "\n", + "\n", + "\n", + "-0.8\n", + "\n", + "\n", + "\n", + "-0.6\n", + "\n", + "\n", + "\n", + "-0.4\n", + "\n", + "\n", + "\n", + "-0.2\n", + "\n", + "\n", + "\n", + "0.0\n", + "\n", + "\n", + "\n", + "0.2\n", + "\n", + "\n", + "\n", + "0.4\n", + "\n", + "\n", + "\n", + "0.6\n", + "\n", + "\n", + "\n", + "0.8\n", + "\n", + "\n", + "\n", + "1.0\n", + "\n", + "\n", + "\n", + "\n", + "-1.0\n", + "\n", + "\n", + "\n", + "-0.8\n", + "\n", + "\n", + "\n", + "-0.6\n", + "\n", + "\n", + "\n", + "-0.4\n", + "\n", + "\n", + "\n", + "-0.2\n", + "\n", + "\n", + "\n", + "0.0\n", + "\n", + "\n", + "\n", + "0.2\n", + "\n", + "\n", + "\n", + "0.4\n", + "\n", + "\n", + "\n", + "0.6\n", + "\n", + "\n", + "\n", + "0.8\n", + "\n", + "\n", + "\n", + "1.0\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "use poisson_sampling::*;\n", + "use plotters::prelude::*;\n", + "use plotters::coord::Shift;\n", + "use plotters::evcxr::SVGWrapper;\n", + "\n", + "// Generate points following Poisson distribution\n", + "//let points = poisson_sampling::poisson_sample_disk(0.1, 30, 1.0);\n", + "let points = poisson_sampling::poisson_sample_box2d(0.05, 30, 1.0, 1.0);\n", + "\n", + "println!( \"Number of samples {}\", points.len());\n", + "\n", + "evcxr_figure((720, 720), |root| {\n", + " root.fill(&WHITE)?;\n", + "\n", + " let mut chart = ChartBuilder::on(&root)\n", + " .caption(\n", + " \"Poisson sampling\",\n", + " (\"Arial\", 30i32),\n", + " )\n", + " .set_label_area_size(LabelAreaPosition::Left, 40i32)\n", + " .set_label_area_size(LabelAreaPosition::Bottom, 40i32)\n", + " .x_label_area_size(40)\n", + " .y_label_area_size(40)\n", + " .build_cartesian_2d(-1f64..1f64, -1f64..1f64)?;\n", + "\n", + " chart\n", + " .configure_mesh()\n", + " .disable_x_mesh()\n", + " .disable_y_mesh()\n", + " .draw()?;\n", + "\n", + " // Finally draw the scatter map\n", + " chart.draw_series(\n", + " points\n", + " .iter()\n", + " .map(|p| Circle::new((p.x, p.y), 2, GREEN.filled()))\n", + " )?;\n", + "\n", + " Ok(())\n", + "})\n" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": { + "vscode": { + "languageId": "rust" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Number of samples 865\n" + ] + }, + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + "Poisson sampling\n", + "\n", + "\n", + "\n", + "-1.0\n", + "\n", + "\n", + "\n", + "-0.8\n", + "\n", + "\n", + "\n", + "-0.6\n", + "\n", + "\n", + "\n", + "-0.4\n", + "\n", + "\n", + "\n", + "-0.2\n", + "\n", + "\n", + "\n", + "0.0\n", + "\n", + "\n", + "\n", + "0.2\n", + "\n", + "\n", + "\n", + "0.4\n", + "\n", + "\n", + "\n", + "0.6\n", + "\n", + "\n", + "\n", + "0.8\n", + "\n", + "\n", + "\n", + "1.0\n", + "\n", + "\n", + "\n", + "\n", + "-1.0\n", + "\n", + "\n", + "\n", + "-0.8\n", + "\n", + "\n", + "\n", + "-0.6\n", + "\n", + "\n", + "\n", + "-0.4\n", + "\n", + "\n", + "\n", + "-0.2\n", + "\n", + "\n", + "\n", + "0.0\n", + "\n", + "\n", + "\n", + "0.2\n", + "\n", + "\n", + "\n", + "0.4\n", + "\n", + "\n", + "\n", + "0.6\n", + "\n", + "\n", + "\n", + "0.8\n", + "\n", + "\n", + "\n", + "1.0\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "
" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "// Generate points following Poisson distribution\n", + "//let points = poisson_sampling::poisson_sample_disk(0.55, 30, 1.0); // around 8 samples\n", + "let points = poisson_sampling::poisson_sample_disk(0.05, 30, 1.0); // around 850 samples\n", + "\n", + "println!( \"Number of samples {}\", points.len());\n", + "\n", + "evcxr_figure((720, 720), |root| {\n", + " root.fill(&WHITE)?;\n", + "\n", + " let mut chart = ChartBuilder::on(&root)\n", + " .caption(\n", + " \"Poisson sampling\",\n", + " (\"Arial\", 30i32),\n", + " )\n", + " .set_label_area_size(LabelAreaPosition::Left, 40i32)\n", + " .set_label_area_size(LabelAreaPosition::Bottom, 40i32)\n", + " .x_label_area_size(40)\n", + " .y_label_area_size(40)\n", + " .build_cartesian_2d(-1f64..1f64, -1f64..1f64)?;\n", + "\n", + " chart\n", + " .configure_mesh()\n", + " .disable_x_mesh()\n", + " .disable_y_mesh()\n", + " .draw()?;\n", + "\n", + " // Finally draw the scatter map\n", + " chart.draw_series(\n", + " points\n", + " .iter()\n", + " .map(|p| Circle::new((p.x, p.y), 2, RED.filled()))\n", + " )?;\n", + "\n", + " Ok(())\n", + "})" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "vscode": { + "languageId": "rust" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "[float2(-0.5348719832393343, -0.33489885827061244), float2(-0.05261731689735866, -0.6259432734150083), float2(0.48150138477198995, -0.45963337638515933), float2(-0.12845744812013066, 0.5154937553050833), float2(0.4276802203570223, 0.13429410703728362), float2(-0.6655129837620624, 0.2530108389214989), float2(0.43596640222298666, 0.7675002636053708), float2(0.9462987844156331, 0.2561557023927691)]\n" + ] + }, + { + "data": { + "text/plain": [ + "()" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "if points.len() < 100 {\n", + " println!(\"{:?}\", points)\n", + "}" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Rust", + "language": "rust", + "name": "rust" + }, + "language_info": { + "codemirror_mode": "rust", + "file_extension": ".rs", + "mimetype": "text/rust", + "name": "Rust", + "pygment_lexer": "rust", + "version": "" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/rustfmt.toml b/rustfmt.toml new file mode 100644 index 0000000..4f5259b --- /dev/null +++ b/rustfmt.toml @@ -0,0 +1,4 @@ +fn_args_layout = "Compressed" +max_width = 120 +use_small_heuristics = "Default" +# trailing_comma = "Never" \ No newline at end of file diff --git a/src/lib.rs b/src/lib.rs new file mode 100644 index 0000000..78141b2 --- /dev/null +++ b/src/lib.rs @@ -0,0 +1,243 @@ +#![allow(dead_code)] + +// Importer le module rand pour la génération de nombres aléatoires +use rand::*; +use std::f64::consts::PI; +use std::fmt; + +// Point structure +#[derive(Default, Clone, Copy)] +pub struct Point { + pub x: f64, + pub y: f64, +} + +impl Point { + fn new() -> Self { + Default::default() + } + + fn sqr_dist(self: Self, point: Point) -> f64 { + (point.x - self.x).powi(2) + (point.y - self.y).powi(2) + } +} + +impl fmt::Debug for Point { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + // HLSL style 'copy-paste'-able + write!(f, "float2({}, {})", self.x, self.y) + } +} + +#[derive(Default)] +pub enum DomainType { + #[default] + Invalid, + Box2D { + w: f64, + h: f64, + }, + Disk { + radius: f64, + }, +} + +// Domain structure which hold a grid and the points on it +#[derive(Default)] +pub struct Domain { + pub grid: Vec>, + pub data: Vec, + pub cell_size: f64, + pub width: f64, + pub height: f64, + pub grid_width: usize, + pub grid_height: usize, + pub domain_type: DomainType, +} + +impl Domain { + pub fn box2d(width: f64, height: f64, sample_dist: f64) -> Self { + const N: u8 = 2; + + let cell_size = sample_dist / (N as f64).sqrt(); + let ncells_width = (width / cell_size).ceil() as usize; + let ncells_height = (height / cell_size).ceil() as usize; + + Domain { + cell_size: cell_size, + width: width, + height: height, + grid_width: ncells_width, + grid_height: ncells_height, + grid: vec![vec![usize::MAX; ncells_width]; ncells_height], + domain_type: DomainType::Box2D { w: width, h: height }, + ..Default::default() + } + } + + pub fn disk(radius: f64, sample_dist: f64) -> Self { + let mut domain = Domain { + domain_type: DomainType::Disk { radius: radius }, + ..Default::default() + }; + + let n = domain.get_dimension(); + let cell_size = sample_dist / (n as f64).sqrt(); + let ncells = (2.0 * radius / cell_size).ceil() as usize; + + domain.cell_size = cell_size; + domain.grid_width = ncells; + domain.grid_height = ncells; + domain.grid = vec![vec![usize::MAX; ncells]; ncells]; + + domain + } + + fn get_dimension(self: &Self) -> u8 { + #![allow(unused_variables)] + match self.domain_type { + DomainType::Box2D { w, h } => 2, + DomainType::Disk { radius } => 2, + DomainType::Invalid => 0, + } + } + + fn extent(self: &Self) -> Point { + match self.domain_type { + DomainType::Box2D { w, h } => Point { x: w, y: h }, + DomainType::Disk { radius } => Point { + x: 2.0 * radius, + y: 2.0 * radius, + }, + DomainType::Invalid => Point::default(), + } + } + + fn is_within_boundaries(self: &Self, point: Point) -> bool { + match self.domain_type { + DomainType::Box2D { w, h } => point.x.abs() < w / 2.0 && point.y.abs() < h / 2.0, + DomainType::Disk { radius } => (point.x.powi(2) + point.y.powi(2)) < radius.powi(2), + DomainType::Invalid => false, + } + } + + pub fn gen_random_point_within_domain(self: &Self) -> Point { + match self.domain_type { + DomainType::Box2D { w, h } => Point { + x: rand::random::() * w - w / 2.0, + y: rand::random::() * h - h / 2.0, + }, + DomainType::Disk { radius } => { + let sin_cos_theta = (rand::random::() * 2.0 * PI).sin_cos(); + let r = rand::random::() * radius; + Point { + x: r * sin_cos_theta.0, + y: r * sin_cos_theta.1, + } + } + DomainType::Invalid => Point::default(), + } + } + + pub fn insert_point(self: &mut Self, point: Point) { + if !self.is_within_boundaries(point) { + return; + } + + let extent = self.extent(); + let x_index = ((point.x + extent.x / 2.0) / self.cell_size).floor() as usize; + let y_index = ((point.y + extent.y / 2.0) / self.cell_size).floor() as usize; + //assert!(self.grid[y_index][x_index] == usize::MAX); + self.grid[y_index][x_index] = self.data.len(); + self.data.push(point); + } + + pub fn is_valid_point(self: &Self, point: Point, radius: f64) -> bool { + if !self.is_within_boundaries(point) { + return false; + } + let extent = self.extent(); + let x_index = ((point.x + extent.x / 2.0) / self.cell_size).floor() as usize; + let y_index = ((point.y + extent.y / 2.0) / self.cell_size).floor() as usize; + + let i0 = 0.max(x_index as i64 - 1) as usize; + let i1 = (self.grid_width - 1).min(x_index + 1); + let j0 = 0.max(y_index as i64 - 1) as usize; + let j1 = (self.grid_height - 1).min(y_index + 1); + + // Rejects the point if too close to the others + for j in j0..=j1 { + for i in i0..=i1 { + if self.grid[j][i] != usize::MAX && point.sqr_dist(self.data[self.grid[j][i]]) < radius.powi(2) { + return false; + } + } + } + + return true; + } +} + +fn poisson_sample_internal(sample_dist: f64, k: usize, domain: &mut Domain) -> Vec { + // The final list of points + let mut points = Vec::::new(); + // The current "active" list of points + let mut active = Vec::::new(); + + // Initial random point p0 + let p0 = domain.gen_random_point_within_domain(); + + domain.insert_point(p0); + points.push(p0); + active.push(p0); + + let mut rng = rand::thread_rng(); + + while active.len() > 0 { + // Pick a point 'p' from our active list + let idx = rng.gen_range(0..active.len()); + let p_active = active[idx]; + + // Try up to 'k' times to find a point that satifies: + // - is at a distance between r and 2r from p + // - is at a distance > r from nearby points + let mut i = k; + while i > 0 { + let theta = rng.gen_range(0.0..2.0 * PI); + + let r = rng.gen_range(sample_dist..2.0 * sample_dist); + + let p = Point { + x: p_active.x + r * theta.cos(), + y: p_active.y + r * theta.sin(), + }; + + // If we succeed in finding a point, add to grid and list + if domain.is_valid_point(p, sample_dist) { + points.push(p); + domain.insert_point(p); + active.push(p); + break; + } + + i -= 1; + } + + // Otherwise, remove point 'p' from the active list + if i == 0 { + active.remove(idx); + } + } + + points +} + +pub fn poisson_sample_box2d(sample_dist: f64, k: usize, width: f64, height: f64) -> Vec { + let mut domain = Domain::box2d(width, height, sample_dist); + poisson_sample_internal(sample_dist, k, &mut domain) +} + +pub fn poisson_sample_disk(sample_dist: f64, k: usize, domain_radius: f64) -> Vec { + let mut domain = Domain::disk(domain_radius, sample_dist); + poisson_sample_internal(sample_dist, k, &mut domain) +} diff --git a/src/main.rs b/src/main.rs new file mode 100644 index 0000000..d41491f --- /dev/null +++ b/src/main.rs @@ -0,0 +1,9 @@ +use poisson_sampling::*; + +fn main() { + // Tester la fonction avec un lambda de 10 et 100 points à générer + let sample = poisson_sample_box2d(0.1, 30, 1.0, 1.0); + + let sample_disk = poisson_sample_disk(0.1, 30, 1.0); + println!("Samples for disk shape {:?}", sample_disk); +}