[orx-noise] Add Rseq and Hammersley sequences

orx-noise/src/commonMain/kotlin/hammersley/Hammersley.kt

@@ -0,0 +1,75 @@
+package org.openrndr.extra.noise.hammersley
+
+import org.openrndr.math.Vector2
+import org.openrndr.math.Vector3
+import org.openrndr.math.Vector4
+
+/**
+ * Computes a 2D Hammersley point based on the given index and total number of samples.
+ *
+ * @param i The index of the sample, typically in the range [0, n).
+ * @param n The total number of samples.
+ * @return A 2D point as a `Vector2` within the unit square [0, 1] x [0, 1].
+ */
+fun hammersley2D(i: Int, n: Int): Vector2 {
+    return Vector2(i.toDouble() / n, radicalInverseBase2(i.toUInt()))
+}
+
+/**
+ * Computes a 3D point in the Hammersley sequence based on the given index and total number of samples.
+ *
+ * @param i The index of the sample, typically in the range [0, n).
+ * @param n The total number of samples.
+ * @return A 3D point as a `Vector3` within the unit cube [0, 1] x [0, 1] x [0, 1].
+ */
+fun hammersley3D(i: Int, n: Int): Vector3 {
+    return Vector3(i.toDouble() / n, radicalInverseBase2(i.toUInt()), radicalInverse(3, i))
+}
+
+/**
+ * Computes a 4D Hammersley point based on the given index and total number of samples.
+ *
+ * @param i The index of the sample, typically in the range [0, n).
+ * @param n The total number of samples.
+ * @return A 4D point as a `Vector4` where each component lies within the range [0, 1].
+ */
+fun hammersley4D(i: Int, n: Int): Vector4 {
+    return Vector4(i.toDouble() / n, radicalInverseBase2(i.toUInt()), radicalInverse(3, i), radicalInverse(5, i))
+}
+
+/**
+ * Computes the radical inverse of a given unsigned integer `i` in base 2.
+ *
+ * @param i The input unsigned integer for which the radical inverse in base 2 is computed.
+ * @return The radical inverse value of the input as a `Double`, mapped to the range [0, 1).
+ */
+fun radicalInverseBase2(i: UInt): Double {
+    var bits = i
+    bits = ((bits shl 16) or (bits shr 16))
+    bits = ((bits and 0x55555555u) shl 1) or ((bits and 0xAAAAAAAAu) shr 1)
+    bits = ((bits and 0x33333333u) shl 2) or ((bits and 0xCCCCCCCCu) shr 2)
+    bits = ((bits and 0x0F0F0F0Fu) shl 4) or ((bits and 0xF0F0F0F0u) shr 4)
+    bits = ((bits and 0x00FF00FFu) shl 8) or ((bits and 0xFF00FF00u) shr 8)
+    return bits.toDouble() * 2.3283064365386963e-10
+}
+
+/**
+ * Computes the radical inverse of an integer `i` in a given base.
+ * This method is often used in quasi-random sequence generation for sampling.
+ *
+ * @param base The base in which to compute the radical inverse. Must be greater than 1.
+ * @param i The integer for which the radical inverse is calculated. Must be non-negative.
+ * @return The radical inverse value as a `Double`, within the range [0, 1).
+ */
+fun radicalInverse(base: Int, i: Int): Double {
+    var v = 0.0
+    var denom = 1.0
+    var n = i
+    while (n > 0) {
+        denom *= base
+        val remainder = n.mod(base)
+        n /= base
+        v += remainder / denom
+    }
+    return v
+}

orx-noise/src/commonMain/kotlin/rseq/Rseq.kt

@@ -0,0 +1,67 @@
+package org.openrndr.extra.noise.rsequence
+
+import org.openrndr.math.Vector2
+import org.openrndr.math.Vector3
+import org.openrndr.math.Vector4
+
+private const val g1 = 1.618033988749895
+private const val a11 = 1.0 / g1
+
+private const val g2 = 1.324717957244746
+private const val a21 = 1.0 / g2
+private const val a22 = 1.0 / (g2 * g2)
+
+private const val g3 = 1.2207440846057596
+private const val a31 = 1.0 / g3
+private const val a32 = 1.0 / (g3 * g3)
+private const val a33 = 1.0 / (g3 * g3 * g3)
+
+private const val g4 = 1.1673039782614187
+private const val a41 = 1.0 / g4
+private const val a42 = 1.0 / (g4 * g4)
+private const val a43 = 1.0 / (g4 * g4 * g4)
+private const val a44 = 1.0 / (g4 * g4 * g4 * g4)
+
+/**
+ * Computes the R1 low-discrepancy quasirandom sequence value for a given index as described by Martin Roberts.
+ *
+ * @param n The index for which the R1 sequence value is to be calculated.
+ * @return The R1 sequence value as a Double, providing a low-discrepancy quasirandom number.
+ */
+fun rSeq1D(n: Int): Double = (0.5 + a11 * n).mod(1.0)
+
+/**
+ * Computes the R2 low-discrepancy quasirandom sequence value for a given index as described by Martin Roberts.
+ *
+ * @param n The index for which the R2 sequence value is to be calculated.
+ * @return The R2 sequence value as a [Vector2], providing a low-discrepancy quasirandom number.
+ */
+fun rSeq2D(n: Int): Vector2 = Vector2(
+    (0.5 + a21 * n).mod(1.0),
+    (0.5 + a22 * n).mod(1.0)
+)
+
+/**
+ * Computes the R3 low-discrepancy quasirandom sequence value for a given index as described by Martin Roberts.
+ *
+ * @param n The index for which the R3 sequence value is to be calculated.
+ * @return The R3 sequence value as a [Vector3], providing a low-discrepancy quasirandom number.
+ */
+fun rSeq3D(n: Int): Vector3 = Vector3(
+    (0.5 + a31 * n).mod(1.0),
+    (0.5 + a32 * n).mod(1.0),
+    (0.5 + a33 * n).mod(1.0)
+)
+
+/**
+ * Computes the R4 low-discrepancy quasirandom sequence value for a given index as described by Martin Roberts.
+ *
+ * @param n The index for which the R4 sequence value is to be calculated.
+ * @return The R4 sequence value as a [Vector4], providing a low-discrepancy quasirandom number.
+ */
+fun rSeq4D(n: Int): Vector4 = Vector4(
+    (0.5 + a41 * n).mod(1.0),
+    (0.5 + a42 * n).mod(1.0),
+    (0.5 + a43 * n).mod(1.0),
+    (0.5 + a44 * n).mod(1.0)
+)