How to open data from World Ocean Atlas in R
#################################################### #################################################### # OPEN data from World Atlas NOAA # http://www.nodc.noaa.gov/cgi-bin/OC5/WOA05/woa05.pl #################################################### ### Tim Bover from WOA09 (how to obtain the file from internet, and download to your computer): # - Go to http://www.nodc.noaa.gov/OC5/SELECT/woaselect/woaselect.html # - select the oceanographic variable you are interested in # - select geographic area (or leave whole world), type of data # (probably climatological mean in your case), time period and depth # (depth is not that important. The data files contain all standard # depths available) # - click on 'Show Figure' # - on the resultant page, on the left, above the figure, you will see # 'Data Access: ASCII ArcGIS' # - click on the ASCII link (or left mouse button and 'save target as'). # readLines(file.choose()) # Let the user choose the local file to read (thank's to Luke Miller) # setwd("C:\\Users\\Marta Rufino\\world") h(readLines("data\\WOA05\\temperature.gz")) #this may take some time... # note the path indicates the folder wher you have the file. # You can get the path, by clicking the left mouse button and properties, # using windows explorer, for example or changing R working directory to the folder where you have your file # Note also that the file name will be different from the current one (it is generated automatly by WOA) # Whatever way you have to obtain the path, note that you need to change the bars: instead of \ put \\ or / T=read.table("data\\wfig1323085620.23584.csv.gz", sep=",", col.names=c("lat","lon", paste("t",c(0,10,20,30,50,75,100,125,150,200,250,300,400,500,600,700,800,900,1000,1100,1200,1300,1400,1500,1750,2000,2500,3000,3500,4000,4500,5000,5500), sep=""))) dim(T)# to check if everything is ok... summary(T) ## If we want to work with the data as a spatial object: require(sp)# this package is required for working on most spatial objects. T.spa=T coordinates(T.spa)=~lon+lat gridded(T.spa)=TRUE #proj4string(T.spa)=CRS("+proj=longlat +datum=WGS84") # if we wanted to define the projection... which is not exactly this one for all but summary(T.spa) #image(T.spa, "t0") #for plotting temperature, at zero depth, for example. Note we can plot the 35 depth levels. spplot(T.spa, c("t0","t10"), col.regions=rev(heat.colors(20))) #spplot(T.spa, c("t0","t10"), col.regions=rev(rainbow(100)), at=seq(30,40,l=100),main="Temperature")
This map represents temperature at 0 and 10 m depth ;)
