ST_DWithin fue más rápido en mi prueba que ST_Intersects. Esto es sorprendente, sobre todo porque se supone que el algoritmo de geometría preparada entra en acción en casos como éste. Creo que existe la posibilidad de que esto sea mucho más rápido de lo que he mostrado aquí.
Hice más pruebas y dos cosas casi duplicaron la velocidad. En primer lugar, he probado en un ordenador más nuevo, pero todavía un portátil bastante ordinario, tal vez a excepción de SATA3 ssd -discos.
Entonces, la consulta que aparece a continuación tardó 18 segundos en lugar de 62 segundos en el antiguo portátil. Luego descubrí que estaba totalmente equivocado cuando escribí que el índice en la tabla de puntos no era necesario. Con ese índice en su lugar ST_Intersects se comportó como se esperaba y las cosas se volvieron muy rápidas. Aumenté el número de puntos en la tabla de puntos a 1 millón de puntos y la consulta:
CREATE TABLE points_ct AS
SELECT imported_ct.gid as ct_gid, t.gid as point_id
FROM imported_ct , t WHERE ST_Intersects(imported_ct.geom , t.geom);
se ejecuta en 72 segundos. Como hay 1249 polígonos, se hacen 1249000000 pruebas en 72 segundos. Eso hace unas 17000000 pruebas por segundo. O bien, probar casi 14000 puntos contra todos los polígonos por segundo.
A partir de esta prueba tus 400000000 puntos a probar deberían tardar unas 8 horas sin ningún problema con la distribución de la carga a varios núcleos. PostGIS nunca deja de impresionarme :-)
En primer lugar, para visualizar el resultado puedes añadir la geometría de puntos a la tabla resultante, abrirla en QGIS por ejemplo y darle estilo con valores únicos en el campo imported_ct.
En segundo lugar, sí, puedes obtener también los puntos que caen fuera de cualquier polígono utilizando una unión a la derecha (o a la izquierda) como esta:
CREATE TABLE points_ct AS
SELECT imported_ct.gid as ct_gid, t.gid as point_id
FROM imported_ct right join t ON ST_Intersects(imported_ct.the_geom , t.geom);
Hice algunas pruebas para verificar si parece posible PostGIS.
Primero algo que no entiendo. Tienes dos puntos por fila. ¿Están siempre los dos puntos en el mismo polígono? Entonces basta con hacer los cálculos sobre uno de los puntos. Si pueden estar en dos polígonos diferentes necesitarás una forma de conectar una fila de puntos a dos polígonos.
Por las pruebas parece que se puede hacer, pero puede que necesites alguna solución creativa para repartir la carga entre más de un núcleo de la cpu.
Lo he probado en un portátil de 4 años con cpu centrino de doble núcleo (unos 2,2GHz creo) , 2GB de RAM. Si tienes 48 BG de RAM supongo que también tienes mucha más potencia de cpu.
Lo que hice fue crear una tabla de puntos aleatorios con 100000 puntos así:
CREATE TABLE t AS
WITH r AS
(SELECT ST_Extent(the_geom)::geometry ext FROM imported_ct)
SELECT ST_Point(x,y) AS geom FROM
(SELECT GENERATE_SERIES(1,100000)) s,
(SELECT ST_Xmin(ext)+(random()*(ST_Xmax(ext)-ST_Xmin(ext))) x, ST_Ymin(ext)+(random()*(ST_Ymax(ext)-ST_Ymin(ext))) y FROM r
) f;
Entonces añadiendo un gid como:
ALTER TABLE t ADD COLUMN GID SERIAL;
Entonces, corriendo :
CREATE TABLE points_ct AS
SELECT imported_ct.gid as ct_gid, t.gid as point_id FROM imported_ct , t WHERE ST_Dwithin(imported_ct.the_geom , t.geom,0);
tarda unos 62 segundos (compárese con el resultado de ArcGIS con la misma cantidad de puntos). El resultado es una tabla que conecta los puntos de mi tabla t con el gid de la tabla con tramo censal.
Con esa velocidad harás 200 millones de puntos en unas 34 horas. Así que, si es suficiente con comprobar uno de los puntos, mi viejo portátil puede hacerlo con un solo núcleo.
Pero si tiene que comprobar ambos puntos puede ser más difícil.
Entonces puedes distribuir manualmente la carga a más de un núcleo iniciando varias sesiones contra la base de datos y ejecutando diferentes consultas.
En mi ejemplo con 50000 puntos y dos núcleos de cpu he probado :
CREATE TABLE t1 as
SELECT imported_ct.gid as ct_gid, t.gid as point_id FROM imported_ct , t WHERE t.gid >50000 and ST_Dwithin(imported_ct.the_geom , t.geom,0);
en una sesión db al mismo tiempo que se ejecuta:
CREATE TABLE t2 as
SELECT imported_ct.gid as ct_gid, t.gid as point_id FROM imported_ct , t WHERE t.gid <=50000 and ST_Dwithin(imported_ct.the_geom , t.geom,0);
en otra sesión db.
Eso tomó alrededor de 36 segundos por lo que es un poco más lento que el primer ejemplo, probablemente en función de la escritura del disco al mismo tiempo. Pero como los dos núcleos están trabajando al mismo tiempo, no me llevó más de 36 segundos.
Para la unión de la tabla t1 y t2 se intentó:
CREATE TABLE t3 AS
SELECT * FROM t1
UNION ALL
SELECT * FROM t2;
utilizando alrededor de medio segundo.
Así que, con un hardware más fresco y distribuyendo la carga entre muchos núcleos, esto debería ser absolutamente posible aunque el mundo real sea más lento que el caso de prueba.
Cabe destacar que el ejemplo es de Linux (Ubuntu). Usar Windows será otra historia. Pero tengo todas las demás aplicaciones diarias en funcionamiento por lo que el portátil está bastante cargado de antes. Así que podría simular el caso de Windows bastante bien tal vez, sin abrir nada más que pgadmin.
