algorithm gis geospatial arcgis heatmap

algorithm - Algoritmo de generación de mapa de altura?



gis geospatial (9)

Kriging es uno de los métodos de peso pesado para hacer esto, particularmente dentro del campo de GIS. Tiene varias buenas propiedades matemáticas, la desventaja es que puede ser lenta dependiendo de su variograma .

Si quieres algo más simple, hay muchas rutinas de interpolación que manejan esto bastante bien. Si puede obtener una copia de Recetas numéricas , el Capítulo 3 está dedicado a explicar muchas variantes para la interpolación, e incluye ejemplos de códigos y descripciones de sus propiedades funcionales.

Estaba buscando en Internet y no pude encontrar un algoritmo perfecto para este problema en particular:

Nuestro cliente tiene un conjunto de puntos y datos de peso junto con cada punto, como se puede demostrar con esta imagen:

puntos ponderados http://chakrit.net/files/stackoverflow/so_heightmap_points.png

De los cuales, tenemos un programa SIG que podría generar un "mapa de alturas" o una especie de datos de terreno de estos puntos y sus valores de peso, pero como tenemos cerca de mil puntos de datos y que estos cambiarán con el tiempo, nos gustaría crea nuestras propias herramientas para autogenerar estos mapas de altura.

Hasta ahora, he intentado calcular el peso de cada píxel desde su distancia hasta el punto de datos más cercano con Sqrt((x1 - x2) ^ 2 + (y1 - y2) ^ 2) y aplicando factor de distancia y peso a los datos color de punto para producir el color de degradado resultante para ese píxel en particular:

resultado del mapa de alturas http://chakrit.net/files/stackoverflow/so_heightmap_result.png

Puede ver que todavía hay problemas con cierta configuración de puntos de datos y el algoritmo a veces produce una imagen bastante poligonal cuando hay muchos puntos de datos. El resultado ideal debería verse más como una elipsis y menos como un polígono.

Aquí hay una imagen de ejemplo del artículo de wikipedia sobre gradiente de ascenso que demuestra el resultado que quiero:

montañas http://chakrit.net/files/stackoverflow/so_gradient_descent.png

El algoritmo de ascenso de gradiente no me interesa. Lo que me interesa es el algoritmo para calcular la función original en esa imagen en primer lugar, proporcionó puntos de datos con pesos.

No he tomado ninguna clase en matemáticas topológicas, pero puedo hacer algo de cálculo. Creo que me puede estar perdiendo algo y estoy perdido en lo que debería escribir en el cuadro de búsqueda de Google.

Necesito algunos consejos.

¡Gracias!


es el algoritmo para calcular la función original en esa imagen en primer lugar, proporcionó puntos de datos con pesos.

Es posible. Si comienzas con puntos individuales, siempre terminarás con círculos, pero si ponderas los puntos de datos y lo tienes en cuenta, puedes aplastar los círculos en óvalos como en la imagen.

La razón por la que termina con polígonos es que está utilizando una función discreta en su cálculo: primero busca el color más cercano y luego determina el color.

En su lugar, debe buscar en los algoritmos de gradiente que asigna un color para un punto en función de la distancia y el peso de los tres puntos de datos que encierran ese punto en un triángulo.

Algoritmo de gradiente

Depende de lo que estás tratando de mostrar. Un algoritmo simplista sería:

Para cada píxel:

  • Encuentra los tres puntos que forman el triángulo más pequeño que rodea este píxel
  • Establezca este punto en el color (sistema de color HSV) que se ve afectado tanto por el peso como por la distancia a cada punto de datos:

    pixel.color = datapoint[1].weight * distance(pixel, datapoint[1]) * datapoint[1].color + datapoint[2].weight * distance(pixel, datapoint[2]) * datapoint[2].color + datapoint[3].weight * distance(pixel, datapoint[3]) * datapoint[3].color

Estoy usando + aquí, pero necesita determinar el algoritmo de "promediado" adecuado para su aplicación.

-Adán


Estás buscando algo que Blender llama " metaballs " ( artículo de Wikipedia con enlaces , ejemplo ). Piénsalo de esta manera:

Tus objetos son conos que sobresalen del suelo. Son todas parábolas y el peso dice cuán lejos sobresalen del suelo. Alternativamente, haga que todos tengan la misma altura y ajuste la "llanura" de la parábola en consecuencia, de modo que un gran peso hace que el cono sea muy ancho, mientras que un peso bajo lo hace más afilado. Tal vez incluso ambos hasta cierto punto.

Sugiero que implemente esto y vea cómo se ve.

A continuación, debe colgar un paño o una sábana de goma sobre el resultado. La tela se estirará en una cierta cantidad y generalmente caerá debido a la gravedad. Los conos lo siguen.

Siempre que esté cerca del centro de un cono, la coordenada Z es solo la posición en la superficie del cono. Al salir del centro del cono, la gravedad comienza a disminuir y la influencia de otros conos aumenta.


Este problema no es tan fácil como parece en la superficie. Su problema es que ambos lados del borde de dos regiones deben tener la misma altura, lo que quiere decir que la altura en un píxel determinado está determinada por más que un vecino más cercano.

Si lo entiendo correctamente, necesitas al menos dos algoritmos (y una tercera pieza de jerga).

Para hacer esto correctamente, debes dividir el plano en una teselación de Voronoi .

Probablemente va a querer usar un árbol kd para ayudarlo a encontrar el vecino más cercano. En lugar de tomar O (n ^ 2), esto lo reducirá a O (n log (n)) (el beneficio adicional es que su fase de generación de la región Voronoi será lo suficientemente rápida en desarrollo para trabajar en la fase de cálculo de altura).

Ahora que tiene un mapa en 2-D que indexa cada punto a su vecino i más cercano, necesita caminar a través de cada punto x, y en el mapa y calcular su altura.

Para hacer eso para un punto dado x, y, primero toma su vecino i más cercano y pégalo en una lista, luego recoge todas las regiones contiguas en el diagrama de Voronoi. Una manera fácil es usar el relleno de inundación para encontrar todos los puntos de la región, luego mirar alrededor de la frontera y recoger las otras identidades.

Usando esta lista de todos los vecinos más cercanos, ¡ahora tiene la oportunidad de interpolar correctamente! (Ver otras respuestas para esquemas de interpolación).


La interpolación de superficies parece ser un problema difícil y matemático. Otra forma más económica de hacerlo es:

For each pixel:
For each point:
pixel.addWeight(weight(point, pixel))

def addWeight(w):
totalweight += w
numberofweights += 1
weight = totalweight / numberofweights

Ejemplo de función de peso:

def weight(point, pixel):
return point.weight * 1/(1 + sqrt((point.x - pixel.x)^2 + (point.y - pixel.y)^2))

Es un enfoque de fuerza bruta, pero es simple.


Lo que estás buscando es interpolación de superficie.

Algunos productos existen para hacer esto (aquí hay uno )

La función resultante / spline / otra construcción matemática puede ser interrogada a la resolución requerida para suministrar el mapa de altura.

Su función de interpolación

Sqrt((x1 - x2) ^ 2 + (y1 - y2) ^ 2)

Es similar a los métodos de distancia inversa ponderada , excepto que está aplicando un filtro arbitrario y descartando muchos de los otros puntos de datos.

La mayoría de estas técnicas se basan en un número razonable de muestras y un comportamiento "similar al terreno" que sustenta los valores.

Sugiero usar el peso como la muestra de altura y probar el Método de Shepard en el segundo enlace (no filtrar ningún pixel para comenzar) tomando la proporción de una contribución de puntos de muestra al valor de altura total en un punto de interpolación que puede mezclar los colores de las muestras en esas proporciones para también colorear el punto. Use la intensidad (aproximadamente la escala de grises en el espacio RGB simple) para mostrar la altura o agregar líneas de contorno en negro como lo hace la imagen de ejemplo.


Ha solicitado información sobre algoritmos para la interpolación bidimensional de datos irregulares, que es un área bastante compleja. Como dice que tiene ArcGIS, le recomiendo encarecidamente que interpole automáticamente en ArcGIS utilizando sus características integradas para cálculos automáticos. Estoy seguro de que será mucho más fácil que escribir su propio algoritmo de interpolación. He hecho algo de automatización de ArcGIS, es bastante sencillo.

Si escribe su propio código de interpolación, le aconsejo que no lo haga, lo primero es elegir el algoritmo apropiado, ya que hay varios, cada uno con sus ventajas y desventajas. Aquí hay algunos consejos extraídos de la ayuda para la herramienta de interpolación excelente Surfer (que BTW también se puede automatizar con bastante facilidad). Hay más algoritmos, estos son solo los que he probado.

  • Kriging es uno de los métodos más flexibles y es útil para la creación de grillas de casi cualquier tipo de conjunto de datos. Con la mayoría de los conjuntos de datos, Kriging con el variograma lineal predeterminado es bastante efectivo. En general, recomendamos este método con mucha frecuencia. Kriging es el método de grillado predeterminado porque genera un buen mapa para la mayoría de los conjuntos de datos. Para conjuntos de datos más grandes, Kriging puede ser bastante lento. Kriging puede extrapolar valores de grilla más allá del rango Z de sus datos.
  • La ponderación de distancia inversa es rápida, pero tiene la tendencia a generar patrones de contornos concéntricos "de ojo de buey" alrededor de los puntos de datos. La distancia inversa a una potencia no extrapola valores Z más allá del rango de datos. Un algoritmo de ponderación de distancia inversa simple es fácil de implementar, pero se ralentizará.
  • La triangulación con interpolación lineal es rápida. Cuando utiliza conjuntos de datos pequeños, la triangulación con interpolación lineal genera distintas caras triangulares entre los puntos de datos. La triangulación con interpolación lineal no extrapola valores Z más allá del rango de datos.
  • El Método de Shephard es similar a la Distancia Inversa a una Potencia, pero no tiende a generar patrones de "ojo de buey", especialmente cuando se utiliza un factor de suavizado. El Método de Shepard puede extrapolar valores más allá del rango Z de sus datos.

Para implementar los algoritmos: puede probar Google, o seguir los enlaces en algunas de las otras respuestas. Hay algunos paquetes de GIS de código abierto que incluyen interpolación, por lo que tal vez pueda extraer los algoritmos de ellos si le gusta escalar a través de C ++. O este libro de David Watson es aparentemente considerado un clásico, aunque es una lectura complicada y el código de muestra es espagueti básico !! Pero, por lo que escuché, es el mejor disponible. Si alguien más en sabe mejor, por favor corrígeme ya que tampoco lo puedo creer.


Implementé algo así en Winamp AVS hace un tiempo. Utiliza un enfoque de tipo "metaballs" para calcular la distancia cuadrada inversa (para evitar el sqrt de velocidad) de cada punto de datos, limitándolo (por ejemplo, a 1.0) y sumando esas distancias para cada punto en la cuadrícula 2D. Esto le dará un mapa de color / altura suavemente variable.

Si desea ver el código, está en el preajuste "Glowy" de mi paquete J10 AVS .

EDITAR: Solo mirándolo agregué otro tipo de jazz para que se vea más bonito, la parte más importante es:

d1=s/(sqr(px1-rx)+sqr(py1-ry)); d2=s/(sqr(px2-rx)+sqr(py2-ry)); d3=s/(sqr(px3-rx)+sqr(py3-ry)); d4=s/(sqr(px4-rx)+sqr(py4-ry)); d5=s/(sqr(px5-rx)+sqr(py5-ry)); d6=s/(sqr(px6-rx)+sqr(py6-ry)); d=d1+d2+d3+d4+d5+d6;

Lo cual toma la suma de los 6 puntos. Todo lo demás hecho a los valores de salida rojo, verde y azul es para que se vea más bonito. 6 puntos no es mucho, pero tenga en cuenta que estaba intentando hacer esta ejecución en tiempo real en una cuadrícula de 320x200 en una máquina de 400MHz cuando era nueva (lo que ocurre a ~ 20 fps). :)

Reemplace las líneas roja =, verde = y azul = ... con rojo = d; etc ... para ver a qué me refiero. Toda la belleza desaparece y se queda con una imagen en escala de grises de manchas que varían suavemente alrededor de los puntos de datos.

Otra edición: Olvidé decir que "s" es el peso compartido para todos los puntos, al cambiarlo por cada uno se obtienen pesos individuales para cada punto, por ej. D1 = 2 / (...) y d2 = 1 / (... .) le daría a d1 el doble de altura en su centro que d2. Es posible que también desee limitar la expresión en la parte inferior con algo como d1 = 2 / max (..., 1.0) para suavizar las partes superiores de los puntos para que no lleguen al infinito en el medio. :)

Perdón por el desorden de la respuesta ... Pensé que publicar el ejemplo de código sería lo suficientemente bueno, pero durante la inspección mi código es confuso y difícil de leer. :(


Sé que esta es una pregunta bastante antigua, pero me encontré con ella mientras intentaba resolver un problema similar.

Hay un proyecto de código abierto llamado Surfit que implementa exactamente este tipo de funcionalidad.