
Contents
In this tutorial we will use geospatial data from GIS.NY.GOV and New York City’s Open Data portal. See my list of geospatial data sources for further sources of urban data.
Download a digital elevation model (DEM)
and digital surface model (DSM)
for Governor’s Island, New York, NY, USA from
GIS.NY.GOV.
NYC commissioned a airborne lidar survey in 2017 and created
a 1 ft resolution bare earth DEM and a DSM from the lidar data.
The DEM represents the shape of the bare ground
with no trees, buildings, or other structures,
while the DSM represents the shape of the earth’s surface
and everything on it.
Large rasters such as these are often divided
into tiles for easier downloading.
For Governor’s Island download DEM and DSM tiles
using an FTP client:
be_NYC_020.tif and
hh_NYC_020.tif.
Alternatively
to download with your browser,
use these Google Drive links:
be_NYC_020.tif
and
hh_NYC_020.tif.
Also download and extract either
Tidally Coordinated Shoreline
or Shoreline
data from NYC Open Data.
For Shoreline
in the web mapping interface
select Export, Download, and then Original
to download Shoreline.zip.
Note:
If you choose Shapefile instead of Original
for Shoreline
you will download a reprojected version in a
geographic CRS in Latitude Longitude.
Start GRASS GIS.
In the startup screen, select New location.
In the Define new GRASS Location window,
first browse to set your GRASS GIS database directory
where the location will be created
and then set a name for the location.
Hit Next to continue.
Select Read projection and datum terms from a georeferenced data file.
Browse to select NYC_2017_LiDAR_Low_Tide_Shoreline.shp.
Hit Next and then Finish to continue.
Accept the Import data to import the data into the new location.
In the startup screen, select the new location
and the PERMANENT mapset and start the GRASS session.
Print information about the projection in the console with
g.proj
with the -p flag.
The map projection should be
NAD 1983 StatePlane New York Long Island FIPS 3104 Feet
in U.S. Surveyor’s Foot.
g.proj -p
Note:
Alternatively you can create a new GRASS location
for New York Long Island State Plane Feet
using the EPSG code 2263.
In the Map Display click the Select vector feature button.
Click on the shore of Governor’s Island to select it.
In the Select features dialog, click Create a new map.
Rename the map to shoreline using
g.rename.
Alternatively query the vector map to find
the category for Governor’s Island’s shore
and then use
v.extract
to extract it by category value.
v.extract input=NYC_2017_LiDAR_Low_Tide_Shoreline cats=184 output=shoreline
Or use an SQL query to select the category value.
v.extract input=NYC_2017_LiDAR_Low_Tide_Shoreline where="Cat = '184'" output=shoreline
Set the computational region and mask to the extracted shoreline for Governor’s Island.
g.region vector=shoreline res=1
r.mask vector=shoreline
Then import the DEM and DSM with the module
r.in.gdal
with the -r flag to limit the import the region.
r.in.gdal -r input=be_NYC_020.tif output=elevation_2017
r.in.gdal -r input=hh_NYC_020.tif output=surface_2017
Note: When you want to import data in a different projection into a location, use r.import to import with automatic reprojection.
Set a color table for each raster map using r.colors. Add a legend with d.legend.
r.colors -e map=elevation_2017 color=viridis
r.colors -e map=surface_2017 color=viridis
| Digital elevation model | Digital surface model |
|---|---|
![]() |
![]() |
For practice, download and import geospatial data for Central Park in New York City. Download the Landcover 2017 geotiff and Open Space shapefile from the NYC Open Data portal. For the Open Space data, in the Export menu select Download Original to download the PARK.zip archive containing PARK.shp in New York Long Island State Plane in U.S. Surveyor’s Feet. (Other export options may be in Lat Long instead.) Check https://orthos.dhses.ny.gov/ to determine which DEM tiles to download for Central Park. Then download the NYC Topobathymetric 2017 DEM (be_NYC_025.tif & be_NYC_029.tif) and NYC Topobathymetric 2017 DSM (hh_NYC_025.tif & hh_NYC_029.tif) tiles from GIS.NY.GOV. Optionally download orthoimagery tiles from 2006, 2008, 2012, and 2016 for Central Park at https://orthos.dhses.ny.gov/.
Create a new location called nyspf_central_park
from PARK.shp.
Then start GRASS GIS.
In the startup screen, select New location.
In the Define new GRASS Location window,
first browse to set your GRASS GIS database directory
where the location will be created
and then set a name for the location.
Hit Next to continue.
Select Read projection and datum terms from a georeferenced data file.
Browse to select PARK.shp.
Hit Next and then Finish to continue.
Accept the Import data to import the data into the new location.
In the startup screen, select the new location
and the PERMANENT mapset and start the GRASS session.
Extract Central Park from the PARK vector map
either by selecting the feature
or running v.extract.
In the Map Display click the Select vector feature button.
Click on the polygon bounding Central Park to select it.
In the Select features dialog, click Create a new map.
Then rename the map to central_park using
g.rename.
Alternatively use
v.extract
with an SQL query of park names in the attribute table.
v.extract input=PARK where="PARK_NAME = 'Central Park'" output=central_park
Then delete the PARK vector map using
g.remove
with flag -f to force the removal of the map.
g.remove type=vector name=PARK
Set your computational region to the new central_park vector map
at 1 foot resolution and the save the region.
This will limit the scope of all raster operations
to our study area in central park.
g.region vector=central_park res=1 save=central_park
Or zoom to an extent containing Central Park, set the computational region extent from the display, and save the region.
g.region res=1 save=central_park
Import the tiles of digital elevation and digital surface models using
r.in.gdal
with the -r flag to limit the import the region.
r.in.gdal -r input=be_NYC_025.tif output=elevatione_2017_a
r.in.gdal -r input=be_NYC_029.tif output=elevation_2017_b
r.in.gdal -r input=hh_NYC_025.tif output=surface_2017_a
r.in.gdal -r input=hh_NYC_029.tif output=surface_2017_b
Note that the color ramps for the pairs of tiles do not match. Each pair of tiles needs to be patched together.
| Elevation Raster Tiles | Surface Raster Tiles |
|---|---|
![]() |
![]() |
Patch the two tiles of the digital elevation and digital surface models together with the module r.patch.
r.patch input=elevation_2017_a,elevation_2017_b output=elevation_2017
r.patch input=surface_2017_a,surface_2017_b output=surface_2017
Remove both original pairs of raster tiles with
g.remove
with flag -f.
g.remove -f type=raster name=elevation_2017_a,elevation_2017_b,surface_2017_a,surface_2017_b
Set a color table for each raster map using
r.colors.
Add a legend with
d.legend.
Then Save .png images using the save display to file button.
r.colors -e map=elevation_2017 color=viridis
r.colors -e map=surface_2017 color=viridis
d.legend raster=elevation_2017 at=75,96,4,6 font=Lato-Regular fontsize=12
d.legend raster=surface_2017 at=75,96,4,6 font=Lato-Regular fontsize=12
| Digital Elevation Model, Central Park | Digital Surface Model, Central Park |
|---|---|
![]() |
![]() |
Optionally set a mask to the central_park vector map.
r.mask vector=central_park
| Digital Elevation Model, Central Park | Digital Surface Model, Central Park |
|---|---|
![]() |
![]() |