add treecrop blog

Author: PaulMelloy

blog/2025/2503_aRcgis.Rmd

@@ -0,0 +1,171 @@
+---
+title: "2503_aRcgis"
+output: html_document
+---
+
+# ArcGIS in R  
+
+ArcGIS is a set of powerful propriety geospatial tools developed by ESRI. 
+It is widely used in the industry and academia. 
+In this blog, I will explore the `arcgis` R package by getting data from
+the Australian tree crop map database, hosted on ArcGIS.
+
+### Austlralian tree crop map
+
+The [Australian tree crop map](https://www.une.edu.au/research/research-centres-institutes/applied-agricultural-remote-sensing-centre/collaborative-r-and-d-opportunities/industry-applications-and-maps/australian-tree-crops) 
+is a database of tree crops in Australia.
+The map is the work of the [Applied Agricultural Remote Sensing Centre](https://www.une.edu.au/research/research-centres-institutes/applied-agricultural-remote-sensing-centre) 
+at the University of New England.  
+
+The data is hosted on ArcGIS and can be accessed at the following link:
+https://www.arcgis.com/home/item.html?id=12bdb57eec7548dd90409603645cfd5e  
+
+
+### Installation  
+
+Full documentation can be found [at the following link on installation and setup.](https://developers.arcgis.com/r-bridge/installation/)
+This is hosted on R-universe.
+
+`arcgis` is a meta-package for `arcgisutils`, `arcgislayers`, `arcgisgeocode` 
+and `arcgisplaces`  
+
+```{r eval = FALSE}
+install.packages("arcgis", repos = c("https://r-arcgis.r-universe.dev", "https://cloud.r-project.org"))
+```
+
+Next lets load the library into the R environment.  
+
+```{r}
+library(arcgis)
+library(sf) # also needed for working with spatial data
+```
+
+### Query features  
+
+The URL to query point features from the Australian Tree crop is as follows:
+https://services5.arcgis.com/3foZbDxfCo9kcPwP/ArcGIS/rest/services/AARSC_treecrop_point/FeatureServer/0  
+
+```{r}
+atcm_pts <- arc_open("https://services5.arcgis.com/3foZbDxfCo9kcPwP/ArcGIS/rest/services/AARSC_treecrop_point/FeatureServer/0")
+atcm_pts
+```
+
+Here we see this feature server contains point geometry, with a coordinate reference
+system 7844 and has Query capabilities.  
+
+We can obtain more metadata information from the feature server using the `arc_meta` function.
+
+```{r}
+list_fields(atcm_pts)
+```
+
+
+### Selecting data  
+
+The package contains This returns a `sf` simple features object.  
+
+We can query the fields in the `sf` object using the `list_fields` function.  
+
+The where argument provides a feature to filter the data based on SQL syntax.
+
+```{r}
+arc_select(atcm_pts, where = c("commodity = 'avocado'"))
+```
+
+
+#### Query data by spatial reference
+
+From a 10km radius around a point.
+I will choose the centre of Alstonville in Northern New South Wales.  
+
+```{r}
+# Define a point (e.g., coordinates for a location)
+point_coords <- c(153.4404,-28.8421)  # Example coordinates (longitude,latitude)
+
+# Create a spatial points data frame
+# We will use the crs 7844 specified in the ATCM metadata
+alston <- st_sfc(st_point(c(point_coords)), crs = 7844)
+
+# Buffer the point to create a polygon (e.g., 10 kilometers)
+buffer_alston <- st_buffer(alston, 10000)
+
+# transform the data to WGS84 so we can plot it on a leaflet map
+buffer_alston <- st_transform(buffer_alston, 4326)
+
+
+library(leaflet)
+
+leaflet(buffer_alston) |>
+  addProviderTiles(provider = "Esri.WorldImagery")|>
+  addMarkers(lng = point_coords[1], lat = point_coords[2], popup = "Alstonville")|>
+  addPolygons(color = "green", weight = 2, fillOpacity = 0.2)
+```
+
+
+Now lets query and download only the data from inside the buffer polygon.
+
+```{r}
+alston_farms <- 
+  arc_select(atcm_pts, 
+             filter_geom = buffer_alston,
+             predicate = "intersects")
+
+# transform the data to WGS84 so we can plot it on a leaflet map
+alston_farms <- st_transform(alston_farms, 4326)
+
+leaflet(alston_farms) |>
+  addProviderTiles(provider = "Esri.WorldImagery")|>
+  addCircleMarkers(color = "orange", 
+                   weight = 1, 
+                   fillOpacity = 0.8,
+                   radius = 4,
+                   popup = alston_farms$commodity)|>
+  addMarkers(lng = point_coords[1], lat = point_coords[2], 
+             popup = "Alstonville")
+  
+```
+
+Lets do this again, however instead of point data, lets get the polygon data from
+the following link:
+https://services5.arcgis.com/3foZbDxfCo9kcPwP/ArcGIS/rest/services/AARSC_treecrop_polygon/FeatureServer/0
+
+
+```{r}
+atcm_poly <- arc_open("https://services5.arcgis.com/3foZbDxfCo9kcPwP/ArcGIS/rest/services/AARSC_treecrop_polygon/FeatureServer/0")
+list_fields(atcm_poly)
+```
+
+
+Now lets query and download only the data from inside the buffer polygon.
+
+```{r}
+# lets get the data from inside a bounding box of the buffer polygon
+buffer_box <- 
+  st_as_sfc(st_bbox(buffer_alston))
+
+# request only the data in the box
+alston_farms <- 
+  arc_select(atcm_poly, 
+             filter_geom = buffer_box,
+             predicate = "intersects")
+
+# transform the data to WGS84 so we can plot it on a leaflet map
+alston_farms <- st_transform(alston_farms, 4326)
+
+leaflet(alston_farms) |>
+  addProviderTiles(provider = "Esri.WorldImagery")|>
+  addPolygons(color = "orange", 
+                   weight = 1, 
+                   fillOpacity = 0.8,
+                   popup = alston_farms$commodity)|>
+  addPolygons(data = buffer_box, color = "green", weight = 2, fillOpacity = 0.2)
+  
+```
+
+
+
+
+
+
+
+