c++ - 具有特定最小距离的特定正方形(而非单位正方形)中的二维泊松盘采样

Question

有什么方法可以修改泊松圆盘点生成器发现 here .我需要使用 textfile.txt 中的点坐标生成新的泊松点以改善分布。下面是单位正方形泊松盘采样的c++代码。

泊松生成器.h:

#include <vector>
#include <random>
#include <stdint.h>
#include <time.h>


namespace PoissoGenerator
{


class DefaultPRNG
{
public:
    DefaultPRNG()
        : m_Gen(std::random_device()())
        , m_Dis(0.0f, 1.f)
    {
        // prepare PRNG
        m_Gen.seed(time(nullptr));
    }

    explicit DefaultPRNG(unsigned short seed)
        : m_Gen(seed)
        , m_Dis(0.0f, 1.f)
    {

    }

    double RandomDouble()
    {
        return static_cast <double>(m_Dis(m_Gen));
    }

    int RandomInt(int Max)
    {
        std::uniform_int_distribution<> DisInt(0, Max);
        return DisInt(m_Gen);
    }

private:
    std::mt19937 m_Gen;
    std::uniform_real_distribution<double> m_Dis;
};



struct sPoint
{
    sPoint()
    : x(0)
    , y(0)
    , m_valid(false){}
    sPoint(double X, double Y)
        : x(X)
        , y(Y)
        , m_valid(true){}
    double x;
    double y;
    bool m_valid;
    //
    bool IsInRectangle() const
    {
        return x >= 0 && y >= 0 && x <= 1 && y <= 1;
    }
    //
    bool IsInCircle() const
    {
        double fx = x - 0.5f;
        double fy = y - 0.5f;
        return (fx*fx + fy*fy) <= 0.25f;

    }


};

struct sGridPoint
{
    sGridPoint(int X, int Y)
    : x(X)
    , y(Y)
    {}
    int x;
    int y;
};

double GetDistance(const sPoint& P1, const sPoint& P2)
{
    return sqrt((P1.x - P2.x)*(P1.x - P2.x) + (P1.y - P2.y)*(P1.y - P2.y));
}

sGridPoint ImageToGrid(const sPoint& P, double CellSize)
{
    return sGridPoint((int)(P.x / CellSize), (int)(P.y / CellSize));
}

struct  sGrid
{
    sGrid(int W, int H, double CellSize)
    : m_W(W)
    , m_H(H)
    , m_CellSize(CellSize)
    {
        m_Grid.resize((m_H));

        for (auto i = m_Grid.begin(); i != m_Grid.end(); i++){ i->resize(m_W); }
    }

    void Insert(const sPoint& P)
    {
        sGridPoint G = ImageToGrid(P, m_CellSize);
        m_Grid[G.x][G.y] = P;
    }
    bool IsInNeighbourhood(sPoint Point, double MinDist, double CellSize)
    {
        sGridPoint G = ImageToGrid(Point, CellSize);

        //number of adjacent cell to look for neighbour points
        const int D = 5;

        // Scan the neighbourhood of the Point in the grid 
        for (int i = G.x - D; i < G.x + D; i++)
        {
            for (int j = G.y - D; j < G.y + D; j++)
            {
                if (i >= 0 && i < m_W && j >= 0 && j < m_H)
                {
                    sPoint P = m_Grid[i][j];
                    if (P.m_valid && GetDistance(P, Point) < MinDist){ return true; }
                }
            }
        }

        return false;
    }

private:
    int m_H;
    int m_W;
    double m_CellSize;

    std::vector< std::vector< sPoint> > m_Grid;
};

template <typename PRNG>
sPoint PopRandom(std::vector<sPoint>& Points, PRNG& Generator)
{
    const int Idx = Generator.RandomInt(Points.size() - 1);
    const sPoint P = Points[Idx];
    Points.erase(Points.begin() + Idx);
    return P;
}

template <typename PRNG>
sPoint GenerateRandomPointAround(const sPoint& P, double MinDist, PRNG& Generator)
{

    // Start with non-uniform distribution
    double R1 = Generator.RandomDouble();
    double R2 = Generator.RandomDouble();

    // radius should be between MinDist and 2 * MinDist
    double Radius = MinDist * (R1 + 1.0f);

    //random angle
    double Angle = 2 * 3.141592653589f * R2;

    // the new point is generated around the point (x, y)
    double X = P.x + Radius * cos(Angle);
    double Y = P.y + Radius * sin(Angle);

    return sPoint(X, Y);
}

// Return a vector of generated points
// NewPointsCount - refer to bridson-siggraph07-poissondisk.pdf
// for details (the value 'k')
// Circle - 'true' to fill a circle, 'false' to fill a rectangle
// MinDist - minimal distance estimator, use negative value for default

template <typename PRNG = DefaultPRNG>
std::vector<sPoint> GeneratePoissonPoints(rsize_t NumPoints, PRNG& Generator, int NewPointsCount = 30,
    bool Circle = true, double MinDist = -1.0f)
{
    if (MinDist < 0.0f)
    {
        MinDist = sqrt(double(NumPoints)) / double(NumPoints);
    }

    std::vector <sPoint> SamplePoints;
    std::vector <sPoint> ProcessList;

    // create the grid
    double CellSize = MinDist / sqrt(2.0f);

    int GridW = (int)(ceil)(1.0f / CellSize);
    int GridH = (int)(ceil)(1.0f / CellSize);

    sGrid Grid(GridW, GridH, CellSize);



    sPoint FirstPoint;
    do
    {
        FirstPoint = sPoint(Generator.RandomDouble(), Generator.RandomDouble());

    } while (!(Circle ? FirstPoint.IsInCircle() : FirstPoint.IsInRectangle()));

    //Update containers

    ProcessList.push_back(FirstPoint);
    SamplePoints.push_back(FirstPoint);
    Grid.Insert(FirstPoint);


    // generate new points for each point in the queue

    while (!ProcessList.empty() && SamplePoints.size() < NumPoints)
    {
             #if POISSON_PROGRESS_INDICATOR
        // a progress indicator, kind of 
        if (SamplePoints.size() % 100 == 0) std::cout << ".";
             #endif // POISSON_PROGRESS_INDICATOR

        sPoint Point = PopRandom<PRNG>(ProcessList, Generator);

        for (int i = 0; i < NewPointsCount; i++)
        {
            sPoint NewPoint = GenerateRandomPointAround(Point, MinDist, Generator);
            bool Fits = Circle ? NewPoint.IsInCircle() : NewPoint.IsInRectangle();

            if (Fits && !Grid.IsInNeighbourhood(NewPoint, MinDist, CellSize))
            {
                ProcessList.push_back(NewPoint);
                SamplePoints.push_back(NewPoint);
                Grid.Insert(NewPoint);
                continue;

            }
        }


    }

         #if POISSON_PROGRESS_INDICATOR
         std::cout << std::endl << std::endl;
        #endif // POISSON_PROGRESS_INDICATOR

    return SamplePoints;

}

主要程序是:

泊松.cpp

 #include "stdafx.h"
 #include <vector>
 #include <iostream>
 #include <fstream>
 #include <memory.h>


 #define POISSON_PROGRESS_INDICATOR 1
 #include "PoissonGenerator.h"


  const int   NumPoints = 20000;    // minimal number of points to generate


int main()
{

 PoissonGenerator::DefaultPRNG PRNG;

 const auto Points =
 PoissonGenerator::GeneratePoissonPoints(NumPoints,PRNG);


std::ofstream File("Poisson.txt", std::ios::out);

File << "NumPoints = " << Points.size() << std::endl;
for (const auto& p : Points)
{
File << "       " << p.x << "       " << p.y << std::endl;
}

system("PAUSE");
return 0;
  }

Answer 1

假设空间中有一个点 [0,1] x [0,1] , 形式为 std::pair<double, double> , 但渴望空间中的点 [x,y] x [w,z] .

函数对象

struct ProjectTo {
    double x, y, w, z;
    std::pair<double, double> operator(std::pair<double, double> in)
    {
        return std::make_pair(in.first * (y - x) + x, in.second * (z - w) + w);
    }
};

会将这样的输入点转换为所需的输出点。

进一步假设您有一个 std::vector<std::pair<double, double>> points , 全部取自输入分布。

std::copy(points.begin(), points.end(), points.begin(), ProjectTo{ x, y, w, z });

现在您在输出空间中有了一个点 vector 。