[R-sig-Geo] Rasters average

[Tomislav Hengl]

Try also using RSAGA. I usually process 100 millions of pixels without a 
problem. See:

 > rsaga.get.usage(lib="geostatistics_grid", module=5)
module name: Statistics for Grids

Here is an example:

[http://spatial-analyst.net/scripts/getMOD12C1.R]

HTH,

T. Hengl
http://www.wewur.wur.nl/popups/vcard.aspx?id=HENGL001


Op 9-2-2011 11:05, luca candeloro schreef:
> Hi all,
> I'm trying to use MODIS 1day raster (.tif file) to obtain a, say, weekly
> raster which is the average of the seven 1 day raster.
> Using the following code (library raster) I can't obtain the results due to
> computational limits.
>
> #Input files (same dimension: [1310,1391]
>
> Map.raster1=raster("C:/MOD11A2/Confronto/1d001.tif")
> Map.raster2=raster("C:/MOD11A2/Confronto/1d002.tif")
> Map.raster3=raster("C:/MOD11A2/Confronto/1d003.tif")
> Map.raster4=raster("C:/MOD11A2/Confronto/1d004.tif")
> Map.raster5=raster("C:/MOD11A2/Confronto/1d005.tif")
> Map.raster6=raster("C:/MOD11A2/Confronto/1d006.tif")
> Map.raster7=raster("C:/MOD11A2/Confronto/1d007.tif")
>
>
> for (r in 1:1310) {
> for (c in 1:1391) {
> if(Map.raster1[r,c]==0) v1=NA else v1= Map.raster1[r,c]
> if(Map.raster2[r,c]==0) v2=NA else v2= Map.raster2[r,c]
> if(Map.raster3[r,c]==0) v3=NA else v3= Map.raster3[r,c]
> if(Map.raster4[r,c]==0) v4=NA else v4= Map.raster4[r,c]
> if(Map.raster5[r,c]==0) v5=NA else v5= Map.raster5[r,c]
> if(Map.raster6[r,c]==0) v6=NA else v6= Map.raster6[r,c]
> if(Map.raster7[r,c]==0) v7=NA else v7= Map.raster7[r,c]
> media=mean(c(v1,v2,v3,v4,v5,v6,v7),na.rm=TRUE)
> if (is.na(media)) media=0
> Map.raster[r,c]=media
> }
> }
>
> Is there an easiest way to average the values of "k" raster?
> Thanks,
> Luca.
>
> 	[[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