Predicting floods caused by the rise of the sea level is a critical task for preventing large scale catastrophes. Such predictions can potentially be made using a forecast of the sea level and a detailed model of the terrain. However, since available terrain datasets can easily exceed the size of the main memory of a standard computer, I/O (rather than internal computation time) can often become the bottleneck when computing such predictions. Thus to perform predictions efficiently we need an I/O-efficient approach, which minimizes the transfer of data blocks between main memory and disk. Given a terrain raster T and a sea-level forecast raster S of N cells each, we examine the problem of computing the water level of the induced flood for each cell in T. We introduce an I/O-efficient algorithm for this problem that uses O((N/B) logM/B (X/B)) I/Os after O((N/B) logM/B (N/B)) I/Os of preprocessing, where X is the number of local minima in T, and M and B are the size of main memory and data block, respectively. When X < M (which holds in practice) our algorithm requires optimal O(N/B) I/Os after preprocessing. We have implemented our algorithm and put considerable effort into engineering it. We present experiments that illustrate the efficiency and practicality of the algorithm, which is so efficient that work is underway to incorporate our results in the forecast services of the Danish Meteorological Institute.

Let T be a terrain, and let P be a set of points (locations) on its surface. An important problem in Geographic Information Science (GIS) is computing the visibility index of a point p on P, that is, the number of points in P that are visible from p. The total visibility-index problem asks for computing the visibility index of every point in P. Most applications of this problem involve 2-dimensional terrains represented by a grid of n × n square cells, where each cell is associated with an elevation value, and P consists of the center-points of these cells. Current approaches for computing the total visibility-index on such a terrain take at least quadratic time with respect to the number of the terrain cells. While finding a subquadratic solution to this 2D total visibility-index problem is an open problem, surprisingly, no subquadratic solution has been proposed for the one-dimensional (1D) version of the problem; in the 1D problem, the terrain is an x-monotone polyline, and P is the set of the polyline vertices. We present an O(n log2 n) algorithm that solves the 1D total visibility-index problem in the RAM model. Our algorithm is based on a geometric dualization technique, which reduces the problem into a set of instances of the red-blue line segment intersection counting problem. We also present a parallel version of this algorithm, which requires O(log2 n) time and O(n log2 n) work in the CREW PRAM model. We implement a naive O(n2) approach and three variations of our algorithm: one employing an existing red-blue line segment intersection algorithm and two new approaches that perform the intersection counting by leveraging features specific to our problem. We present experimental results for both serial and parallel implementations on large synthetic and real-world datasets, using two distinct hardware platforms. Results show that all variants of our algorithm outperform the naive approach by several orders of magnitude on large datasets. Furthermore, we show that our new intersection counting implementations achieve more than 8 times speedup over the existing red-blue line segment intersection algorithm. Our parallel implementation is able to process a terrain of 224 vertices in under 1 minute using 16 cores, achieving more than 7 times speedup over serial execution.

Let P be a set of points in ℝd, and let M be a function that maps any subset of P to a positive real number. We examine the problem of computing the exact mean and variance of M when a subset of points in P is selected according to a well-defined random distribution. We consider two distributions; in the first distribution (which we call the Bernoulli distribution), each point p ∊ P is included in the random subset independently, with probability π(p). In the second distribution (the fixed-size distribution), a subset of exactly s points is selected uniformly at random among all possible subsets of s points in P. This problem is a crucial part of modern ecological analyses; each point in P represents a species in d-dimensional trait space, and the goal is to compute the statistics of a geometric measure on this trait space, when subsets of species are selected under random processes. We present efficient exact algorithms for computing the mean and variance of several geometric measures when point sets are selected under one of the described random distributions. More specifically, we provide algorithms for the following measures: the bounding box volume, the convex hull volume, the mean pairwise distance (MPD), the squared Euclidean distance from the centroid, and the diameter of the minimum enclosing disk. We also describe an efficient (1 – ∊)-approximation algorithm for computing the mean and variance of the mean pairwise distance. We implemented three of our algorithms: an algorithm that computes the exact mean volume of the 2D bounding box in the Bernoulli distribution, an algorithm that computes the exact mean and variance of the MPD for d-dimensional point sets in the fixed-size distribution, and an (1 – ∊)-approximation algorithm for the same measure. We conducted experiments where we compared the performance of our implementations with a standard heuristic approach used in ecological applications. We show that our implementations can provide major speedups compared to the standard approach, and they produce results of higher precision, especially for the calculation of the variance. We also compared the implementation of our exact MPD algorithm with the corresponding (1 – ∊)-approximation method; we show that the approximation method performs faster in certain cases, while also providing high-precision approximations. We thus demonstrate that, as an alternative to the exact algorithm, this method can also be used as a reliable tool for ecological analysis.

In many scientific fields, it is common to encounter raster data sets consisting of categorical data, such as soil type or land usage of a terrain. A problem that arises in the presence of such data is the following: given a raster G of n cells storing categorical data, compute for every cell c in G the number of different categories appearing within a window centered at c. The window can either be a (2r+1)x(2r+1) square or a disk of radius r for a positive integer parameter r. We call this the categorical richness problem, and we present two algorithms for it: one for square windows that runs in O(n) time and one for circular windows that runs in O((1+K/r)n) time, where K is the number of distinct categories in G. The algorithms are not only efficient in theory, but also in practice: our experiments show that our algorithms can handle raster data of hundreds of millions of cells.

In many scientific applications it is required to reconstruct a raster dataset many times, each time using a different resolution. This leads to the following problem; let $\mathcal{G}$ be a raster of $\sqrt{N}$ x $\sqrt{N}$ cells. We want to compute for every integer 2 $\leq \mu \leq \sqrt{N}$ a raster $\mathcal{G}_\mu$ of [ $\sqrt{N}/\mu$ ] x [ $\sqrt{N}/\mu$ ] cells where each cell of $\mathcal{G}_\mu$ stores the average of the values of μ x μ cells of $\mathcal{G}$ . Here we consider the case where $\mathcal{G}$ is so large that it does not fit in the main memory of the computer. We present a novel algorithm that solves this problem in O(scan(N)) data block transfers from/to the external memory, and in θ(N) CPU operations; here scan(N) is the number of block transfers that are needed to read the entire dataset from the external memory. Unlike previous results on this problem, our algorithm achieves this optimal performance without making any assumptions on the size of the main memory of the computer. Moreover, this algorithm is cache-oblivious; its performance does not depend on the data block size and the main memory size. We have implemented the new algorithm and we evaluate its performance on datasets of various sizes; we show that it clearly outperforms previous approaches on this problem. In this way, we provide solid evidence that non-trivial cache-oblivious algorithms can be implemented so that they perform efficiently in practice.

In many GIS applications it is important to study the characteristics of a raster data set at multiple resolutions. Often this is done by generating several coarser resolution rasters from a fine resolution raster. In this paper we describe efficient algorithms for different variants of this problem.

The natural way of modeling water flow on a triangulated terrain is to make the fundamental assumption that water follows the direction of steepest descent (dsd). However, computing watersheds and other flow-related structures according to the dsd model in an exact manner is difficult: the dsd model implies that water does not necessarily follow terrain edges, which makes designing exact algorithms difficult and causes robustness problems when implementing them. As a result, existing software implementations for computing watersheds are inexact: they either assume a simplified flow model or they perform computations using inexact arithmetic, which leads to inexact and sometimes inconsistent results. We perform a detailed study of various issues concerning the exact or approximate computation of watersheds according to the dsd model. Our main contributions are the following.

Computing watersheds on triangulated terrain models in a robust manner is a difficult task as it is sensitive to noise that appears in the elevation values of the input. This is amplified by the existence of many very small watersheds (corresponding to spurious minima) that obscure the overall hydrological structure of the terrain. In the present work we perform an experimental evaluation of various algorithms that may help alleviate these problems:

Flow-related structures on terrains are defined in terms of paths of steepest descent (or ascent). A steepest descent path on a polyhedral terrain T with n vertices can have Θ(n2) complexity. The watershed of a point p—the set of points on T whose paths of steepest descent reach p—can have complexity Θ(n3). We present a technique for tracing a collection of n paths of steepest descent on T implicitly in O(n log n) time. We then derive O(n log n) time algorithms for: (i) computing for each local minimum p of T the triangles contained in the watershed of p and (ii) computing the surface network graph of T. We also present an O(n2) time algorithm that computes the watershed area for each local minimum of T.

We study the complexity of the visibility map of terrains whose triangles are fat, not too steep and have roughly the same size. It is known that the complexity of the visibility map of such a terrain with n triangles is θ(n2) in the worst case. We prove that if the elevations of the vertices of the terrain are subject to uniform noise which is proportional to the edge lengths, then the worst-case expected (smoothed) complexity is only θ(n). This provides an explanation why visibility maps of superlinear complexity are unlikely to be encountered in practice.