Part 1 探索矢量
一、在GEE中使用几何工具创建要素
二、加载现有的特征和特征集合在地球
var roi = table;
//var tiger = ee.FeatureCollection('TIGER/2010/Blocks');
// Add the new feature collection to the map, but do not display.
Map.addLayer(roi, { 'color': 'black' }, 'roi', true);
Map.centerObject(roi, 10);
在Google Earth Engine (GEE) 中，为矢量数据定义样式时，可以通过一个样式字典（JavaScript中的对象字面量）来指定各种样式属性。以下是一些主要的样式参数，可用于控制矢量数据（如Feature或FeatureCollection）的显示方式：
-
color: 定义要素（如线或点）的颜色。颜色可以是一个十六进制的字符串（例如'#ff0000'表示红色），也可以是一个包含RGBA值的数组。 -
fillColor: 仅对多边形要素有效，用于定义多边形内部的填充颜色。颜色的格式同color属性。 -
strokeWidth: 定义线要素的宽度或多边形边界的宽度，单位是像素。 -
strokeColor: 定义线要素或多边形边界的颜色。颜色的格式同color属性。 -
pointSize: 定义点要素的大小，单位是像素。 -
pointShape: 定义点要素的形状。可选值包括'circle'、'square'、'cross'、'x'等。 -
opacity: 定义要素的不透明度，取值范围从0（完全透明）到1（完全不透明）。 -
lineType: 定义线的样式。可选值有'solid'（实线）、'dot'（点线）、'dash'（短划线）、'dash-dot'（点划线）等。 -
lineDash: 自定义线型的数组，表示交替的绘制线段和空白（如[5, 5]表示绘制长度为5的线段，然后空白长度为5，然后重复）。 -
scale: 用于定义要素的显示比例。主要影响点要素和线宽。
这些样式参数使得矢量数据的视觉表示更加灵活和多样化。通过调整这些参数，可以根据具体的应用场景和个人偏好来定制矢量图层的外观。
请注意，不是所有参数都适用于每种类型的矢量要素。例如，fillColor主要用于多边形，而pointSize和pointShape专用于点要素。根据你的数据类型（点、线、多边形）和展示需求选择合适的样式属性。
三、上传矢量到GEE
四、按属性过滤要素集合
4.1 Filter by Geometry of Another F eature
4.2 Filter by F eature (Attribute) Properties
Part 2 栅格/矢量转换
一、栅格转多边形
1.1
// Load raster (elevation) and vector (colombia) datasets.
var elevation = ee.Image('USGS/GMTED2010').rename('elevation');
var colombia = ee.FeatureCollection( 'FAO/GAUL_SIMPLIFIED_500m/2015/level0')
.filter(ee.Filter.equals('ADM0_NAME', 'Colombia'));
print('数据信息', colombia);
// Display elevation image.
Map.addLayer(elevation, {
min: 0,
max: 4000
}, 'Elevation');
// Add the vector to the map.
Map.addLayer(colombia, {
color: 'green', // 线和点的颜色
fillColor: '00000000', // 填充颜色，这里是透明的
strokeWidth: 2 ,// 线的宽度
strokeColor:'red',
lineType:'dash',
opacity:0.4
}, 'colombia', true);
// Initialize image with zeros and define elevation zones.
var zones = ee.Image(0)
.where(elevation.gt(100), 1)
.where(elevation.gt(200), 2)
.where(elevation.gt(500), 3);
print('zones',zones);
// Display zones image.
Map.addLayer(zones, {
min: 0,
max: 3
}, 'zones');
// Mask pixels below sea level (<= 0 m) to retain only land areas.
// Name the band with values 0-3 as 'zone'.
zones = zones.updateMask(elevation.gt(0)).rename('zone');
print('Elevation zones',zones);
Map.addLayer(zones, {
min: 0,
max: 3,
palette: ['white', 'yellow', 'lime', 'green'],
opacity: 0.7
}, 'Elevation zones');
//获取数据集的投影信息
var projection = elevation.projection();
print('projection',projection);
var scale = elevation.projection().nominalScale();
print('elevation.projection().nominalScale()',scale);
//使用reduceToVectors进行矢量化
var elevationVector = zones.reduceToVectors({
geometry: colombia.geometry(),
crs: projection,
scale: 1000,
geometryType: 'polygon',
eightConnected: false,
labelProperty: 'zone',
bestEffort: true,
maxPixels: 1e13,
tileScale: 3 // In case of error.
});
print(elevationVector.limit(10));
//绘制并添加矢量化结果到地图:
var elevationDrawn = elevationVector.draw({
color: 'black',
strokeWidth: 1
});
Map.addLayer(elevationDrawn, {}, 'Elevation zone polygon');
在这个reduceToVectors调用中，各参数的含义如下：
-
geometry: 指定了要进行矢量化操作的区域。这里使用colombia.geometry()表示操作仅限于哥伦比亚国界内的区域。colombia.geometry()从前面定义的哥伦比亚国界的FeatureCollection中提取出几何形状。 -
crs: 指定了输出矢量的坐标参考系统（Coordinate Reference System）。projection变量应该是从某个ee.Image对象（在这个例子中是elevation图像）获取的投影信息，确保矢量化的结果与原始图像在同一投影下进行，以便保持地理精度。 -
scale: 定义了矢量化时使用的空间分辨率（单位是米）。这个参数设置为1000，意味着每个矢量多边形将代表大约1000米×1000米的区域。 -
geometryType: 指定了输出几何的类型。这里设置为'polygon'，表示矢量化的结果将是多边形。这是处理栅格数据时常见的选择，特别是当你想要基于某些条件（如高度区间）划分区域时。 -
eightConnected: 控制矢量化算法的连接性。设置为false意味着使用四连通算法（只考虑像素的上下左右邻域），而true会使用八连通算法（考虑像素的八个邻域，即包括对角线方向）。这个选择影响矢量化的精细度和形状。 -
labelProperty: 为结果矢量图层的每个要素指定一个属性标签。这里用'zone'属性来标记每个多边形，这个属性应该是在zones数据中已定义，用于表示每个区域的分类或标识。 -
bestEffort: 当设置为true时，GEE会尝试在给定的资源限制下（如maxPixels）尽可能地完成矢量化任务。这可能会导致分辨率的自动调整，以便在内存和计算时间的限制内完成操作。 -
maxPixels: 定义了操作可以使用的最大像素数量。这个限制有助于避免因处理大量数据而导致的内存溢出或长时间计算。在这个例子中，设置为1e13（即10万亿像素），是一个相对宽松的限制，旨在确保即使是大范围的矢量化操作也能完成。 -
tileScale: 用于优化处理的参数。增加tileScale的值可以减小服务器处理每个瓦片的区域大小，有助于减轻内存压力，特别是在处理大范围或复杂数据时。这里的设置为3，意味着将使用更小的瓦片来分布计算任务。
1.2 .focalMode() 方法
对栅格数据集应用一个平滑操作。
使用 .focalMode() 方法，这是一种邻域分析技术，用于每个像素计算其邻域（即周围像素）中最常见的值（众数），并以此作为该像素的新值。
这种方法通常用于数据平滑或噪声减少，特别是在分类地图上，可以帮助消除“椒盐”噪声（即随机散布的、孤立的误分类像素）。
举例：
.focalMode(4, 'square');
-
4: 这个参数指定了操作的半径为4个像素。这意味着对于每个像素，其周围4个像素的范围内的邻域将被用来计算众数。较大的半径值将考虑更宽的邻域，可能导致更平滑的结果，但也可能导致更多的细节丢失。
-
'square': 这个字符串参数指定了邻域的形状。在这个例子中，
'square'意味着使用正方形的邻域。GEE 也支持圆形邻域（'circle'），其中邻域的形状会有所不同，这影响了哪些周围像素被包含在计算中。
//续上面代码
var zonesSmooth = zones.focalMode(4, 'square'); zonesSmooth = zonesSmooth.reproject(projection.atScale(scale)); print('zonesSmooth',zonesSmooth); Map.addLayer(zonesSmooth, { min: 1, max: 3, palette: ['yellow', 'lime', 'green'], opacity: 0.7 }, 'Elevation zones (smooth)' ); var elevationVectorSmooth = zonesSmooth.reduceToVectors({ geometry: colombia.geometry(), crs: projection, scale: scale, geometryType: 'polygon', eightConnected: false, labelProperty: 'zone', bestEffort: true, maxPixels: 1e13, tileScale: 3 }); print('elevationVectorSmooth',elevationVectorSmooth); var smoothDrawn = elevationVectorSmooth.draw({ color: 'black', strokeWidth: 1 }); Map.addLayer(smoothDrawn, {}, 'Elevation zone polygon (smooth)');
二、栅格转点
2.1 每一个像素转点 .sample
//=========================栅格转点==================================
var geometry = ee.Geometry.Polygon([
[-89.553, -0.929],
[-89.436, -0.929],
[-89.436, -0.866],
[-89.553, -0.866],
[-89.553, -0.929]
]);
// To zoom into the area, un-comment and run below 取消注释并运行
Map.centerObject(geometry,12);
Map.addLayer(geometry, {}, 'Areas to extract points');
var elevationSamples = elevation.sample({
region: geometry,
projection: projection,
scale: scale,
geometries: true,
});
Map.addLayer(elevationSamples, {}, 'Points extracted');
// Add three properties to the output table:
// 'Elevation', 'Longitude', and 'Latitude'.
elevationSamples = elevationSamples.map(function(feature) {
var geom = feature.geometry().coordinates();
return ee.Feature(null, {
'Elevation': ee.Number(feature.get('elevation')),
'Long': ee.Number(geom.get(0)),
'Lat': ee.Number(geom.get(1))
});
});
// Export as CSV.
Export.table.toDrive({
collection: elevationSamples,
description: 'extracted_points',
fileFormat: 'CSV'
});
2.2 每个区域门选取部分样本点 .stratifiedSample
代码示例：
var sample = zones.stratifiedSample({
numPoints: 100, // 每个类别抽取100个点
classBand: 'landcover', // 使用名为'landcover'的波段进行分层
region: regionGeometry, // 指定抽样区域
scale: 30, // 指定空间分辨率为30米
seed: 0, // 指定随机种子
geometries: true // 返回样本点的地理位置信息
});
调用 stratifiedSample 方法时，可以指定多个参数，包括但不限于：
numPoints：每个分层中要抽取的样本点数量。classBand：指定用于分层的栅格图层的波段名称，或是矢量数据的属性字段，根据这一信息进行分层。region：指定抽样的地理区域，通常是一个ee.Geometry对象。scale：指定分析时的空间分辨率。seed：用于随机抽样过程的种子，确保结果可复现。geometries：一个布尔值，指定是否返回样本的地理位置信息（作为点的几何数据）。
var elevationSamplesStratified = zones.stratifiedSample({
numPoints: 10,
classBand: 'zone',
region: geometry,
scale: scale,
projection: projection,
geometries: true
});
Map.addLayer(elevationSamplesStratified, {}, 'Stratified samples');
三、一个更复杂的例子
我们将使用两个全局数据集，一个用于表示栅格格式，另一个用于表示矢量:
The Global Forest Change (GFC) :描述2001年至今全球树木覆盖及其变化的栅格数据集。
The World Database on Protected Areas (WDPA) :全球保护区矢量数据库。
3.1 GFC
lossyear band: 森林损失的分类栅格(1-20，对应2001-2020年期间的森林砍伐)，0表示没有变化.
3.2 WDPA
NAME ':每个保护区的名称
' WDPA_PID ':每个保护区的唯一数字ID
3.3 ee.Image().byte().paint() 创建轮廓图像
// 假设有一个名为myFeatures的ee.FeatureCollection
var paintedImage = ee.Image().byte().paint({
featureCollection: myFeatures,
color: 1, // 所有特征绘制为值1（通常为白色）
width: 2 // 线宽为2像素
});
// 将绘制的图像添加到地图上
Map.addLayer(paintedImage, {palette: ['000000', 'FFFFFF']}, 'Features Outline');
在Google Earth Engine (GEE) 中，ee.Image().byte().paint()函数是一种强大的工具，用于在图像上绘制特征（如线或多边形）。这个方法通常用于将矢量数据（如ee.Feature或ee.FeatureCollection）的轮廓或边界绘制到一个图像上，从而可以在GEE的地图上可视化这些边界或轮廓。下面是paint方法的参数详解：
-
featureCollection: 需要被绘制到图像上的特征集合（ee.FeatureCollection）或单个特征（ee.Feature）。这是你想要可视化的轮廓或边界数据。 -
color: 用于绘制特征的颜色。这个参数可以是一个常量，表示所有特征使用相同的颜色绘制；也可以是一个字符串，指定一个属性名称，根据每个特征的该属性值来动态设置颜色。 -
width: 线条宽度，仅对线特征或多边形特征的边界有效。这个参数定义了绘制线条的像素宽度。
实际上，你需要先调用ee.Image().byte()来创建一个空的单波段字节图像（值范围在0到255之间），然后使用.paint()方法在其上绘制特征。这样做的好处是可以生成一个轻量级的图像，适用于将矢量数据覆盖到栅格图层上。
3.4 地图的底图样式
Map.setOptions('ROADMAP'); // 设置为街道地图视图
Map.setOptions('SATELLITE'); // 设置为卫星图像视图
Map.setOptions('HYBRID'); // 设置为卫星和街道混合视图
Map.setOptions('TERRAIN'); // 设置为地形视图
在Google Earth Engine (GEE) 的JavaScript API中，Map.setOptions()方法用于设置地图的显示选项，包括地图的底图类型。'SATELLITE'是这些选项之一，提供了卫星图像的视图。除了卫星图像外，还有其他几种预设的底图类型可以选择，常见的包括：
'ROADMAP': 显示谷歌地图的标准街道地图视图。'SATELLITE': 显示谷歌地图的卫星图像。'HYBRID': 显示卫星视图叠加了主要街道和地名的图层，结合了'SATELLITE'和'ROADMAP'的特点。'TERRAIN': 显示地形信息的地图，包括高度阴影和自然地貌标识（如山脉、河流等）。
3.5
//=========================== Read input data. =====================================
// Note: these datasets are periodically updated.
// Consider searching the Data Catalog for newer versions.
var gfc = ee.Image('UMD/hansen/global_forest_change_2020_v1_8');
var wdpa = ee.FeatureCollection('WCMC/WDPA/current/polygons');
// Print assets to show available layers and properties.
print(gfc);
print(wdpa.limit(10)); // Show first 10 records.
// Add the vector to the map.
Map.addLayer(wdpa, {
color: 'green', // 线和点的颜色
fillColor: '00000000', // 填充颜色，这里是透明的
strokeWidth: 2 ,// 线的宽度
strokeColor:'red',
lineType:'dash',
opacity:0.4
}, 'wdpa', true); //=============================section 2============================================
// Display deforestation.
var deforestation = gfc.select('lossyear');
Map.addLayer(deforestation, {
min: 1,
max: 20,
palette: ['yellow', 'orange', 'red']
}, 'Deforestation raster');
// Display WDPA data.
//var protectedArea = wdpa.filter(ee.Filter.equals('NAME', 'La Paya'));
var protectedArea = wdpa.filter(ee.Filter.equals('NAME', 'Fanjingshan'));
// Display protected area as an outline (see F5.3 for paint()).
var protectedAreaOutline = ee.Image().byte().paint({
featureCollection: protectedArea,
color: 1,
width: 3
});
Map.addLayer(protectedAreaOutline, {
palette: 'white'
}, 'Protected area');
// Set up map display.
Map.centerObject(protectedArea,11);
Map.setOptions('SATELLITE');//设置地图底图样式,地图的底图样式被设置为卫星图像
//===========================section 3 Convert from a deforestation raster to vector. ==========
var deforestationVector = deforestation.reduceToVectors({
scale: deforestation.projection().nominalScale(),
geometry: protectedArea.geometry(),
labelProperty: 'lossyear', // Label polygons with a change year.
maxPixels: 1e13
});
// Count the number of individual change events
print('Number of change events:', deforestationVector.size());
// Display deforestation polygons. Color outline by change year.
var deforestationVectorOutline = ee.Image().byte().paint({
featureCollection: deforestationVector,
color: 'lossyear',
width: 1
});
Map.addLayer(deforestationVectorOutline, {
palette: ['yellow', 'orange', 'red'],
min: 1,
max: 20
}, 'Deforestation vector');
//绘制20年每年的一个毁林时间数量
//==========Plot of the number of deforestation events in La Paya for the years 2001–2020 ========
var chart = ui.Chart.feature
.histogram({
features: deforestationVector,
property: 'lossyear'
})
.setOptions({
hAxis: {
title: 'Year'
},
vAxis: {
title: 'Number of deforestation events'
},
legend: {
position: 'none'
}
});
print(chart);
//================section 5 计算每个多边形的质心 Generate deforestation point locations. ==================================
var deforestationCentroids =
deforestationVector.map(function(feat) {
return feat.centroid();
});
Map.addLayer(deforestationCentroids, {
color: 'darkblue'
}, 'Deforestation centroids');
//================section 6 筛选毁林面积大于3公顷的图斑====================================
// 给毁林变化矢量添加一个“面积”的字段，并计算每个多边形的面积
deforestationVector =
deforestationVector.map(function(feat) {
return feat.set('area', feat.geometry().area({
maxError: 10
}).divide(10000)); // Convert m^2 to hectare.
});
//代码分解
// Filter the deforestation FeatureCollection for only large-scale (>10 ha) changes
var deforestationLarge =
deforestationVector.filter(ee.Filter.gt(
'area', 3));
// Display deforestation area outline by year.
var deforestationLargeOutline = ee.Image().byte().paint({
featureCollection: deforestationLarge,
color: 'lossyear',
width: 1
});
Map.addLayer(deforestationLargeOutline, {
palette: ['yellow', 'orange', 'red'],
min: 1,
max: 20
}, 'Deforestation (>3 ha)');
-
对
deforestationVector中的每个特征应用函数:deforestationVector = deforestationVector.map(function(feat) { ... });deforestationVector是一个ee.FeatureCollection对象，包含了多个地理特征（可能是表示森林砍伐区域的多边形）。.map()方法对deforestationVector中的每个特征执行给定的函数。这里，函数接受一个特征feat作为输入。
-
计算每个特征的面积:
var area = feat.geometry().area({ maxError: 10 }).divide(10000);feat.geometry().area({maxError: 10})计算特征的面积（单位为平方米）。maxError参数用于指定计算面积时的最大允许误差，单位为米。较小的maxError值可以提高面积计算的精度，但可能增加计算成本。.divide(10000)将面积从平方米转换为公顷（1公顷=10,000平方米）。这是因为在很多环境和农业应用中，公顷是更常用的面积单位。
-
更新特征属性:
return feat.set('area', area);.set('area', area)方法用于在特征的属性中添加（或更新）一个名为'area'的属性，值为计算得到的面积（单位为公顷）。
四、栅格属性 至 矢量字段
有时我们想从栅格中提取信息，以便包含在转换后的矢量数据集中。
// Load required datasets.
var gfc = ee.Image('UMD/hansen/global_forest_change_2020_v1_8');
var wdpa = ee.FeatureCollection('WCMC/WDPA/current/polygons');
// Display deforestation.
var deforestation = gfc.select('lossyear');
Map.addLayer(deforestation, {
min: 1,
max: 20,
palette: ['yellow', 'orange', 'red']
}, 'Deforestation raster');
// Select protected areas in the Colombian Amazon.
var amazonianProtectedAreas = [
'Cordillera de los Picachos', 'La Paya', 'Nukak',
'Serrania de Chiribiquete',
'Sierra de la Macarena', 'Tinigua'
];
var wdpaSubset = wdpa.filter(ee.Filter.inList('NAME', amazonianProtectedAreas));
// Display protected areas as an outline.
var protectedAreasOutline = ee.Image().byte().paint({
featureCollection: wdpaSubset,
color: 1,
width: 1
});
Map.addLayer(protectedAreasOutline, {
palette: 'white'
}, 'Amazonian protected areas');
// Set up map display.
Map.centerObject(wdpaSubset);
Map.setOptions('SATELLITE');
var scale = deforestation.projection().nominalScale();
// Use 'reduceRegions' to sum together pixel areas in each protected area.
wdpaSubset = deforestation.gte(1)
.multiply(ee.Image.pixelArea().divide(10000)).reduceRegions
({
collection: wdpaSubset,
reducer: ee.Reducer.sum().setOutputs([
'deforestation_area']),
scale: scale
});
print(wdpaSubset); // Note the new 'deforestation_area' property.
// Normalize by area.
wdpaSubset = wdpaSubset.map(
function(feat) {
return feat.set('deforestation_rate',
ee.Number(feat.get('deforestation_area'))
.divide(feat.area().divide(10000)) // m2 to ha
.divide(20) // number of years
.multiply(100)); // to percentage points
});
// Print to identify rates of change per protected area.
// Which has the fastest rate of loss?
print(wdpaSubset.reduceColumns({
reducer: ee.Reducer.toList().repeat(2),
selectors: ['NAME', 'deforestation_rate']
}));
代码解释：
- 表达式
deforestation.gte(1)创建了一个二值掩膜，其中砍伐水平大于或等于 1 的像素被标记为 1（真），其余被标记为 0（假）。这一步骤旨在隔离感兴趣的区域——即被砍伐的区域。 ee.Image.pixelArea().divide(10000)计算每个像素的面积，单位为平方公里（假设原始像素面积单位为平方米，因此除以 10,000）。然后，将此面积计算乘以二值砍伐掩膜，得到的图像中每个像素的值代表该像素中的砍伐面积。
在保护区内聚合：然后使用 reduceRegions 函数在指定区域内聚合这些像素区域值。函数应用如下：
collection: wdpaSubset指定要执行聚合的区域集合（在此上下文中为保护区）。此集合中的每个区域都将计算其边界内砍伐区域的总和。reducer: ee.Reducer.sum().setOutputs(['deforestation_area'])定义聚合方法。这里，它对每个区域内的砍伐像素面积进行求和。集合中每个区域的输出属性将命名为deforestation_area。scale: scale设置执行计算的空间分辨率。scale的值应该与输入数据的分辨率相匹配，以获得准确的结果。-
-
归一化砍伐面积：通过
wdpaSubset.map函数，代码对wdpaSubset集合中的每个保护区（feat）进行遍历，计算其砍伐率。砍伐率是基于保护区的砍伐面积与总面积的比例，并考虑了一定时间跨度（此处为20年）。 -
计算步骤：
- 使用
ee.Number(feat.get('deforestation_area'))获取每个保护区的砍伐面积。 - 通过
feat.area().divide(10000)计算保护区的总面积，单位从平方米转换为公顷（1公顷=10,000平方米）。 - 通过
.divide(20)将砍伐面积分摊到20年中，得到年均砍伐面积。 - 通过
.multiply(100)将结果转换为百分比形式，以便更直观地理解砍伐率。
- 使用
-
输出：通过使用
print语句和wdpaSubset.reduceColumns函数，代码计算并打印出每个保护区的名字和对应的森林砍伐率。ee.Reducer.toList().repeat(2)确保将保护区的名称和砍伐率作为两列数据输出。.repeat(2)：这个方法用于指定简化器操作的重复次数 -
识别砍伐率最高的保护区：通过比较输出列表中的
deforestation_rate，可以识别出砍伐率最高的保护区。这一步骤虽然没有直接在代码中体现，但可以通过输出结果的分析来完成。
-
五、矢量转栅格
5.1 多边形生成掩膜 Polygons to a Mask
.reduceToImage（）用法
用法
reduceToImage() 方法通常用于将地理特征集合（如一系列具有地理属性的区域或点）中的数值属性转换成一个栅格图像，其中图像中的每个像素值代表了对应地理区域中特定属性的值。
参数
-
属性列表 (
properties)：这是一个字符串数组，指定了要从每个特征中提取哪些属性以生成图像。这些属性应该是数值类型，因为它们需要转换成图像中的像素值。 -
简化器 (
reducer)：这是指定如何从每个特征的属性列表中生成单个像素值的Reducer对象。例如，如果多个特征覆盖了图像的同一个像素位置，简化器将决定如何合并这些特征的属性值以产生最终的像素值。常用的简化器有ee.Reducer.first()（取第一个特征的值）、ee.Reducer.mean()（计算平均值）、ee.Reducer.sum()（求和）等。
// 假设有一个FeatureCollection，每个特征都有一个名为'height'的属性
var featureCollection = ee.FeatureCollection([...]);
// 使用.reduceToImage()将'height'属性转换成图像
var heightImage = featureCollection.reduceToImage({
properties: ['height'], // 要转换的属性列表
reducer: ee.Reducer.first() // 使用的简化器
});
// Load required datasets.
var gfc = ee.Image('UMD/hansen/global_forest_change_2020_v1_8');
var wdpa = ee.FeatureCollection('WCMC/WDPA/current/polygons');
print('wdpa', wdpa.limit(2));
// Get deforestation.
var deforestation = gfc.select('lossyear');
// Generate a new property called 'protected' to apply to the output mask.
var wdpa = wdpa.map(function(feat) {
return feat.set('protected', 1);
});
print('wdpa1', wdpa.limit(2));
// Rasterize using the new property.
// unmask() sets areas outside protected area polygons to 0.
var wdpaMask = wdpa.reduceToImage(['protected'],
ee.Reducer.first())
.unmask();
print('wdpaMask', wdpaMask);
// Center on Colombia.
Map.setCenter(-75, 3, 6);
// Display on map.
Map.addLayer(wdpaMask, {
min: 0,
max: 1
}, 'Protected areas (mask)');
// Produce an image with unique ID of protected areas.
//reduceToImage还可以保留输入多边形的数值属性
var wdpaId = wdpa.reduceToImage(['WDPAID'],
ee.Reducer.first());
print('wdpaMask-wdpaId', wdpaMask);
Map.addLayer(wdpaId, {
min: 1,
max: 100000
}, 'Protected area ID');
代码解释：
-
reduceToImage方法：这个方法将一个FeatureCollection（如保护区的集合）转换为一个Image，方法中指定的属性值（此例中为 'protected'）被用作像素值。这里使用的ee.Reducer.first()确保如果有多个特征重叠，只选择第一个特征的属性值。 -
.unmask()：该方法用于处理图像中的无数据（masked）像素。在这个上下文中，它可能被用来确保图像在没有数据的地方不会被掩膜，但具体行为取决于unmask()是否接收了参数（在这里没有给出参数，所以它将简单地显示所有没有数据的像素）。
5.2 一个更复杂的例子
// Load required datasets.
var gfc = ee.Image('UMD/hansen/global_forest_change_2020_v1_8');
var wdpa = ee.FeatureCollection('WCMC/WDPA/current/polygons');
// Select a single protected area.
var protectedArea = wdpa.filter(ee.Filter.equals('NAME', 'La Paya'));
var protectedAreaOutline = ee.Image().byte().paint({
featureCollection: protectedArea,
color: 'red',
width: 1
});
Map.addLayer( protectedAreaOutline,{palette: 'red'},'protectedArea');
// Maximum distance in meters is set in the brackets.
var distance = protectedArea.distance(1000000);
Map.addLayer(distance, {
min: 0,
max: 20000,
palette: ['white', 'grey', 'black'],
opacity: 0.6
}, 'Distance');
Map.centerObject(protectedArea);
// Produce a raster of inside/outside the protected area.
var protectedAreaRaster = protectedArea.map(function(feat)
{
return feat.set('protected', 1);
}).reduceToImage(['protected'], ee.Reducer.first());
Map.addLayer(distance.updateMask(protectedAreaRaster), {
min: 0,
max: 20000
}, 'Distance inside protected area');
Map.addLayer(distance.updateMask(protectedAreaRaster.unmask
()
.not()), {
min: 0,
max: 20000
}, 'Distance outside protected area');
var distanceZones = ee.Image(0)
.where(distance.gt(0), 1)
.where(distance.gt(1000), 2)
.where(distance.gt(3000), 3)
.updateMask(distance.lte(5000));
Map.addLayer(distanceZones, {}, 'Distance zones');
var deforestation = gfc.select('loss');
var deforestation1km = deforestation.updateMask(distanceZones.eq(1));
var deforestation3km = deforestation.updateMask(distanceZones.lte(2));
var deforestation5km = deforestation.updateMask(distanceZones.lte(3));
Map.addLayer(deforestation1km, {
min: 0,
max: 1
}, 'Deforestation within a 1km buffer');
Map.addLayer(deforestation3km, {
min: 0,
max: 1,
opacity: 0.5
}, 'Deforestation within a 3km buffer');
Map.addLayer(deforestation5km, {
min: 0,
max: 1,
opacity: 0.5
}, 'Deforestation within a 5km buffer');
var deforestation1kmOutside = deforestation1km
.updateMask(protectedAreaRaster.unmask().not());
// Get the value of each pixel in square meters
// and divide by 10000 to convert to hectares.
var deforestation1kmOutsideArea =
deforestation1kmOutside.eq(1)
.multiply(ee.Image.pixelArea()).divide(10000);
// We need to set a larger geometry than the protected area
// for the geometry parameter in reduceRegion().
var deforestationEstimate = deforestation1kmOutsideArea
.reduceRegion({
reducer: ee.Reducer.sum(),
geometry: protectedArea.geometry().buffer(1000),
scale: deforestation.projection().nominalScale()
});
print('Deforestation within a 1km buffer outside the protected area (ha)', deforestationEstimate);
代码解释：
var distance = protectedArea.distance(1000000);
-
.distance()方法：这个方法计算从protectedArea到每个像素位置的距离。它返回一个Image，其中的每个像素值表示该像素到最近的protectedArea边界的距离。 -
1000000参数：这个数字参数（以米为单位）是.distance()方法的最大距离阈值。这意味着计算将限制在距离保护区边界最多 1,000,000 米（或 1000 公里）的范围内。距离超过这个阈值的像素将被赋予该阈值作为距离。这个参数有助于优化计算过程，防止在远离感兴趣区域的地方进行不必要的计算。 -
1 // Produce a raster of inside/outside the protected area. 2 var protectedAreaRaster = protectedArea.map(function(feat) 3 { 4 return feat.set('protected', 1); 5 }).reduceToImage(['protected'], ee.Reducer.first()); 6 7 Map.addLayer(distance.updateMask(protectedAreaRaster), { 8 min: 0, 9 max: 20000 10 }, 'Distance inside protected area'); 11 12 Map.addLayer(distance.updateMask(protectedAreaRaster.unmask 13 () 14 .not()), { 15 min: 0, 16 max: 20000 17 }, 'Distance outside protected area');这段代码使用 Google Earth Engine (GEE) 创建了两个栅格（raster）图层：一个表示在保护区内部的距离，另一个表示在保护区外部的距离，并将这两个图层添加到地图上。
-
创建保护区内外的栅格图层
-
生成表示保护区的栅格图像：首先，使用
protectedArea.map函数遍历保护区集合中的每个特征（feat），给每个特征设置一个名为'protected'的属性，值为 1。这意味着所有保护区特征都标记为“受保护”。然后，通过.reduceToImage(['protected'], ee.Reducer.first())方法将这些特征集合转换为一个栅格图像（protectedAreaRaster），在这个图像中，保护区内的像素值为 1，而保护区外的像素值为无数据状态。 -
在保护区内显示距离：
distance.updateMask(protectedAreaRaster)方法用于更新distance图像的掩膜。这个掩膜基于protectedAreaRaster图像，只有当protectedAreaRaster的像素值为 1（即在保护区内）时，distance图像的对应像素才会被保留。结果是只有保护区内的距离被显示。这个处理后的图像通过Map.addLayer方法添加到地图上，图层名为 'Distance inside protected area'，显示的是保护区内各点到保护区边界的距离。 -
在保护区外显示距离：通过
protectedAreaRaster.unmask().not()创建一个新的掩膜，这个掩膜在保护区外的地方为真（即保护区外的像素值变为 1）。然后，使用distance.updateMask方法应用这个掩膜，仅显示保护区外的距离。这个处理后的图像也通过Map.addLayer方法添加到地图上，图层名为 'Distance outside protected area'，显示的是保护区外各点到保护区最近边界的距离。 -
-
inside outside
1 var distanceZones = ee.Image(0)
2 .where(distance.gt(0), 1)
3 .where(distance.gt(1000), 2)
4 .where(distance.gt(3000), 3)
5 .updateMask(distance.lte(5000));
6
7
8 Map.addLayer(distanceZones, {}, 'Distance zones');
-
初始化图像：
ee.Image(0)创建了一个初始值全为 0 的图像。这个图像将作为基础，通过后续操作修改其值来表示不同的距离分区。 -
创建距离分区：
.where(distance.gt(0), 1)：这个操作检查之前计算的distance图像中的每个像素值。如果像素值大于 0（意味着位于保护区边界之外的任意位置），则将该像素的值设置为 1。这个步骤定义了最内层的距离分区。.where(distance.gt(1000), 2)：然后，对于距离大于 1000 米的区域，将像素值更新为 2。这创建了第二个距离分区。.where(distance.gt(3000), 3)：接着，对于距离大于 3000 米的区域，将像素值更新为 3。这定义了第三个距离分区。
-
应用掩膜以限制分区的显示范围：
.updateMask(distance.lte(5000))这行代码应用了一个掩膜，仅显示距离小于或等于 5000 米的区域。距离大于 5000 米的区域将不会在最终的图层中显示。这个步骤限制了距离分区图层的显示范围，确保只有距离在 5000 米以内的分区被展示。 -
添加图层到地图：
Map.addLayer(distanceZones, {}, 'Distance zones')这行代码将处理后的distanceZones图像作为一个新图层添加到地图上，图层名称为 'Distance zones'。这个图层展示了不同的距离分区，分别由不同的像素值（1、2、3）表示，覆盖了从保护区边界向外至多 5000 米的区域。
- Distance zones
-
1 var deforestation = gfc.select('loss'); 2 var deforestation1km = deforestation.updateMask(distanceZones.eq(1)); 3 var deforestation3km = deforestation.updateMask(distanceZones.lte(2)); 4 var deforestation5km = deforestation.updateMask(distanceZones.lte(3)); 5 6 7 8 Map.addLayer(deforestation1km, { 9 min: 0, 10 max: 1 11 }, 'Deforestation within a 1km buffer'); 12 13 Map.addLayer(deforestation3km, { 14 min: 0, 15 max: 1, 16 opacity: 0.5 17 }, 'Deforestation within a 3km buffer'); 18 19 Map.addLayer(deforestation5km, { 20 min: 0, 21 max: 1, 22 opacity: 0.5 23 }, 'Deforestation within a 5km buffer'); 24 25 26 var deforestation1kmOutside = deforestation1km 27 .updateMask(protectedAreaRaster.unmask().not()); 28 29 // Get the value of each pixel in square meters 30 // and divide by 10000 to convert to hectares. 31 var deforestation1kmOutsideArea = 32 deforestation1kmOutside.eq(1) 33 .multiply(ee.Image.pixelArea()).divide(10000); 34 35 36 // We need to set a larger geometry than the protected area 37 // for the geometry parameter in reduceRegion(). 38 var deforestationEstimate = deforestation1kmOutsideArea 39 .reduceRegion({ 40 reducer: ee.Reducer.sum(), 41 geometry: protectedArea.geometry().buffer(1000), 42 scale: deforestation.projection().nominalScale() 43 }); 44 print('Deforestation within a 1km buffer outside the protected area (ha)', deforestationEstimate);分析不同距离范围内的森林砍伐
-
选择砍伐数据：
var deforestation = gfc.select('loss');从 Global Forest Change (GFC) 数据集中选择了森林砍伐（'loss'）图层。 -
应用距离分区掩膜：使用
distanceZones图层中定义的距离分区（前面步骤中创建），通过.updateMask()方法为不同的距离范围创建砍伐图层：deforestation1km：只显示距离保护区边界1公里内的砍伐情况。deforestation3km：显示距离保护区边界3公里内的砍伐情况。deforestation5km：显示距离保护区边界5公里内的砍伐情况。
-
在地图上添加图层：使用
Map.addLayer方法将这些砍伐图层添加到地图上，通过不同的名称标识每个图层（'Deforestation within a 1km buffer' 等），并设置透明度以便图层重叠时能够区分。
计算保护区外1公里缓冲区内的砍伐面积
-
保护区外1公里砍伐情况：首先，
deforestation1kmOutside通过.updateMask(protectedAreaRaster.unmask().not())从deforestation1km中排除保护区内的砍伐，仅保留保护区外1公里范围内的砍伐数据。 -
计算砍伐面积：
deforestation1kmOutsideArea使用.eq(1)确定砍伐像素，然后乘以每个像素的面积（.multiply(ee.Image.pixelArea())），并除以10000将结果从平方米转换为公顷。 -
估算总砍伐面积：
deforestationEstimate使用.reduceRegion()方法计算deforestation1kmOutsideArea图像的总砍伐面积（以公顷为单位）。为了确保覆盖整个感兴趣的区域，使用了保护区加上1公里缓冲区的几何形状作为geometry参数，并将缩放比例（scale）设置为砍伐数据的原始分辨率。
输出结果
最后，使用
print语句输出保护区外1公里缓冲区内的总砍伐面积（单位为公顷），提供了对保护区周边砍伐影响的量化评估。 -
-
标签：C22,Map,area,ee,矢量,像素,deforestation,栅格,var From： https://www.cnblogs.com/bltstop/p/18094957