Skip to content

Instantly share code, notes, and snippets.

@dannguyen
Last active May 5, 2016 00:17
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save dannguyen/51ffaf546b16e1a214aa85bb0c623e1a to your computer and use it in GitHub Desktop.
Save dannguyen/51ffaf546b16e1a214aa85bb0c623e1a to your computer and use it in GitHub Desktop.
Notes on geocoding and GIS in Python 3.x, OGR, osgeo, pyproj and all that jazz
# had problems getting Basemap and ogr2ogr (which depends on updated libgdal)
# working at the same time
# in the end it kind of worked but don't know how to zoom the basemap in
# need to look at library that just handles shapefile reading
%matplotlib
from shutil import unpack_archive
from pathlib import Path
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
import requests
URL = 'http://apps.sfgov.org/datafiles/view.php?file=sfgis/bayarea_zipcodes.zip'
DDIR = Path('apps.sfgov.org')
DDIR.mkdir(exist_ok=True, parents=True)
ZIPNAME = Path(URL).name
SHPNAME = DDIR.joinpath('bayarea_zipcodes_4326')
resp = requests.get(URL)
with open(ZIPNAME, 'wb') as f:
f.write(resp.content)
unpack_archive(ZIPNAME, extract_dir=str(DDIR))
# must convert to lat/lng
# ESRI:102643
# https://epsg.io/102643
# OK, let's do this
# ogr2ogr -t_srs EPSG:4326 bayarea_zipcodes_4326.shp bayarea_zipcodes.shp
# Got this error:
# ERROR 4: Unable to open EPSG support file gcs.csv.
# Try setting the GDAL_DATA environment variable to point to the
# directory containing EPSG csv files.
# Failed to process SRS definition: EPSG:4326
# But had to manually add EXPORT GDAL_PATH= {{gdal-config --datadir}}
#
# https://github.com/ContinuumIO/anaconda-issues/issues/221
fig, ax = plt.subplots()
m = Basemap(ax=ax)
m.readshapefile(str(SHPNAME), 'stuffwhatisthis')

Nice tutorial here on using basemap:

Visualization: Mapping Global Earthquake Activity: This project introduces the Basemap library, which can be used to create maps and plot geographical datasets.

http://introtopython.org/visualization_earthquakes.html


osgeo

  • How to open a projection file?
  • How to reproject?

https://pcjericks.github.io/py-gdalogr-cookbook/projection.html

http://gis.stackexchange.com/questions/145015/is-it-possible-to-look-at-the-contents-of-shapefile-using-python-without-an-arcm/145029

This kind of works?

from pathlib import Path
from pyproj import Proj, transform
import shapefile

SHP_FOLDER = Path('nycb2010_16a')
SHP_FILENAME = SHP_FOLDER.joinpath('nycb2010.shp')
sf = shapefile.Reader(str(SHP_FILENAME))

shapelib

This is supposedly useful for converting a shapefile from one projection to another.

This was a good readme, kind of...not sure what parts of the code are germane to the task: https://glenbambrick.com/2016/01/24/reproject-shapefile/

pyshp

I installed it via pip.

https://pypi.python.org/pypi/pyshp


On pyproj:

If you do need to deal with projections programmatically you basically have one choice: the PROJ4 library. It is one of the few free libraries, if not the only library period, that comprehensively deals with re-projecting goespatial data. Fortunately it has bindings for just about every language out there and is incorporated into many libraries including OGR. There is a Python project called pyproj which provides python bindings.

http://geospatialpython.com/2011/09/map-projections.html

from shutil import unpack_archive
from pathlib import Path
import requests
URL = 'http://www2.census.gov/geo/tiger/GENZ2015/shp/cb_2015_us_cd114_20m.zip'
DDIR = Path('census-shp')
DDIR.mkdir(exist_ok=True, parents=True)
ZIPNAME = Path(URL).name
SHPNAME = Path(URL).stem + '.shp'
resp = requests.get(URL)
with open(ZIPNAME, 'wb') as f:
f.write(resp.content)
unpack_archive(ZIPNAME, extract_dir=str(DDIR))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment