[R-sig-Geo] DEMs

[Tomislav Hengl]

> -----Original Message-----
> From: Aben Woha [mailto:abenwoha at gmail.com]
> Sent: Tuesday, September 15, 2009 7:52 AM
> To: Tomislav Hengl
> Cc: r-sig-geo at stat.math.ethz.ch
> Subject: Re: [R-sig-Geo] DEMs
> 
> Hi Tomislav,
> 
> Thanks for the detailed response.  I tried to use the globedem.zip from your website, but then
> when I call overlay() with my data, there are many NAs produced in the output: it seems the DEM
> covers  X= -180 to 180 & Y= -65 to 65, my dataset covers X= -180 to 180 & Y= -90 to 90
> So when I call krige, it complains:
> Error: dimensions do not match: locations 129600 and data 47029
> 
> 129600 = (number of sample points) * 2
> 47029 = (number of sample points) - (number of NAs in overlay dem data)
> 
> Any suggestions on how to handle this?  Should I subset the data?  Find a different dem with full
> global coverage?

You should mask the NA values prior to variogram fitting using e.g.:

> sel <- !is.na(point.ov$SOLAR)
> plot(variogram(SOLAR~1, point.ov[sel,]))

SRTM DEM is available only for 56°S to 60°N. The 10 minutes global DEM (GTOPO?) at is probably your
best choice then (it is 2160 x 900 pixels):

http://biogeo.berkeley.edu/worldclim1_4/grid/cur/alt_10m_bil.zip  (metadata missing :( )

> 
> Also when you say a global dataset <5km resolution, is that as measured at the equator?

Yes. The ETOPO1 (10 km) DEM is 3600 x 1300 pixels, which is already on the edge of what R can
handle. There is also the GDEM (ASTER DEM, which has global coverage), but nobody yet produced a
single layer (e.g. at resolution of 1 km).

> 
> The current variable I am working with is solar radiation.  I have plans to work on several other
> variables as well.
> 

But do you need to produce maps of solar radiation for artics/antartica, sea etc? Besides, plenty of
global solar radiation maps already exist e.g.:

1. http://www.meteonorm.com/pages/en/downloads/maps.php
2. http://re.jrc.ec.europa.eu/pvgis/countries/europe.htm 
...

> 
> On Mon, Sep 14, 2009 at 1:58 AM, Tomislav Hengl <hengl at spatial-analyst.net> wrote:
> 
> 
> 
> 	Hi Abe,
> 
> 	Did you try using the ETOPO1? This (0.1 degree) DEM can be well loaded to R:
> 
> 	> URL <- "http://spatial-analyst.net/worldmaps/"
> 	# list of maps:
> 	> map.list <- c("globedem", "chlo08", "dcoast")
> 	# download the zipped maps one by one:
> 	> for(i in 1:length(map.list)) {
> 	>  download.file(paste(URL, map.list[i], ".zip", sep=""),
> 	+   destfile=paste(getwd(), "/", map.list[i], ".zip", sep=""))
> 	>  unzip(zipfile=paste(map.list[i], ".zip", sep=""), exdir=getwd())
> 	>  unlink(paste(map.list[i], ".zip", sep=""))
> 	> }
> 	# read maps to R:
> 	> worldmaps <- readGDAL(paste(map.list[i], ".asc", sep=""))
> 	# fix the layer name:
> 	> names(worldmaps)[1] <- map.list[1]
> 	> for(i in map.list[-1]) {
> 	> worldmaps at data[map.list[i]] <- readGDAL(paste(map.list[i], ".asc",
> 	sep=""))$band1
> 	> }
> 	> proj4string(meuse.grid) <- CRS("+proj=longlat +ellps=WGS84")
> 
> 	I do not think that it will be easy to generate kriging with any global
> 	dataset that is <5 km resolution. You could try using the kriging in SAGA
> 	(it can handle maps up to 2 GB in size; it is exp(2) times faster than R).
> 	Here is an example:
> 
> 	http://spatial-analyst.net/wiki/index.php?title=Interpolation_of_ISRIC-
> WISE_international_soil_profile_data
> 
> 	T. Hengl
> 
> 	PS: Which variable are you interpolating?
> 
> 
> 	> Hello All,
> 	>
> 	> My questions are regarding the best way to obtain DEMs for the entire
> 	> globe.  I am trying to krige a dataset of global scope.  The dataset has
> 	> points at all latitude/longitude intersects and so there are 64,800 sample
> 	> points.  The problem that I am having is that the DEMs I have been using
> 	> are
> 	> too large (memory) to be used as R variables.  I also have had trouble
> 	> finding a way to download DEMs for the entire globe.
> 	>
> 	> 1.  Given the size of my sample dataset: What is the maximum resolution of
> 	> DEM that I should be using?  Do I only need values at prediction points or
> 	> will higher res yield higher accuracy?
> 	>
> 	> 2.  Is there a way to programmatically access DEMs using GDAL or R in
> 	> general?  ie. internally or through integration with McIDAS, NetCDF, etc.
> 	>
> 	> 3. Do I need to break the DEMs and sample dataset into smaller tiles to
> 	> avoid the memory issues?  If so how do I do that in R / GDAL?
> 	>
> 	> Thanks,
> 	> Abe
> 	>
> 
> 	>       [[alternative HTML version deleted]]
> 	>
> 	> _______________________________________________
> 	> R-sig-Geo mailing list
> 	> R-sig-Geo at stat.math.ethz.ch
> 	> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 	>
> 
> 
> 
>