Nuevo PostGIS Geometría y Aclaración de Multilínea

Question

Tengo problemas para averiguar si he cargado los datos correctos para PgRouting.

Usando el script bash incluido en PostGIS 2.0, cargué los datos Tiger2010 para el estado de California en EE. UU.

1. La tabla de aristas contiene la columna the_geom, con el tipo de datos geometry. Siguiendo el ejemplo de Underdark, parece que necesito las carreteras en formato Multilinestring para comenzar a generar rutas. ¿Se puede convertir los datos en la tabla de aristas a tipo multilínea? La consulta de la tabla está abajo.
2. Modifiqué mis datos de ubicación geocodificados al tipo de datos geography. La función shortest_path de PgRouting necesita los datos en formato entero. ¿Cómo puedo convertir el tipo point de geography en un entero que shortest_path pueda usar?

Gracias por tu paciencia.

CREATE TABLE tiger.edges
(
  gid integer NOT NULL DEFAULT nextval('edges_gid_seq'::regclass),
  statefp character varying(2),
  countyfp character varying(3),
  tlid numeric(10,0),
  tfidl numeric(10,0),
  tfidr numeric(10,0),
  mtfcc character varying(5),
  fullname character varying(100),
  smid character varying(22),
  lfromadd character varying(12),
  ltoadd character varying(12),
  rfromadd character varying(12),
  rtoadd character varying(12),
  zipl character varying(5),
  zipr character varying(5),
  featcat character varying(1),
  hydroflg character varying(1),
  railflg character varying(1),
  roadflg character varying(1),
  olfflg character varying(1),
  passflg character varying(1),
  divroad character varying(1),
  exttyp character varying(1),
  ttyp character varying(1),
  deckedroad character varying(1),
  artpath character varying(1),
  persist character varying(1),
  gcseflg character varying(1),
  offsetl character varying(1),
  offsetr character varying(1),
  tnidf numeric(10,0),
  tnidt numeric(10,0),
  the_geom geometry,
  CONSTRAINT edges_pkey PRIMARY KEY (gid),
  CONSTRAINT enforce_dims_the_geom CHECK (st_ndims(the_geom) = 2),
  CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'MULTILINESTRING'::text OR the_geom IS NULL),
  CONSTRAINT enforce_srid_the_geom CHECK (st_srid(the_geom) = 4269)
)
WITH (
  OIDS=FALSE
);
ALTER TABLE tiger.edges OWNER TO postgres;
GRANT ALL ON TABLE tiger.edges TO postgres;
GRANT ALL ON TABLE tiger.edges TO gis_group;

-- Índice: tiger.idx_edges_tlid

-- DROP INDEX tiger.idx_edges_tlid;

CREATE INDEX idx_edges_tlid
  ON tiger.edges
  USING btree
  (tlid);

-- Índice: tiger.idx_tiger_edges_countyfp

-- DROP INDEX tiger.idx_tiger_edges_countyfp;

CREATE INDEX idx_tiger_edges_countyfp
  ON tiger.edges
  USING btree
  (countyfp);

-- Índice: tiger.idx_tiger_edges_tfidl

-- DROP INDEX tiger.idx_tiger_edges_tfidl;

CREATE INDEX idx_tiger_edges_tfidl
  ON tiger.edges
  USING btree
  (tfidl);

-- Índice: tiger.idx_tiger_edges_tfidr

-- DROP INDEX tiger.idx_tiger_edges_tfidr;

CREATE INDEX idx_tiger_edges_tfidr
  ON tiger.edges
  USING btree
  (tfidr);

-- Índice: tiger.tiger_edges_the_geom_gist

-- DROP INDEX tiger.tiger_edges_the_geom_gist;

CREATE INDEX tiger_edges_the_geom_gist
  ON tiger.edges
  USING gist
  (the_geom);

Answer 1

Creo que puede haber tenido un problema similar. El assign_vertex_id() de pgRouting en la versión 1.05 está un poco desactualizado. Para ver si se solucionó en una versión posterior, descargué la fuente más reciente del repositorio Git (https://github.com/pgRouting/pgrouting). El problema del prefijo "ST_" se había resuelto en routing_topology.sql.

Desafortunadamente, no se puede ejecutar assign_vertex_id() directamente en geometrías MULTILINESTRING porque ST_StartPoint() y ST_EndPoint() devolverán valores NULL. No creo que esto siempre haya sido cierto, pero hasta donde yo sé, es un problema actual.

SELECT ST_StartPoint('MULTILINESTRING((-29 -27,-30 -29.7,-36 -31,-45 -33),(-45.2 -33.2,-46 -32))'::geometry) IS NULL AS is_null;

 is_null
---------
 t
(1 row)

En mi caso particular, pude modificar mis datos y eso pareció solucionar el problema.

Answer 2

A continuación se presenta una versión modificada del script "routing_topology.sql". Proporciona una solución para dos de los problemas mencionados anteriormente:

• Tus tablas ya no tienen que estar en el esquema "public".
• Si utilizas PostGIS 2.0 y las funciones ST_StartPoint() y ST_EndPoint(), devolverán NULL si tu geometría es de tipo "MULTILINESTRING". Eso ahora está solucionado, por lo que incluso con un punto de inicio y fin MULTILINESTRING se devolverán correctamente.

Reemplaza el contenido de tu archivo "routing_topology.sql" con el código a continuación y ejecuta el script nuevamente.

-----------------------------------------------------------------------
-- Esta función no debe ser utilizada directamente. En su lugar, usa assign_vertex_id
-- 
-- Inserta un punto en una tabla de vértices temporal y devuelve un id
--  de un nuevo punto o un punto existente. La tolerancia es la distancia mínima
--  entre puntos existentes y el nuevo punto para crear un nuevo punto.
--
-- Últimos cambios: 16.04.2008
-- Autor: Christian Gonzalez
-----------------------------------------------------------------------
CREATE OR REPLACE FUNCTION point_to_id(p geometry, tolerance double precision)
RETURNS BIGINT 
AS 
$$ 

DECLARE
    _r record; 
    _id bigint; 
    _srid integer;

BEGIN

    BEGIN   
        _srid := Find_SRID(user::text, 'vertices_tmp', 'the_geom');
    EXCEPTION
        WHEN SQLSTATE 'P0001' THEN
            _srid := Find_SRID('public', 'vertices_tmp', 'the_geom');
    END;

    SELECT

        ST_Distance(the_geom, p) AS d, id, the_geom

    INTO _r FROM vertices_tmp WHERE

        ST_DWithin(the_geom, p, tolerance)

    ORDER BY d LIMIT 1; IF FOUND THEN

        _id:= _r.id;

    ELSE

        INSERT INTO vertices_tmp(the_geom) VALUES (p); _id:=lastval();

    END IF;

    RETURN _id;

END; $$ LANGUAGE 'plpgsql' VOLATILE STRICT; 

-----------------------------------------------------------------------
-- Llena la columna source y target_id para todas las líneas. Todos los extremos de línea
--  con una distancia menor que la tolerancia, se les asigna el mismo id
--
-- Últimos cambios: 16.04.2008
-- Autor: Christian Gonzalez
-----------------------------------------------------------------------

CREATE OR REPLACE FUNCTION assign_vertex_id(geom_table varchar, tolerance double precision, geo_cname varchar, gid_cname varchar)
RETURNS VARCHAR AS
$$
DECLARE
    _r record;
    source_id int;
    target_id int;
    srid integer;
BEGIN

    BEGIN
    DROP TABLE vertices_tmp;
    EXCEPTION 
    WHEN UNDEFINED_TABLE THEN
    END;

    EXECUTE 'CREATE TABLE vertices_tmp (id serial)';

    BEGIN   
        srid := Find_SRID(user::text, quote_ident(geom_table), quote_ident(geo_cname));
    EXCEPTION
        WHEN SQLSTATE 'P0001' THEN
            srid := Find_SRID('public', quote_ident(geom_table), quote_ident(geo_cname));
    END;

    EXECUTE 'SELECT addGeometryColumn(''vertices_tmp'', ''the_geom'', '||srid||', ''POINT'', 2)';
    CREATE INDEX vertices_tmp_idx ON vertices_tmp USING GIST (the_geom);

    FOR _r IN EXECUTE 'SELECT ' || quote_ident(gid_cname) || ' AS id,'
        || ' ST_StartPoint(ST_GeometryN('|| quote_ident(geo_cname) ||', 1)) AS source,'
            || ' ST_EndPoint(ST_GeometryN('|| quote_ident(geo_cname) ||', 1)) as target'
        || ' FROM ' || quote_ident(geom_table) || ' WHERE ' || quote_ident(geo_cname) || ' IS NOT NULL '
    LOOP

        source_id := point_to_id(_r.source, tolerance);
    target_id := point_to_id(_r.target, tolerance);

    EXECUTE 'update ' || quote_ident(geom_table) || 
        ' SET source = ' || source_id || 
        ', target = ' || target_id || 
        ' WHERE ' || quote_ident(gid_cname) || ' =  ' || _r.id;
    END LOOP;

    RETURN 'OK';

END;
$$
LANGUAGE 'plpgsql' VOLATILE STRICT;

Answer 3

Tus geometrías ya son MULTILINESTRING. De lo contrario, no cumplirían con la siguiente restricción:

CONSTRAINT enforce_geotype_the_geom VERIFICA (geometrytype(the_geom) = 'MULTILINESTRING'::text O the_geom ES NULL),

Para el enrutamiento con shortest_path() necesitas IDs de nodo no IDs de borde. Para generar IDs de nodo, utiliza la función assign_vertex_id().

assign_vertex_id() completará ciertas columnas en tu tabla de bordes, pero primero debes crearlas. Desde http://www.pgrouting.org/docs/foss4g2008/ch06.html:

ALTER TABLE ways AGREGAR COLUMNA source integer;
ALTER TABLE ways AGREGAR COLUMNA target integer;
SELECT assign_vertex_id('ways', 0.00001, 'the_geom', 'gid');