[R-sig-Geo] extracting oceanographic data

Michael Sumner mdsumner at utas.edu.au
Thu Jun 26 17:21:34 CEST 2008


LOUZAO wrote:
> Hello Mike,
>
> Thank you for your answer but I have some questions for you:
>>
>> HDF and netCDF can be difficult - you need to compile rgdal from 
>> against FWTools
> Do you mean installing the FWTools software opening the HDF or netCDF 
> files and using later rgdal from R? Could you give me an example, please?

Yes, sorry I mean using readGDAL to read HDF/netCDF directly can be 
complicated since the default binary build won't have those drivers.  I 
would use gdal_translate to convert them to GeoTIFF first. GeoTIFF is 
easily the best raster format overall, but it depends on your needs.

Try this at the command line once FWTools is installed:

gdal_translate file.hdf out.tif

Then in R:

d <- readGDAL("out.tif")  ## with appropriate working directory

It all depends on exactly which files you mean, as it's likely these 
files have subdatasets - you can interrogate these names using gdalinfo, 
e.g.

gdalinfo file.hdf

It's well worth becoming familiar with FWTools and the GDAL command line 
utilities, and also the classes and methods in the sp package. rgdal 
package uses these classes and methods with the underlying GDAL 
framework for file I/O and projections, amongst other things.

That mkLookup function expects an xyz image list (see ?image ), and you 
can convert the results of a readGDAL call to that with 
as.image.SpatialGridDataFrame in the sp package.

I hope some of that helps, I can't spend time on the rest right now.

Cheers, Mike.


>> , or from source to have those work with readGDAL.
> Do you mean those sources of oceanographic information which work with 
> readGDAL? Do you know which website work with that? I use to download 
> information from http://poet.jpl.nasa.gov/, 
> http://oceancolor.gsfc.nasa.gov/cgi/level3.pl,etc.
>> I would just use gdal_translate to convert data to GeoTIFF 
>> (installing FWTools is trivial). 
> gdal_translate is a function in FWTools? GeoTIFF is a georeferenced 
> file type? Can you load them in R? Which function can I use?
>> Those formats can offer further challenges, just in interpreting the 
>> structure and content, but it depends on the actual data you want to 
>> use.
> I would like to use SST, CHL, bathymetry, dynamic height, wind, etc. 
> Until now, I have been downloading data from the website previously 
> cited as ESRI ASC which are loaded into R without problem using 
> import.asc() function.
>>
>> A question to ask here is what is the service providing the data in 
>> the first place? Ideally you can read directly from the server, but 
>> more likely you have a local cache of files. (If you have these as 
>> netCDF/HDF it's probably best to extract what you need to GeoTIFF 
>> using gdal_translate. )  The kftools package provides a direct read 
>> for SST, but it's a bit complicated. 
> I haven't found the kftools package!Do you know where I can find it? 
> Also could you send me an example of how is possible to download 
> directly from a ftp server SST, chlorophyll, wind data, etc.?
>> I had success with reading directly from the service behind the NASA 
>> POET site, but that wasn't useful for many datasets. I haven't 
>> explored the variety of services for some time, apparently there is a 
>> NOAA service that covers most needs.
> Could you explain me a little bit how can I download directly from 
> POET for example? Which is the command line to connect with the ftp?
>>
>> For getting a range of pixel values around each point I would use 
>> overlay() in sp package, and work on creating circle polygons around 
>> your locations.
> Yes, I have tried two different functions. One is polycirc2() from the 
> pgirmess package. I have tried the following code (locs are the lat 
> and long of the locations) and it works (because there is no error) 
> but I cannot visualize anything and also there is no information on 
> the help about the units it is using.
>
> polygon(polycirc2(radius = 1, center = c(x=locs$Long,y=locs$Lat)))
>
> Also I have used the grid.circle() function from the grid package but 
> the output circles seem of the same size that the original one.
>
> Do you know how could I efficiently create circles around my points 
> and later overlay with oceanographic information?
>> You didn't mention the temporal aspect that you need when overlaying 
>> track data with time series of oceanographic data, but this can be 
>> very simple when using rgdal since all your time series layers can be 
>> stored in one SpatialGridDataFrame dataset as separate attributes 
>> (depending on the size) and you simply find the interval in which 
>> your locations fall as well as the pixel for each location. 
>> Determining a spatial range around these is more complicated but not 
>> greatly so. The power of general indexing in R, along with the 
>> spatial extensions provided in sp is very helpful here.
> Yes, this would be very useful. For example, imaging that I have a 
> track corresponding to two weeks and I would like to extract weekly 
> chlorophyll data. First, I would like to match track data with the 
> corresponding chlorophyll week and later extract the data from a 
> determining spatial range around locations of 10 km for example. Could 
> you send me an example of this?
>>
>> How are you determining an appropriate range to query around your 
>> track locations? The tripEstimation package is a suite of functions 
>> written for estimating position from raw archival tags or a remote 
>> service, and I mention that as the outputs there can be very flexible 
>> for spatial/temporal overlays when location precision is important 
>> (but may be overkill here, and it's not user friendly).
>>
>> The trip package extends the foundation provided by sp for track 
>> data, so that also may be useful for dealing with sets of tracks.
> I have been working with the trip package a little bit and thanks for 
> the tripEstimation information.
>>
>> I hope some of that is helpful - sorry I'm not addressing your 
>> questions specifically. If you post more detail about the actual 
>> datasets you are trying to deal with I may be able to help more.
> Thank you for you answer and hope you could help me a little bit with 
> me doubts. It seems to me that it would be more useful to use Matlab 
> instead of R for oceanographic data.
>
> Thanks,
> Maite
>>
>> Cheers, Mike.
>>
>>
>> LOUZAO wrote:
>>> Dear list,
>>>
>>> I'm trying to extract different oceanographic data (chlorophyll, sea 
>>> surface temperature, etc.) for doing habitat modelling. The input 
>>> data are locations (lat, lon) of tracking data of marine top 
>>> predators and I would like to extract oceanographic data from a 
>>> given area around each specific location, from which compute  the 
>>> median, maximum and minimum values of the corresponding data (sst, 
>>> chl, etc.). For doing this,
>>>
>>> 1. I need to import the oceanographic data into R. Which format 
>>> (ascii, hdf, netcdf,etc.) would you recommend? I have been working 
>>> with ascii, importing files with import.asc().
>>> 2. For defining the area around each location, I could define a 
>>> circle (radius of 1° in decimal degree) or a square (a total heigth 
>>> of 1° in decimal degree) centered in each location.
>>> In the grid package, the circleGrob() or the viewport() fucntions 
>>> provide this approach. However, I have applied the circleGrob but I 
>>> do not reach to obtain a circle around each location (when I plot it 
>>> I only obtain the same configuration than when plotting the 
>>> locations). I think that I'm missing something related with the 
>>> units. In this case, x and y are longitud and latitud in decimal 
>>> degress.
>>>
>>> circle<-circleGrob(x=locs$Long, y=locs$Lat, r=2, 
>>> default.units="native", name=NULL,
>>>            gp=gpar(), vp=NULL)
>>>
>>> I have found other two functions that I do not know if someone has 
>>> used for doing the same.
>>>
>>> - symbolsInPolys() in maptools
>>>
>>> - spsample() in sp
>>>
>>> 3. The last step would be to overlay the oceanographic layer and the 
>>> object of circles (or squares) to extract the median, maximum and 
>>> minimum values of the corresponding data (sst, chl, etc.)
>>>
>>> Has anybody worked on this problem? Could you please point to any 
>>> useful documents?
>>> Thanks in advance,
>>>
>>> Maite Louzao
>>>
>>> ------------------------------------------------------------------------ 
>>>
>>>
>>>
>>> No virus found in this incoming message.
>>> Checked by AVG. Version: 8.0.100 / Virus Database: 270.4.1/1515 - 
>>> Release Date: 6/23/2008 7:16 PM
>>>   
>>
>>
>>
>>
>
>
> ------------------------------------------------------------------------
>
>
> No virus found in this incoming message.
> Checked by AVG. 
> Version: 8.0.101 / Virus Database: 270.4.1/1519 - Release Date: 6/25/2008 4:13 PM
>




More information about the R-sig-Geo mailing list