[R-sig-Geo] kml generation question

[Roger Bivand]

On Mon, 13 Jul 2009, Barry Rowlingson wrote:

> On Mon, Jul 13, 2009 at 3:26 PM, Matt Oliver<moliver at udel.edu> wrote:
>> Dear r-sig-geo,
>>
>> I am trying to generate a kml for an image I am producing using image().
>> After reading about kml generation, I'm unsure if I sould be using
>> kmlOverlay() or writeOGR, or some other function.
>>
>> This is the data
>>
>> http://www.ceoe.udel.edu/cms/moliver/20071003.276.0237.n17.nc
>>
>> require(ncdf)
>> require(fields)
>> f <- open.ncdf("20071003.276.0237.n17.nc")
>>
>> lon <- get.var.ncdf(f, "lon")
>> lat <- get.var.ncdf(f, "lat")
>> mcsst <- get.var.ncdf(f, "mcsst")
>> plot(diff(lon)) #####notice decreasing
>> par(mar=c(0, 0, 0, 0))
>> par(bty="n")
>>
>> image(lon, lat, mcsst, col = tim.colors(64)) #####image I want to make kml
>> for
>>
>> I don't seem to be able to make a proper "Spatial" object because of the
>> unequally spaced geographic coordinates.This seems to be necessary to to
>> proceed with a kml generation.
>>
>> I'm probably missing something simple so any help or example code would be
>> wonderful
>>
>
> KML does this kind of image overlay by taking an image file (PNG,
> jpeg etc) and it's bounding box in lat-long. It can't cope with
> irregular grids.
>
> My first thought as a workaround was to turn every pixel into a
> rectangular polygon. But then I got the data and saw we were dealing
> with 1445626 pixels and that google earth would probably crawl...
>
> So I think you'll need to resample your data onto a regular lat-long
> grid, then save it as an image, and then write the relevant KML with
> the bounds. The KML is quite simple:
>
> http://code.google.com/apis/kml/documentation/kml_tut.html#ground_overlays
>
> but the resampling might not be. I think I've written some resampling
> code somewhere, it just works out where in the old array a bunch of
> regularly spaced points are that will make the new array are, and
> samples. There may be R code for doing interpolation - maybe the
> rimage or RImageJ packages...
>
> But no, I don't think anything existing can do it [waits for Roger to
> prove him wrong again...].

Maybe ... ?

Something like (not yet ideal, some stripey artefacts caused by the cell 
steppings interfering with each other):

require(ncdf)
f <- open.ncdf("20071003.276.0237.n17.nc")
lon <- get.var.ncdf(f, "lon")
lat <- get.var.ncdf(f, "lat")
mcsst <- get.var.ncdf(f, "mcsst")
require(fields)
image(lon, lat, mcsst, col = tim.colors(64))
obj <- list(x=lon, y=lat, z=mcsst)
library(sp)
crds <- expand.grid(lon, lat)
SPDF <- SpatialPointsDataFrame(SpatialPoints(crds, 
proj4string=CRS("+proj=longlat +datum=WGS84")), 
data=data.frame(mcsst=c(mcsst)))
summary(SPDF)
bb <- bbox(SPDF)
grd <- GridTopology(bb[,1], c(diff(bb[1,])/1200, diff(bb[2,])/1200),
   c(1200,1200))
SG <- SpatialGrid(grd, CRS("+proj=longlat +datum=WGS84"))
o <- overlay(SG, SPDF)
oo <- tapply(SPDF$mcsst, o, mean, na.rm=TRUE)
df <- data.frame(mcsst=rep(NA, 1200*1200))
df$mcsst[as.integer(names(oo))] <- oo
summary(df)
SGDF <- SpatialGridDataFrame(grd, CRS("+proj=longlat +datum=WGS84"),
   data=df)
x11()
image(SGDF, col = tim.colors(64), axes=TRUE)


from there follow ?kmlOverlay, using the example as a template, with 
image() instead of plot(). At least it's something.

Roger


>
> Barry
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no