前言
對(duì)森林高度進(jìn)行持續(xù)的大尺度監(jiān)測(cè)，對(duì)于估算森林碳排放技竟、分析森林退化情況以及量化森林恢復(fù)措施的效果具有重要意義。星載激光雷達(dá)GEDI（Global Ecosystem Dynamics Investigation）通過(guò)發(fā)射激光脈沖俗壹，記錄植被和地表反射回波的能量恰梢。利用森林冠層與地表回波的時(shí)間差，可以計(jì)算出森林冠層的高度酿秸。本文將介紹如何下載和讀取GEDI_01B數(shù)據(jù)灭翔，并講解使用GEDI數(shù)據(jù)提取冠層高度的相關(guān)理論知識(shí)。
數(shù)據(jù)下載與讀取
如圖所示辣苏，GEDI提供了不同級(jí)別的產(chǎn)品數(shù)據(jù)肝箱，其中01_B是經(jīng)過(guò)地理校正的全波形數(shù)據(jù)哄褒，通常用于更靈活的處理和應(yīng)用。此外煌张，GEDI提供的樹(shù)高和生物量產(chǎn)品數(shù)據(jù)也具有很高的應(yīng)用價(jià)值呐赡。
1、GEDI01_B數(shù)據(jù)下載
下載地址：https://search.earthdata.nasa.gov/
先注冊(cè)賬號(hào)
在搜索窗口骏融，搜索GEDI链嘀，并填上相應(yīng)的篩選條件，最后進(jìn)行下載
2绎谦、數(shù)據(jù)讀取
在Python中使用h5py或pyGEDI這兩個(gè)庫(kù)可以實(shí)現(xiàn)01_B數(shù)據(jù)的讀取管闷，但是h5py更應(yīng)用性更高，pyGEDI可以讀取窃肠，并展示高程與回波能量的圖包个，但是我怎么也沒(méi)找到地理信息數(shù)據(jù)的讀取與導(dǎo)出函數(shù)，單純用pyGEDI應(yīng)用很受限冤留。
2.1 pyGEDI
在使用pyGEDI需要安裝gdal庫(kù)碧囊，但是即使安裝好，你運(yùn)行也會(huì)顯示沒(méi)有g(shù)dal模塊纤怒，這時(shí)你需要在安裝pyGEDI包的位置打開(kāi)int文件糯而。打開(kāi)文件后在開(kāi)頭加入代碼：from osgeo import gdal，就可以了泊窘。因?yàn)檫@個(gè)文件中是直接import gdal熄驼，這是無(wú)法引入的。
from pyGEDI import *
import matplotlib.pyplot as plt
fileh5_1B=r'C:\Users\YSY\Desktop\pyGEDI-master\pyGEDI-master\notebook\data\GEDI01_B_2019108080338_O01964_T05337_02_003_01_sub.h5'
h5_1B=getH5(fileh5_1B)
print(h5_1B)
shot_number=19640521100108408
rxwaveform,elevation = waveForm(shot_number,h5_1B)
beam=getBeam(shot_number,h5_1B)
print(beam)
fig = plt.subplots(figsize=(7,7))
plt.plot(rxwaveform,elevation, 'green')
plt.xlabel("Waveform Amplitude")
plt.ylabel("Elevation (m)")
plt.show()
2.2 h5py
import h5py
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
# d讀取h5格式的波形數(shù)據(jù)
fileh5_1B = r'C:\Users\YSY\Desktop\pyGEDI-master\pyGEDI-master\notebook\data\GEDI01_B_2019108080338_O01964_T05337_02_003_01_sub.h5'
with h5py.File(fileh5_1B, 'r') as h5_file:
# 選擇激光束
beam_group = h5_file['BEAM0101']
# 選擇光斑
shot_number = 19640521100108408
shots = beam_group['shot_number'][:]# 可以提前打印看光斑號(hào)
# 找到先對(duì)應(yīng)的光斑
shot_index = np.where(shots == shot_number)[0][0]
# 提取地理位置和回波能量
rx_sample_count = beam_group['rx_sample_count'][shot_index]
rx_start_index = beam_group['rx_sample_start_index'][shot_index]
rx_waveform = beam_group['rxwaveform'][rx_start_index:rx_start_index + rx_sample_count]
latitudes = beam_group['geolocation']['latitude_bin0'][shot_index]
longitudes = beam_group['geolocation']['longitude_bin0'][shot_index]
elevation_bin0 = beam_group['geolocation']['elevation_bin0'][shot_index]
elevation_lastbin = beam_group['geolocation']['elevation_lastbin'][shot_index]
# Generate elevation range
elevation_range = np.linspace(elevation_lastbin, elevation_bin0, rx_sample_count)
# 將地理位置和回波能量導(dǎo)出
data = {
'Latitude': [latitudes] * len(rx_waveform),
'Longitude': [longitudes] * len(rx_waveform),
'Elevation (m)': elevation_range,
'Waveform Amplitude': rx_waveform
}
df = pd.DataFrame(data)
df.to_csv('gedi_waveform_data.csv', index=False)
# Plot waveform amplitude vs. elevation
plt.figure(figsize=(7, 7))
plt.plot(rx_waveform, elevation_range, color='green')
plt.xlabel("Waveform Amplitude")
plt.ylabel("Elevation (m)")
plt.title(f"Waveform Amplitude vs Elevation for Shot {shot_number}")
plt.gca().invert_yaxis() # Invert y-axis for correct elevation view
plt.show()
這幾天我翻閱資料有以下幾點(diǎn)需要在注意：
- 01_B數(shù)據(jù)需要根據(jù)“stable_return_flag”和“degrade”的值進(jìn)行質(zhì)量過(guò)濾烘豹，通常會(huì)刪除這些值為0的行瓜贾，以確保數(shù)據(jù)質(zhì)量；
- GEDI光斑的地理位置信息通常具有15至20米的誤差携悯，因此在應(yīng)用中需考慮這種空間偏差祭芦；
- BEAM0101、BEAM0110憔鬼、BEAM1000龟劲、BEAM1011是全功率脈沖（Full power beam），具有更強(qiáng)的穿透能力轴或，由兩臺(tái)激光器生成昌跌，并快速偏轉(zhuǎn)一次形成。而B(niǎo)EAM0000照雁、BEAM0001避矢、BEAM0010、BEAM0011是覆蓋光束（Coverage beam），由一臺(tái)激光器分成兩束审胸，并以1.5毫弧度（mrad）的角度偏轉(zhuǎn)亥宿，生成四條光束。
結(jié)果如下圖所示砂沛，高度變化相對(duì)正常烫扼，每個(gè)bin之間的差距約為15厘米，這與GEDI記錄時(shí)間間隔為1ns相對(duì)應(yīng)碍庵。
3映企、冠層高度提取
在我撰寫(xiě)的第一篇推文《利用GEDI數(shù)據(jù)反演森林冠層高度》中，已經(jīng)介紹了如何在GEE平臺(tái)上使用GEDI數(shù)據(jù)與被動(dòng)光學(xué)影像（如Sentinel-2）融合静浴，實(shí)現(xiàn)大范圍高分辨率的森林冠層高度制圖堰氓。那時(shí)我還不熟悉Markdown的使用，代碼是以圖片形式上傳的苹享，也缺少理論的詳細(xì)解釋双絮，當(dāng)時(shí)心情比較緊張。提取冠層高度時(shí)得问，使用的是GEDI_02A數(shù)據(jù)囤攀，關(guān)鍵在于我們應(yīng)使用相對(duì)高度中的哪一個(gè)作為森林冠層高度值。根據(jù)《Mapping global forest canopy height through integration of GEDI and Landsat data》一文宫纬，參考基于機(jī)載點(diǎn)云計(jì)算的冠層高度焚挠，使用RH95高度作為冠層高度會(huì)更加準(zhǔn)確。
那么漓骚，為什么不直接使用提取的RH100蝌衔，而是要將高度乘以一個(gè)百分比呢？
我認(rèn)為主要有兩個(gè)原因：第一是脈沖寬度引起的誤差蝌蹂，即使激光打到平地噩斟，回波能量的波形仍會(huì)顯示大約四米的高度；第二是坡度引起的誤差叉信，坡度會(huì)導(dǎo)致土壤層的回波寬度增加亩冬，從而影響測(cè)量精度艘希。
因此硼身，為了減少這些誤差，我們通常會(huì)將波形高度乘以一個(gè)小于1的比值（如0.95覆享，RH95）佳遂，作為冠層高度。注意一點(diǎn)使用的波形數(shù)據(jù)對(duì)應(yīng)位置的森林林木狀態(tài)應(yīng)該未落葉撒顿，否則高度會(huì)被低估丑罪。
