植被覆盖度获取
植被覆盖度(Fractional Vegetation Cover,FVC),是指植被(包括叶、茎、枝)在地面的垂直投影面积占统计区总面积的百分比,范围在 [0,1] 之间。FVC 是刻画地表植被覆盖的重要参数,能够直观的反映一个地区绿的程度,是反应植被生长状态的重要指标,在植被变化、生态环境研究、水土保持、城市宜居等方面问题研究中起到重要作用。本案例以 Landsat-8 数据为例,计算贵州省区域的 FVC 指数。
初始化环境
import aie aie.Authenticate() aie.Initialize()
Landsat-8 数据检索
指定区域、时间、云量检索 Landsat-8 ,并对数据进行去云处理。
def removeLandsatCloud(image): cloudShadowBitMask = (1 << 4) cloudsBitMask = (1 << 3) qa = image.select('QA_PIXEL') mask = qa.bitwiseAnd(aie.Image(cloudShadowBitMask)).eq(aie.Image(0)).And(qa.bitwiseAnd(aie.Image(cloudsBitMask)).eq(aie.Image(0))) return image.updateMask(mask)
feature_collection = aie.FeatureCollection('China_Province') .filter(aie.Filter.eq('province', '贵州省')) geometry = feature_collection.geometry() dataset = aie.ImageCollection('LANDSAT_LC08_C02_T1_L2') .filterBounds(geometry) .filterDate('2021-05-01', '2021-10-31') .filter(aie.Filter.lte('eo:cloud_cover', 30.0)) print(dataset.size().getInfo())
dataset = dataset.map(removeLandsatCloud) image = dataset.median()
裁剪影像
image = image.clip(geometry)
计算 NDVI 指数
ndvi = image.normalizedDifference(['SR_B5', 'SR_B4']).rename(['NDVI']) ndvi_vis = { 'min': -0.2, 'max': 0.6, 'palette': ['#d7191c', '#fdae61', '#ffffc0', '#a6d96a', '#1a9641'] } map = aie.Map( center=ndvi.getCenter(), height=800, zoom=6 ) map.addLayer( ndvi, ndvi_vis, 'NDVI', bounds=ndvi.getBounds() ) map
定义植被覆盖度算法
使用像元二分模型法进行 FVC 估算。 利用 aie.Reducer.histogram 实现输入影像的直方图统计。通过 numpy 调用数组运算,计算生长季的 NDVI 像元百分比统计中 5% 位置 NDVI 值作为土壤部分 NDVIsoil 、95% 位置的 NDVI 值作为植被部分 NDVIveg ,并通过 FVC = (NDVI - NDVIsoil)/ (NDVIveg - NDVIsoil ) 计算 FCV ,得出 FVC 。
import numpy as np import pandas as pd
def calculateFVC(image, scale): histogram = image.reduceRegion(aie.Reducer.histogram(2000), None, scale) histogram_info = histogram.getInfo() # print(histogram_info) bucketKey = histogram_info['NDVI_range'] bucketValue = histogram_info['NDVI_counts'] key = np.array(bucketValue) accSum = np.cumsum(key) # print(accSum[20]) # print(accSum[-1]) accPercent = accSum / accSum[-1] p5 = np.searchsorted(accPercent, 0.5) min_ndvi = bucketKey[p5 + 1] # print(min_ndvi) p95 = np.searchsorted(accPercent, 0.95) max_ndvi = bucketKey[p95] # print(max_ndvi) higher_ndvi_mask = image.gt(aie.Image(max_ndvi)) lower_ndvi_mask = image.lt(aie.Image(min_ndvi)) middle_ndvi_mask = aie.Image(1).subtract(higher_ndvi_mask).subtract(lower_ndvi_mask) tmp = image.subtract(aie.Image(min_ndvi)).divide(aie.Image(max_ndvi).subtract(aie.Image(min_ndvi))) FVC = aie.Image(1).multiply(higher_ndvi_mask).add(aie.Image(0).multiply(lower_ndvi_mask)).add(tmp.multiply(middle_ndvi_mask)) return FVC
数据可视化
FVC = calculateFVC(ndvi, 1000) vis_params = { 'min': 0, 'max': 1, 'palette': [ '#a1a1a1', '#008000' ] } map.addLayer( FVC, vis_params, 'fvc', bounds=ndvi.getBounds() ) map
导出数据
task = aie.Export.image.toAsset(FVC, 'FVC_export_result', 100) task.start()
后记
AIE进行遥感云计算的时候还是很方便,可能刚刚出来,很多地方还是需要完善,这个案例里面,我导出数据以后要到ArcGIS里面再出来一下下。接下来,我利用自然间断法分成了五类,然后再统计这五类的面积,这ArcGIS操作都很简单了,这里就不多说,还有就是阿里云的小哥哥特别有耐心,特别负责任,计算也很强。
本案例主要引用AIE官方案例。