In envima/GEDItools:

To do: * automatic download function of Sentinel data with matching dates

getting started

#####
###load packages
#####
library(rGEDI)
library(sf)
library(mapview)
library(GEDItools)
library(raster)

#####
###settings
#####

data_path <- "C:/Users/Alice/Uni/Projekte/GEDI/data/rGEDI_2B/"
clip_dir = "rGEDI_2B_clipped" # for clip of whole study area
sen_inpath <- "C:/Users/Alice/Uni/Projekte/GEDI/data/Sentinel/"

##university forest
# ymin<- 50.831518
# ymax<- 50.848372
# xmin<- 8.661830
# xmax<- 8.694962

##max sentinel overlap ##134 gedi orbits
# xmin <- 7.560548 
# xmax <- 9.140458 
# ymin <- 50.45527 
# ymax <- 51.4511

##larger MR area 20 tracks (16.7.2020)
# ymax <- 50.848933
# ymin <- 50.740749
# xmax <- 8.820696
# xmin <- 8.516494

# ##corine_area > 100
# ymax <- 51.027524
# ymin <- 50.747085
# xmax <- 8.256862
# xmin <- 8.972561


# 50.924463, 8.460928 50.742007, 8.835955 33 tracks stand:16.7.2020
ymax <- 50.924463
ymin <- 50.742007
xmax <- 8.835955
xmin <- 8.460928

Preprocessing

Finding and downloading data for study area

#####
###finding the data
#####
gLevel2B <- gedifinder(product="GEDI02_B",ymax, xmin, ymin, xmax, version="001")

#####
###downloading data
#####
gediDownload(filepath=gLevel2B, outdir=data_path)

Workflow for single files

exmpl_file <- "GEDI02_B_2019175071755_O03004_T04466_02_001_01.h5"

Reading, clipping and getting biophysical variables

#####
###reading data
#####
gedilevel2b<-readLevel2B(level2Bpath = file.path(data_path, exmpl_file))

#####
###clipping
#####
file_nm_raw <- substr(exmpl_file, 1,(nchar(exmpl_file)-3)) #filename without .h5 ending
level2b_clip_bb <- clipLevel2B(gedilevel2b, xmin, xmax, ymin, ymax,
                                   output=paste0(data_path, "../", clip_dir, "/", file_nm_raw, "_clip.h5"))

level2BVPM<-getLevel2BVPM(level2b_clip_bb)
# head(level2BVPM[,c("beam","shot_number","pai","fhd_normal","omega","pgap_theta","cover")])

Preparing and plotting

# Converting shot_number as "integer64" to "character"
level2BVPM$shot_number<-paste0(level2BVPM$shot_number)

# Converting GEDI Vegetation Profile Biophysical Variables as data.table to sf object
level2BVPM_sf <- st_as_sf(level2BVPM, coords = c("longitude_lastbin", "latitude_lastbin"), crs = 4326)


# removing -9999 values
level2BVPM_sf <- level2BVPM_sf[(level2BVPM_sf$pai != -9999),]


#plotting
mapview(level2BVPM_sf, zcol = "pai")@map

Workflow for several files for

Functions for that matter can be found in GEDItools Package (https://github.com/envima/GEDItools)

Reading, clipping and getting biophysical variables

#####
###clipping data and saving new hdf5 files with gedi_clip()
#####

gedi_clip(inpath = data_path, xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, clip_dir = clip_dir)

#####
###getting biophysical Variables and merging all files into one dataframe with date to distinguish between the files
#####
bv_df <- gedi_merge_bv(inpath = paste0(data_path, "../", clip_dir)) #package problem...wenn ich funktion einzeln einlade kein problem

#####
###Converting GEDI Vegetation Profile Biophysical Variables as data.table to sf object
##### 
bv_df_sf <- st_as_sf(bv_df, coords = c("longitude_lastbin", "latitude_lastbin"), crs = 4326)

#####
###removing -9999 values
#####
bv_df_sf <- bv_df_sf[(bv_df_sf$pai != -9999),]

#####
###write out .rds
#####
saveRDS(bv_df_sf, file = paste0(data_path, "../rGEDI_2B_clipped/bv_df_sf_big_200717.rds"))
st_write(bv_df_sf, dsn = paste0(data_path, "../rGEDI_2B_clipped/bv_df_sf_big_200717.json"), driver = "GeoJSON")
#readRDS()

mapview(bv_df_sf, zcol = "pai")@map
hist(bv_df_sf$pai)

corine and extract

corine <- raster(paste0("C:/Users/Alice/Uni/Projekte/GEDI/data/CORINE/",
                        "83684d24c50f069b613e0dc8e12529b893dc172f/u2018_clc2018_v2020_20u1_raster100m/",
                        "u2018_clc2018_v2020_20u1_raster100m/DATA/U2018_CLC2018_V2020_20u1.tif"))
bv_sf_prj <- st_transform(bv_df_sf, "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")

extracted <- extract(corine, bv_sf_prj, buffer = 12.5, df = T)
a <- as.numeric(extracted)
bv_sf_prj$cor <- a

#####
###write out .rds
#####
saveRDS(bv_sf_prj, file = paste0(data_path, "../rGEDI_2B_clipped/bv_df_sf_big_200717_cor.rds"))
# st_write(bv_sf_prj, dsn = paste0(data_path, "../rGEDI_2B_clipped/bv_df_sf_big_200717_cor.json"), driver = "GeoJSON", overwrite = T)

table(bv_sf_prj$cor)

crop and stack sentinel files

At this point temporally best fitting Sentinel scenes are downloaded manually.

# crop and stack all relavant sentinel scenes
sen_scenes <- list.dirs(sen_inpath, recursive = F)
for(i in sen_scenes){
  # i <- sen_scenes[1]
  sen_clip_stack(inpath = i, xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, write_out = T)
}

matching with Sentinel files

gedi_dates <- unique(bv_df_sf$date)
date_match_lst <- lapply(setNames(gedi_dates, gedi_dates), function(j){
  match <- sen_match_date(sen_inpath = sen_inpath, date = j)
  return(match)
})