markdown/access_on_gee.md

@@ -39,6 +39,18 @@ var image = ee.Image("projects/soilgrids-isric/layer_name");
 image = ee.Image("projects/soilgrids-isric/layer_name");
 ```
 
+### Export SoilGrids from Google Earth Engine to One Drive
+
+The following script exports selected SoilGrids soil property 
+layers via Google Earth Engine to Google Drive, 
+allowing users to specify the region, soil property variable, 
+resolution, and whether to export each depth layer as a 
+separate file or all layers in one file
+
+https://code.earthengine.google.com/08f74bc377b253523cd4cffe4268f66d
+
 ### Example script in Google Earth Engine code editor
 
+The following script displays a sand content layer from SoilGrids with customized visualization settings
+
 https://code.earthengine.google.com/9e11d149a1d28d10619fa5353b895ed2

markdown/webdav_from_R.md

@@ -4,14 +4,13 @@
 # [WebDAV](https://files.isric.org/soilgrids/data/recent/): direct access to the maps in VRT format.
 
 ### R
-**This code is tested for GDAL version 2.4.x. When running with GDAL version 3.x.x there may be problems. We are working on providing updates**
 
 #### Load libraries, select boundary box and cell size
 Initially you need to load the following libraries and select the bounding box of the area you are interested in:
 
 ```R
-library(rgdal)
-library(gdalUtils)
+library(terra)
+library(gdalUtilities)
 bb=c(-337500.000,1242500.000,152500.000,527500.000) # Example bounding box (homolosine) for Ghana
 igh='+proj=igh +lat_0=0 +lon_0=0 +datum=WGS84 +units=m +no_defs' # proj string for Homolosine projection
 
@@ -26,8 +25,7 @@ gdal_translate(paste0(sg_url,'ocs/ocs_0-30cm_mean.vrt'),
     "./crop_roi_igh_r.tif",
     tr=c(250,250),
     projwin=bb,
-    projwin_srs =igh,
-    verbose=TRUE)
+    projwin_srs =igh)
 ```
 
 #### To a geotiff in a different projection
@@ -41,8 +39,7 @@ gdal_translate(paste0(sg_url,'ocs/ocs_0-30cm_mean.vrt'),
     "./crop_roi_igh_r.vrt",
     of="VRT",tr=c(250,250),
     projwin=bb,
-    projwin_srs =igh,
-    verbose=TRUE)
+    projwin_srs =igh)
 ```
 
 ##### To a VRT in, for example, LatLong
@@ -51,9 +48,10 @@ The following command will generate a VRT in the projection of your choice:
 ```R
 gdalwarp("./crop_roi_igh_r.vrt",
     "./crop_roi_ll_r.vrt", 
-    s_src=igh, 
+    s_srs=igh, 
     t_srs="EPSG:4326", 
-    of="VRT")
+    of="VRT",
+    overwrite = TRUE)
 ```
 
 ##### To a final Geotiff
@@ -69,13 +67,14 @@ gdal_translate("./crop_roi_ll_r.vrt",
 Finally you can read any of the generated VRTs or GeoTiffs in R for further analysis
 
 ```R
-r=readGDAL("./crop_roi_igh_r.vrt") # Or any other of the files produce above
+rst=rast("./crop_roi_igh_r.vrt") # Or any other of the files produce above
+
+plot(rst) # Plot the file
 ```
-​
 
-# To download a global geotiff in Homolosine projection
+#### To download a global geotiff in Homolosine projection
+
 ```R
 gdal_translate(paste0(sg_url,'ocs/ocs_0-30cm_mean.vrt'),
-    "./crop_roi_igh_r.tif",
-    verbose=TRUE)
+    "./crop_roi_igh_r.tif")
 ```

markdown/webdav_from_R_terra.md

+<a href="https://www.isric.org" rel="isric.org"> <img src="https://www.isric.org/themes/custom/basic/logo.svg" height="130"/>       <a href="https://soilgrids.org" rel="soilgrids.org"> <img src="https://www.isric.org/sites/default/files/styles/gallery_big_image_900x700/public/SoilGrids_banner_web.png" height="130"/>
+
+# [WebDAV](https://files.isric.org/soilgrids/data/recent/): direct access to the maps using terra.
+
+### R
+
+By following the four steps outlined below, you will be able to access and download your desired SoilGrids layer using the terra package in R. Step 1 involves specifying the output file path and defining necessary variables such as resolution, variable of interest, depth, and quantile to download your selected SoilGrids layer. In Step 2 you will define a region of interest (roi)\*. If you require global data then skip Step 2 and 3b. Steps 3 and 4 utilise the variables defined in Step 1 (& optionally Step 2) to download your SoilGrids layer. Steps 3 & 4 can be executed by the user without any need for modifications.
+
+*\*Please note that downloading global data may take quite some time. SoilGrids is intended for global applications, however in some cases it is appropriate to use a roi, i.e., continental scale analyses.*
+
+#### Step 1 - Load library & select layer of interest
+
+To begin, load the "terra" library, define the required resolution, output location and output filename, and specify the SoilGrids' layer of interest. You can find a list of available layers at [**https://files.isric.org/soilgrids/latest/data/**](https://files.isric.org/soilgrids/latest/data/).
+
+SoilGrids native resolution is 250m but aggregated products for the **mean only** of selected soil properties are available at 1000m and 5000m as .tif files , check the availability here: [**https://files.isric.org/soilgrids/latest/data_aggregated/**](https://files.isric.org/soilgrids/latest/data_aggregated/)
+
+If your variable of interest is a **soil property** then this should be defined in **Step 1a.** If your variable of interest is a soil class (reference soil group) from the **World Reference Base** (WRB) then this should be defined in **Step 1b**. Please note that SoilGrids and WRB use different coordinate-reference systems, take care that the correct one is assigned ("ESRI:54052" for soil properties; 'EPSG:4326' for WRB reference groups).
+
+``` r
+library(terra)
+
+resolution = "250m" # or "1000m"" or "5000m"
+
+output_folder = "output/folder/location/" # filepath to output folder on your computer
+
+output_filename = "output.tif" # choose a name for your file
+```
+
+##### **1a Define layer of interest for soil properties**
+
+``` r
+voi = "nitrogen" # variable of interest
+depth = "5-15cm" # 0-5cm, 5-15cm, 15-30cm, 30-60cm, 60-100cm or 100-200cm
+quantile = "mean" # mean, Q0.05, Q0.5, Q0.95 or uncertainty
+
+voi_layer = paste(voi,depth,quantile, sep="_") # layer of interest
+
+rst_crs = "ESRI:54052" # ESRI code of the interrupted Goode Homolosine projection used by SoilGrids
+```
+
+##### **1b Define layer of interest for soil class**
+
+*Only available at 250m resolution, please make sure resolution is 250m in Step 1.*
+
+``` r
+voi = "wrb" # this should not change
+soil_class = "Acrisols" # choose soil class or "MostProbable"
+
+voi_layer = soil_class
+
+rst_crs = "EPSG:4326" # EPSG code for the WRB layers
+```
+
+#### Step 2 - Define Extent
+
+Next, you will define the extent of the region of interest (ROI). You have two options for this:
+
+1.  **Step 2a:** If you are already familiar with the minimum and maximum x (longitude) and y (latitude) coordinates, along with the Coordinate Reference System (CRS) of your ROI, you can manually input these values.
+
+2.  **Step 2b:** Alternatively, you can extract the extent from a raster or vector file that represents your ROI.
+
+##### 2a Get bounding box from known extent for region of interest
+
+``` r
+in_crs = "EPSG:4326" # CRS of the bounding box. use EPSG code or WKT 
+
+xmin = 0 #bounding box xmin
+xmax = 2 #bounding box xmax
+ymin = 40 #bounding box ymin
+ymax = 42 #bounding box ymax
+
+bb = ext(c(xmin,xmax,ymin,ymax))
+```
+
+##### 2b Get bounding box from raster or vector file of region of interest
+
+``` r
+
+roi_path = "roi/file/path/roi_filename.gpkg" #set the file path to your roi file (vector or raster format)
+
+roi = vect(roi_path) # use vect(filename) for vector data or rast(filename) for raster
+
+bb = ext(roi)
+
+in_crs = crs(roi)
+```
+
+#### Step 3 - Download SoilGrids layer
+
+In this step the SoilGrids layer is loaded in terra. In step 3b the 'window' function is used to select the region of interest. Skip step 3b if a global coverage is required.
+
+##### 3a Get global SoilGrids layer of interest
+
+``` r
+if (resolution == "250m") { 
+rstFile = paste0("/vsicurl/https://files.isric.org/soilgrids/latest/data/", voi, '/', voi_layer,'.vrt')
+} else if (resolution == "1000m") {
+rstFile = paste0("/vsicurl/https://files.isric.org/soilgrids/latest/data_aggregated/1000m/", voi, '/', voi_layer,'_1000.tif')
+} else {
+rstFile = paste0("/vsicurl/https://files.isric.org/soilgrids/latest/data_aggregated/5000m/", voi, '/', voi_layer,'_5000.tif')
+}
+
+rst = rast(rstFile)
+crs(rst) = rst_crs # assign crs as the ESRI code rather than proj string
+
+plot(rst) # view layer of interest for the whole world (may take a minute or ten)
+```
+
+##### 3b Get SoilGrids layer of interest for given ROI
+
+``` r
+bb_proj = project(bb,from=crs(in_crs),to=crs(rst_crs)) #project bb to same crs as raster layer
+
+window(rst) = NULL # remove any existing window
+window(rst) = bb_proj # get just roi
+
+plot(rst) # quick visual check
+```
+
+#### Step 4 - Save file
+
+Optionally, re-project the layer to the CRS of the ROI (un-comment the line if required), if not, go ahead and save your file.
+
+``` r
+# project to initial crs of roi                    
+#rst = project(rst, in_crs) #optional
+
+# save  your file                  
+writeRaster(rst, paste0(output_folder, output_filename), overwrite = TRUE)
+```
\ No newline at end of file