簡(jiǎn)介全國(guó)生態(tài)功能區(qū)劃修編版矢量數(shù)據(jù)資源包專為地理信息、生態(tài)評(píng)價(jià)、國(guó)土空間規(guī)劃、環(huán)境管理領(lǐng)域的科研人員和工程師準(zhǔn)備幫助快速獲得標(biāo)準(zhǔn)化的全國(guó)生態(tài)功能分區(qū)邊界及屬性信息。壓縮包共含七個(gè)文件核心為shp矢量圖層并配套dbf屬性表、prj投影坐標(biāo)、shx、sbn、sbx索引文件和xml元數(shù)據(jù)合計(jì)七點(diǎn)七二兆字節(jié)結(jié)構(gòu)緊湊可被ArcGIS、QGIS等主流軟件直接調(diào)用。該資源經(jīng)由作者整理分享當(dāng)前已有三百一十三人學(xué)習(xí)或下載。在實(shí)際使用中用戶能夠直接開展疊加分析、屬性查詢和專題圖繪制用于生態(tài)紅線評(píng)估、環(huán)境承載力測(cè)算、區(qū)域開發(fā)適宜性評(píng)價(jià)等任務(wù)數(shù)據(jù)邊界與分類代碼源自修編版具備較好權(quán)威性免去自行找圖和矢量化流程大幅節(jié)省前期處理時(shí)間。對(duì)于需要全國(guó)尺度生態(tài)區(qū)劃數(shù)據(jù)支撐的研究、規(guī)劃與教學(xué)場(chǎng)景這是一份可直接入庫(kù)使用的基礎(chǔ)數(shù)據(jù)有助于快速搭建分析底圖、開展空間統(tǒng)計(jì)與決策輔助。1. 全國(guó)生態(tài)功能區(qū)劃 2015修編版shp tif 一份底圖生態(tài)評(píng)估與選線避讓的起點(diǎn)在哪里做環(huán)境影響識(shí)別、區(qū)域規(guī)劃里的生態(tài)紅線對(duì)接或者把一條擬建公路沿線的生態(tài)敏感類型拉出來分析手里沒有一份全國(guó)生態(tài)功能區(qū)劃 2015修編版底圖報(bào)告里很難給出讓人信服的空間依據(jù)。這份數(shù)據(jù)和常見路網(wǎng)、行政區(qū)劃不一樣發(fā)布形態(tài)一般是 shp 矢量邊界配合 tif 柵格shp 用來查功能區(qū)名稱、等級(jí)和代碼tif 像元值則可以直接參與柵格計(jì)算、面積統(tǒng)計(jì)以及和 DEM、土地利用這類數(shù)據(jù)做疊置分析。兩套格式看著互補(bǔ)實(shí)際用起來卻往往因?yàn)樽鴺?biāo)系不統(tǒng)一、NoData 沒處理、柵格和面邊界對(duì)不齊在一開始就把人卡住。下面按我實(shí)際操作的順序把讀入、預(yù)處理、裁剪、轉(zhuǎn)換和避坑的完整流程梳理一遍照著走能少花很多前期摸索的時(shí)間。2. 數(shù)據(jù)讀入與預(yù)處理先搞清 shp 字段和 tif 值域再做坐標(biāo)系對(duì)齊2.1 屬性表里查生態(tài)功能區(qū)常見字段、區(qū)劃代碼與名稱的對(duì)應(yīng)不同渠道拿到的數(shù)據(jù)屬性字段命名經(jīng)常不統(tǒng)一我用過的版本里既有叫 code、name、type 的也有直接用拼音縮寫如 xzq行政區(qū)、stgnq生態(tài)功能區(qū)的。不管字段長(zhǎng)什么樣建議進(jìn) ArcMap 后用 Python 窗口先把字段名全部打出來再確認(rèn)區(qū)劃代碼到底存在哪個(gè)字段里。# ArcMap Python 窗口 import arcpy fc rD:\eco\eco_function_area.shp fields [f.name for f in arcpy.ListFields(fc)] print(fields) with arcpy.da.SearchCursor(fc, fields[:6]) as cur: for i, row in enumerate(cur): if i 6: break print([str(v)[:30] for v in row])這段代碼用 arcpy.ListFields 把屬性字段全部列出來再用 SearchCursor 讀前 6 條記錄的字段值。這樣能快速判斷 shp 里到底有沒有“區(qū)劃名稱”“生態(tài)功能類型”“所屬省份”這類信息避免做圖例時(shí)到處找字段。打印時(shí)把字符串截到 30 位只是為了讓輸出整齊不影響源數(shù)據(jù)本身?？炊侄魏笾攸c(diǎn)看兩個(gè)關(guān)鍵項(xiàng)區(qū)劃代碼和功能類型。2015 修編版的區(qū)劃體系是“生態(tài)功能區(qū)—生態(tài)功能亞區(qū)—生態(tài)功能小區(qū)”三個(gè)層級(jí)但屬性表不一定每層都列全。大多數(shù)發(fā)布版本會(huì)有一個(gè)漢字字段寫一級(jí)功能類型比如水源涵養(yǎng)、水土保持、防風(fēng)固沙、生物多樣性保護(hù)等而 tif 柵格里的像元值往往是針對(duì)這套類型做了數(shù)值映射。所以先把代碼和名稱的對(duì)照關(guān)系記下來后面做柵格轉(zhuǎn)面、按類型統(tǒng)計(jì)面積全靠這層映射不出錯(cuò)。空間參考也建議在圖層屬性里看一眼。有些渠道下載的華東、華北分幅數(shù)據(jù)自帶投影坐標(biāo)有些直接是 WGS84 地理坐標(biāo)。拿北京地區(qū)舉例如果 shp 是 CGCS2000 下的高斯投影帶tif 卻是 WGS84兩者疊加會(huì)整體錯(cuò)位幾米到幾十米必須先統(tǒng)一。常見做法是全部轉(zhuǎn)到 CGCS2000或按你最終制圖要求轉(zhuǎn)到所在區(qū)域的高斯分帶然后矢量做投影、柵格做重采樣一步到位。2.2 柵格 tif 的像元值、NoData 和空值處理工具怎么選打開 tif 后第一件事不是看顏色而是看屬性里的像元大小、波段數(shù)和 NoData 值。生態(tài)功能區(qū)劃 tif 如果是單波段整數(shù)型像元值通常直接對(duì)應(yīng)區(qū)劃代碼如果包含多波段通常第一波段是區(qū)劃類型后續(xù)波段可能存放生態(tài)敏感性等級(jí)等輔助信息。不少發(fā)布版本為了壓縮體積把 NoData 設(shè)成 -9999 或 0這兩個(gè)值一旦參與面積統(tǒng)計(jì)結(jié)果會(huì)多出一大片“未知類型”區(qū)域。提示tif 加載后整幅圖是同一個(gè)色多半不是數(shù)據(jù)壞了而是符號(hào)化用了連續(xù)拉伸不適合整數(shù)類型的分類數(shù)據(jù)。改成“唯一值”分類符號(hào)化區(qū)劃界線立刻就能顯示出來。處理 NoData我常用的工具是柵格計(jì)算器里的 SetNull 和 IsNull也可以在 ArcToolbox 里找“條件函數(shù)”完成。下面表達(dá)式等價(jià)于一個(gè)條件賦值# ArcGIS 柵格計(jì)算器表達(dá)式 out SetNull(IsNull(eco_tif) | (eco_tif -9999), eco_tif)這句的意思是當(dāng)像元為 NoData或像元值等于 -9999 時(shí)把該位置設(shè)為空其余像元保留原值。IsNull 先生成一個(gè) 0/1 掩膜SetNull 再根據(jù)條件把對(duì)應(yīng)像元置空。之所以必須處理 -9999是因?yàn)楹芏嘬浖憱鸥駮r(shí)把“數(shù)據(jù)缺失”寫成了 -9999 而不是標(biāo)準(zhǔn) NoData不轉(zhuǎn)換的話后面轉(zhuǎn)面、統(tǒng)計(jì)會(huì)把 -9999 當(dāng)成一個(gè)真實(shí)的區(qū)劃類型生成整片無效圖斑。用 QGIS 或 GDAL 命令行等價(jià)操作是gdal_calc.py -A eco_2015.tif --outfileeco_clean.tif \ --calcA*(A0) --NoDataValue0這條命令把小于等于 0 的像元全部寫成 0 并標(biāo)記為 NoData。執(zhí)行前先看一眼直方圖確認(rèn) tif 的值是連續(xù)的短代碼如 101、102、201還是大數(shù)字編號(hào)如 1000001避免條件判斷寫錯(cuò)把有效類型一起抹掉。2.3 坐標(biāo)系與像元對(duì)齊投影轉(zhuǎn)換、重采樣方法和 Snap Rastershp 和 tif 都加載后如果兩者“看起來”重疊但相交面積統(tǒng)計(jì)出來卻少了一截多半是投影基準(zhǔn)和像元對(duì)齊的問題。常見做法是先把 tif 投影到和 shp 同一個(gè)坐標(biāo)系隨后在環(huán)境設(shè)置里把 Snap Raster 設(shè)為生態(tài)區(qū)劃 tifCell Size 也固定成它的像元大小。投影柵格時(shí)重采樣方法必須選“最近鄰”NEAREST。生態(tài)功能區(qū)劃是離散類型數(shù)據(jù)雙線性或三次卷積會(huì)對(duì)邊界像元做插值比如說 101 和 102 之間插出 101.5這個(gè)值轉(zhuǎn)面后就成了“未知類型”。而最近鄰法只取原始值不會(huì)生成新類別。矢量 shp 的投影用 Project 工具目標(biāo)坐標(biāo)系選 CGCS2000 的地區(qū)分帶或 Albers 等面積投影都可以關(guān)鍵是 shp 和 tif 最終必須在同一套坐標(biāo)系里。環(huán)境設(shè)置往往是新手最容易忽略的一步。在 ArcMap 的“環(huán)境設(shè)置”中把處理范圍設(shè)為生態(tài)區(qū)劃 shp 的范圍捕捉柵格設(shè)為生態(tài)區(qū)劃 tif像元大小填 tif 的像元尺寸。這樣后續(xù)任何柵格運(yùn)算的輸出都和源 tif 網(wǎng)格嚴(yán)格對(duì)齊不會(huì)出現(xiàn)錯(cuò)半格的情況。如果輸出偏差半個(gè)像元邊界會(huì)呈現(xiàn)鋸齒狀轉(zhuǎn)面后也很容易產(chǎn)生細(xì)碎窄條多邊形。對(duì)齊完成后可以用下面幾句檢查兩個(gè)數(shù)據(jù)的范圍是否一致import arcpy eco_tif rD:\eco_2015.tif eco_shp rD:\eco_function_area.shp r1 arcpy.Describe(eco_tif).extent r2 arcpy.Describe(eco_shp).extent print(r1.XMin, r1.YMin, r1.XMax, r1.YMax) print(r2.XMin, r2.YMin, r2.XMax, r2.YMax)前后兩行范圍值基本重合說明后續(xù)裁剪、提取都能正常工作差別明顯時(shí)以 shp 面范圍為準(zhǔn)做裁剪不要反過來用 tif 范圍去約束 shp否則容易把邊界切掉一圈。3. 用生態(tài)區(qū)劃 shp 裁剪 DEM 柵格Clip 與掩膜提取的差異和參數(shù)選擇3.1 Clip 和 Extract by Mask 在輸出邊界上的真實(shí)差異很多人搜“arcmap 中依靠面圖層裁剪 dem 柵格 tif 文件”和“依靠面圖層掩碼提取”分不清兩個(gè)工具到底有何區(qū)別。ArcGIS 的柵格 Clip數(shù)據(jù)管理可以按矩形范圍裁剪也可以勾選“使用輸入要素裁剪幾何”選項(xiàng)使輸出范圍貼合面要素邊界。Extract by Mask Spatial Analyst則是把掩膜面柵格化后保留掩膜內(nèi)部的像元掩膜外一律寫成 NoData。兩者從結(jié)果看經(jīng)常很接近但有三個(gè)實(shí)際差異NoData 策略不同。Clip 在裁剪范圍內(nèi)的 NoData 原樣保留Extract by Mask 會(huì)把掩膜邊界外全部設(shè)為 NoData輸出的有效像元范圍看起來更“干凈”。表達(dá)式能力不同。Extract by Mask 配合柵格計(jì)算器可以做帶條件的提取比如只保留生態(tài)類型等于“水源涵養(yǎng)”的像元Clip 沒有表達(dá)式選項(xiàng)只能按幾何切。性能不同。數(shù)據(jù)量大的時(shí)候Clip 更快因?yàn)樗蛔鲅谀ぶ胤诸愔苯影捶秶袎KExtract by Mask 多了掩膜柵格化這一步速度略慢。所以我的選擇策略是只做幾何范圍切割用 Clip簡(jiǎn)單直接后續(xù)還要按類型疊加統(tǒng)計(jì)用 Extract by Mask 更順手。但這里有個(gè)高頻坑之前的工程在環(huán)境設(shè)置里掛了“分析掩膜”或“捕捉柵格”ArcMap 會(huì)默默繼承這些設(shè)置導(dǎo)致明明選的全國(guó)范圍輸出還是別處的矩形塊。跑之前先打開環(huán)境設(shè)置把不相關(guān)的掩膜清掉。3.2 跟著操作ArcMap 中用生態(tài)區(qū)劃面裁剪 DEM 的掩膜提取步驟下面這套是我在 ArcMap 里最常用的流程前提是數(shù)據(jù)已經(jīng)按第 2 章做完坐標(biāo)系和 NoData 處理。我之前拿這套方案把全國(guó)生態(tài)功能區(qū)劃里的“水源涵養(yǎng)”區(qū)單獨(dú)提取出來再接上 30 米 DEM 算每個(gè)子區(qū)域的平均坡度。import arcpy arcpy.env.workspace rD:\eco_work arcpy.env.extent rD:\eco_work\eco_function_area.shp arcpy.env.snapRaster rD:\eco_work\eco_tif.tif arcpy.env.cellSize rD:\eco_work\eco_tif.tif arcpy.env.mask rD:\eco_work\eco_function_area.shp out arcpy.sa.ExtractByMask(rD:\dem_30m.tif, rD:\eco_work\eco_function_area.shp) out.save(rD:\eco_work\dem_eco_clip.tif)這段代碼把處理范圍、捕捉柵格、像元大小和分析掩膜統(tǒng)一設(shè)到生態(tài)功能區(qū) shp 上。ExtractByMask 的第二個(gè)參數(shù)可以是矢量面工具內(nèi)部會(huì)先把它柵格化。注意 Spatial Analyst 需要啟用擴(kuò)展模塊否則會(huì)直接報(bào)“工具不可用”這也是個(gè)高頻入門問題。環(huán)境參數(shù)里最有玄機(jī)的是 snapRaster。它保證輸出像元和生態(tài)區(qū) tif 的網(wǎng)格完全重合否則會(huì)把柵格重新對(duì)齊到當(dāng)前坐標(biāo)系默認(rèn)網(wǎng)格邊界像元錯(cuò)開半個(gè)像元。對(duì)于生態(tài)功能區(qū)的面積統(tǒng)計(jì)這半個(gè)像元在邊界上累計(jì)起來能差出幾百上千平方米。在環(huán)評(píng)報(bào)告里“面積對(duì)不上”往往會(huì)被評(píng)審盯住所以這個(gè)參數(shù)我每次都會(huì)顯式設(shè)置。如果手頭沒有 DEM這一步也可以直接裁生態(tài)區(qū) tif 本身邏輯完全一樣。很多人下載了全國(guó) tif 后先按省份邊界裁成省圖再按市、縣 shp 繼續(xù)切分這樣后續(xù)每次加載不用扛著全國(guó)范圍的大文件效率明顯更好。3.3 驗(yàn)證裁剪結(jié)果像元數(shù)、唯一值和邊界貼合度裁剪完成后不要直接拿去制圖先做三個(gè)快速檢查確認(rèn)輸出像元大小沒被重采樣改掉用唯一值符號(hào)化確認(rèn)值集合是原 tif 的子集疊加 shp 看邊界是否貼合。對(duì)應(yīng)命令行檢查import arcpy from arcpy.sa import * r Raster(rD:\eco_work\dem_eco_clip.tif) print(r.width, r.height, r.cellSize) uc arcpy.UniqueValues(r) print(uc[:20])這段代碼打印輸出柵格的寬度、高度、像元大小以及前 20 個(gè)唯一值?？吹?-9999 或 0 混在列表里就回頭看掩膜和 NoData 是否處理干凈。更直觀的辦法是把生態(tài)區(qū)劃 shp 和裁出來的 tif 疊加在 ArcMap 里沿著邊界放大到 1:5000 左右肉眼就能看出有沒有錯(cuò)開。不要嫌這一步啰嗦數(shù)據(jù)問題大多出在源數(shù)據(jù)本身把驗(yàn)證步驟固定下來每次跑都少踩很多坑。4. 生態(tài)區(qū)劃數(shù)據(jù)處理排查與避坑五個(gè)高頻問題的現(xiàn)象、原因、解決4.1 現(xiàn)象tif 加載后整幅圖灰黑一片看不到任何區(qū)劃圖斑原因出在符號(hào)化。生態(tài)功能區(qū)劃是離散分類柵格很多軟件默認(rèn)按連續(xù)拉伸顯示把所有像元值壓到一個(gè)灰度區(qū)間看起來就是一塊灰。有些版本還把 NoData 區(qū)域和有效區(qū)域混在一起有效圖斑占比又低符號(hào)化后幾乎看不見。解決方法是先把符號(hào)化方式從“拉伸”改成“唯一值”并給不同功能類型分配不同色帶。如果改成唯一值后仍是一片黑用第 2 章的 SetNull 表達(dá)式把 NoData 和 -9999 轉(zhuǎn)成空再重新符號(hào)化。絕大多數(shù)“數(shù)據(jù)壞了”的錯(cuò)覺到這里就好了。4.2 現(xiàn)象shp 與 tif 看似重疊邊界卻總有一圈錯(cuò)位現(xiàn)象是疊加后大輪廓基本對(duì)得上但放大看 shp 邊界和 tif 的像元邊界差出幾個(gè)像元像一圈貼邊。原因通常有兩個(gè)一個(gè)坐標(biāo)系不是同一套基準(zhǔn)比如 WGS84 和 CGCS2000 之間的差異在高精度要求下會(huì)表現(xiàn)出來另一個(gè)是 tif 經(jīng)過某次重采樣像元網(wǎng)格和 shp 的投影網(wǎng)格不再對(duì)齊。解決方法是把兩者統(tǒng)一投影到 CGCS2000 的同一分帶并在環(huán)境設(shè)置里把 Snap Raster 指向 tif。檢查辦法是在 ArcMap 里把 shp 和 tif 都打開用“按屬性選擇”選中任意一個(gè)圖斑放大后對(duì)比邊界。如果錯(cuò)位均勻且都在 10 米以內(nèi)通常就是投影差異如果錯(cuò)位方向不一致還要懷疑 shp 本身是否做過局部編輯這種只能對(duì)照原始發(fā)布數(shù)據(jù)重新下載。4.3 現(xiàn)象裁剪輸出邊緣鋸齒明顯柵格轉(zhuǎn)面后出現(xiàn)大量細(xì)碎多邊形裁剪之后邊緣呈明顯階梯狀轉(zhuǎn)成 shp 后沿邊界出現(xiàn)一串長(zhǎng)條狀多邊形這種問題多半是環(huán)境設(shè)置里沒有指定 Snap Raster 和 Cell Size系統(tǒng)按當(dāng)前數(shù)據(jù)框的默認(rèn)網(wǎng)格重新對(duì)齊了一次。生態(tài)區(qū)劃 tif 的像元原本是正南正北方向排列的重對(duì)齊后像元邊界就歪了。解決方法是重新執(zhí)行裁剪在環(huán)境設(shè)置里把捕捉柵格設(shè)為生態(tài)區(qū)劃 tifCell Size 填 tif 的原始像元尺寸。轉(zhuǎn)面時(shí)如果碎多邊形已經(jīng)產(chǎn)生會(huì)用“消除”工具合并但最穩(wěn)妥的辦法還是在裁剪階段就把網(wǎng)格對(duì)齊磨刀不誤砍柴工。4.4 現(xiàn)象柵格轉(zhuǎn)面后屬性表里的代碼全是 -9999類型字段為空柵格轉(zhuǎn)面工具Raster to Polygon在沒有先處理 NoData 的情況下會(huì)把 -9999 或 0 也當(dāng)成一類像元轉(zhuǎn)出來的面屬性里出現(xiàn)大量“-9999”記錄真正的區(qū)劃代碼和名稱反而查不到。原因是 NoData 在轉(zhuǎn)面時(shí)被當(dāng)作普通類值參與建面。解決方法是先按第 2 章做法把 -9999 設(shè)為 NoData或者用柵格計(jì)算器把非目標(biāo)范圍的像元賦空值再執(zhí)行轉(zhuǎn)面。另外一個(gè)常見習(xí)慣是直接把 tif 用“按掩膜提取”先裁到研究區(qū)范圍內(nèi)邊角就不會(huì)有 -9999 混進(jìn)來。轉(zhuǎn)面前順手跑一次唯一值統(tǒng)計(jì)看到列表里只有目標(biāo)代碼再轉(zhuǎn)。4.5 現(xiàn)象用漁網(wǎng)工具分割 shp 后各部分面積總和比原始面小用創(chuàng)建漁網(wǎng)工具把生態(tài)功能區(qū)分成網(wǎng)格再逐格裁剪最后匯總面積卻發(fā)現(xiàn)少了。原因是漁網(wǎng)是獨(dú)立生成的規(guī)則格網(wǎng)它的邊界不可能和原始面邊界完全重合。裁剪時(shí)漁網(wǎng)格與面邊界重疊的窄條被丟掉共同產(chǎn)生的邊緣損失累計(jì)后相當(dāng)可觀。解決方法是不要用漁網(wǎng)直接裁剪原始面而應(yīng)該用“相交”工具讓原始面與漁網(wǎng)線做空間相交相交后每個(gè)格網(wǎng)部分保留原面屬性面積由相交結(jié)果重新計(jì)算。這樣所有面積都遵循原始邊界不會(huì)出現(xiàn)丟邊。硬要用裁剪也是一條路但必須接受邊緣損失并在匯總時(shí)單獨(dú)說明。5. 讓生態(tài)區(qū)劃柵格快速?gòu)?fù)用柵格轉(zhuǎn)面、歸一化處理和導(dǎo)出 WKT、3D Tiles生態(tài)區(qū)劃 shp 和 tif 的價(jià)值往往不在全國(guó)圖本身而在把它轉(zhuǎn)成可復(fù)用的工程數(shù)據(jù)。我常用的三條路徑是柵格轉(zhuǎn)面、歸一化處理和轉(zhuǎn)成通用文本格式。柵格轉(zhuǎn)面用來把 tif 的離散類型變成矢量圖斑方便和其他面圖層做空間連接歸一化則解決多源柵格之間的數(shù)值尺度問題。歸一化處理我通常用柵格計(jì)算器表達(dá)式是(A - min) / (max - min)把生態(tài)區(qū)劃代碼壓到 0 到 1 之間。這對(duì)類型型數(shù)據(jù)而言只適合做顯示真正做疊加分析時(shí)還是保留原始代碼更可靠別為了統(tǒng)一尺度把分類語(yǔ)義丟掉。導(dǎo)出 WKT 可以用 GDAL 的 ogr2ogrogr2ogr -f CSV eco_wkt.csv eco_function_area.shp -dialect sqlite \ -sql SELECT code, name, ST_AsText(geometry) AS geom FROM eco_function_area這條命令把 shp 的每個(gè)面要素轉(zhuǎn)成一行 WKT 文本適合做數(shù)據(jù)庫(kù)入庫(kù)或程序?qū)?。如果目?biāo)是 3D Tiles常見做法是先用轉(zhuǎn)面工具把 tif 轉(zhuǎn)成 shp再通過其他工具鏈做瓦片化這一步并非 ArcMap 自帶需要單獨(dú)搭建。最后說一句我的個(gè)人習(xí)慣每次拿到生態(tài)區(qū)劃數(shù)據(jù)第一件事永遠(yuǎn)是復(fù)制一份原始文件然后在副本上進(jìn)坐標(biāo)系、NoData 和投影轉(zhuǎn)換絕不直接動(dòng)下載源文件。柵格數(shù)據(jù)破壞性操作沒有后悔藥寧可多占一點(diǎn)硬盤也別在原始底圖上反復(fù)折騰。這套流程幫我避開了很多沒法回退的局面希望也能幫到你。本文還有配套的精品資源點(diǎn)擊獲取