And here is the corrected version of the function, sorry....


asc2ras_off <- function(filename, offset=6) {
 splitasc <- function(s) {
s <- trim(s)
spl <- unlist(strsplit(s, ''))
pos <- which(spl==' ')[1]
first <- substr(s, 1, (pos-1))
second <- substr(s, (pos+1), nchar(s))
return(trim(c(first, second)))
}
 filename <- trim(filename)
    if (!file.exists(filename)) { stop(filename, " does not exist") }
lines <- readLines(filename, n=offset)
ini <- lapply(lines, splitasc)
ini <- matrix(unlist(ini), ncol=2, byrow=TRUE)
ini[,1] = toupper(ini[,1])
 nodataval <- -Inf
for (i in 1:length(ini[,1])) {
if (ini[i,1] == "NCOLS") {nc <- as.integer(ini[i,2])
} else if (ini[i,1] == "NROWS") {nr <- as.integer(ini[i,2])
} else if (ini[i,1] == "XLLCORNER") {xn <- as.numeric(ini[i,2])
} else if (ini[i,1] == "YLLCORNER") {yn <- as.numeric(ini[i,2])
} else if (ini[i,1] == "CELLSIZE") {d <- as.numeric(ini[i,2])
} else if (ini[i,1] == "NODATA_VALUE") { try (nodataval <-
as.numeric(ini[i,2]), silent=TRUE)
}
    }
xx <- xn + nc * d
yx <- yn + nr * d
x <- raster(ncols=nc, nrows=nr, xmn=xn, ymn=yn, xmx=xx, ymx=yx, crs='')

v <- as.numeric( scan(filename, skip=offset, what='character', quiet=TRUE) )
v[v == nodataval ] <- NA
setValues(x, v)
}


On Fri, Oct 21, 2011 at 10:18 AM, Robert J. Hijmans <r.hijmans@gmail.com>wrote:

> Matthias,
>
> To finish my message to slipped away:
>
> If you are done batting, you can use the function below, slight adjustment
> of what 'raster' does, adjusting the offset value (and library(raster) ).
> In future versions of raster you can do
>
> r <- raster('file.asc', offset=6)
>
> Robert
>
>
> On Fri, Oct 21, 2011 at 10:15 AM, Robert J. Hijmans <r.hijmans@gmail.com>wrote:
>
>> I
>>
>>
>>
>> asc2ras_off <- function(filename, offset=6) {
>>  splitasc <- function(s) {
>>  s <- trim(s)
>> spl <- unlist(strsplit(s, ''))
>> pos <- which(spl==' ')[1]
>>  first <- substr(s, 1, (pos-1))
>> second <- substr(s, (pos+1), nchar(s))
>>  return(trim(c(first, second)))
>> }
>>  filename <- trim(filename)
>>     if (!file.exists(filename)) { stop(filename, " does not exist") }
>> lines <- readLines(filename, n=offset)
>>  ini <- lapply(lines, splitasc)
>> ini <- matrix(unlist(ini), ncol=2, byrow=TRUE)
>>  ini[,1] = toupper(ini[,1])
>>  nodataval <- -Inf
>>  for (i in 1:length(ini[,1])) {
>> if (ini[i,1] == "NCOLS") {nc <- as.integer(ini[i,2])
>>  } else if (ini[i,1] == "NROWS") {nr <- as.integer(ini[i,2])
>> } else if (ini[i,1] == "XLLCORNER") {xn <- as.numeric(ini[i,2])
>>  } else if (ini[i,1] == "YLLCORNER") {yn <- as.numeric(ini[i,2])
>> } else if (ini[i,1] == "CELLSIZE") {d <- as.numeric(ini[i,2])
>>  } else if (ini[i,1] == "NODATA_VALUE") { try (nodataval <-
>> as.numeric(ini[i,2]), silent=TRUE)
>> }
>>     }
>> xx <- xn + nc * d
>> yx <- yn + nr * d
>> x <- raster(ncols=nc, nrows=nr, xmn=xn, ymn=yn, xmx=xx, ymx=yx, crs='')
>>
>> v <- as.numeric( scan(filename, skip=offset, what='character', quiet=TRUE)
>> )
>> v[v == x@file@nodatavalue ] <- NA
>>  setValues(x, v)
>> }
>>
>>
>>
>> On Fri, Oct 21, 2011 at 6:08 AM, Barry Rowlingson <
>> b.rowlingson@lancaster.ac.uk> wrote:
>>
>>> On Fri, Oct 21, 2011 at 12:29 PM, Matthias Boer <M.Boer@uws.edu.au>
>>> wrote:
>>> > Dear list members
>>> >
>>> > I have a very large number of ascii map layers representing
>>> interpolated daily weather observations for Australia.
>>> > Each file is an ascii grid with a header consisting of 15 lines (6
>>> lines for the usual header of Arcraster ascii grids plus 9 lines of
>>> additional info). I want to import the maps one by one, to compute monthly
>>> and annual values, overlay those with veg type maps, etc. I am having
>>> trouble with the spatial data formats.
>>> >
>>> > I usually use 'read.asciigrid' to import ascii grids in Arcraster
>>> format but this fails here due to the header format. I can read the actual
>>> grid data using read.delim and create the GridTopology seperately but then
>>> struggle to create a SpatialGridDataFrame or other spatial data format. I
>>> assume there is a simple way to do this, but I haven't managed to work it
>>> out. Some help from the list would be great.
>>>
>>>  GDAL doesn't seem to mind if extra headers come after an initial
>>> NCOLS header. Do your extra bits come before that? That messes up
>>> GDAL's detection. You could bodge this by adding the right NCOLS
>>> header (a duplicate of the NCOLS in the file) at the start of the
>>> file.
>>>
>>>  I also tried opening a file connection, reading the header lines I
>>> stuck in, then tried reading the raster from the connection - both
>>> raster and readGDAL failed, but I would have been surprised if it had
>>> worked.
>>>
>>>  Others have suggested some general text editing tools - I would also
>>> suggest appropriate use of the ClueBat on whoever created these things
>>> - open geospatial data sets with rich metadata options are hardly
>>> rarities!
>>>
>>> Barry
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo@r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>
>>
>

	[[alternative HTML version deleted]]

