[R-sig-Geo] reproject coordinates from meters to decimal degrees

Thu Jan 5 12:40:31 CET 2017

On Thu, 5 Jan 2017 at 22:24 Roger Bivand <Roger.Bivand at nhh.no> wrote:

> On Thu, 5 Jan 2017, Michael Sumner wrote:
>
> > On Thu, 5 Jan 2017 at 19:19 Nikolioudakis Nikolaos <
> > nikolaos.nikolioudakis at imr.no> wrote:
> >
> > Thank you for the prompt response Mike. It would be great if you could do
> > so and I can later post the solution on the list in case someone needs it
> > in the future
> >
> >
> >
> >
> > Ugh, it's not what I thought, raster is just wrong on what this file has
> in
> > it, but there's no point laying blame, NetCDF is so general and so
> commonly
> > abused that software everywhere has endless heuristics to figure out what
> > to assume.
> >
> > Basically your data ship with longitude/latitudes, and you need to figure
> > out what that is in the native polar stereographic projection as a
> regular
> > grid for use in raster().  You can't easily reproject that raster as it
> is,
> > and I'd advise against it since it's lossy, and doesn't scale well to
> long
> > time series, you could use nn techniques to transform to the longlat
> values
> > in a brute force way, but that's a horrible hack - I think it's worth
> > figuring this one out.
> >
> > I've put my startings here: https://github.com/mdsumner/ranch/issues/3
> >
> > If I can I'll figure it out and let you know.  Feel free to ping on that
> > issue, but also share more info you find from the sources -  it's a
> > repository I made to collect NetCDF adventures - it  might be worthwhile
> .
> > . .
> >
> > Similar data (from raw binary but using the same-mostly projection) is
> > wrapped by us here:
> >
> https://github.com/AustralianAntarcticDivision/raadtools/blob/master/R/ice.R
> > This is also built into raster itself (but only to read those particular
> > binary files, not NetCDF ones).
> >
> > I kind of agree with Roger about a vignette for rgdal, but really the
> > information belongs at the source of the data (and inside the NetCDF file
> > itself), and I'd be more interested if rgdal was easier to work with as a
> > collaborator. (Most users would think it belongs in raster, since that
> > hides rgdal for the most part - neither package is easy to contribute
> to).
> >
> > I agree we need a better community resource for these file-bogeys, but
> > it's something we share with the GDAL community and many others so it's
> > bigger than any particular software user group.
>
> Yes, I agree. Does the use of +proj=obtran - oblique transformations -
> come into the same general area of file-bogeys which could be documented
> in the same place?
>
>
Yes it does, I just don't seem to have the capacity to collect them in a
sensible way. I have a little here and there:

https://github.com/mdsumner/angstroms

https://github.com/mdsumner/roc

https://github.com/AustralianAntarcticDivision/raadtools

Perhaps Edzer would like us to collect file-bogeys as part of
http://r-spatial.org/ for future developments?

If there was a single place that people would actually look at I'd get
behind it, I hit these kinds of examples regularly.

Cheers, Mike.

Roger
>
> >
> > (For this problem generally, we just collect our fixes in our local tool
> > kit (raadtools) for files that will be of long-term use, and because it
> > provides the rest of the workflows we need on top of huge file
> collections
> > on local shares. I don't see any general solutions that are better than
> > that on the horizon, though incremental improvements to raster, GDAL and
> > NetCDF have seen a lot of our workarounds updated).
> >
> > Cheers, Mike.
> >
> >
> >
> > Cheers,
> >
> > Nikos
> >
> >
> >
> > *From:* Michael Sumner [mailto:mdsumner at gmail.com]
> > *Sent:* Thursday, January 05, 2017 9:01 AM
> > *To:* Nikolioudakis Nikolaos <nikolaos.nikolioudakis at imr.no>;
> > r-sig-geo at r-project.org
> > *Subject:* Re: [R-sig-Geo] reproject coordinates from meters to decimal
> > degrees
> >
> >
> >
> > This data is in polar stereographic so it's best to ignore the lonlats
> and
> > apply the right affine transform (extent and resolution) directly. I.e.
> no
> > transformation or hacky workarounds required. We do this routinely from
> > similar files, happy to explain off-list.
> >
> > For annoying reasons this is poorly understood on both the GIS and NetCDF
> > (modelling) sides.
> >
> > Cheers, Mike
> >
> >
> >
> > On Thu, 5 Jan 2017, 18:43 Nikolioudakis Nikolaos <
> > nikolaos.nikolioudakis at imr.no> wrote:
> >
> > Dear members of the list,
> > My question might be trivial but please bear with me as I am quite new in
> > GIS and R so my terminology and/or expressions might not be the right
> ones.
> > I have a set of coordinates in lat/lon for which I need to extract values
> > from a netCDF file. The file is downloaded from
> http://marine.copernicus.eu/
> > (needs registration to reach data from there). For this reason you can
> find
> > in the following link a .nc subset (single raster) from a multiband .nc I
> > am using.
> > .nc file:
> >
> > https://dl.dropboxusercontent.com/u/2156982/subset.nc
> >
> > The coordinates of my locations (in lon/lat) are in the coords.csv file
> can
> > be found here:
> >
> > https://dl.dropboxusercontent.com/u/2156982/coords.csv
> >
> > The raster's projection is (according to the data source EPSG: 3411
> (NSIDC
> > Polar Stereographic North)) so the obvious thing to do is reproject my
> > points to the same projection as the raster's in order to overplot them
> > (and subsequently extract values from the raster based on the points).
> > Reprojecting the rasters would also be an option, but in my specific
> case I
> > have multiband (many dates) netCDFs and also many different variables in
> > each netCDF so reprojecting the coordinates if possible and performing
> the
> > extraction from the epsg 3411 rasters would be ideal for me. My problem
> is
> > that when I reproject my coordinates I get them in meters and not in
> > decimal degrees as the raster is.
> >
> >
> >
> > So the question is how can I obtain my reprojected coordinates in decimal
> > degrees instead of meters in order to be able to perform the extraction
> of
> > values from the raster (I have done this for another dataset that both
> the
> > locations and the rasters were in the same projection and units). I
> > understand that my problem is probably related with how to declare that
> the
> > output units of my reprojected coordinates should be in dd instead of
> > meters but I am not sure on the syntax of the proj4 string and on the web
> > the inverse procedure is usually described (dd to meters). The r code I
> am
> > using is the following:
> >
> > library(raster)
> >
> > library(rgdal)
> >
> >
> >
> > #setwd to the location of the files
> >
> > subset.r = raster('subset.nc', varname='chla')
> >
> > #print(subset.r)
> >
> >
> >
> > locs <- read.table('coords.csv',header = T,sep = ';')
> >
> > coordinates(locs) <- ~ LON+LAT
> >
> > # par(mfrow=c(1,2))
> >
> > # plot(subset.r)
> >
> > # plot(locs at coords)
> >
> >
> >
> > #set the projection to EPSG: 4326 for the locs
> >
> > wgs84 <- CRS('+init=epsg:4326')
> >
> > proj4string(locs) <- wgs84
> >
> > #reproject locs to polar stereographic (EPSG 3411)
> >
> > locs.repr <- spTransform(locs, crs(subset.r))
> >
> > #compareCRS(locs.repr, subset.r)
> >
> > plot(subset.r)
> >
> > plot(locs.repr at coords<mailto:locs.repr at coords>)
> >
> > # extent(locs.repr)
> >
> > # extent(subset.r)
> >
> >
> >
> > Any help would be greatly appreciated as well as any comments on things I
> > might be forgetting to provide so as to make my problem better
> understood.
> >
> > Best regards,
> > Nikos
> > __________________________
> > Nikos Nikolioudakis
> > PostDoc Researcher
> > Pelagic Fish Research Group
> > Institute of Marine Research (IMR)
> > P.O. Box 1870 - Nordnes
> > 5817 Bergen, Norway
> > E-mail:nikolaos.nikolioudakis at imr.no
> > Tel:0047489981 <0047489981>
> >
> > [cid:image001.jpg at 01D0B00F.4BA30FC0]
> > __________________________
> >
> >        [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo at r-project.org
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
> >
>
> --
> Roger Bivand
> Department of Economics, Norwegian School of Economics,
> Helleveien 30, N-5045 Bergen, Norway.
> voice: +47 55 95 93 55 <+47%2055%2095%2093%2055>; e-mail:
> Roger.Bivand at nhh.no
> http://orcid.org/0000-0003-2392-6140
> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
> http://depsy.org/person/434412
>
-- 
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

	[[alternative HTML version deleted]]