8

私は一対の点を持っており、これらの2つの点によって決定される既知のrの円を見つけたいと思います。これをシミュレーションで使用し、境界がある可能性のあるスペース(たとえばxy-200、200のボックス)を使用します。

半径の二乗

(x-x1)^2 + (y-y1)^2 = r^2
(x-x2)^2 + (y-y2)^2 = r^2

ここで、この非線形連立方程式を解いて、2つの潜在的な円の中心を取得したいと思います。パッケージを使ってみBBました。これが私の微妙な試みであり、1つのポイントしか与えられません。私が得たいのは、両方の可能な点です。正しい方向へのポインターは、最初の可能な機会に無料のビールと出会うでしょう。

library(BB)
known.pair <- structure(c(-46.9531139599816, -62.1874917150412, 25.9011462171242, 
16.7441676243879), .Dim = c(2L, 2L), .Dimnames = list(NULL, c("x", 
"y")))

getPoints <- function(ps, r, tr) {
    # get parameters
     x <- ps[1]
     y <- ps[2]

     # known coordinates of two points
     x1 <- tr[1, 1]
     y1 <- tr[1, 2]
     x2 <- tr[2, 1]
     y2 <- tr[2, 2]

     out <- rep(NA, 2)
     out[1] <- (x-x1)^2 + (y-y1)^2 - r^2
     out[2] <- (x-x2)^2 + (y-y2)^2 - r^2
     out
}

slvd <- BBsolve(par = c(0, 0),
                fn = getPoints,
                method = "L-BFGS-B",
                tr = known.pair,
                r = 40
                )

グラフィカルにこれを次のコードで確認できますが、いくつかの追加パッケージが必要になります。

library(sp)
library(rgeos)
plot(0,0, xlim = c(-200, 200), ylim = c(-200, 200), type = "n", asp = 1)
points(known.pair)
found.pt <- SpatialPoints(matrix(slvd$par, nrow = 1))
plot(gBuffer(found.pt, width = 40), add = T)

ここに画像の説明を入力してください

補遺

貴重なコメントとコードをありがとうございました。回答をコードで補完した投稿者による回答のタイミングを提供します。

    test replications elapsed relative user.self sys.self user.child sys.child
4   alex          100    0.00       NA      0.00        0         NA        NA
2  dason          100    0.01       NA      0.02        0         NA        NA
3   josh          100    0.01       NA      0.02        0         NA        NA
1 roland          100    0.15       NA      0.14        0         NA        NA
4

7 に答える 7

9

次のコードは、目的の 2 つの円の中心にあるポイントを取得します。これをコメントアップしたり、結果をオブジェクトに変換したりする時間は今のところありませんSpatial*が、これで良いスタートが切れるはずです。

まず、ポイント名を紹介する ASCII アートの図を示します。kKは既知の点、Bは を通る水平線上の点kC1C2は求める円の中心です。

        C2





                            K


                  k----------------------B






                                       C1

今コード:

# Example inputs
r <- 40
known.pair <- structure(c(-46.9531139599816, -62.1874917150412, 
25.9011462171242, 16.7441676243879), .Dim = c(2L, 2L), 
.Dimnames = list(NULL, c("x", "y")))

## Distance and angle (/_KkB) between the two known points
d1 <- sqrt(sum(diff(known.pair)^2))
theta1 <- atan(do.call("/", as.list(rev(diff(known.pair)))))

## Calculate magnitude of /_KkC1 and /_KkC2
theta2 <- acos((d1/2)/r)

## Find center of one circle (using /_BkC1)
dx1 <- cos(theta1 + theta2)*r
dy1 <- sin(theta1 + theta2)*r
p1 <- known.pair[2,] + c(dx1, dy1)

## Find center of other circle (using /_BkC2)
dx2 <- cos(theta1 - theta2)*r
dy2 <- sin(theta1 - theta2)*r
p2 <- known.pair[2,] + c(dx2, dy2)

## Showing that it worked
library(sp)
library(rgeos)
plot(0,0, xlim = c(-200, 200), ylim = c(-200, 200), type = "n", asp = 1)
points(known.pair)
found.pt <- SpatialPoints(matrix(slvd$par, nrow = 1))
points(p1[1], p1[2], col="blue", pch=16)
points(p2[1], p2[2], col="green", pch=16)

ここに画像の説明を入力

于 2012-09-04T15:04:28.763 に答える
4

これは、他の誰もが言及している、それを解決するための基本的な幾何学的方法です。得られた二次方程式の根を取得するために polyroot を使用しますが、二次方程式を直接使用することも簡単にできます。

# x is a vector containing the two x coordinates
# y is a vector containing the two y coordinates
# R is a scalar for the desired radius
findCenter <- function(x, y, R){
    dy <- diff(y)
    dx <- diff(x)
    # The radius needs to be at least as large as half the distance
    # between the two points of interest
    minrad <- (1/2)*sqrt(dx^2 + dy^2)
    if(R < minrad){
        stop("Specified radius can't be achieved with this data")
    }

    # I used a parametric equation to create the line going through
    # the mean of the two points that is perpendicular to the line
    # connecting the two points
    # 
    # f(k) = ((x1+x2)/2, (y1+y2)/2) + k*(y2-y1, x1-x2)
    # That is the vector equation for our line.  Then we can
    # for any given value of k calculate the radius of the circle
    # since we have the center and a value for a point on the
    # edge of the circle.  Squaring the radius, subtracting R^2,
    # and equating to 0 gives us the value of t to get a circle
    # with the desired radius.  The following are the coefficients
    # we get from doing that
    A <- (dy^2 + dx^2)
    B <- 0
    C <- (1/4)*(dx^2 + dy^2) - R^2

    # We could just solve the quadratic equation but eh... polyroot is good enough
    k <- as.numeric(polyroot(c(C, B, A)))

    # Now we just plug our solution in to get the centers
    # of the circles that meet our specifications
    mn <- c(mean(x), mean(y))
    ans <- rbind(mn + k[1]*c(dy, -dx),
                 mn + k[2]*c(dy, -dx))

    colnames(ans) = c("x", "y")

    ans
}

findCenter(c(-2, 0), c(1, 1), 3)
于 2012-09-04T15:53:33.797 に答える
4

@PhilHのソリューションに従って、Rで三角法を使用するだけです:

radius=40

半径上に元の点を描きます

plot(known.pair,xlim=100*c(-1,1),ylim=100*c(-1,1),asp=1,pch=c("a","b"),cex=0.8)

cの中点を見つけますab(これはde2 つの円の中心の中点でもあります)

AB.bisect=known.pair[2,,drop=F]/2+known.pair[1,,drop=F]/2
C=AB.bisect
points(AB.bisect,pch="c",cex=0.5)

弦の長さと角度を見つけるab

AB.vector=known.pair[2,,drop=F]-known.pair[1,,drop=F]
AB.len=sqrt(sum(AB.vector^2))
AB.angle=atan2(AB.vector[2],AB.vector[1])
names(AB.angle)<-NULL

cから2 つの円の中心までの線の長さと角度を計算します

CD.len=sqrt(diff(c(AB.len/2,radius)^2))
CD.angle=AB.angle-pi/2

d2 つの中心の位置とe垂線abと長さを計算してプロットします。

center1=C+CD.len*c(x=cos(CD.angle),y=sin(CD.angle))
center2=C-CD.len*c(x=cos(CD.angle),y=sin(CD.angle))
points(center1[1],center1[2],col="blue",cex=0.8,pch="d")
points(center2[1],center2[2],col="blue",cex=0.8,pch="e")

ショー:

ここに画像の説明を入力

于 2012-09-04T18:11:46.490 に答える
3

数値方程式を解く必要はありません。数式だけ:

  1. 点AとBの両方が円上にあるため、それぞれから特定の中心までの距離が半径rであることがわかります。
  2. 底辺にある2つの既知の点と、円の中心にある3番目の点の弦で二等辺三角形を形成します。
  3. AとBの中間で三角形を二等分し、直角三角形を作成します。
  4. http://mathworld.wolfram.com/IsoscelesTriangle.htmlは、ベースの長さと半径の観点から高さを示します。
  5. ポイントから各方向に計算された高さの距離については、ABコードの法線(このSO回答を参照)に従ってください。
于 2012-09-04T14:46:44.387 に答える
1

取り除くべき警告がいくつかありますが、これで作業を開始できます。パフォーマンスの問題が発生する可能性があるため、基本的なジオメトリで完全に解決する方が適切なアプローチになる可能性があります。

known.pair <- structure(c(-46.9531139599816, -62.1874917150412, 25.9011462171242, 
                          16.7441676243879), .Dim = c(2L, 2L), .Dimnames = list(NULL, c("x", 
                                                                                        "y")))

findCenter <- function(p,r) {
  yplus <- function(y) {
    ((p[1,1]+sqrt(r^2-(y-p[1,2])^2)-p[2,1])^2+(y-p[2,2])^2-r^2)^2
  }


 yp <- optimize(yplus,interval=c(min(p[,2]-r),max(p[,2]+r)))$minimum
 xp <- p[1,1]+sqrt(r^2-(yp-p[1,2])^2)  
 cp <- c(xp,yp)
 names(cp)<-c("x","y") 

 yminus <- function(y) {
   ((p[1,1]-sqrt(r^2-(y-p[1,2])^2)-p[2,1])^2+(y-p[2,2])^2-r^2)^2
 }


 ym <- optimize(yminus,interval=c(min(p[,2]-r),max(p[,2]+r)))$minimum
 xm <- p[1,1]-sqrt(r^2-(ym-p[1,2])^2)  
 cm <- c(xm,ym)
 names(cm)<-c("x","y")


 list(c1=cp,c2=cm)
}

cent <- findCenter(known.pair,40)
于 2012-09-04T14:22:38.380 に答える
1

これが答えの骨です。後で時間があれば、それらを肉付けします。これは、単語に沿って描画すれば簡単に理解できるはずです。申し訳ありませんが、このコンピューターには絵を描くための適切なソフトウェアがありません。

点が同一である (無限解) か、離れすぎて選択した半径の同じ円上にない (解がない) 退化したケースは脇に置いておきます。

点 とおよび 2 つの円の未知の中心点Xとにラベルを付けます。の垂直二等分線上にあり、この行を と呼びます。この段階では、との場所の詳細がすべてわかっていないことは重要ではありません。Yc1c2c1c2XYc1c2c1c2

というわけで、直線の方程式を考えてみましょうc1c2XYそれは(この点を と呼びます)の中間点を通過しZ、 の負の逆数に等しい勾配を持ちXYます。これで、次の方程式が得られましたc1c2(または、これらの骨に肉があればそうなるでしょう)。

ここで、ある点から直線とその垂直二等分線と円の中心点の交点までの三角形を作成します (たとえばXZc1)。どこにあるかはまだ正確にc1はわかりませんが、ジオメトリのスケッチを止めることはありませんでした。XZ既知の 2 つの辺の長さ (および)を持つ直角三角形があるXc1ので、簡単に見つけることができますZc1。他の三角形と円の中心に対してプロセスを繰り返します。

もちろん、このアプローチは OP の最初のアプローチとはかなり異なり、魅力的ではないかもしれません。

于 2012-09-04T14:18:22.547 に答える
0

残念ながら私はそれを描くことができないので、いくつかの基本的な幾何学を知っていることを願っています.

垂直二等分線は、A と B の両方を横切る円のすべての中点を結ぶ線です。

これで、AB と r の中心ができたので、点 A、AB の中心、未知の円の中点で直角三角形を描くことができます。

次に、ピタゴラスの定理を使用して、AB の中点から円の中点までの距離を取得し、基本的な sin/cos の組み合わせを使用して、ここから円の位置を計算することは難しくありません。

于 2012-09-04T14:31:21.960 に答える