Below are projects I have worked on with students which resulted in publications and presentations.

We describe an efficient algorithm for
computing flood maps caused by overlapping sealevel rise (SLR) and
storm surges, based on the method used by the second New York City
Panel on Climate Change (NPCC2). Presently, areas at risk of storm
flooding are given by FEMA's baseflood elevation (BFE) maps which
depict the flooded areas as well as surge heights. Given a grid
digital elevation model (DEM), a BFE map and an arbitrary
amount x of SLR, our algorithm computes what part of the
terrain gets flooded if a sealevel rise of x overlaps the
storm surges in the BFE map.
As an example, we present results for Lincoln County in Maine, using a
2m resolution raster DEM obtained from LiDAR data, totaling 919
million points. Our algorithm overlaps FEMA's 1%annualchance storm
BFE map with an arbitrary SLR in less than 100 seconds on a laptop,
and produces results that are very close to the ones obtained using
NOAA's protocol in ArcGIS, which runs for hours.
Papers and presentations:


The viewshed of a point v on a
grid terrain T,
viewshed_{T}(v), is the set of grid points
in T that are visible from v. We describe a novel
algorithm for computing viewshed_{T}(v) using a
multiresolution approach: Given a parameter k that
represents the block size, we create a grid T' which is a
smaller, lowerresolution version of T, such that each
point in T' corresponds to a block of
√ k by
√ k
points in T. The key of our approach is using T' to
speed up the computation of viewshed_{T}(v) while not
introducing approximation. We describe a suite of experimental
results showing the performance of our algorithm in practice and a
speedup of more than an order of magnitude compared to previous
algorithms.
Papers and presentations:


We consider the problem of building a quadtree
subdivision for a set of edges in the place. Our approach is to
first build a quadtree on the vertices corresponding to the
endpoints of the edges, and then compute the intersections
between the edges and the cells in the subdivision. We
introduce the Kquadtree: for any
k>=1, a kquadtree is a linear compressed
quadtree that has O(n/k) cells
with O(k) points each. We show how to build a
kquadtree cacheefficienctly. As an application, we
consider one of the basic problems in spatial data
structures, map overlay, and implement it via
kquadtrees. We report experimental results.
Papers and presentations:


We describe an algorithm to compute the
viewshed with linear interpolation and without using any
approximation. The algorithm is based on computing and
merging horizons, and in practice it is faster than our
previous viewshed algorithms, while being more
accurate. We prove that the complexity of horizons on a
grid of n points is O(n), improving on the
general O(n \alpha (n)) bound on triangulated
terrains.
Papers and presentations:


We describe the design and engineering
of three new algorithms for computing viewsheds on massive
grid terrains. The first two algorithms sweep the terrain
by rotating a ray around the viewpoint, while maintaining
the terrain profile along the ray. The third algorithm
sweeps the terrain centrifugally, growing a starshaped
region around the viewpoint, while maintaining the
approximate visible horizon of the terrain within the
swept region. On a terrain of n grid points,
these algorithms run in O(sort(n)) I/O's in the
I/Omodel if Agarwal and Vitter. This algorithm runs in
O(n) time and O(scan(n)) I/O's and is cacheoblivious.
Papers and presentations:


Given a terrain T and a point v, the
viewshed of v is the set of points in T that are visible from
v. We consider the problem of computing the viewshed of a
point on a very large grid terrain in external memory. We
describe an algorithm for this problem in the cacheaware and
cacheoblivious models, together with an implementation and an
experimental evaluation. Our algorithm is based on Van
Kreveld's O(n lg n) algorithm for this problem in internal
memory and uses O(sort(n)) I/Os, where sort(n) is the
complexity of sorting n items of data in the I/Omodel of
Agarwal and Vitter.
r.viewshed: We ported the algorithm to GRASS, where it was adopted as one of their core modules and extended to use linear interpolation by their developer group. Papers and presentations:

 Our module Terracost computes leastcost surfaces on massive grid terrains and is motivated by cost analysis in GIS. Computing a leastcost surface corresponds to the problem of computing multisource shortest paths (MSSP) using an implicit graph given by grid topology. From a memory point of view, MSSP has a larger cache footprint and poorer spatial locality compared to SSSP. When the standard MSSP algorithm runs in external memory, it starts to thrash. Terracost is based on an I/Oefficient algorithm where the idea is to partition the graph in subgraphs with small boundaries, and replace each subgraph by a complete graph on its boundary vertices. The resulting graph has fewer vertices, preserves shortest paths, and the computation can be performed I/Oefficiently. Terracost is available as r.terracost in GRASS (extension library). Papers and presentations:
A tarball of the code is available here: r.terracost1.1 (works with GRASS 6.x)  r.terracost1.0 (works with GRASS 5.3). 
While most terrain data in GIS is represented in grid form, grids store a lot of redundancy in areas that are flat and could be represented with fewer points. Simplifying a grid into a TIN (triangular irregular networks) is promising because TINs have the potential of being more space efficient than grids, especially for very large terrains. In this project we took a first step at simplifying a grid terrain into a TIN within a specified error threshold using an approach that scales to terrains that are arbitrarily larger than main memory. Our module, r.refine can work standalone or as an addon in GRASS open source GIS. Papers and presentations:
