# 7. Vector Analysis

# Project set up
# Only if needed:
# install.packages(c("sf", "ggplot2", "dplyr"))

# Load required libraries
library(sf)        # For spatial vector data -> package used for anything related to vector data
library(dplyr)     # For data manipulation
library(ggplot2)   # For plotting

# Step 1: Read the wetlands and roads shapefiles
# Set your data directory, this is another method of getting data from a file path
# Best if you want to access many files from the same folder
data_dir <- "C:/Users/schul/Documents/Berufliches/Hiwi EOCap4Africa/Exercices Finalised/Geodata/7 Vector Analysis"

wetlands <- st_read(file.path(data_dir, "Wetlands.shp"))
roads <- st_read(file.path(data_dir, "Roads.shp"))
aoi <- st_read(file.path(data_dir, "Clip.shp"))

# Step 2: Reproject to a projected CRS suitable for Botswana (e.g. UTM Zone 35S)
wetlands_proj <- st_transform(wetlands, crs = 32735)  # EPSG:32735 = WGS 84 / UTM zone 35S
roads_proj <- st_transform(roads, crs = 32735)
aoi_proj <- st_transform(aoi, crs = 32735)

# Step 3: Create a 500 meter buffer around wetlands
wetlands_buffer <- st_buffer(wetlands_proj, dist = 500)

# Step 4: Clip roads to the Aoi, we do not need the entire road network here
# This is generally good practise as large datasets can slow you down
roads_clipped <- st_intersection(roads_proj, st_union(aoi_proj))

# Step 5: Find intersections between roads and buffered wetlands
intersections <- st_intersection(roads_clipped, wetlands_buffer)

# Optional: Save outputs (shapefiles)
st_write(wetlands_buffer, "output/wetlands_buffer.shp", delete_dsn = TRUE)
st_write(roads_clipped, "output/roads_clipped.shp", delete_dsn = TRUE)
st_write(intersections, "output/road_wetland_intersections.shp", delete_dsn = TRUE)

# Step 6: Plot the result
ggplot() +
  geom_sf(data = wetlands_buffer, fill = "lightblue", color = "blue", alpha = 0.4) + # alpha indicated the transparency
  geom_sf(data = roads_clipped, color = "gray60") +
  geom_sf(data = intersections, color = "red", size = 1) + #s size refers to line thickness
  labs(title = "Wetlands Buffer and Road Intersections in Botswana",
       subtitle = "Buffered wetlands (500m), roads, and intersections",
       caption = "Data source: User-provided shapefiles") +
  theme_minimal()
