grokkingstuff/thesis
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Stealing Kenohori's thesis to use as a template for my thesis report

Browse Source
This commit is contained in:
Vishakh Kumar 2025-07-14 08:17:30 +04:00
parent dd6ef5922a
commit 878300738b
313 changed files with 14900 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

20
work_report/.gitignore vendored Normal file
View File

@ -0,0 +1,20 @@
thesis.toc
thesis.synctex.gz
thesis.sublime-workspace
thesis.sublime-project
thesis.pdf
thesis.lot
thesis.log
thesis.fls
thesis.fdb_latexmk
thesis.brf
thesis.blg
thesis.bbl
thesis.aux
.DS_Store
propositions.synctex.gz
propositions.pdf
propositions.log
propositions.fls
propositions.fdb_latexmk
propositions.aux

117
work_report/LICENSE Normal file
View File

@ -0,0 +1,117 @@
CC0 1.0 Universal
Statement of Purpose
The laws of most jurisdictions throughout the world automatically confer
exclusive Copyright and Related Rights (defined below) upon the creator and
subsequent owner(s) (each and all, an "owner") of an original work of
authorship and/or a database (each, a "Work").
Certain owners wish to permanently relinquish those rights to a Work for the
purpose of contributing to a commons of creative, cultural and scientific
works ("Commons") that the public can reliably and without fear of later
claims of infringement build upon, modify, incorporate in other works, reuse
and redistribute as freely as possible in any form whatsoever and for any
purposes, including without limitation commercial purposes. These owners may
contribute to the Commons to promote the ideal of a free culture and the
further production of creative, cultural and scientific works, or to gain
reputation or greater distribution for their Work in part through the use and
efforts of others.
For these and/or other purposes and motivations, and without any expectation
of additional consideration or compensation, the person associating CC0 with a
Work (the "Affirmer"), to the extent that he or she is an owner of Copyright
and Related Rights in the Work, voluntarily elects to apply CC0 to the Work
and publicly distribute the Work under its terms, with knowledge of his or her
Copyright and Related Rights in the Work and the meaning and intended legal
effect of CC0 on those rights.
1. Copyright and Related Rights. A Work made available under CC0 may be
protected by copyright and related or neighboring rights ("Copyright and
Related Rights"). Copyright and Related Rights include, but are not limited
to, the following:
i. the right to reproduce, adapt, distribute, perform, display, communicate,
and translate a Work;
ii. moral rights retained by the original author(s) and/or performer(s);
iii. publicity and privacy rights pertaining to a person's image or likeness
depicted in a Work;
iv. rights protecting against unfair competition in regards to a Work,
subject to the limitations in paragraph 4(a), below;
v. rights protecting the extraction, dissemination, use and reuse of data in
a Work;
vi. database rights (such as those arising under Directive 96/9/EC of the
European Parliament and of the Council of 11 March 1996 on the legal
protection of databases, and under any national implementation thereof,
including any amended or successor version of such directive); and
vii. other similar, equivalent or corresponding rights throughout the world
based on applicable law or treaty, and any national implementations thereof.
2. Waiver. To the greatest extent permitted by, but not in contravention of,
applicable law, Affirmer hereby overtly, fully, permanently, irrevocably and
unconditionally waives, abandons, and surrenders all of Affirmer's Copyright
and Related Rights and associated claims and causes of action, whether now
known or unknown (including existing as well as future claims and causes of
action), in the Work (i) in all territories worldwide, (ii) for the maximum
duration provided by applicable law or treaty (including future time
extensions), (iii) in any current or future medium and for any number of
copies, and (iv) for any purpose whatsoever, including without limitation
commercial, advertising or promotional purposes (the "Waiver"). Affirmer makes
the Waiver for the benefit of each member of the public at large and to the
detriment of Affirmer's heirs and successors, fully intending that such Waiver
shall not be subject to revocation, rescission, cancellation, termination, or
any other legal or equitable action to disrupt the quiet enjoyment of the Work
by the public as contemplated by Affirmer's express Statement of Purpose.
3. Public License Fallback. Should any part of the Waiver for any reason be
judged legally invalid or ineffective under applicable law, then the Waiver
shall be preserved to the maximum extent permitted taking into account
Affirmer's express Statement of Purpose. In addition, to the extent the Waiver
is so judged Affirmer hereby grants to each affected person a royalty-free,
non transferable, non sublicensable, non exclusive, irrevocable and
unconditional license to exercise Affirmer's Copyright and Related Rights in
the Work (i) in all territories worldwide, (ii) for the maximum duration
provided by applicable law or treaty (including future time extensions), (iii)
in any current or future medium and for any number of copies, and (iv) for any
purpose whatsoever, including without limitation commercial, advertising or
promotional purposes (the "License"). The License shall be deemed effective as
of the date CC0 was applied by Affirmer to the Work. Should any part of the
License for any reason be judged legally invalid or ineffective under
applicable law, such partial invalidity or ineffectiveness shall not
invalidate the remainder of the License, and in such case Affirmer hereby
affirms that he or she will not (i) exercise any of his or her remaining
Copyright and Related Rights in the Work or (ii) assert any associated claims
and causes of action with respect to the Work, in either case contrary to
Affirmer's express Statement of Purpose.
4. Limitations and Disclaimers.
a. No trademark or patent rights held by Affirmer are waived, abandoned,
surrendered, licensed or otherwise affected by this document.
b. Affirmer offers the Work as-is and makes no representations or warranties
of any kind concerning the Work, express, implied, statutory or otherwise,
including without limitation warranties of title, merchantability, fitness
for a particular purpose, non infringement, or the absence of latent or
other defects, accuracy, or the present or absence of errors, whether or not
discoverable, all to the greatest extent permissible under applicable law.
c. Affirmer disclaims responsibility for clearing rights of other persons
that may apply to the Work or any use thereof, including without limitation
any person's Copyright and Related Rights in the Work. Further, Affirmer
disclaims responsibility for obtaining any necessary consents, permissions
or other rights required for any use of the Work.
d. Affirmer understands and acknowledges that Creative Commons is not a
party to this document and has no duty or obligation with respect to this
CC0 or use of the Work.
For more information, please see
<http://creativecommons.org/publicdomain/zero/1.0/>

6
work_report/README.md Normal file
View File

@ -0,0 +1,6 @@
# My PhD thesis
Some related blog posts:
* [How to self-publish your (LaTeX) thesis](https://3d.bk.tudelft.nl/ken/en/2016/03/20/self-publishing-your-latex-thesis.html)
* [What is my LaTeX workflow](https://3d.bk.tudelft.nl/ken/en/2016/04/02/my-latex-thesis-workflow.html)
* [How my 1.5 column layout works](https://3d.bk.tudelft.nl/ken/en/2016/04/17/a-1.5-column-layout-in-latex.html)

616
work_report/cleaning.tex Normal file
View File

@ -0,0 +1,616 @@
% !TEX root = thesis.tex
\chapter{Processing real-world datasets into clean geometric models}
\label{ch:cleaning}
The representations described in \refpa{pa:representation} and the operations described in \refpa{pa:operations} work well on perfectly \emph{valid} data.
Objects are assumed not to overlap each other, to be properly closed with no degenerate geometries, and to have perfectly planar faces which are consistently oriented, among many other validity criteria.
Unfortunately, as \refse{se:repair-motivation} explains, GIS processes often fail to clearly specify which criteria are expected and real-world data often fails to meet them, with consequences ranging from the innocuous to a complete inability to use a desired tool, including instances where software gives erroneous results unbeknownst to the user.
As GIS datasets can be rather intricate and expensive to acquire, they are not easily replaceable and must therefore be \emph{repaired}, \ie\ they must be processed into geometric models that conform to certain validity specifications, so as to make them fit for use.
In the context of this thesis, clean geometric models are important as they are the base for the higher-dimensional models using the representations and operations described in previous chapters.
The following sections thus describe the background and a particular solution used in this thesis to obtain valid polygons and planar partitions in \refse{se:pprepair}, and valid polyhedra and 3D space partitions in \refse{se:3drepair}.
The chapter concludes with a generalisation of the definition of validity in arbitrary dimensions in \refse{se:ndrepair}, which can be used for both higher-dimensional objects and space partitions.
\refse{se:pprepair} is largely based on the papers:
\begin{itemize}
\paperpfgpprepair%
\papercgeoprepair%
\end{itemize}
\section{Motivation}
\label{se:repair-motivation}
The representations described in \refpa{pa:representation} are each able to effectively represent a particular class of objects.
For example, most data structures are intended for 3D space partitions whose 3D objects have surfaces that form 2-manifolds, and $n$D combinatorial maps are capable of representing subdivisions of orientable quasi-manifolds, which within this thesis are generally further limited to having only linear geometries.
Similarly, the operations described in \refpa{pa:operations} can only return good results when the input fulfils certain requirements.
For instance, the extrusion operation from \refch{ch:extrusion} requires the input data to form a space partition, the incremental construction algorithm from \refch{ch:incremental-construction} requires the ridges of the facets of an object to form matching pairs (\ie\ a quasi-manifold), and the linking approach from \refch{ch:linking-lods} will only form a valid 4D cell complex with 4-cells if a matching scheme that preserves all topological relationships between cells can be found or is provided.
All these operations are also based on the assumption that the input cells themselves are valid, and are being completely bounded by valid lower-dimensional cells (\ie\ facets, ridges, etc.) so as to form a valid cell complex.
The representations and operations presented in this thesis are not special in this sense.
Validity assumptions are widely used in all software, especially when complex data is used as input, and are used to make many tasks more manageable, such as to interpret a dataset and load it into a particular data structure for internal usage, as well as for further operations that may be performed using this structure.
However, GIS datasets are more complex than most other data\footnote{For instance, unlike other types of datasets, GIS objects often cannot be stored directly as a plain list of tuples, but are instead decomposed into primitives of a certain shape (\refse{se:spatial-modelling}), sometimes recursively, and these have to be defined in terms of its dimension and structure, topological relations, geometry and attributes, sometimes including rich semantics as well.}, and GIS software thus tends to make more assumptions than most other software.
Moreover, though making sure that these assumptions are true is highly desirable, testing for every possible invalid configuration in a spatial dataset is cumbersome and often unnecessary for a particular task at hand, while testing for only some invalid configurations depending on what needs to be done can easily become intractable and can result in a large number of redundant tests.
Running these validation tests can also be difficult and computationally expensive, as many tests take longer to execute than some of the common tasks that a GIS is used for (\eg\ visualisation of a dataset, simple statistical analyses, or checking the attributes of some objects).
At the same time, the datasets found in the GIS world are very diverse and their properties vary significantly.
They can be generated using a large variety of GIS, CAD and 3D modelling software based on different processes, complying to different specifications and stored in different formats, each of which follows its own internal logic and structure.
Most importantly, GIS datasets are created for different purposes or meant for general-purpose applications (\eg\ CityGML \citep{CityGML2}).
As such, a dataset's specifications are often only vaguely defined, are defined only with regard for a particular application, or are not conformed with in practice.
\footnotetext[141]{\url{http://www.esri.com/software/arcgis/}}
\footnotetext[142]{\url{http://www.qgis.org}}
\footnotetext[143]{\url{http://www.safe.com/fme/}}
\footnotetext[144]{\url{http://grass.osgeo.org/}}
It is thus perhaps unsurprising that invalid GIS datasets are prevalent \citep[Ch.~7]{Panigrahi14} and a major source of problems for those who work with them.
As shown in \reffig{fig:partially}, invalid datasets can be interpreted inconsistently in different software, leading to inconsistent or erroneous results.
They can also make it impossible to perform a certain operation, either due to a failing precondition check or due to software crashes.
% \footnote{Admittedly based on anecdotal evidence. However, it seems to be a consensus opinion.}.
A partial solution to this issue lies in \emph{validating} these datasets, \ie\ identifying the problematic objects in the data so that they can be discarded or (manually) fixed.
There are a variety of checklists and (semi) automatic tools for this purpose, such as those provided by SAFE\footnote{\url{http://blog.safe.com/2014/11/data-quality-checklist/}} and ESRI\footnote{\url{http://www.esri.com/software/arcgis/extensions/arcgis-data-reviewer/~/media/Files/Pdfs/library/fliers/pdfs/arcgis-data-reviewer-checks.pdf}}.
\begin{figure}[tbp]
\centering
\subfloat[]{
\includegraphics[height=0.35\linewidth]{figs/part-arcgis}
\label{fig:part-arcgis}}
\quad
\subfloat[]{
\includegraphics[height=0.35\linewidth]{figs/part-grass}
\label{fig:part-grass}}
\caption[Different interpretations of a polygon]{Different interpretations of a polygon with a hole that is partly outside its outer boundary ($p_3$ in \reffig{fig:p}).
(a) ArcGIS\protect\footnotemark\ considers the overlapping region as a hole, but the non-overlapping part of the hole as a new polygon (QGIS\protect\footnotemark\ and FME\protect\footnotemark\ do this as well).
(b) GRASS\protect\footnotemark\ removes the overlapping part from the polygon, becoming a new polygon with a different shape. No warning is shown in any software.}
\label{fig:partially}
\end{figure}
For example, it is possible to incorporate a set of formal preconditions for every operation, as is done in the design by contract software engineering pattern \citep{Meyer86} or the Eiffel programming language \citep{ISO/IEC25436:2007}.
However, simply discarding problematic (subsets of) datasets is not always feasible, as GIS datasets can be expensive to acquire and thus irreplaceable in practice, and manually fixing errors can be an extremely time-consuming process.
In fact, according to \citet{McKenney98}, users of 3D CAD models for finite element analysis---which has similar requirements as certain computations in GIS, such as well-shaped and non-overlapping mesh elements---spend up to 70\% of their time fixing the input CAD models.
While similar figures for GIS are to the best of my knowledge not available, it is worth noting that CAD software tends to produce better quality models than GIS software\footnote{There are many reasons for this. For instance, CAD software makes wider use of topological data structures, and also has topology-aware and smart interactive editing tools (\eg\ snapping to guide lines and nearby objects), which help to avoid problems where objects seem to be valid but have small errors, such as sliver polygons and shells that are not properly closed.}.
A more complete solution therefore lies in using methods to automatically \emph{repair} a dataset, \ie\ to make it conform to a particular set of validity criteria, enabling the full use of many more datasets than would otherwise be possible.
As this requires a rather complex defensive programming approach, testing for various types of partially overlapping cascading errors and fixing them accordingly, it is best performed separately and called as needed rather than integrated into every operation of a GIS.\@
The following sections describe such independent repair methods as were used in this thesis, which involve the creation of valid (multi)polygons and planar partitions from 2D GIS datasets (\refse{se:pprepair}), and the creation of valid polyhedra and 3D space partitions from 3D BIM datasets (\refse{se:3drepair}).
These were then used as input for the different experiments described in \refpa{pa:operations}.
\section{Creating valid (multi)polygons and planar partitions}
\label{se:pprepair}
\subsection{What is a valid (multi)polygon or planar partition?}
In most GIS file formats and the software that reads and writes them, polygons and multipolygons are defined in a manner that is consistent with the definitions in the Simple Features Specification \citep{SimpleFeatures1,ISO19125-1:2006}---an implementation of the ISO 19107 standard \citep{ISO19107:2005}.
The specification states that: `\emph{A \texttt{Polygon} is a planar \texttt{Surface} defined by 1 exterior boundary and 0 or more interior boundaries. Each interior boundary defines a hole in the \texttt{Polygon}}'.
Each of these boundaries is described as a \texttt{LinearRing} (\cf\ \reffigp{fig:sfs}).
According to the specification, an outer ring should be oriented \emph{anticlockwise} when viewed from a predefined \emph{top} direction, which is generally (but not necessarily) the viewing direction in 2D or \emph{outwards} when the polygon specifies part of the boundary of a polyhedron.
Inner rings should be oppositely oriented, \ie\ generally \emph{clockwise} when viewed from the top direction.
The Simple Features Specification provides several \textbf{validity rules for polygons}, the most relevant of which are described below with examples of invalid polygons provided in \reffig{fig:p}.
The rules can be summarised as follows:
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/unitpolygons}
\caption[Several invalid polygons]{Several invalid polygons, with their outer boundaries shown in black and their inner boundaries in grey. The orange areas represent one possible interpretation of the interior of the polygons. Polygon $p_{12}$ has an exterior and an interior boundary with the same geometry.}
\label{fig:p}
}
\begin{itemize}
\item each ring defining the exterior and interior boundaries is \emph{simple}, \ie\ non-self-intersecting ($p_{1}$ and $p_{10}$);
\item each ring is closed ($p_{11}$), \ie\ its first and its last points should be the same;
\item the rings of a polygon do not cross ($p_{3}$, $p_{7}$, $p_{8}$ and $p_{12}$), but they may intersect at one tangent point;
\item a polygon does not have cut lines, spikes or punctures ($p_{5}$ and $p_{6}$);
\item the interior of every polygon is a connected point set ($p_{4}$);
\item each interior ring creates a new area that is disconnected from the exterior ($p_{2}$ and $p_{9}$).
\end{itemize}
Similarly, the specification provides a definition and some \textbf{validity rules for multipolygons}.
A \texttt{MultiPolygon} is defined as a \texttt{MultiSurface} forming an aggregation of \texttt{Polygons}, which also follows certain validity criteria, which can be summarised as follows:
\begin{itemize}
\item the interiors of its polygons do not overlap, \ie\ their point set intersection should be empty;
\item the boundaries of its polygons may only touch at a finite number of points;
\item a multipolygon does not have cut lines, spikes or punctures;
\item the interior of a multipolygon with more than one polygon is \emph{not} a connected point set.
\end{itemize}
Intuitively, a \emph{planar partition} is a set of polygons that form a subdivision of a region of the plane.
Planar partitions are thus commonly used to model concepts where objects are expected not to overlap, such as land cover, cadastral parcels, or the administrative boundaries of a given country.
Despite being a very frequently used representation in GIS, planar partitions are not explicitly defined in the main GIS standards.
Within the classes in the ISO 19107 standard \citep[\S{}6.6]{ISO19107:2005}, a planar partition could be considered as a \texttt{GM\_CompositeSurface}, defined in the standard as `\emph{a collection of oriented surfaces that join in pairs on common boundary curves and which, when considered as a whole, form a single surface}'.
By following this definition, overlaps between polygons are explicitly forbidden, as a \texttt{GM\_Complex} (a parent of \texttt{GM\_CompositeSurface}) is defined as `\emph{a set of primitive geometric objects (in a common coordinate system) whose interiors are disjoint}'.
However, a \texttt{GM\_CompositeSurface} explicitly allows gaps between the surfaces, as these would simply result in inner rings within the overarching single surface.
An alternative definition could be created based on the ISO 19123 standard \citep[\S{}6.8]{ISO19123:2007}---a standard focusing on coverages of various types.
According to the standard, a planar partition can be considered as a type of \texttt{CV\_DiscreteSurfaceCoverage} where `\emph{the surfaces that constitute the domain of a coverage are mutually exclusive and exhaustively partition the extent of the coverage}'.
Overlapping polygons are disallowed by them being `\emph{mutually exclusive}' and gaps are disallowed by the surfaces `\emph{exhaustively partitioning}' the extent.
However, the standard states these conditions as something that occurs `\emph{in most cases}', whereas in a planar partition it should be considered as a strict prerequisite.
In a \textbf{valid planar partition}, there should thus be no overlapping polygons, and no gaps between them either unless these gaps are considered to be outside of the region.
These two conditions are covered by the ISO 19107 standard in a different context, when it lists some possible inconsistencies of `spaghetti' datasets represented as a \texttt{GM\_Complex}, stating that `\emph{slivers and gaps are multiple lines that should represent the same geometry, but do not coincide, leaving areas of overlap between two surface boundaries (slivers), and gaps between them}' \citep[\S{}6.2.2.6]{ISO19107:2005}.
\subsection{Commonly used validation and repair methods}
Starting from an arrangement of line segments that are meant to define the boundary of a (multi)polygon, there are various rules that can be used to define its interior and exterior.
\citet{Foley95} discusses two commonly used sets of rules in vector-based graphic software, which are shown in \reffig{fig:polygoninterpretation}\footnote{\citet{Foley95} considers only polygons, but the rules as explained here also cover holes in polygons and multipolygons.}.
\marginpar{
\captionsetup{type=figure}
\centering
\subfloat[]{\includegraphics[width=\marginparwidth]{figs/oddeven}} \\
\subfloat[]{\includegraphics[width=\marginparwidth]{figs/nonzerowinding}}
\caption[Rules used to interpret the interior of a polygon]{(a) According to the \emph{odd-even rule}, a polygon's interior is the region(s) that can be accessed by passing through an odd number of edges from its exterior. (b) According to the \emph{non-zero winding rule}, a polygon's interior is the region(s) around which the boundaries of a polygon do a non-zero number of revolutions.}
\label{fig:polygoninterpretation}
}
In practice, GIS users often repair invalid polygons manually.
Among the few documented automatic solutions, it is possible to use a `buffer-by-0' operation \citep{Ramsey10} or PostGIS 2.0's \texttt{ST\_MakeValid} function\footnote{\url{http://postgis.org/documentation/manual-svn/ST_MakeValid.html}}.
The validation of planar partitions is usually performed using a checklist of individual tests that together ensure its validity.
For instance, \citet{Plumer97} specify that a valid planar partition consists of: no dangling edges, no zero-length edges, planarity, no holes, no self-intersections, no overlaps, and having a connected graph.
However, it is worth noting that without the use of a supporting structure, some of these tests can be problematic or computationally expensive.
For instance, checking whether any possible pair of polygons overlap can have quadratic behaviour even when heuristics to speed up the process are used \citep{Badawy99,Kirkpatrick00}, and robustness issues are significant in polygon intersection tests \citep{Hoffmann88}.
Finding the potential gaps in a planar partition is also a problem, as it can require computing the union of the entire set of polygons \citep{Margalit89,Rivero00}.
The most common method used to repair a planar partition is based on the assumption that polygons \emph{approximately} match each other at their common boundaries.
If the adjacent polygons are within a certain distance \emph{threshold} of each other along their common boundaries (\reffig{subfig:threshold-gap}), and all parts further apart than this threshold are known not to be common boundaries (\reffig{subfig:threshold-overlap}), it is possible to \emph{snap} together the polygons using this threshold, while in theory keeping the rest untouched.
This method of planar partition repair is available in many GIS packages, including ArcGIS, FME, GRASS and Radius Topology\footnote{\url{http://www.1spatial.com/software/radius_topology/}}.
\marginpar{
\captionsetup{type=figure}
\centering
\subfloat[]{\includegraphics[scale=0.36]{figs/threshold-gap}\label{subfig:threshold-gap}} \\
\subfloat[]{\includegraphics[scale=0.36]{figs/threshold-overlap}\label{subfig:threshold-overlap}}
\caption[Defining a snapping threshold]{Defining a snapping threshold taking into account: (a) gaps and (b) overlaps.}
\label{fig:threshold}
}
The threshold value for certain input dataset(s) is then usually manually determined, either by trial and error, or by analysing certain properties of the datasets involved (\eg\ point spacing, precision, or map scale).
However, it is often hard to find an optimal threshold for certain datasets, and sometimes impossible as such a threshold does not even exist (\eg\ because point spacing in some places might be smaller than the width of the gaps and overlaps present).
\subsection{Repair using a constrained triangulation}
The method developed to repair polygons and planar partitions uses a constrained triangulation of the input polygons as a base structure.
Constrained triangulations have distinct properties that make them useful as a base for a repair algorithm.
The can be built efficiently with a variety of approaches \citep{Guibas85,Clarkson92}, can easily be made numerically robust\footnote{Based on one robust geometric predicate that tests whether three successive points are collinear, have a clockwise or a anticlockwise orientation (\eg\ \citet{Shewchuk96a}) and the computation of new vertices at the intersections of line segments.
A constrained Delaunay triangulation only requires an additional predicate that determines whether a point lies inside, on or outside the circle defined by three other points.}, can be used for quick traversal and point location \citep{Mucke99}, and have fast and robust implementations in several libraries, such as CGAL \citep{Boissonnat02}, Triangle\footnote{\url{https://www.cs.cmu.edu/~quake/triangle.html}} \citep{Shewchuk97} and GTS\footnote{\url{http://gts.sourceforge.net}}.
The method used to \textbf{repair individual (multi)polygons} exploits these properties and consists of three broad steps, which are shown in \reffig{fig:prepair-workflow} and described as follows:
\begin{enumerate}
\item construction of the constrained triangulation of the line segments in the input, processing outer and inner rings identically and including an extra edge that connects the first and last vertices of a ring when these are not the same (\reffig{subfig:prepair-workflow-3});
\item labelling of each triangle as either \emph{outside} or \emph{inside}, which is based on an extension of the odd-even rule that supports overlapping lines by adding or removing (parts of) constraints in the triangulation (\reffig{subfig:prepair-workflow-4}), taking only edges that are \emph{constraints} into account;
\item reconstruction of the interior areas as a repaired multipolygon\footnote{The output might only be representable as a multipolygon even if the input was a polygon} (\reffig{subfig:prepair-workflow-5}).
\end{enumerate}
\begin{figure}[b]
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-workflow-1}\label{subfig:prepair-workflow-1}}
% \subfloat[]{\includegraphics[width=0.2\linewidth]{figs/prepair-workflow-2}\label{subfig:prepair-workflow-2}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-workflow-3}\label{subfig:prepair-workflow-3}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-workflow-4}\label{subfig:prepair-workflow-4}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-workflow-5}\label{subfig:prepair-workflow-5}}
\caption[Steps to repair a (multi)polygon using a constrained triangulation]{Steps to repair a (multi)polygon using a constrained triangulation: (a) input data, (b) triangulation, (c) labelling, and (d) reconstruction.}
\label{fig:prepair-workflow}
\end{figure}
This method is remarkably efficient, and its implementation based on CGAL classes is able to process large polygons quickly.
As an example, \reffig{fig:prepair} shows the process on the largest polygon in the CORINE land cover dataset\footnote{\url{http://www.eea.europa.eu/publications/COR0-landcover}}, which consists of almost 1 189 903 vertices and 7 672 holes, which is processed in under a second.
\begin{figure*}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-1}\label{subfig:prepair-1}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-2}\label{subfig:prepair-2}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-3}\label{subfig:prepair-3}}
\subfloat[]{\includegraphics[width=0.25\linewidth]{figs/prepair-4}\label{subfig:prepair-4}}
\caption[Processing the largest polygon in the CORINE land cover dataset]{Processing the largest polygon in the CORINE land cover dataset: (a) outer and inner boundaries, (b) triangulation, (c) labelling, (d) reconstruction.}
\label{fig:prepair}
\end{figure*}
The method to \textbf{repair a planar partition} then uses polygons which are known to be valid based on the previous method.
It consists of four main steps, shown in \reffig{fig:pprepair-workflow}, and is as follows:
\begin{enumerate}
\item the constrained triangulation of the input segments forming the (now valid) polygons is constructed;
\item each triangle is labelled with the labels of the polygons inside which it is located, such that problems are detected by identifying triangles having no or multiple labels;
\item problems are fixed by re-labelling triangles according to customisable criteria, such that each triangle has exactly one label;
\item the polygons are reconstructed from the triangulation.
\end{enumerate}
\begin{figure}
\subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-a}\label{subfig:prepair-workflow-a}}
\subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-b}\label{subfig:prepair-workflow-b}}
\subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-c}\label{subfig:prepair-workflow-c}}
% \subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-d}\label{subfig:prepair-workflow-d}}
\subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-e}\label{subfig:prepair-workflow-e}}
\subfloat[]{\includegraphics[width=0.2\linewidth]{figs/pprepair-workflow-f}\label{subfig:prepair-workflow-f}}
\caption[Steps to repair a planar partition using a constrained triangulation]{Steps to repair a planar partition using a constrained triangulation: (a) input data, (b) triangulation, (c) labelling, (d) relabelling problematic triangles, (e) reconstruction.}
\label{fig:pprepair-workflow}
\end{figure}
Various local repair methods can thus be defined, all of which are based on choosing how to relabel a triangle which no label or with multiple labels.
\reffig{fig:prepair-method} shows a few examples of such methods.
\begin{figure*}
\subfloat[]{\includegraphics[height=22mm]{figs/pprepair-method-a}\label{subfig:prepair-method-a}}
\quad
\subfloat[]{\includegraphics[height=22mm]{figs/pprepair-method-b}\label{subfig:prepair-method-b}}
\quad
\subfloat[]{\includegraphics[height=22mm]{figs/pprepair-method-d}\label{subfig:prepair-method-d}}
\quad
\subfloat[]{\includegraphics[height=22mm]{figs/pprepair-method-c}\label{subfig:prepair-method-c}}
\caption[Various repair methods based on relabelling (sets of) triangles]{Various repair methods can be defined based on relabelling triangles or sets of connected triangles.
For instance, based on (a) the input data, it is possible to relabel: (b) an invalid triangle based using its longest boundary with a neighbour, (c) an invalid region of connected triangles using its longest boundary with a neighbour, or (d) an invalid region of connected triangles using a random neighbour.}
\label{fig:prepair-method}
\end{figure*}
Within this thesis, planar partitions are obtained by repairing invalid regions using the longest boundary with a neighbour (\reffig{subfig:prepair-method-d}) if possible, as this tends to produce a cartographically more pleasing result, and a random neighbour (\reffig{subfig:prepair-method-c}) otherwise.
This method is able to process large datasets quickly, such as the one shown in \reffig{fig:16tiles}, which consists of 16 tiles of the CORINE land cover dataset in a 4$\times$4 configuration.
It has 63 868 polygons with a total of 6 622 133 vertices and was processed in 4 minutes 47 seconds.
By comparison, both ArcGIS and GRASS are unable to repair this dataset by snapping geometries, while FME repairs it in 15 minutes 48 seconds\footnote{On tests where ArcGIS and GRASS were able to process the data, they were slower than FME.}, also using snapping.
Note that apart from the fact that snapping is slower, it does not guarantee a valid result.
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/16tiles}
\caption[A planar partition of 16 tiles of the CORINE land cover dataset]{A planar partition made from 16 tiles of the CORINE land cover dataset.}
\label{fig:16tiles}
}
\section{Creating valid polyhedra and 3D space partitions}
\label{se:3drepair}
\subsection{Motivation: IFC input data}
As discussed in \refse{ss:formats}, the data models used in 3D GIS have often favoured an approach where objects that are volumetric in reality are modelled as a set of their (visible) surfaces, which can be easily captured after the objects are built.
In addition, the volumetric objects that are small, elongated or thin are sometimes not modelled as volumes, but they are instead respectively modelled as points, curves or surfaces.
For instance, in 3D models encoded in CityGML \citep{CityGML2}, all objects are modelled as surfaces that do not necessarily form closed volumes, and thin objects (\eg\ roof overhangs, windows and doors) are often modelled as single surfaces.
The volumetric models that are often used in CAD and BIM, such as those in IFC format \citep{ISO16739:2013}, follow a different approach, which is shown in \reffig{fig:volumetric}.
In it, almost all real-world volumetric objects, \ie\ the objects with a non-zero volume, are modelled as volumes as well.
This approach is more expensive in terms of space and makes certain computations more difficult, such as obtaining the volume of a room that is only represented implicitly by a set of (volumetric) walls, doors and windows around it.
However, the approach is ultimately a more powerful representation that is closer to reality, enabling more complex operations, such as the structural analysis of a building, and eliminates the ambiguities inherent in interpreting volumetric real-world objects that have been modelled as points, curves or surfaces.
\begin{figure}[tb]
\centering
\includegraphics[width=\linewidth]{figs/volumetric}
\caption[Surface-based models vs.\ volumetric models]{A surface-based model (left) consists of a set of semantically labelled surfaces, which are shown here in different colours. A volumetric model (centre+right) instead consists of a set of semantically labelled \emph{volumes}.}
\label{fig:volumetric}
\end{figure}
Although the volumes in such a model do not generally fit together perfectly, they can nevertheless be better processed\footnote{Akin to how a set of polygons can be better processed into a planar partition compared to a rough line drawing with undershoots and overshoots (\ie\ a typical spaghetti dataset).} into a 3D space partition consisting of a set of non-overlapping 3D objects, where each of the volumes---and optionally also their vertices, edges and faces---is labelled with appropriate semantic information\footnote{As shown in \refse{se:duality}, these lower-dimensional cells can store important information about the relationships between the objects.}.
The result is thus a \emph{space-partitioning 3D model}, which is analogous to the 2D planar partitions often used to model coverages in GIS.\@
Such a model can be stored in a computer using the topological data structures for 3D (\refse{ss:data-structures}) or $n$D (\refse{ss:cell-complexes}) cell complexes.
This makes a volumetric model \emph{a better base for a higher-dimensional representation}, as it fully exploits the properties of higher-dimensional data structures and their operations.
For instance, taking the examples of this thesis, it can be extruded into higher dimensions (\refch{ch:extrusion}) or multiple such 3D models can be linked into a single 4D model (\refch{ch:linking-lods}).
A space-partitioning model is also better able to take advantage of the complex 3D models that are already created during the design and construction processes of a building, such as the faces that form a room, storey or building.
This avoids the need to extract or abstract appropriate surfaces for their recreation in a GIS model, such as was done by \citet{Donkers13}.
Since architectural and BIM models have a strong emphasis on representing individual 3D elements that need to be designed, manufactured and put together, such as those shown in \reffig{fig:fzk}, as well as the fact that these elements generally do not overlap in reality (except as part of hierarchies), they are ideally suited to be used in a space-partitioning model.
\begin{figure}[tbp]
\centering
\includegraphics[width=\linewidth]{figs/ifc}
\caption[An IFC model of the FZK-house]{An IFC model of the FZK-house\footnotemark\ contains a complex representation of the structural elements of the building. Note how the volumes in this model fit together, having no overlaps even in places that are normally not visible (\eg\ the beam ends which are embedded in the walls).}
\label{fig:fzk}
\end{figure}
\footnotetext{\url{http://iai-typo3.iai.fzk.de/www-extern/index.php?id=1174&L=1}}
\subsection{What is a valid solid or 3D space partition?}
The ISO 19107 standard \citep[\S{}6.3.18]{ISO19107:2005} defines 3D objects with 3D holes that are known as \emph{solids}, which are specified based on a boundary representation scheme.
As shown in \reffig{fig:iso-brep}, the standard thus defines a \texttt{GM\_Solid} with a \emph{boundary} operation returning a \texttt{GM\_SolidBoundary}, which is a `\emph{sequence of sets of \texttt{GM\_Surfaces} that limit the extent of [the] \texttt{GM\_Solid}}'.
Each of these sets of surfaces describes one of the boundaries of the \texttt{GM\_Solid} as a \texttt{GM\_Shell}, corresponding to either the outer boundary for the solid\footnote{In some cases, there might not be an outer boundary of a solid, such as in non-Euclidean spaces or in the representation of unbounded solids. However, there is nearly always an outer boundary in the context of geographic information.} or one of its holes.
\begin{figure}[tbp]
\centering
\includegraphics[width=\linewidth]{figs/iso-brep}
\caption[The boundary representation scheme in ISO 19107]{The ISO 19107 standard \citep[\S{}6.3.2]{ISO19107:2005} is able to specify the boundaries of \texttt{GM\_Curve}, \texttt{GM\_Surface} and \texttt{GM\_Solid} as subclasses of \texttt{GM\_Boundary}, respectively a \texttt{GM\_CurveBoundary} linked to a pair of \texttt{GM\_Point} (the end-points of a line segment), \texttt{GM\_SurfaceBoundary} linked to a set of instances of \texttt{GM\_Ring}, and a \texttt{GM\_SolidBoundary} linked to set of instances of \texttt{GM\_Shell}.}
\label{fig:iso-brep}
\end{figure}
A \texttt{GM\_Shell} \citep[\S{}6.3.8]{ISO19107:2005} thus represents `\emph{a single connected component of a \texttt{GM\_SolidBoundary}}'.
It is known to be \emph{simple}, and consists of a set of oriented instances of \texttt{GM\_Surface} composed of instances of \texttt{GM\_SurfacePatch}, which intuitively form a cellular subdivision of the surface and themselves have a \texttt{GM\_SurfaceBoundary}.
A \texttt{GM\_SurfaceBoundary} represents an area potentially with any number of holes, each of which is stored as a reference to a \texttt{GM\_Ring}.
A \texttt{GM\_Ring} \citep[\S{}6.3.6]{ISO19107:2005} is additionally defined as being \emph{simple}.
\texttt{GM\_Object} \citep[\S{}6.2.2]{ISO19107:2005}, a parent class to all the classes previously mentioned, defines every object as a point set and provides the definition of simple as a `\emph{GM\_Object [that] has no interior point of self-intersection or self-tangency.
In mathematical formalisms, this means that every point in the interior of the object must have a metric neighborhood whose intersection with the object is isomorphic to an $n$-sphere, where $n$ is the dimension of this \texttt{GM\_Object}}'.
As discussed by \citet{Ledoux13}, this implies that shells are effectively 2-manifolds.
Rings are similarly 1-manifolds.
It is important to note that even though each \texttt{GM\_Ring} and \texttt{GM\_Shell} is individually simple, the boundary of the \texttt{GM\_Surface} or \texttt{GM\_Solid} that they together describe does not need to be simple.
A common example would involve an inner ring/shell tangent to the outer ring/shell containing it.
Arguably, the standard does appear to explicitly forbid intersections between the interior of rings or shells as \texttt{GM\_Complex} is a parent class of \texttt{GM\_SurfaceBoundary} and \texttt{GM\_SolidBoundary} and this class requires its composing primitives to be `\emph{geometrically disjoint}'.
However, this interpretation is problematic as it would arguably also forbid inner rings being inside their containing outer ring.
Alternatively, it is possible to consider that the standard does not specify \emph{any} restrictions regarding the interactions between rings of a surface or between shells of a solid.
As the standard explicitly states that `\emph{implementations may enforce stronger restrictions on the interaction of boundary elements}', it might be the responsibility of other implementing standards to place appropriate restrictions.
Although the GML standard \citep{GML3.2.1} implementing ISO 19107 does not specify such restrictions, it is possible to use those defined in the Simple Features Specification in 2D (\refse{se:pprepair}) and define analogous ones in 3D \citep{Ledoux13}.
One possible formulation of these could be as follows:
\begin{itemize}
\item the shells of a solid do not cross, but the shells on the boundary of a solid may intersect only at a vertex or edge;
\item the interior of every solid is a connected point set;
\item each interior shell creates a new volume that is disconnected from the exterior.
\end{itemize}
Intuitively, a 3D space partition is a subdivision of a region of 3D space into non-overlapping solids.
However, just as with planar partitions, 3D space partitions are usually not strictly defined.
Following the same logic as with planar partitions in \refse{se:pprepair}, a 3D space partition can be considered as an ISO 19107 \texttt{GM\_CompositeSolid} \citep[\S{}6.6.13]{ISO19107:2005}, which is defined in the standard as a `\emph{a set of solids that join in pairs on common boundary surfaces to form a single solid}'.
While overlapping solids are explicitly forbidden by a \texttt{GM\_CompositeSolid} inheriting from \texttt{GM\_Complex} in which `\emph{[primitive] interiors are disjoint}', gaps between the solids are explicitly allowed.
An alternative definition could also be created based on the ISO 19123 standard by considering a 3D space partition as a type of \texttt{CV\_DiscreteSolidCoverage} \citep[\S{}6.10]{ISO19123:2007}, which states that `\emph{generally, the solids that constitute the domain of a coverage are mutually exclusive and exhaustively partition the extent of the coverage}'.
While overlaps and gaps are respectively eliminated by the `\emph{mutually exclusive}' and `\emph{exhaustively partition}' conditions, the word `\emph{generally}' implies that these are not always enforced.
\subsection{Common problems for 3D objects in an IFC file}
An IFC model in theory consists of a set of 3D objects that mostly do not overlap.
However, there are a few common problems that cause some of these 3D objects to be invalid or prevents them from forming a clean space partition.
The most significant of these are the following:
\begin{figure*}[tbp]
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-diff-w}\label{subfig:ioh-diff-w}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-diff-wo}\label{subfig:ioh-diff-wo}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-diff-o}\label{subfig:ioh-diff-o}}
\caption[Openings in the IfcOpenHouse dataset]{In the (a) IfcOpenHouse dataset\protect\footnotemark, the openings where the windows and door fit are not explicitly carved out from the wall volumes.
Instead, the dataset consists of (b) walls with simple shapes and (c) the openings themselves as large boxes.
A Boolean point set difference thus needs to be computed in order to obtain the final explicit geometries.}
\label{fig:ioh-diff}
\end{figure*}
\footnotetext{\url{http://blog.ifcopenshell.org/2012/11/say-hi-to-ifcopenhouse.html}}
\begin{description}
\item[Implicit geometries] Objects are often represented using sweeps, intersections of half-spaces and Boolean set operations, such as in the case of openings as shown in \reffig{fig:ioh-diff}.
These need to be converted to explicit (boundary representation) objects that can be made to fit with other objects.
As these need to be discretised, they will often not perfectly match the shape of the original objects, sometimes creating problems in their interaction with other objects.
\item[Local coordinate systems] Objects are defined using a local coordinate system, which might differ per object.
As the parameters of these coordinate systems are stored (and computed) in a finite computer representation, this will cause the objects not to fit together perfectly.
In addition, applying a transformation to embed all objects into a unique coordinate system using computed arithmetic will cause additional problems.
The result is that objects that visually appear to fit together do not actually do so, having small gaps and overlaps between them as shown in \reffig{fig:ioh-small}.
\item[Hidden intentional overlaps] Not caring about hidden object intersections is common practice in most 3D modelling approaches.
By not caring about these intersections, it is possible to ease and speed up the modelling process by using simpler volumes (\eg\ boxes or rectangles) than would otherwise be required and hiding their undesirable parts behind or inside other objects, as is shown in \reffig{fig:ioh-large}.
\end{description}
\begin{figure*}[tbp]
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-1}\label{subfig:ioh-1}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-2}\label{subfig:ioh-2}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-3}\label{subfig:ioh-3}}
\caption[The volumes in the IfcOpenHouse dataset do not fit together]{The volumes in the (a) IfcOpenHouse dataset do not perfectly match each other at their common boundaries.
(b) The two volumes forming the roof actually do not touch, causing the house not to be closed.
(c) The left wall and the foundation volumes have a small overlapping portion (\ie\ the wall sinks inside the foundation).}
\label{fig:ioh-small}
\end{figure*}
\begin{figure*}[tbp]
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-4}\label{subfig:ioh-4}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-5}\label{subfig:ioh-5}}
\enspace%
\subfloat[]{\includegraphics[width=0.325\linewidth]{figs/ioh-6}\label{subfig:ioh-6}}
\caption[A large overlap in the IfcOpenHouse dataset]{(a) The right wall and the right roof volumes in the IfcOpenHouse dataset visually appear to match each other, but the (b) right wall and (c) right roof are modelled as parallelepipeds, so they have a non-empty intersection (a triangular pyramid).}
\label{fig:ioh-large}
\end{figure*}
\subsection{Commonly used validation and repair methods}
In order for software to be able to process invalid polyhedra, various automatic validation and repair methods have been developed, detecting and/or fixing some of the possible invalid configurations that can exist.
The possible invalid configurations are many and partly overlap or cascade (from lower to higher dimensions), but a possible list can be generated by systematically exploring those that occur at the ring, polygon, shell and solid levels\footnote{In theory, errors can occur at the point and edge level, such as a lack or excess of coordinates, but on a typical GIS input file these would be generally considered as syntactic rather than geometric errors.}.
For instance, \citet{Wagner13} tests for 13 types of invalid configurations while \citet{Ledoux13} tests for 26, notably including: consecutive points in a ring with the same coordinates, rings that are not closed or self-intersect, polygons with intersecting rings or with an inner ring outside the outer ring, shells that are not closed or are not 2-manifold, and solids with intersecting shells.
\textbf{Repairing individual rings and polygons} can be basically done using the processes outlined in \refse{se:pprepair}, even if these are embedded in 3D rather than 2D.
For this, rings and polygons can be first projected onto a certain plane (\eg\ the best-fitting one, or one obtained by disregarding a coordinate of its vertices in a manner that does not create new degenerate shapes) or a restricted triangulation \citep[Ch.~13]{Cheng12} could be used as a basis for a repair procedure.
\textbf{Repairing shells and solids} are problems that are also partly related to those discussed previously, as
a shell can be seen as a planar partition that wraps around an object, \ie\ it is watertight, while the constraints that define how the inner and outer shells of a valid solid should interact are similar to those defining the interaction of inner and outer rings within a valid polygon.
However, the methods that have been presented previously are much less applicable to shells and solids.
Instead, a shell can often be seen as a mesh that \emph{should} be closed, and it is thus possible to repair individual shells with the procedures used for surface reconstruction and mesh repair.
There are various good surveys of the methods that can be used to repair various problems in polygonal meshes, such as \citet{Ju09} and \citet{Attene13}.
Some of these problems and their respective methods are summarised below.
However, it is important to notice that many methods are intended for meshes representing the boundary of a single `smooth' 2-manifold, making them not applicable to the great majority of BIM and GIS models where perpendicular angles are common (\eg\ those between walls and floors/ceilings).
\citet{Rossignac99a} propose a method to make polygonal meshes combinatorial manifolds, determining an appropriate order for a set of edges around a face or for a set of faces bounding a volume such that non-manifolds can be stored using a manifold data structure (\reffigp{fig:nonmanifold-012}).
\citet{Gueziec01} treat the problem from a geometric point of view, converting non-manifold edges into thin volumes by cutting and joining the mesh around them.
\citet{Attene09} propose a method to make tetrahedral meshes\footnote{assuming a tetrahedron-based data structure with adjacency relationships to other tetrahedra} manifold---something that can be used also for solids by computing their constrained tetrahedralisation.
In a similar manner as in GIS, small gaps in a mesh (causing a shell not to enclose any space) can be naively repaired by snapping vertices \citep{Rock92}, but this requires an error-prone threshold and can lead to topological errors.
Iteratively snapping boundary edges together works better \citep{Sheng95}, as it is possible to start from a single corresponding pair of edges, and then iteratively `zip' together corresponding edges that are adjacent to these.
\citet{Barequet95} follows a related approach, matching certain edges and triangulating the remaining gaps.
\citet{Turk94} shows how overlapping triangular meshes\footnote{As overlaps are embedded in 3D, they need to be defined from a given point of view.} can be fixed by clipping all but one of a set of overlapping triangles, retriangulating them afterwards.
Holes in a mesh, which tend to be bigger than gaps and reflect missing parts of a surface (\eg\ the bottom of a house, the sides of a terraced house, or a surface that is hidden from a typical point of view), require different methods.
In many cases, a hole is close to planar and can thus be simply projected to 2D and triangulated \citep{Bohn92}.
However, in other cases the holes are far from planar, and it is thus necessary to use more complex methods.
For instance, \citet{Levy03} fills holes while attempting to minimise a certain objective function representing the energy needed to fill it, \citet{Wang07} uses moving least squares fitting of the points around it\footnote{Proposed by \citet{Lancaster81}, it is a widely used surface reconstruction method based on interpolating a set of points.} and \citet{Podolak05} does so by subdividing space into regions deemed to be completely in or out of the shell.
\citet{Nooruddin03}, \citet{Bischoff05} and \citet{Hetroy11} use voxel-based methods to attempt to determine the interior and exterior of a shell and thus close a mesh.
Specifically in the context of buildings, which have different characteristics from many other meshes such as sharp corners and orthogonal surfaces, \citet{Bogdahn10} and \citet{Alam13} propose two methods that create a smooth 2-manifold, but do not guarantee their results.
\citet{Zhao14} attempts to solve gaps, holes and overlaps in a building mesh simultaneously by using a constrained tetrahedralisation \citep{Si05} of a set of faces, progressively carving away tetrahedra that are deemed as not belonging to its interior based on a set of rules.
The possible intersections between the faces are explicitly computed by the constrained tetrahedralisation, so that the starting tetrahedra form a 3D space subdivision.
Finally, the problem of creating a valid 3D space partition is loosely related to mesh simplification and more closely to the \emph{topological reconstruction} of a 3D model, which sometimes deals with computation of topological relationships from imperfect datasets.
Generally, the latter methods work by snapping the geometries that lie within a threshold.
For instance, \citet{Horna06} snaps together generalised map darts that lie within a threshold $\varepsilon$, connecting them by the appropriate involution $\alpha$, which is based on the steps of the reconstruction process.
\subsection{Repair using snapping and Boolean set operations}
Despite the fact that the methods presented above fix many of the problems in individual 3D objects, these methods are not always sufficient to create a valid 3D space subdivision.
For instance, the methods that are based on snapping are unable to create a space partition when adjacent geometries are farther than the snapping threshold or their overlapping regions have a width larger than the threshold.
Similarly, most methods to repair meshes do not guarantee that they form properly enclosed spaces and cannot deal well with solids with inner shells.
A different technique that solves many of these issues was thus developed for this thesis.
It starts from a set of separate 3D objects that \emph{approximately} form a 3D space partition, such as those that are commonly found in IFC building models.
The method starts by handling every object separately, obtaining a valid interpretation of each that is stored in an exact representation.
It then snaps objects together in order to remove small gaps and overlaps at their common boundaries.
Finally, it uses Boolean set operations to remove larger overlaps and in order guarantee that a 3D space partition is obtained.
Unlike most repair methods used in GIS, it tries to avoid triangulating every face of an object but still manages to obtain perfectly planar faces (in memory).
These steps required are described in detail as follows:
\begin{enumerate}
\item
\textbf{Rough individual polyhedra}
One or more `rough' polyhedral representations of every object are first extracted from the input model.
These rough representations should be combinatorially valid quasi-manifolds, being composed of patches that join in pairs at their common boundaries.
However, at this stage the polygonal patches can have geometric issues, such as not being planar and might intersect geometrically.
In order to obtain the rough polyhedral representations, every object is first individually parsed, converting implicit geometries into explicit geometries using Boolean set operations (\reffig{fig:ioh-diff}) and triangulating curved surfaces.
A transformation is then applied to convert every object's coordinates to a global coordinate system.
Finally, the faces of every object are extracted, object by object, and used to incrementally construct a set of polyhedra.
This is done by starting from a given face and attempting to add adjacent faces until a closed shell is formed, which is repeated until all faces of an object are processed or the remaining faces cannot form any closed shells.
\item
\textbf{Clean individual polyhedra}
Clean individual polyhedra are then created from the rough polyhedra.
In order to do this, the best fitting plane to the vertices of each face is first computed using linear least squares (\reffig{fig:faceplanes}) and is stored as a plane equation.
Considering these planes as \emph{constraints}, all the vertices of every polyhedron are moved to an exact intersection of as many as possible of its incident face planes using a greedy algorithm unless this would result in a too large shift (as defined by a threshold).
If some plane constraints could not be met, the corresponding faces are triangulated, thus becoming perfectly planar.
The planes of these new triangular faces are computed.
\begin{figure}[tbp]
\subfloat[]{\includegraphics[width=0.45\linewidth]{figs/faceplanes-1}}
\quad
\subfloat[]{\includegraphics[width=0.45\linewidth]{figs/faceplanes-2}}
\caption[Computing the best fitting plane of every face]{The best fitting plane to the vertices of every face is computed, here showing those belonging to (a) the floor, back right wall, back left wall and back eave, and (b) the front left wall, front right wall and front eave.}
\label{fig:faceplanes}
\end{figure}
\item
\textbf{Snapping vertices together}
The vertices of the polyhedra that lie within a threshold are snapped together, which removes most of the small gaps and overlaps in the model (\reffig{fig:snapping}).
This applies to vertices belonging to the same or to different polyhedra.
The snapped vertices are now considered \emph{immovable}.
Iterating through all the faces of the polyhedra, if a face has at least three immovable non-collinear vertices, the plane passing through these vertices is computed and it is considered as \emph{fixed}.
When there are more than three non-coplanar vertices, the face is triangulated and the faces with three immovable vertices are also considered as \emph{fixed}.
\marginpar{
\captionsetup{type=figure}
\centering
\subfloat[]{\includegraphics[width=\marginparwidth]{figs/snapping-1}}\\
\subfloat[]{\includegraphics[width=\marginparwidth]{figs/snapping-2}}
\caption[Vertex snapping]{(a) The vertices of different polyhedra that lie within a threshold are snapped together, thus (b) removing a small gap.
A small overlap works in the same manner.}
\label{fig:snapping}
}
\item
\textbf{Snapping vertices to fixed planes}
The vertices that are still considered movable are then snapped to nearby fixed planes, if any, eliminating certain other small gaps and overlaps that do not have vertices in common (\eg\ the steps in front of the IfcOpenHouse shown in \reffig{fig:ioh-steps}, which are actually not touching the house's foundation).
For this, iterating through every movable vertex, if it is incident to three or more faces with non-coplanar fixed planes and their intersection lies within a threshold of the vertex's current position, the vertex is moved to the intersection of three of these planes and considered as immovable.
If it has more than three incident faces with non-coplanar fixed planes, the faces of the planes that were not used, and are thus now not perfectly planar, are triangulated.
This step can be repeated a given number of times, increasing the number of immovable vertices.
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/ioh-steps}
\caption[The steps of the IfcOpenHouse]{The steps of the IfcOpenHouse do not actually touch the house's foundation.
Note that it also does not have common vertices with the foundation, so this gap cannot be closed by vertex-to-vertex snapping, but it can be closed by snapping its vertices on the left to the foundation's right plane.}
\label{fig:ioh-steps}
}
\item
\textbf{Fixing the remaining vertices}
The remaining movable vertices are fixed to their incident faces' fixed planes or to their current location.
For this, iterating through every movable vertex, the same procedure as the step above is followed.
However, if a vertex has less than three incident faces with non-coplanar fixed planes, the vertex is fixed to the position on the intersection of its incident faces' fixed planes that is closest to the vertex's current position.
These moved vertices are also considered as fixed and their incident faces' planes are recomputed if necessary.
\item
\textbf{Creating individual Nef polyhedra}
A Nef polyhedron is created from every polyhedral representation using the precomputed planes for each face.
Note that the exact representations of each plane (in the form of plane equations) are thus kept in this process.
\item
\textbf{Add the structural types}
The Nef polyhedra representing structural types (\eg\ walls, slabs and beams) are incrementally added to a model, making sure that a new Nef polyhedron does not intersect the previously added ones.
This is done using a Boolean set difference with a Nef polyhedron containing all previously added polyhedra, which is then regularised.
\item
\textbf{Remove the opening types}
The Nef polyhedra representing openings are carved out from the structural types by a Boolean set difference whose result is then regularised.
They are also carved out from the Nef polyhedron representing the entire model.
\item
\textbf{Add the fixture types}
The Nef polyhedra representing fixtures (\eg\ windows, doors, railings and frames) are incrementally added to the model in the same manner as the structural types.
The fixtures will often fit into the openings carved out in the previous step.
\end{enumerate}
The output of these steps is thus a list of Nef polyhedra representing each object using exact arithmetic.
As these are regularised, they are known to contain their boundary.
Because of this, the common faces of a pair of adjacent polyhedra are easy to obtain through the computation of their Boolean set intersection, which can be implemented quickly.
These faces can therefore be used to easily compute a topological representation of the 3D space subdivision.
This 3D repair method was implemented with the help of several libraries: IfcOpenShell is used to parse the IFC file, Open CASCADE is used to triangulate implicit representations and to transform the objects' coordinates to a global coordinate system, CGAL \texttt{Polyhedron\_3} is used to construct a half-edge representation of every polyhedron, and CGAL \texttt{Nef\_polyhedron\_3} is used to store every Nef polyhedron and to perform the Boolean set operations between them.
\section{Dimension-independent validity criteria}
\label{se:ndrepair}
The previous sections have expanded on the criteria that define what is a valid object or space partition of objects in 2D and 3D.
They also described some methods that can be used to make real-world data comply with these criteria, enabling the data to be used for more applications.
As this thesis aims at utilising real-world higher-dimensional data, it is important to also consider what criteria can be used to define validity in higher-dimensional data.
The standards for geographic information in 2D and 3D described previously (Simple Features \citep{SimpleFeatures1}, GML \citep{GML3.3} and ISO 19107 \citep{ISO19107:2005}) are in theory limited to 2D and 3D.
Concretely, the ISO 19107 standard explicitly states that `\emph{this International Standard is restricted to at most three dimensions}'.
However, all of these standards are easily extensible to higher dimensions.
This would mostly involve the addition of new classes and corresponding definitions.
However, the standards do contain minor hard-coded assumptions that are only valid for the 2D and 3D cases, such as how ISO 19107 and GML consider orientable curves and surfaces, but not orientable solids (\reffigp{fig:gml}).
This section therefore defines higher-dimensional objects in a manner that is (mostly) harmonious with the standards used in the GIS world.
An $n$-cell can be thus represented by the set of $(n-1)$-cells in its (outer) boundary, using a similar mechanism as how other boundaries are represented in the ISO 19107 standard, which was shown previously in \reffig{fig:iso-brep}.
\begin{figure}
\centering
\includegraphics[width=0.7\linewidth]{figs/iso-cell-brep}
\caption[A dimension-independent cell harmonised with ISO 19107]{A dimension-independent definition of a cell in a harmonised manner with other classes in the ISO 19107 standard \citep{ISO19107:2005}.}
\label{fig:iso-cell-brep}
\end{figure}
Following the terminology used in the standard and as shown in \reffig{fig:iso-cell-brep}, such an extension of the would mainly entail a \texttt{GM\_OrientableGeometricPrimitive} with a \texttt{dimension} attribute, which would set to $n$.
This class would be analogous to \texttt{GM\_OrientableCurve} for \texttt{dimension} 1, \texttt{GM\_OrientableSurface} for \texttt{dimension} 2 and a newly created \texttt{GM\_OrientableSolid} for \texttt{dimension} 3, which would be a subclass of \texttt{GM\_OrientablePrimitive}.
The \texttt{GM\_OrientableGeometricPrimitive} would be bounded by a \texttt{GM\_GeometricPrimitiveBoundary}, which would be linked to aggregations of $(n-1)$-dimensional instances of a newly created \texttt{GM\_Cell} (with their \texttt{dimension} attribute set to $n-1$).
This \texttt{GM\_Cell} class would be analogous to \texttt{GM\_Point} for \texttt{dimension} 0, \texttt{GM\_Curve} for \texttt{dimension} 1, \texttt{GM\_Ring} for \texttt{dimension} 2, and \texttt{GM\_Shell} for \texttt{dimension} 3.
Each of the instances of \texttt{GM\_Cell} bounding a \texttt{GM\_OrientablePrimitive} would represent either the outer boundary of the geometric primitive (if any), or one of any number of inner boundaries representing $n$-dimensional holes.
This extension of the standard would seem to follow most in the spirit of ISO 19107.
However, other alternative extensions could be considered.
As \texttt{GM\_Curve}, \texttt{GM\_Surface}, and \texttt{GM\_Solid} would essentially be special cases of \texttt{GM\_Cell}, all of the former could be seen as redundant and eliminated.
However, the standard already contains many specialisations that are somewhat redundant but that cover common use cases in geographic information, such as \texttt{GM\_Triangle} and \texttt{GM\_Tin}.
Another possibility would be considering \texttt{GM\_Curve}, \texttt{GM\_Surface}, and \texttt{GM\_Solid} as subclasses of \texttt{GM\_Cell} or substituting the abstract \texttt{GM\_Primitive} for a non-abstract \texttt{GM\_Cell}, but this would involve a major change in the standard and seems to run counter to the preferred use of abstract top classes in the standard.
The definition of an \texttt{GM\_OrientableGeometricPrimitive} as explained above also lends itself to the definition of sets of disjoint cells (akin to the \texttt{Multi\ldots}\ classes in the standard) and cell complexes (akin to the \texttt{Composite\ldots}\ classes in the standard), which could also be handled in the same manner as in the ISO 19107 standard.
As shown in \reffig{fig:iso-composites}, the standard already defines composite curves, surfaces and solids, which are equivalent to 1-, 2- and 3-dimensional cell complexes.
A \texttt{GM\_CompositeCurve} is `\emph{a list of orientable curves (\texttt{GM\_OrientableCurve}) agreeing in orientation in a manner such that each curve (except the first) begins where the previous one ends.}', a \texttt{GM\_CompositeSurface} is `\emph{a collection of oriented surfaces that join in pairs on common boundary curves}', and a \texttt{GM\_CompositeSolid} is `\emph{a set of solids that join in pairs on common boundary surfaces}'.
\begin{figure}
\centering
\includegraphics[width=\linewidth]{figs/iso-composites}
\caption[Cell complexes in the ISO 19107 standard]{The cell complexes of dimension 1, 2 and 3 are respectively defined in the ISO 19107 standard \citep[\S{}6.6.3]{ISO19107:2005} as the classes \texttt{GM\_CompositeCurve}, \texttt{GM\_CompositeSurface} and \texttt{GM\_CompositeSolid}.}
\label{fig:iso-composites}
\end{figure}
A similarly defined \texttt{GM\_CompositeGeometricPrimitive}, shown in \reffig{fig:iso-compositecell} which should contain the \texttt{dimension} as a parameter, would thus be equivalent to a representation of a space partition of any dimension that allows objects with holes.
It could be defined as `a set of $n$-dimensional orientable geometric primitives (\texttt{GM\_OrientableGeometricPrimitive}) that join in pairs on common $(n-1)$-dimensional boundary geometric primitives'.
Note that this implies that the primitives combinatorially form an $n$-quasi-manifold, although geometrically they might not do so due to the presence of holes.
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/iso-compositecell}
\caption[An $n$D space subdivision harmonised with ISO 19107]{A definition of an $n$-dimensional space subdivision in a harmonised manner with other classes in the ISO 19107 standard \citep{ISO19107:2005}.}
\label{fig:iso-compositecell}
}
Following the validity criteria previously described in \refse{se:pprepair} and \refse{se:3drepair}, it is possible to define additional validity criteria for an $n$-dimensional geometric primitive, which would serve to specify the conditions upon which its bounding cells may interact.
These would be as follows:
\begin{itemize}
\item the bounding $n$-cells of an $n$-dimensional geometric primitive do not cross, but they might intersect only at a cell of dimension $n-1$ or lower;
\item the interior of every geometric primitive is a connected point set;
\item each interior $n$-cell creates a new point set in $\mathbb{R}^n$ that is disconnected from the exterior.
\end{itemize}
Meanwhile, an $n$-dimensional space subdivision should consist of a set of $n$-dimensional geometric primitives that are mutually exclusive and exhaustively partition an extent, itself a well-defined subset of $\mathbb{R}^n$.
As with the definitions of a planar partition and 3D space subdivision, this implies that there should be no overlapping primitives, and no gaps between them unless these gaps are considered to be outside the extent.

272
work_report/conclusions.tex Normal file
View File

@ -0,0 +1,272 @@
% !TEX root = thesis.tex
\chapter{Conclusions and future work}
\label{ch:conclusions}
% Introduction / about topic?
There are a great number of possible representations for spatial information, each of which with its own benefits and drawbacks and only some of which have been discussed in this thesis.
From an engineering point of view, choosing an appropriate representation for a given system or task is a fundamental issue, as this choice cascades down to almost every engineering decision and indirectly affects a GIS program's every functionality.
Among other aspects, it affects the type of objects that can be efficiently stored and the operations that can be easily performed on them, and it also has important computational consequences in terms of both memory and processing time.
The 2D representations that dominate the GIS world work well for 2D datasets and problems that are essentially two-dimensional, but many problems arise when they are adapted to model 3D, spatiotemporal and multi-scale geographic information.
Most research in GIS is devoted to improving these adaptations, as well as to developing new methods that build on them to solve problems both old and new.
This thesis pushes GIS research in a different direction, starting from the assumption that many issues in GIS can probably be better solved by using a new, fundamentally different modelling approach---modelling both spatial and non-spatial characteristics as dimensions in the geometric sense, thus using higher-dimensional representations to create, manipulate and visualise geographic information.
Accordingly, this thesis' main research objective was to \textbf{realise the fundamental aspects of a higher-dimensional Geographic Information System}, and therefore focusing on the development of higher-dimensional representations and methods for GIS.\@
This concluding chapter starts with a short outlook on higher-dimensional GIS in \refse{se:isitworthit}, describing concisely
% some of the more enticing possibilities that can be opened up by higher-dimensional models in a GIS, and
when and where it makes sense to use higher-dimensional models.
Afterwards, \refse{se:conclusions} describes in detail the lessons learned by pursuing this thesis' research objective.
Finally, \refse{se:contributions} lists the main contributions of this thesis, and \refse{se:futurework} discusses the topics that I think would be most useful for future research on this topic.
\section{An outlook on higher-dimensional GIS}
\label{se:isitworthit}
As this thesis has shown, there are many potential advantages to the use of higher-dimensional modelling in GIS.\@
This approach provides \emph{a simple and consistent way to store geometry, attributes and topological relationships between objects of any dimension}.
This generic technique can be easily extended to handle other non-spatial characteristics, enabling better data management and more powerful operations.
At the same time, higher-dimensional representations are undoubtedly memory-intensive and often hard to work with, both due to their level of abstraction and the unintuitiveness of working with dimensions higher than three.
Admittedly, \emph{current GIS use cases do not easily justify the higher-dimensional approach}.
% Envisioning things
However, a far stronger case for higher-dimensional representations emerges when considering what sort of \emph{new tools could be developed using this type of representations}, both in GIS and in related fields.
For instance, 3D modelling software could consider time-varying topology, much as \citet{Dalstein15} do for 2D vector drawings.
4D topology (as 3D+time) could then be used to automatically generate smooth transitions for animations.
Similarly, CAD tools could use 4D topology to model buildings at different scales and timeframes automatically, keeping track of all relationships between objects and providing immediate user input during interactive editing.
For example, a program could display how different changes affect construction time, the size of the model at predefined LODs and when certain safety constraints were violated.
At a lower-level, 4D geometric modelling operators could be developed to operate directly on 4D primitives, such as splitting and merging 4-cells.
These could also be used intuitively in an interactive environment, allowing for instance to change a building's configuration by adding and removing walls while always ensuring that the representation remains a valid 4D space partition.
Considering that current GIS are very often used to manage large heterogeneous datasets and keep them up to date, a future system using a higher-dimensional underlying representation could be used to enforce certain validity constraints at the data structure level, such as avoiding 4D intersections or preserving a certain degree of continuity at LOD transitions.
% Requirements, perfect data
\subsection*{Is the higher-dimensional approach worthwhile?}
In short, yes, but only given certain conditions.
\emph{Higher-dimensional modelling is advantageous, but only when the added functionality that will be built with them---as compared to the simpler and more compact 2D/3D models---justifies it}, such as when queries across space and time can be implemented as higher-dimensional geometric/topological operations.
As discussed in \refse{ss:dimensions}, higher-dimensional modelling makes sense only when the characteristics depicted as dimensions are \emph{parametrisable and independent from each other}, as is the case for space/time/scale, and where objects occur along a dimension as intervals, not as discrete points (which can be easily stored as attributes).
Based on current hardware, software and the typical GIS datasets, there is another important practical requirement: the manageable number of dimensions is limited to 6--8.
Finally, it is also worth noting that higher-dimensional models are not incompatible with standard 2D/3D data structures and methods---the best tool for the job can be chosen depending on the need at hand and both approaches can be combined.
\section{Lessons learned}
\label{se:conclusions}
\begin{description}
% Chapter 2 + 3
\item[Current representations are not suitable in higher dimensions]
The mathematical foundations of spatial data modelling (\refch{ch:modelling-mathematics}) are defined in a dimension-independent manner, including all the basic tenets of geometry and topology.
This dimensional independence is also true for all spatial data models or representation schemes (\refch{ch:modelling-background})---at least when they are analysed at a high level---but is generally not preserved when they are implemented into more concrete data structures.
Most 2D and 3D data structures in GIS thus only encode a few chosen geometric and topological properties (\eg\ the coordinates of each point and the adjacencies between polygons), which are often defined with a formulation that is different per dimension.
This modelling approach can make for custom structures that are compact and efficient when used exactly as intended (\ie\ for a particular class of objects of a given dimension), but it can limit functionality or introduce inefficiencies when the data structures are adapted to be used under a different set of circumstances.
Case in point, the typical data structures of 2D GIS are frequently used with minimal changes for 3D, spatiotemporal and multi-scale GIS.\@
This results various problems, such as inefficient representations (\eg\ due to duplicate elements), an inability to represent common 3D objects (\eg\ those with a non-2-manifold boundary), difficulties in expressing 3D topological relationships (\eg\ adjacencies between solids), and the widespread availability of invalid datasets (\eg\ 3D models that do not formally enclose any space), among others.
% Chapter 4 + 5
\item[Higher-dimensional modelling as a solution]
An alternative to the use of ad hoc adaptations to 2D data structures is to model both spatial and non-spatial characteristics as dimensions in the geometric sense (\refch{ch:nd-modelling}).
While this approach can be memory-intensive, it provides a generic solution that can be applied to the representation of $n$D space, time, scale and any other parametrisable characteristics.
A tuple of $n$ parametrisable spatial and non-spatial characteristics can thus define a coordinate system in $\mathbb{R}^n$, and lower-dimensional 0D--3D objects existing across these characteristics can thus be modelled as higher-dimensional 0D--$n$D objects embedded in higher-dimensional space.
As GIS objects are usually non-overlapping\footnote{At least in theory.}, they should form an $n$D space partition and can thus be represented using $n$D topological data structures, reducing the total number of elements and ensuring that it is easy to navigate between the objects.
However, even when the objects do not form a space partition, the objects themselves can be partitioned by using an intermediate representation where the original objects are transformed into a set of non-overlapping regions, such that each of these regions represents \emph{a set of the original objects} \citep{Rossignac89}.
$n$D space partitions are also ideal in a practical sense, as they simplify many of the operations that can be defined with them, including simple point-in-polytope queries and all of the constructions methods presented in \refchs{ch:extrusion}--\ref{ch:linking-lods}.
$n$D point clouds, while outside the scope of this thesis, are a good complement to $n$D space partitions, as they are close to data as it is acquired and are relatively easy to store and manipulate.
\item[The most promising higher-dimensional representations]
Spatial data structures often consist of two aspects: 1.\ a combinatorial part, which consists of a set of primitives and some topological relationships between them, and 2.\ an embedding that links these primitives to their geometry and attributes \citep{Lienhardt94}.
As argued in this thesis, the knowledge of higher-dimensional topological relationships in a data structure is the main aspect that differentiates higher-dimensional data structures from 2D/3D ones.
Some data structures that are frequently used in 2D/3D GIS have straightforward extensions to higher dimensions, such as $n$D \emph{rasters} and \emph{hierarchies of trees}.
These use similar structures and algorithms as their 2D/3D counterparts and are therefore easy to understand and implement.
Other data structures implement the models of an $n$D simplicial complex or an $n$D cell complex.
As the simplices in a simplicial complex have a known number of adjacent simplices and bounding facets, they are most efficiently stored using simplex-based data structures.
Meanwhile, cell complexes can be easily stored using incidence graphs and related structures.
However, Nef polyhedra \citep{Bieri88} are probably the most promising representation for an $n$D cell complex, as they provide a good base to develop Boolean set operations, enabling a wide range of geometric operations.
Ordered topological models such as the cell-tuple \citep{Brisson93} and generalised/combinatorial maps \citep{Lienhardt94} also deserve a special mention.
By combining the strong algebra and easy navigation of a simplicial complex with the easy representation of a cell complex, they provide the most important benefits of both.
They are rather memory-intensive, but it is important to note that can still be more compact than a non-topological approach (\refch{ch:operations-background}).
Nevertheless, non-topological higher-dimensional representations do have a clear role to play as exchange formats, much as is the case for those based around Simple Features in 2D \citep{SimpleFeatures1}, and CityGML \citep{CityGML2} and IFC\footnote{\url{http://www.buildingsmart-tech.org/specifications/ifc-releases}} in 3D.
\item[Three construction methods for higher-dimensional objects]
Creating computer representations of higher-dimensional objects can be complex.
Common construction methods used in 2D and 3D, such as directly manipulating combinatorial primitives, or using primitive-level construction operations (\eg\ Euler operators \citep{Mantyla88}), rely on our intuition of 2D/3D geometry, and thus do not work well in higher dimensions.
It is therefore all too easy to create invalid objects, which then cannot be easily interpreted or fixed---a problem that is already exceedingly apparent in most 3D datasets.
As an alternative to the use of simple operations on combinatorial primitives, this thesis thus proposed three higher-level methods, all of which are relatively easy to use and attempt to create valid output.
% Chapter 6
\item[Method I.\ constructing objects using $n$D extrusion]
\hspace{15mm}
Extrusion as used in GIS has a natural extension into a dimension-independent formulation (\refch{ch:extrusion}).
Starting from an $(n-1)$-dimensional space partition as an $(n-1)$-dimensional cell complex and a set of intervals per cell, it is possible to extrude them to create an $n$-dimensional cell complex.
It is the easiest method to load existing 2D or 3D data into a higher-dimensional structure, representing a set of cells that exist along a given dimension, such as a length of time or a range of scales.
It is also easy to guarantee that the output cell complex is valid and can be used as a base for further operations, such as dimension-independent generalisation algorithms.
The extrusion algorithm developed in this thesis works on the basis of a generalised map representation of the cell complex and is relatively fast, with a worst case complexity of $O(ndr)$ in the main algorithm, where $n$ is the extrusion dimension, $d$ is the total number of darts in the input map and $r$ is the total number of intervals in the input, but offers better complexity in practice.
It is also memory-efficient, as only three layers of darts (of the size of the input cell complex) need to be kept in memory at the same time.
% Chapter 7
\item[Method II.\ constructing $n$D objects incrementally]
Based on the Jordan-Brouwer separation theorem \citep{Lebesgue11,Brouwer11}, it is known that an $i$-cell can be described based on a set of its bounding $(i-1)$-cells (\refch{ch:incremental-construction}).
Since individual $(i-1)$-cells are easier to describe than the $i$-cell, this can be used to subdivide a complex representation problem into a set of simpler, more intuitive ones.
This method can be incrementally applied to construct cell complexes of any dimension, starting from a set of vertices in $\mathbb{R}^n$ defined by a $n$-tuple of their coordinates, and continuing with cells of increasing dimension---creating edges from vertices, faces from vertices or edges, volumes from faces and so on.
The incremental construction algorithm developed in this thesis solves this problem in a practical setting by computing the topological relationships connecting the bounding $(i-1)$-cells.
It uses indices on the lexicographically smallest vertex of every cell per dimension, as well as an added index using the lexicographically smallest vertex of the ridges around the bounding facets of the cell that is being built.
It generates an $i$-cell in $O(d^{2})$ in the worst case, with $d$ the total number of darts in the cell.
However, it fares markedly better in real-world datasets, as cells do not generally share the same lexicographically smallest vertex.
By checking all matching ridges within a cell's facets, the algorithm can optionally verify that the cell being constructed forms a combinatorially valid quasi-manifold, avoiding the construction of invalid configurations.
% Chapter 8
\item[Method III.\ linking 3D models at different LODs into a 4D model]
As an example high-level higher-dimensional object construction method, a 4D model can be constructed from a series of different 3D models at different LODs (\refch{ch:linking-lods}).
The method presented in this thesis consists of three steps: identifying corresponding elements in different LODs, deciding how these should be connected according to a linking scheme, and finally linking relevant 3-cells into 4-cells.
Different linking schemes yield 4D models having different properties, such as objects that suddenly appear and disappear, gradually change in size or morph into different objects along the fourth dimension.
By modelling the LOD as a dimension, the correspondences between equivalent objects across LODs become geometric primitives, making it possible to perform geometric operations with them (\eg\ extracting an intermediate LOD for visualisation purposes) or to attach attributes to them (\eg\ general semantics or the meaning of these correspondences), just as is done to other geometric primitives.
These topological relationships and correspondences can then be used for multiple applications, such as updating and maintaining series of 3D models at different LODs, or testing the consistency of multi-LOD models (\eg\ by using the validity checks in \citet{Groger11}).
% Chapter 9
\item[Extracting 2D/3D subsets from an $n$D model]
The process to obtain a lower-dimensional subset of a higher-dimensional dataset can be regarded as a function that maps a subset of $\mathbb{R}^n$ to a subset of $\mathbb{R}^m$, $m < n$, which is obtained by cutting through the dataset in a geometrically meaningful way (\refch{ch:slicing}).
Broadly, this process consists of two steps: (i) selecting a subset of the objects in the model and (ii) projecting this subset to a lower dimension.
Both of these steps can vary substantially.
Selecting a subset of the objects can be as simple as obtaining those within a axis-aligned bounding box, or can be as complex as a Boolean set intersection operation, such as for the computation of cross-sections.
Meanwhile, there are a wide variety of transformations that apply different projections with different properties, such as the $n$-dimensional to ($n-1$)-dimensional orthographic and perspective projections derived in this thesis and the $\mathbb{R}^n$ to $S^{n-1}$ spherical projection used in its cover.
% Chapter 10 + Appendix A
\item[Methods to create valid objects and space partitions in 2D and 3D]
Most algorithms described in computational geometry and GIS assume that their input datasets are flawless and they are processable using real numbers.
However, invalid datasets are widespread in GIS (\refch{ch:cleaning}), and they are represented and processed using limited-precision arithmetic (\refap{ch:implementation}).
In order to cope with 2D invalid datasets, this thesis further developed methods to create valid polygons and planar partitions using a constrained triangulation of the input.
These were based on the work done in \citet{ArroyoOhori10}, improving the reconstruction algorithm, fixing edge cases, implementing an odd-even constraint counting mechanism and improving the quality of the implementation.
Similarly, a method to repair 3D objects and space subdivisions was developed by snapping together lower-dimensional primitives and removing overlaps using Boolean set operations on Nef polyhedra \citep{Bieri88,Hachenberger06}.
These methods were used in this thesis in order to use real-world datasets in practice, such as when applying the construction algorithms.
% There is no systemic reason to tolerate invalid data
\end{description}
\section{Contributions}
\label{se:contributions}
The main contribution of this thesis is the realisation of the fundamental aspects of a higher-dimensional Geographic Information System.
By approaching this problem in a practical manner, many of the technical issues of its development were investigated, including an analysis of its possible internal (in-memory) and external (exchange format) representations, the development of basic algorithms for object construction and visualisation, and the development of GIS data repair tools for 2D and 3D datasets.
By taking a model that was previously only described at a conceptual level \citep{vanOosterom10} and realising it, it is now possible to more fully evaluate the consequences of this higher-dimensional approach.
In more concrete terms, there were several smaller contributions that were necessary to be able to achieve this realisation.
The most significant ones are:
\begin{description}
\item[Survey and analysis of higher-dimensional models and structures] I conducted a survey of all of the main data models and data structures used in GIS, geometric modelling and related fields, considering 2D, 3D and $n$D data structures.
These were analysed in terms of their feasibility for the higher-dimensional modelling of geographic information, including how they could handle different geometry classes, topology and attributes, either in their current form or through modifications, as well as their ease of implementation in practice.
\item[Three construction methods] I developed three easy-to-use construction methods for objects of any dimension.
The two methods lower-level methods---extrusion and incremental construction---were implemented using CGAL Combinatorial Maps and tested using real-world datasets.
The third method has been tested with a few synthetic datasets, but more work is necessary to fully realise it and automate it.
All of the implementations were made available publicly under an open source licence.
\item[Higher-dimensional real-world models] As part of this thesis, I created higher-dimensional models from real-world 2D and 3D datasets.
To the best of my knowledge, these are the only datasets consisting of realistic higher-dimensional objects in a GIS setting.
\item[Combinatorial map reversal] As part of the development of the incremental construction operation, I developed additional functions that were added to CGAL Combinatorial Maps after being approved by the CGAL Editorial Board.
These functions involved reversing the orientation of a combinatorial map of any dimension.
\item[Simple formulation of $n$D to ($n-1$)D projections] I developed intuitive formulations of $n$-dimensional to ($n-1$)-dimensional orthographic and perspective projections.
While other formulations exist, my formulations based on normal vectors are in my opinion the easiest to understand and manipulate, at least for a GIS audience.
\item[Repair methods tools] I developed methods to automatically repair 2D polygons and planar partitions, as well as 3D polyhedra and space partitions.
The 2D methods were also released publicly under an open source licence as \texttt{prepair} and \texttt{pprepair}.
The 3D methods will also be released publicly after further improvements are made, together with more thorough testing and the addition of basic documentation.
\end{description}
\section{Future work}
\label{se:futurework}
As this thesis challenges many of the assumptions underpinning current GIS, there are many potential lines of research that can be formulated for higher-dimensional GIS and higher-dimensional modelling in general.
Most algorithms for 2D/3D GIS have a dimension-dependent formulation and would result in open problems in an $n$D context.
However, while these are worthy of attention, I would like to focus on what are still significant gaps in knowledge for the implementation of a higher-dimensional GIS and that could not be solved in this thesis' timeframe.
While some of these research topics are not within the main subject matter of GIS research, they are what I consider to be the key steps for a more complete implementation of a working system.
They are:
\begin{description}
\item[Low-level linking algorithms] The linking schemes from \refch{ch:linking-lods} have only been described in terms of high-level algorithms.
These should be further developed by finding adequate low-level algorithms to identify matching elements using customisable constraints.
For instance, it is reasonable to attempt to minimise a certain distance function between two models (\eg\ Earth mover's distance), but it is also important to do so in a manner that preserves the topological relationships between the objects.
The special treatment of holes of different dimensions should also be investigated, together with adequate methods to ensure that holes are linked correctly between themselves and to other primitives.
\item[High-level construction algorithms] The three object construction methods described in this thesis operate mostly on lower-dimensional primitives and consequently cover only a few use cases.
There is a need to develop intuitive methods that operate directly on higher-dimensional primitives, which should preferably be usable in an interactive environment.
Note however that this does not mean that the end user would be viewing the higher-dimensional model directly, as it could be shown in simplified form or as a 2D or 3D representative subset.
\item[$n$D constrained triangulator] There are high-quality robust implementations of constrained Delaunay triangulations in 2D \citep{Shewchuk96} and 3D \citep{Si05}, as well as good descriptions of a constrained Delaunay triangulation in $n$D \citep{Shewchuk07}.
However, in order to realise a higher-dimensional GIS based on a simplicial complex model, it is necessary to have a robust $n$D constrained triangulator which should be preferably Delaunay.
\item[Hyperspherical projective geometry kernel] I deemed Nef polyhedra one of the most promising models for a higher-dimensional GIS, but so far they have only been implemented in 2D and 3D.\@
In order to implement $n$-dimensional Nef polyhedra, it is necessary to develop a hyperspherical projective geometry kernel.
This kernel would then be used to compute the local $(n-1)$-dimensional pyramids around every vertex.
While the projective mathematics for this are relatively simple, it is a complex engineering problem to implement it robustly.
\item[$n$D Boolean set operations] Using the above mentioned hyperspherical projective kernel or another method, it is necessary to develop robust algorithms to compute $n$D Boolean set operations.
For instance, an $n$D Nef polyhedra implementation using recursive boundary definitions could compute these operations at the local pyramid level.
$n$D Boolean set operations would be an excellent base for most geometric operations in a higher-dimensional GIS.\@
\item[$n$D to/from 2D and 3D projective kernel] Many operations that are required for object manipulation in $n$D are actually well-defined operations applied to 2D and 3D objects that are merely embedded in $n$D.
A robust kernel that could handle on-the-fly conversions of $n$D geometries to 2D and 3D without loss of precision would enable many of these operations to be applied.
Relatedly, CGAL currently has simple kernels that apply 3D to 2D orthographic projections using the planes defined by the $xy$, $yz$ or $zx$ planes.
These are useful as they can be wrapped around the basic 2D kernels (\ie\ $\mathbb{R}^2$ with floating-point, interval or exact representations) using the traits programming paradigm available throughout CGAL, and so they can be used to apply 2D operations to 2D objects that are embedded in 3D.
However, these kernels do not handle the conversions from 2D back to 3D automatically nor are easily extensible to higher dimensions.
\item[Visualisation of higher-dimensional geographic information] \refch{ch:slicing} described how data selection and projection methods can work in arbitrary dimensions, while \refse{se:ndmath} described the $n$D mathematics behind the basic manipulation operations for higher-dimensional objects.
However, there are significant issues to tackle in order to create a useful visualiser for higher-dimensional datasets.
Namely, there are significant hurdles in user interaction, dealing with large datasets, computing higher-dimensional cross-sections and the definition of useful visual cues in higher dimensions.
The simple implementation used for the cover of this thesis make several hard-coded assumptions for the dataset that was used and is far from optimal, but it will be published together with the rest of the open source tools developed for this thesis once it is improved to handle more general input in the form of $n$D linear cell complexes.
\item[2D/3D/$n$D repair methods with quality guarantees] The data repair methods described in \refch{ch:cleaning} are able to recover from most simple invalid 2D and 3D configurations.
However, they are not easily extensible to higher dimensions due to the lack of a robust $n$D constrained triangulator (see above), they do not provide quality guarantees in their output, and can be numerically unstable when there are many nearly coplanar planes.
In order to develop robust systems, it is highly desirable to be able to specify a robustness criterion (\eg\ a minimum distance between vertices or ensuring that the geometries are not collapsed/flipped in a floating-point representation), which is guaranteed by a data repair algorithm.
Additional geometric constraints could also be implemented, such as guaranteeing that coplanar planes stay coplanar after repair.
\item[Higher-dimensional modification operations] By necessity, this thesis focused on operations for object creation in order to generate initial higher-dimensional datasets.
Now that it is possible to generate them, it is important to think of intuitive higher-dimensional object modification operations.
For instance, how can a 4-cell split operation be intuitively defined, or how can collapsed geometries (\eg\ from the extruded and collapsed models in \refse{se:extrusion-generalisation}) be processed to remove combinatorial elements.
Ideally, these operations should also be intuitively usable in an interactive environment, such as a geometric modeller.
\item[Real-world 4D spatiotemporal datasets] Using timestamped 3D volumetric datasets, it should be possible to create true 4D datasets using spatiotemporal information.
However, obtaining reasonably clean volumetric datasets is nearly impossible at this point.
Every dataset that was found during this project had only surfaces embedded in 3D, had severe validity problems up to the point that it would require substantial manual work to fix, or was missing the temporal information.
Exporting the temporal information that is present in some closed commercial formats is also an issue, as doing do in a na{\"\i}ve way generally means losing the links to the timestamps' corresponding geometries.
\end{description}

BIN
work_report/cover-back.pdf Normal file
View File

Binary file not shown.

BIN
work_report/cover-front.pdf Normal file
View File

Binary file not shown.

51
work_report/cv.tex Normal file
View File

@ -0,0 +1,51 @@
% !TEX root = thesis.tex
\chapter{Curriculum vitae}
\noindent
Gustavo Adolfo \href{http://ken.mx}{\textbf{Ken Arroyo Ohori}} was born on June 12, 1985 in Mexico City.
He graduated in 2003 from the Bicultural High School (`\emph{Preparatoria Bicultural}') programme and in 2007 from the BSc in Computer Science and Technology programme (`\emph{ITC01: Ingenier\'\i{}a en Tecnolog\'\i{}as Computacionales}') at the Mexico City campus of the Monterrey Institute of Technology and Higher Education (`\emph{Instituto Tecnol\'ogico y de Estudios Superiores de Monterrey}').
\marginpar{\includegraphics[width=\marginparwidth]{figs/me}}
Arriving to the Delft University of Technology in 2008, Ken graduated from the MSc in Geomatics programme in 2010 with his thesis on the \emph{Validation and automatic repair of planar partitions using a constrained triangulation}, which later evolved into some of the data repair ideas exposed in \refch{ch:cleaning} of this thesis.
In 2011, Ken started his PhD project on the \emph{Higher-dimensional modelling of geographic information}, supervised by \href{https://3d.bk.tudelft.nl/jstoter/}{Jantien Stoter} and \href{http://tudelft.nl/hledoux}{Hugo Ledoux}, working under the umbrella of the project \emph{5D Data Modelling: Full Integration of 2D/3D Space, Time and Scale Dimensions} of the Dutch Technology Foundation (\emph{STW}).
In October 2012, he visited \href{http://liris.cnrs.fr/guillaume.damiand/}{Guillaume Damiand}, collaborating on the incremental construction method that later became \refch{ch:incremental-construction} of this thesis.
In 2016, Ken started working as a postdoc in the same \emph{5D Data Modelling} project, focusing on 4D visualisation and in the integrated modelling of space and time.
\clearpage
\section*{Publications}
{\small
\begin{itemize}
\papermethodsxvoxelisation%
\paperijgisroeland%
\paperudmvobj%
\paperisprsnd%
\paperijgind\mparshift{-3.3\baselineskip}\marginpar{\raggedleft{}In \refch{ch:linking-lods}}%
\paperijgisextrusion\mparshift{-3.3\baselineskip}\marginpar{\raggedleft{}In \refch{ch:extrusion}}%
\paperijgisndstructures\mparshift{-3.3\baselineskip}\marginpar{\raggedleft{}In \refse{se:data-structures} \& \refse{se:nd-modelling-conclusions}}%
\papercgeoprepair\mparshift{-2.3\baselineskip}\marginpar{\raggedleft{}In \refse{se:pprepair}}%
\papericaaincrementalconstruction\mparshift{-6.3\baselineskip}\marginpar{\raggedleft{}In \refch{ch:incremental-construction}}%
\paperacmsigspatialextrusion%
\clearpage%
\papericcsand%
\papergeoadvancesnd\mparshift{-5.2\baselineskip}\marginpar{In \refse{se:duality}}%
\paperagileslicing%
\paperpfgpprepair\mparshift{-3.2\baselineskip}\marginpar{In \refse{se:pprepair}}%
\paperosgisrepair%
\papertdgeoinfond%
\paperagileprepair%
\paperostravaedgematching%
\end{itemize}
}
\clearpage
\null%
\marginpar{
\includegraphics[width=\marginparwidth]{figs/3dgeoinfo}\\
\vspace*{11.5cm}
\includegraphics[width=\marginparwidth]{figs/barcode}
}

278
work_report/dictionary.tex Normal file
View File

@ -0,0 +1,278 @@
% !TEX root = thesis.tex
\chapter{A short dictionary of dimension-based GIS terms}
\label{ch:dictionary}
\begin{multicols}{2}
\begin{description}
\footnotesize
\item[ambient space]
the space in which objects are embedded.
\item[area]
1.\ measure of the 2D extent of an object;
2.\ 2D combinatorial element;
{\color{gray} 3.\ 2D geometric object}.
\item[ball]
1.\ ($\rightarrow$\ same as 3-\emph{ball}) topological definition of the space within a 2-sphere (\ie\ a sphere), often defined as the space within a certain distance (the radius) from a point in $\mathbb{R}^3$ (the centre);
2.\ ($\rightarrow$\ often $n$-\emph{ball}) topological definition of the space within an $(n-1)$-sphere, often defined as the space within a certain distance (the radius) from a point in $\mathbb{R}^n$ (the centre);
3.\ geometric definition of a perfectly round filled 3D object.
\item[body]
$n$-dimensional combinatorial element in an $n$D context.
\item[box]
1.\ cuboid;
2.\ orthotope of any dimension.
\item[cavity]
3D hole.
\item[cell]
1.\ ($\rightarrow$\ sometimes $n$-\emph{cell}) $n$D combinatorial element;
{\color{gray} 2.\ 3D combinatorial element}.
\item[cell complex]
topological space formed by a set of cells glued along common faces, all faces of a cell in the complex should be in the complex.
\item[circle]
1.\ topological definition of the space at a certain distance (the radius) from a point in $\mathbb{R}^2$ (the centre);
2.\ geometric definition of a perfectly round hollow 2D object.
\item[closed]
1.\ topological definition of an object that includes its boundary;
2.\ geometric definition of an object in $n$D, usually defined using boundary representation, that encloses an $n$D subspace.
\item[congruent]
having the same shape and size.
\item[cube]
1.\ ($\rightarrow$\ same as $3$-\emph{cube}) 3D polyhedron with 6 square facets;
2.\ ($\rightarrow$\ usually $n$-\emph{cube}) $n$-orthotope with identical $(n-1)$-cube facets, akin to a square in 2D and a cube in 3D.
\item[cuboid]
1.\ box-shaped 3D polyhedron with 6 rectangular facets, its opposite facets are congruent and parallel;
{\color{gray} 2.\ parallelotope}.
\item[curve]
1.\ 1D geometric object, often with a non-linear geometry;
2.\ 1D combinatorial element;
3.\ 1-manifold.
\item[cut-line]
degenerate part of a 2D shape forming a curve and protruding inwards from its boundary.
\item[disk]
topological definition of the space within a circle, often defined as the space within a certain distance (the radius) from a point in $\mathbb{R}^2$ (the centre).
\item[edge]
1D combinatorial element.
\item[face]
1.\ 2D combinatorial element;
2.\ ($\rightarrow$\ sometimes \emph{$n$-face of $X$}) $n$D combinatorial element on the boundary of $X$;
{\color{gray} 3.\ ($\rightarrow$\ usually \emph{face of $X$}) $(n-1)$D combinatorial element on the boundary of an $n$D combinatorial element $X$}.
\item[facet]
1.\ $(n-1)$D combinatorial element in an $n$D context;
2.\ ($\rightarrow$\ often \emph{facet of $X$}) $(n-1)$D combinatorial element on the boundary of an $n$D combinatorial element $X$;
{\color{gray} 3.\ 2D combinatorial element}.
\item[flat]
unbounded geometric object with linear geometry of any dimension, such as a point, line or plane.
\item[hole]
1.\ 2D void region;
2.\ $n$D void region.
\item[hyperball]
higher-dimensional ball.
\item[hypercell]
1.\ higher-dimensional combinatorial element;
2.\ 4D combinatorial element.
\item[hypercube]
1.\ higher-dimensional cube;
2.\ 4-cube.
\item[hyperplane]
1.\ $(n-1)$D linear subspace in an $n$D context;
{\color{gray} 2.\ plane in a higher-dimensional context.}
\item[hyperrectangle]
higher-dimensional orthotope.
\item[hypersphere]
higher-dimensional sphere.
\item[hypersurface]
$(n-1)$D subspace in an $n$D context, often curved;
\item[interval]
topological definition of the space between two points in $\mathbb{R}$ (the endpoints).
\item[Lebesgue measure]
measure of the $n$D extent of an object, akin to length in 1D, area in 2D and volume in 3D.
\item[length]
measure of the 1D extent of an object.
\item[line segment]
1D geometric object.
\item[manifold]
($\rightarrow$\ sometimes $n$-\emph{manifold}) topological space that resembles $\mathbb{R}^n$ at every point.
\item[manifold with boundary]
($\rightarrow$\ sometimes $n$-\emph{manifold with boundary}) topological space that resembles $\mathbb{R}^n$ at every point in its interior.
\item[node]
0D combinatorial element.
\item[open]
1.\ topological definition of an object that does not include its boundary;
2.\ geometric definition of an object in $n$D, usually defined using boundary representation, that does not enclose any $n$D subspace.
\item[orthotope]
($\rightarrow$\ sometimes $n$-\emph{orthotope}) $n$-polytope with congruent parallel facets forming right angles to each other, akin to a rectangle in 2D and a cuboid in 3D.
\item[parallelepiped]
3D polyhedron bounded by 6 parallelogram-shaped 2D faces.
\item[parallelogram]
polygon bounded by 4 parallel and have the same shape.
\item[parallelotope]
a polytope whose opposite facets are parallel and have the same shape, akin to a parallelogram in 2D and a parallelepiped in 3D.
\item[peak]
1.\ $(n-3)$D combinatorial element in an $n$D context;
2.\ ($\rightarrow$\ \emph{peak of $X$}) $(n-3)$D combinatorial element on the boundary of an $n$D combinatorial element $X$;
\item[plane]
2D unbounded linear subspace.
\item[point]
1.\ 0D geometric element;
2.\ topological definition of the space at one location.
\item[polychoron]
4D geometric object with linear geometry.
\item[polygon]
2D geometric object with linear geometry, possibly with holes.
\item[polygonal curve]
curve formed by a sequence of line segments joined at their endpoints.
\item[polyhedron]
1.\ 3D geometric object with linear geometry, possibly with holes;
{\color{gray} 2.\ an $n$D geometric object with linear geometry}.
\item[polyline]
curve formed by a sequence of line segments joined at their endpoints.
\item[polyteron]
5D geometric object with linear geometry.
\item[polytope]
$n$D geometric object with linear geometry.
\item[puncture]
0D (point) hole, often formed from a degenerate polyline or polygon.
\item[ridge]
1.\ $(n-2)$D combinatorial element in an $n$D context;
2.\ ($\rightarrow$\ \emph{ridge of $X$}) $(n-2)$D combinatorial element on the boundary of an $n$D combinatorial element $X$.
\item[line string]
curve formed by joining a sequence of vertices by straight line segments, represented by these vertices.
\item[linear ring]
2D geometric object with linear geometry represented as a sequence of vertices.
\item[ring]
2D combinatorial element represented by its 1D boundary.
\item[subfacet]
1.\ $(n-2)$D combinatorial element in an $n$D context;
2.\ ($\rightarrow$\ \emph{subfacet of $X$}) $(n-2)$D combinatorial element on the boundary of an $n$D combinatorial element $X$.
\item[shell]
3D geometric object with linear geometry without 3D holes, usually defined as the volume enclosed by its 2D boundary, often represented as a set of 2D faces.
\item[simplex]
1.\ ($\rightarrow$\ sometimes \emph{$n$-simplex}) combinatorial element with $n+1$ vertices and $n+1$ facets;
2.\ ($\rightarrow$\ sometimes \emph{$n$-simplex}) geometric object based on $n+1$ affinely independent vertices, akin to a point in 0D, a line segment in 1D, a triangle in 2D, or a tetrahedron in 3D.
\item[simplicial complex]
1.\ ($\rightarrow$\ sometimes \emph{abstract simplicial complex}) topological space formed by a set of simplices glued along common faces, all faces of a simplex in the complex should be in the complex;
2.\ ($\rightarrow$\ sometimes \emph{geometric simplicial complex}) abstract simplicial complex where simplices are embedded into Euclidean space, their interiors should not intersect geometrically.
\item[solid]
1.\ 3D geometric object with linear geometry, possibly with 3D holes, often represented as a set of shells;
2.\ object defined based on solid modelling, often as opposed one based on boundary representation;
3.\ an object containing its interior (as opposed to \emph{hollow}).
\item[spike]
degenerate part of a 2D shape forming a curve, often protruding outwards from its boundary.
\item[sphere]
1.\ ($\rightarrow$\ same as $2$-\emph{sphere}) topological definition of the space at a certain distance (the radius) from a point in $\mathbb{R}^3$ (the centre);
2.\ ($\rightarrow$\ often $n$-\emph{sphere}) topological definition of the space at a certain distance (the radius) from a point in $\mathbb{R}^{n+1}$ (the centre);
3.\ geometric definition of a perfectly round hollow 3D object.
\item[surface]
1.\ 2-manifold;
2.\ 2D combinatorial element.
\item[tesseract]
4-cube: polychoron with 8 cubical facets.
\item[vertex]
0D combinatorial element.
\item[volume]
1.\ measure of the 3D extent of an object;
2.\ 3D combinatorial element.
\item[wire]
2D combinatorial element represented by its 1D boundary.
\end{description}
\end{multicols}
% \begin{table}
% \scriptsize
% \begin{tabular}{ccccccc}
% \toprule
% \multicolumn{7}{c}{dimension} \\
% $n$ & 0 & 1 & 2 & 3 & 4 \\
% \midrule
% \multicolumn{7}{c}{topological} \\
% $n$-simplex & vertex & edge & face & volume \\
% $n$-cell & vertex & edge & face & volume \\
% $n$-sphere & pair of points & circle & sphere \\
% $n$-ball & & interval & disk & ball \\
% open $n$-ball & & open interval & open disk & open ball \\
% $n$-manifold & discrete space & curve & surface \\
% \midrule
% \multicolumn{7}{c}{geometric} \\
% $n$-simplex & point & line segment & triangle & tetrahedron & pentachoron \\
% $n$-cube & point & square & cube & tesseract \\
% $n$-orthotope & point & line segment & rectangle & cuboid/box \\
% $n$-parallelotope & point & line segment & parallelogram & parallelepiped \\
% prismatic $n$-polytope & & line segment & rectangle & prism & prismatic polychoron \\
% $n$-polytope & point & line segment & polygon & polyhedron & polychoron \\
% \midrule
% \multicolumn{7}{c}{other} \\
% $n$-th Lebesgue measure & length & area & volume \\
% & & collinear & coplanar \\
% \bottomrule
% \end{tabular}
% \end{table}

6541
work_report/docs.bib Normal file
View File

File diff suppressed because it is too large Load Diff

693
work_report/extrusion.tex Normal file
View File

@ -0,0 +1,693 @@
% !TEX root = thesis.tex
\chapter{Extrusion}
\label{ch:extrusion}
Extrusion is a commonly used technique in GIS to construct simple 3D models.
Starting from a planar partition of polygons and a height interval associated to each of them, it generates a set of space-partitioning box-shaped polyhedra by considering that each polygon exists all along its related interval.
For instance, a set of building footprints and associated heights is extruded into a set of simple prismatic buildings.
Based on the fundamental operations on generalised and combinatorial maps discussed in \refch{ch:operations-background}, this chapter presents a generalisation of this technique to higher dimensions: an $(n-1)$-dimensional cell complex and a set of associated intervals per $(n-1)$-cell is thus transformed into an $n$-dimensional cell complex.
The chapter starts by presenting some background on 2D-to-3D extrusion and an intuitive description of dimension-independent extrusion in \refse{se:extrusion-background}.
Afterwards, \refse{se:extrusion-algos} describes an dimension-independent extrusion algorithm.
% that can create $n$-dimensional cell complexes from $(n-1)$-dimensional space partitions.
% by assigning one or more intervals to each $(n-1)$-cell
\refse{se:extrusion-implementation} describes how this algorithm was implemented based on the implementation of combinatorial maps in CGAL.\@
\refse{se:extrusion-experiments} summarises experiments using this implementation, which created consistent objects in up to 6D by combining publicly available GIS datasets.
\refse{se:extrusion-generalisation} describes how extrusion can be used not only directly, but as a base to create dimension-independent generalisation operations.
Finally, \refse{se:extrusion-conclusions} concludes the chapter with the main findings and the possibilities to use extrusion for higher-dimensional datasets.
Most of this chapter is based on the paper:
\begin{itemize}
\paperijgisextrusion%
\end{itemize}
\section{Background}
\label{se:extrusion-background}
Extrusion is most commonly used in GIS\footnote{This contrasts with the use of the term `extrusion' in other fields.
For instance, in geometric modelling, extrusion is a well-known operation in which all of the objects in the model are `dragged' along a predefined path given by a curve, more akin to an actual (physical) extrusion process, but whose output topology can be computed without any geometric computations.} in order to create simple 3D city models with box-shaped buildings.
In such a process, a user takes a set of building footprints, optionally subdivided into parts of similar height, and extrudes them using intervals that extend from the ground level to the building height at the location of each polygon.
This process thus creates simple 3D representations of each building, where a building is represented by one or more box-shaped polyhedra.
\reffig{fig:3dtop10nl} shows an example of such a model over a triangulated terrain (\ie\ a 2D simplicial complex embedded in 3D).
\begin{figure*}
\centering
\includegraphics[width=\linewidth]{figs/3dtop10nl}
\caption[A view of the 3D TOP10NL dataset]{A view of the 3D TOP10NL dataset\protect\footnotemark\ in the area of Valkenburg, the Netherlands.}
\label{fig:3dtop10nl}
\end{figure*}
\footnotetext{\url{https://www.pdok.nl/nl/producten/pdok-downloads/basis-registratie-topografie/top10nl-3d}}
Such an extrusion process is appealing largely because of its simplicity.
All the information that is required for the previous example is notoriously easy to obtain: building footprints are widely available---including in open datasets such as OpenStreetMap\footnote{\url{http://www.openstreetmap.org/}}---or can be easily obtained from satellite imagery, while the height information can be acquired with techniques like airborne laser scanning or photogrammetry.
Not coincidentally, this type of 3D models are commonly generated in practice and often appear in standards such as CityGML \citep{CityGML2}, where it is defined as the level of detail (LOD) 1.
Moreover, it is easy to ensure that extruded models are generated in a topologically consistent manner~\citep{Ledoux11a}.
If the input forms a planar partition---or is processed to form a planar partition, \eg\ using the method described in \refse{se:pprepair}---and the extrusion direction is orthogonal to the planar partition, \emph{the output polyhedra are guaranteed not to overlap}.
While not every possible 3D shape can be generated using extrusion---the top and bottom faces will always be horizontal, and the side faces connecting them will always be vertical---, a remarkable number of shapes can be constructed using this type of method \citep{Ferrucci93}.
This is especially true when simple extensions to the method are considered, such as multiple intervals per 2-cell, possibly overlapping 2-cells in the input, or an arbitrary extrusion direction.
Viewed at a more fundamental level, shown in \reffig{fig:extrusion}, extrusion `lifts' a 2D cell complex to form a 3D cell complex.
In order to form a valid cell complex as defined in \refse{ss:algebraic-topology}, cells must not overlap, but instead they should form a structure of cells of dimension from 0 to 3, where $i$-cells ($\forall i > 0$) are bounded by sets of $(i-1)$-cells.
The extruded cells must therefore take into account the incidence and adjacency relationships between cells, splitting cells of every dimension at the places where the common boundaries of higher-dimensional cells start and end.
\begin{figure}[tbp]
\centering
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/extrusion-1}
\label{subfig:extrusion-1}}
% \quad
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/extrusion-2}
\label{subfig:extrusion-2}}
\caption[Extrusion from 2D to 3D]{(a) A 2D cell complex, where each 2-cell is associated with a height, is extruded to generate (b) a 3D cell complex. Note how the extrusion of vertex $v$ in (a) causes the generation of several vertices and edges.
In particular, the vertex $v^\prime$ in (b), which is at the height of the top face of the back left polyhedron.
This vertex needs to be used in the representation of the all 4 edges, 5 faces and 3 volumes that are incident to it, not only in those that are part of the back left polyhedron. This ensures that the cells in (b) are pairwise disjoint and have the expected correct topology.}
\label{fig:extrusion}
\end{figure}
Extrusion, in the sense of its GIS definition as exposed in the previous paragraph, has a natural extension to higher dimensions.
A dimension-independent extrusion operation thus `lifts' an $(n-1)$-dimensional cell complex, to an $n$-dimensional cell complex by assigning one or more intervals to each $(n-1)$-cell, \ie\ a subset of 1D Euclidean space\footnote{This is equivalent to one or more edges in \citet{Lienhardt04} or 1-dimensional polyhedra in \citet{Ferrucci93}.}, along which this cell is defined.
For instance, just as a polygon representing the footprint of a building is commonly extruded using an interval $(0, height)$, a polyhedron representing a building can be extruded along the time interval $(construction, destruction)$.
As \reffig{fig:Frauenkirche} shows, extrusion using multiple unconnected intervals is also a possible use case.
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/Frauenkirche}
\caption[The Frauenkirche in Dresden]{The Frauenkirche in Dresden, Germany, was originally built in 1743, destroyed during the bombing of Dresden in 1945, and reconstructed in 2005.
A simple 4D model of the building could be made by extruding a 3D model of it along the interval $(1743, 1945) \cup (2005, \infty)$.
Photograph by David M\"uller in Wikimedia Commons.
}
\label{fig:Frauenkirche}
}
In the related literature, there are two methods that can be used for a similar purpose.
The Cartesian product, a more general operation defined for generalised maps in \citet{Lienhardt04}, generates a combinatorial structure equivalent to the extrusion of all the objects along a single interval by computing the Cartesian product of the original (unextruded) cell complex with an edge.
However, since such an operation is limited to a single interval (which is the same for all the input objects), the output cannot be directly applied for modelling real-world datasets.
A related possibility, presented by \citet{Ferrucci93}, is to first pre-process all the intervals (for all cells), splitting them into fragments at the other intervals' endpoints so that each fragment intersects another one if and only if they have the same endpoints, and then computing the output purely combinatorially for each possible combination of fragment and for each input cell.
While this allows us to model real-world objects, the process of splitting the input intervals might greatly increase the number of output cells since two intersecting intervals are split even if their corresponding cells are far away from each other.
The method presented in this chapter follows a similar approach, but instead process the input intervals for each input cell separately based on the incidence relationships with other cells, thus generating fewer total intervals, a smaller number of cells and a considerably smaller combinatorial structure.
\section{A dimension-independent extrusion algorithm}
\label{se:extrusion-algos}
The dimension-independent extrusion algorithm requires two input arguments: an $(n-1)$-dimensional space partition of $(n-1)$-polytopes embedded into $(n-1)$-dimensional space, stored as an $(n-1)$-dimensional generalised map $G$; and a map of extrusion intervals $\rho$\footnote{I originally referred to these as \emph{ranges} and thus used variations of $R$ and $\rho$ (rho).
While a kind anonymous reviewer correctly suggested that it is more precise to refer to them as \emph{intervals}, $\rho$ is more distinctive than iota ($\iota$) and upsilon ($\upsilon$), and so it remains as $\rho$ here.} that links each $(n-1)$-cell $c$ in $G$ to a set $R$ of 1-dimensional intervals, where every interval $r$ in $R$ is represented as a pair of values $(r_{\min}, r_{\max})$ where $c$ is extruded along the $n$-th dimension.
The intervals in $\rho$ for the cells of lower dimension do not need to be given, since they can be computed by the algorithm based on their incidence relationships to the $(n-1)$-cells.
Note that, as \reffig{fig:extrusion-intervals} shows, multiple intervals for the same cell are possible.
\begin{figure}[tbp]
\centering
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/intervals-1}
\label{subfig:extrusion-intervals-1}}
% \quad
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/intervals-2}
\label{subfig:extrusion-intervals-2}}
\caption[The extrusion intervals of the lower-dimensional cells]{All the cells are given extrusion intervals based on the intervals for the $(n-1)$-cells.
Note how it is possible to have multiple intervals per cell both in the input, \eg\ the purple square on the right is extruded along two intervals, and in the output, \eg\ the two edges $v^{(0,1)}$, $v^{(1.5,2)}$ and $v^{(2,3)}$ are the result of the extrusion of a single vertex $v$ whose extrusion interval was not directly given.}
\label{fig:extrusion-intervals}
\end{figure}
There are two cases in which a cell has multiple intervals.
One of these is that for an $(n-1)$-dimensional cell, multiple intervals may be explicitly provided in the input (\eg\ the example previously given in \reffig{fig:Frauenkirche}).
The other case is that for a lower-dimensional cell, multiple extrusion intervals may be passed to it by several adjacent higher-dimensional cells.
The result of the extrusion algorithm is an $n$-dimensional generalised map $G^{\prime}$ representing a $n$-dimensional cell complex containing a set of prismatic $n$-polytopes, the $n$-dimensional analogue of a set of prisms, which also form a $n$-dimensional space partition.
It is important to note that the output map $G^\prime$ creates entirely new structures, \ie\ it does not reuse the darts or embeddings of $G$, since the $(n-1)$-simplices in $G$ are similar but not identical to the $n$-simplices in the base of $G^\prime$.
These differ in terms of the highest-dimensional involution $\alpha_{n-1}$ and---importantly for an implementation---in the total number of involutions per dart.
As \reffigs{fig:trianglecomplex} and~\ref{fig:blockscomplex} show, the cells in the new $n$-dimensional cell complex in $G^\prime$ have a direct relation to and can be expressed in terms of the $(n-1)$-cells in $G$ and their extrusion intervals in $\rho$.
This property is used in order to define the cells of the output cell complex, which are equivalent to the embeddings in $G^\prime$.
These consist of:
\begin{itemize}
\item
`Base' and `top' $(n-1)$-cells (\ie\ facets), which are constructed from every $(n-1)$-cell $c$ in $G$ at the $r_{\min}$ and $r_{\max}$ end points of every interval in $\rho(c)$.
\item
A series of prismatic facets linking corresponding $(n-2)$-cells (\ie\ ridges) $c$ of the above mentioned $r_{\min}$ and $r_{\max}$ facets for every interval, such that every one of these facets corresponds to the extrusion of the ridge along an interval in $\rho(c)$.
\end{itemize}
\begin{figure}[tb]
\centering
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/trianglecomplex-1}
\label{subfig:trianglecomplex-1}}
\quad
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/trianglecomplex-2}
\label{subfig:trianglecomplex-2}}
\caption[Extruding the embeddings]{Extruding the embeddings of an $i$-cell $c$ along a single interval $r = (a,b)$ such that $a,b \in \mathbb{R}, a \neq b$, generates the embeddings of three cells: two $i$-cells $c^{a}$ and $c^{b}$, and an $(i+1)$-cell $c^{r} = c^{(a,b)}$ lying between them.}
\label{fig:trianglecomplex}
\end{figure}
\begin{figure}[tb]
\centering
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/blockscomplex-1}
\label{subfig:blockscomplex-1}}
\quad
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/blockscomplex-2}
\label{subfig:blockscomplex-2}}
\caption[Extruding a 2D cell complex]{Extruding a 2D cell complex using an interval that is defined per 2-cell.
The facet $f_0$ is to be extruded along the interval $r_0 = (a,c)$ and $f_1$ along $r_1 = (a,b)$.
Note how the vertices and edges incident to multiple facets in the input are extruded along the intervals of these facets, generating a series of cells connecting the end points of their intervals.}
\label{fig:blockscomplex}
\end{figure}
For the combinatorial structure, the algorithm takes advantage of the fact that the darts in the generalised map $G$ can be extruded largely independently based on querying the combinatorial structure of $G$ and simple 1D geometric queries along the $n$-th dimension.
Intuitively, the extrusion of a single $(n-1)$-simplex in $G$ consists of layers of $n$-simplices that are `stacked'
% \footnote{Not in the sense of a stack as a data structure, but as a pile of simplices arranged vertically along the $n$-th dimension.}
so as to form one part of the prism-shaped output.
Each new layer corresponds to a new $n$-simplex that shares all but one of the nodes with the $n$-simplex below it.
The extrusion algorithm is thus divided into three parts that are performed sequentially and explained in detail in the following sections: (1) propagating the input intervals to all the cells in the input cell complex, (2) generating the new embeddings (\ie\ the attributes and geometry) for each input cell, and (3) generating the combinatorial structure (\ie\ the darts and involutions) and linking each dart to its correct embeddings for every dimension.
\subsection{Propagating the extrusion intervals to all cells}
\label{ss:propagation}
While the extrusion intervals are defined only for the $(n-1)$-cells, the cells of every dimension need to be extruded, and therefore the extrusion intervals need to be propagated from the $(n-1)$-cells to the cells of lower dimension.
This is done recursively in decreasing dimension, using the incidence relationships between every $i$-cell and the $(i-1)$-cells on its boundary to pass the intervals of the former to the latter, using the same map of extrusion intervals $\rho$.
Because incidence is a transitive relation, an interval attached to a particular lower dimensional cell thus indicates that it is incident to an $(n-1)$-cell that will be extruded along that interval.
As $(i-1)$-cells can be on the boundary of multiple $i$-cells, several intersecting intervals can be passed to the map of extrusion intervals of a lower dimensional cell.
In order to generate non-intersecting cells, these need to be split into a set of non-intersecting intervals as shown in \reffig{fig:intersectingintervals}, each of which will represent the incidence of this cell to an equal set of higher-dimensional cells along the entirety of the interval.
\marginpar{
\captionsetup{type=figure}
\centering
\includegraphics[width=\marginparwidth]{figs/intersectingintervals}
\caption[Propagating the extrusion intervals]{The extrusion intervals from the cells incident to the middle vertex of \reffig{fig:extrusion-intervals} in (a), are passed on to it, resulting in a new set of non-intersecting intervals in (b).}
\label{fig:intersectingintervals}
}
A sketch of a simple process to do this is shown in \refalgo{algo:propagateintervals}, which assumes that the sets of intervals are kept sorted.
Note that this is essentially a one dimensional analogue of other operations used in GIS, \eg\ vector map overlays \citep[\S{}2.3]{deBerg08} and time-composites in spatio-temporal modelling \citep{Langran88}.
% \begin{figure}
\begin{algorithm}[tb]
\small
\SetKwInOut{Input}{Input}
\SetKwInOut{Output}{Output}
\DontPrintSemicolon\Input{generalised map $G$ of the input cell complex,\\
map $\rho$ with the extrusion intervals of the $(n-1)$-cells in $G$}
\Output{map $\rho$ with the extrusion intervals of all the cells in $G$}
\For{$i \leftarrow n-1$ \KwTo$0$}{
\ForEach{$i$-cell $c$ in $G$}{
\ForEach{$(i-1)$-cell $b$ on the boundary of $c$}{
\ForEach{extrusion interval $r = (r_{\min}, r_{\max}) \in \rho(c)$}{
Find the intervals $R^\prime \subseteq \rho(b)$ whose interiors intersect with $r$ \;
\If{$r_{\min}$ is in the interior of an interval $r^\prime = (r_{\min}^\prime, r_{\max}^\prime) \in R^\prime$}{
Remove $r^\prime$ from $\rho(b)$ \;
Add $(r_{\min}^\prime, r_{\min})$ to $\rho(b)$ \;
Add $(r_{\min}, r_{\max}^\prime)$ to $\rho(b)$ \;
}
\If{$r_{\min}$ is outside all intervals in $R^\prime$}{
Add a new interval in $R^\prime$ from $r_{\min}$ to the minimum of the lowest interval in $R^\prime$ \;
}
\If{$r_{\max}$ is in the interior of an interval $r^\prime = (r_{\min}^\prime, r_{\max}^\prime) \in R^\prime$}{
Remove $r^\prime$ from $\rho(b)$ \;
Add $(r_{\min}^\prime, r_{\max})$ to $\rho(b)$ \;
Add $(r_{\max}, r_{\max}^\prime)$ to $\rho(b)$ \;
}
\If{$r_{\max}$ is outside all intervals in $R^\prime$}{
Add a new interval in $R^\prime$ from the maximum of the highest interval in $R^\prime$ to $r_{\max}$ \;
}
\ForEach{empty interval $(r_{\min}^\prime, r_{\max}^\prime)$ between consecutive intervals in $R^\prime$}{
Add $(r_{\min}^\prime, r_{\max}^\prime)$ to $\rho(b)$ \;
}
}
}
}
}
\caption{\textsc{PropagateRanges}}\label{algo:propagateintervals}
\end{algorithm}
% \end{figure}
Considering what happens when a single new interval is passed to a cell in the incremental process shown in \refalgo{algo:propagateintervals}, the total number of intervals of that cell can increase from $j$ to $2j+1$ in the worst case, which happens when the new interval starts before and ends after all other intervals and these are all not contiguous (\ie\ they are all separated by empty intervals).
However, when taking into account $r$ intervals passed to the same cell, the total number of intervals for the cell can only increase to $2r-1$.
As each extrusion interval is bounded by two values, $r$ intervals lead to at most $2r$ boundary values and thus the number of intervals cannot be higher than $2r-1$.
Assuming that a new interval can be added to a cell in $O(\log r)$ time, \eg\ by maintaining a sorted list or an augmented red-black tree with the endpoints of the intervals, the total time to construct the ordered set of all intervals for a cell $c$ is $O(r \log r)$, where $r$ is the number of intervals that are passed onto $c$.
The overall computational complexity of this step depends on the incidence relationships between the cells in the complex.
\subsection{Generating the new embeddings}
\label{ss:embeddings}
Intuitively, when a single $i$-cell is extruded along a single interval $r = (r_{\min}, r_{\max})$, three new cells are generated: two $i$-cells that correspond to the end points of the interval $r_{\min}$ and $r_{\max}$, and a prismatic $(i+1)$-cell along all of $r$, \ie\ lying between the two $i$-cells.
For instance, extruding a polygon results in `top' and `bottom' polygonal faces and a polyhedron that lies between the two polygons.
Analogously, considering the \emph{embedding} structures that store the necessary attributes for each cell, such that an $i$-dimensional embedding ($i$-embedding) structure keeps the attributes of an $i$-cell, \emph{extruding an $i$-dimensional embedding results in two $i$-embeddings and one $(i+1)$-embedding}.
In the case of linear geometries, extruding a point entails the creation of two additional points, one with an appended $r_{\min}$ coordinate, and one with an $r_{\max}$ one.
With multiple intervals, this same procedure can then be used by applying it per embedding and per interval.
Notice that some embeddings are shared by multiple intervals (when the minimum of an interval is equal to the maximum of another), but these only need to be created once.
The extrusion algorithm for the embeddings receives the set of input embeddings $E$ and the map of extruded intervals $\rho$ which, being the output of \textsc{PropagateRanges} (\refalgo{algo:propagateintervals}), now contains a set of non-intersecting intervals for the cells of every dimension, and it returns an entirely new set of extruded embeddings $E^{\prime}$, as well a function $ex(e, v) \rightarrow E^\prime$ linking an input embedding $e \in E$ and an interval or end point of an interval $v$ to an extruded (output) embedding $e^{\prime} \in E^{\prime}$.
For instance, given an interval $r = (r_{\min}, r_{\max}) \in \rho(e)$, $v$ can be $r$, $r_{\min}$ or $r_{\max}$, reflecting the fact that extruding an $i$-embedding results in two $i$-embeddings (respectively for $r_{\min}$ and $r_{\max}$) and one $(i+1)$-embedding (for $r$).
The general procedure to generate the new embeddings and their relation to the old embeddings and the intervals is shown in \refalgo{algo:embeddingsextrusion}.
Note that a practical implementation of this algorithm has to deal with the desired attributes for each of the extruded cells rather than simply making a copy of the input ones (lines~\ref{ln:copyattributemin},~\ref{ln:coordinatemin} and~\ref{ln:copyattributemax}) and appending one more coordinate to the 0-embeddings (lines~\ref{ln:coordinatemin} and~\ref{ln:coordinatemax}).
% \begin{figure*}[tb]
\begin{algorithm}[tb]
\small
\SetKwInOut{Input}{Input}
\SetKwInOut{Output}{Output}
\DontPrintSemicolon\Input{set $E$ of the embeddings in the input cell complex,\\
map $\rho$ of the extrusion intervals for all cells of the input complex}
\Output{set $E^\prime$ of the embeddings for the output cell complex,\\
map $ex$ that links an input embedding to its extruded embeddings}
\ForEach{$e \in E$}{
\ForEach{$r = (r_{\min}, r_{\max}) \in \rho(e)$}{
\If{$ex(e, r_{\min}) = \emptyset$}{
$ex(e, r_{\min}) \leftarrow e$ \;
\label{ln:copyattributemin}
\If{$e.dimension = 0$}{
Append $r_{\min}$ to the coordinates of $ex(e, r_{\min})$ \;
\label{ln:coordinatemin}
}
}
$ex(e, r) \leftarrow e$ \;
\label{ln:copyattributeinterval}
$ex(e, r).dimension \leftarrow ex(e, r).dimension+1$ \;
\If{$ex(e, r_{\max}) = \emptyset$}{
$ex(e, r_{\max}) \leftarrow e$ \;
\label{ln:copyattributemax}
\If{$e.dimension = 0$}{
Append $r_{\max}$ to the coordinates of $ex(e, r_{\max})$ \;
\label{ln:coordinatemax}
}
}
Put $ex(e, r_{\min}), ex(e, r)$ and $ex(e, r_{\max})$ in $E^{\prime}$ \;
}
}
\caption{\textsc{EmbeddingsExtrusion}}\label{algo:embeddingsextrusion}
\end{algorithm}
% \end{figure*}
Since this part of the algorithm iterates through all the cells in the input cell complex, and it generates at most three new embeddings per interval for each cell, the computational complexity can be $O(rn)$ per cell, where $r$ is the total number of intervals and $n$ is the dimension, \ie\ as long as the map $ex$ can be queried and the new embeddings can be created in time that is linear on the dimension (which in practice is a relatively small constant).
\subsection{Generating the new combinatorial structure and linking it to its correct embeddings}
\label{ss:combinatorial}
The extrusion of a \emph{single} dart in the input map $G$ along a single interval $r = (r_{\min}, r_{\max})$ generates a series of connected darts in $G^\prime$ connected by a sequence of involutions $\alpha_{n-1}, \alpha_{n-2}, \ldots, \alpha_{1}, \alpha_{0}, \alpha_{1}, \ldots, \alpha_{n-2}, \alpha_{n-1}$.
Intuitively, these are equivalent to stacked $n$-simplices that together form a prism\footnote{This is always true combinatorially, but it might not be true geometrically if the $(n-1)$-simplex is not embedded so as to lie in the interior of the cell, as shown in \refse{ss:ordered-topological-models}.}.
\reffig{fig:triangle} shows a simple extrusion of the cell complex from \reffig{fig:trianglecomplex} including all its darts.
\begin{figure}[tbp]
\centering
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/triangle-1}
\label{subfig:triangle-1}}
\quad
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/triangle-2}
\label{subfig:triangle-2}}
\caption[The darts in the cell complexes in \reffig{fig:trianglecomplex}]{The darts in the cell complexes in \reffig{fig:trianglecomplex}.
Note that one can obtain a dart's representation as a simplex by considering additional nodes at its corresponding facet and one in the interior of the volume.}
\label{fig:triangle}
\end{figure}
\marginpar{
\captionsetup{type=figure}
\centering
\setcounter{subfigure}{6}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-6}\label{subfig:stack-6}}\\
\vspace{-10mm}
\setcounter{subfigure}{5}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-5}\label{subfig:stack-5}}\\
\vspace{-8mm}
\setcounter{subfigure}{4}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-4}\label{subfig:stack-4}}\\
\vspace{-6mm}
\setcounter{subfigure}{3}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-3}\label{subfig:stack-3}}\\
\vspace{-5mm}
\setcounter{subfigure}{2}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-2}\label{subfig:stack-2}}\\
\vspace{-6mm}
\setcounter{subfigure}{1}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-1}\label{subfig:stack-1}}\\
\vspace{-3mm}
\setcounter{subfigure}{0}\subfloat[]{\includegraphics[width=0.7\marginparwidth]{figs/stack-0}\label{subfig:stack-0}}\\
\caption[The stack of simplices in an extruded dart]{The extrusion of a single dart (blue) results in a stack of simplices that together form a prism.
In order to visualise it more intuitively, the procedure here is shown from bottom (a) to top (g).
}
\label{fig:stack}
}
Since two darts connected by an $\alpha_i$ involution share all but the $i$-th node, and the $i$-th node in a generalised map means that a dart belongs to a certain $i$-cell\footnote{In this specific case, where cells have linear geometries.}, this series of darts represents a succession of simplices that progressively change from the cells at the `base' (at $r_{\min}$) to those at the `side' (at the interval $(r_{\min}, r_{\max})$) to those at the `top' (at $r_{\max}$).
This can be more clearly seen in \reffig{fig:stack}.
As the pattern that is followed in each stack of darts in $G^\prime$ is the same for all the darts in $G$, assuming that their extrusion intervals are the same, they can also be expressed in terms of the dart in $G$ being extruded and their position in the stack.
The darts at a certain level in the stack are thus called a \emph{layer}, and various functions will map an input dart in $G$ to the equivalent dart at a specific layer in $G^\prime$.
An intuitive justification for this is that if one considers the darts in the base and top facets in the extrusion of a cell $e$ along a single interval $r = (r_{\min}, r_{\max})$, these belong to different cells of every dimension except for their $n$-cell---which is $ex(e, r)$.
Considering a function $e_i: D \rightarrow E$ that links a dart $d \in D$ to its $i$-embedding in $E$, for every dimension $i < n$, the darts of the base facet belong to the extruded $ex\left(e_i\left(d\right), r_{\min}\right)$ cells, while the top darts belong to the extruded $ex\left(e_i\left(d\right), r_{\max}\right)$ cells.
The darts of the faces on the sides, which connect corresponding ridges of the base and top, belong instead to a mixture of cells in $ex\left(e_i\left(d\right), r_{\min}\right)$, $ex\left(e_i\left(d\right), r_{\max}\right)$ and $ex\left(e_i\left(d\right), r\right)$.
The darts closer to the base face belong to cells in $ex\left(e_i\left(d\right), r_{\min}\right)$ and $ex\left(e_i\left(d\right), r\right)$, while those closer to the top belong to cells in $ex\left(e_i\left(d\right), r_{\max}\right)$ and $ex\left(e_i\left(d\right), r\right)$.
This is natural when one considers that $\alpha_i$-linked darts differ only in one of their cells (\ie\ the $i$-cells), all other cells being the same for both darts.
This means that there must be a natural progression of layers of darts: starting from the base, they change their cells one by one from $ex\left(e_i\left(d\right), r_{\min}\right)$ to $ex\left(e_i\left(d\right), r\right)$ from the highest dimension down, and then change their cells one by one from $ex\left(e_i\left(d\right), r\right)$ to $ex\left(e_i\left(d\right), r_{\max}\right)$ from the lowest dimension up until reaching the top.
The algorithm to generate the extruded combinatorial structure therefore works by generating layers of darts that follow a pattern based on the darts in the input cell complex, starting from those for the base face, moving on to the $2n-2$ layers for the side faces, and finishing with the top face.
The output darts and the involutions between them are expressed in terms of the input generalised map $G$.
For this, a function $cur: G \rightarrow G^\prime$ maps a dart in $G$ to the corresponding dart in the currently being generated layer of the output map $G^\prime$.
Note that this means that the algorithm needs to keep track of only two layers of darts at any given time, one in $G$ and one in $G^\prime$, each of which has at most the number of darts in the input map.
This bounds the memory usage of this part of the algorithm, which is therefore on the order of $O(d)$, with $d$ the number of darts in the input space partition.
As shown by \citet{Ferrucci93}, when multiple intervals are involved this procedure can simply be repeated for all intervals, assuming that these have all been subdivided so as not to intersect one another.
However, it is possible to greatly reduce the number of darts generated by skipping the creation of some of the darts.
This is possible because the extrusion intervals have been independently propagated to each cell of every dimension.
Based on the algorithm, the lower-dimensional cells in the complex have received all the intervals from their incident higher-dimensional cells so that the intervals in the lower-dimensional cells contain all the endpoints of the intervals of their incidences.
However, the same is not true in the opposite direction: the intervals for the higher-dimensional cells have not been subdivided so as to contain all the ones of their incident lower-dimensional cells.
Nevertheless, as \reffig{fig:blocks} shows, even when an extrusion interval that would be used by the pattern described above is not in a cell, it is possible to map it to a \emph{bigger} interval that contains it.
If one considers the darts that would be generated using the above mentioned approach for all the non-intersecting intervals in the cell complex, but at the same time mapping the non-existent extrusion intervals to their immediately bigger containing ones, this can result in many darts that are equivalent (\ie\ they have nodes at the same cells).
These darts are those that can be skipped during the extrusion process.
\begin{figure*}[tbp]
\centering
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/blocks-1}
\label{subfig:blocks-1}}
\quad
\subfloat[]{
\includegraphics[width=0.45\linewidth]{figs/blocks-2}
\label{subfig:blocks-2}}
\caption[The darts in the cell complexes in \reffig{fig:blockscomplex}]{The darts in the cell complexes in \reffig{fig:blockscomplex}.
Note that the darts of the right cube in (b), when viewed as individual stacks, are similar to the ones in \reffigs{fig:triangle} and~\ref{fig:stack}.
The darts of the left box are however different: those on the left all involve an extrusion along the interval $(a,c)$ due to the fact that there is no vertex, edge or facet extruded to $b$.
On the other hand, those in the right do have a vertex and an edge at $b$, but not a facet and still belong to the same volume.}
\label{fig:blocks}
\end{figure*}
Therefore, in order to generate the darts for all intervals of all cells, this procedure is repeated for all the intervals in $\rho$, after making them non-intersecting using the same procedure delineated in \refse{ss:propagation}, and skipping all the darts that would be equivalent to those that have already been created.
This is done using a sweep-hyperplane-like algorithm that generates up to $n$ layers of darts at the events at the beginning or end of an interval, creating darts only when the sweep-hyperplane\footnote{\ie\ a $(n-1)$-dimensional shape that is unbounded along $n-1$ dimensions in $n$D space, \eg\ a line in $\mathbb{R}^2$ or a plane in $\mathbb{R}^3$.} passes by the beginning or end of their extrusion intervals (\ie\ not while it is in their interior).
The darts are linked to their appropriate embedding according to the same pattern.
The complete procedure to generate the combinatorial structure is presented in \refalgo{algo:gmapextrusion}, which uses the \refalgo{algo:gmaplayerbase} when the sweep hyperplane passes by the beginning of an interval and the similar \refalgo{algo:gmaplayertop} when the sweep hyperplane by the end of one.
If an input dart is denoted as $d$, its corresponding dart in the currently being generated layer of the output is denoted as $cur(d)$, the dart linked to $d$ by an $i$-involution as $\alpha_{i}(d)$, and the $i$-embedding of $d$ as $e_{i}(d)$.
% \begin{figure}
\begin{algorithm}[tb]
\small
\SetKwInOut{Input}{Input}
\SetKwInOut{Output}{Output}
\DontPrintSemicolon\Input{generalised map $G$ of the input cell complex,\\
set $E$ of the embeddings in the input cell complex,\\
set $E^\prime$ of the embeddings for the output cell complex,\\
map $\rho$ of the extrusion intervals for all cells of the input complex,\\
map $ex$ that links an input embedding to its extruded embeddings}
\Output{generalised map $G^\prime$ of the output cell complex}
Compute an ordered set of non-intersecting intervals $r_{all}$ using all the intervals for all cells in $\rho$ \;
Consider a sweep plane that passes through all the intervals $r_{all}$ in increasing order along dimension $n$ \;
\If{the sweep plane passes by the beginning of an interval $r = (r_{\min}, r_{\max}) \in r_{all}$}{
\For{$i \leftarrow n$ \KwTo$0$}{
$\textsc{GMapLayerBegin}(G, G^{\prime}, i, E, E^{\prime}, \rho, r, ex)$ \;
}
}
\If{the sweep plane passes by the end of an interval $r = (r_{\min}, r_{\max}) \in r_{all}$}{
\For{$i \leftarrow 0$ \KwTo$n$}{
$\textsc{GMapLayerEnd}(G, G^{\prime}, i, E, E^{\prime}, \rho, r, ex)$ \;
}
}
\caption{\textsc{GMapExtrusion}}\label{algo:gmapextrusion}
\end{algorithm}
% \end{figure}
% \begin{figure}[tb]
\begin{algorithm}[tb]
\small
\SetKwInOut{Input}{Input}
\SetKwInOut{Output}{Output}
\DontPrintSemicolon\Input{generalised map $G$ of the input cell complex,\\
generalised map $G^\prime$ of the output cell complex,\\
dimension $i$ of the current layer,\\
% map $last$ that links a dart in $G$ to a dart in the most recent layer of $G^\prime$,\\
set $E$ of the embeddings in the input cell complex,\\
set $E^\prime$ of the embeddings of the output cell complex,\\
map $\rho$ of the extrusion intervals of all cells of the input complex,\\
current interval $r$,\\
map $ex$ that links an input embedding to its extruded embeddings}
\Output{generalised map $G^\prime$ of the output cell complex}
\ForEach{dart $d$ in $G$}{
% \If{$r \subseteq \rho(e_{n-1}(d))$}{
% \If{$\exists r^\prime \mid r \subseteq \rho(e_{n-1}(d), r^\prime)$}{
\If{$\exists r^\prime \in \rho(e_{n-1}(d)) \mid r \subseteq r^\prime$}{
% \If{$\rho(e_{i}(d), r_{\min}) \neq \emptyset$}{
\If{$\exists r^{\prime\prime} \in \rho(e_{i}(d)) \mid r_{\min}^{\prime\prime} = r_{\min}$}{
$last \leftarrow cur(d)$ \;
$cur(d) \leftarrow$ new dart \;
Put $cur(d)$ in $G^{\prime}$ \;
$\alpha_{i+1}(cur(d)) \leftarrow last$ \;
$\alpha_{i+1}(last) \leftarrow cur(d)$ \;
}
}
}
\ForEach{dart $d$ in $G$}{
\If{$\exists r^\prime \in \rho(e_{n-1}(d)) \mid r \subseteq r^\prime$}{
\If{$\exists r^{\prime\prime} \in \rho(e_{i}(d)) \mid r_{\min}^{\prime\prime} = r_{\min}$}{
\For{$inv \leftarrow 0$ \KwTo$i-1$}{
$\alpha_{inv}^{\prime}(cur(d)) \leftarrow cur(\alpha_{inv}(d))$ \;
$\alpha_{inv}^{\prime}(cur(\alpha_{inv}(d))) \leftarrow cur(d)$ \;
}
\For{$inv \leftarrow i+2$ \KwTo$n$}{
$\alpha_{inv}^{\prime}(cur(d)) \leftarrow cur(\alpha_{inv-1}(d))$ \;
$\alpha_{inv}^{\prime}(cur(\alpha_{inv-1}(d))) \leftarrow cur(d)$ \;
}
\For{$emb \leftarrow 0$ \KwTo$i$}{
$e_{emb}^{\prime}(cur(d)) \leftarrow ex(e_{emb}(d), r_{\min})$ \;
}
\For{$emb \leftarrow i+1$ \KwTo$n$}{
$e_{emb}^{\prime}(cur(d)) \leftarrow ex(e_{emb-1}(d), r)$ \;
}
}
}
}
\caption{\textsc{GMapLayerBegin}}\label{algo:gmaplayerbase}
\end{algorithm}
% \end{figure}
% \begin{figure}[tb]
\begin{algorithm}[tb]
\small
\SetKwInOut{Input}{Input}
\SetKwInOut{Output}{Output}
\DontPrintSemicolon\Input{generalised map $G$ of the input cell complex,\\
generalised map $G^\prime$ of the output cell complex,\\
dimension $i$ of the current layer,\\
% map $last$ that links a dart in $G$ to a dart in the most recent layer of $G^\prime$,\\
set $E$ of the embeddings in the input cell complex,\\
set $E^\prime$ of the embeddings of the output cell complex,\\
map $\rho$ of the extrusion intervals of all cells of the input complex,\\
current interval $r$,\\
map $ex$ that links an input embedding to its extruded embeddings}
\Output{generalised map $G^\prime$ of the output cell complex}
\ForEach{dart $d$ in $G$}{
\label{ln:newlayerbegin}
% \If{$\exists r^\prime \mid r \subseteq \rho(e_{n-1}(d), r^\prime)$}{
\If{$\exists r^\prime \in \rho(e_{n-1}(d)) \mid r \subseteq r^\prime$}{
% \If{$\rho(e_{i}(d), r_{\max}) \neq \emptyset$}{
\If{$\exists r^{\prime\prime} \in \rho(e_{i}(d)) \mid r_{\max}^{\prime\prime} = r_{\max}$}{
$last \leftarrow cur(d)$ \;
$cur(d) \leftarrow$ new dart \;
Put $cur(d)$ in $G^{\prime}$ \;
$\alpha_{i}(cur(d)) \leftarrow last$ \;
$\alpha_{i}(last) \leftarrow cur(d)$ \;
}
}
}
\label{ln:newlayerend}
\ForEach{dart $d$ in $G$}{
\If{$\exists r^\prime \in \rho(e_{n-1}(d)) \mid r \subseteq r^\prime$}{
\If{$\exists r^{\prime\prime} \in \rho(e_{i}(d)) \mid r_{\max}^{\prime\prime} = r_{\max}$}{
\For{$inv \leftarrow 0$ \KwTo$i-1$}{
\label{ln:setinvolutionsbegin}
$\alpha_{inv}^{\prime}(cur(d)) \leftarrow cur(\alpha_{inv}(d))$ \;
$\alpha_{inv}^{\prime}(cur(\alpha_{inv}(d))) \leftarrow cur(d)$ \;
}
\For{$inv \leftarrow i+2$ \KwTo$n$}{
$\alpha_{inv}^{\prime}(cur(d)) \leftarrow cur(\alpha_{inv-1}(d))$ \;
$\alpha_{inv}^{\prime}(cur(\alpha_{inv-1}(d))) \leftarrow cur(d)$ \;
}
\label{ln:setinvolutionsend}
\For{$emb \leftarrow 0$ \KwTo$i$}{
\label{ln:setembeddingsbegin}
$e_{emb}^{\prime}(cur(d)) \leftarrow ex(e_{emb}(d), r_{\max})$ \;
}
\For{$emb \leftarrow i+1$ \KwTo$n$}{
$e_{emb}^{\prime}(cur(d)) \leftarrow ex(e_{emb-1}(d), r)$ \;
}
\label{ln:setembeddingsend}
}
}
}
\caption{\textsc{GMapLayerEnd}}\label{algo:gmaplayertop}
\end{algorithm}
% \end{figure}
For the latter two algorithms, lines~\ref{ln:newlayerbegin}--\ref{ln:newlayerend} show the generation of a new layer of darts and its linking to the previous one, lines~\ref{ln:setinvolutionsbegin}--\ref{ln:setinvolutionsend} show how the darts within the layer are linked based on the pattern of the input map, and lines~\ref{ln:setembeddingsbegin}--\ref{ln:setembeddingsend} show how the darts within the layer are related to their correct embeddings, which are also patterned after the embeddings in the input.
Notice that this method generates layers of darts in a grid-like fashion, each layer containing at most the number of darts in the input map, and calling \textsc{GMapLayerBegin} and \textsc{GMapLayerEnd} to create at most $2n$ layers per non-intersecting interval.
The time complexity of computing the set of non-intersecting intervals in \refalgo{algo:propagateintervals} is $O(r \log r)$ as before, while the number of darts in the output map is bounded by $O(ndr)$, where $n$ is the extrusion dimension, $d$ is the total number of darts in the input map and $r$ is the total number of intervals in the input.
\section{Implementation}
\label{se:extrusion-implementation}
The extrusion algorithm has been implemented in C++11 and made available under the open source MIT licence at \url{https://github.com/kenohori/lcc-tools}.
It requires and builds upon the CGAL packages Combinatorial Maps and Linear Cell Complex, among others.
The first package provides data structures and algorithms to store and to efficiently iterate over the darts of a combinatorial map, and the second links the 0-embeddings to other CGAL types in order to store the geometry of a model.
As the algorithm is described on the basis of generalised maps while the CGAL packages instead use combinatorial maps, all operations are done through a small wrapper that converts operations for the former to operations for the latter.
Since a combinatorial maps dart is equivalent to two generalised map darts (see \refse{ss:ordered-topological-models}), some redundant operations are ignored in the process, such as requests to create darts that already exist or to connect darts that are already connected.
These redundant operations are instead used as assertions to verify the validity of the map during its creation.
Both these packages and the implementation of the extrusion algorithm make heavy use of the traits programming technique \citep{Myers95} and recursive templates (TMP or template meta-programming) in order to produce efficient code.
These techniques are described in \refse{se:traits}.
One shortcoming of the current prototype implementation is that it uses the C++ type \texttt{std::map} in order to link cells to their extrusion intervals, which offers only logarithmic time access \citep[\S{}23.4]{ISO14882:2015}\footnote{See \citet{Austern00} for the consequences of this.}, rather than the constant time that would be possible by integrating this into the templated structures.
This would involve storing the set of extrusion intervals of a cell directly in the data structure that is used for its attributes.
In order to input and output data, as well as to visualise the results, the implementation uses the OGR Simple Feature Library\footnote{\url{http://gdal.org}} to read standard GIS data formats.
It is also able to output Wavefront OBJ\footnote{\url{http://www.martinreddy.net/gfx/3d/OBJ.spec}} files, in which faces or darts can be directly exported for visualisation.
In the latter case, darts are exported as triangles with vertices at the two 0-cells of a dart, as well as an added vertex at the centroid of its 2-cell.
\section{Experiments}
\label{se:extrusion-experiments}
The algorithm has been tested by extruding various 2D datasets to higher dimensions.
For this, several free and open datasets in the area of Delft are used, matching the geometries in some with the attributes present in others so as to obtain new attributes and appropriate extrusion intervals for each geometry.
The tests were performed on a Mac OS X computer with a 2.7 GHz Intel Core 2 Duo processor and 12 GB of RAM.\@
The main characteristics of the datasets tested are shown in \reftab{table:tests}.
\begin{table*}[tbp]
\caption[Extrusion test results]{Characteristics of the tested datasets. The cell counts represent the total number of cells for all dimensions.}
\centering
\small
\begin{tabular}{lccccc}
\toprule
& \multicolumn{2}{c}{Input} & \multicolumn{2}{c}{Output} & \\
\cmidrule(lr){2-3}
\cmidrule(lr){4-5}
Test & darts & cells & darts & cells & time \\
\midrule
One GBKN building to 5D & 14 & 29 & 80 640 & 783 & 3 s \\
370 buildings to 4D & 3 268 & 6 065 & 123 184 & 42 552 & 12 s \\
TOP10NL to 3D & 30 098 & 48 562 & 181 640 & 148 102 & 2 m 39 s\\
\bottomrule
\end{tabular}
\label{table:tests}
\end{table*}
Note however that since the implementation uses an \texttt{std::map} to access the extrusion intervals of each cell, the running time is dominated by the large number of times this query is performed (several times per interval and per cell).
The times provided in \reftab{table:tests} are therefore not indicative of the theoretical complexity of the algorithm, which is instead dominated by the generation and linking of the darts in the extruded dataset (\refse{ss:combinatorial}).
The generation of the non-intersecting intervals for all cells (\refse{ss:propagation}) and the generation of the extruded embeddings (\refse{ss:embeddings}) always ran in under a second, even for very large datasets.
\marginpar{
\captionsetup{type=figure}
\centering
\subfloat[]{
\includegraphics[width=\marginparwidth]{figs/aula-1}
\label{subfig:aula-1}}
\\
\subfloat[]{
\includegraphics[width=\marginparwidth]{figs/aula-2}
\label{subfig:aula-2}}
\caption[Extruding the footprint of the Aula Congress Centre in Delft to 3D]{Extruding the (a) footprint of the Aula Congress Centre in Delft to (b) a 3D representation showing individual darts. }
\label{fig:aula-ex}
}
One test involved the Aula Congress Centre in Delft, a building represented by a single polygon with 14 vertices extracted from the GBKN dataset\footnote{\url{http://www.gbkn.nl}}, as shown in \reffig{fig:aula-ex}.
It was extruded from 2D up to 5D using some manually added attributes: its height, construction date and a specified level of detail for the model.
\reffig{fig:aula4d} shows the extruded model at an intermediate stage in this process, in 4D.
The end result was a generalised map with 80 640 darts, 112 vertices, 280 edges, 260 facets, 110 volumes, 20 4-cells and 1 5-cell.
It was generated in 3 seconds using 15 MB of RAM.\@
\begin{figure*}[b]
\centering
\subfloat[]{
\includegraphics[width=0.33\linewidth]{figs/aula4d-ends}
\label{subfig:aula4d-ends}}
% \quad
\subfloat[]{
\includegraphics[width=0.33\linewidth]{figs/aula4d-sides}
\label{subfig:aula4d-sides}}
% \quad
\subfloat[]{
\includegraphics[width=0.33\linewidth]{figs/aula4d-top-bottom}
\label{subfig:aula4d-top-bottom}}
\caption[Extruding the footprint of the Aula Congress Centre in Delft to 4D]{The footprint of the Aula Congress Centre (\reffig{fig:aula-ex}) extruded twice to create a 4D model.
It is shown here in three parts for clarity: (a) the faces in the two end volumes, (b) the lateral faces connecting corresponding vertical edges, and (c) the top and bottom faces connecting corresponding horizontal edges.
Not shown here are the 16 different polyhedra that bound the polychoron and which are all bounded by different sets of these faces.}
\label{fig:aula4d}
\end{figure*}
Another test, shown in 2D and 3D in \reffig{fig:campus}, involved a previously generated dataset of the campus of the Delft University of Technology (see \citet{Ledoux11a}), consisting of 2 749 vertices, 2 946 edges, and 370 facets.
\begin{figure}[tb]
\centering
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/campus}
\label{subfig:campus}}
% \quad
\subfloat[]{
\includegraphics[width=0.5\linewidth]{figs/campusextruded}
\label{subfig:campusextruded}}
\\
\subfloat[]{
\includegraphics[width=\linewidth]{figs/campusdarts}
\label{subfig:campusdarts}}
\caption{Extruding a dataset of the campus of the Delft University of Technology}
\label{fig:campus}
\end{figure}
This dataset, covering 2.3 km$^{2}$, was originally built from the GBKN dataset by manually forming footprint polygons from the lines in the dataset.
Building heights for each polygon were obtained from the AHN dataset\footnote{\url{http://www.ahn.nl}} (airborne laser altimetry), while building dates were obtained from the BAG dataset\footnote{\url{https://www.kadaster.nl/bag}}.
This dataset, including the added attributes, is available together with the source code of the program.
It was therefore extruded from its original 2D representation to 4D using building heights for the third dimension and dates for the fourth dimension.
The result was a generalised map with 123 184 darts, 8 613 vertices, 17 919 edges, 12 471 facets, 3 310 volumes and 239 4-cells.
It was generated in 12 seconds using 46 MB of RAM.\@
One more test used 1 836 buildings from Delft from the TOP10NL dataset\footnote{\url{http://www.kadaster.nl/web/artikel/producten/TOP10NL.htm}}, which was extruded to 3D using intervals from building dates also obtained from the BAG dataset.
The result was a generalised map with 181 640 darts, 46 464 vertices, 71 826 edges, 27 975 facets and 1 837 volumes.
An algorithm was used to verify that the constructed datasets conform to the definition of a combinatorial map\footnote{\ie\ checking that the links between darts correctly form partial permutations or involutions, and all the darts in the orbit of an $i$-cell are linked to its correct $i$-attribute and vice versa.}.
Additional tests were made to ensure that all the darts of an $i$-cell correctly point to it in their $i$-embeddings, and to verify that all the darts of an $n$-cell are linked to point embeddings within the extrusion interval given for the $n$-cell, among other tests.
The extruded datasets were also inspected visually in 2D and 3D by exporting 2-cells as polygons and verifying that they form a valid cell complex, \ie\ that 2-cells intersect only at their common boundaries, forming 1-cells that are also in the complex.
In 2D and 3D, individual darts were also exported as triangles to visually verify that they form a valid generalised map (as shown in \reffigs{subfig:aula-2} and~\ref{subfig:campusdarts}).
\section{Extrusion-based generalisation operations}
\label{se:extrusion-generalisation}
Apart from using extrusion to directly create finished higher-dimensional datasets, it is also interesting to explore the possibility to use it as a base for further operations.
For instance, it can be used as the basis of a simple set of dimension-independent generalisation operations.
Starting from a detailed $n$-dimensional cell complex, it is possible to first extrude the complex to create an $(n+1)$-dimensional cell complex, then to obtain the desired geometry by applying certain transformations to the vertices of the `top' or `bottom' facet of the extruded model, even if this procedure results in degenerate cells, and finally to remove the degenerate cells from the model.
The transformations that can be applied to the model include all those that can be applied through transformation matrices, such as those shown in \reffigp{fig:aula-transformations} (translation, rotation and scale), as well as others such as reflections and shears.
Perhaps more importantly, it is possible to gradually reduce the complexity of the model by collapsing cells of any dimension, which is achieved by simply moving all vertices of a cell to the same location, as is shown in \reffig{fig:aula-collapse}.
\marginpar{
\captionsetup{type=figure}
\centering
\subfloat[]{
\includegraphics[width=\marginparwidth]{figs/aula-collapse}
\label{subfig:aula-collapse}}
\\
\subfloat[]{
\includegraphics[width=\marginparwidth]{figs/aula-edge-collapse}
\label{subfig:aula-edge-collapse}}
\caption[Collapsing cells]{Cells can be collapsed by moving all of their vertices to the same location. Shown here are the collapse of: (a) a face and (b) an edge (in red).}
\label{fig:aula-collapse}
}
After such transformations, it is possible to easily detect and remove degenerate cells, as they can be identified because they have all their vertices at the same location.
In the case of a topological model, such as a generalised or combinatorial map, the topological relationships around the collapsed cells will have to be recomputed.
However, this is possible to achieve either by finding the changed incidences and adjacencies, or by regenerating all topological relationships using the method described in \refch{ch:incremental-construction}.
For the latter, it is simply necessary to ignore all degenerate cells (as defined by having all their vertices at the same embedding in $\mathbb{R}^n$).
\section{Conclusions}
\label{se:extrusion-conclusions}
Extrusion in the GIS sense has a natural extension into a dimension-independent formulation.
It can be used to load existing 2D or 3D data into a higher-dimensional structure, either directly or as a base for further (dimension-independent) operations, such as the generalisation operations described in \refse{se:extrusion-generalisation}.
The conceptual simplicity of extrusion makes it particularly suitable for higher-dimensional data, as the complex problem of its generation can be reduced to the definition of extrusion intervals per $n$-cell in an $n$-dimensional cell complex.
By ensuring that the input $n$-cells do not overlap---which is a general precondition in a geometric cell complex---, and the extrusion intervals per cell also do not overlap, it is easy to guarantee that the output cell complex forms a valid $(n+1)$-dimensional space partition.
The dimension-independent extrusion formulation has been realised into an algorithm, presented in this chapter, which supports multiple intervals per input object and it is also memory efficient---only three layers of darts of the same size as the input map need to be kept in main memory at the same time.
It is also relatively fast, with a worst case complexity of $O(ndr)$ in the main algorithm, where $n$ is the extrusion dimension\footnote{In practice a small constant.}, $d$ is the total number of darts in the input map and $r$ is the total number of intervals in the input, but offers better complexity in practice.
Moreover, the algorithm has been implemented using CGAL Combinatorial Maps and the source code has been made publicly available under a permissive licence.
It has also been tested with publicly available datasets commonly used in GIS.\@
Apart from being applied cell by cell, it can also be combined with other operations, such as the generalisation operations described here or building certain $n$D cells incrementally, where more complex models are required only for a subset of the objects in a model.
While the implementation presented here performs well enough for typical GIS datasets, it could be made faster by integrating the extrusion intervals into the embeddings of the cells, allowing for constant time access to the extrusion intervals.
This would make it possible to support much larger datasets and make the implementation match the theoretical complexity of the algorithm.
In order to improve memory usage, another possibility would be implementing the algorithm in a more compact representation of a simplicial complex, \eg\ \citet{Boissonnat12}.
In the future, it would be interesting to consider other types of models that can be generated using extrusion-like algorithms, such as sweeps and $n$-dimensional cells of revolution.
These are both simple extensions either using either non-linear embeddings or discretised approximations with small cells with linear geometries.
A more complex extension would involve starting from a \emph{series of cell complexes} that each form a planar partition, such as series of 3D models of the same region at different levels of detail.
While most of the algorithm described here should work directly on this latter example, managing the differing boundaries of the cells along the $n$-th dimension could be challenging---possibly requiring the computation of all intersections in an arrangement of $(n-1)$-cells.

BIN
work_report/figs/0-cell-double-star.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/0-cell-embedding.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 70 KiB

BIN
work_report/figs/0-cell-star.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/0-cell.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/0-simplex.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/1-cell-boundary.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/1-cell-star.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/1-simplex-boundary.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/1-simplex-star.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/16tiles.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 109 KiB

BIN
work_report/figs/2-cell-boundary.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2-simplex-adjacency.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2-simplex-boundary.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2-simplex-vertices.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2d+scale.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 128 KiB

BIN
work_report/figs/2d+time.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 93 KiB

BIN
work_report/figs/2dcc-alphas.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2dcc-betas.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2dcc-gmap.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2dcc.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2tetra-before-3sew.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2tetra-map.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/2tetra.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/3-sew.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/3dgeoinfo.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/3dtop10nl.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 744 KiB

BIN
work_report/figs/4DslicingHouses.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/4dhouse.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 286 KiB

BIN
work_report/figs/Frauenkirche.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 458 KiB

BIN
work_report/figs/SlicingPrism2D.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/SlicingPrism2DResult.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/agg2.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/aluminum.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 339 KiB

BIN
work_report/figs/annulus.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/atkis.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 373 KiB

BIN
work_report/figs/aula-1.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 37 KiB

BIN
work_report/figs/aula-2.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 109 KiB

BIN
work_report/figs/aula-collapse.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 80 KiB

BIN
work_report/figs/aula-edge-collapse.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 83 KiB

BIN
work_report/figs/aula-extruded.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 67 KiB

BIN
work_report/figs/aula-moved.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 71 KiB

BIN
work_report/figs/aula-rotated.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 70 KiB

BIN
work_report/figs/aula-scaled.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 68 KiB

BIN
work_report/figs/aula.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 60 KiB

BIN
work_report/figs/aula4d-ends.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 96 KiB

BIN
work_report/figs/aula4d-sides.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 122 KiB

BIN
work_report/figs/aula4d-top-bottom.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 147 KiB

BIN
work_report/figs/axes-es.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/axes-nl.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/axes.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/barcode.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/bk.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 201 KiB

BIN
work_report/figs/bksep.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 94 KiB

BIN
work_report/figs/blocks-1.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/blocks-2.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/blockscomplex-1.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 42 KiB

BIN
work_report/figs/blockscomplex-2.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 50 KiB

BIN
work_report/figs/boolean-difference.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 61 KiB

BIN
work_report/figs/boolean-intersection.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 44 KiB

BIN
work_report/figs/boolean-union.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 89 KiB

BIN
work_report/figs/boolean.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 85 KiB

BIN
work_report/figs/brep.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/bridge-1.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/bridge-2.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/campus.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 90 KiB

BIN
work_report/figs/campusdarts.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 380 KiB

BIN
work_report/figs/campusextruded.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 104 KiB

BIN
work_report/figs/cellcomplex-cells.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cellcomplex-polygons.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cellscampus.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 168 KiB

BIN
work_report/figs/citygml-filip.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 118 KiB

BIN
work_report/figs/citygml.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cmap-3d-vertices.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 96 KiB

BIN
work_report/figs/cmap-3d.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 188 KiB

BIN
work_report/figs/cmap-free.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cmaps-orientation-1.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cmaps-orientation-2.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/cmaps-simplices.pdf Normal file
View File

Binary file not shown.

BIN
work_report/figs/conic-circle.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 90 KiB

BIN
work_report/figs/conic-ellipse.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 92 KiB

BIN
work_report/figs/conic-hyperbola.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 95 KiB

BIN
work_report/figs/conic-parabola.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 94 KiB

BIN
work_report/figs/construction-faces.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 69 KiB

BIN
work_report/figs/construction-vertices.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 65 KiB

BIN
work_report/figs/construction-volumes.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 75 KiB

47
work_report/figs/construction.cpp Normal file
View File

@ -0,0 +1,47 @@
// 2: Left of first house (yellow)
faceLists.push_back(...);
faceLists.back().push_back(vertices[0]);
faceLists.back().push_back(vertices[2]);
faceLists.back().push_back(vertices[6]);
faceLists.back().push_back(vertices[4]);
//21: Back right vertical edge (blue)
faceLists.push_back(...);
faceLists.back().push_back(vertices[3]);
faceLists.back().push_back(vertices[7]);
faceLists.back().push_back(vertices[15]);
faceLists.back().push_back(vertices[11]);
for (int i = 0; i < 35; ++i)
faces.push_back(
builder.get_facet_from_vertices(
faceLists[i].begin(),
faceLists[i].end(), false).first);
// 1: Right house (green)
volumeLists.push_back(...);
volumeLists.back().push_back(faces[6]);
volumeLists.back().push_back(faces[7]);
volumeLists.back().push_back(faces[8]);
volumeLists.back().push_back(faces[9]);
volumeLists.back().push_back(faces[10]);
volumeLists.back().push_back(faces[11]);
volumeLists.back().push_back(faces[12]);
volumeLists.back().push_back(faces[13]);
volumeLists.back().push_back(faces[14]);
// 11: Roof from left to right tip (red)
volumeLists.push_back(...);
volumeLists.back().push_back(faces[31]);
volumeLists.back().push_back(faces[32]);
volumeLists.back().push_back(faces[33]);
volumeLists.back().push_back(faces[34]);
volumeLists.back().push_back(faces[5]);
for (int i = 0; i < 12; ++i)
volumes.push_back(builder.get_cell<3>(
volumeLists[i].begin(),
volumeLists[i].end()).first);
LCC::Dart_handle house4d = builder.get_cell<4>(
volumes.begin(), volumes.end()).first;

BIN
work_report/figs/corine.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 155 KiB

BIN
work_report/figs/csg.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 50 KiB

BIN
work_report/figs/cube.jpg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 87 KiB

Some files were not shown because too many files have changed in this diff Show More