[R-sig-Geo] Impoting multiple asci grids into R?

Edzer Pebesma edzer.pebesma at uni-muenster.de
Mon Oct 31 09:58:37 CET 2011


So, raster::aggregate does spatial aggregation, and calc does temporal 
aggregation when the stack reflect time. Do you see possibilities to 
extend the behaviour of raster::aggregate to generalize 
stats::aggregate, and e.g. in the temporal domain to extend/reuse 
zoo::aggregate? Besides a grouping predicate ("by"), zoo::aggregate 
allows a grouping function that works on time, giving you e.g. 
aggregations to monthly, quarterly and yearly values without effort.

I tried this approach for the classes in package spacetime, some of 
which generalize the raster time-stack but lack its powerful disk 
caching. For spatio-temporal objects this was documented here:
http://cran.r-project.org/web/packages/spacetime/vignettes/sto.pdf
and for spatial (sp) objects here:
http://cran.r-project.org/web/packages/sp/vignettes/over.pdf

Any comments welcome. Best regards,
--
Edzer

On 10/31/2011 04:31 AM, Robert J. Hijmans wrote:
> Thanks Tom,
>
> "aggregate" does spatial aggregation of cell values, whereas I think Swen
> asked for means of cell values across layers (temporal aggregation, in this
> case). I would do:
>
> library(raster)
> asc_files<- list.files(pattern = glob2rx('EVAPTR_dem30m_*01.sum'))
>
> # with multiple files, use stack rather than brick
> EtRR<- stack(asc_files)
>
> # mean values
> # EtRR.sum<- mean(EtRR)
> # or (currently faster)
> EtRR.sum<- calc(EtRR, mean)
>
> #plot(EtRR.sum)
>
> # write to geotiff (?). For speed's sake, avoid ascii files if possible
> EtRR.sum<- writeRaster(EtRR.sum, "EtRR_sum.tif")
>
> Robert
>
>
> On Sun, Oct 30, 2011 at 3:21 AM, Tomislav Hengl
> <hengl at spatial-analyst.net>wrote:
>
>>
>> Carsten,
>>
>> I highly recommend using the Robert Hijmans' raster package [
>> http://CRAN.R-project.org/**package=raster<http://CRAN.R-project.org/package=raster>]
>> to run such operation on a stack of grids (raster package is really a
>> brilliant implementation of raster processing as your final code looks
>> fairly short, and it does not uses so much RAM; for more info see [
>> http://cran.r-project.org/**web/packages/raster/vignettes/**Raster.pdf]<http://cran.r-project.org/web/packages/raster/vignettes/Raster.pdf%5D>
>> ).
>>
>> This is an example of how to run your stats on grid nodes in 4 lines:
>>
>>> library(raster); library(rgdal)
>> # list all files of interest:
>>> asc_files<- list.files(pattern = glob2rx('EVAPTR_dem30m_*01.**sum'))
>> # Read to R:
>>> EtRR<- brick(lapply(asc_files, raster))
>> # derive a sum:
>>> EtRR.sum<- aggregate(EtRR, fun=mean)
>>> image(EtRR.sum)
>> # write back to ASCII format:
>>> writeRaster(EtRR.sum, "EtRR_sum.asc")
>>
>>
>> Otherwise, SAGA GIS is also efficient in processing large grids:
>>
>>> library(RSAGA)
>>> rsaga.get.usage("**geostatistics_grid", 4)
>> library path:   C:/PROGRA~1/R/R-212~1.2/**library/RSAGA/saga_vc/modules
>> library name:   geostatistics_grid
>> module name :   Statistics for Grids
>> Usage: 4 -GRIDS<str>  [-MEAN<str>] [-MIN<str>] [-MAX<str>] [-VAR<str>]
>> [-STDDEV<str>] [-STDDEVLO<str>] [-STDDEVHI<str>]
>>   -GRIDS:<str>           Grids
>>         Grid list (input)
>>   -MEAN:<str>            Arithmetic Mean
>>         Grid (optional output)
>>   -MIN:<str>             Minimum
>>         Grid (optional output)
>>   -MAX:<str>             Maximum
>>         Grid (optional output)
>>   -VAR:<str>             Variance
>>         Grid (optional output)
>>   -STDDEV:<str>          Standard Deviation
>>         Grid (optional output)
>>   -STDDEVLO:<str>        Mean less Standard Deviation
>>         Grid (optional output)
>>   -STDDEVHI:<str>        Mean plus Standard Deviation
>>         Grid (optional output)
>>
>>
>> see also: http://spatial-analyst.net/**wiki/index.php?title=List_of_**
>> SAGA_modules<http://spatial-analyst.net/wiki/index.php?title=List_of_SAGA_modules>
>>
>> cheers,
>>
>> T. Hengl
>> http://www.wewur.wur.nl/**popups/vcard.aspx?id=HENGL001<http://www.wewur.wur.nl/popups/vcard.aspx?id=HENGL001>
>>
>>
>>
>> On 28-10-2011 14:49, Carsten Neumann wrote:
>>
>>>    Am 27.10.2011 17:12, schrieb Swen Meyer:
>>>
>>>> Dear All,
>>>> I like to import a series of asci grids into R. It is a time series of
>>>> 30 years and I like to calculate the mean monthly sum value for each
>>>> grid node for every month of the year. This is way I´m doing this at the
>>>> moment.
>>>>
>>>> #read asci files
>>>> ####Month01
>>>> EtRR.asci01<-  read.asciigrid("EVAPTR_dem30m_**7101.sum",
>>>> as.image=FALSE,
>>>> plot.image=TRUE)
>>>> EtRR.asci01 at data$EVAPTR7101<- EtRR.asci01 at data$EVAPTR_**dem30m_7101.sum
>>>> #1
>>>> EVAPTR7201<- read.asciigrid("EVAPTR_dem30m_**7201.sum", as.image=FALSE,
>>>> plot.image=TRUE)
>>>> EtRR.asci01 at data$EVAPTR7201<- EVAPTR7201 at data$EVAPTR_dem30m_**7201.sum
>>>> #2
>>>> EVAPTR7301<- read.asciigrid("EVAPTR_dem30m_**7301.sum", as.image=FALSE,
>>>> plot.image=TRUE)
>>>> EtRR.asci01 at data$EVAPTR7301<- EVAPTR7301 at data$EVAPTR_dem30m_**7301.sum
>>>> #3
>>>> .
>>>> .
>>>> .
>>>> #29
>>>> EVAPTR0001<- read.asciigrid("EVAPTR_dem30m_**0001.sum", as.image=FALSE,
>>>> plot.image=TRUE)
>>>> EtRR.asci01 at data$EVAPTR0001<- EVAPTR0001 at data$EVAPTR_dem30m_**0001.sum
>>>> #30
>>>>
>>>> # claculating the mean sum ofer 30 years
>>>> EtRR.asci01 at data$sum<- ((EtRR.asci01[[1]] + EtRR.asci01[[2]] +......+
>>>> EtRR.asci01[[30]]) /30)
>>>>
>>>> My script is working well, but it is pretty time and workspace
>>>> consuming. I tried to use wildcards like this:
>>>>
>>>> txt_files = list.files(pattern = 'EVAPTR_dem30m_*01.sum')
>>>>
>>>> But the result is then:
>>>>
>>>> character(0)
>>>>
>>>> Does anyone maybe have an idea how to import multiple asci-grids in a
>>>> faster a more elegant way?
>>>>
>>>> Thank you in advance,
>>>>
>>>> Swen
>>>>
>>>>
>>>>         [[alternative HTML version deleted]]
>>>>
>>>>
>>>>
>>>> ______________________________**_________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at r-project.org
>>>> https://stat.ethz.ch/mailman/**listinfo/r-sig-geo<https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>>
>>> Dear Sven,
>>>
>>> First: wildcards with * do not work. Just use pattern = 'EVAPTR_dem30m'
>>> to get all files in a folder within a list.
>>> Second: use raster package and then function stack() to read all ascii
>>> grids without accessing memory.
>>>
>>> So you can simply do a layer stack with your list within two r command
>>> lines and without memory handling problems ;-)
>>>
>>> With Best regards
>>>
>>> Carsten
>>>
>>>
>>>
>>>
>>> ______________________________**_________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at r-project.org
>>> https://stat.ethz.ch/mailman/**listinfo/r-sig-geo<https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>
>>
>> ______________________________**_________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/**listinfo/r-sig-geo<https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>
>
> 	[[alternative HTML version deleted]]
>
>
>
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
8333081, Fax: +49 251 8339763  http://ifgi.uni-muenster.de
http://www.52north.org/geostatistics      e.pebesma at wwu.de



More information about the R-sig-Geo mailing list