I use a PostGIS database with temperature data. I want to use a spatial interpolation method but not the supported IDW method. I need to use kriging, but it is not part of the supported interpolation methods of PostGIS.

I found kriging interpolation method written in R. Now I want to know how can I combine the R code with the data? The problem is that I've never worked with R before.


To clarify: I build 4 climate station with Arduino Megas. They send their climate data to PostGIS over wifi/gsm. The received data should now be spatially interpolated and saved to a new table, so I can use in GeoServer for visualisation. Only the interpolated data is shown on a map. So this all should be done on server side.

  • You can also use GRASS v.krige, look at gis.stackexchange.com/questions/1041/… – Tomek Nov 7 '14 at 10:27
  • If you want to do this on the server side you just need PL/R and some custom scripts, but you'd have trouble plotting variograms and all that. Not easy. Much easier to pull the data out of the DB and work with R on the client side. – Spacedman Nov 7 '14 at 11:54
  • Do you want to interpolate from just four measurement locations? That only gives you six points (4*3/2) for working out the covariance for the Kriging equations. It'll be rubbish. – Spacedman Nov 9 '14 at 20:01

I would write an R script that worked as a client, but will run on the database server. This will save the complication of trying to hook into PostGIS's backend and using PL/R (as I said in comments).

The script will look something like this (which is practically pseudo-code here):

 > con = dbConnect(PG,"localhost","weather")  # connect to local DB
 > m = dbReadTable(con, "readings") # read the `readings` table
 > surf = autoKrige(rain ~ x + y, data=m) # smooth rain over x and y
 > writeRaster(surf, "rain.tif") # save surface as a geoTIFF

ie connect to your database, read the readings, do the kriging (here using autoKrige from the automap package, save as a geoTIFF using the raster package.

Then somehow get that geoTIFF into geoserver.

Now, how do you trigger running this script on the server? Either set up a cron job to do it every ten minutes, or some kind of hook in PostGIS to run it as a system job every time new data is added.

Note that all this is vague and untested but is an outline and should at least give you an idea of the complexity involved. I've probably got some of the function names wrong etc.

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.