在地球引擎/GIS中提取样带数据



我正在努力研究如何在谷歌地球引擎中沿着样带提取数据——最好是以定义的间隔提取。使用Jupyter notebook中的PythonAPI,我绘制了一些数据,缓冲了一个点(以定义感兴趣的区域(,并绘制了一个样带。我不知道我是否应该使用ee方法来提取沿lineString约束的数据(我认为它不是形状?(,或者我是否朝着错误的方向前进,并且应该将缓冲区域导出为GeoTIFF以处理QGIS/ArcGIS中的样带。

import ee
ee.Initialize()
def maskS2clouds(image):
qa = image.select('QA60')
# Bits 10 and 11 are clouds and cirrus, respectively.
cloudBitMask = 1 << 10
cirrusBitMask = 1 << 11
# Both flags should be set to zero, indicating clear conditions.
mask = qa.bitwiseAnd(cloudBitMask).eq(0) 
.And(qa.bitwiseAnd(cirrusBitMask).eq(0))
return image.updateMask(mask).divide(10000)
centre = ee.Geometry.Point(-115.435272, 35.584001,)
centreBuffer = centre.buffer(5000)
dataset = ee.ImageCollection('COPERNICUS/S2_SR') 
.filterDate('2020-07-01', '2020-07-31') 
.filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE',20)) 
.map(maskS2clouds)
visualization = {
'min': 0.0,
'max': 0.3,
'bands': ['B4', 'B3', 'B2'], # Red, Green, BLue.  10m resolution.
}
Map.setCenter(-115.435272, 35.584001, 12)
# Map = geemap.Map(center=[35.584001, -115.435272], zoom=14)
Map.addLayer(dataset.mean(), visualization, 'RGB')
Map.addLayer(centre,
{'color': 'black'},
'Geometry [black]: point');
Map.addLayer(centreBuffer,
{'color': 'red'},
'Result [red]: point.buffer')
transect = ee.Geometry.LineString([[-115.4, 35.584001], [-115.45, 35.584001]]);
Map.addLayer(transect, {'color': 'Green'}, 'transect');
Map

您可能会发现吴秋生的这些资源很有用。

  • https://www.youtube.com/watch?v=0TNXSs6fwrg
  • https://geemap.org/common/#geemap.common.extract_transect
  • https://geemap.org/notebooks/73_transect/

他的Common moduleAPI具有extract_transect(image, line, reducer='mean', n_segments=100, dist_interval=None, scale=None, crs=None, crsTransform=None, tileScale=1.0, to_pandas=False, **kwargs)功能

这将从图像中提取横断面。

从上面的第三个链接来看,Jupyter笔记本是:

import ee
import geemap
from bqplot import pyplot as plt
Map = geemap.Map()
Map.add_basemap("TERRAIN")
Map
image = ee.Image('USGS/SRTMGL1_003')
vis_params = {
'min': 0,
'max': 4000,
'palette': ['006633', 'E5FFCC', '662A00', 'D8D8D8', 'F5F5F5']}
Map.addLayer(image, vis_params, 'SRTM DEM', True, 0.5)
# Use the drawing tool to draw any line on the map.
line = Map.user_roi
if line is None:
line = ee.Geometry.LineString([[-120.223279, 36.314849], [-118.926969, 36.712192], [-117.202217, 36.756215]])
Map.addLayer(line, {}, "ROI")   
Map.centerObject(line)
line.getInfo()
reducer='mean'  # Any ee.Reducer, e.g., mean, median, min, max, stdDev
transect = geemap.extract_transect(image, line, n_segments=100, reducer=reducer, to_pandas=True)
transect
fig = plt.figure()
plt.plot(transect['distance'], transect[reducer])
plt.xlabel('Distance')
plt.ylabel("Elevation")
plt.show()

JavaScript版本是:


/***
* Reduces image values along the given line string geometry using given reducer.
* 
* Samples image values using image native scale, or opt_scale
*/
function reduceImageProfile(image, line, reducer, scale, crs) {
var length = line.length();
var distances = ee.List.sequence(0, length, scale)
var lines = line.cutLines(distances, ee.Number(scale).divide(5)).geometries();
lines = lines.zip(distances).map(function(l) { 
l = ee.List(l)

var geom = ee.Geometry(l.get(0))
var distance = ee.Number(l.get(1))

geom = ee.Geometry.LineString(geom.coordinates())

return ee.Feature(geom, {distance: distance})
})
lines = ee.FeatureCollection(lines)
// reduce image for every segment
var values = image.reduceRegions( {
collection: ee.FeatureCollection(lines), 
reducer: reducer, 
scale: scale, 
crs: crs
})

return values
}

// Define a line across the Olympic Peninsula, USA.
// Import a digital surface model and add latitude and longitude bands.
var elevImg = ee.Image('JAXA/ALOS/AW3D30/V2_2').select('AVE_DSM');
var profile = reduceImageProfile(elevImg, transect, ee.Reducer.mean(), 100)
print(ui.Chart.feature.byFeature(profile, 'distance', ['mean']))

最新更新