Plotting GPS Coordinates on a Polar Projection Using rgdal

I would definitely defer to Michael and the other true "Geo" folks on this
list**, but if you're ok with base graphics (and a much less technical
interface) you could accomplish what you want using the `oce` package:

library(oce)
data(coastlineWorld)
circ <- list(longitude=seq(-180, 180), latitude=-rep(60, 361)) # fake
circumnavigation
## use mapPlot to make a polar stereographic plot of antarctica
mapPlot(coastlineWorld, projection="+proj=stere +lat_0=-90",
latitudelim=c(-90, -45), longitudelim = c(-180, 180), col='grey')
## now add the circumnavigation with mapPoints
mapPoints(circ)

Hope that helps!
Clark

**as you can see mapPlot is using proj/gdal interface under the hood, but
being oceanographers and only partially understanding the details of the
transforms, datums, projections, etc, we tried to write an interface that
will allow us to make pretty maps without getting too far into the details.


Message: 4

Date: Mon, 23 Oct 2017 22:53:22 +0000
From: "Mortimer, Linsey Anne" <l.mortimer.14 at aberdeen.ac.uk>
To: Michael Sumner <mdsumner at gmail.com>
Cc: "r-sig-geo at r-project.org" <r-sig-geo at r-project.org>

Subject: Re: [R-sig-Geo] Plotting GPS Coordinates on a Polar
        Projection Using rgdal

Message-ID:

        <
VI1PR0101MB23173329C1E17653F09ECC8086460 at VI1PR0101MB2317.eurprd01.prod.exchangelabs.com

Content-Type: text/plain; charset="Windows-1252"

This was a big help, thank you. The coordinates are plotted and look good
except some of my GPS coordinates are outside of the plot - I should have
check this beforehand. Is there a way to produce a plot using spbabel
similar to this:?

library(ggplot2)
world <- map_data("world")
worldmap <- ggplot(world, aes(x=long, y=lat, group=group)) +

? geom_path() +

? scale_y_continuous(breaks=(-2:2) * 30) +
? scale_x_continuous(breaks=(-4:4) * 45)

worldmap + coord_map("ortho", orientation=c(-90, 0, 0))


or would it be better to simply use ggplot for this?


Many thanks,

Linsey


From: Michael Sumner <mdsumner at gmail.com>
Sent: 23 October 2017 21:35
To: Mortimer, Linsey Anne
Cc: r-sig-geo at r-project.org
Subject: Re: [R-sig-Geo] Plotting GPS Coordinates on a Polar Projection
Using rgdal

?

On Tue, 24 Oct 2017 at 03:42 Mortimer, Linsey Anne <
l.mortimer.14 at aberdeen.ac.uk> wrote:
  Hello,

I am new to R and attempting to plot a map of Antarctica with the GPS
coordinates from a circumpolar navigation of the Southern Ocean. I have
searched the archives and rgdal-related manuals/blogs/papers to find
variations of the code below. I have two issues:  firstly, removing the
vertical line at -90 degrees on the plot of Antarctica. Secondly, adding
the coordinates to this plot using points(). I have tried transforming the
coordinates of the points to fit the projection but every time the output
is a single  point plotted at the end of the vertical line on the
Antarctica map, a separate plot of points in a world map projection
(stretched) rather than polar projection, or simply an error about coercing
an S4 class to a vector. Any advice on how to correct the code  or a point
in the direction of a guide to this kind of mapping and spatial analysis
would be of great help.








You can remove the dummy seam at the south pole by decomposing to raw
coordinates and reconstructing with spbabel:?

data("wrld_simpl", package = "maptools")
crs <- sp::proj4string(wrld_simpl)
antarctica0 <- wrld_simpl[wrld_simpl$NAME == "Antarctica", ]
library(spbabel)
library(dplyr)

ant <- spbabel::sptable(antarctica0) %>%?

? dplyr::filter(y_ > -90) %>%?
? spbabel::sp(crs = crs)

pr <- "+proj=laea +lat_0=-90 +ellps=WGS84 +datum=WGS84 +no_defs
+towgs84=0,0,0"



antarctica <- sp::spTransform(ant, sp::CRS(pr))




That makes sense when plotted in polar form, but will wrap improperly in
native longitude/latitude form.?


?  An idea of what my dataset looks like (10774579 obs. of 9 variables):

head(GPS)

?..GPS_fix? ? GPS_date? GPS_time GPS_milliseconds? Latitude Longitude
GPS_status

1? ? ? ? ? 0? 2016-12-24? 14:26:03? ? ? ? ? ? ? 878 -44.76846? 35.36859?
? ? ? ? 2

2? ? ? ? ? 1? 2016-12-24? 14:26:04? ? ? ? ? ? ? 309 -44.76846? 35.36859?
? ? ? ? 2

3? ? ? ? ? 2? 2016-12-24? 14:26:04? ? ? ? ? ? ? 878 -44.76849? 35.36866?
? ? ? ? 1

4? ? ? ? ? 3? 2016-12-24? 14:26:05? ? ? ? ? ? ? 309 -44.76849? 35.36866?
? ? ? ? 2

5? ? ? ? ? 4? 2016-12-24? 14:26:05? ? ? ? ? ? ? 878 -44.76852? 35.36873?
? ? ? ? 1

6? ? ? ? ? 5? 2016-12-24? 14:26:06? ? ? ? ? ? ? 305 -44.76852? 35.36873?
? ? ? ? 2

GPS_filename? ? Distance

1? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw 0.000000000

2? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw 0.000000000

3? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw? ? ? ? ? NA

4? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw 0.003509715

5? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw? ? ? ? ? NA

6? D:\\Leg1\\0 CT to Marion\\12500\\ACE-D20161224-T142556.raw 0.007045702




The code I?m using:

install.packages(c("sp", "raster", ?crs?))

library(maptools)

gpclibPermit()

data(wrld_simpl)

antarctica <- wrld_simpl[wrld_simpl$NAME == "Antarctica", ]

library(rgdal)

pr <- "+proj=laea +lat_0=-90 +ellps=WGS84 +datum=WGS84 +no_defs

+towgs84=0,0,0" #define projection

antarctica.laea <- spTransform(antarctica, CRS(pr))

plot(antarctica.laea)

## Now I attempt to make the GPS table into a spatial object and plot

the coordinates

attach(GPS)

coordinates(GPS) <- c("Longitude","Latitude")

points(Longitude, Latitude, pch=19, col="red", cex=0.5) #This only

plots a singular point at the end of vertical line on the map






Here you are mixing old attach() idiom with more modern sp idiom to
assign the coordinates, it's probably best to not use attach at all.?

coordinates(GPS) <- c("Longitude","Latitude")



## despite the name of the coordinates, we don't?
## otherwise know the coordinate system in use,?
## so assign is here:?

proj4string(GPS) <- sp::CRS("+init=epsg:4326")


## now we can transform to the projection in use
gps.ant <- sp::spTransform(GPS, sp::CRS(pr))

?

You might find the sf package easier to use now, though there are dozens
of tracking-specific packages with some other improvements. I find
generally that ggplot2 provides the best overall approach to plotting, but
it requires specific handling of the  coordinate systems for each data set,
as done here.?


The trip package has some methods to generate SpatialLinesDataFrame (and
other high level classes), but to do this by removing track-information to
get to GIS-friendly form. (It's also really painful to use and is in
desperate need of a user-friendly update.)  There's no real stand out good
approach for this general problem, so feel free to ask more specifics.?
There are a huge number of disparate tracking-packages listed on CRAN so
you might explore them for your needs:?

https://cran.r-project.org/web/views/SpatioTemporal.html

(sere under "Moving objects, trajectories")


Cheers, Mike.?

  > ## Also tried to transform the coordinates but receive an error

proj4string(GPS) <- CRS("+proj=laea +lat_0=-90 +ellps=WGS84 +datum=WGS84

+no_defs +towgs84=0,0,0")

gps.ant <- spTransform(GPS, CRS(antarctica.laea))

Error in CRS(antarctica.laea) :

? no method for coercing this S4 class to a vector

In addition: Warning message:

In is.na(projargs) :  is.na() applied to non-(list or vector) of type
'S4'

sessionInfo()

R version 3.4.2 (2017-09-28)

Platform: x86_64-w64-mingw32/x64 (64-bit)

Running under: Windows >= 8 x64 (build 9200)



Matrix products: default



locale:

[1] LC_COLLATE=English_United Kingdom.1252? LC_CTYPE=English_United
Kingdom.1252

[3] LC_MONETARY=English_United Kingdom.1252 LC_NUMERIC=C

[5] LC_TIME=English_United Kingdom.1252



attached base packages:

[1] stats? ? ?graphics? grDevices utils? ? ?datasets? methods? ?base



other attached packages:

[1] maptools_0.9-2 rgdal_1.2-13? ?sp_1.2-5

loaded via a namespace (and not attached):

[1] compiler_3.4.2? tools_3.4.2? ? ?foreign_0.8-69? grid_3.4.2? ? ?
lattice_0.20-35



Many thanks,



Linsey Mortimer



? ? ? ? [[alternative HTML version deleted]]