Apply correction to AGBD estimates for selected L4A shots

Apply correction to AGBD estimates for selected L4A shots#

GEDI L4A Version 2 (V002)[Dubayah et al., 2021] files, with granule production version 01, are affected by an issue that results in incorrect estimates of the aboveground biomass density (AGBD) for a subset of GEDI shots. The issue only affects the shots with selected algorithm setting group 10. There was an error in some of the relative height (RH) metrics used as predictor variables in the AGBD models for these shots. AGBD estimates based on other algorithm setting groups are not affected. The Version 1 (V001) GEDI L4A files [Dubayah et al., 2021] are also not affected by this issue. This issue is specific to GEDI L4A only and does not affect GEDI L2 products.

The last two sections of the filenames indicate the granule production and version number of the file, e.g., the bold section of the GEDI04_A_2020207182449_O09168_03_T03028_02_002_01_V002.h5. The GEDI L4A includes estimates and parameters for seven algorithm setting groups; please refer to the user guide and algorithm theoretical basis document (ATBD) for details. The selected algorithm setting group for each shot is provided the variable selected_algorithm within the root group.

Within the root group /BEAMXXXX/, the following variables are affected for the shots with selected_algorithm=10: agbd,agbd_pi_lower,agbd_pi_upper, agbd_se, agbd_t, agbd_t_se, xvar.

Within the /BEAMXXXX/agbd_prediction/ subgroup, the following variables are affected for all shots: agbd_a10, agbd_pi_lower_a10, agbd_pi_upper_a10, agbd_se_a10, agbd_t_a10, agbd_t_pi_lower_a10, agbd_t_pi_upper_a10, agbd_t_se_a10, xvar_a10.

Luckily, this issue only affects algorithm setting group 10 (a10), which is a special case because it uses the same settings as the algorithm setting 5 (a05) but uses a higher mode to calculate the RH metrics when the lowest detected mode is likely noise. This means we can use the xvar data provided for a05 to apply correction to the AGBD estimates for a10.

Learning Objectives

  • While the above issue is being fixed with a new granule production version (02_V002) [Dubayah et al., 2022], we can apply correction to affected files (version 01_V002) [Dubayah et al., 2021] to correct AGBD estimates. This tutorial shows how to apply such corrections to selected variables.

Requirements#

Additional prerequisites are provided here.

import shutil
from os import path
import h5py
import pandas as pd
import numpy as np
from glob import glob

# input directory containing GEDI L4A *01_V002.h5 files
indir = 'full_orbits'

# output directory to save the corrected GEDI L4A files
outdir = 'corrected'

# loop over all h5 files at the input directory
for infile in glob(path.join(indir, 'GEDI04_A*01_V002.h5')):
    
    name, ext = path.splitext(path.basename(infile))
    corrfilename = "{name}_corr{ext}".format(name=name, ext=ext)
    outfile = path.join(outdir, corrfilename)
    
    # copy the input file to output location
    shutil.copy(infile, outfile)
    
    # read the L4A file at output location
    hf_l4a = h5py.File(outfile, 'r+')

    # get the model_parameters
    model_data=hf_l4a['ANCILLARY']['model_data']
    predict_stratum = model_data['predict_stratum'].astype('U13')
    bias_correction_name = model_data['bias_correction_name'].astype('U13')
    bias_correction_value = model_data['bias_correction_value']
    x_transform = model_data['x_transform'].astype('U13')
    y_transform = model_data['y_transform'].astype('U13')
    npar = model_data['npar']
    par = model_data['par']
    rh_index = model_data['rh_index']

    # storing model parameters to a pandas dataframe
    df_l4a_model = pd.DataFrame(list(zip(predict_stratum, bias_correction_name, bias_correction_value,
                                        x_transform,y_transform, npar, par, rh_index)), 
                         columns=['predict_stratum', 'bias_correction_name', 'bias_correction_value',
                                        'x_transform', 'y_transform', 'npar', 'par', 'rh_index'])
    df_l4a_model.set_index('predict_stratum', inplace=True)
    
    for beams in list(hf_l4a.keys()):
        if beams.startswith('BEAM0110'):
            beam = hf_l4a.get(beams)

            # get root-level variables
            agbd = beam['agbd'] 
            agbd_t = beam['agbd_t'] 
            xvar = beam['xvar'] 
            selected_algorithm = beam['selected_algorithm'][:]
            predict_stratum = beam['predict_stratum'][:].astype('U13') 

            # get variables within agbd_prediction group
            predictor_offset = beam['agbd_prediction'].attrs['predictor_offset'] 
            xvar_a5 = beam['agbd_prediction']['xvar_a5'][:]
            xvar_a10 = beam['agbd_prediction']['xvar_a10'] 
            agbd_a10 = beam['agbd_prediction']['agbd_a10']
            agbd_t_a10 = beam['agbd_prediction']['agbd_t_a10'] 

            # get variables within geolocation group
            elev_lowestmode_a10 = beam['geolocation']['elev_lowestmode_a10'][:] 
            elev_lowestmode_a5 = beam['geolocation']['elev_lowestmode_a5'][:] 
            
            # backtransform xvar_a10 from xvar_a5
            xvar_a10_0 = np.square(xvar_a5) - predictor_offset
            xvar_a10_1 = xvar_a10_0 - (elev_lowestmode_a10[:, None] - elev_lowestmode_a5[:, None])
            xvar_a10[:] =  np.sqrt(xvar_a10_1 + predictor_offset) # sqrt x transformation
            
            # compute agbd for a10
            for j, ps in enumerate(predict_stratum):
                if ps in df_l4a_model.index:
                    # get model for the predict stratum ps
                    shot_model = df_l4a_model.loc[ps]
                    # model intercept
                    a_t = shot_model['par'][0] 
                    # apply model coefficients to xvars
                    for i in range(shot_model['npar']-1):
                        a_t += shot_model['par'][i+1] * xvar_a10[j,i]
                    
                    # apply correction to agbd_t_a10
                    agbd_t_a10[j] = a_t
                    # apply correction to agbd_t
                    agbd_a10[j] = np.square(a_t) * shot_model['bias_correction_value']
            
            # apply correction to root variables where selected_algorithm is 10
            xvar[:][selected_algorithm==10] = xvar_a10[:][selected_algorithm==10]
            agbd_t[selected_algorithm==10] = agbd_t_a10[selected_algorithm==10]
            agbd[selected_algorithm==10] = agbd_a10[selected_algorithm==10]
            
#     close hdf file
    hf_l4a.close()

References#

[1] (1,2)

R.O. Dubayah, J. Armston, J.R. Kellner, L. Duncanson, S.P. Healey, P.L. Patterson, S. Hancock, H. Tang, J. Bruening, M.A. Hofton, J.B. Blair, and S.B. Luthcke. GEDI L4A Footprint Level Aboveground Biomass Density, Version 2. 2021. doi:10.3334/ORNLDAAC/1986.

[2]

R.O. Dubayah, J. Armston, J.R. Kellner, L. Duncanson, S.P. Healey, P.L. Patterson, S. Hancock, H. Tang, J. Bruening, M.A. Hofton, J.B. Blair, and S.B. Luthcke. GEDI L4A Footprint Level Aboveground Biomass Density, Version 2.1. 2022. doi:10.3334/ORNLDAAC/2056.

[3]

R.O. Dubayah, J. Armston, J.R. Kellner, L. Duncanson, S.P. Healey, P.L. Patterson, S. Hancock, H. Tang, M.A. Hofton, J.B. Blair, and S.B. Luthcke. GEDI L4A Footprint Level Aboveground Biomass Density, Version 1. 2021. doi:10.3334/ORNLDAAC/1907.