Anisotropic shortest paths across the Strait of Gibraltar

Javier Fernández-López and Klaus Schliep

October 08, 2026

Overview

This vignette shows the core rWind workflow: obtain wind data, represent them as modern spatial rasters, translate wind speed and direction into anisotropic movement conductance, and calculate directional least-cost paths.

The example crosses the Strait of Gibraltar. Because wind-assisted movement is direction dependent, the least-cost path from Spain to Morocco does not necessarily match the path in the opposite direction.

Obtain wind data

wind.dl() automatically selects the appropriate data source. Current dates are requested from PacIOOS, whereas older dates are retrieved from the NOAA/NCEI GFS archive. GFS data are available at three-hour UTC intervals.

The download call used for this example is:

w <- wind.dl(
  2015, 2, 12, 12,
  lon1 = -7, lon2 = -4,
  lat1 = 34.5, lat2 = 37.5
)

To make the vignette reproducible without a network connection, it uses the result saved with the package.

w <- readRDS("w.rds")
head(w)
#>                   time  lat  lon ugrd10m vgrd10m      dir     speed
#> 43 2015-02-12 12:00:00 37.5 -7.0    0.74   -0.30 112.0679 0.7984986
#> 44 2015-02-12 12:00:00 37.5 -6.5   -0.90   -0.27 253.3008 0.9396276
#> 45 2015-02-12 12:00:00 37.5 -6.0   -2.20   -0.97 246.2069 2.4043503
#> 46 2015-02-12 12:00:00 37.5 -5.5   -2.29   -0.72 252.5463 2.4005207
#> 47 2015-02-12 12:00:00 37.5 -5.0   -1.42    0.02 270.8069 1.4201408
#> 48 2015-02-12 12:00:00 37.5 -4.5   -0.26    0.18 304.6952 0.3162278

Create a SpatRaster

wind2raster() converts the tabular result into a two-layer terra::SpatRaster (Hijmans et al. 2026). The layers are named direction and speed, and the coordinate reference system is WGS 84 (EPSG:4326).

wind_layer <- wind2raster(w)
wind_layer
#> class       : SpatRaster
#> size        : 7, 7, 2  (nrow, ncol, nlyr)
#> resolution  : 0.5, 0.5  (x, y)
#> extent      : -7.25, -3.75, 34.25, 37.75  (xmin, xmax, ymin, ymax)
#> coord. ref. : lon/lat WGS 84 (EPSG:4326)
#> source(s)   : memory
#> names       :  direction,    speed
#> min values  :  37.568591, 0.316228
#> max values  : 304.695153, 6.931695

Individual layers use standard terra syntax:

wind_layer[["speed"]]
#> class       : SpatRaster
#> size        : 7, 7, 1  (nrow, ncol, nlyr)
#> resolution  : 0.5, 0.5  (x, y)
#> extent      : -7.25, -3.75, 34.25, 37.75  (xmin, xmax, ymin, ymax)
#> coord. ref. : lon/lat WGS 84 (EPSG:4326)
#> source(s)   : memory
#> name        :    speed
#> min value   : 0.316228
#> max value   : 6.931695

Build directional conductance

For passive movement, flow.dispersion() applies a restricted version of the formulation of Felicísimo et al. (2008). Deviations of 90 degrees or more from the flow are treated as barriers. The resulting TransitionLayer stores directional conductance and can be used by gdistance (van Etten 2017). Internally, rWind converts only the direction layer to the legacy raster representation required by gdistance; the public spatial representation remains a SpatRaster.

conductance <- flow.dispersion(wind_layer, type = "passive")
conductance
#> class      : TransitionLayer 
#> dimensions : 7, 7, 49  (nrow, ncol, ncell)
#> resolution : 0.5, 0.5  (x, y)
#> extent     : -7.25, -3.75, 34.25, 37.75  (xmin, xmax, ymin, ymax)
#> crs        : +proj=longlat +datum=WGS84 +no_defs 
#> values      : conductance 
#> matrix class: dgCMatrix

Every transition cost is the mean of the local costs at the origin and destination cells. Diagonal transitions are multiplied by sqrt(2). This corrects diagonal length on an abstract square grid, but not the changing physical distances of longitude/latitude cells.

Active movement

The optional active model is deliberately general rather than mechanistic. It uses the signed component of flow parallel to movement. Let z = speed * cos(angle) / speed.scale: supportive flow (z >= 0) has cost 1 / (1 + z), while opposing flow has cost 1 - z. Perpendicular flow is neutral and movement against the flow remains possible.

active_conductance <- flow.dispersion(
  wind_layer,
  type = "active",
  speed.scale = median(w$speed[w$speed > 0], na.rm = TRUE)
)
active_conductance
#> class      : TransitionLayer 
#> dimensions : 7, 7, 49  (nrow, ncol, ncell)
#> resolution : 0.5, 0.5  (x, y)
#> extent     : -7.25, -3.75, 34.25, 37.75  (xmin, xmax, ymin, ymax)
#> crs        : +proj=longlat +datum=WGS84 +no_defs 
#> values      : conductance 
#> matrix class: dgCMatrix

If speed.scale is omitted, rWind uses the median positive speed of the raster. Supply one fixed value when comparing dates or regions. Active costs are relative connectivity indices: they are not estimates of travel time, energy, or an organism’s intrinsic movement speed.

Calculate paths in both directions

spain <- c(-5.5, 37)
morocco <- c(-5.5, 35)

spain_to_morocco <- shortestPath(
  conductance, spain, morocco,
  output = "SpatialLines"
)
morocco_to_spain <- shortestPath(
  conductance, morocco, spain,
  output = "SpatialLines"
)

The following map uses terra for the raster and path representation. Wind arrows are drawn from the downloaded U and V components; their lengths are scaled for display.

plot(
  wind_layer[["speed"]],
  col = hcl.colors(20, "YlOrRd", rev = TRUE),
  main = "Wind-assisted paths across the Strait of Gibraltar",
  xlab = "Longitude", ylab = "Latitude"
)

lines(vect(spain_to_morocco), col = "#C62828", lwd = 3, lty = 2)
lines(vect(morocco_to_spain), col = "#1565C0", lwd = 3, lty = 2)
points(spain[1], spain[2], pch = 19, cex = 1.3, col = "#C62828")
points(morocco[1], morocco[2], pch = 19, cex = 1.3, col = "#1565C0")

arrow_scale <- 0.12 / max(w$speed, na.rm = TRUE)
arrows(
  w$lon, w$lat,
  w$lon + w$ugrd10m * arrow_scale,
  w$lat + w$vgrd10m * arrow_scale,
  length = 0.035, col = "grey20"
)

text(spain[1], spain[2] + 0.18, "Spain", col = "#C62828", font = 2)
text(morocco[1], morocco[2] - 0.18, "Morocco", col = "#1565C0", font = 2)
legend(
  "topleft",
  legend = c("Spain to Morocco", "Morocco to Spain"),
  col = c("#C62828", "#1565C0"),
  lwd = 3, lty = 2, cex = 0.8, bg = "white"
)

The two paths differ because the conductance matrix is asymmetric: travelling with the wind has a different cost from travelling against it. For larger areas or long time series, both raster size and the number of possible cell transitions increase rapidly, so start with the smallest spatial extent and resolution needed for the biological question.

References

Felicísimo, Ángel M., Jesús Muñoz, and Jacob González-Solís. 2008. “Ocean Surface Winds Drive Dynamics of Transoceanic Aerial Movements.” PLoS ONE 3 (8): e2928. https://doi.org/10.1371/journal.pone.0002928.
Hijmans, Robert J., Andrew Brown, A. Márcia Barbosa, Emanuele Cordano, and Krzysztof Dyba. 2026. Terra: Spatial Data Analysis. https://doi.org/10.32614/CRAN.package.terra.
van Etten, Jacob. 2017. “R Package Gdistance: Distances and Routes on Geographical Grids.” Journal of Statistical Software 76 (13): 1–21. https://doi.org/10.18637/jss.v076.i13.