Suponiendo una serie de puntos en el espacio 2d que no se intersecan, ¿cuál es un método eficiente para determinar el área del polígono resultante?

Como nota al margen, esto no es tarea y no estoy buscando código. Estoy buscando una descripción que pueda usar para implementar mi propio método. Tengo mis ideas sobre cómo extraer una secuencia de triángulos de la lista de puntos, pero sé que hay un montón de casos extremos con respecto a polígonos convexos y cóncavos que probablemente no captaré.

Aquí está el método estándar , AFAIK. Básicamente, suma los productos cruzados alrededor de cada vértice. Mucho más simple que la triangulación.

Código de Python, dado un polígono representado como una lista de coordenadas de vértice (x, y), envolviendo implícitamente desde el último vértice hasta el primero:

def area(p):
    return 0.5 * abs(sum(x0*y1 - x1*y0
                         for ((x0, y0), (x1, y1)) in segments(p)))

def segments(p):
    return zip(p, p[1:] + [p[0]])

David Lehavi comenta: Vale la pena mencionar por qué funciona este algoritmo: es una aplicación del teorema de Green para las funciones −y y x; exactamente como funciona un planímetro . Más específicamente:

Fórmula anterior =
integral_over_perimeter(-y dx + x dy) =
integral_over_area((-(-dy)/dy+dx/dx) dy dx) =
2 Area

El producto cruzado es un clásico.

Si tiene un trillón de tales cálculos por hacer, pruebe la siguiente versión optimizada que requiere la mitad de multiplicaciones menos:

area = 0;
for( i = 0; i < N; i += 2 )
   area += x[i+1]*(y[i+2]-y[i]) + y[i+1]*(x[i]-x[i+2]);
area /= 2;

Utilizo un subíndice de matriz para mayor claridad. Es más eficaz utilizar punteros. Aunque los buenos compiladores lo harán por usted.

Se asume que el polígono está "cerrado", lo que significa que copia el primer punto como punto con subíndice N. También asume que el polígono tiene un número par de puntos. Agregue una copia adicional del primer punto si N no es par.

El algoritmo se obtiene desenrollando y combinando dos iteraciones sucesivas del algoritmo clásico de productos cruzados.

No estoy tan seguro de cómo se comparan los dos algoritmos con respecto a la precisión numérica. Mi impresión es que el algoritmo anterior es mejor que el clásico porque la multiplicación tiende a restaurar la pérdida de precisión de la resta. Cuando está limitado a usar flotadores, como con GPU, esto puede marcar una diferencia significativa.

EDITAR: "Área de triángulos y polígonos 2D y 3D" describe un método aún más eficiente

// "close" polygon
x[N] = x[0];
x[N+1] = x[1];
y[N] = y[0];
y[N+1] = y[1];

// compute area
area = 0;
for( size_t i = 1; i <= N; ++i )
  area += x[i]*( y[i+1] - y[i-1] );
area /= 2;

Esta página muestra que la fórmula

se puede simplificar a:

Si escribe algunos términos y los agrupa de acuerdo con factores comunes de xi, la igualdad no es difícil de ver.

La suma final es más eficiente ya que solo requiere nmultiplicaciones en lugar de 2n.

def area(x, y):
    return abs(sum(x[i] * (y[i + 1] - y[i - 1]) for i in xrange(-1, len(x) - 1))) / 2.0

Aprendí esta simplificación de Joe Kington, aquí .

Si tiene NumPy, esta versión es más rápida (para todas las matrices, excepto las muy pequeñas):

def area_np(x, y):        
    x = np.asanyarray(x)
    y = np.asanyarray(y)
    n = len(x)
    shift_up = np.arange(-n+1, 1)
    shift_down = np.arange(-1, n-1)    
    return (x * (y.take(shift_up) - y.take(shift_down))).sum() / 2.0

Un conjunto de puntos sin ninguna otra restricción no necesariamente define de forma única un polígono.

Entonces, primero debe decidir qué polígono construir desde estos puntos, ¿quizás el casco convexo? http://en.wikipedia.org/wiki/Convex_hull

Luego triangula y calcula el área. http://www.mathopenref.com/polygonirregulararea.html

Para expandir el triangular y sumar áreas de triángulo, funcionan si tiene un polígono convexo O si elige un punto que no genera líneas a todos los demás puntos que se cruzan con el polígono.

Para un polígono general que no se interseca, debe sumar el producto cruzado de los vectores (punto de referencia, punto a), (punto de referencia, punto b) donde ayb están "próximos" entre sí.

Suponiendo que tiene una lista de puntos que definen el polígono en orden (el orden son los puntos i e i + 1 que forman una línea del polígono):

Suma (producto cruzado ((punto 0, punto i), (punto 0, punto i + 1)) para i = 1 an - 1.

Toma la magnitud de ese producto cruzado y tienes el área de superficie.

Esto manejará polígonos cóncavos sin tener que preocuparse por elegir un buen punto de referencia; cualesquiera tres puntos que generen un triángulo que no esté dentro del polígono tendrán un producto cruzado que apunte en la dirección opuesta a cualquier triángulo que esté dentro del polígono, por lo que las áreas se suman correctamente.

Para calcular el área del polígono

http://community.topcoder.com/tc?module=Static&d1=tutorials&d2=geometry1#polygon_area

int cross(vct a,vct b,vct c)
{
    vct ab,bc;
    ab=b-a;
    bc=c-b;
    return ab.x*bc.y-ab.y*bc.x;
}    
double area(vct p[],int n)
{ 
    int ar=0;
    for(i=1;i+1<n;i++)
    {
        vct a=p[i]-p[0];
        vct b=p[i+1]-p[0];
        area+=cross(a,b);
    }
    return abs(area/2.0);
}    

O haz una integral de contorno. El teorema de Stokes le permite expresar una integral de área como una integral de contorno. Una pequeña cuadratura de Gauss y Bob es tu tío.

solución independiente del idioma:

DADO: un polígono SIEMPRE puede estar compuesto por n-2 triángulos que no se superponen (n = número de puntos O lados). 1 triángulo = polígono de 3 lados = 1 triángulo; 1 cuadrado = polígono de 4 lados = 2 triángulos; etc ad nauseam QED

por lo tanto, un polígono se puede reducir "cortando" triángulos y el área total será la suma de las áreas de estos triángulos. Pruébelo con un papel y unas tijeras, es mejor si puede visualizar el proceso antes de seguir.

Si toma 3 puntos consecutivos en una ruta de polígonos y crea un triángulo con estos puntos, tendrá uno y solo uno de los tres escenarios posibles:

1. El triángulo resultante está completamente dentro del polígono original.
2. El triángulo resultante está totalmente fuera del polígono original.
3. El triángulo resultante está parcialmente contenido en el polígono original.

solo nos interesan los casos que caen en la primera opción (totalmente contenido).

cada vez que encontramos uno de estos, lo cortamos, calculamos su área (fácil, no explicaré la fórmula aquí) y hacemos un nuevo polígono con un lado menos (equivalente al polígono con este triángulo cortado). hasta que solo nos quede un triángulo.

cómo implementar esto programáticamente:

Cree una matriz de puntos (consecutivos) que representen la ruta ALREDEDOR del polígono. comience en el punto 0. ejecute la matriz formando triángulos (uno a la vez) desde los puntos x, x + 1 y x + 2. transforma cada triángulo de una forma a un área e interseca con el área creada a partir de un polígono. SI la intersección resultante es idéntica al triángulo original, entonces dicho triángulo está totalmente contenido en un polígono y se puede cortar. elimine x + 1 de la matriz y comience de nuevo desde x = 0. de lo contrario (si el triángulo está fuera [parcial o completamente] del polígono), muévase al siguiente punto x + 1 en la matriz.

Además, si está buscando integrarse con el mapeo y está comenzando desde geopuntos, debe convertir de geopuntos a puntos de pantalla PRIMERO. esto requiere decidir un modelo y una fórmula para la forma de la tierra (aunque tendemos a pensar en la tierra como una esfera, en realidad es un ovoide irregular (forma de huevo), con abolladuras). hay muchos modelos por ahí, para más información wiki. una cuestión importante es si considerará o no que el área es un plano o una curva. en general, las áreas "pequeñas", donde los puntos están separados hasta algunos kilómetros, no generarán un error significativo si se consideran planas y no convexas.

Una forma de hacerlo sería descomponer el polígono en triángulos , calcular el área de los triángulos y tomar la suma como el área del polígono.

1. Establezca un punto base (el punto más convexo). Este será su punto de pivote de los triángulos.
2. Calcule el punto más a la izquierda (arbitrario), que no sea su punto base.
3. Calcula el segundo punto más a la izquierda para completar tu triángulo.
4. Guarde esta área triangulada.
5. Desplaza un punto hacia la derecha en cada iteración.
6. Suma las áreas trianguladas

Mejor que sumar triángulos es sumar trapezoides en el espacio cartesiano:

area = 0;
for (i = 0; i < n; i++) {
  i1 = (i + 1) % n;
  area += (vertex[i].y + vertex[i1].y) * (vertex[i1].x - vertex[i].x) / 2.0;
}

La implementación de la fórmula de Shoelace se podría realizar en Numpy. Asumiendo estos vértices:

import numpy as np
x = np.arange(0,1,0.001)
y = np.sqrt(1-x**2)

Podemos definir la siguiente función para encontrar el área:

def PolyArea(x,y):
    return 0.5*np.abs(np.dot(x,np.roll(y,1))-np.dot(y,np.roll(x,1)))

Y obteniendo resultados:

print PolyArea(x,y)
# 0.26353377782163534

Evitar el bucle hace que esta función sea ~ 50 veces más rápida que PolygonArea:

%timeit PolyArea(x,y)
# 10000 loops, best of 3: 42 µs per loop
%timeit PolygonArea(zip(x,y))
# 100 loops, best of 3: 2.09 ms per loop

Nota: había escrito esta respuesta para otra pregunta , solo menciono esto aquí para tener una lista completa de soluciones.

Mi inclinación sería simplemente comenzar a cortar triángulos. No veo cómo otra cosa podría evitar ser terriblemente peluda.

Tome tres puntos secuenciales que componen el polígono. Asegúrese de que el ángulo sea inferior a 180. Ahora tiene un nuevo triángulo que no debería ser un problema para calcular, elimine el punto medio de la lista de puntos del polígono. Repite hasta que solo te queden tres puntos.

C forma de hacer eso:

float areaForPoly(const int numVerts, const Point *verts)
{
    Point v2;
    float area = 0.0f;

    for (int i = 0; i<numVerts; i++){
        v2 = verts[(i + 1) % numVerts];
        area += verts[i].x*v2.y - verts[i].y*v2.x;
    }

    return area / 2.0f;
}

Código Python

Como se describe aquí: http://www.wikihow.com/Calculate-the-Area-of-a-Polygon

Con pandas

import pandas as pd

df = pd.DataFrame({'x': [10, 20, 20, 30, 20, 10, 0], 'y': [-10, -10, -10, 0, 10, 30, 20]})
df = df.append(df.loc[0])

first_product = (df['x'].shift(1) * df['y']).fillna(0).sum()
second_product = (df['y'].shift(1) * df['x']).fillna(0).sum()

(first_product - second_product) / 2
600

Voy a dar algunas funciones simples para calcular el área de un polígono 2d. Esto funciona tanto para polígonos convexos como cóncavos. simplemente dividimos el polígono en muchos sub-triángulos.

//don't forget to include cmath for abs function
struct Point{
  double x;
  double y;
}
// cross_product
double cp(Point a, Point b){ //returns cross product
  return a.x*b.y-a.y*b.x;
}

double area(Point * vertices, int n){  //n is number of sides
  double sum=0.0;
  for(i=0; i<n; i++){
    sum+=cp(vertices[i], vertices[(i+1)%n]); //%n is for last triangle
  }
  return abs(sum)/2.0;
}