Flat Earth Mathematics in the R Language

I am soon embarking on a trip around the world to speak at the 2018 World Water Congress in Tokyo, visit family and friends in the Netherlands and some tourism in San Francisco. One of the workshops I am involved with is about Fake News in the water industry.

One of the craziest fake news trends is the flat earth conspiracy. Some of you might ask whether my trip actually around the world or am I travelling across a flat disk?

This article discusses flat earth mathematics, and how to convert and visualise map projections in R. This article uses the code I published earlier in articles about creating flight maps and Pacific island hopping. You can view the code below or download the full version from GitHub.

The Flat Earth

YouTube contains thousands of videos from people that claim the earth is flat and that science as we know it is a “globalist conspiracy”. While their claims are hilarious, people that believe our planet is flat are often passionate and even conduct simple experiments and engage in naive flat earth mathematics.

Adherents to the idea that the world is flat often propose Gleason’s Map (shown above) as their version of the correct representation of our planet. This map is not a map of the flat earth but a polar Azimuthal equidistant projection of the globe. Even though the map itself states that it is a projection of the world, flat earth believers nevertheless see it as a literal representation.

Projecting the spherical earth on a flat surface is a complex task which will always require compromise as it is impossible to truthfully draw the surface of the globe on a piece of paper. The video below provides a practical introduction to map projections that show how maps are stretched to display them on a flat surface.

We can recreate Gleason’s map with ggplot, which incorporates the mapproj package to show maps in various projections. The Azimuthal equidistant method is popularised in the flag of the United Nations. Antarctica is in this projection displayed as a ring around the world. Flat earth evangelists believe that the South Pole a ring of ice that prevents us from proceeding beyond the disc. The edge of this disc cannot be shown because this method projects the south pole at an infinite distance from the centre of the map.

library(tidyverse)
world <- map_data("world")
worldmap <- ggplot(world) +
    geom_path(aes(x = long, y = lat, group = group)) +
    theme(axis.text = element_blank(),
          axis.ticks = element_blank()) +
    labs(x = "", y = "")
worldmap + 
    coord_map("azequidistant", orientation = c(90, 0, 270))

The coord_map function allows you to change projections. The orientation argument defines the centre point of the projection. The most common location is the North Pole (90, 0). The third value gives the clockwise rotation in degrees. Gleason’s map uses a rotation of 270 degrees.

Polar azimuthal equidistant projection with ggmap
Polar azimuthal equidistant projection with ggmap.

My Round-the-World Itinerary

My trip will take me from Australia to Japan, from where I go to the Netherlands. The last two legs will bring me back to Australia via San Francisco.

The itinerary is stored in a data frame, and the ggmap package geocodes the longitude and latitude of each of the locations on my trip. As the earth is a sphere, an intermediate point needs to be added for trips that pass the dateline, as I explained in my article about flight maps.

The itinerary is visualised using the same method as above but centring on Antarctica. The geosphere package helps to estimate the total travel distance, which is approximately 38,995 km, slightly less than a trip around the equator. This distance is the great circle distance, which is the shortest distance between two points on a sphere, adjusted for a spheroid (a flattened sphere).

The flight paths on this map are curved because of the inevitable distortions when projecting a sphere on a flat surface.

## Define itinerary
library(ggmap)
airports <- tibble(city = c("Melbourne", "Tokyo", "Amsterdam", "San Francisco", "Melbourne"))

## Find coordinates
## Where you receive an OVER_QUERY_LIMIT error, repeat the code
itinerary <- geocode(airports$city) %>%
    mutate(location = airports$city)

## Split travel past dateline
dl <- which(diff(itinerary$lon) > 180)
dr <- ifelse(itinerary$lon[dl] < 0, -180, 180)
dateline <- tibble(lon = c(dr, -dr),
                   lat = rep(mean(itinerary$lat[dl:(dl + 1)]), 2),
                   location = "dateline")
itinerary <- rbind(itinerary[1:dl, ], dateline, itinerary[(dl + 1):nrow(itinerary), ]) itinerary ## Visualise worldmap + geom_point(data = itinerary, aes(lon, lat), colour = "red", size = 4) + geom_path(data = itinerary, aes(lon, lat), colour = "red", size = 1) + coord_map("azequidistant", orientation = c(-90, 0, 270)) ## Great Circle Distance library(geosphere) sapply(1:(nrow(itinerary) - 1), function(l) distVincentyEllipsoid(itinerary[l, 1:2], itinerary[(l + 1), 1:2]) / 1000) %>%
    sum()
Rount the world trip
Round the world trip in polar Azimuthal equidistant projection.

Flat Earth Mathematics

If the Gleason map was an actual map of the flat earth, then the flight paths on the map would show as straight lines.

The mapproj package contains the mapproject function that calculates the projected coordinates based on longitude and latitude. The output of this function is a grid with limits from -\pi to \pi. The first part of the code converts the longitude and latitude from the world data frame to projected coordinates.

A line from the lon/lat 0,0 to the north pole has a projected distance of \pi/2, which in the spherical world is \pi / 2 * 6378.137 = 10018.75 km. In other words, we need to multiply the Euclidean distances with the radius of the Earth to derive the Gleason map coordinates.

This last code snippet converts the world map to flat earth coordinates, calculates the Euclidean distance between the points on the itinerary and multiplies this with the Earth’s diameter.

This last code snippet shows why the Gleason map is not a map of a flat earth. On this map, the shortest distance to fly from Sydney to Santiago de Chili is about 25,000 km, more than twice the actual distance. The actual travel time is about 14 hours, which implies that passenger jets break the sound barrier. This mainly is a problem for journeys along lines of lattitude in the Southern hemisphere as the distortion in this projection increases with the distance from the centre of the map.

This map looks like the first one, but the coordinate system is now Euclidean instead of longitudes and lattitudes, as indicated by the square grid.

Flat earth mathematics: Sydney to Santiago de Cili
Sydney to Santiago de Chili on a flat earth map.

This article proves that the Gleason map is not a representation of a flat earth because aeroplanes would have to break the sound barrier to fly these distances in the time it takes to travel. Whichever way you project the globe on a flat map will lead to inevitable distortions. The Gleason map itself mentions that it is a projection.

However, these facts will not dissuade people from believing in a flat earth, I am after all an engineer and thus part of the globalist science conspiracy.

You can view the complete code on my GitHub repository.

## Flat Eart Mathematics
library(mapproj)
flatearth <- mapproject(world$long, world$lat,
           "azequidistant", orientation = c(90, 0, 270))
world <- mutate(world, x = flatearth$x, y = flatearth$y)

## Australia - South America
airports <- tibble(city = c("Sydney", "Santiago de Chili"))
itinerary <- geocode(airports$city) %>%
    mutate(location = airports$city)
itinerary
coords <- mapproject(itinerary$lon, itinerary$lat, "azequidistant",
                     orientation = c(90, 0, 270))
coords <- tibble(x = coords$x, y = coords$y)
sum(sqrt(diff(coords$x)^2 + diff(coords$y)^2) * 6378.137)

ggplot(world) +
    geom_path(aes(x, y, group = group)) +
    theme(axis.text = element_blank(),
          axis.ticks = element_blank()) +
    labs(x = "", y = "") + 
    geom_point(data = coords, aes(x, y), colour = "red", size = 4) +
    geom_path(data = coords, aes(x, y), colour = "red", size = 1)

Leave a Reply

This site uses Akismet to reduce spam. Learn how your comment data is processed.