removing SOA; padding Jittered + LP sampler

CMakeLists.txt

@@ -272,7 +272,6 @@ add_library(
   include/sampler/Random.h
   include/sampler/RandomPermutation.h
   include/sampler/Sampler.h
-  include/sampler/SOA.h
   include/sampler/Sobol.h
   include/sampler/Sudoku.h
   include/sampler/xi.h
@@ -294,7 +293,6 @@ add_library(
   src/sampler/OABush.cpp
   src/sampler/OACMJND.cpp
   src/sampler/Random.cpp
-  src/sampler/SOA.cpp
   src/sampler/Sobol.cpp
   src/sampler/Sudoku.cpp
   src/sampler/XiSequence.cpp

include/sampler/GrayCode.h

@@ -3,7 +3,6 @@
 */
 #pragma once
 
-#include <memory>
 #include <pcg32.h>
 #include <sampler/Sampler.h>
 #include <vector>
@@ -18,7 +17,7 @@
 class GrayCode : public TSamplerDim<2>
 {
 public:
-    GrayCode(unsigned n = 3);
+    GrayCode(unsigned n = 2);
 
     void sample(float[], unsigned i) override;
 
@@ -31,17 +30,7 @@ class GrayCode : public TSamplerDim<2>
     {
         return N;
     }
-    int setNumSamples(unsigned num) override
-    {
-        int log2N = std::round(std::log2(num));
-        if (log2N & 1)
-            log2N++;        // Only allow even powers of 2.
-        N     = 1 << log2N; // Make N a power of 4, if not.
-        log2n = log2N / 2;
-        n     = 1 << log2n;
-        regenerate();
-        return N;
-    }
+    int setNumSamples(unsigned num) override;
 
     bool randomized() const override
     {

include/sampler/Jittered.h

@@ -7,35 +7,41 @@
 #include <sampler/Sampler.h>
 
 /// Encapsulate a 2D stratified or "jittered" point set.
-class Jittered : public TSamplerDim<2>
+class Jittered : public TSamplerMinMaxDim<1, 1024>
 {
 public:
     Jittered(unsigned resX, unsigned resY, float jitter = 1.0f);
 
-    void reset() override;
     void sample(float[], unsigned i) override;
 
-    std::string name() const override
+    unsigned dimensions() const override
     {
-        return "Jittered";
+        return m_numDimensions;
+    }
+    void setDimensions(unsigned n) override
+    {
+        m_numDimensions = n;
     }
 
     int numSamples() const override
     {
-        return m_resX * m_resY;
+        return m_numSamples;
     }
     int setNumSamples(unsigned n) override
     {
+        if (n == m_numSamples)
+            return m_numSamples;
         int sqrtVal = (n == 0) ? 1 : (int)(std::sqrt((float)n) + 0.5f);
-        m_resX = m_resY = sqrtVal;
-        reset();
-        return m_resX * m_resY;
+        setNumSamples(sqrtVal, sqrtVal);
+        return m_numSamples;
     }
     void setNumSamples(unsigned x, unsigned y)
     {
-        m_resX = x;
-        m_resY = y;
-        reset();
+        m_resX       = x == 0 ? 1 : x;
+        m_resY       = y == 0 ? 1 : y;
+        m_numSamples = m_resX * m_resY;
+        m_xScale     = 1.0f / m_resX;
+        m_yScale     = 1.0f / m_resY;
     }
 
     bool randomized() const override
@@ -44,8 +50,9 @@ class Jittered : public TSamplerDim<2>
     }
     void setRandomized(bool r = true) override
     {
-        if ((m_randomize = r))
-            m_seed++;
+        m_randomize = r;
+        m_rand.seed(m_permutation);
+        m_permutation = r ? m_rand.nextUInt() : 0;
     }
 
     /// Gets/sets how much the points are jittered
@@ -58,12 +65,18 @@ class Jittered : public TSamplerDim<2>
         return m_maxJit = j;
     }
 
+    std::string name() const override
+    {
+        return "Jittered";
+    }
+
 private:
-    unsigned m_resX, m_resY;
+    unsigned m_resX, m_resY, m_numSamples, m_numDimensions;
     float    m_maxJit;
-
-    float    m_xScale, m_yScale;
     bool     m_randomize = true;
     pcg32    m_rand;
-    unsigned m_seed = 13;
+    unsigned m_seed        = 13;
+    unsigned m_permutation = 13;
+
+    float m_xScale, m_yScale;
 };

include/sampler/LP.h

@@ -24,7 +24,7 @@
     in P. L'Ecuyer and A. Owen (eds.),
     Monte Carlo and Quasi-Monte Carlo Methods 2008, Springer-Verlag, 2009.
 */
-class LarcherPillichshammerGK : public TSamplerMinMaxDim<1, 3>
+class LarcherPillichshammerGK : public TSamplerMinMaxDim<1, 1024>
 {
 public:
     LarcherPillichshammerGK(unsigned dimensions = 2, unsigned numSamples = 64, bool randomize = false);

include/sampler/MultiJittered.h

@@ -20,7 +20,7 @@
 class MultiJittered : public TSamplerDim<2>
 {
 public:
-    MultiJittered(unsigned, unsigned, bool randomize = true, float jitter = 0.0f);
+    MultiJittered(unsigned x, unsigned y, bool randomize = true, float jitter = 0.0f);
     ~MultiJittered() override;
     void clear();
 
@@ -225,9 +225,7 @@ class CorrelatedMultiJitteredInPlace : public TSamplerMinMaxDim<1, 1024>
     }
     float setJitter(float j = 1.0f) override
     {
-        m_maxJit = j;
-        reset();
-        return m_maxJit;
+        return m_maxJit = j;
     }
 
     std::string name() const override
@@ -238,9 +236,9 @@ class CorrelatedMultiJitteredInPlace : public TSamplerMinMaxDim<1, 1024>
 protected:
     unsigned m_resX, m_resY, m_numSamples, m_numDimensions;
     float    m_maxJit;
-    bool     m_randomize;
+    bool     m_randomize = true;
     pcg32    m_rand;
-    unsigned m_seed = 13;
-    unsigned m_permutation;
+    unsigned m_seed        = 13;
+    unsigned m_permutation = 13;
     unsigned m_decorrelate;
 };

src/app.cpp

@@ -117,7 +117,7 @@ SampleViewer::SampleViewer()
     m_samplers.emplace_back(new Hammersley<Halton>(m_num_dimensions, 1));
     m_samplers.emplace_back(new Hammersley<HaltonZaremba>(m_num_dimensions, 1));
     m_samplers.emplace_back(new LarcherPillichshammerGK(3, 1, false));
-    m_samplers.emplace_back(new GrayCode(4));
+    m_samplers.emplace_back(new GrayCode(1));
     m_samplers.emplace_back(new XiSequence(1));
     m_samplers.emplace_back(new CSVFile());
 
@@ -602,7 +602,7 @@ void SampleViewer::draw_editor()
     if (big_header(ICON_FA_SLIDERS_H "  Sampler settings"))
     // =========================================================
     {
-        if (ImGui::BeginCombo("##Sampler combo", m_samplers[m_sampler]->name().c_str()))
+        if (ImGui::BeginCombo("##Sampler combo", m_samplers[m_sampler]->name().c_str(), ImGuiComboFlags_HeightLargest))
         {
             for (int n = 0; n < (int)m_samplers.size(); n++)
             {

src/sampler/GrayCode.cpp

@@ -4,6 +4,7 @@
 
 #include <bitcount.h>
 #include <sampler/GrayCode.h>
+#include <sampler/Misc.h>
 
 GrayCode::GrayCode(unsigned n)
 {
@@ -17,9 +18,11 @@ void GrayCode::regenerate()
     m_samples.clear();
     m_samples.resize(N);
 
-    double res = 1.0 / N;
+    double offset = m_randomize ? 0. : 0.5;
+    double res    = 1.0 / N;
+    // Sample in first stratum of each row and column
     for (unsigned V = 0; V < n; V++)
-        m_samples[V].x = m_samples[V * n].y = res * (V * n); // Sample in first stratum of each row and column
+        m_samples[V].x = m_samples[V * n].y = res * (V * n + offset);
 
     uint32_t u = 0;
     // Iterate through strata up the column
@@ -39,8 +42,29 @@ void GrayCode::regenerate()
 
         // Iterate through columns
         for (uint32_t V = 0; V < n; V++)
-            m_samples[U * n + V].x = m_samples[V * n + U].y = res * (V * n + u);
+            m_samples[U * n + V].x = m_samples[V * n + U].y = res * (V * n + u + offset);
     }
+
+    if (m_randomize)
+    {
+        // perform additional xor scrambling
+        auto s1 = m_randomize ? m_rand.nextUInt() : 0;
+        auto s2 = m_randomize ? m_rand.nextUInt() : 0;
+        for (unsigned i = 0; i < N; ++i)
+            m_samples[i] = Point{randomDigitScramble(m_samples[i].x, s1), randomDigitScramble(m_samples[i].y, s2)};
+    }
+}
+
+int GrayCode::setNumSamples(unsigned num)
+{
+    int log2N = std::round(std::log2(num));
+    if (log2N & 1)
+        log2N++;        // Only allow even powers of 2.
+    N     = 1 << log2N; // Make N a power of 4, if not.
+    log2n = log2N / 2;
+    n     = 1 << log2n;
+    regenerate();
+    return N;
 }
 
 void GrayCode::sample(float r[], unsigned i)

src/sampler/Jittered.cpp

@@ -3,29 +3,34 @@
 */
 
 #include <sampler/Jittered.h>
+#include <sampler/Misc.h> // for permute
 
-Jittered::Jittered(unsigned x, unsigned y, float jitter) : m_resX(x), m_resY(y), m_maxJit(jitter)
+Jittered::Jittered(unsigned x, unsigned y, float jitter) : m_maxJit(jitter)
 {
-    reset();
-}
-
-void Jittered::reset()
-{
-    if (m_resX == 0)
-        m_resX = 1;
-    if (m_resY == 0)
-        m_resY = 1;
-    m_xScale = 1.0f / m_resX;
-    m_yScale = 1.0f / m_resY;
-    // m_rand.seed(m_seed);
+    setNumSamples(x, y);
 }
 
 void Jittered::sample(float r[], unsigned i)
 {
+    i %= m_numSamples;
+
     if (i == 0)
         m_rand.seed(m_seed);
-    float jx = 0.5f + int(m_randomize) * m_maxJit * (m_rand.nextFloat() - 0.5f);
-    float jy = 0.5f + int(m_randomize) * m_maxJit * (m_rand.nextFloat() - 0.5f);
-    r[0]     = ((i % m_resX) + jx) * m_xScale;
-    r[1]     = ((i / m_resY) + jy) * m_yScale;
+
+    for (unsigned d = 0; d < dimensions(); d += 2)
+    {
+        int s = permute(i, m_numSamples, m_permutation * 0x51633e2d * (d + 1));
+
+        // horizontal and vertical indices of the stratum in the jittered grid
+        int x = s % m_resX;
+        int y = s / m_resX;
+
+        // jitter in the d and d+1 dimensions
+        float jx = 0.5f + m_randomize * m_maxJit * (m_rand.nextFloat() - 0.5f);
+        float jy = 0.5f + m_randomize * m_maxJit * (m_rand.nextFloat() - 0.5f);
+
+        r[d] = (x + jx) * m_xScale;
+        if (d + 1 < dimensions())
+            r[d + 1] = (y + jy) * m_yScale;
+    }
 }