Supersonic Routes provides a
quick end-to-end introduction to using the
and is the place to start. This vignette provides advice on more
advanced use, explaining details that the introduction skates over.
Much of this vignette is optional or for occasional use, but the advice on saving and reading the cache is likely to be essential for a speedy workflow.
#the libraries needed for the vignette are library(himach) library(dplyr, quietly = TRUE, warn.conflicts = FALSE) library(ggplot2) library(sf) # and we'll load a full set of test data <- hm_get_test("coast") NZ_coast <- hm_get_test("buffer") NZ_buffer30 <- hm_get_test("nofly") NZ_Buller_buffer40 <- hm_get_test("grid") NZ_grid <- hm_get_test("route")NZ_routes
himach uses caching to speed things up. Legs are cached
route_cache and arrival-departure links to airports are
star_cache (STAR is short for standard arrival
route, its counterpart being SID the standard instrument departure).
You will want to save and load the cache (meaning the combination of
star_cache) as part of your
standard workflow. Quite where you save it is up to you, but a set of
routes is sensitive to (a) the route grid on which it is calculated (b)
the list of aircraft used, and their performance. The saving function
hm_save_cache forces you to refer to these two datasets,
and uses metadata from them in the file name for the cache.
If you change either of these, then you can use
hm_clean_cache() to empty the cache. You will also note
that if you run
find_routes and the map has changed, or
findToCToD and the map or aircraft have changed, then the
cache will be cleared automatically.
For the vignette, we save to a temporary directory. You really don’t want to do this in practice ;-)
hm_clean_cache() #start without cache # need to load some of the built-in data for this example <- make_aircraft(warn = FALSE) aircraft <- make_airports(crs = crs_Pacific) airports #> Using default airport data: airportr::airport. options("quiet"= 2) # for a little reporting # how long does it take with an empty cache? system.time( <- find_route(aircraft[1, ], routes make_AP2("NZAA", "NZDN", airports), fat_map = NZ_buffer30, route_grid = NZ_grid, ap_loc = airports) )#> Route:-NZAA<>NZDN---- #> Leg: NZAA<>NZDN Aircraft: SST M2.2 #> Cut envelope from lattice: 0.3 #> Calculated phase changes #> Done recursion #> Checking Shortcuts #> user system elapsed #> 1.294 0.033 1.332 # test saving of cache to a disposable directory <- tempdir() tmp_dir # for convenience, hm_save_cache gives the full name, including path <- hm_save_cache("test_v", NZ_grid, aircraft, path = tmp_dir) full_filename #empty cache - just to demonstrate the re-loading # this isn't part of your normal workflow! hm_clean_cache() # but normally a session will begin with loading a cache like this hm_load_cache(full_filename) # how long does it take with a cache? system.time( <- find_route(aircraft[1, ], routes make_AP2("NZAA", "NZDN", airports), fat_map = NZ_buffer30, route_grid = NZ_grid, ap_loc = airports) )#> Route:-NZAA<>NZDN---- #> user system elapsed #> 0.037 0.001 0.037 # if you want to see a map # map_routes(NZ_coast, routes, crs_Pacific, fat_map = NZ_buffer30, simplify_km = 2)
The cache just works invisibly in the background - you will notice it
speeds up finding of routes no end: in that example, from 1.5s (user) to
0.04s (user) on my machine. In particular, it helps with refuelling,
route_cache quickly remembers the routes from
major hub airports to the main refuelling points, so they don’t need to
be calculated again.
Incidentally, if you add a new refuelling point, then the cache
remains valid because only legs are cached, not
routes. With a new refuelling point,
will check both old legs and new (to the new refuelling points), gaining
where the legs are cached, before selecting the best combination of legs
to make the route.
There may be times when you have built a large cache, of several
thousand routes for a handful of aircraft. If you want to make a change
to the performance of one aircraft, say, you either need to
hm_clean_cache() and start everything again, or do a bit of
housekeeping. The same is true if an airport location needs correcting:
ideally we delete all references to just that airport.
If you’re not comfortable housekeeping, then just let
himach run overnight to redo the whole thing.
Housekeeping involves loading the cache manually, with
load. It contains two environments
star_cache (segments to and from airports). The items
in the cache are named with strings which capture a number of run
parameters, including the aircraft ID. So we identify items for
rm them. Then save the cache manually with
# for this example, add a second route to the cache <- find_route(aircraft[2, ], routes make_AP2("NZAA", "NZDN", airports), fat_map = NZ_buffer30, route_grid = NZ_grid, ap_loc = airports) #> Route:-NZAA<>NZDN---- #> Leg: NZAA<>NZDN Aircraft: SST M1.6 8.8 #> Cut envelope from lattice: 0.2 #> Calculated phase changes #> Done recursion #> Checking Shortcuts # save the cache, which has NZAA-NZDN for 2 aircraft now hm_save_cache("test_v", NZ_grid, aircraft, path = tmp_dir) #now do housekeeping load(full_filename) # filename from the previous chunk ls(route_cache) # show the contents, just for information #>  "M16_9-NZAA-NZDN--TRUE-TRUE-NA-TRUE-100-12" #>  "M22-NZAA-NZDN--TRUE-TRUE-NA-TRUE-100-12" # we want to delete instances of aircraft with ID that includes 'M22' <- ls(route_cache, pattern="M22") |> as.list() z length(route_cache) # before deletion #>  2 do.call(rm, z, envir = route_cache) # delete the M22 items length(route_cache) #after deletion, 1 less #>  1 # then repeat for star_cache <- ls(star_cache, pattern="M22") |> as.list() z length(star_cache) #>  4 do.call(rm, z, envir = star_cache) length(star_cache) #>  2 # then save the result (you might want to change the filename, or backup the old cache beforehand) save("route_cache", "star_cache", file = full_filename)
It is not unusual for parts of the airspace to be closed, or be
considered unsafe for flying.
himach allows regions to be
marked as ‘avoid’. They will not feature in the grid, so routes will
avoid them, with one exception: an arrival or departure airport can be
inside a no-fly zone, as long as at least one connection point to the
grid is outside. So they might more precisely be called ‘no-overfly’
A no-fly zone is prepared in the same way as a map of land. If specific countries are to be avoided, this is where having a country name in the geographic data comes in handy.
One essential item is the
avoid attribute of the no-fly
zone. This is used to distinguish sets of legs with different, or no,
no-fly zone. Set
attr(your_avoid_map, "avoid") <- "your summary of that avoid map"
which will (a) remind you what was used (b) tell
recalculate all legs that have not already been calculated with that
avoid. If you were to add an avoid area for North
Korean airspace, say, then in reality North Atlantic routes are not
affected, but currently
himach plays safe and assumes that
In this example, no offence is intended to the citizens of Buller District of New Zealand; it is a convenient example for showing how routes are forced to change when airspace is unavailable.
# using your own shp file # NZ_Buller <- sf::read_sf("...../territorial-authority-2020-clipped-generalised.shp") %>% # filter(TA2020_V_1 == "Buller District") # NZ_Buller_u <- sf::st_union(sf::st_simplify(NZ_Buller, dTolerance = 1000)) # NZ_Buller_buffer50 <- sf::st_union(sf::st_buffer(NZ_Buller_u, 50 * 1000)) # attr(NZ_Buller_buffer50, "avoid") <- "Buller+50km" # the quicker version, using a built-in no fly zone # this uses data as in the previous code chunk <- make_aircraft(warn = FALSE) aircraft <- make_airports(crs = crs_Pacific) airports #> Using default airport data: airportr::airport. # run the same route, but with the avoid region options("quiet"= 2) #just the progress bar <- aircraft[c(1, 4), ]$id ac <- find_routes(ac, routes data.frame(ADEP = "NZAA", ADES = "NZDN"), aircraft, airports,fat_map = NZ_buffer30, route_grid = NZ_grid, cf_subsonic = aircraft[3,], avoid = NZ_Buller_buffer40) #> Route:-NZAA<>NZDN---- #> Leg: NZAA<>NZDN Aircraft: SST M2.2 #> Cut envelope from lattice: 0.3 #> Calculated phase changes #> Done recursion #> Checking Shortcuts #> Adding subsonic, without range bounds. #> Leg: NZAA<>NZDN Aircraft: 777-300ER #> Calculated phase changes #> Done recursion #> Checking Shortcuts #> #> Route:-NZAA<>NZDN---- #> Too far for one leg. #> Adding subsonic, without range bounds. #> #this shows versions of the legs with and without no-fly # ls(route_cache, pattern = "NZCH", envir = .hm_cache) # create route summary <- summarise_routes(routes, airports) rtes # draw a basic map map_routes(NZ_coast, routes, crs_Pacific, fat_map = NZ_buffer30, avoid_map = NZ_Buller_buffer40, simplify_km = 2)
map_routes(NZ_coast, routes, show_route = "aircraft", crs = crs_Pacific, fat_map = NZ_buffer30, avoid_map = NZ_Buller_buffer40, simplify_km = 2)
After a call to
find_routes, the output can have
NA entries in some columns for some routes. There are two
reasons for this:
refuel = xxx, then you will find other entries for the same
routeID(eg “EGLL<>KSFO”) but with different
fullRouteID(eg “EGLL<>PANC<>KSFO”) showing a good route including refuelling.
So these appear when the specified route is not possible.
Above 60 or 70 (North or South), the approximations used by the
st_buffer function begin to show signs of exceeding their
limits. In particular, if you’re adding a 50km coastal buffer, for
example, there are separations between Canadian islands which are just
under 100km. Borden and Ellef Ringnes are examples. A buffer generated
st_buffer shows the strait between them as open water,
where it should be closed.
This can lead to over-optimistic routings: supersonic where they should not be.
The solution is to use the links from
sf to the
s2 package which come in more recent versions of the
sf package. This does require you to use quite a high value
max_cells parameter of
<- s2::s2_data_countries(c("Greenland", "Canada", "Iceland")) gr <- s2::s2_buffer_cells(gr, distance = 50000, max_cells = 20000) %>% gr_buffer_s2 st_as_sfc() <- ggplot(st_transform(gr_buffer_s2, crs_Atlantic)) + geom_sf(fill = "grey40") + m_s2 geom_sf(data = st_transform(st_as_sfc(gr), crs_Atlantic)) sf_use_s2(FALSE) # to be sure <- gr %>% gr_transf st_as_sfc() %>% st_transform(crs_Atlantic) <- gr_transf %>% gr_t_buffer st_buffer(dist = 50000) <- ggplot(gr_t_buffer) + geom_sf(fill = "grey40") + geom_sf(data = gr_transf) m_old ::plot_grid(m_old, m_s2, labels = c("bad", "good"), cowplotncol = 1)
In fact, the problem of finding too many apparently over-ocean routes is broader than this. The other main contributor to this is missing islands from the map. See the comments in the first vignette.
An example of this is in the same place. Some maps omit small islands (well, larger ones like Killniq down to tiny ones like Goodwin Island) at the mouth of the Hudson Strait. This affects the apparent width of the opening. Given the islands, and a 50km buffer, the Strait is not open as the next example shows.
This uses a non-CRAN, but public package of hi-resolution maps,
rnaturalearthhires. If you don’t want to load this package,
just note the results shown in the figure.
::sf_use_s2(TRUE) sf<- sf::st_as_sf(rnaturalearthhires::countries10) %>% hires filter(NAME %in% c("Greenland", "Canada", "Iceland")) <- s2::s2_buffer_cells(hires, distance = 50000, max_cells = 20000) %>% hires_buffer_s2 st_as_sfc() <- ggplot(st_transform(hires_buffer_s2, crs_Atlantic)) + m_hires geom_sf(fill = "grey40") + geom_sf(data = st_transform(hires, crs_Atlantic)) ::plot_grid(m_s2, m_hires, labels = c("good", "better"), cowplotncol = 1)
There are a number of place in the vignettes, eg making an airport
dataset, where we have shown the use of a parameter to specify a
coordinate reference system.
himach has recently
transitioned to using spherical geometry directly using the
s2 package, both directly and through the
s2 was available in
was a constant need to align the coordinate reference systems of objects
before combining them.
Now, in theory, all geometrical operations use spherical geometry, so
a coordinate reference system should only be needed when you plot a map.
At that point, the coordinate reference system is saying how to move
from spherical coordinates to a flat projection. Four basic projections
crs_South which you can use in
map_routes to get the right map for your particular set of
routes. You can create others as shown in the vignette.
We will remove remaining references to coordinate reference systems
during route creation in later versions of