programing

구면에서 n개의 점을 고르게 분포

linuxpc 2023. 7. 15. 09:45
반응형

구면에서 n개의 점을 고르게 분포

저는 N점(아마도 20점 미만)에 대해 구 주위의 위치를 제공할 수 있는 알고리즘이 필요합니다."완벽함"은 필요하지 않지만, 저는 그것이 필요하기 때문에 그들 중 어떤 것도 함께 뭉치지 않습니다.

  • 이 질문은 좋은 코드를 제공했지만, 이것이 100% 무작위화된 것처럼 보여서 이 유니폼을 만드는 방법을 찾을 수 없었습니다.
  • 이 블로그 게시물은 구체의 포인트 수를 입력할 수 있는 두 가지 방법이 있었지만, Saff와 Kuijlars 알고리즘은 정확히 내가 기록할 수 있는 psuedocode이고, 내가 찾은 코드 예제는 "node[k]"를 포함하고 있어 설명할 수 없고 그 가능성을 망쳤습니다.두 번째 블로그의 예는 골든 섹션 스파이럴이었는데, 일정한 반경을 정의할 수 있는 명확한 방법이 없이 이상하고 뭉친 결과를 저에게 주었습니다.
  • 질문에서 나온 알고리즘은 작동할 수 있을 것 같지만, 저는 그 페이지에 있는 것들을 통합해서 코드를 만들 수는 없습니다.

제가 만난 몇 가지 다른 질문 스레드는 무작위 균일 분포에 대해 언급했는데, 이는 제가 걱정하지 않는 복잡성 수준을 더합니다.너무 바보 같은 질문이라 죄송합니다만, 저는 제가 정말 열심히 보았지만 여전히 부족하다는 것을 보여주고 싶었습니다.

그래서, 제가 찾고 있는 것은 N개의 점들을 구면 또는 데카르트 좌표로 반환하는 단위 구 주위에 균등하게 분포시키는 간단한 의사 코드입니다.만약 그것이 약간의 무작위화(별 주위의 행성들이 적절하게 퍼져 있지만 여유가 있다고 생각하면)로도 분포할 수 있다면 훨씬 더 좋습니다.

Fibonacci알고리즘은 이에 매우 유용합니다.그것은 빠르고 한눈에 사람의 눈을 쉽게 속일 수 있는 결과를 제공합니다.처리가 완료된 예제를 볼 수 있으며, 시간이 지남에 따라 점이 추가됨에 따라 결과가 표시됩니다.여기 @gman에 의해 만들어진 또 다른 훌륭한 상호작용 사례가 있습니다.여기 파이썬의 간단한 구현이 있습니다.

import math


def fibonacci_sphere(samples=1000):

    points = []
    phi = math.pi * (math.sqrt(5.) - 1.)  # golden angle in radians

    for i in range(samples):
        y = 1 - (i / float(samples - 1)) * 2  # y goes from 1 to -1
        radius = math.sqrt(1 - y * y)  # radius at y

        theta = phi * i  # golden angle increment

        x = math.cos(theta) * radius
        z = math.sin(theta) * radius

        points.append((x, y, z))

    return points

1000개의 샘플은 다음을 제공합니다.

enter image description here

황금나선법

당신은 골든 스파이럴 방법을 사용할 수 없다고 말했고 그것은 정말 정말 좋기 때문에 유감입니다.저는 당신에게 그것에 대한 완전한 이해를 주고 싶습니다. 그래서 아마도 당신은 이것이 "뭉개지는" 것을 막는 방법을 이해할 수 있을 것입니다.

그래서 여기 거의 정확한 격자를 만드는 빠르고, 무작위적이지 않은 방법이 있습니다. 위에서 논의했듯이, 어떤 격자도 완벽하지 않을 것이지만, 이것은 충분할 것입니다.BendWavy.org 같은 다른 방법과 비교되지만, 이것은 보기 좋고 예쁜 것뿐만 아니라 한도 내의 짝수 간격에 대한 보장을 가지고 있습니다.

프라이머: 장치 디스크의 해바라기 나선형

이 알고리즘을 이해하기 위해 먼저 2D 해바라기 나선 알고리즘을 살펴보도록 하겠습니다.은 가장 이는가불합한숫자가황리금이사것근다입실니거한에는라율비장▁라는 사실에 근거한 것입니다.(1 + sqrt(5))/2그리고 "중심에 서서, 전체 회전의 황금 비율을 돌린 다음 그 방향으로 다른 점을 방출"하는 접근법에 의해 점을 방출하면, 점의 수가 점점 더 많아질수록, 그럼에도 불구하고 점이 (Note 1.)정렬되는 잘 정의된 '바'를 갖는 것을 거부하는 나선형을 자연스럽게 구성합니다.

디스크의 짝수 간격 알고리즘은 다음과 같습니다.

from numpy import pi, cos, sin, sqrt, arange
import matplotlib.pyplot as pp

num_pts = 100
indices = arange(0, num_pts, dtype=float) + 0.5

r = sqrt(indices/num_pts)
theta = pi * (1 + 5**0.5) * indices

pp.scatter(r*cos(theta), r*sin(theta))
pp.show()

(n=100 및 n=1000)과 같은 결과가 생성됩니다.

enter image description here

방사형으로 점 간격 지정

가장 이상한 것은 공식입니다.r = sqrt(indices / num_pts)내가 어떻게 그것에 도달했을까요?(Note 2.)

여기서 제곱근을 사용하는 이유는 디스크 주변에 균등한 공간이 있기를 원하기 때문입니다.이는 N의 한계에서 작은 영역 R ∈(r, r + dr), θ(θ, θ + dr)이 영역에 비례하는 점들의 수를 포함하기를 원한다는 것과 같습니다. 는 r r rdr입니다.만약 우리가 여기서 무작위 변수에 대해 이야기하고 있다고 가정한다면, 이것은 (R, Δ)에 대한 공동 확률 밀도가 어떤 상수 c에 대해 cr일 뿐이라는 간단한 해석을 가지고 있습니다.그러면 장치 디스크에서 정규화하면 c = 1/dll이 강제됩니다.

이제 한 가지 요령을 소개하겠습니다.는 역CDF를 샘플링하는 것으로 알려진 확률 이론에서 비롯됩니다. 확률 밀도 f(z)를 갖는 무작위 변수를 생성하고 다음과 같이 임의 변수 U ~ 균일(0, 1)을 갖는다고 가정합니다.random()대부분의 프로그래밍 언어에서.이걸 어떻게 하나요?

  1. 먼저 밀도를 누적 분포 함수 또는 CDF로 변환합니다. 이를 F(z)라고 합니다.CDF는 도함수 f(z)를 사용하면 0에서 1로 단조롭게 증가합니다.
  2. 그런 다음 CDF의 역함수-1 F(z)를 계산합니다.
  3. Z = F-1(U)가 목표 밀도에 따라 분포되어 있음을 알 수 있습니다.(Note 3).

이제 황금색 나선형 트릭은 π에 대해 점들을 잘 고른 패턴으로 공간을 만들어내자. 단위 디스크에 대해 우리는 F(r) = r2 남게 된다. 그래서 역함수는 F-1(u) = u이고1/2 따라서 우리는 다음과 같은 극좌표로 디스크에 임의의 점들을 생성할 것입니다.r = sqrt(random()); theta = 2 * pi * random().

이제 우리는 이 역함수를 무작위로 샘플링하는 대신 균일하게 샘플링합니다. 그리고 균일한 샘플링의 좋은 점은 큰 N의 한계에서 어떻게 점들이 퍼져있는지에 대한 우리의 결과가 무작위로 샘플링한 것처럼 행동한다는 것입니다.이 조합이 비결입니다.대신에random()우리는 사용합니다.(arange(0, num_pts, dtype=float) + 0.5)/num_pts예를 그서래, 말자면점, 만우가리 10점, 은고싶하면그, 들약하다.r = 0.05, 0.15, 0.25, ... 0.95동일한 면적 간격을 얻기 위해 균일하게 r을 샘플링하고, 출력에서 끔찍한 "바"를 피하기 위해 해바라기 증분을 사용합니다.

지금 구 위에서 해바라기를 하고 있습니다.

우리가 구에 점을 찍기 위해 해야 할 변화는 구면 좌표에 대한 극좌표를 바꾸는 것을 포함합니다.물론 방사 좌표는 여기에 포함되지 않습니다. 왜냐하면 우리는 단위 구에 있기 때문입니다.여기서 좀 더 일관성을 유지하기 위해, 물리학자로서 훈련을 받았음에도 불구하고, 저는 수학자들의 좌표를 사용할 것입니다. 여기서 0 ≤ ≤ ≤ ≤ ≤ 2 ≤는 극에서 내려오는 위도이고, 0 ≤ ≤ ≤ 경도입니다.따라서 위와 다른 점은 기본적으로 변수 r을 φ로 대체하고 있다는 것입니다.

rdrdθ이었던 우리의 면적 요소는 이제 그다지 복잡하지 않은 죄(φ)dθ가 됩니다.따라서 균일한 간격을 위한 우리의 접합 밀도는 sin(θ)/4θ입니다.π를 적분하면, 우리는 f(π) = sin(π)/2, 즉 F(π) = (1 - cos(π))/2를 찾습니다.이를 뒤집으면 균일한 무작위 변수가 코스(1 - 2 u)처럼 보이지만 무작위로 샘플링하는 대신 균일하게 샘플링하므로 대신 φk = 아코스(1 - 2(k + 0.5)/N)를 사용합니다.나머지 알고리즘은 이것을 x, y, z 좌표에 투영합니다.

from numpy import pi, cos, sin, arccos, arange
import mpl_toolkits.mplot3d
import matplotlib.pyplot as pp

num_pts = 1000
indices = arange(0, num_pts, dtype=float) + 0.5

phi = arccos(1 - 2*indices/num_pts)
theta = pi * (1 + 5**0.5) * indices

x, y, z = cos(theta) * sin(phi), sin(theta) * sin(phi), cos(phi);

pp.figure().add_subplot(111, projection='3d').scatter(x, y, z);
pp.show()

다시 n=100 및 n=1000의 경우 결과는 다음과 같습니다.

추가 연구

저는 마틴 로버츠의 블로그에 큰 소리를 내고 싶었습니다.위에서 각 인덱스에 0.5를 추가하여 인덱스의 오프셋을 만들었습니다.이것은 단지 시각적으로 저에게 매력적이었지만 오프셋의 선택은 매우 중요하고 간격에 걸쳐 일정하지 않으며 올바르게 선택하면 포장의 정확도가 8% 향상될 수 있다는 것으로 나타났습니다.또한 2 R 시퀀스가 구체를 덮을 수 있는 방법이 있어야 하며, 이것이 또한 좋은 고른 덮개를 생성했는지 여부를 보는 것은 흥미로울 것입니다. 아마도 그대로이지만 대각선 또는 그 이상으로 절단된 단위 사각형의 절반에서만 가져와 원을 얻기 위해 주변에 늘어져야 합니다.

메모들

  1. 이러한 "bar"는 숫자에 대한 합리적인 근사치에 의해 형성되며, 숫자에 대한 최상의 합리적 근사치는 연속 분수 표현에서 나옵니다.z + 1/(n_1 + 1/(n_2 + 1/(n_3 + ...)))z이고 는정이고입니다.n_1, n_2, n_3, ...정수의 입니다.

    def continued_fraction(r):
        while r != 0:
            n = floor(r)
            yield n
            r = 1/(r - n)
    

    분수 부분 이후로1/(...)항상 0과 1 사이에 있으며, 연속 분수의 큰 정수를 사용하면 특히 우수한 합리적인 근사치를 얻을 수 있습니다. 것으로 것"이 "과 2 것으로 것"보다 "1과 100과 2 사이의 어떤 것으로 나눈 것", "1과 100과 101 사이의 어떤 것으로 나눈 것", "1과 2 사이의 것", "1과 2 사이의 것"보다 좋습니다.그러므로 가장 비합리적인 숫자는 다음과 같은 것입니다.1 + 1/(1 + 1/(1 + ...))그리고 특별히 좋은 합리적 근사치는 없으며, π = 1 + 1/π를 곱하여 황금비 공식을 얻을 수 있습니다.

  2. NumPy에 -- 있기 때문에 NumPy는 "벡터화"됩니다.sqrt(array)는 다른 가 쓸 수 과 같습니다.map(sqrt, array)로 그서이은구성소별로요것래▁so▁a▁is로별▁component-.sqrt 를 사용하여 입니다. 이는 모든 됩니다.스칼라로 나눗셈하거나 스칼라로 덧셈하는 경우에도 마찬가지입니다. 모든 구성 요소에 병렬로 적용됩니다.

  3. 이것이 결과라는 것을 알면 그 증거는 간단합니다.만약 당신z < z < z + dz의 확률을 묻는다면, 이것-1 z < F(U) < z + dz, 단조롭게 증가하는 함수임을 주목하면서 F를 세 개의 에 모두 적용하여 F(z) < U(z + dz)를 찾는 과 같습니다.그리고 U가 균일하기 때문에 이 확률은 약속한 대로 f(z) dz입니다.

이를 구체의 포장점이라고 하며, 일반적인 (알려진) 완벽한 솔루션은 없습니다.하지만 불완전한 해결책들이 많이 있습니다.가장 인기 있는 세 가지는 다음과 같습니다.

  1. 시뮬레이션을 만듭니다.각 점을 구에 구속된 전자로 처리한 다음 특정 단계에 대해 시뮬레이션을 실행합니다.전자의 반발은 자연스럽게 점들이 서로로부터 가능한 한 멀리 떨어져 있는 더 안정적인 상태로 시스템을 기울이게 됩니다.
  2. 하이퍼큐브 거부.이 환상적인 방법은 실제로 매우 간단합니다. 구체를 둘러싸고 있는 입방체 내부의 점(그보다 훨씬 더 많은 점)을 균일하게 선택한 다음, 구체 외부의 점을 거부합니다.나머지 점을 벡터로 처리하고 정규화합니다.이것들은 당신의 "샘플"입니다 - 선택하세요n어떤 방법(욕심, 탐욕 등)을 사용하는 것.
  3. 나선형 근사.구면을 기준으로 완화곡선을 추적하고 완화곡선 주위에 점을 고르게 분포시킵니다.관련된 수학 때문에, 이것들은 시뮬레이션보다 이해하기가 더 복잡하지만 훨씬 더 빠릅니다(그리고 아마도 코드를 덜 포함할 것입니다).가장 인기 있는 것은 샤프 등인 것 같습니다.

이 문제에 대한 자세한 내용은 여기에서 확인할 수 있습니다.

이 예제 코드에서는 node[k]k번째 노드일 뿐입니다. N 와 "N"을 하고 있습니다.node[k]k번째(0 ~ N-1)입니다.그것이 당신을 혼란스럽게 하는 전부라면, 지금 그것을 사용할 수 있기를 바랍니다.

말로 하면, (다시말해서,서말해다▁(▁()k는 코드 조각이 시작되기 전에 정의되고 점 목록을 포함하는 크기 N의 배열입니다.

또는 여기에 있는 다른 답변(및 Python 사용)을 기반으로 합니다.

> cat ll.py
from math import asin
nx = 4; ny = 5
for x in range(nx):
    lon = 360 * ((x+0.5) / nx)
    for y in range(ny):                                                         
        midpt = (y+0.5) / ny                                                    
        lat = 180 * asin(2*((y+0.5)/ny-0.5))                                    
        print lon,lat                                                           
> python2.7 ll.py                                                      
45.0 -166.91313924                                                              
45.0 -74.0730322921                                                             
45.0 0.0                                                                        
45.0 74.0730322921                                                              
45.0 166.91313924                                                               
135.0 -166.91313924                                                             
135.0 -74.0730322921                                                            
135.0 0.0                                                                       
135.0 74.0730322921                                                             
135.0 166.91313924                                                              
225.0 -166.91313924                                                             
225.0 -74.0730322921                                                            
225.0 0.0                                                                       
225.0 74.0730322921                                                             
225.0 166.91313924
315.0 -166.91313924
315.0 -74.0730322921
315.0 0.0
315.0 74.0730322921
315.0 166.91313924

이 그래프를 표시하면 각 점이 거의 동일한 총 공간 영역에 위치하도록 극 근처의 수직 간격이 더 커집니다(극 근처에는 "수평" 공간이 더 적기 때문에 "수직" 공간이 더 많이 제공됨).

이것은 모든 점이 이웃과의 거리가 거의 같은 것과 같지는 않지만(이것이 링크가 말하는 것이라고 생각합니다), 원하는 것에 충분할 수 있으며 단순히 균일한 lat/lon 그리드를 만드는 것을 개선할 수 있습니다.

당신이 찾고 있는 것을 구형 덮개라고 합니다.구형 피복 문제는 매우 어렵고 소수의 점을 제외하고는 해결책을 알 수 없습니다.한 가지 확실하게 알려진 것은 구에 n개의 점이 주어지면 항상 두 개의 거리가 존재한다는 것입니다.d = (4-csc^2(\pi n/6(n-2)))^(1/2)아니면 더 가까이.

구에 균일하게 분포된 점을 생성하기 위한 확률론적 방법을 원한다면 가우스 분포에 의해 공간에서 점을 균일하게 생성하는 것이 쉽습니다(다른 언어의 코드를 찾는 것은 어렵지 않습니다).그래서 3차원 공간에서, 당신은 다음과 같은 것이 필요합니다.

Random r = new Random();
double[] p = { r.nextGaussian(), r.nextGaussian(), r.nextGaussian() };

그런 다음 원점으로부터의 거리를 정규화하여 구에 점을 투영합니다.

double norm = Math.sqrt( (p[0])^2 + (p[1])^2 + (p[2])^2 ); 
double[] sphereRandomPoint = { p[0]/norm, p[1]/norm, p[2]/norm };

n차원의 가우스 분포는 구면 대칭이므로 구면으로의 투영은 균일합니다.

물론 균일하게 생성된 점 집합에서 두 점 사이의 거리가 아래로 제한된다는 보장은 없으므로 거부를 사용하여 사용자가 가질 수 있는 모든 조건을 적용할 수 있습니다. 아마도 전체 집합을 생성한 다음 필요한 경우 전체 집합을 거부하는 것이 가장 좋습니다.(또는 "초기 거부"를 사용하여 지금까지 생성한 전체 컬렉션을 거부할 수 있습니다. 일부 포인트만 유지하고 다른 포인트는 삭제하지 마십시오.다음 공식을 사용할 수 있습니다.d위에서 주어진 값에서 약간의 슬랙을 뺀 값으로 포인트 집합을 기각할 포인트 사이의 최소 거리를 결정합니다.두 가지 거리를 계산하고 두 가지 거리를 선택해야 하며, 기각 확률은 느슨함에 따라 달라질 것입니다. 방법을 설명하기가 어려우므로 시뮬레이션을 실행하여 관련 통계를 파악하십시오.

이 답변은 이 답변에서 잘 설명된 동일한 '이론'을 기반으로 합니다.

이 답변을 다음과 같이 추가합니다.
'균일성'에 맞는 다른 옵션은 '스폿 온'(또는 분명히 그렇지 않음)을 필요로 합니다.(원래 요청에서 특별히 원하는 분포와 같은 동작을 행성으로 가져오기 위해서는 k개의 균일하게 생성된 점의 유한 목록에서 임의로 거부합니다(K개의 항목에 인덱스 카운트를 다시 씁니다).
는 '을으로두값 모두에서 은 'N'만하도록 했습니다. (에서는 무엇이가 될가 되지 않을 . 를 들어, '점 - fun입니다. (N의 낮은 카운트에서는 무엇이 중요하거나 중요하지 않은지 알기가 매우 까다롭습니다(예: '5'점을 원함 -- 재미있게 보내세요).
'groking'. 이 제공되는 에 있습니다. --complete 하그이없구방것어매게다니있아나구가제실즉함래습에이이와한현시능행니는되공다옵께과션습다가는렵법미는을지은우킹로이하분션을다른옵▁----▁it,fur▁without가아구있나실ation하즉에래니▁that시▁any행제한능▁looks▁options습gro와이현

N이 20일 때:

enter image description here
80에서 N: enter image description here


여기 실행 준비가 된 python3 코드가 있습니다. 여기서 에뮬레이션은 다른 사람들이 찾은 동일한 소스입니다. "http://web.archive.org/web/20120421191837/http ://www.cgafaq.info/wiki/Evenly_distributed_points_on_sphere ". ('main'으로 실행될 때 실행되는 플롯은 http://www.scipy.org/Cookbook/Matplotlib/mplot3D 에서 가져온 것입니다.)

from math import cos, sin, pi, sqrt

def GetPointsEquiAngularlyDistancedOnSphere(numberOfPoints=45):
    """ each point you get will be of form 'x, y, z'; in cartesian coordinates
        eg. the 'l2 distance' from the origion [0., 0., 0.] for each point will be 1.0 
        ------------
        converted from:  http://web.archive.org/web/20120421191837/http://www.cgafaq.info/wiki/Evenly_distributed_points_on_sphere ) 
    """
    dlong = pi*(3.0-sqrt(5.0))  # ~2.39996323 
    dz   =  2.0/numberOfPoints
    long =  0.0
    z    =  1.0 - dz/2.0
    ptsOnSphere =[]
    for k in range( 0, numberOfPoints): 
        r    = sqrt(1.0-z*z)
        ptNew = (cos(long)*r, sin(long)*r, z)
        ptsOnSphere.append( ptNew )
        z    = z - dz
        long = long + dlong
    return ptsOnSphere

if __name__ == '__main__':                
    ptsOnSphere = GetPointsEquiAngularlyDistancedOnSphere( 80)    

    #toggle True/False to print them
    if( True ):    
        for pt in ptsOnSphere:  print( pt)

    #toggle True/False to plot them
    if(True):
        from numpy import *
        import pylab as p
        import mpl_toolkits.mplot3d.axes3d as p3

        fig=p.figure()
        ax = p3.Axes3D(fig)

        x_s=[];y_s=[]; z_s=[]

        for pt in ptsOnSphere:
            x_s.append( pt[0]); y_s.append( pt[1]); z_s.append( pt[2])

        ax.scatter3D( array( x_s), array( y_s), array( z_s) )                
        ax.set_xlabel('X'); ax.set_ylabel('Y'); ax.set_zlabel('Z')
        p.show()
        #end

낮은 카운트(2, 5, 7, 13 등의 N)에서 테스트되며 '좋음'으로 작동하는 것 같습니다.

시도:

function sphere ( N:float,k:int):Vector3 {
    var inc =  Mathf.PI  * (3 - Mathf.Sqrt(5));
    var off = 2 / N;
    var y = k * off - 1 + (off / 2);
    var r = Mathf.Sqrt(1 - y*y);
    var phi = k * inc;
    return Vector3((Mathf.Cos(phi)*r), y, Mathf.Sin(phi)*r); 
};

위의 함수는 N개의 루프 합계와 k개의 루프 전류 반복으로 실행되어야 합니다.

해바라기 씨 패턴에 기반을 두고 있습니다. 해바라기 씨는 반구형으로 구부러져 있고 다시 구형으로 되어 있습니다.

여기 사진이 있습니다. 제가 카메라를 공의 중간쯤에 둔 것을 제외하고는요.

enter image description here

Healpix는 밀접하게 관련된 문제(구면을 동일한 면적 픽셀로 픽셀화)를 해결합니다.

http://healpix.sourceforge.net/

아마 과잉 살상일 수도 있지만, 그것을 보고 난 후에 여러분은 다른 멋진 특성들이 여러분에게 흥미롭다는 것을 알게 될 것입니다.이것은 단순히 포인트 클라우드를 출력하는 기능 이상의 것입니다.

다시 찾으려다가 착륙했어요 "힐픽스"라는 이름이 구체를 떠올리게 하는 건 아니지만...

편집: OP가 의도한 질문에 답하지 않고, 사람들이 어떻게든 유용하다고 느낄 경우를 대비해 여기에 남깁니다.

우리는 확률의 곱셈 법칙을 무한 시멘탈과 결합하여 사용합니다.이렇게 하면 원하는 결과를 얻을 수 있는 두 줄의 코드가 생성됩니다.

longitude: φ = uniform([0,2pi))
azimuth:   θ = -arcsin(1 - 2*uniform([0,1]))

(다음 좌표계에서 정의)

enter image description here

사용자의 언어는 일반적으로 균일한 난수 기본값을 사용합니다.를 들어 에서는 예에들어서니다사용수있습할를을 사용할 수 .random.random()범위 내의 번호를 반환합니다[0,1) 이 숫 를 에 자 여 내 임 숫 의 수 있 니 습 다 얻 를 을 자 의 의 범 위 곱 하 ▁range ▁in 니 ▁you ▁number 다 있 ▁a 습 ▁random ▁this이 ▁the ▁to ▁get ▁by ▁multiply ▁can ▁number 수 자 숫[0,k)따라서 파이썬에서는uniform([0,2pi))을 의미합니다.random.random()*2*math.pi.


증명

지금은 θ을 일률적으로 할당할 수 없습니다. 그렇지 않으면 극지방에 뭉치게 될 것입니다.구면 쐐기의 표면적에 비례하는 확률을 할당하고자 합니다(이 다이어그램의 θ는 실제로 θ입니다).

enter image description here

적도에서 각 변위 dθ는 dθ*r의 변위를 초래합니다.임의의 방위각 µ에서 변위는 얼마입니까?은 음, z축로부의반지은름터으▁well은반▁from.r*sin(θ)은 "접서쐐교와의차는하사 "미그기같다 "래다니과습 "호는음길이 "의"▁of"접사 "▁intersectcl▁so같다 "입니다.dφ * r*sin(θ)따라서 우리는 남극에서 북극까지의 슬라이스 면적을 통합하여 샘플링할 면적의 누적 분포를 계산합니다.

enter image description here stuff(여기서)dφ*r)

이제 CDF의 역방향을 샘플로 얻으려고 합니다. http://en.wikipedia.org/wiki/Inverse_transform_sampling

먼저 거의 CDF를 최대값으로 나누어 정규화합니다.이것은 dφ와 r을 상쇄시키는 부작용이 있습니다.

azimuthalCDF: cumProb = (sin(θ)+1)/2 from -pi/2 to pi/2

inverseCDF: θ = -sin^(-1)(1 - 2*cumProb)

따라서:

let x by a random float in range [0,1]
θ = -arcsin(1-2*x)

소수의 점으로 시뮬레이션을 실행할 수 있습니다.

from random import random,randint
r = 10
n = 20
best_closest_d = 0
best_points = []
points = [(r,0,0) for i in range(n)]
for simulation in range(10000):
    x = random()*r
    y = random()*r
    z = r-(x**2+y**2)**0.5
    if randint(0,1):
        x = -x
    if randint(0,1):
        y = -y
    if randint(0,1):
        z = -z
    closest_dist = (2*r)**2
    closest_index = None
    for i in range(n):
        for j in range(n):
            if i==j:
                continue
            p1,p2 = points[i],points[j]
            x1,y1,z1 = p1
            x2,y2,z2 = p2
            d = (x1-x2)**2+(y1-y2)**2+(z1-z2)**2
            if d < closest_dist:
                closest_dist = d
                closest_index = i
    if simulation % 100 == 0:
        print simulation,closest_dist
    if closest_dist > best_closest_d:
        best_closest_d = closest_dist
        best_points = points[:]
    points[closest_index]=(x,y,z)


print best_points
>>> best_points
[(9.921692138442777, -9.930808529773849, 4.037839326088124),
 (5.141893371460546, 1.7274947332807744, -4.575674650522637),
 (-4.917695758662436, -1.090127967097737, -4.9629263893193745),
 (3.6164803265540666, 7.004158551438312, -2.1172868271109184),
 (-9.550655088997003, -9.580386054762917, 3.5277052594769422),
 (-0.062238110294250415, 6.803105171979587, 3.1966101417463655),
 (-9.600996012203195, 9.488067284474834, -3.498242301168819),
 (-8.601522086624803, 4.519484132245867, -0.2834204048792728),
 (-1.1198210500791472, -2.2916581379035694, 7.44937337008726),
 (7.981831370440529, 8.539378431788634, 1.6889099589074377),
 (0.513546008372332, -2.974333486904779, -6.981657873262494),
 (-4.13615438946178, -6.707488383678717, 2.1197605651446807),
 (2.2859494919024326, -8.14336582650039, 1.5418694699275672),
 (-7.241410895247996, 9.907335206038226, 2.271647103735541),
 (-9.433349952523232, -7.999106443463781, -2.3682575660694347),
 (3.704772125650199, 1.0526567864085812, 6.148581714099761),
 (-3.5710511242327048, 5.512552040316693, -3.4318468250897647),
 (-7.483466337225052, -1.506434920354559, 2.36641535124918),
 (7.73363824231576, -8.460241422163824, -1.4623228616326003),
 (10, 0, 0)]

fnord의 답변에 따르면 범위가 추가된 Unity3D 버전이 있습니다.

코드:

// golden angle in radians
static float Phi = Mathf.PI * ( 3f - Mathf.Sqrt( 5f ) );
static float Pi2 = Mathf.PI * 2;

public static Vector3 Point( float radius , int index , int total , float min = 0f, float max = 1f , float angleStartDeg = 0f, float angleRangeDeg = 360 )
{
    // y goes from min (-) to max (+)
    var y = ( ( index / ( total - 1f ) ) * ( max - min ) + min ) * 2f - 1f;

    // golden angle increment
    var theta = Phi * index ; 
        
    if( angleStartDeg != 0 || angleRangeDeg != 360 )
    {
        theta = ( theta % ( Pi2 ) ) ;
        theta = theta < 0 ? theta + Pi2 : theta ;
            
        var a1 = angleStartDeg * Mathf.Deg2Rad;
        var a2 = angleRangeDeg * Mathf.Deg2Rad;
            
        theta = theta * a2 / Pi2 + a1;
    }

    // https://stackoverflow.com/a/26127012/2496170
    
    // radius at y
    var rY = Mathf.Sqrt( 1 - y * y ); 
    
    var x = Mathf.Cos( theta ) * rY;
    var z = Mathf.Sin( theta ) * rY;

    return  new Vector3( x, y, z ) * radius;
}

요지: https://gist.github.com/nukadelic/7449f0872f708065bc1afeb19df666f7/edit

미리 보기:

MP4

의 사의가큰장요고보오십시려해인을지의 가장 큰 두 을 생각해 .N,한다면N==20 큰 두 이 다면가큰두요가은인지그렇장▁are▁then은▁largest요인▁the▁factors가▁two입니다.{5,4}더 일반적으로 또는일, 다적으로반보.{a,b}합니다.계산한다.

dlat  = 180/(a+1)
dlong = 360/(b+1})

첫 번째 요점을 다음에 제시합니다.{90-dlat/2,(dlong/2)-180}에 대한 당신의 두 번째.{90-dlat/2,(3*dlong/2)-180}당신의 3번째 에서의{90-dlat/2,(5*dlong/2)-180}당신이 세계를 한 바퀴 돌 때까지, 그 때쯤이면 당신은 대략.{75,150} 옆에 {90-3*dlat/2,(dlong/2)-180}.

분명히 저는 이것을 구면의 표면에서 +/-를 N/S 또는 E/W로 변환하는 일반적인 규칙으로 작업하고 있습니다.그리고 분명히 이것은 여러분에게 완전한 비랜덤 분포를 제공하지만, 그것은 균일하고 점들은 함께 뭉치지 않습니다.

랜덤성을 어느 정도 추가하려면 정규 분포를 따르는 2개(필요에 따라 평균 0 및 표준 편차 {dlat/3, dlong/3})를 생성하여 균일하게 분포된 점에 추가할 수 있습니다.

또는... 20개의 점을 배치하려면 정이십면체 면의 중심을 계산합니다.12개의 점에 대해 정이십면체의 꼭짓점을 찾습니다.30개의 점의 경우, 정이십면체 모서리의 중간점입니다.당신은 사면체, 입방체, 십이면체, 팔면체에서 같은 일을 할 수 있습니다: 한 세트의 점들은 꼭짓점에 있고, 다른 점들은 면의 중심에 있고, 또 다른 점들은 모서리의 중심에 있습니다.하지만 그것들은 섞일 수 없습니다.

# create uniform spiral grid
numOfPoints = varargin[0]
vxyz = zeros((numOfPoints,3),dtype=float)
sq0 = 0.00033333333**2
sq2 = 0.9999998**2
sumsq = 2*sq0 + sq2
vxyz[numOfPoints -1] = array([(sqrt(sq0/sumsq)), 
                              (sqrt(sq0/sumsq)), 
                              (-sqrt(sq2/sumsq))])
vxyz[0] = -vxyz[numOfPoints -1] 
phi2 = sqrt(5)*0.5 + 2.5
rootCnt = sqrt(numOfPoints)
prevLongitude = 0
for index in arange(1, (numOfPoints -1), 1, dtype=float):
  zInc = (2*index)/(numOfPoints) -1
  radius = sqrt(1-zInc**2)

  longitude = phi2/(rootCnt*radius)
  longitude = longitude + prevLongitude
  while (longitude > 2*pi): 
    longitude = longitude - 2*pi

  prevLongitude = longitude
  if (longitude > pi):
    longitude = longitude - 2*pi

  latitude = arccos(zInc) - pi/2
  vxyz[index] = array([ (cos(latitude) * cos(longitude)) ,
                        (cos(latitude) * sin(longitude)), 
                        sin(latitude)])

@robert king 정말 좋은 해결책이지만 몇 가지 엉성한 버그가 있습니다.그래도 도움이 많이 됐으니, 나른한 건 신경 쓰지 마세요.:) 정리된 버전이 있습니다.

from math import pi, asin, sin, degrees
halfpi, twopi = .5 * pi, 2 * pi
sphere_area = lambda R=1.0: 4 * pi * R ** 2

lat_dist = lambda lat, R=1.0: R*(1-sin(lat))

#A = 2*pi*R^2(1-sin(lat))
def sphere_latarea(lat, R=1.0):
    if -halfpi > lat or lat > halfpi:
        raise ValueError("lat must be between -halfpi and halfpi")
    return 2 * pi * R ** 2 * (1-sin(lat))

sphere_lonarea = lambda lon, R=1.0: \
        4 * pi * R ** 2 * lon / twopi

#A = 2*pi*R^2 |sin(lat1)-sin(lat2)| |lon1-lon2|/360
#    = (pi/180)R^2 |sin(lat1)-sin(lat2)| |lon1-lon2|
sphere_rectarea = lambda lat0, lat1, lon0, lon1, R=1.0: \
        (sphere_latarea(lat0, R)-sphere_latarea(lat1, R)) * (lon1-lon0) / twopi


def test_sphere(n_lats=10, n_lons=19, radius=540.0):
    total_area = 0.0
    for i_lons in range(n_lons):
        lon0 = twopi * float(i_lons) / n_lons
        lon1 = twopi * float(i_lons+1) / n_lons
        for i_lats in range(n_lats):
            lat0 = asin(2 * float(i_lats) / n_lats - 1)
            lat1 = asin(2 * float(i_lats+1)/n_lats - 1)
            area = sphere_rectarea(lat0, lat1, lon0, lon1, radius)
            print("{:} {:}: {:9.4f} to  {:9.4f}, {:9.4f} to  {:9.4f} => area {:10.4f}"
                    .format(i_lats, i_lons
                    , degrees(lat0), degrees(lat1)
                    , degrees(lon0), degrees(lon1)
                    , area))
            total_area += area
    print("total_area = {:10.4f} (difference of {:10.4f})"
            .format(total_area, abs(total_area) - sphere_area(radius)))

test_sphere()

Jonathan Kogan은 이 문제를 반정확하게 해결하는 동시에 속도에 영향을 주지 않는 새로운 방법(2016, pdf)을 제안했습니다.

import math
def spherical_coordinate(x, y):
    return [
        math.cos(x) * math.cos(y),
        math.sin(x) * math.cos(y),
        math.sin(y)
    ]

def NX(n, x):
    pts = []
    start = (-1. + 1. / (n - 1.))
    increment = (2. - 2. / (n - 1.)) / (n - 1.)
    for j in xrange(0, n):
        s = start + j * increment
        pts.append(
            spherical_cordinate(
                s * x, math.pi / 2. *
                math.copysign(1, s) *
                (1. - math.sqrt(1. - abs(s)))
            )
        )
    return pts

def generate_points(n):
    return NX(n, 0.1 + 1.2 * n)

이것은 효과가 있고 매우 간단합니다.원하는 만큼의 포인트:

    private function moveTweets():void {


        var newScale:Number=Scale(meshes.length,50,500,6,2);
        trace("new scale:"+newScale);


        var l:Number=this.meshes.length;
        var tweetMeshInstance:TweetMesh;
        var destx:Number;
        var desty:Number;
        var destz:Number;
        for (var i:Number=0;i<this.meshes.length;i++){

            tweetMeshInstance=meshes[i];

            var phi:Number = Math.acos( -1 + ( 2 * i ) / l );
            var theta:Number = Math.sqrt( l * Math.PI ) * phi;

            tweetMeshInstance.origX = (sphereRadius+5) * Math.cos( theta ) * Math.sin( phi );
            tweetMeshInstance.origY= (sphereRadius+5) * Math.sin( theta ) * Math.sin( phi );
            tweetMeshInstance.origZ = (sphereRadius+5) * Math.cos( phi );

            destx=sphereRadius * Math.cos( theta ) * Math.sin( phi );
            desty=sphereRadius * Math.sin( theta ) * Math.sin( phi );
            destz=sphereRadius * Math.cos( phi );

            tweetMeshInstance.lookAt(new Vector3D());


            TweenMax.to(tweetMeshInstance, 1, {scaleX:newScale,scaleY:newScale,x:destx,y:desty,z:destz,onUpdate:onLookAtTween, onUpdateParams:[tweetMeshInstance]});

        }

    }
    private function onLookAtTween(theMesh:TweetMesh):void {
        theMesh.lookAt(new Vector3D());
    }

언급URL : https://stackoverflow.com/questions/9600801/evenly-distributing-n-points-on-a-sphere

반응형