Parallel Processing for Geometric and Related Operations in Geospatial Analysis: A Survey
1. Introduction
Space and time are two significant aspects of data objects from domains such as geographical information science, social science, neuroscience, epidemiology, meteorology, transportation, and criminology. These data are collected using a variety of sensors and devices such as weather-capturing remote sensors, ground and ocean sensors, traffic cameras and sensors, medical imaging devices, handheld devices, social media, and simulations.
This study focuses on parallel algorithms used for efficient geospatial analysis using large volumes of data. Spatial and spatiotemporal data models and operations are the building blocks of spatial analysis [1]. Traditional relational database systems cannot efficiently handle this kind of data due to the high compute-intensive operations, high volume of data, and the different relationships among data objects. To mitigate these limitations, efficient data models and parallel algorithms that can harness high-performance computing has been introduced.
There are different data models used in geospatial analysis such as vectors and attributes, point clouds, raster and satellite imagery, census data, cell phone data, drawn images, and social media data [1], [2]. These data models have different characteristics and there are different operations performed on them to exploit useful information.
This survey discusses parallel algorithms used for geometric and related operations on geospatial vector data. Vector data has location information of points or collections of points in the form of lines and polygons. Geometric operations (intersection, union, and difference), spatial join, and overlay are key basic operations that are frequently used and important for many other operations. To handle large volumes of such data and to perform these operations efficiently parallel processing has been introduced by harnessing computing from multi-core, heterogenous, and distributed parallel hardware platforms.
In this work, we propose an extensive literature survey on parallel algorithms used for geometric and related operations in geospatial vector data. This study proposes to survey stateof-the-art parallel and scalable algorithms which can harness the computing of different hardware platforms handling large volumes of data. Geospatial data comes in mainly two types: spatial data and spatiotemporal data. Both have location information of an object(s) associated with other features information. Additionally, spatiotemporal data has a time stamp at each of those locations. Efficient data structures are used in both serial and parallel environments to improve the performance and scalability of the overall algorithms.
Data structures that are based on space division have shown improved efficiency with spatial data. They have advantages such as space efficiency, high query speeds, and amendable for parallelization [3]. Spatiotemporal data can be used with both SQL and NoSQL systems [4], [5]. Specialized data structures for this kind of data are few and have not been parallelized in all of the parallel computing platforms due to their complexities.
Uncertainty is another challenge when using location data. The locations are dependent on the GPS signals and loss of signals can lead to inefficiencies and errors. Especially, range queries performed on spatiotemporal data need models to provide correct results [4], [5]. This survey is beneficial to geographic information science communities and GIS data-related parallel computing researchers to learn novel algorithms and their usage.
A. Geospatial Data
Geospatial data refers to information that describes objects, events, or features with a specific location on or near the Earth's surface. It combines location information, which includes coordinates on the Earth, attribute information that describes the characteristics of the object or event, and temporal information that indicates the time or duration of the location and attributes. The location can be either static, such as the position of equipment, an earthquake event, children in poverty, or dynamic, like a moving vehicle, pedestrian, or the spread of an infectious disease [2].
This data type encompasses large sets of spatial information obtained from various sources in different formats. It can include data from censuses, satellite imagery, weather reports, cell phones, drawings, and social media. Geospatial data becomes most valuable when it can be discovered, shared, analyzed, and combined with traditional business data [2].
Common data types:
Geospatial data can be categorized into two major categories: spatial data and spatiotemporal data. Spatial data refers to a static location on the earth near the earth's surface. Spatiotemporal data carries time data along with spatial data adding a dynamic nature to the data. Geospatial data can also incorporate supplementary information such as elevation, color data, and many more. There are various data formats used to store these types of data.
-
Vector and attributes: Location descriptions in the forms of points, polylines, and polygons.
-
Raster and satellite imagery: High-resolution images taken by orbiting satellites.
-
Point clouds: Collections of points in higher dimensional space encoding additional data such as elevation, and color.
-
Census data: Census data connected to geographic locations or regions.
-
Cell phone data: Location and other data from cell phone call routes.
-
Drawn images: CAD images of man-made buildings and structures.
-
Social media data: Social media posts and videos that are used to identify emerging trends.
The rest of the paper is organized as follows. Section II discusses the serial polygon clipping and overlay algorithms. Section III presents parallel polygon clipping and overlay algorithms. Section IV presents a discussion on serial and parallel spatiotemporal data structures and algorithms used with polygonal data. Section V presents open research questions.
2. Serial Polygon Clipping and Overlay Algorithms
There are well-known polygon clipping sequential and parallel algorithms which are used in Computer Graphics and GIS domains. Maillot's algorithm can only clip a polygon against a rectangle, but not against polygons [6]. Sutherland-Hodgman, Weiler-Atherton, Vatti's, and Greiner-Hormann algorithms can clip a concave polygon against another concave polygon [7]-[10]. But they have limited capability when convex and selfintersections polygons are used. Vattiβs and Greiner-Hormann (GH) algorithms stand out since they can clip arbitrary polygons [9], [10].
There are other serial polygon clipping algorithms in the literature they are listed below. Liang et al. 1983 [11], Skala et al. 1993 [12], Hoffman et al. 1993 [13], Skala et al. 1994 [14], Karinithi et al. 1995s [15], Ng et al. 1998 [16], Zhang et al. 2002 [17], Peng et al. 2005 [18], Liu et al. 2017 [19], Eberly et al. 2008 [20], Qatawneh et al. 2009 [21], Martinez et al. 2009 [22], Feng et al. 2010 [23], Xie et al. 2010 [24], Sakala et al. 2012 [25], Foster et al. 2019 [26]. Evangeline et al. 2014 [27], Yang et al. 2016 [28],
Vatti's clipping uses a sweep-line approach. Its run time is sensitive to output size. GH algorithm brute forces allto-all-edge pairs to find the intersections and label them to identify the clipping regions. This runs in quadratic time. GH algorithm has a simpler way to represent polygons than Vatti's and its time complexity is not output sensitive. However, the GH algorithm is not able to handle degenerate cases, but Vatti's is able to do so. Perturbation has been used to solve this limitation of the GH algorithm, but it can produce ambiguous results based on the perturbation direction. Foster et al. [29] present a solution to overcome this limitation in the GH algorithm with the use of more advanced labeling and guarantees unique results.
There is a multicore Vatti's algorithm implementation presented in [30]. Parallel many-core and multicore implementations for GH algorithm are presented in [31]. However, this version is unable to handle degenerate cases. Parallel polygon clipping algorithm research lacks a simple algorithm that can handle complex situations such as degenerate cases in polygons. In the current literature, the existing systems are either single node sequential algorithms [29], heterogeneous algorithms which can not handle degenerate boundaries [31], or compute cluster based (CPU only) using MapReduce/MPI that do not harness GPUs [32].
3. Parallel Polygon Clipping and Overlay Algorithms
A. All-to-all Edge Pair Based Polygon Clipping Algorithms
This is a naive approach to clip two arbitrary polygons. The major phases in this approach are: finding intersections, labeling them, and tracing the resulting polygon where the dominating phase is intersection finding [11, 15, 39]. In the intersection-finding phase, all edge pairs of the subject and clipping polygons are tested against each other to find the intersections. This step runs in O(n2) time, where n is the number of edges in the input polygons. Labeling can vary depending on the algorithm and the tracing phase uses the labels the evaluate the edges of the input polygons to consider them in the result (see Fig. 3). The rules defined in the tracing phase yield the geometric operations: intersection, union, and set difference [11, 15, 38, 39].
Puri et al. [31] present two key findings. 1) the GreinerHormann algorithm has been successfully parallelized, which can perform polygon clipping in O(logn) time using O(n+k) processors. 2) the practical, parallel GIS system called MPIGIS that uses R-tree for efficient indexing for polygon overlay achieving a 44X speedup compared to state-of-art ArcGIS while processing 600K polygons in two GIS layers within 19 seconds on a 32-node NERSC's CARVER cluster.
Spark is a popular cloud computing paradigm for processing large-scale data due to its low latency and scalability. Zhao at al. [33] present a parallel processing algorithm based on the Spark paradigm that takes into account the shape complexity of polygons. The algorithm uses data partitioning, distributed spatial indexing, minimum boundary rectangular filtering, and other optimization techniques to improve overlay analysis speed while maintaining high speed and parallel efficiency. The major steps of the algorithm involve preliminary filtering using MBRs of the polygons, Hibert grid-based polygon distribution to partition data, R-tree indexing of the partitions, and final target polygon filtering. The experimental results demonstrate the high efficiency of this work.
B. Plane-sweep Based Polygon Clipping Algorithms
Serial plane-sweep approach can clip a pair of arbitrary polygons in O((π +π) log π) time. This approach exhibits better time complexity than the previous all-to-all edge comparison method. But parallelizing the phases in this approach is more challenging and requires parallel data structures [38]. [38] presents a multi-core implementation of this approach. [34] presents a GPU parallelized method for line segment intersection. But there is no GPU parallel algorithm based on this approach.
Puri et al [30] presents an algorithm that can perform polygon clipping in O(logn) time using (n+k+k') processors. The algorithm is the first output-sensitive CREW PRAM algorithm of its kind, which is more efficient than the previous algorithm developed by Karinthi, Srinivas, and Almasi. The new algorithm can handle self-intersecting polygons and is practical to use on multicores, achieving a 30x speedup for real-world datasets. It can perform typical clipping operations such as intersection, union, and difference.
Paudel et al. [34] present a parallel algorithm that can perform geometric intersection using OpenMP and OpenACC. The algorithm uses a compiler directive-based approach for parallelization, making it simpler to parallelize sequential code. This is the first work to demonstrate effective parallelization of plane sweep on GPUs. When using an Nvidia Tesla P100 GPU, the implementation achieves approximately 40 times the speed of the sequential CGAL library for the line segment intersection problem with 40K and 80K data sets.
C. Raster Based Polygon Clipping Algorithms
Zhou et al [35] introduced boundary algebra filling to rasterize input polygons and detect intersections in raster cells. Parallel strategies have been used based on workload decomposition and task scheduling to balance the workload across different parallel processes, resulting in effective acceleration of the polygon intersection process. This method is tested on datasets with over one million overlapping polygons, with a significant reduction in execution time and good scalability demonstrated. Reducing total runtime through task scheduling was found to be especially effective with a lower number of processes.
Gao et al. [36] propose a GPU-based rasterization algorithm to handle arbitrary closed polygons in Boolean operations. The algorithm uses efficient raster and vector data structures in the GPU and CPU respectively. The algorithm is based on a fast GPU-based contour extraction method, and a novel traversing strategy to ensure an error-free calculation of the intersection points. The algorithm is evaluated and shows higher performance than existing algorithms for processing polygons with a large number of vertices. However, rasterization and pair-wise clipping are two major limitations in this work.
Fan et al. [37] compare two polygon overlay analysis algorithms for handling overlap among polygons with different numbers of vertices. Traditional vector polygon clipping algorithms are less efficient when dealing with large quantities of vertices, while the RaPC algorithm, based on a discretization processing method and rasterization, is more efficient. The paper presents two types of parallel polygon overlay analysis algorithms and shows that the RaPC algorithm is advantageous for processing large data sets, while the Vatti algorithm is better for small data sets. The fine-grained parallelization method of the RaPC algorithm could provide a potential solution for the polygon overlay problem. However, due to the rasterization errors, the overall overlay process lacks accuracy.
D. Optimizations for Polygon Clipping and Overlay
Liu et al. [38] introduce a new geospatial filter called PolySketch-based CMBR (PSCMBR) filter. This filter uses a two-step process to detect workload and refine sketches. It speeds up line segment intersections (LSI) reporting tasks and reduces computational workload in PNP tests. The paper also integrates these new filters into the hierarchical filter and refinement system to solve geometric intersection problems. The new filters have been implemented on GPU using CUDA, resulting in a processing rate of 61 million/sec on average, making them more efficient than existing filters.
E. Parallel Data Structures Useful in Spatial Data Analysis
B-Tree performs better than the state-of-the-art, including a GPU log-structured merge tree and a GPU sorted array [39]. It supports concurrent queries and updates. Point and range queries are significantly faster, and B-Tree insertions are also faster unless insertions come in batches of more than roughly 100k. Caching the upper levels of the tree enables lookup throughput that exceeds the DRAM bandwidth of the GPU. Performance on a GPU is limited by contention, and the design choices made in this implementation overcome this issue.
Efficient R-tree construction and search algorithms on GPUs are in high demand due to the power of GPUs as computing platforms. However, this is a challenging problem due to the non-linear tree topology of the data structure and the unconventional architecture of modern GPUs. Prasad et al. [40] propose a space-efficient data structure design and a bottomup construction algorithm, resulting in a 226-fold speedup in parallel construction of an R-tree on a GPU compared to onecore execution on a CPU. This work also presents innovative R-tree search algorithms that can overcome GPU's limitations, resulting in a speedup of 91-fold to 180-fold on an R-tree with 16384 base objects for query sizes ranging from 2k to 16k.
Chaudhary et al. [41] demonstrate a novel method for implementing octree manipulations on a shared memory architecture, which is efficient for performing independent computations on each node of the octree. Octrees are hierarchical structures used to model 3-D objects, which are much easier to manipulate than other representations. The proposed method distributes operations among multiple processing elements (PEs) resulting in a uniform load balance, allowing for a significant speedup. This approach includes various manipulations such as complement, union, intersection, volume, and centroid calculations, which are implemented on the Sequent Balance multiprocessor. Experimental results confirm the feasibility of this approach and provide insights for evaluating other architectures.
4. Spatiotemporal Data
A. Data Structures
The G-tree data structure combines features of grids and Btrees. It's designed for organizing multidimensional data and is efficient in dynamic data spaces with a high frequency of insertions and deletions. The partitions are spatially ordered to perform operations like insertion, retrieval, deletion, and range queries efficiently. This work compares the G-tree with other trees such as the BD tree, zkdb tree, and KDB tree, and presents a discussion of the advantages of the G-tree over other data structures. Simulated bucket utilization rates for the G-tree are also reported [43].
The MP-tree is an extension of conventional spatial data structures using the techniques of the Persistent Search Tree. Teraoka et al. [44] describe the structure and operations of the MP-tree and demonstrates its superiority through computer simulations. For multi-version range searching, the MP-tree also performs better than the alternatives when the query time interval is short.
Ni et al. [45] propose a new approximation method that uses a single continuous polynomial to approximate a sequence of movement functions. The resulting index is called the PA-tree which yields much finer approximation quality than MBRs and has significantly lower construction costs. Experimental results show that PA-trees with clustered or non-clustered indices outperform MBR-based methods in spatiotemporal range queries by reducing I/O by 20%-60% and 300%-400%, respectively.
The TPR-tree of Tao et al. [46] is used for predictive spatiotemporal queries. However, the authors have discovered that the TPR-tree does not consider all the factors that affect query performance. They propose the TPR*-tree, a new index structure that considers the unique features of dynamic objects. This work presents cost models to determine the optimal performance achievable by any data-partition spatiotemporal access method. The experimental results show that the TPR*- tree is nearly-optimal and outperforms the TPR-tree in all conditions.
The TPR-tree is a data structure originally used for diskoriented environments. Nitkiewicz et al. [47] investigated its behavior in main memory. Their experiments verify the common belief that memory access time is a bottleneck for inmemory indexes due to its CPU-heavy and the limited number of computations that can improve its performance. This work introduces three improvement methods for update algorithms: lazy delete, bottom-up delete, and new penalty insert.
The velocity distribution of objects in an index tree determines the search space for a query. Nguyen et al. [48] proposed a velocity partitioning (VP) technique to speed up query processing using moving object indexes. VP identifies the dominant velocity axes using PCA and k-means clustering and creates a moving object index for each axis and maintains objects whose DVA is closest to their current direction. This technique reduces the search space expansion with time from a quadratic function to a near-linear function of maximum speed. The VP technique can be applied to various moving object index structures and consistently improves their performance, as validated by extensive experiments.
The TPR-tree uses VBR to save space by representing different positions of moving objects with one regional rectangle. This means that during a search, only VBR-based rectangles that overlap with a given query range are followed toward candidate leaf nodes. However, the TPR-tree suffers from dead space due to constant expansions of VBR-based rectangles, especially when stored in flash storage, which has expensive update costs. To address this issue, a new caching scheme based on Bloom filters is proposed, which can efficiently control the frequency of updates on a leaf node and improve the performance of indexing moving objects in flash storage [49].
B. Scalable Solutions
1) Database Based Approach:
Time-parameterized queries (TP queries) are queries that retrieve the actual result at the time of the query, the validity period of the result, and the change that causes the expiration of the result. This paper proposes a general framework that covers time-parameterized window queries, k-nearest neighbors, and spatial joins, reducing each of these to nearest neighbor search with defined distance functions. This allows for the application and extension of well-known branch and bound techniques using R-trees and their extensions for dynamic objects [50].
Raptopoulou et al. [51] propose an efficient method for processing nearest-neighbor queries in moving-object databases. A future query is made to find the set of objects that satisfy the query within a given time interval, which is difficult due to the continuously changing positions of both the query and data objects. The proposed approach issues only one query per time interval and uses a time-parameterized R-tree structure to index the objects. The performance evaluation shows significant improvements in CPU and I/O time compared to existing techniques based on either sampling or repetitive application of time-parameterized queries.
Xiong et al. [52] propose a new index called LUGrid that is designed to minimize the cost of object updates in spatiotemporal databases. LUGrid has two important features, lazy insertion and lazy deletion that reduce the I/Os required for updates and avoid deleting single obsolete entries immediately. It adapts to object distributions through cell splitting and merging. Theoretical analysis and experiments show that LUGrid outperforms previous work by up to eight times for intensive updates while still providing similar search performance.
Robert et al. [53] focuse on analyzing changes in spatial patterns over time. The authors propose five new movement events, including displacement, convergence, divergence, fragmentation, and concentration. They also introduce size and directional measures for these events in different time periods. The metrics are applied to a case study of a wildfire in northwestern Montana, showcasing the utility of the approach.
2) Grid Based Approach:
Xia et al. [54] propose a parallel Grid-Search-based Support Vector Machine (GS-SVM) optimization algorithm on Spark to predict passenger hotspots. The algorithm employs the grid search approach to optimize the radial basis function of the support vector machine and the cross-validation method to find the global optimal parameter combination. The proposed algorithm is successfully applied to predict passenger hotspots and outperforms several state-ofthe-art algorithms, including ARIMA, SVR, LSTM, and CNN, with a lower MAPE value of at least 78.4%.
3) Grid with B-tree Approach:
Jung et al. [55] propose a new index that combines TPR*-tree and hash map to quickly update and query target objects in real-time. When there are many target objects collected by multiple sensors, the update performance of the TPR*-tree needs improvement. The paper introduces a non-leaf insertion policy that includes target objects and MBRs to reduce the time spent on MBR expansion in TPR*-tree. The proposed hybrid index is tested through experiments and is proven to perform well in updating and querying target objects in real time.
4) Indexing Based Approach:
Teng et al. [56] propose a real-time contact tracing system to deal with the challenges of processing enormous volumes of data. The system is based on spatial proximity queries with temporal constraints, and it uses location data for dynamic indexing of moving objects. The system is optimized for parallelism and provides an efficient contact evaluation mechanism to ensure only spatially and temporally valid contacts are kept. The experiments show that the system can achieve sub-second level response for largescale contact tracing with a significant performance boost over CPU-based approaches.
Prabhakar et al. [42] propose a solution that combines Query Indexing and Velocity Constrained Indexing (VCI) to efficiently and scalably evaluate multiple continuous queries on moving objects. Query Indexing uses incremental evaluation, reverses the role of queries and data, and exploits the relative locations of objects and queries. VCI takes advantage of the maximum possible speed of objects to delay the costly operation of updating an index to reflect the movement of objects. The proposed solution outperforms traditional approaches by almost two orders of magnitude in experimental evaluations.
5) Cloud Computing Based Systems:
To support the increasing demand for data storage, Big Data requires horizontally scalable data stores. This need led to the development of various NoSQL databases, each with distinct query and persistence methodologies. One of the most massive data types being gathered is spatiotemporal data, including location data. To achieve effective query and transformation semantics, a unique spatiotemporal index structure is described, which uses the horizontal scalability of NoSQL databases. Performance data obtained from testing with Accumulo is presented [57].
6) Other Scalable Systems:
Large amounts of data on taxi trajectories can be quickly collected using advanced positioning and communication technologies. However, handling such large amounts of data is still challenging for data users. To address this issue, a software system called XStar has been developed, which has a scalable index and storage structure. This allows for more efficient storage and access to raw data. The system has been tested using real taxi trajectory data and has proven to be effective in facilitating the processing and analysis of such data. XStar has already been downloaded over 4000 times since its release in January 2019 and new analytical functions are currently being developed [58].
C. Uncertainty in Spatiotemporal Data
Saltenis et al. [61] propose a new indexing technique based on R-tree. This technique supports efficient querying of the current and future positions of moving objects in 1D, 2D, and 3D spaces. With updated algorithms, the index can manage a dynamic data set, where objects can come and go, and changes occur in the expected positions of the objects. Furthermore, the paper reports a thorough performance analysis of the technique.
Lou et al. [62] present a new map-matching algorithm named ST-Matching, which is designed for GPS trajectories with low-sampling rates. ST-Matching combines spatial geometry and topology of the road network with the speed and temporal constraints of the trajectories to create a candidate graph. The best matching path sequence is then identified from this graph. The algorithm is compared to the incremental algorithm and AFD-based global map-matching algorithm on synthetic and real datasets. Results indicate that ST-Matching outperforms the incremental algorithm in terms of accuracy for low-sampling trajectories and also improves accuracy and running time compared to the AFD-based global algorithm.
Spatial query problems with uncertain spatiotemporal data are complex, especially when a time component is involved. There are currently no efficient solutions available for managing such data, despite its wide range of applications. Bernecker et al [63] propose general models for spatiotemporal uncertain data that have the potential to enable the efficient processing of various queries. The main challenge is to develop new algorithms based on these models, which will unlock their potential. The paper also includes examples of interesting spatiotemporal queries on uncertain data.
5. Open Research Problems
Plane-sweep sequential polygon clipping algorithms run in O((n+k)logn) time. There are multi-core implementations of such algorithms. But there are no efficient GPU algorithms in the literature. Parallel plane-sweep clipping algorithms require parallel data structures and it is a challenging task.
Efficient parallel data structures are the key to achieving better run times. However, parallel data structure and query algorithms choices for spatiotemporal data are limited. Current literature heavily depends on no-SQL databases and cloud computing solutions for handling such data over scalable systems.
Uncertainty is a key aspect when using spatial and spatiotemporal data. The literature lacks integration of uncertainty with parallel algorithms.