library(terra)
## terra 1.7.78
library(stringr)

aggregate_byfile <- function(filename_vegheight, raster_target, outdir){
  
  # Load the two raster files
  r1 <- rast(filename_vegheight)  # obtained from   https://libdrive.ethz.ch/index.php/s/cO8or7iOe5dT2Rt?path=%2F
  
  # Aggregate r1 to match r2's resolution
  # The factor is determined by the ratio of resolutions
  fact_x <- (xres(raster_target) / xres(r1))
  fact_y <- (yres(raster_target) / yres(r1))
  
  # Aggregate using mean (can be changed to other functions like max, min, sum)
  r1_aggregated <- aggregate(r1, fact = c(fact_x, fact_y), fun = mean, na.rm = TRUE)
  
  # # Optional: Align extents using resample
  # r1_resampled <- resample(r1_aggregated, r2, method = "bilinear")  # Options: "bilinear", "near", etc.
  
  # create output file name and write to file
  outfilnam <- paste0(outdir, str_remove(basename(filename_vegheight), ".tif"), "_15arcsec.nc")
  message(paste("Writing to file", outfilnam, "..."))
  writeCDF(r1_aggregated, outfilnam, overwrite = TRUE)
}

raster_target <- rast("~/data/archive/gti_marthews_2015/data/ga2.nc")   # The target resolution raster

# small file
aggregate_byfile(
  filename_vegheight = "~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E006_Map.tif",
  raster_target,
  "~/Downloads/"
)
## Writing to file ~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E006_Map_15arcsec.nc ...
# bigger file
aggregate_byfile(
  filename_vegheight = "~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E009_Map.tif",
  raster_target,
  "~/Downloads/"
)
## Writing to file ~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E009_Map_15arcsec.nc ...
# test
rasta <- rast("~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E009_Map_15arcsec.nc")
plot(rasta)

Combine multiple rasters into a mosaic.

raster_files <- c(
  "~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E006_Map_15arcsec.nc",
  "~/Downloads/ETH_GlobalCanopyHeight_10m_2020_N00E009_Map_15arcsec.nc"
  )

# Load rasters as a SpatRaster object
rasters <- lapply(raster_files, rast)

# Merge rasters into a single mosaic
mosaic_raster <- do.call(merge, rasters)

# Save the merged raster as a GeoTIFF or NetCDF file
writeCDF(mosaic_raster, "mosaic_output.nc", overwrite = TRUE)

# Plot the mosaic raster
plot(mosaic_raster)