Как найти пересечение кривой и окружности? - PullRequest
Новая
       26

Как найти пересечение кривой и окружности?

0 голосов
/ 29 апреля 2018

У меня есть кривая, полученная из эмпирических данных, и я могу получить ее разумную модель. Мне нужно определить точку (x, y), где кривая пересекает круг известного центра и радиуса. Следующий код иллюстрирует вопрос. 

x <- c(0.05, 0.20, 0.35, 0.50, 0.65, 0.80, 0.95, 
   1.10, 1.25, 1.40, 1.55, 1.70, 1.85, 2.00, 
   2.15, 2.30, 2.45, 2.60, 2.75, 2.90, 3.05)

y <- c(1.52, 1.44, 1.38, 1.31, 1.23, 1.15, 1.06,
   0.96, 0.86, 0.76, 0.68, 0.61, 0.54, 0.47, 
   0.41, 0.36, 0.32, 0.29, 0.27, 0.26, 0.26)

fit <- loess(y ~ x, control = loess.control(surface = "direct"))
newx <- data.frame(x = seq(0, 3, 0.01))
fitline <- predict(fit, newdata = newx)
est <- data.frame(newx, fitline)

plot(x, y, type = "o",lwd = 2)
lines(est, col = "blue", lwd = 2)

library(plotrix)
draw.circle(x = 3, y = 0, radius = 2, nv = 1000, lty = 1, lwd = 1)

enter image description here

Ответы [ 2 ]

0 голосов
/ 29 апреля 2018

Найти пересечение с помощью функций из пакета sf очень просто.

Рассчитать значения круга (на основе этого ответа и , как сделано @ Onyambu ) 

circ <- function(xc = 0, yc = 0, r = 1, n = 100){
  v <- seq(0, 2 * pi, len = n)
  cbind(x = xc + r * cos(v),
        y = yc + r * sin(v))
}

m <- circ(xc = 3, yc = 0, r = 2)

Преобразование прогнозируемых значений и значений окружности в "простые элементы" (LINESTRING) и определение их пересечения (a POINT): 

library(sf)
int <- st_intersection(st_linestring(as.matrix(est)),
                       st_linestring(m))
int
# POINT (1.2091 0.8886608)

Добавьте пересечение к вашему участку: 

plot(x, y, type = "o", lwd = 2)
lines(est, col = "blue", lwd = 2)
lines(m)
points(int[1], int[2], col = "red", pch = 19)

enter image description here

0 голосов
/ 29 апреля 2018

Чтобы получить точку пересечения, мы можем использовать функцию optim в r, чтобы сделать так: 

circle=function(x){
  if(4<(x-3)^2) return(NA)# Ensure it is limited within the radius
  sqrt(4-(x-3)^2)
}
fun=function(x)predict(fit,data.frame(x=x))  
g=function(x)(circle(x)-fun(x))# We need to set this to zero. Ie solve this
sol1=optimise(function(x)abs(g(x)),c(1,5))$min
 [1] 1.208466

Таким образом, две функции должны оценить одно и то же значение в x=1.208466 ..

Чтобы сделать его еще более точным, вы можете использовать функцию optim: 

sol2= optim(1,function(x)abs(g(x)),g,method="Brent",upper=5,lower=1)$par
 [1] 1.208473

Теперь вы можете оценить: 

circle(sol1)
[1] 0.889047
fun(sol1)
        1 
0.8890654 
circle(sol2)
[1] 0.889061
fun(sol2)
       1 
0.889061

Из вышесказанного вы можете сказать, что решение 2 очень близко ..

Построение этой точки на графике будет сложной задачей, поскольку функция draw.circle рисует круги пропорционально zxes. Таким образом, каждый раз меняется в зависимости от того, насколько велика область графика.

Если бы вы написали свою собственную функцию круга: 

circleplot=function(x,y,r){
  theta=seq(0,2*pi,length.out = 150)
  cbind(x+r*cos(theta),y+r*sin(theta))
}

Тогда вы можете сделать: 

plot(x, y, type = "o",lwd = 2)
lines(est, col = "blue", lwd = 2)
lines(circleplot(3,0,2))
abline(v=sol2,col=2) 
points(sol2,fun(sol2),col=2,pch=16)

enter image description here

...