6

I'm working with data on sea surface currents from ftp://podaac-ftp.jpl.nasa.gov/allData/oscar/preview/L4/oscar_third_deg The data has the latitudes and longitudes saved in the following format:

latitude  Size:481
units: degrees-north
long_name: Latitude
longitude  Size:1201
units: degrees-east
long_name: Longitude

head(lon)
[1] 20.00000 20.33333 20.66667 21.00000 21.33333 21.66667
tail(lon)
[1] 418.3333 418.6667 419.0000 419.3333 419.6667 420.0000

How can I convert this into latitude-longitude units in R?

The lon object goes from 20 to 420 whereas I want it to be from -180 to 180.

I have seen Convert North and East coordinates to longitude and latitude in R at Stack Overflow but I am unable to adapt it to my setting.

x <- nc_open("oscar_vel1992.nc")
# lon
lon <- ncvar_get(x,"longitude")
#lat
lat <- ncvar_get(x,"latitude")
z <-as.numeric(char2dms(lat))
6

The data goes from 20 to 420, which is including some areas more than once because the metadata says:

NC_GLOBAL#NOTE2=Longitude extends from 20 E to 420 E to avoid a break in major ocean basins. Data repeats in overlap region.

So we need to split the data into the 20 to 360 section and the 360 to 380 section and put that 360 to 380 section on the left of the 20 to 360 section. We don't need 380 to 420 because that's overlap.

So...

library(raster)

define two extents to cut the raster with. The "Y" extent just has to be bigger than the raster so I go up to 90:

e1 = extent(c(xmin=20,xmax=360,ymin=-90,ymax=90))
e2 = extent(c(xmin=360, xmax=360+20, ymin=-90, ymax=90))

Load the raster - I'm loading just the "vm" layer, you'll have to repeat over the others:

d = raster("./oscar_vel7001.nc",varname="vm")

Now split the bits that we want using the extents:

d1 = crop(d,e1)
d2 = crop(d,e2)

Next shift the bit on the right to fill the space from 0 to 20:

xmin(d2)=xmin(d2)-360
xmax(d2)=xmax(d2)-360

merge will paste the two parts together, giving us a raster from 0 to 360, and rotate will shift it to be -180 to 180:

dm = rotate(merge(d1,d2))
plot(d)
plot(dm)

shows the original shifted overlapping map and the corrected -180 to +180 map:

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.