16 votos

Intersección de una trama con un polígono de usar PostGIS - artefacto de error

Estoy usando PostGIS2.0 para hacer algo de trama/polígono de las intersecciones. Estoy teniendo dificultades para entender qué operación debo usar, y lo que la forma más rápida de realizar esto es. Mi problema es el siguiente:

  • Tengo un polígono y una trama
  • Quiero encontrar a todos los píxeles que están dentro del polígono, y obtener la suma del valor de píxel
  • Y (actualizado problema): estoy recibiendo masiva de los valores de los píxeles que no existen en el original de la trama, cuando realizo la consulta

Estoy teniendo dificultades para entender si debo utilizar ST_Intersects() o ST_Intersection(). Yo tampoco sé cuál es el mejor enfoque para sumar mi píxeles. Aquí es el primer acercamiento que he probado (#1):

SELECT 
    r.rast 
FROM
    raster as r, 
    polygon as p
WHERE
    ST_Intersects(r.rast, p.geom)

Esto devuelve una lista de rast valores, que no estoy seguro de qué hacer con ellos. Traté de calcular el resumen de las estadísticas de uso de ST_SummaryStats() pero no estoy seguro de si este es el promedio ponderado de la suma de todos los píxeles que se encuentran dentro del polígono.

SELECT  
        (result).count,
        (result).sum    
FROM (  
         SELECT 
            ST_SummaryStats(r.rast) As result
         FROM
            raster As r, 
            polygon As p
         WHERE
            ST_Intersects(r.rast, p.geom)    
    ) As tmp

El otro enfoque que he probado (#2) utiliza ST_Intersection():

 SELECT
        (gv).geom,         
        (gv).val
 FROM 
 (
    SELECT 
        ST_Intersection(r.rast, p.geom) AS gv
    FROM 
        raster as r, 
        polygon as p           
    WHERE 
        ST_Intersects(r.rast, p.geom)

      ) as foo;

Esto devuelve una lista de geometrías que analizo más, pero supongo que esto es menos eficiente.

Estoy claro en que es el más rápido, el orden de operación también. Siempre debo elegir raster, polygon o polygon, raster, o convertir el polígono en una trama por lo que es raster, raster ?

EDIT: he actualizado enfoque #2 con algunos detalles de R.K.'s de enlace.

El uso de enfoque #2, me he dado cuenta el siguiente error en los resultados, que es parte de la razón por la que no entiendo la salida. Aquí está la imagen original de mi trama, y un contorno del polígono que está siendo utilizado a lo cruzan, se superponen en la parte superior:

enter image description here

Y aquí está el resultado de la intersección con PostGIS:

enter image description here

El problema con la consecuencia es que no hay valores de 21474836 de ser devuelto, que no están en el original de la trama. Yo no tengo ni idea de por qué está sucediendo esto. Yo sospecho que es relativa a un pequeño número en algún lugar (dividiendo por casi 0), pero devuelve el resultado es incorrecto.

7voto

shsteimer Puntos 8749

He encontrado un tutorial en la intersección del vector de polígonos con una gran trama de cobertura utilizando PostGIS WKT Trama. Podría tener las respuestas que estás buscando.

El tutorial utiliza dos conjuntos de datos, un punto archivo de forma que se búfer para producir polígonos y una serie de 13 de elevación SRTM rásteres. Había un montón de pasos, pero la consulta que se utiliza para cruzar el raster y vector se veía así:

 CREATE TABLE caribou_srtm_inter AS
 SELECT id, 
        (gv).geom AS the_geom, 
        (gv).val
 FROM (SELECT id, 
              ST_Intersection(rast, the_geom) AS gv
       FROM srtm_tiled,
            cariboupoint_buffers_wgs
       WHERE ST_Intersects(rast, the_geom)
      ) foo;

Los valores se resumen de la siguiente manera:

 CREATE TABLE result01 AS
 SELECT id, 
        sum(ST_Area(ST_Transform(the_geom, 32198)) * val) / 
        sum(ST_Area(ST_Transform(the_geom, 32198))) AS meanelev
 FROM caribou_srtm_inter
 GROUP BY id
 ORDER BY id;

Yo realmente no sé lo suficiente de PostGIS para explicar esto, pero seguro que suena como lo que estaban tratando de lograr. El tutorial debe arrojar luz sobre los pasos intermedios. Buena suerte :)

3voto

aalaap Puntos 216

Con respecto al punto 2 en la pregunta original - varios de los postgis 2.0 versiones en desarrollo utiliza una versión de la librería GDAL que el reparto de los flotadores a enteros. Si la trama tiene valores de coma flotante en ella, y que estaba utilizando una versión de GDAL inferior 1.9.0, o una versión de PostGIS 2.0 preliminar que no llaman GDALFPolygonize(), entonces usted puede ser que se produzca este error. Entradas en la PostGIS y GDAL los gestores de fallos fueron presentadas y cerrado. Este error se activa en el momento de la pregunta original.

En términos de rendimiento, encontrará que el uso de ST_Intersects(raster, geom) es mucho, mucho más rápido que usando ST_Intersects(geom, raster). La primera versión, se rasteriza la geometría y hace un mapa de bits de espacio de intersección. La segunda versión vectorizes la geometría) y tiene un vector de espacio de intersección, que puede ser mucho más caro el proceso.

0voto

jczaplew Puntos 312

Yo también tenía extraño con los problemas de ST_SummaryStats con ST_Clip. La consulta de los datos de manera diferente, me dijo que el valor min de mi trama era de 32, y, a continuación, max 300, sin embargo, ST_SummaryStats regresaba -32700 para los valores de píxel en mi polígono de destino.

Terminé de hacking en torno a la cuestión así:

WITH first AS (
   SELECT id, (ST_Intersection(geom, rast)).val
   FROM raster_table
   INNER JOIN vector_table ON ST_Intersects(rast, geom)
)
SELECT id, COUNT(val), SUM(val), AVG(val), stddev(val), MIN(val), MAX(val)
FROM first
GROUP BY id

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X