Вывести центр окружности из части окружности

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

Алгоритм был бы великолепен (особенно если бы он был объяснен так, чтобы я мог его понять), но, поскольку это, вероятно, известная территория, указатели на ответы также были бы хороши. Если все остальное терпит неудачу, некоторые лучшие условия поиска будут оценены.

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

(Я свободно владею языками c/java/kotlin, поэтому решения в чем-то подобном были бы для меня самыми простыми, но я рад работать с любым языком или псевдокодом).

заранее спасибо


Спасибо Ренату за его решение - оно работает. Я немного обижен, что не понимаю алгоритма, но ссылка есть и пока можно жить с незнанием.

Для тех, кому нужен код, переведенный на kotlin:

class Circle(val x: Int, val y: Int, val radius: Int) {
    companion object {
        /**
         * Find the best fit circle for a collection of points.
         * With many thanks to Renat:
         * https://stackoverflow.com/questions/74493481/deduce-center-of-circle-from-portion-of-circumference
         */
        fun fit(points: List<Point>): Circle {
            val sx = points.sumOf { it.x.toDouble() }
            val sy = points.sumOf { it.y.toDouble() }
            val sx2 = points.sumOf { it.x.toDouble() * it.x }
            val sy2 = points.sumOf { it.y.toDouble() * it.y }
            val sxy = points.sumOf { it.x.toDouble() * it.y }
            val sx3 = points.sumOf { it.x.toDouble() * it.x * it.x }
            val sy3 = points.sumOf { it.y.toDouble() * it.y * it.y }
            val sx2y = points.sumOf { it.x.toDouble() * it.x * it.y }
            val sxy2 = points.sumOf { it.x.toDouble() * it.y * it.y }

            val n = points.size;
            val s11 = n * sxy - sx * sy
            val s20 = n * sx2 - sx * sx
            val s02 = n * sy2 - sy * sy
            val s30 = n * sx3 - sx2 * sx
            val s03 = n * sy3 - sy * sy2
            val s21 = n * sx2y - sx2 * sy
            val s12 = n * sxy2 - sx * sy2
            val a = ((s30 + s12)*s02 - (s03 + s21)*s11) / (2*( s20 * s02 - s11 * s11));
            val b = ((s03 + s21)*s20 - (s30 + s12)*s11) / (2*( s20 * s02 - s11 * s11));
            val c = (sx2 + sy2 - 2*a*sx - 2*b*sy) / n;
            val radius = Math.sqrt(c + a*a + b*b);
            return Circle(a.roundToInt(), b.roundToInt(), radius.roundToInt())
        }
    }
}
Стоит ли изучать PHP в 2026-2027 годах?
Стоит ли изучать PHP в 2026-2027 годах?
Привет всем, сегодня я хочу высказать свои соображения по поводу вопроса, который я уже много раз получал в своем сообществе: "Стоит ли изучать PHP в...
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
В JavaScript одним из самых запутанных понятий является поведение ключевого слова "this" в стрелочной и обычной функциях.
Приемы CSS-макетирования - floats и Flexbox
Приемы CSS-макетирования - floats и Flexbox
Здравствуйте, друзья-студенты! Готовы совершенствовать свои навыки веб-дизайна? Сегодня в нашем путешествии мы рассмотрим приемы CSS-верстки - в...
Тестирование функциональных ngrx-эффектов в Angular 16 с помощью Jest
В системе управления состояниями ngrx, совместимой с Angular 16, появились функциональные эффекты. Это здорово и делает код определенно легче для...
Концепция локализации и ее применение в приложениях React ⚡️
Концепция локализации и ее применение в приложениях React ⚡️
Локализация - это процесс адаптации приложения к различным языкам и культурным требованиям. Это позволяет пользователям получить опыт, соответствующий...
Пользовательский скаляр GraphQL
Пользовательский скаляр GraphQL
Листовые узлы системы типов GraphQL называются скалярами. Достигнув скалярного типа, невозможно спуститься дальше по иерархии типов. Скалярный тип...
2
0
74
1
Перейти к ответу Данный вопрос помечен как решенный

Ответы 1

Ответ принят как подходящий

Итак, это задача Подгонка круга (или эллипса). Которая разрешима за линейное время.

Для приведенного ниже кода я использовал ответ на этот вопрос из Math.stackexchange , который ссылается на коники регрессии, квадрики. Линейные и явные регрессии, круговая, сферическая, стр. 12–13, в которых используется метод наименьших квадратов для минимизации суммы.

(нажмите «Выполнить фрагмент кода» и рисуйте мышью)

var canvas = document.createElement('canvas');
document.body.appendChild(canvas);

document.body.style.margin = 0;
canvas.style.position = 'fixed';

// get canvas 2D context and set him correct size
var ctx = canvas.getContext('2d');
resize();

let points = [];
let unclicked = false;

window.addEventListener('resize', resize);
document.addEventListener('mousemove', draw);
document.addEventListener('mousedown', setPosition);

function getApproxCircle(p) {
  // Régressions coniques, quadriques, régressions linéaires et apparentées, circulaire, sphérique
  // pp 12-13
    let Sx = 0 
    let Sy = 0 
    let Sx2 = 0 
    let Sy2 = 0 
    let Sxy = 0 
    let Sx3 = 0 
    let Sy3 = 0 
    let Sx2y = 0
    let Sxy2 = 0
    for(point of p ){
        const [x,y] = [point.x,point.y]
        const [x2,y2] = [x*x,y*y]
        const [x3,y3] = [x2*x,y2*y]
        const xy    = x*y
        const [x2y,xy2] = [xy*x,xy*y]

        Sx += x
        Sy += y
        Sx2 += x2
        Sy2 += y2
        Sxy += xy
        Sx3 += x3
        Sy3 += y3
        Sx2y += x2y
        Sxy2 += xy2
    }

    const n = points.length;
    const s11 = n * Sxy - Sx * Sy
    const s20 = n * Sx2 - Sx * Sx
    const s02 = n * Sy2 - Sy * Sy
    const s30 = n * Sx3 - Sx2 * Sx
    const s03 = n * Sy3 - Sy * Sy2
    const s21 = n * Sx2y - Sx2 * Sy
    const s12 = n * Sxy2 - Sx * Sy2
    const a = ((s30 + s12)*s02 - (s03 + s21)*s11)
            / (2*( s20 * s02 - s11 * s11));
    const b = ((s03 + s21)*s20 - (s30 + s12)*s11)
            / (2*( s20 * s02 - s11 * s11));
    const c = (Sx2 + Sy2 - 2*a*Sx - 2*b*Sy) / n;
    const radius = Math.sqrt(c + a*a + b*b);
        
    return {
        center: {
            x: a,
          y: b
        },
        radius: radius
    };
}
function setPosition(e) {
  if (unclicked && e.buttons === 1) {
    unclicked = false;
        points = [];
  }
  points.push({ x: e.clientX, y: e.clientY });
}

function resize() {
  ctx.canvas.width = window.innerWidth;
  ctx.canvas.height = window.innerHeight;
}

function clear() {
    ctx.clearRect(0, 0, window.innerWidth, window.innerHeight);
}

function draw(e) {
  if (e.buttons !== 1) { 
    unclicked = true;
    return;
  }
  
  points.push({ x: e.clientX, y: e.clientY });
    clear();  

  ctx.lineWidth = 5;
  ctx.lineCap = 'round';
  ctx.strokeStyle = '#c0392b';
  
  //point of drawn arc
  ctx.beginPath(); 
  ctx.moveTo(points[0].x, points[0].y);
  for(i=1; i<points.length; i++) {
    ctx.lineTo(points[i].x, points[i].y); 
  }
  ctx.stroke(); 
  
  let approxCircle = getApproxCircle(points);
  
  //circle
  ctx.lineWidth = 2;
  ctx.strokeStyle = '#3039cb77';
    ctx.beginPath();
  ctx.arc(approxCircle.center.x, approxCircle.center.y, approxCircle.radius, 0, 2 * Math.PI, false);
  ctx.stroke();
  
  // center of the circle
    ctx.beginPath();
  ctx.arc(approxCircle.center.x, approxCircle.center.y, 1, 0, 2 * Math.PI, false);
  ctx.fillStyle = 'green';
  ctx.fill();
  ctx.lineWidth = 5;
  ctx.strokeStyle = '#003300';
  ctx.stroke();
}

А если нужно нарисовать окружность по 3 точкам, то есть формула из Википедии:

в Котлине (ссылка на песочницу):

data class Point(
        val x: Float,
        val y: Float
) {
    fun lenSquare() = x*x + y*y
}

fun main() {
    val A = Point(0f, 0f)
    val B = Point(0f, 2f)
    val C = Point(2f, 0f)
    
    val D = 2 * (A.x * (B.y - C.y)
                + B.x * (C.y - A.y)
                + C.x * (A.y - B.y))

    var lenASqr = A.lenSquare()
    var lenBSqr = B.lenSquare()
    var lenCSqr = C.lenSquare()
    val Ux = (lenASqr * (B.y - C.y)
                + lenBSqr * (C.y - A.y)
                + lenCSqr * (A.y - B.y)
             ) / D;
    val Uy = (lenASqr * (C.x - B.x)
                + lenBSqr * (A.x - C.x)
                + lenCSqr * (B.x - A.x)
             ) / D;
    val U = Point(Ux, Uy)

    println("A: $A")
    println("B: $B")
    println("C: $C")
    println("Circumcircle center: $U")
}

Спасибо за это ... Я как бы следую матричной алгебре, но не с какой-либо уверенностью. Я не уверен, что смогу эффективно закодировать решение на моем целевом языке (котлин). Я умею читать python — я просто забыл спросить «без матриц, пожалуйста». Базовая теория алгебры имеет для меня смысл, но, поскольку дуги, которые у меня есть, только приближают линию, мне нужно будет выполнять вычисления с несколькими тройками точек и усреднять результаты, поэтому она должна быть эффективной.

DrPhill 18.11.2022 20:35

@DrPhill: я не вижу матриц в опубликованном коде!? Усреднение нескольких решений — не лучший метод. Существуют более глобальные подходы, такие как минимизация по методу наименьших квадратов.

Yves Daoust 19.11.2022 16:09

Я только что обновил свой ответ. Предыдущая версия с матрицами на самом деле не подходила для вопроса.

Renat 19.11.2022 16:10

Спасибо Ренат. Подходящий код круга - это именно то, что я хотел... не зная, что просить. Хотел бы я дважды проголосовать за ваши усилия.

DrPhill 19.11.2022 20:31

Нет проблем, для меня тоже был интересный опыт, я выучил этот новый термин подгонка-проблема, и немного Kotlin, и немного работы с Canvas в JS.

Renat 19.11.2022 20:32

Прохладный. Это заставляет меня задаться вопросом, есть ли лучший алгоритм для эллипса .......

DrPhill 19.11.2022 21:33

@DrPhill, нашел пакет python lsq-ellipse с исходным кодом на Github, github.com/bdhammel/least-squares-ellipse-fitting/blob/…, он использует матричные операции из библиотеки numpy, если необходимо, чтобы перевести его на другой язык, который потребуется для поиска альтернативной библиотеки/реализации для этих операций

Renat 20.11.2022 01:26

Спасибо @Renat, я посмотрю на это, когда у меня будет остальная часть структуры. Обнаружение эллипсов в настоящее время выходит за рамки, но я уверен, что он проберется обратно.....

DrPhill 20.11.2022 09:15

Другие вопросы по теме