A reproducible R GIS study of the Quabbin Reservoir, Massachusetts, created between 1938 and 1946 by damming and flooding the Swift River Valley. Four towns — Dana, Enfield, Greenwich, and Prescott — were disincorporated on 28 April 1938 and about 2,500 residents were relocated; the reservoir supplies metropolitan Boston, roughly 105 km (65 miles) to the east.
It is a multi-layer study: several spatial layers that, read together, describe how the reservoir was sited and what it replaced — a valley whose terrain forms a natural basin, four towns that had been losing population for decades, and a present-day map in which their land has been absorbed by the surrounding towns.
| # | Figure | What it shows | Source |
|---|---|---|---|
| 1 | 01_locator.png |
The valley in shaded relief, with a Massachusetts "you are here" chip (Quabbin → Boston) | AWS Terrain Tiles + TIGER |
| 2 | 02_dem_hillshade.png |
The Swift River Valley basin (shaded relief) | AWS Terrain Tiles (elevatr) |
| 3 | 03_reservoir_towns.png |
The four former towns over the present reservoir | DEM + town points |
| 4 | 04_watershed.png |
The reservoir within its drainage (regional context) | USGS WBD HUC-10 |
| 5 | 05_erasure.png |
The former town land, now divided among surrounding towns | TIGER municipalities |
| 6 | 06_town_lifelines.png |
Each town's span, charter to the 1938 disincorporation | town records |
| 7 | 07_population_decline.png |
Decennial census population, 1900–1920 | US Census 1920 (Number of Inhabitants) |
| 8 | 08_hero.png |
Terrain, reservoir, and the former town sites together | all of the above |
| 9 | 09_floodfill.png |
The reservoir filling in equal-area stages to the 530-ft full pool (schematic synthetic basin) | DEM + modeled depth |
| 10 | 10_aqueduct.png |
The aqueduct route, ~105 km east to Boston | hand-placed coordinates |
| 11 | 11_crosssection.png |
West-east valley cross-section at the 530-ft pool | DEM |
| 12 | 12_losses.png |
"By the numbers": what was removed and what is supplied | MWRA / DCR / histories |
| 13 | 13_terrain3d.png |
3D view of the valley and reservoir | DEM (persp) |
| 14 | 14_dana_lidar.png |
Dana Common in 1 m LiDAR (surviving town site) | USGS 3DEP LiDAR |
| 15 | 15_prescott_lidar.png |
The Prescott Peninsula in 1 m LiDAR | USGS 3DEP LiDAR |
| 16 | 16_roads.png |
The valley's 1890s road network, with the reservoir overlaid | USGS Belchertown 1893 + Barre 1894 |
| 17 | 24_prescott_survey.png |
The whole Prescott Peninsula: 1893 survey vs. the LiDAR imprints surviving today | MassGIS LiDAR (3-tile mosaic) + 1893 quad |
| 18 | 24_enfield_survey.png |
Enfield (Winsor Dam): surviving roads on the dry south vs. the drowned center | MassGIS 1 m LiDAR + 1893 quad |
| 19 | 24_dana_survey.png |
Dana Common: the surviving common/ridge vs. the drowned village | MassGIS 1 m LiDAR + 1893 quad |
| 20 | 25_prescott_xref.png |
Ground-truth: the 1893 road network extracted and cross-referenced against the LiDAR traces | MassGIS LiDAR + 1893 quad |
Plus a schematic reservoir-filling animation (quabbin_floodfill.gif) and an
interactive LiDAR imprint explorer in map/ (mobile-first Leaflet):
a "LiDAR relief" layer covering the whole reservoir in bare-earth LiDAR —
pan and zoom across every acre the reservoir spared to hunt the relict streets,
house-lot outlines and cellar-hole pits of the drowned villages still imprinted
in the ground (crisp 1 m tiles at the surviving village sites). Toggle the
auto-traced lines and the historical ground-truth surveys (1890s and 1940s),
raise the pool over the valley, and follow the aqueduct east to Boston.
R is the only requirement. On Ubuntu the heavy spatial stack installs as binaries (no source compiles):
sudo apt-get install -y r-base-core r-cran-sf r-cran-terra r-cran-raster \
r-cran-ggplot2 r-cran-dplyr libcurl4-openssl-dev libxml2-dev libjpeg-dev
Rscript -e 'install.packages(c("elevatr","osmdata","tigris","ggspatial","ggnewscale","ggrepel","patchwork"))'Then, from the study folder:
Rscript run_all.R(run_all.R locates its own folder, so it also works from anywhere —
e.g. Rscript path/to/run_all.R.)
Downloads are cached under data/cache/ (git-ignored), so the first
run takes a few minutes and every run after that is ~90 seconds — except the
LiDAR stages (14–16), which mosaic the Prescott Peninsula and tile the whole
reservoir, and are cache-guarded: their first build takes several minutes,
then later runs skip any area whose output/24_*_survey.png (or
25_prescott_xref.png) already exists, and stage 16 skips entirely while
map/data/reservoir_ghost.json and every overlay it lists are present
(delete a figure — or reservoir_ghost.json — to force a rebuild). The
twenty figures and the GIF land in output/; the web-map
GeoJSON, the historical map overlays, and the LiDAR relief/imprint overlays land in
map/data/. (The GIF needs
ImageMagick — apt-get install imagemagick; without it the pipeline still
produces the panel 09.)
The interactive map is static files — serve map/ over HTTP, e.g.
python3 -m http.server --directory map, then open localhost:8000.
quabbin/
├── run_all.R one-command reproduction
├── R/
│ ├── 00_setup.R packages, CRS (EPSG:26986), area of interest, palettes
│ ├── 01_fetch_data.R elevatr DEM · OSM/DEM reservoir · TIGER towns · USGS watershed
│ ├── 02_build_layers.R reproject, hillshade, carve the reservoir, assemble layers
│ ├── 03_maps.R the six spatial figures (shaded relief + vector overlays)
│ ├── 04_population.R the two displacement charts (real census + lifelines)
│ ├── 05_floodfill.R schematic reservoir-filling frames + GIF + panel + stage GeoJSON
│ ├── 06_aqueduct.R the aqueduct-to-Boston map + infrastructure GeoJSON
│ ├── 07_export_web.R export towns/reservoir/watershed GeoJSON for the web map
│ ├── 08_profile.R west-east valley cross-section
│ ├── 09_losses.R the "by the numbers" figure
│ ├── 10_terrain3d.R 3D terrain view (base-R persp)
│ ├── 11_preflood.R historical map overlays (1890s + 1940s) for the web map
│ ├── 12_lidar.R USGS 3DEP LiDAR of Dana Common & the Prescott Peninsula
│ ├── 13_roads.R the valley's 1890s road network, reservoir overlaid
│ ├── 14_imprints.R LiDAR imprint survey (MassGIS) + per-town survey figures
│ ├── 15_xref.R cross-reference: extracted 1893 roads vs the LiDAR traces
│ ├── 16_reservoir.R full-reservoir "LiDAR relief" coverage for the explorer
│ └── lidar_utils.R shared helpers for the LiDAR stages 14–16 (LRM, hillshade, line tracing, MassGIS export)
├── data/
│ ├── drowned_towns.csv the four towns: location, county, charter & end dates
│ └── town_population.csv real US Census counts 1900–1920 + peaks + 1938 dissolution
├── map/ interactive imprint explorer (index.html + vendored Leaflet)
│ └── data/ GeoJSON + historical overlays (1890s/1940s) + LiDAR relief/imprint overlays
└── output/ the rendered figures + GIF (committed)
Every network fetch in 01_fetch_data.R is wrapped so one unreachable service
never breaks the run — it degrades to a documented fallback instead.
- Elevation — AWS Terrain Tiles (SRTM/USGS, public domain) via
elevatr::get_elev_raster(z = 11), reprojected to NAD83 / Massachusetts Mainland (EPSG:26986) and hillshaded withterra. - Reservoir — derived from MassGIS 1 m LiDAR (downsampled to 10 m): the
largest contiguous area at/below the 530-ft full pool. The LiDAR resolves the
terrain far better than the coarse regional DEM, but at 10 m even it still leaks
~2 sq km of below-dam ground into the largest patch (the dams don't perfectly
seal), so the footprint is then clipped at the Winsor Dam / Goodnough Dike
line — the dams are the reservoir's true southern boundary — and the water no
longer spills south past the dams onto Belchertown/Ware (an earlier artifact
that had leaked into the gallery). If MassGIS is unreachable it falls back to the coarse-DEM carve
(largest polygon below the pool), and an OpenStreetMap shoreline path
(
RESERVOIR_METHOD <- "osm") is also available but opt-in (Overpass rate-limits cloud IPs). - Watershed — USGS Watershed Boundary Dataset HUC-10 units (public domain), queried live from the National Map ArcGIS service. These are deliberately shown as regional context, not as the catchment boundary: the dissolved HUC-10s are several times larger than the ~120 sq mi DCR-defined Quabbin watershed (a MassGIS layer that was not reachable here). The map says so.
- Modern municipalities — US Census TIGER county subdivisions (2021,
public domain) via
tigris. - The four towns — there is no clean, freely downloadable historic GIS boundary for towns abolished in 1938, so they are plotted as labeled points at their historic centers (coordinates from the towns' Wikipedia pages).
- Population — real decennial U.S. Census counts for all four towns, read straight from the 1920 Number of Inhabitants, Massachusetts bulletin (Table 2, public domain), via OCR + a hand check of the scanned page: for 1900 / 1910 / 1920, Dana 790 → 736 → 599, Enfield 1,036 → 874 → 790, Greenwich 491 → 452 → 309, Prescott 380 → 320 → 236. Every town was already shrinking before the reservoir. Earlier peaks (Enfield ~1,100 in 1850, Prescott ~750 in 1830) and the 1938 dissolution (~2,500 displaced in all) are shown as annotated context, not plotted as census points. All four disincorporated 28 April 1938.
- Reservoir-filling animation —
05_floodfill.R, explicitly schematic. Modern DEMs (AWS Terrain and LiDAR alike) capture only today's water surface (~530 ft), not the drowned valley floor, so a true elevation-based fill is impossible. Instead the basin is modeled synthetically: depth grows with distance from the shoreline (the broad main basin deepest, the narrow arms shallow), scaled to the reservoir's surveyed ~150 ft maximum depth; the pool then rises over that synthetic bed in equal-area stages, filling the deep central channel first and spreading to the arms (≈7 % → 100 %). Frames become a GIF (ImageMagick), a small-multiples panel, and per-stage GeoJSON for the map slider. It is a schematic of how the basin filled (1939–1946), not surveyed bathymetry — captioned as such. - Aqueduct & dams — the route (Quabbin → Wachusett → Boston) and the dams (Winsor Dam, Goodnough Dike) are hand-placed from known coordinates and labeled schematic; they convey the ~105 km eastward course of the water, not a surveyed alignment.
- Cross-section (
08_profile.R) — a west-east transect of the DEM at the reservoir's widest point. The DEM retains sub-pool relief of the drowned valley, but its underwater values are approximate, so the figure cites the surveyed maximum depth (~150 ft) rather than asserting DEM depths. - By the numbers (
09_losses.R) — documented quantities (surface area, volume, shoreline, displaced residents, relocated graves, buildings razed, people supplied) compiled from MWRA, Massachusetts DCR, and regional histories; several vary by source and are shown with ranges. - 3D view (
10_terrain3d.R) — base-Rpersp()(no GPU required), with the reservoir drawn as a flat pool over the relief. - Historical map overlays (
11_preflood.R) — two eras of USGS topographic survey (Historical Topographic Map Collection, public domain), each a mosaic of the quadrangles that cover the valley (one 15' sheet only spans half of it): the 1890s (Belchertown 1893 + Barre 1894, 1:62,500) and the more-detailed 1940s (Winsor Dam + Quabbin Reservoir 1944, 1:31,680, mapped as the valley was being taken). Each sheet is clipped to its neatline (dropping the collar and forming a clean seam), the sheets are tonally harmonized so the join is seamless, reprojected to EPSG:4326, cropped to the reservoir window and exported as a JPEG- bounds. A
histmaps.jsonmanifest drives the explorer's year selector.
- bounds. A
- LiDAR of the surviving sites (
12_lidar.R) — 1 m bare-earth LiDAR from the USGS 3DEP dynamic elevation service (public domain) for Dana Common and the Prescott Peninsula, the two areas above the full pool. A denoised low-sun hillshade carries the terrain; a local relief model (elevation minus its local mean) then flags genuine depressions deeper than ~0.5 m — cellar holes and road cuts — in red, with a slope mask so natural gullies on the flanks aren't false- flagged. Paired with the 1893 quad as static then-&-now figures (the explorer's relief now comes from the full-reservoir LiDAR relief layer,16_reservoir.R). Submerged areas are not shown: the buildings were demolished and LiDAR cannot penetrate water, and the conservative threshold under-flags rather than fills the frame with canopy noise. - The valley's road network (
13_roads.R) — renders the real Swift River Valley roads from the georeferenced 1890s USGS survey (the Belchertown 1893 + Barre 1894 mosaic, so the whole reservoir is covered), with the present reservoir overlaid in blue so you can see which roads and villages drowned. The four town centres are marked and the routes out to the surviving neighbours are labelled at the frame edges. - The LiDAR imprint survey (
14_imprints.R) — the heart of the relict-landscape work. For the dry land that survives in each town it pulls MassGIS 1 m bare-earth LiDAR (the 2013–2021 statewide DEM ImageServer, public domain; native NAD83 / MA metres, cleaner than the 3DEP seamless used in12_lidar.R) and renders a composite relief — an 8-direction hillshade emphasised by a Local Relief Model — that makes faint linear features read. It then auto-traces the relict network: elongated negative relief (sunken roadbeds, cart paths) and elongated positive relief (banks, stone walls), keeping only long, straight connected components on gentle ground so blobs and slope noise are dropped. Every panel is ground-truthed against the 1893 quad. Findings, honestly: Prescott is the richest (a clear road plus many lineations); Enfield's dry south (toward Winsor Dam) keeps road traces; Dana's common/ridge keeps a road cut; Greenwich's center is entirely under water (a static "what drowned" panel, no imprints). The auto-trace is a candidate finder — it cannot perfectly separate man-made lines from natural slope features, so the relief itself is the primary evidence and the traces are an explorable overlay. This is, as far as we found, a gap in the public record: the lost-town history is well documented (J.R. Greene's Atlas of the Quabbin Valley; the Swift River Valley Historical Society) and LiDAR for New England relict landscapes is proven (UConn's Ouimet Lab; Johnson & Ouimet 2014), but not combined into a LiDAR imprint survey of the four towns. The Prescott Peninsula runs ~12 km — beyond the server's single-export cap — so it is mosaicked from three ~2 m strips into one DEM for full-length coverage; the stage is cache-guarded (skips any area whose survey figure — and, for web areas, trace overlay — already exists), so only the first build pays the render cost. Exported as static survey figures (output/24_*) and as auto-trace web overlays (transparent over water) with bounds inmap/data/imprints.json; the explorer's relief itself comes from the full-reservoir LiDAR relief layer (16_reservoir.R). - Interactive imprint explorer — Leaflet (vendored locally, no CDN dependency), mobile-first: a full-screen map, bottom-sheet Layers control, big touch targets, a full-width pool slider, and zoom to z19. It serves the full-reservoir LiDAR relief (above) as its primary layer, the auto-traced lines, "jump to" buttons (whole reservoir / Prescott Center / Dana Common), the historical map overlays (1890s/1940s, year-selectable), the drowned-town popups, the flood stages, and the aqueduct, over Esri World Hillshade + CARTO label tiles. The basemap tiles need internet; every layer the study itself produces is served locally.
- Full-reservoir "LiDAR relief" (
16_reservoir.R) — the explorer's headline layer. It tiles the entire Quabbin land area in MassGIS bare-earth LiDAR and renders each tile as a fine Local Relief Model (elevation minus its local mean over ~13 cells) — the rendering that makes the drowned villages' street plans, house-lot outlines and cellar-hole pits read directly in the bare ground — transparent over water. Broad coverage is ~2 m with a fixed contrast span so tiles match seamlessly; the surviving village sites (Prescott Center, Dana Common) get crisp ~1 m tiles. Turn on LiDAR relief in the explorer and zoom in to hunt footprints anywhere the reservoir spared. Cache-guarded twice: the LiDAR downloads (the slow part), and the whole stage — re-runs skip the re-export whilemap/data/reservoir_ghost.jsonand every overlay it lists exist (delete the manifest to force a rebuild). Exports web overlays +map/data/reservoir_ghost.json. - Ground-truth cross-reference (
15_xref.R) — closes the loop on the imprints: it extracts the 1893 road network straight from the quad (dark, low-saturation linework, morphologically closed, kept only where elongated → roads, not text), then classifies the LiDAR road traces by proximity — within 14 m of a mapped 1893 road = a confirmed surviving roadbed (green), the rest unverified (orange). Demonstrated on the mid-peninsula (the old Prescott village area). Honest limit: auto-extraction of the scanned linework catches the main roads, not the full network, so "unverified" is not "newly discovered"; the explorer's 1893 fade overlay remains the fuller visual ground-truth. Rendersoutput/25_prescott_xref.png.
R · sf · terra · elevatr · osmdata · tigris · ggplot2 ·
ggnewscale · ggrepel · ggspatial · patchwork · ImageMagick · Leaflet · GDAL 3.8
The study exists in two places:
mapzimus/quabbin— the standalone repository and source of truth. Changes are made here first.- The
quabbin/folder of the portfolio repository, a copy that serves the live project page and the interactive explorer at maxwellhowegis.com. After a change lands in the standalone repo, copy it into the portfolio folder to publish it.
The code (the R/ scripts, run_all.R, and the explorer in map/) is
released under the MIT License. The underlying elevation, LiDAR,
census, and historical-map data are U.S. public-domain sources, credited in
Data, methods, and honest caveats above.
Part of an ongoing series of multi-layer GIS studies of geography-shaped American places. Data is open; figures are reproducible from the scripts above.

