Note
Go to the end to download the full example code.
Changing the Projection system for a model#
In this example we show how to change the projection system for a Polaris model.
Image credit to: https://ubc-library-rc.github.io/map-projections/content/diffs-geo-proj.html
from pathlib import Path
from tempfile import gettempdir
from polaris import Polaris
from polaris.network.network import Network
from polaris.project.project_restorer import create_db_from_csv
from polaris.utils.database.db_utils import commit_and_close
from polaris.utils.database.standard_database import DatabaseType
sphinx_gallery_thumbnail_path = ‘../../examples/editing_models/projections.jpg’
model_path = Path("/tmp/Grid")
model = Polaris.from_dir(model_path)
model.upgrade()
old_net = model.network
WARNING
BE CAREFUL WHEN CHOOSING A PROJECTION SYSTEM. It must be projected in meters for Polaris to function properly
reprojected_folder = Path(gettempdir()) / "new_projection"
# Let's say we choose the projection system 5069
old_net.ie.dump(reprojected_folder, target_crs=5069)
new_net_file = Path(gettempdir()) / "Grid-Supply_reprojected.sqlite"
create_db_from_csv(new_net_file, reprojected_folder, DatabaseType.Supply)
We also need to update the length of a few elements in the network to be consistent with the new projection
with commit_and_close(new_net_file, spatial=True) as conn:
for table in ("Link", "Road_Connectors", "Transit_Links", "Transit_Bike", "Transit_Walk"):
conn.execute(f"update {table} set length=round(ST_Length(geo), 8)")
Compute the overlay between the zone and the block group geometries
new_net = Network.from_file(new_net_file)
Let’s compare the distances now
new_lengths = new_net.tables.get("Link")[["link", "length"]].set_index("link")
new_lengths.rename(columns={"length": "new_length"}, inplace=True)
old_lengths = old_net.tables.get("Link")[["link", "length"]].set_index("link")
old_lengths.rename(columns={"length": "old_length"}, inplace=True)
old_lengths.join(new_lengths).plot.scatter("old_length", "new_length")
<Axes: xlabel='old_length', ylabel='new_length'>
new_net.close()
old_net.close()
Total running time of the script: (0 minutes 2.587 seconds)