Speed up R sripts – external program (RSAGA)
What if all the features we have tested/written do not allow us to get satisfactory execution times for the intended task? We can use other programs that have functions we are interested in and have a library implemented in R.
We will show you how to do this using old familiar Zonal Statistics, which we have been using for some time. To process our data, we will use another open source software, namely SAGA GIS (http://www.saga-gis.org).
For this we need the RSAGA library:
install.packages(“RSAGA”)
library(RSAGA)The library requires the program SAGA-GIS installed on our computer. By default, it looks for access to it under the paths:
C:/Program Files/SAGA-GIS
C:/SAGA-GIS
C:/OSGeo4W64/apps/saga
C:/OSGeo4W64/apps/saga-ltrIf you install the library in a different location, it will be searched everywhere on the disk. To speed up the whole process, we can specify the path to the SAGA-GIS environment itself:
my_saga_env = rsaga.env("D:/programs/saga-7.0.0_x64/")Let’s verify that we created the environment correctly by displaying the version of our SAGA+GIS:
rsaga.get.version(env = my_rsaga_env)
"7.0.0"SAGA- GIS allows data processing through its GUI (saga_gui.exe) and command line (saga_cmd.exe). The RSAGA library uses the command line. You can read how to create commands in the documentation for each version, available here.
In short, the tools are grouped into libraries. For each tool, we can define parameters. In our case, we will use the tool “Grid Statistics for Polygons” which is located in the library “shapes_grid”. In R, we will run saga_cmd using the rsaga.geoprocessor function:
rsaga.geoprocessor(lib = "shapes_grid",module = 2,param = list(GRIDS="d:/GIS_in_R/B08.tif",POLYGONS="d:/GIS_in_R/grid_10x100.shp",RESULT="d:/GIS_in_R/grid_100x100_result.shp"),env = my_rsaga_env), where:
lib – is the library with the tool we are interested in
module – is a sequence number of our tool
param – is a list of parameters of the tool
To see what our command should look like, it is best to look at the documentation. At the end of each tool description there is an example of a command:
http://www.saga-gis.org/saga_tool_doc/7.0.0/shapes_grid_2.html
After running the command, you will see a vector file with associated values of raster statistics in the location specified in RESULT. Let’s load these:
library(sf)
shp = read_sf("d:/GIS_in_R/grid_100x100_result.shp")And let’s check the statistics:
shp
Simple feature collection with 1000 features and 9 fields
geometry type: POLYGON
dimension: XY
bbox: xmin: 405440 ymin: 5690759 xmax: 504280 ymax: 5799481
epsg (SRID): 32634
proj4string: +proj=utm +zone=34 +datum=WGS84 +units=m +no_defs
# A tibble: 1,000 x 10
FID B08..CELLS. B08..MIN. B08..MAX. B08..RANGE. B08..SUM. B08..MEAN. B08..VARIAN B08..STDDEV geometry
‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹dbl› ‹POLYGON [m]›
1 0 4 2157 2219 62 8776 2194 601 24.5 ((405460 5690769, 405460 5690768~
2 1 4 1758 2468 710 8434 2108. 67641. 260. ((416440 5690769, 416440 5690768~
3 2 4 2217 2313 96 9042 2260. 1173. 34.2 ((427420 5690769, 427420 5690768~
4 3 4 2561 2735 174 10478 2620. 5017. 70.8 ((438400 5690769, 438400 5690768~
5 4 4 2770 2861 91 11274 2818. 1197. 34.6 ((449380 5690769, 449380 5690768~
6 5 4 2222 2367 145 9184 2296 3534. 59.4 ((460360 5690769, 460360 5690768~
7 6 4 2568 2869 301 10939 2735. 12941. 114. ((471340 5690769, 471340 5690768~
8 7 4 1774 1902 128 7360 1840 3252. 57.0 ((482320 5690769, 482320 5690768~
9 8 4 1970 2113 143 8161 2040. 4542. 67.4 ((493300 5690769, 493300 5690768~
10 9 4 1916 2092 176 7974 1994. 4463. 66.8 ((504280 5690769, 504280 5690768~
# ... with 990 more rowsThe RSAGA library has some SAGA-GIS tools implemented directly in R that can be used without detailed knowledge of the documentation. One such function is rsaga.contour for generating contours. We leave the familiarization with these functions to you.
In the next entry, we will show you how to use functions from other programs in your code that do not yet have a library in R.