diff --git a/scripts_v2/phase_3_compute.py b/scripts_v2/phase_3_compute.py index b678552..f21a1e8 100644 --- a/scripts_v2/phase_3_compute.py +++ b/scripts_v2/phase_3_compute.py @@ -69,6 +69,77 @@ def closeConnRedadegDB(): # ============================================================================== +def getPKfromRouting(start, distance): + + sql_routage = """ +WITH t AS ( +SELECT * +FROM pgr_drivingDistance('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr +WHERE SOURCE IS NOT NULL AND id > 0', +"""+str(start)+""", """+str(distance)+""") +) +SELECT node, edge, round(agg_cost) FROM t ORDER BY seq DESC LIMIT 1;""" + + #print(sql_routage) + + cursor = db_redadeg_pg_conn.cursor() + cursor.execute(sql_routage) + + # on récupère l'id du nœud de fin + data = cursor.fetchone() + node_end = data[0] + edge = data[1] + distance = data[2] + cursor.close() + + return([node_end, edge, distance]) + + +# ============================================================================== + +def getPgrNodeInfos(node_id): + + # cette fonction va chercher les infos où il faut pour le PK + cursor = db_redadeg_pg_conn.cursor() + + # géométrie… + sql_get_from_pgr_node = """ +SELECT + the_geom, + TRUNC(ST_X(the_geom)::numeric,1) AS x, + TRUNC(ST_Y(the_geom)::numeric,1) AS y, + TRUNC(ST_X(ST_Transform(the_geom,4326)::geometry(Point, 4326))::numeric,8) AS long, + TRUNC(ST_Y(ST_Transform(the_geom,4326)::geometry(Point, 4326))::numeric,8) AS lat +FROM phase_3_troncons_pgr_vertices_pgr v WHERE id = """+ str(node_id) +""";""" + #print(sql_get_from_pgr_node) + + cursor.execute(sql_get_from_pgr_node) + data = cursor.fetchone() + the_geom = data[0] + x = data[1] + y = data[2] + long = data[3] + lat = data[4] + + cursor.close() + + return([the_geom,x,y,long,lat]) + + +# ============================================================================== + +def getLongueurParcourue(node_start, node_end): + + # cette fonction sert à retourner la longueur parcourue entre 2 nœuds + sql_get_longueur = "SELECT round(max(agg_cost))::integer " + sql_get_longueur += "FROM pgr_dijkstra('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr " + sql_get_longueur += " WHERE SOURCE IS NOT NULL', " + str(node_start) + ", " + str(node_end) + ")" + + cursor = db_redadeg_pg_conn.cursor() + cursor.execute(sql_get_longueur) + data = cursor.fetchone()[0] + + return (data) # ============================================================================== @@ -152,8 +223,11 @@ try: # ------------------------------------------------------ print(" Récupération des infos du secteur") + + # PK de départ de la Redadeg (à cause de l'avant course) + pk_start_redadeg = config.get('redadeg', 'pk_start') - sql_get_infos_secteur = "SELECT pk_start, pk_stop, node_start, node_stop, longueur " + sql_get_infos_secteur = "SELECT pk_start, node_start, node_stop, longueur, longueur_km_redadeg " sql_get_infos_secteur += "FROM secteur " sql_get_infos_secteur += "WHERE id = "+secteur+" ;" @@ -161,69 +235,136 @@ try: infos_secteur = db_redadeg_cursor.fetchone() secteur_pk_start = infos_secteur[0] - secteur_pk_stop = infos_secteur[1] - secteur_node_start = infos_secteur[2] - secteur_node_stop = infos_secteur[3] - secteur_longueur = infos_secteur[4] + secteur_node_start = infos_secteur[1] + secteur_node_stop = infos_secteur[2] + secteur_longueur = infos_secteur[3] + secteur_longueur_km_redadeg = infos_secteur[4] # pour test #secteur_longueur = 10000 + # longueur de découpage des tronçons de la phase 2 + longueur_densification = config.get('redadeg', 'longueur_densification') + print(" fait") print("") # on détermine le nb théorique de PK pour ce secteur - secteur_nb_pk = secteur_pk_stop - secteur_pk_start - # et ainsi la longueur réelle entre chaque PK - longueur_decoupage = int(secteur_longueur / secteur_nb_pk) + secteur_nb_pk = int(secteur_longueur / secteur_longueur_km_redadeg) + # et la longueur réelle demandée au calculateur + longueur_decoupage = int(secteur_longueur_km_redadeg) - print(" " + str(secteur_nb_pk) + " KM redadeg de " + str(longueur_decoupage) + " m vont être créés") - print(" pour une longeur réelle de " + '{:,}'.format(secteur_longueur).replace(',', ' ') + " m") + print(" " + str(secteur_nb_pk) + " KM redadeg de " + str(secteur_longueur_km_redadeg) + " m vont être créés") + print(" pour une longeur réelle de " + str(secteur_longueur) + " m") print("") - # ------------------------------------------------------ - print(" Suppression des PK du secteur") + # cette variable pour stocker la requête SQL de création des PK + # et la première requête sera de supprimer les données du secteur + sql_insert_pks = "DELETE FROM phase_3_pk WHERE secteur_id = "+secteur+" ;\n" - sql_delete_pks = f"DELETE FROM phase_3_pk WHERE secteur_id = {secteur} ;" - db_redadeg_cursor.execute(sql_delete_pks) + # ------------------------------------------------------ + print(" Calcul du 1er PK du secteur") + + # on a les infos -> on calcule la route qui va du 1er nœud de départ et qui fait la distance demandée + # pour récupérer l'id du noeud de fin qui va devenir notre PK + node_zero = secteur_node_start + node_zero_data = getPgrNodeInfos(node_zero) + + sql_insert_pks += "INSERT INTO phase_3_pk (secteur_id, pk_id, the_geom, pk_x, pk_y, pk_long, pk_lat, length_theorical, length_real) VALUES (" + sql_insert_pks += secteur + ", " + str(secteur_pk_start) + sql_insert_pks += ",'" + node_zero_data[0] + "'" + sql_insert_pks += "," + str(node_zero_data[1]) + "," + str(node_zero_data[2]) + sql_insert_pks += "," + str(node_zero_data[3]) + "," + str(node_zero_data[4]) + sql_insert_pks += f",{longueur_decoupage},0" + sql_insert_pks += ");\n" + + print(" nœud du PK " + str(secteur_pk_start) + " : " + str(node_zero)) + print("") + + # ------------------------------------------------------ + print(" Calcul des autres PK") + print("") + + # maintenant on peut itérer jusqu'à la fin du secteur + node_x = node_zero + + # en sa basant sur la longueur des PK posés et la longueur totale du secteur + longueur_parcourue = getLongueurParcourue(node_zero, node_x) + if longueur_parcourue is None: longueur_parcourue = 0 + longueur_restante = secteur_longueur - longueur_parcourue + + # un compteur pour la boucle + i = 1 + # début de l'id des PK qui commence avec le PK de début du secteur + pk_id = secteur_pk_start + + # tant que la distance restante est supérieure à la distance de découpage + # on boucle + while longueur_restante >= longueur_decoupage: + + # incrément du compteur de la boucle + i += 1 + # incrément du compteur de PK + pk_id += 1 + + # on va trouver le prochain PK + pk_data = getPKfromRouting(node_x , longueur_decoupage) + node_x = pk_data[0] + previous_pk_edge = pk_data[1] + longueur_km_redadeg = pk_data[2] + longueur_parcourue = getLongueurParcourue(node_zero,node_x) + longueur_restante = secteur_longueur - longueur_parcourue + + #print(" nouveau nœud : " + str(node_x)) + #print(" previous_pk_edge : "+ str(previous_pk_edge)) + + # on sort une infos pour suivre si le traitement bosse + if (i <= 5) or (i % 10 == 0) or (i >= secteur_nb_pk - 5): + print(" PK " + str(pk_id)) + print(" id du nœud : " + str(node_x)) + print(" " + str(longueur_parcourue) + " m jusqu'à maintenant") + print(" " + str(longueur_restante) + " m restant jusqu'à la fin du secteur") + + + # ici on construit la requête avec les données du PK + node_x_data = getPgrNodeInfos(node_x) + + # on fait une requête SQL d'insert de ce PK + sql_insert_pks += "INSERT INTO phase_3_pk (secteur_id, pk_id, the_geom, pk_x, pk_y, pk_long, pk_lat, length_theorical, length_real) VALUES (" + sql_insert_pks += secteur + "," + str(pk_id) + sql_insert_pks += ",'" + node_x_data[0] + "'" + sql_insert_pks += "," + str(node_x_data[1]) + "," + str(node_x_data[2]) + sql_insert_pks += "," + str(node_x_data[3]) + "," + str(node_x_data[4]) + sql_insert_pks += f",{longueur_decoupage},{longueur_km_redadeg}" + sql_insert_pks += ");\n" + + # on met en négatif l'info de routage du précédent tronçon afin de l'écarter du prochain calcul de routage + sql_neutralisation = "UPDATE phase_3_troncons_pgr SET id = -ABS("+str(previous_pk_edge)+") WHERE id = "+str(previous_pk_edge)+" ;" + #print(sql_neutralisation) + db_redadeg_cursor.execute(sql_neutralisation) + + print("") + print(" Fin de la boucle") + print("") + + print(" RAZ de la neutralisation des infos de routage pour la boucle") + sql_reset_neutralisation = "UPDATE phase_3_troncons_pgr SET id = -1*id WHERE id < 0 ;" + db_redadeg_cursor.execute(sql_reset_neutralisation) print(" fait") print("") - - # ------------------------------------------------------ - print(" Création des nouveaux PK du secteur tous les " + str(longueur_decoupage) + " m") - - sql_generate_pks = f""" -WITH linemeasure AS ( - WITH line AS ( - -- on récupère un itinéraire calculée par pgRouting - SELECT ST_Union(the_geom) AS the_geom - FROM pgr_dijkstra('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr', - {secteur_node_start},{secteur_node_stop}) a - JOIN phase_3_troncons_pgr b ON a.edge = b.id - ) -SELECT - generate_series(0, (ST_Length(line.the_geom))::int, {longueur_decoupage}) AS i, - ST_AddMeasure(the_geom,0,ST_length(the_geom)) AS the_geom -FROM line -) -INSERT INTO phase_3_pk (pk_id,secteur_id,length_real,length_total,the_geom) -SELECT - ROW_NUMBER() OVER() + ({secteur_pk_start}-1) - ,{secteur} - ,{longueur_decoupage} - ,i - ,(ST_Dump(ST_GeometryN(ST_LocateAlong(the_geom, i), 1))).geom AS the_geom -FROM linemeasure ; -""" - - db_redadeg_cursor.execute(sql_generate_pks) + print(" Écriture des PK dans la couche phase_3_pk") + db_redadeg_cursor.execute(sql_insert_pks) print(" fait") print("") - - - + print(" sauvegarde du dernier PK calculé pour ce secteur") + # on est sorti de la boucle alors on va écrire en base l'id du dernier PK calculé + sql_update_pk_end = "UPDATE secteur SET pk_stop = " + str(pk_id) + sql_update_pk_end += "WHERE id = " + secteur + " ;" + db_redadeg_cursor.execute(sql_update_pk_end) + print(" fait") + print("") db_redadeg_cursor.close() diff --git a/scripts_v2/phase_3_prepare.py b/scripts_v2/phase_3_prepare.py index 7207f64..a9e6414 100644 --- a/scripts_v2/phase_3_prepare.py +++ b/scripts_v2/phase_3_prepare.py @@ -159,14 +159,31 @@ try: # ------------------------------------------------------ - print(" Chargement des tronçons de la phase 2") + print(" Chargement de tronçons découpés tous les "+longueur_densification+" m depuis la couche des tronçons phase 2") + + # on charge, pour le secteur concerné des tronçons courts découpés tous les x mètres + # (densification avec ST_LineSubstring ) - sql_insert = f""" -INSERT INTO phase_3_troncons_pgr (secteur_id, path_seq, osm_id, highway, type, oneway, ref, name_fr, name_br, the_geom) + sql_insert = """ +INSERT INTO phase_3_troncons_pgr (secteur_id, osm_id, highway, type, oneway, ref, name_fr, name_br, the_geom) + SELECT + secteur_id, osm_id, highway, type, oneway, ref, name_fr, name_br, + ST_LineSubstring(the_geom, """+longueur_densification+"""*n/length, + CASE + WHEN """+longueur_densification+"""*(n+1) < length THEN """+longueur_densification+"""*(n+1)/length + ELSE 1 + END) As the_geom +FROM + ( SELECT - secteur_id, path_seq, osm_id, highway, type, oneway, ref, name_fr, name_br, the_geom - FROM phase_2_trace_pgr -WHERE secteur_id = {secteur} ;""" + secteur_id, osm_id, highway, type, oneway, ref, name_fr, name_br, + ST_Length(the_geom) AS length, + the_geom + FROM phase_2_trace_troncons + WHERE secteur_id = """+secteur+""" + ) AS t +CROSS JOIN generate_series(0,10000) AS n +WHERE n*"""+longueur_densification+"""/length < 1;""" db_redadeg_cursor.execute(sql_insert) @@ -177,12 +194,20 @@ WHERE secteur_id = {secteur} ;""" print(" Calcul des coûts...") # calcul des attributs de support du calcul pour PGR - sql_update_costs = f""" + sql_update_costs = """ UPDATE phase_3_troncons_pgr SET -cost = trunc(st_length(the_geom)::numeric,2), -reverse_cost = trunc(st_length(the_geom)::numeric,2) -WHERE secteur_id = {secteur} ;""" +cost = + CASE + WHEN trunc(st_length(the_geom)::numeric,2) = """ + str( float(longueur_densification) - 0.01 ) + """ THEN """ + longueur_densification + """ + ELSE trunc(st_length(the_geom)::numeric,2) + END, +reverse_cost = + CASE + WHEN trunc(st_length(the_geom)::numeric,2) = """ + str( float(longueur_densification) - 0.01 ) + """ THEN """ + longueur_densification + """ + ELSE trunc(st_length(the_geom)::numeric,2) + END +WHERE secteur_id = """ + secteur + """ ;""" db_redadeg_cursor.execute(sql_update_costs) diff --git a/scripts_v2/sql/create_tables.sql b/scripts_v2/sql/create_tables.sql index 33dd841..d52b594 100644 --- a/scripts_v2/sql/create_tables.sql +++ b/scripts_v2/sql/create_tables.sql @@ -507,15 +507,31 @@ ALTER TABLE phase_2_tdb OWNER TO redadeg; ========================================================================== */ +-- cette table contient les infos de découpage des tronçons +DROP TABLE IF EXISTS phase_3_secteurs CASCADE ; +CREATE TABLE phase_3_secteurs +( + secteur_id integer, + nom_br text, + nom_fr text, + longueur_km_redadeg integer, + node_start integer, + node_stop integer, + pk_start integer, + pk_stop integer, + CONSTRAINT phase_3_secteurs_pkey PRIMARY KEY (secteur_id) +); +ALTER TABLE phase_3_secteurs OWNER TO redadeg; + -- la table qui va accueillir une couche support de calcul itinéraire phase 3 +-- à savoir les tronçons phase 2 découpés tous les x mètres DROP TABLE IF EXISTS phase_3_troncons_pgr CASCADE ; CREATE TABLE phase_3_troncons_pgr ( secteur_id integer, -- info de routage id serial, - path_seq bigint, source bigint, target bigint, cost double precision, @@ -536,7 +552,6 @@ CREATE INDEX phase_3_troncons_pgr_geom_idx ON phase_3_troncons_pgr USING gist(th ALTER TABLE phase_3_troncons_pgr OWNER to redadeg; -/* DROP TABLE IF EXISTS phase_3_troncons CASCADE ; CREATE TABLE phase_3_troncons ( @@ -560,7 +575,7 @@ CREATE VIEW phase_3_troncons_4326 AS ST_Transform(the_geom,4326) AS the_geom FROM phase_3_troncons ; ALTER TABLE phase_3_troncons_4326 OWNER TO redadeg; -*/ + DROP TABLE IF EXISTS phase_3_trace_secteurs CASCADE ; @@ -602,7 +617,6 @@ CREATE TABLE phase_3_pk pk_lat numeric(10,8), length_real integer, length_theorical integer, - length_total integer, secteur_id integer, municipality_admincode text, municipality_postcode text, @@ -617,7 +631,7 @@ CREATE TABLE phase_3_pk way_name_br text, the_geom geometry, CONSTRAINT phase_3_pk_pkey PRIMARY KEY (pk_id), - --CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text), + CONSTRAINT enforce_geotype_the_geom CHECK (geometrytype(the_geom) = 'POINT'::text), CONSTRAINT enforce_srid_the_geom CHECK (st_srid(the_geom) = 2154) ) ; CREATE INDEX phase_3_pk_geom_idx ON phase_3_pk USING gist(the_geom); diff --git a/scripts_v2/sql/phase_3_brouillon.sql b/scripts_v2/sql/phase_3_brouillon.sql index df6a5ea..6bbdfe5 100644 --- a/scripts_v2/sql/phase_3_brouillon.sql +++ b/scripts_v2/sql/phase_3_brouillon.sql @@ -134,14 +134,6 @@ FROM pgr_drivingDistance( 107, 300); -WITH t AS ( -SELECT * -FROM pgr_drivingDistance('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr -WHERE SOURCE IS NOT NULL AND id > 0', -1,915) -) -SELECT node, edge, round(agg_cost) FROM t ORDER BY seq DESC LIMIT 1; - -- RAZ de la topologie pgRouting diff --git a/scripts_v2/sql/phase_3_linear_referencing.sql b/scripts_v2/sql/phase_3_linear_referencing.sql deleted file mode 100644 index 28f0c66..0000000 --- a/scripts_v2/sql/phase_3_linear_referencing.sql +++ /dev/null @@ -1,115 +0,0 @@ - - -/* - -https://postgis.net/docs/ST_AddMeasure.html -http://postgis.net/docs/ST_LocateAlong.html -https://qastack.fr/gis/115341/point-sampling-along-a-pole-wrapping-coastline-with-postgis -https://qastack.fr/gis/88196/how-can-i-transform-polylines-into-points-every-n-metres-in-postgis -https://www.ibm.com/docs/en/informix-servers/12.10?topic=functions-st-locatealong-function - -*/ - - --- RAZ de la topologie pgRouting -TRUNCATE TABLE phase_3_troncons_pgr; -ALTER SEQUENCE phase_3_troncons_pgr_id_seq RESTART WITH 1; -VACUUM phase_3_troncons_pgr; - -TRUNCATE TABLE phase_3_troncons_pgr_vertices_pgr; -ALTER SEQUENCE phase_3_troncons_pgr_vertices_pgr_id_seq RESTART WITH 1; -VACUUM phase_3_troncons_pgr_vertices_pgr; - --- RAZ de la couche de PK -TRUNCATE TABLE phase_3_pk ; -VACUUM phase_3_pk; - - - --- on supprime ce qui concerne le secteur -DELETE FROM phase_3_troncons_pgr WHERE secteur_id = 100 ; -INSERT INTO phase_3_troncons_pgr (secteur_id, path_seq, osm_id, highway, type, oneway, ref, name_fr, name_br, the_geom) - SELECT - secteur_id, path_seq, osm_id, highway, type, oneway, ref, name_fr, name_br, the_geom - FROM phase_2_trace_pgr -WHERE secteur_id = 100 ; - - --- calcul des coûts (longueur) -UPDATE phase_3_troncons_pgr -SET -cost = trunc(st_length(the_geom)::numeric,2), -reverse_cost = trunc(st_length(the_geom)::numeric,2) -WHERE secteur_id = 100 ; - --- calcul de la topologie -SELECT pgr_createTopology('phase_3_troncons_pgr', 0.001, rows_where:='secteur_id=100', clean:=false); - - - --- calcul d'un point placé à 920 m sur une ligne de 1000 m -WITH linemeasure AS ( - WITH line AS ( - -- on récupère une ligne de 1000 m calculée par pgRouting - SELECT ST_Union(the_geom) AS the_geom - FROM pgr_drivingDistance('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr - WHERE SOURCE IS NOT NULL AND id > 0', - 9,1000) a - JOIN phase_3_troncons_pgr b ON a.edge = b.id - ) -SELECT - ST_AddMeasure(the_geom,0,ST_length(the_geom)) AS the_geom -FROM line -) -SELECT - ST_LocateAlong(the_geom,920) AS the_geom -FROM linemeasure; - - --- calcul de points tous les 200 m --- depuis une ligne calculée par pgRouting -WITH linemeasure AS ( - WITH line AS ( - -- on récupère une ligne de 1000 m calculée par pgRouting - SELECT ST_Union(the_geom) AS the_geom - FROM pgr_drivingDistance('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr - WHERE SOURCE IS NOT NULL AND id > 0', - 9,2000) a - JOIN phase_3_troncons_pgr b ON a.edge = b.id - ) -SELECT - generate_series(0, (ST_Length(line.the_geom))::int, 200) AS i, - ST_AddMeasure(the_geom,0,ST_length(the_geom)) AS the_geom -FROM line -) -SELECT - i - ,(ST_Dump(ST_GeometryN(ST_LocateAlong(the_geom, i), 1))).geom AS geom -FROM linemeasure; - - --- calcul de points tous les 920 m --- depuis une ligne calculée par pgRouting -WITH linemeasure AS ( - WITH line AS ( - -- on récupère un itinéraire calculée par pgRouting - SELECT ST_Union(the_geom) AS the_geom - FROM pgr_dijkstra('SELECT id, source, target, cost, reverse_cost FROM phase_3_troncons_pgr', - 9,1506) a - JOIN phase_3_troncons_pgr b ON a.edge = b.id - ) -SELECT - generate_series(0, (ST_Length(line.the_geom))::int, 920) AS i, - ST_AddMeasure(the_geom,0,ST_length(the_geom)) AS the_geom -FROM line -) -INSERT INTO phase_3_pk (pk_id,length_real,length_total,the_geom) -SELECT - ROW_NUMBER() OVER() + 11 AS pk - ,920 - ,i AS longueur_cumulee - ,(ST_Dump(ST_GeometryN(ST_LocateAlong(the_geom, i), 1))).geom AS the_geom -FROM linemeasure ; - - -