случайным образом генерировать точки в определенных сферических координатах в R - PullRequest
2 голосов
/ 07 апреля 2019

У меня есть конкретные координаты x, y, z.Я хочу сгенерировать случайные точки в сфере, в которой x указан как центр, а x2 - из другого фрейма данных в качестве границы радиуса (поэтому расстояние от x до x2 будет равно длине радиуса сферы).

Я видел много дискуссий о том, как сделать это соответствующим образом математически (распределить точки случайным образом, чтобы избежать кластеризации), и смог собрать самые простые примеры здесь и здесь для примера кода R.

Я также нашел этот [сферный пакет R] (https://cran.r -project.org / web / packages /phereplot /phereplot.pdf ), который можетбудет проще, но мне трудно понять, как его применить.

Все это хорошие начальные точки, но, используя приведенный ниже пример кода, я не уверен, как применить его к определенным начальным точкам / сферическим координатам?

set.seed(101)
n <- 50
theta <- runif(n,0,2*pi)
u <- runif(n,-1,1)
x <- sqrt(1-u^2)*cos(theta)
y <- sqrt(1-u^2)*sin(theta)
z <- u

Используя только один набор / строку координат x, y, z из моего фрейма данных:

x = -0.0684486861
y= 0.0125857380
z= 0.0201056441

x2= -0.0684486861
y2 = 0.0125857380
z2= -0.0228805516


Я хочу, чтобы x, y, z был центром сферы ирасстояние до x2, y2, z2 будет радиусом длины / края сферы.Затем создайте случайные точки внутри сферы с x, y, z в качестве центра.

В конце концов, я пытаюсь сделать это с 100 сферами, чтобы сравнить, все ли точки во втором наборе координат движутся в одинаковых углах / направлениях в пространстве.

Спасибо за руководство.

1 Ответ

1 голос
/ 07 апреля 2019

Хорошо, давайте разделим задачу на несколько подзадач.

Во-первых, генерирует точки, равномерно распределенные по сфере (или по объему, или по поверхности), с центром в (0,0,0) и заданным радиусом. После http://mathworld.wolfram.com/SpherePointPicking.html, и совсем близко к показанному вами коду,

rsphere <- function(n, r = 1.0, surface_only = FALSE) {
    phi       <- runif(n, 0.0, 2.0 * pi)
    cos_theta <- runif(n, -1.0, 1.0)
    sin_theta <- sqrt((1.0-cos_theta)*(1.0+cos_theta))
    radius <- r
    if (surface_only == FALSE) {
        radius <- r * runif(n, 0.0, 1.0)^(1.0/3.0)
    }

    x <- radius * sin_theta * cos(phi)
    y <- radius * sin_theta * sin(phi)
    z <- radius * cos_theta

    cbind(x, y, z)
}

set.seed(312345)
sphere_points <- rsphere(10000)

Вторая проблема - переместить эти точки в центр в точке X

rsphere <- function(n, r = 1.0, surface_only = FALSE, center=cbind(Xx, Xy, Xz)) {
    ....
    cbind(x+center[1], y+center[2], z+center[3])
}

Третья проблема - вычислить радиус заданного центра в (Xx, Xy, Xz) и точке поверхности (Yx, Yy, Yz))

radius <- sqrt((Xx-Yx)**2+(Xy-Yy)**2+(Xz-Yz)**2)

Объедините их всех вместе для полного удовлетворения. Хорошо, теперь, когда вы указали значения для центра и радиуса, давайте соберем все вместе

rsphere <- function(n, r = 1.0, surface_only = FALSE, center=cbind(0.0, 0.0, 0.0)) {
    phi       <- runif(n, 0.0, 2.0 * pi)
    cos_theta <- runif(n, -1.0, 1.0)
    sin_theta <- sqrt((1.0-cos_theta)*(1.0+cos_theta))
    radius <- r
    if (surface_only == FALSE) {
        radius <- r * runif(n, 0.0, 1.0)^(1.0/3.0)
    }

    x <- radius * sin_theta * cos(phi)
    y <- radius * sin_theta * sin(phi)
    z <- radius * cos_theta

    # if radius is fixed, we could check it
    # rr = sqrt(x^2+y^2+z^2)
    # print(rr)

    cbind(x+center[1], y+center[2], z+center[3])
}

x1 = -0.0684486861
y1 = 0.0125857380
z1 = 0.0201056441

x2 = -0.0684486861
y2 = 0.0125857380
z2 = -0.0228805516

R = sqrt((x2-x1)^2 + (y2-y1)^2 + (z2-z1)^2)
print(R)

set.seed(32345)
sphere_points <- rsphere(100000, R, FALSE, cbind(x1, y1, z1))

Как это выглядит?

UPDATE

Сгенерировал 10 точек на поверхности и в объеме и напечатал их, радиус = 2 выглядит нормально для меня

# 10 points uniform on surface, supposed to have fixed radius
sphere_points <- rsphere(10, 2, TRUE, cbind(x1, y1, z1))
for (k in 1:10) {
    rr <- sqrt((sphere_points[k,1]-x1)^2+(sphere_points[k,2]-y1)^2+(sphere_points[k,3]-z1)^2)
    print(rr)
}

# 10 points uniform in the sphere, supposed to have varying radius
sphere_points <- rsphere(10, 2, FALSE, cbind(x1, y1, z1))
for (k in 1:10) {
    rr <- sqrt((sphere_points[k,1]-x1)^2+(sphere_points[k,2]-y1)^2+(sphere_points[k,3]-z1)^2)
    print(rr)
}

получил

[1] 2
[1] 2
[1] 2
[1] 2
[1] 2
[1] 2
[1] 2
[1] 2
[1] 2
[1] 2

и

[1] 1.32571
[1] 1.505066
[1] 1.255023
[1] 1.82773
[1] 1.219957
[1] 1.641258
[1] 1.881937
[1] 1.083975
[1] 0.4745712
[1] 1.900066

что вы получили?

...