??本文介紹在谷歌地球引擎（Google Earth Engine，GEE）中，按照給定的地表分類數(shù)據(jù)，對(duì)每一種不同的地物類型，分別加以全球范圍內(nèi)隨機(jī)抽樣點(diǎn)自動(dòng)批量選取的方法。
??本文是谷歌地球引擎（Google Earth Engine，GEE）系列教學(xué)文章的第22篇，更多GEE文章請(qǐng)參考專欄：GEE學(xué)習(xí)與應(yīng)用。
??首先，來具體明確一下本文的需求。我們現(xiàn)在有一個(gè)MODIS的地表分類數(shù)據(jù)產(chǎn)品MCD12Q1，其LC_Type3波段可以提供全球陸地范圍內(nèi)，每一年的地表植被分類數(shù)據(jù)，空間分辨率為500 m。如下圖所示，就是其對(duì)于地表植被類型的分類；紫色框內(nèi)的8種地表類型，就是8種不同的植被類型。
??我們希望實(shí)現(xiàn)的是，在全球陸地范圍內(nèi)，針對(duì)上述8種植被類型，對(duì)于每一種植被類型，分別篩選3000個(gè)左右的隨機(jī)采樣點(diǎn)；因?yàn)橐还彩?code>8種植被類型，所以相當(dāng)于我們最終希望得到的全部采樣點(diǎn)個(gè)數(shù)大約為8 * 3000 = 24000。其中，每一種植被類型對(duì)應(yīng)的3000個(gè)左右的采樣點(diǎn)，都單獨(dú)作為一個(gè)FeatureCollection，然后保存到Assets中。
??明確了需求，即可開始代碼的撰寫。本文用到的代碼如下。
var modis_type = ee.Image("projects/ee-ucanuse6/assets/vegetation_type"),
rectangle =
/* color: #98ff00 */
/* shown: false */
/* displayProperties: [
{
"type": "rectangle"
}
] */
ee.Geometry.Polygon(
[[[-164.97163926414936, 74.9796307938761],
[-164.97163926414936, -57.6624255131542],
[176.74711073585064, -57.6624255131542],
[176.74711073585064, 74.9796307938761]]], null, false);
// Map.addLayer(rectangle);
var landcoverTypes = [1, 2, 3, 4, 5, 6, 7, 8];
var modis_type_vis = {
min: 1.0,
max: 17.0,
palette: [
'05450a', '086a10', '54a708', '78d203', '009900', 'c6b044', 'dcd159',
'dade48', 'fbff13', 'b6ff05', '27ff87', 'c24f44', 'a5a5a5', 'ff6d4c',
'69fff8', 'f9ffa4', '1c0dff'
],
};
Map.addLayer(modis_type, modis_type_vis, "MODIS_Type");
landcoverTypes.forEach(function(type) {
var type_mask = modis_type.eq(type);
var type_image = modis_type.updateMask(type_mask);
var type_point = type_image.sample({
scale: 500,
region: rectangle,
numPixels: 140000,
seed: 7,
dropNulls: true,
geometries: true
});
print(type_point);
Map.addLayer(type_point, {}, "Point");
Export.table.toAsset({
collection: type_point,
description: 'point_' + type,
assetId: 'point_' + type,
});
});
??其中，// Map.addLayer(rectangle);這句代碼以上的內(nèi)容，是獲取我這里處理好的MODIS地表分類數(shù)據(jù)，以及我手動(dòng)繪制的一個(gè)包含了全球主要陸地部分的Polygon。這一部分的代碼，大家復(fù)制之后，可以轉(zhuǎn)換成Imports的格式，如下圖紅色框內(nèi)所示。
??代碼接下來的部分的含義如下。
??首先，var landcoverTypes = [1, 2, 3, 4, 5, 6, 7, 8];定義了一個(gè)包含地物類型的數(shù)組landcoverTypes，這里的值1到8分別對(duì)應(yīng)著我們所需的8種土地植被覆蓋類別。
??隨后，var modis_type_vis = { ... };這部分代碼定義了MODIS地物類型數(shù)據(jù)的可視化參數(shù)，其中包含了最小值、最大值和調(diào)色板信息，用于根據(jù)類型值將地物類型數(shù)據(jù)映射到不同顏色，從而可以更好地可視化MODIS地表分類數(shù)據(jù)產(chǎn)品。隨后，Map.addLayer(modis_type, modis_type_vis, "MODIS_Type");將MODIS地表分類數(shù)據(jù)產(chǎn)品添加到地圖中，并使用上述定義的可視化參數(shù)進(jìn)行顯示。
??接下來，landcoverTypes.forEach(function(type) { ... });表示對(duì)landcoverTypes數(shù)組中的每個(gè)值進(jìn)行循環(huán)，也就是對(duì)8種地表類型數(shù)據(jù)分別加以循環(huán)。
??隨后，var type_mask = modis_type.eq(type);表示根據(jù)當(dāng)前循環(huán)的類型值，創(chuàng)建一個(gè)布爾類型的掩膜type_mask，用于篩選出與該類型相符的像素。其次，var type_image = modis_type.updateMask(type_mask);表示使用掩膜type_mask對(duì)地表類型數(shù)據(jù)進(jìn)行屏蔽（掩膜），保留與當(dāng)前類型相符的像素，并將其他像素設(shè)置為null。
??其次，var type_point = type_image.sample({ ... });表示對(duì)掩膜后的地表類型數(shù)據(jù)進(jìn)行抽樣，生成樣本點(diǎn)；在這里，我們使用給定的參數(shù)（比例尺、區(qū)域、像素?cái)?shù)等）指定抽樣的方式。
??最后，使用print(type_point);打印該地表類型的樣本點(diǎn)信息，以查看抽樣結(jié)果；通過Map.addLayer(type_point, {}, "Point");將樣本點(diǎn)圖層添加到地圖中進(jìn)行可視化。
??完成上述操作后，通過Export.table.toAsset({ ... });將樣本點(diǎn)導(dǎo)出，將其存儲(chǔ)為GEE中的Asset。
??運(yùn)行上述代碼，首先我們看一下Console一欄的結(jié)果——也就是print(type_point);這句代碼，打印出的8種地表類型對(duì)應(yīng)的抽樣點(diǎn)信息；如下圖所示。
??在這里，首先可以看到，由于部分地表類型數(shù)據(jù)對(duì)應(yīng)的抽樣點(diǎn)個(gè)數(shù)比較多，所以對(duì)于該類型（如上圖中第一種和第四種植被類型）抽樣點(diǎn)，print(type_point);這句代碼會(huì)報(bào)錯(cuò)，是打印不出來結(jié)果的。
??其次，就是有一點(diǎn)，我這里這個(gè)代碼暫時(shí)不能直接限定每一種植被類型對(duì)應(yīng)的抽樣點(diǎn)個(gè)數(shù)——如上圖所示，有的植被類型抽樣點(diǎn)可以多到print(type_point);這句代碼直接報(bào)錯(cuò)（也就是超過了5000個(gè)），而有的植被類型抽樣點(diǎn)只有248個(gè)。這個(gè)問題，是var type_point = type_image.sample()這句代碼中，.sample()函數(shù)不能很好地限制抽樣點(diǎn)個(gè)數(shù)導(dǎo)致的。我的理解是，由于全球范圍內(nèi)，不同植被類型對(duì)應(yīng)的面積不一樣，也就是不同植被類型對(duì)應(yīng)的MODIS地表分類數(shù)據(jù)產(chǎn)品的像元總數(shù)也就不一樣；但是我們?cè)趫?zhí)行.sample()函數(shù)時(shí)，其抽樣比例是固定的（例如每3000個(gè)像素中，選出1個(gè)作為抽樣點(diǎn)）——那么對(duì)于在全球面積占比大，也就是像元個(gè)數(shù)多的植被類型，其得到的抽樣點(diǎn)個(gè)數(shù)自然就多；反之自然就少。針對(duì)這種情況，如果大家希望對(duì)于每一種植被類型，都獲取比較一致的抽樣點(diǎn)個(gè)數(shù)，那么就可以單獨(dú)對(duì)這種（或者這幾種）植被類型運(yùn)行上述代碼（即修改landcoverTypes中的類型數(shù)值），然后調(diào)整.sample()函數(shù)中的numPixels參數(shù)——這個(gè)參數(shù)越大，對(duì)應(yīng)的抽樣點(diǎn)個(gè)數(shù)就越多；反之，抽樣點(diǎn)個(gè)數(shù)則越少。大家可以結(jié)合實(shí)際情況，自行調(diào)整這個(gè)參數(shù)，直到能夠獲取到合適數(shù)量的抽樣點(diǎn)。
??接下來，我們看一下地圖。如下圖所示，可以看到8種植被類型對(duì)應(yīng)的抽樣點(diǎn)，以及MODIS地表分類數(shù)據(jù)產(chǎn)品，都已經(jīng)顯示出來了。由于抽樣點(diǎn)個(gè)數(shù)比較多，所以看起來密密麻麻的；甚至有的地方都重疊在了一起，看起來仿佛是純黑色的。
??最后，我們看一下Tasks一欄。如下圖所示，可以看到8種植被類型導(dǎo)出的任務(wù)都已經(jīng)生成了。
??我們逐一執(zhí)行這些任務(wù)，即可將不同種植被類型對(duì)應(yīng)的抽樣點(diǎn)合集導(dǎo)出到Assets；如下圖所示。
??其中，我們打開第一種植被類型對(duì)應(yīng)的抽樣點(diǎn)，也就是上圖中的point_1，如下圖所示。
??可以看到，這個(gè)植被類型對(duì)應(yīng)的抽樣點(diǎn)遍布全球，且總數(shù)有12292個(gè)，數(shù)量巨大——怪不得前面的print(type_point);這句代碼，對(duì)這個(gè)植被類型的抽樣點(diǎn)而言會(huì)報(bào)錯(cuò)。
??至此，大功告成。
