# Chapter 4 Spherical Geometries

“*There are too many false conclusions drawn and stupid measurements
made when geographic software, built for projected Cartesian
coordinates in a local setting, is applied at the global scale*”
(Chrisman 2012)

The previous chapter discussed geometries defined on the plane, \(R^2\). This chapter discusses what changes when we consider geometries not on the plane, but on the sphere (\(S^2\)).

Although we learned in chapter 2 that the shape of the Earth is usually approximated by an ellipsoid, none of the libraries shown in green in figure 1.6 provide access to a comprehensive set of functions that compute on an ellipsoid. Only the s2geometry (Dunnington, Pebesma, and Rubak 2021; Veach et al. 2020) library does provide it using a sphere rather than an ellipsoid. However, when compared to using a flat (projected) space we did in the previous chapter, a sphere is a much better approximation to an ellipsoid.

## 4.1 Straight lines

The basic premise of *simple features* of chapter 3
is that geometries are represented by sequences of points *connected
by straight lines*. On \(R^2\) (or any Cartesian space), this is
trivial, but on a sphere straight lines do not exist. The shortest
line connecting two points is an arc of the circle through both
points and the center of the sphere, also called a *great circle
segment*. A consequence is that “the” shortest distance line
connecting two points on opposing sides of the sphere does not exist,
as any great circle segment connecting them has equal length.

## 4.2 Ring direction

Any polygon on the sphere divides the sphere surface in two parts with finite area: the inside and the outside. Using the “counter clockwise rule” as was done for \(R^2\) will not work, because the direction interpretation depends on what is defined as inside. The convention here is to define the inside as the left (or right) side of the polygon boundary when traversing its points in sequence. Reversal of the node order then switches inside and outside.

## 4.3 Full polygon

In addition to empty polygons, one can define the full polygon on a sphere, which comprises its entire surface. This is useful, for instance for computing the oceans as the geometric difference between the full polygon and those of the land mass.

## 4.4 Bounding box, rectangle, and cap

Where in \(R^2\) one can easily define bounding boxes as the range
of the \(x\) and \(y\) coordinates, for ellipsoidal coordinates these
ranges are not of much use when geometries cross the antimeridian
(longitude +/- 180) or one of the poles. The assumption in \(R^2\)
that lower \(x\) values are Westwards of higher ones does not hold
when crossing the antimeridian. An alternative to delineating
an area on a sphere that is more natural is the *bounding cap*,
defined by its center coordinates and a radius. For Antarctica,
as depicted in figure 4.1 (a) and (c), the
bounding box formed by coordinate ranges is

```
library(sf)
suppressPackageStartupMessages(library(maps))
# maps:
par(mfrow = c(2,2))
par(mar = c(1,1.2,1,1))
m = st_as_sf(map(fill=TRUE, plot=FALSE))
a = m[m$ID == "Antarctica", ]
st_bbox(a)
```

```
# xmin ymin xmax ymax
# -180.0 -85.2 179.6 -60.5
```

which clearly does not contain the region (`ymin`

being -90 and `xmax`

180).
Two geometries that do contain the region are the bounding cap:

```
# lng lat angle
# 1 0 -90 29.5
```

and the bounding *rectangle*:

```
# lng_lo lat_lo lng_hi lat_hi
# 1 -180 -90 180 -60.5
```

For an area spanning the antimeridian, here the Fiji island country, the bounding box:

```
# xmin ymin xmax ymax
# -179.9 -21.7 180.2 -12.5
```

seems to span most of the Earth, as opposed to the bounding rectangle:

```
# lng_lo lat_lo lng_hi lat_hi
# 1 175 -21.7 -178 -12.5
```

where a value `lng_lo`

*larger* than `lng_hi`

indicates that the
bounding rectangle spans the antimeridian. This property could not
be inferred from the coordinate ranges.

## 4.5 Validity on the sphere

Many global datasets are given in ellipsoidal coordinates but are prepared in a way that they “work” when interpreted on the \(R^2\) space [-180,180] \(\times\) [-90,90]. This means that:

- geometries crossing the antimeridian (longitude +/- 180) are cut in halves, such that they no longer cross it (but nearly touch each other)
- geometries including a pole, like Antarctica, are cut at +/- 180 and make an excursion through -180,-90 and 180,-90 (both representing the Geographic South Pole)

Figure 4.1 shows two different representation of Antarctica, plotted with ellipsoidal coordinates taken as \(R^2\) (top) and in a Polar Stereographic projection (bottom), without (left) and with (right) an excursion through the Geographic South Pole. In the projections as plotted, polygons (b) and (c) are valid; polygon (a) is not valid as it self-intersects, polygon (d) is not valid because it traverses the same edge to the South Pole twice. On the sphere (\(S^2\)), polygon (a) is valid but (b) is not, for the same reason as (d) is not valid.

```
library(sf)
suppressPackageStartupMessages(library(maps))
# maps:
par(mfrow = c(2,2))
par(mar = c(1,1.2,1,1))
m = st_as_sf(map(fill=TRUE, plot=FALSE))
m = m[m$ID == "Antarctica", ]
plot(st_geometry(m), asp = 2)
title("a (not valid)")
# ne:
library(rnaturalearth)
ne = ne_countries(returnclass = "sf")
ne = ne[ne$region_un == "Antarctica", "region_un"]
plot(st_geometry(ne), asp = 2)
title("b (valid)")
# 3031
m %>%
st_geometry() %>%
st_transform(3031) %>%
plot()
title("c (valid)")
ne %>%
st_geometry() %>%
st_transform(3031) %>%
plot()
title("d (not valid)")
```

## 4.6 Exercises

For the following exercises, use R where possible or relevant.

- How does the GeoJSON format define
“straight” lines between ellipsoidal coordinates (section 3.1.1)?
Using this definition of straight, how would
`LINESTRING(0 85,180 85)`

look like in a polar projection? How could this geometry be modified to have it cross the North Pole? - For a typical polygon on \(S^2\), how can you find out ring direction?
- Are there advantages of using bounding caps over using bounding boxes? If so, list them.
- Why is, for small areas, the orthographic projection centered at the area a good approximation of the geometry as handled on \(S^2\)?
- For
`rnaturalearth::ne_countries(country = "Fiji", returnclass="sf")`

, check whether the geometry is valid on \(R^2\), on an orthographic projection centered on the country, and on \(S^2\). How can the geometry be made valid on`S^2`

? Plot the resulting geometry back on \(R^2\). Compare the centroid of the country, as computed on \(R^2\) and on \(S^2\), and the distance between the two.

### References

Chrisman, Nicholas. 2012. “A Deflationary Approach to Fundamental Principles in GIScience.” In *Francis Harvey (Ed.) Are There Fundamental Principles in Geographic Information Science?*, 42–64. CreateSpace, United States.

Dunnington, Dewey, Edzer Pebesma, and Ege Rubak. 2021. *S2: Spherical Geometry Operators Using the S2 Geometry Library*. https://CRAN.R-project.org/package=s2.

Veach, Eric, Jesse Rosenstock, Eric Engle, Robert Snedegar, Julien Basch, and Tom Manshreck. 2020. “S2 Geometry.” *Website*. https://s2geometry.io/.