在网格中生成位于圆内的坐标答案

【问题标题】：Generate coordinates in grid that lie within a circle在网格中生成位于圆内的坐标
【发布时间】：2017-02-13 05:37:04
【问题描述】：

我找到了this answer，这似乎与这个问题有些相关，但我想知道是否可以一一生成坐标而不会产生额外的~22% (1 - pi / 4) 损失将每个点与圆的半径进行比较（通过计算圆心与该点之间的距离）。

到目前为止，我在 Python 中有以下函数。我通过Gauss' circle problem 知道我最终会得到的坐标个数，但我也想一一生成这些点。

from typing import Iterable
from math import sqrt, floor

def circCoord(sigma: float =1.0, centroid: tuple =(0, 0)) -> Iterable[tuple]:
    r""" Generate all coords within $3\vec{\sigma}$ of the centroid """

    # The number of least iterations is given by Gauss' circle problem:
    # http://mathworld.wolfram.com/GausssCircleProblem.html

    maxiterations = 1 + 4 * floor(3 * sigma) + 4 * sum(\
      floor(sqrt(9 * sigma**2 - i**2)) for i in range(1, floor(3 * sigma) + 1)
    )

    for it in range(maxiterations):
       # `yield` points in image about `centroid` over which we loop

我想要做的是只迭代那些位于像素的 3 * sigma 范围内的像素（在上述函数中位于 centroid）。

此后，我编写了以下示例脚本，证明以下解决方案是准确的。

#! /usr/bin/env python3
# -*- coding: utf-8 -*-


import matplotlib.pyplot as plt
from matplotlib.patches import Ellipse
import numpy as np
import argparse
from typing import List, Tuple
from math import sqrt


def collect(x: int, y: int, sigma: float =3.0) -> List[Tuple[int, int]]:
    """ create a small collection of points in a neighborhood of some point 
    """
    neighborhood = []

    X = int(sigma)
    for i in range(-X, X + 1):
        Y = int(pow(sigma * sigma - i * i, 1/2))
        for j in range(-Y, Y + 1):
            neighborhood.append((x + i, y + j))

    return neighborhood


def plotter(sigma: float =3.0) -> None:
    """ Plot a binary image """    
    arr = np.zeros([sigma * 2 + 1] * 2)

    points = collect(int(sigma), int(sigma), sigma)

    # flip pixel value if it lies inside (or on) the circle
    for p in points:
        arr[p] = 1

    # plot ellipse on top of boxes to show their centroids lie inside
    circ = Ellipse(\
        xy=(int(sigma), int(sigma)), 
        width=2 * sigma,
        height=2 * sigma,
        angle=0.0
    )

    fig = plt.figure(0)
    ax  = fig.add_subplot(111, aspect='equal')
    ax.add_artist(circ)
    circ.set_clip_box(ax.bbox)
    circ.set_alpha(0.2)
    circ.set_facecolor((1, 1, 1))
    ax.set_xlim(-0.5, 2 * sigma + 0.5)
    ax.set_ylim(-0.5, 2 * sigma + 0.5)

    plt.scatter(*zip(*points), marker='.', color='white')

    # now plot the array that's been created
    plt.imshow(-arr, interpolation='none', cmap='gray')
    #plt.colorbar()

    plt.show()


if __name__ == '__main__':
    parser = argparse.ArgumentParser()

    parser.add_argument('-s', '--sigma', type=int, \
      help='Circle about which to collect points'
    )

    args = parser.parse_args()

    plotter(args.sigma)

输出

./circleCheck.py -s 4

是：

【问题讨论】：

标签： python generator

【解决方案1】：

像这样简单的东西（对于原点的圆）怎么样？

X = int(R) # R is the radius
for x in range(-X,X+1):
    Y = int((R*R-x*x)**0.5) # bound for y given x
    for y in range(-Y,Y+1):
        yield (x,y)

当圆不以原点为中心时，这可以很容易地适应一般情况。

【讨论】：

这样看起来更好，因为它不浪费。
我在上面添加了一个小测试脚本，它也可以制作一个很好的视觉“确认”图。

【解决方案2】：

您可能想考虑我的高斯圆问题算法（带有一些 Java 源代码和一个丑陋但方便的插图）：https://stackoverflow.com/a/42373448/5298879

它比您现在正在做的计算其中一个四分之一、中心和轴上的点的速度快 3.4 倍左右，同时只需要多行代码。

您只需想象一个内接正方形，然后只计算该正方形外的八分之一。

【讨论】：