library(oce)
# 1. Today's value at Halifax NS
magneticField(-(63 + 36 / 60), 44 + 39 / 60, Sys.Date())
# 2. World map of declination in year 2025.
# \donttest{
data(coastlineWorld)
par(mar = rep(0.5, 4)) # no axes on whole-world projection
mapPlot(coastlineWorld, projection = "+proj=robin", col = "lightgray")
# Construct matrix holding declination
lon <- seq(-180, 180)
lat <- seq(-90, 90)
dec2025 <- function(lon, lat) {
magneticField(lon, lat, 2025)$declination
}
dec <- outer(lon, lat, dec2025) # hint: outer() is very handy!
# Contour, unlabelled for small increments, labeled for
# larger increments.
mapContour(lon, lat, dec,
col = "blue", levels = seq(-180, -5, 5),
lty = 3, drawlabels = FALSE
)
mapContour(lon, lat, dec, col = "blue", levels = seq(-180, -20, 20))
mapContour(lon, lat, dec,
col = "red", levels = seq(5, 180, 5),
lty = 3, drawlabels = FALSE
)
mapContour(lon, lat, dec, col = "red", levels = seq(20, 180, 20))
mapContour(lon, lat, dec, levels = 180, col = "black", lwd = 2, drawlabels = FALSE)
mapContour(lon, lat, dec, levels = 0, col = "black", lwd = 2)
# }
# 3. Declination differences between versions 13 and 14
lon <- seq(-180, 180)
lat <- seq(-90, 90)
decDiff <- function(lon, lat) {
old <- magneticField(lon, lat, 2025, version = 13)$declination
new <- magneticField(lon, lat, 2025, version = 14)$declination
new - old
}
decDiff <- outer(lon, lat, decDiff)
t.test(decDiff)
Run the code above in your browser using DataLab