2009-02-12 29 views
18

¿Cómo se calcula el área de intersección entre un triángulo (especificado como tres pares (X, Y)) y un círculo (X, Y, R)? He hecho algunas búsquedas en vano. Esto es por trabajo, no por la escuela. :)¿Calcula el área de intersección entre un círculo y un triángulo?

se vería algo como esto en C#:

struct { PointF vert[3]; } Triangle; 
struct { PointF center; float radius; } Circle; 

// returns the area of intersection, e.g.: 
// if the circle contains the triangle, return area of triangle 
// if the triangle contains the circle, return area of circle 
// if partial intersection, figure that out 
// if no intersection, return 0 
double AreaOfIntersection(Triangle t, Circle c) 
{ 
... 
} 

Respuesta

27

Si desea una solución exacta (o al menos tan exacta como usted puede conseguir usando aritmética de punto flotante), entonces esto va a implica mucho trabajo de campo, porque hay tantos casos para considerar.

I cuentan nueve casos diferentes (categorizados en la siguiente figura por el número de vértices del triángulo dentro del círculo, y el número de bordes del triángulo que se cruzan o están contenidos en el círculo):

Nine cases for intersection: 1, 2. no vertices, no edges; 3. no vertices, one edge; 4. no vertices, two edges; 5. no vertices, three edges; 6. one vertex, two edges; 7. one vertex, three edges; 8. two vertices, three edges; 9. three vertices, three edges.

(sin embargo, este tipo de enumeración de los casos geométrica es bien conocido por ser difícil, y no me sorprendería en absoluto si echaba de menos una o dos!)

lo que el enfoque es:

  1. Determine para cada vértice del triángulo si está dentro del círculo. Voy a suponer que sabes cómo hacer eso.

  2. Determine para cada borde del triángulo si se cruza con el círculo. (Escribí un método here, o veo cualquier libro de geometría computacional). Tendrá que calcular el punto o puntos de intersección (si corresponde) para usar en el paso 4.

  3. Determine cuál de los nueve casos tener.

  4. Calcule el área de la intersección. Los casos 1, 2 y 9 son fáciles. En los seis casos restantes dibujé líneas punteadas para mostrar cómo dividir el área de intersección en triángulos y circular segments según los vértices originales del triángulo y en los puntos de intersección que calculó en el paso 2.

Este algoritmo va a ser bastante delicado y propenso a errores que afectan sólo a uno de los casos, por lo que asegúrese de que tiene casos de prueba que cubren los nueve casos (y sugiero permutando los vértices de los triángulos de prueba también). Preste especial atención a los casos en que uno de los vértices del triángulo se encuentra en el borde del círculo.

Si no necesita una solución exacta, rastrillar las cifras y contar los píxeles en la intersección (como lo sugirieron otros encuestados) parece ser un enfoque mucho más sencillo del código y, por consiguiente, menos propenso a errores .

+0

+1 matemática! Parece que la solución exacta se ejecutaría mucho más rápido que un techinque rasterizado también. – Crashworks

+1

Estoy debidamente impresionado por su minuciosidad. –

+0

Tenga en cuenta que la manera más fácil de hacer # 4 y # 5 es tomar el área del círculo y restar los segmentos fuera del triángulo (en lugar de sumar todos los subtrágulos y segmentos dentro de él). Estoy realmente impresionado, Gareth. – Crashworks

1

Suponiendo que estás hablando píxeles enteros, no real, la aplicación ingenua sería colocar a través de cada píxel del triángulo y comprobar la distancia desde el centro del círculo contra su radio.

No es una fórmula linda, o particularmente rápida, pero hace el trabajo bien.

0

Cómo exacta es lo que necesita ser? Si puede aproximar el círculo con formas más simples, puede simplificar el problema. No sería difícil modelar un círculo como un conjunto de triángulos muy estrechos que se reúnen en el centro, por ejemplo.

1

Si usted tiene una GPU a su disposición, usted podría utilizar this técnica para la obtención de un número de píxeles de la intersección ..

0

Si sólo uno de los segmentos de línea del triángulo se cruza con el círculo, la solución pura matemática ISN' demasiado duro. Una vez que sepa cuándo son los dos puntos de intersección, puede usar la fórmula de distancia para encontrar la longitud de la cuerda.

Según these equations:

ϑ = 2 sin⁻¹(0.5 c/r) 
A = 0.5 r² (ϑ - sin(ϑ)) 

donde c es la longitud de la cuerda, r es el radio, θ se convierte en el ángulo a través del centro, y A es el área. Tenga en cuenta que esta solución se rompe si se corta más de la mitad del círculo.

Probablemente no vale la pena el esfuerzo si solo necesita una aproximación, ya que hace varias suposiciones sobre cómo se ve la intersección real.

1

Creo que no debe aproximar el círculo como un conjunto de triángulos, en lugar de eso, puede aproximar su forma con un polígono. El algoritmo ingenuo puede verse como:

  1. Convierte tu círculo en un polígono con el número deseado de vértices.
  2. Calcula la intersección de dos polígonos (círculo convertido y un triángulo).
  3. Calcula el cuadrado de esa intersección.

Puede optimizar este algoritmo combinando los pasos 2 y 3 en una sola función.

Lea este enlaces:
Area of convex polygon
Intersection of convex polygons

0

Mi primer instinto sería transformar todo para que el círculo está centrado en el origen, trans del triángulo a coordenadas polares, y resolver por la intersección (o encompassment) del triángulo con el círculo. Todavía no lo he trabajado en papel, así que eso es solo una corazonada.

+0

Estoy investigando ese enfoque en este momento ... en el caso general, hay una integración bastante fea involucrada. No creo que haya una buena fórmula simple que una computadora pueda calcular. –

+2

Esto se siente como el tipo de cosa que debe haber sido resuelto por algún matemático del siglo XIX, ¡pero desafortunadamente Google Scholar no regresa tan lejos! =) – Crashworks

1

Dado que sus formas son convexas, puede utilizar la estimación de área de Monte Carlo.

Dibuja un cuadro alrededor del círculo y el triángulo.

Elija puntos al azar en la casilla y lleve un recuento de cuántos caen en el círculo, y cuántos caen en el círculo y el triángulo.

área de intersección ≅ Área del círculo * # puntos en el círculo y el triángulo/# puntos en círculo

Detener la elección de puntos cuando la superficie estimada no cambia en más de una cierta cantidad durante un determinado número de rondas o simplemente elija una cantidad fija de puntos según el área de la caja. La estimación del área debe converger bastante rápido a menos que una de tus formas tenga muy poca área.

Nota: Aquí es cómo determinar si un punto está en un triángulo: Barycentric coordinates

2

Llego casi un año y medio tarde, pero pensé que tal vez la gente estaría interesada en code here que escribí y creo que lo hace correctamente. Mire en la función IntersectionArea cerca de la parte inferior. El enfoque general es seleccionar el polígono convexo circunscrito por el círculo, y luego tratar con las pequeñas tapas circulares.

31

Primero, nos recordará cómo encontrar el área de un polígono. Una vez que hayamos hecho esto, el algoritmo para encontrar la intersección entre un polígono y un círculo debería ser fácil de entender.

cómo encontrar el área de un polígono

Veamos el caso de un triángulo, porque toda la lógica esencial aparece allí. Supongamos que tenemos un triángulo con vértices (x1, y1), (x2, y2) y (x3, y3) a medida que avanza alrededor del triángulo hacia la izquierda, como se muestra en la siguiente figura: triangleFigure

Entonces puede calcular el área con la fórmula

A = (x1 y2 + x2 y3 + x3 y1 - x2y1- x3 y2 - x1y3)/2.

ver por qué funciona esta fórmula, vamos a reorganizar lo que es en la forma

A = (x1 y2 - x2 y1)/2 + (x2 y3 - x3 y2)/2 + (y1 x3 - x1y3)/2.

Ahora el primer término es el área siguiente, lo cual es positivo en nuestro caso: enter image description here

Si no está claro que el área de la región verde es de hecho (x1 y2 - y1 x2)/2 , luego lea this.

el segundo término es esta zona, que es de nuevo positivo:

enter image description here

Y la tercera área se muestra en la siguiente figura. Esta vez, el área es negativo

enter image description here

Adición de estos tres arriba obtenemos la siguiente imagen

enter image description here

Vemos que el área verde que estaba fuera del triángulo es cancelada por el área roja , de modo que el área neta es solo el área del triángulo, y esto muestra por qué nuestra fórmula fue verdadera en este caso.

Lo que dije antes fue la explicación intuitiva de por qué la fórmula del área era correcta. Una explicación más rigurosa sería observar que cuando calculamos el área desde un borde, el área que obtenemos es la misma que obtendríamos de la integración r^2dθ/2, de modo que estamos integrando r^2dθ/2 de forma efectiva alrededor del límite del polígono, y por el teorema de Stokes, esto da el mismo resultado que la integración de rdrdθ sobre la región delimitada por el polígono. Dado que la integración de rdrdθ sobre la región delimitada por el polígono da el área, concluimos que nuestro procedimiento debe proporcionar correctamente el área.

Área de la intersección de un círculo con un polígono

Ahora vamos a hablar de cómo encontrar el área de la intersección de un círculo de radio R con un polígono como se muestra en la siguiente figura:

enter image description here

Estamos interesados ​​en encontrar el área de la región verde. Podemos, como en el caso del polígono único, dividir nuestro cálculo para encontrar un área para cada lado del polígono y luego agregar esas áreas.

Nuestra primera área se verá como: enter image description here

La segunda área se verá como enter image description here

Y la tercera área se enter image description here

Una vez más, las dos primeras áreas son positivos en nuestro caso, mientras que el tercero será negativo. Afortunadamente, las cancelaciones funcionarán para que el área neta sea de hecho el área en la que estamos interesados. Veamos.

enter image description here

De hecho, la suma de las áreas que será el área que nos interesa.

Una vez más, podemos dar una explicación más rigurosa de por qué esto funciona. Deje que yo sea la región definida por la intersección y deje que P sea el polígono. Luego, de la discusión anterior, sabemos que queremos computarizar la integral de r^2dθ/2 alrededor del límite de I. Sin embargo, esto es difícil de hacer porque requiere encontrar la intersección.

En cambio, hicimos una integral sobre el polígono. Integramos max (r, R)^2 dθ/2 sobre el límite del polígono. Para ver por qué esto da la respuesta correcta, definamos una función π que toma un punto en coordenadas polares (r, θ) hasta el punto (max (r, R), θ). No debería ser confuso referirse a las funciones de coordenadas de π (r) = max (r, R) y π (θ) = θ. Entonces, lo que hicimos fue integrar π (r)^2 dθ/2 sobre el límite del polígono.

Por otro lado, como π (θ) = θ, esto es lo mismo que integrar π (r)^2 dπ (θ)/2 sobre el límite del polígono.

Ahora haciendo un cambio de variable, encontramos que obtendríamos la misma respuesta si integramos r^2 dθ/2 sobre el límite de π (P), donde π (P) es la imagen de P en π .

Usando el teorema de Stokes nuevamente, sabemos que la integración de r^2 dθ/2 sobre el límite de π (P) nos da el área de π (P). En otras palabras, da la misma respuesta que integrar dxdy sobre π (P).

Usando un cambio de variable nuevamente, sabemos que la integración de dxdy sobre π (P) es lo mismo que integrar Jdxdy sobre P, donde J es el jacobiano de π.

Ahora podemos dividir la integral de Jdxdy en dos regiones: la parte en el círculo y la parte fuera del círculo. Ahora π deja puntos en el círculo solo así que J = 1 allí, por lo que la contribución de esta parte de P es el área de la parte de P que se encuentra en el círculo, es decir, el área de la intersección. La segunda región es la región fuera del círculo. Hay J = 0 ya que π colapsa esta parte hasta el límite del círculo.

Por lo tanto, lo que computamos es de hecho el área de la intersección.

Ahora que estamos relativamente seguros de que sabemos conceptualmente cómo encontrar el área, hablemos más específicamente sobre cómo calcular la contribución de un solo segmento. Comencemos mirando un segmento en lo que llamaré "geometría estándar". Se muestra a continuación.

enter image description here

En la geometría estándar, el borde va horizontalmente de izquierda a derecha. Se describe con tres números: xi, la coordenada x donde comienza el borde, xf, la coordenada x donde termina el borde, y y, la coordenada y del borde.

Ahora vemos que si | y | < R, como en la figura, entonces el borde intersectará el círculo en los puntos (-xint, y) y (xint, y) donde xint = (R^2-y^2)^(1/2). Luego, el área que necesitamos calcular se divide en tres partes que se etiquetan en la figura. Para obtener las áreas de las regiones 1 y 3, podemos usar arctan para obtener los ángulos de los diversos puntos y luego igualar el área a R^2 Δθ/2. Entonces, por ejemplo, estableceríamos θi = atan2 (y, xi) y θl = atan2 (y, -xint). Entonces el área de la región uno es R^2 (θl-θi)/2. Podemos obtener el área de la región 3 de manera similar.

El área de la región 2 es simplemente el área de un triángulo. Sin embargo, debemos tener cuidado con la señal. Queremos que el área muestre ser positiva, por lo que diremos que el área es - (xint - (-xint)) y/2.

Otra cosa a tener en cuenta es que, en general, xi no tiene que ser menor que -xint y xf no tiene que ser mayor que xint.

El otro caso a considerar es | y | > R. Este caso es más simple, porque solo hay una pieza que es similar a la región 1 en la figura.

Ahora que sabemos cómo calcular el área desde un borde en la geometría estándar, lo único que queda por hacer es describir cómo transformar cualquier borde en geometría estándar.

Pero esto solo es un simple cambio de coordenadas. Dado que algunos tienen el vértice inicial vi y el vértice final vf, el nuevo vector unidad x será el vector unitario que apunta de vi a vf. Entonces xi es solo el desplazamiento de vi desde el centro del círculo punteado en x, y xf es solo xi más la distancia entre vi y vf. Mientras tanto, y está dada por el producto de cuña de x con el desplazamiento de vi desde el centro del círculo.

Código

que completa la descripción del algoritmo, ahora es el momento de escribir código. Usaré Java.

En primer lugar, ya que estamos trabajando con los círculos, debemos tener una clase de círculo

public class Circle { 

    final Point2D center; 
    final double radius; 

    public Circle(double x, double y, double radius) { 
     center = new Point2D.Double(x, y); 
     this.radius = radius; 
    } 

    public Circle(Point2D.Double center, double radius) { 
     this(center.getX(), center.getY(), radius); 
    } 

    public Point2D getCenter() { 
     return new Point2D.Double(getCenterX(), getCenterY()); 
    } 

    public double getCenterX() { 
     return center.getX(); 
    } 

    public double getCenterY() { 
     return center.getY(); 
    } 

    public double getRadius() { 
     return radius; 
    } 

} 

Para polígonos, voy a utilizar la clase de Java Shape. Shape s tienen un PathIterator que puedo usar para recorrer los bordes del polígono.

Ahora para el trabajo real. Separará la lógica de iteración a través de los bordes, poniendo los bordes en geometría estándar, etc., a partir de la lógica de computación del área una vez que se hace esto. La razón para esto es que en el futuro puede que desee calcular algo además o además del área y desea poder reutilizar el código teniendo que lidiar con iterar a través de los bordes.

Así que tengo una clase genérica que calcula alguna propiedad de la clase T sobre nuestra intersección de círculo de polígono.

public abstract class CircleShapeIntersectionFinder<T> {

Tiene tres métodos estáticos que solo ayudan a la geometría de cómputo:

private static double[] displacment2D(final double[] initialPoint, final double[] finalPoint) { 
    return new double[]{finalPoint[0] - initialPoint[0], finalPoint[1] - initialPoint[1]}; 
} 

private static double wedgeProduct2D(final double[] firstFactor, final double[] secondFactor) { 
    return firstFactor[0] * secondFactor[1] - firstFactor[1] * secondFactor[0]; 
} 

static private double dotProduct2D(final double[] firstFactor, final double[] secondFactor) { 
    return firstFactor[0] * secondFactor[0] + firstFactor[1] * secondFactor[1]; 
} 

Hay dos campos de instancia, un Circle, que sólo mantiene una copia del círculo, y el currentSquareRadius, que mantiene una copia del radio cuadrado Esto puede parecer extraño, pero la clase que estoy usando está realmente equipada para encontrar las áreas de una colección completa de intersecciones círculo-polígono. Es por eso que me refiero a uno de los círculos como "actual".

private Circle currentCircle; 
private double currentSquareRadius; 

A continuación viene el método para el cálculo de lo que queremos calcular:

public final T computeValue(Circle circle, Shape shape) { 
    initialize(); 
    processCircleShape(circle, shape); 
    return getValue(); 
} 

initialize() y getValue() son abstractas. initialize() establecería la variable que mantiene un total del área en cero, y getValue() simplemente devolvería el área. La definición para processCircleShape es

private void processCircleShape(Circle circle, final Shape cellBoundaryPolygon) { 
    initializeForNewCirclePrivate(circle); 
    if (cellBoundaryPolygon == null) { 
     return; 
    } 
    PathIterator boundaryPathIterator = cellBoundaryPolygon.getPathIterator(null); 
    double[] firstVertex = new double[2]; 
    double[] oldVertex = new double[2]; 
    double[] newVertex = new double[2]; 
    int segmentType = boundaryPathIterator.currentSegment(firstVertex); 
    if (segmentType != PathIterator.SEG_MOVETO) { 
     throw new AssertionError(); 
    } 
    System.arraycopy(firstVertex, 0, newVertex, 0, 2); 
    boundaryPathIterator.next(); 
    System.arraycopy(newVertex, 0, oldVertex, 0, 2); 
    segmentType = boundaryPathIterator.currentSegment(newVertex); 
    while (segmentType != PathIterator.SEG_CLOSE) { 
     processSegment(oldVertex, newVertex); 
     boundaryPathIterator.next(); 
     System.arraycopy(newVertex, 0, oldVertex, 0, 2); 
     segmentType = boundaryPathIterator.currentSegment(newVertex); 
    } 
    processSegment(newVertex, firstVertex); 
} 

Tomemos un segundo para mirar initializeForNewCirclePrivate rápidamente. Este método simplemente establece los campos de instancia y permite que la clase derivada almacene cualquier propiedad del círculo. Su definición es

private void initializeForNewCirclePrivate(Circle circle) { 
    currentCircle = circle; 
    currentSquareRadius = currentCircle.getRadius() * currentCircle.getRadius(); 
    initializeForNewCircle(circle); 
} 

initializeForNewCircle es abstracta y una implementación sería para él para almacenar el radio de los círculos para evitar tener que hacer raíces cuadradas. De todas formas, volvamos al processCircleShape. Después de llamar al initializeForNewCirclePrivate, verificamos si el polígono es null (que interpreto como un polígono vacío) y lo devolvemos si es null. En este caso, nuestra área calculada sería cero. Si el polígono no es null, obtenemos el PathIterator del polígono. El argumento del método getPathIterator al que llamo es una transformación afín que se puede aplicar a la ruta. No quiero aplicar uno, así que acabo de pasar null.

A continuación declaro el double[] s que hará un seguimiento de los vértices. Debo recordar el primer vértice porque el PathIterator solo me da cada vértice una vez, así que tengo que regresar después de que me haya dado el último vértice, y formar un borde con este último vértice y el primer vértice.

El método currentSegment en la línea siguiente coloca el siguiente vértice en su argumento. Devuelve un código que le dice cuando está fuera de los vértices. Esta es la razón por la cual la expresión de control para mi ciclo while es lo que es.

La mayor parte del resto del código de este método no tiene ninguna lógica relacionada con la iteración a través de los vértices.Lo importante es que una vez por iteración del ciclo while llamo a processSegment y luego llamo al processSegment nuevamente al final del método para procesar el borde que conecta el último vértice con el primer vértice. mirada

Vamos en el código para processSegment:

private void processSegment(double[] initialVertex, double[] finalVertex) { 
    double[] segmentDisplacement = displacment2D(initialVertex, finalVertex); 
    if (segmentDisplacement[0] == 0 && segmentDisplacement[1] == 0) { 
     return; 
    } 
    double segmentLength = Math.sqrt(dotProduct2D(segmentDisplacement, segmentDisplacement)); 
    double[] centerToInitialDisplacement = new double[]{initialVertex[0] - getCurrentCircle().getCenterX(), initialVertex[1] - getCurrentCircle().getCenterY()}; 
    final double leftX = dotProduct2D(centerToInitialDisplacement, segmentDisplacement)/segmentLength; 
    final double rightX = leftX + segmentLength; 
    final double y = wedgeProduct2D(segmentDisplacement, centerToInitialDisplacement)/segmentLength; 
    processSegmentStandardGeometry(leftX, rightX, y); 
} 

En este método que implementar los pasos para transformar una ventaja en la geometría estándar como se ha descrito anteriormente. Primero calculo segmentDisplacement, el desplazamiento desde el vértice inicial hasta el vértice final. Esto define el eje x de la geometría estándar. Hago un regreso temprano si este desplazamiento es cero.

Siguiente Calculo la longitud del desplazamiento, porque esto es necesario para obtener el vector x unidad. Una vez que tengo esta información, calculo el desplazamiento desde el centro del círculo hasta el vértice inicial. El producto escalar de esto con segmentDisplacement me da leftX que había estado llamando xi. Entonces rightX, que había estado llamando xf, es solo leftX + segmentLength. Finalmente hago el producto de cuña para obtener y como se describió anteriormente.

Ahora que he transformado el problema en la geometría estándar, será fácil de tratar. Eso es lo que hace el método processSegmentStandardGeometry. Veamos el código

private void processSegmentStandardGeometry(double leftX, double rightX, double y) { 
    if (y * y > getCurrentSquareRadius()) { 
     processNonIntersectingRegion(leftX, rightX, y); 
    } else { 
     final double intersectionX = Math.sqrt(getCurrentSquareRadius() - y * y); 
     if (leftX < -intersectionX) { 
      final double leftRegionRightEndpoint = Math.min(-intersectionX, rightX); 
      processNonIntersectingRegion(leftX, leftRegionRightEndpoint, y); 
     } 
     if (intersectionX < rightX) { 
      final double rightRegionLeftEndpoint = Math.max(intersectionX, leftX); 
      processNonIntersectingRegion(rightRegionLeftEndpoint, rightX, y); 
     } 
     final double middleRegionLeftEndpoint = Math.max(-intersectionX, leftX); 
     final double middleRegionRightEndpoint = Math.min(intersectionX, rightX); 
     final double middleRegionLength = Math.max(middleRegionRightEndpoint - middleRegionLeftEndpoint, 0); 
     processIntersectingRegion(middleRegionLength, y); 
    } 
} 

La primera if distingue los casos en los que y es lo suficientemente pequeño que el borde puede intersectar el círculo. Si y es grande y no hay posibilidad de intersección, entonces invoco el método para manejar ese caso. De lo contrario, manejo el caso donde la intersección es posible.

Si la intersección es posible, calculo la coordenada x de la intersección, intersectionX, y divido el borde en tres porciones, que corresponden a las regiones 1, 2 y 3 de la figura de geometría estándar anterior. Primero manejo la región 1.

Para manejar la región 1, verifico si leftX es de hecho menor que -intersectionX, de lo contrario no habría región 1. Si hay una región 1, entonces necesito saber cuándo termina. Termina al mínimo de rightX y -intersectionX. Después de encontrar estas coordenadas x, me ocupo de esta región sin intersecciones.

que hacer algo similar para manejar la región 3.

para la Región 2, que tengo que hacer algo de lógica para comprobar que leftX y rightX tiene que realmente soporte de alguna región en el medio y -intersectionXintersectionX. Después de encontrar la región, sólo necesito la longitud de la región y y, así que pasar estos dos números a un método abstracto que se ocupa de la región 2.

Ahora veamos el código para processNonIntersectingRegion

private void processNonIntersectingRegion(double leftX, double rightX, double y) { 
    final double initialTheta = Math.atan2(y, leftX); 
    final double finalTheta = Math.atan2(y, rightX); 
    double deltaTheta = finalTheta - initialTheta; 
    if (deltaTheta < -Math.PI) { 
     deltaTheta += 2 * Math.PI; 
    } else if (deltaTheta > Math.PI) { 
     deltaTheta -= 2 * Math.PI; 
    } 
    processNonIntersectingRegion(deltaTheta); 
} 

Simplemente uso atan2 para calcular la diferencia en el ángulo entre leftX y rightX. Luego agrego un código para tratar la discontinuidad en atan2, pero esto probablemente sea innecesario, ya que la discontinuidad ocurre a 180 grados o 0 grados. Luego paso la diferencia de ángulo en un método abstracto.Por último sólo tenemos métodos abstractos y captadores:

protected abstract void initialize(); 

    protected abstract void initializeForNewCircle(Circle circle); 

    protected abstract void processNonIntersectingRegion(double deltaTheta); 

    protected abstract void processIntersectingRegion(double length, double y); 

    protected abstract T getValue(); 

    protected final Circle getCurrentCircle() { 
     return currentCircle; 
    } 

    protected final double getCurrentSquareRadius() { 
     return currentSquareRadius; 
    } 

} 

Ahora vamos a ver la clase que se extiende, CircleAreaFinder

public class CircleAreaFinder extends CircleShapeIntersectionFinder<Double> { 

public static double findAreaOfCircle(Circle circle, Shape shape) { 
    CircleAreaFinder circleAreaFinder = new CircleAreaFinder(); 
    return circleAreaFinder.computeValue(circle, shape); 
} 

double area; 

@Override 
protected void initialize() { 
    area = 0; 
} 

@Override 
protected void processNonIntersectingRegion(double deltaTheta) { 
    area += getCurrentSquareRadius() * deltaTheta/2; 
} 

@Override 
protected void processIntersectingRegion(double length, double y) { 
    area -= length * y/2; 
} 

@Override 
protected Double getValue() { 
    return area; 
} 

@Override 
protected void initializeForNewCircle(Circle circle) { 

} 

}

Tiene un campo area hacer un seguimiento de la zona. initialize establece el área en cero, como se esperaba. Cuando procesamos un borde que no se intersecta, incrementamos el área en R^2 Δθ/2, ya que concluimos que debemos hacerlo más arriba. Para un borde de intersección, disminuimos el área por y*length/2. Esto fue para que los valores negativos para y correspondan a áreas positivas, como decidimos que deberían.

Lo bueno es que si queremos hacer un seguimiento del perímetro, no tenemos que hacer mucho más trabajo. He definido una clase AreaPerimeter:

public class AreaPerimeter { 

    final double area; 
    final double perimeter; 

    public AreaPerimeter(double area, double perimeter) { 
     this.area = area; 
     this.perimeter = perimeter; 
    } 

    public double getArea() { 
     return area; 
    } 

    public double getPerimeter() { 
     return perimeter; 
    } 

} 

y ahora sólo tenemos que ampliar nuestra clase abstracta de nuevo utilizando AreaPerimeter como el tipo.

public class CircleAreaPerimeterFinder extends CircleShapeIntersectionFinder<AreaPerimeter> { 

    public static AreaPerimeter findAreaPerimeterOfCircle(Circle circle, Shape shape) { 
     CircleAreaPerimeterFinder circleAreaPerimeterFinder = new CircleAreaPerimeterFinder(); 
     return circleAreaPerimeterFinder.computeValue(circle, shape); 
    } 

    double perimeter; 
    double radius; 
    CircleAreaFinder circleAreaFinder; 

    @Override 
    protected void initialize() { 
     perimeter = 0; 
     circleAreaFinder = new CircleAreaFinder(); 
    } 

    @Override 
    protected void initializeForNewCircle(Circle circle) { 
     radius = Math.sqrt(getCurrentSquareRadius()); 
    } 

    @Override 
    protected void processNonIntersectingRegion(double deltaTheta) { 
     perimeter += deltaTheta * radius; 
     circleAreaFinder.processNonIntersectingRegion(deltaTheta); 
    } 

    @Override 
    protected void processIntersectingRegion(double length, double y) { 
     perimeter += Math.abs(length); 
     circleAreaFinder.processIntersectingRegion(length, y); 
    } 

    @Override 
    protected AreaPerimeter getValue() { 
     return new AreaPerimeter(circleAreaFinder.getValue(), perimeter); 
    } 

} 

tenemos una variable perimeter para realizar un seguimiento del perímetro, nos recuerda el valor de la radius para evitar tener que llamar Math.sqrt mucho, y que delegar el cálculo del área en nuestro CircleAreaFinder. Podemos ver que las fórmulas para el perímetro son fáciles.

Como referencia aquí es el código completo de CircleShapeIntersectionFinder

private static double[] displacment2D(final double[] initialPoint, final double[] finalPoint) { 
     return new double[]{finalPoint[0] - initialPoint[0], finalPoint[1] - initialPoint[1]}; 
    } 

    private static double wedgeProduct2D(final double[] firstFactor, final double[] secondFactor) { 
     return firstFactor[0] * secondFactor[1] - firstFactor[1] * secondFactor[0]; 
    } 

    static private double dotProduct2D(final double[] firstFactor, final double[] secondFactor) { 
     return firstFactor[0] * secondFactor[0] + firstFactor[1] * secondFactor[1]; 
    } 

    private Circle currentCircle; 
    private double currentSquareRadius; 

    public final T computeValue(Circle circle, Shape shape) { 
     initialize(); 
     processCircleShape(circle, shape); 
     return getValue(); 
    } 

    private void processCircleShape(Circle circle, final Shape cellBoundaryPolygon) { 
     initializeForNewCirclePrivate(circle); 
     if (cellBoundaryPolygon == null) { 
      return; 
     } 
     PathIterator boundaryPathIterator = cellBoundaryPolygon.getPathIterator(null); 
     double[] firstVertex = new double[2]; 
     double[] oldVertex = new double[2]; 
     double[] newVertex = new double[2]; 
     int segmentType = boundaryPathIterator.currentSegment(firstVertex); 
     if (segmentType != PathIterator.SEG_MOVETO) { 
      throw new AssertionError(); 
     } 
     System.arraycopy(firstVertex, 0, newVertex, 0, 2); 
     boundaryPathIterator.next(); 
     System.arraycopy(newVertex, 0, oldVertex, 0, 2); 
     segmentType = boundaryPathIterator.currentSegment(newVertex); 
     while (segmentType != PathIterator.SEG_CLOSE) { 
      processSegment(oldVertex, newVertex); 
      boundaryPathIterator.next(); 
      System.arraycopy(newVertex, 0, oldVertex, 0, 2); 
      segmentType = boundaryPathIterator.currentSegment(newVertex); 
     } 
     processSegment(newVertex, firstVertex); 
    } 

    private void initializeForNewCirclePrivate(Circle circle) { 
     currentCircle = circle; 
     currentSquareRadius = currentCircle.getRadius() * currentCircle.getRadius(); 
     initializeForNewCircle(circle); 
    } 

    private void processSegment(double[] initialVertex, double[] finalVertex) { 
     double[] segmentDisplacement = displacment2D(initialVertex, finalVertex); 
     if (segmentDisplacement[0] == 0 && segmentDisplacement[1] == 0) { 
      return; 
     } 
     double segmentLength = Math.sqrt(dotProduct2D(segmentDisplacement, segmentDisplacement)); 
     double[] centerToInitialDisplacement = new double[]{initialVertex[0] - getCurrentCircle().getCenterX(), initialVertex[1] - getCurrentCircle().getCenterY()}; 
     final double leftX = dotProduct2D(centerToInitialDisplacement, segmentDisplacement)/segmentLength; 
     final double rightX = leftX + segmentLength; 
     final double y = wedgeProduct2D(segmentDisplacement, centerToInitialDisplacement)/segmentLength; 
     processSegmentStandardGeometry(leftX, rightX, y); 
    } 

    private void processSegmentStandardGeometry(double leftX, double rightX, double y) { 
     if (y * y > getCurrentSquareRadius()) { 
      processNonIntersectingRegion(leftX, rightX, y); 
     } else { 
      final double intersectionX = Math.sqrt(getCurrentSquareRadius() - y * y); 
      if (leftX < -intersectionX) { 
       final double leftRegionRightEndpoint = Math.min(-intersectionX, rightX); 
       processNonIntersectingRegion(leftX, leftRegionRightEndpoint, y); 
      } 
      if (intersectionX < rightX) { 
       final double rightRegionLeftEndpoint = Math.max(intersectionX, leftX); 
       processNonIntersectingRegion(rightRegionLeftEndpoint, rightX, y); 
      } 
      final double middleRegionLeftEndpoint = Math.max(-intersectionX, leftX); 
      final double middleRegionRightEndpoint = Math.min(intersectionX, rightX); 
      final double middleRegionLength = Math.max(middleRegionRightEndpoint - middleRegionLeftEndpoint, 0); 
      processIntersectingRegion(middleRegionLength, y); 
     } 
    } 

    private void processNonIntersectingRegion(double leftX, double rightX, double y) { 
     final double initialTheta = Math.atan2(y, leftX); 
     final double finalTheta = Math.atan2(y, rightX); 
     double deltaTheta = finalTheta - initialTheta; 
     if (deltaTheta < -Math.PI) { 
      deltaTheta += 2 * Math.PI; 
     } else if (deltaTheta > Math.PI) { 
      deltaTheta -= 2 * Math.PI; 
     } 
     processNonIntersectingRegion(deltaTheta); 
    } 

    protected abstract void initialize(); 

    protected abstract void initializeForNewCircle(Circle circle); 

    protected abstract void processNonIntersectingRegion(double deltaTheta); 

    protected abstract void processIntersectingRegion(double length, double y); 

    protected abstract T getValue(); 

    protected final Circle getCurrentCircle() { 
     return currentCircle; 
    } 

    protected final double getCurrentSquareRadius() { 
     return currentSquareRadius; 
    } 

De todos modos, esa es mi descripción del algoritmo. Creo que es bueno porque es exacto y no hay muchos casos que verificar.

+5

¡Respuesta intensa!Debería tenerlo por separado en una publicación de blog Creo que – Lodewijk

+2

Creo que la cantidad de tiempo y esfuerzo para poner esta respuesta bien merecen una apreciación. Y aquí está el mío. ¡Gracias! – Alp

+0

Desearía poder votar esto dos veces. ¡Gracias! – Daerst

Cuestiones relacionadas