企业官网建设流程全解析

简介：这份资源是中国1km土地利用遥感监测数据集，覆盖1980至2015年间共45个年份，面向地理学、环境科学、城乡规划与经济学等领域的研究者及政策分析人员，用于分析近四十年土地利用变迁、城市扩张、生态保护与农业动态等议题。压缩包共126个文件，以adf、nit、dat等ArcInfo栅格与属性文件为主，辅以doc说明文档、log日志及xml元数据，整体约21.83MB，按年份分文件组织，便于逐年对比与时间序列分析。目前已有2896人学习下载，具备一定参考热度。数据以1km×1km格网记录耕地、林地、草地、建设用地、水域等类型，可用于识别耕地转建设用地、林地恢复等变化，评估城市化方向与生态政策成效，也可与社会经济数据融合构建土地利用变化模型，为气候变化适应与灾后重建研究提供基础支撑。

1. 拿到「中国1km土地利用数据1980-2015.zip」先别急着解压：这份数据到底能干什么

很多人第一次拿到这份数据，第一反应是解压、打开、看颜色，然后发现一堆 GeoTIFF 和几个说明文件，瞬间不知道从哪下手。我见过太多人卡在这一步：数据是真的，但不知道它能回答什么问题。这份数据本质上是 1980 到 2015 年、空间分辨率 1km 的全国土地利用栅格，每个像元用一个整数编码代表当时的土地类型——耕地、林地、草地、水域、建设用地、未利用地，通常还有二级分类。它的核心价值不是「看」，而是「算变化」：哪块地在哪一年变成了什么，变化速率是多少，空间上集中在哪。适合做国土空间演变分析、生态评估、城市扩张监测、碳储量估算这类工作的人。如果你只是想要一张好看的土地利用图，这份数据有点大材小用；但如果你想做时间序列上的空间分析，它是绕不开的基础底图。

2. 从压缩包到可分析数据：解压、编码对齐与投影检查

2.1 先看清压缩包里有什么，再决定怎么解压

拿到压缩包后不要直接双击解压到桌面。先在命令行里看一眼文件清单，确认年份命名规则和文件格式。常见做法是用unzip -l列出内容，或者用 Python 的zipfile模块读取目录。这一步能帮你判断数据是按年份分文件，还是按大类分文件，以及有没有附带分类编码表。

# 列出压缩包内容，不实际解压 unzip -l "中国1km土地利用数据1980-2015.zip" # 如果文件很多，只看前 30 行 unzip -l "中国1km土地利用数据1980-2015.zip" | head -30

逻辑说明：-l参数只列出文件清单，不写入磁盘，避免解压出一堆文件后才发现命名混乱。参数说明：如果压缩包有密码，这里会提示输入；如果没有，直接输出文件名、大小、日期。重点看文件名里有没有年份、有没有tif或img后缀，以及有没有class或code之类的说明文件。

2.2 解压后第一件事：确认分类编码和投影

解压到一个独立目录，比如landuse_1km/。然后打开说明文件，找到分类编码表。这份数据通常采用二级分类体系，一级类 6 个，二级类 20 多个。编码可能是 1、2、3……也可能是 11、12、21……不同版本编码规则不一样，必须对着说明文件确认。投影一般是 Albers 等积投影或经纬度投影，这直接影响你后面算面积时的单位。

import rasterio # 打开某一年的栅格，检查元数据 with rasterio.open("landuse_1km/1980.tif") as src: print("投影:", src.crs) print("尺寸:", src.width, "x", src.height) print("波段数:", src.count) print("像元大小:", src.res) print("无效值:", src.nodata) # 读取第一个波段的前 5 行 5 列，看看编码值 print(src.read(1)[:5, :5])

逻辑说明：rasterio是读栅格最稳的库，src.crs告诉你投影，src.res告诉你像元大小。如果投影是经纬度，像元大小大约是 0.0083 度左右；如果是 Albers，单位是米，像元大小接近 1000 米。参数说明：src.nodata很关键，很多分析翻车就是因为把无效值当成了某一类土地。读出来的前几行数值应该和编码表对得上，如果出现 0 或 255 这类值，先查说明文件确认是不是无效值。

2.3 用 Python 批量读取多年数据并统一无效值

单年读取只是热身，真正要做的是把 1980 到 2015 年所有年份读进来，统一无效值，再叠加分析。常见做法是写一个循环，把每年数据读成数组，存到一个字典或列表中。注意内存：全国 1km 栅格单年大约几千万像元，30 多年全读进内存可能吃紧，建议按需读取或分块处理。

import rasterio import numpy as np import os years = range(1980, 2016, 5) # 常见是 5 年一期，具体看文件 data = {} for year in years: path = f"landuse_1km/{year}.tif" if not os.path.exists(path): print(f"跳过缺失年份: {year}") continue with rasterio.open(path) as src: arr = src.read(1) nodata = src.nodata # 把无效值统一成 -1，方便后续掩膜 if nodata is not None: arr = np.where(arr == nodata, -1, arr) data[year] = arr print(f"{year} 读取完成，形状 {arr.shape}，唯一值前 10 个: {np.unique(arr)[:10]}")

逻辑说明：np.where把无效值替换成 -1，避免后续统计时把无效值算进某一类。参数说明：years的步长要根据实际文件调整，有的数据集是 5 年一期，有的是逐年。np.unique打印前 10 个唯一值，用来快速检查编码是否和说明文件一致。如果发现某个年份的唯一值数量异常少，可能是文件损坏或读取错误。

3. 土地利用变化分析：转移矩阵、面积统计与空间制图

3.1 用转移矩阵看两期之间到底变了什么

转移矩阵是土地利用变化分析里最核心的一张表。它告诉你从起始年到终止年，每一类土地有多少面积变成了其他类。计算方法是先做交叉制表，再乘以像元面积。注意：如果投影是经纬度，像元面积随纬度变化，不能直接用固定值乘，建议先投影到等积投影再算。

import numpy as np import pandas as pd def transition_matrix(arr_start, arr_end, pixel_area_km2=1.0): # 只保留有效像元 mask = (arr_start > 0) & (arr_end > 0) start = arr_start[mask] end = arr_end[mask] # 交叉制表 df = pd.crosstab(start, end) # 乘以像元面积得到面积矩阵 area_df = df * pixel_area_km2 return area_df # 假设 data 字典里已有 1980 和 2015 的数组 tm = transition_matrix(data[1980], data[2015], pixel_area_km2=1.0) print(tm.head())

逻辑说明：mask同时排除起始和终止的无效值，保证统计口径一致。pd.crosstab生成计数矩阵，乘以像元面积后得到面积矩阵。参数说明：pixel_area_km2在 Albers 投影下约等于 1，在经纬度投影下需要按纬度带计算，不能直接填 1。如果数据是 1km 分辨率，像元面积就是 1 平方公里，但前提是投影正确。

3.2 面积统计与变化速率：别把百分比算错

面积统计看起来简单，但坑很多。最常见的是把像元数直接当面积，忽略了投影单位。另一个坑是变化速率的分母选择：有人用起始年面积，有人用总面积，结果差很多。我一般用起始年该类面积作分母，算年均变化率，这样能反映该类土地自身的增减趋势。

def area_change_rate(tm, class_id): # tm 是转移矩阵，行是起始类，列是终止类 start_area = tm.loc[class_id].sum() # 起始年该类总面积 end_area = tm[class_id].sum() # 终止年该类总面积 if start_area == 0: return None # 假设时间跨度 35 年 years = 35 rate = ((end_area - start_area) / start_area) / years * 100 return { "起始面积": start_area, "终止面积": end_area, "净变化": end_area - start_area, "年均变化率(%)": rate } print(area_change_rate(tm, 1)) # 假设 1 是耕地

逻辑说明：tm.loc[class_id].sum()是起始年该类对所有终止类的面积之和，tm[class_id].sum()是终止年所有起始类变成该类的面积之和。参数说明：years要根据实际年份差调整，1980 到 2015 是 35 年。年均变化率用百分比表示，正数表示增加，负数表示减少。如果起始面积为 0，说明该类在起始年不存在，直接返回 None。

3.3 空间制图：把变化热点画出来

表格能告诉你变了多少，但变在哪要靠图。常见做法是生成变化栅格：把起始年和终止年不一样的像元标记出来，再按变化类型上色。用matplotlib加rasterio就能快速出图，不需要 GIS 桌面软件。

import matplotlib.pyplot as plt import rasterio with rasterio.open("landuse_1km/1980.tif") as src: arr1980 = src.read(1) extent = [src.bounds.left, src.bounds.right, src.bounds.bottom, src.bounds.top] with rasterio.open("landuse_1km/2015.tif") as src: arr2015 = src.read(1) # 变化像元：两期不一致且都有效 change = np.where((arr1980 > 0) & (arr2015 > 0) & (arr1980 != arr2015), 1, 0) plt.figure(figsize=(12, 8)) plt.imshow(change, cmap="Reds", extent=extent) plt.title("1980-2015 土地利用变化像元分布") plt.colorbar(label="是否变化") plt.savefig("change_map.png", dpi=300, bbox_inches="tight")

逻辑说明：change数组只有 0 和 1，1 表示该像元在 35 年间发生了类型变化。extent用地理坐标，保证图上有正确的空间范围。参数说明：cmap="Reds"让变化区域更醒目，dpi=300适合出版级出图。如果图太大，可以先裁剪到研究区再画，避免全国图看不清细节。

4. 避坑与排查：这份数据最容易翻车的 5 个地方

4.1 无效值没处理，统计结果全错

现象：面积统计出来某一类大得离谱，或者转移矩阵里出现 0 类到 0 类的大量转移。原因：栅格里的无效值（如 -9999、255、0）没被掩膜，被当成了某一类土地。解决：读数据后第一件事就是查nodata，用np.where统一替换成 -1，再做任何统计。如果说明文件没写无效值，用np.unique看有没有异常大的负数或异常小的正数。

4.2 投影没统一，面积算出来差几倍

现象：同样的像元数，算出来的面积和官方公布数据对不上。原因：经纬度投影下像元面积随纬度变化，直接用 1km×1km 乘会高估或低估。解决：先用rasterio或gdal把数据重投影到 Albers 等积投影，再算面积。如果不想重投影，至少按纬度带分段计算像元面积。

4.3 年份命名不连续，循环读取漏文件

现象：跑完循环发现少了几个年份，但代码没报错。原因：文件名不是严格的1980.tif、1985.tif，可能带前缀或后缀，或者某些年份缺失。解决：先用os.listdir打印所有文件名，确认命名规则，再用os.path.exists逐个检查。缺失年份要么跳过，要么用相邻年份插值，但要在报告里说明。

4.4 分类编码版本不一致，两期数据对不上

现象：转移矩阵里出现大量「未知类」或编码对不上的情况。原因：不同年份的数据可能用了不同的分类体系，比如早期是 6 类，后期是 20 多类。解决：先读说明文件，确认每一年的编码表。如果编码不一致，需要做重分类，把二级类合并到一级类再分析。重分类映射表要单独存成 CSV，方便复查。

4.5 内存不够，读全国数据直接崩

现象：循环读到一半程序被系统杀掉，或者报MemoryError。原因：全国 1km 栅格单年几千万像元，30 多年全读进内存超过普通机器上限。解决：不要一次性读所有年份。可以按年份逐对分析，或者用rasterio的窗口读取，分块处理。如果一定要全读，至少用float32而不是float64，能省一半内存。

5. 进阶技巧：用窗口读取和并行加速处理全国多年数据

全国 1km 数据做多年分析，最耗时的不是算法，是 I/O 和内存。我一般用两种策略：窗口读取和并行。窗口读取是把大栅格切成小块，逐块处理，内存占用可控；并行是按年份或按窗口分给多个进程，利用多核。下面是一个窗口读取的示例，配合concurrent.futures做并行。

import rasterio from rasterio.windows import Window import numpy as np from concurrent.futures import ProcessPoolExecutor def process_window(args): path, row_off, col_off, height, width = args with rasterio.open(path) as src: window = Window(col_off, row_off, width, height) arr = src.read(1, window=window) # 这里做你的分析，比如统计某一类面积 count = np.sum(arr == 1) # 假设 1 是耕地 return count def parallel_process(path, block_size=2048): with rasterio.open(path) as src: height, width = src.height, src.width tasks = [] for row_off in range(0, height, block_size): for col_off in range(0, width, block_size): h = min(block_size, height - row_off) w = min(block_size, width - col_off) tasks.append((path, row_off, col_off, h, w)) with ProcessPoolExecutor(max_workers=4) as executor: results = list(executor.map(process_window, tasks)) return sum(results) total = parallel_process("landuse_1km/2015.tif") print("耕地像元总数:", total)

逻辑说明：Window定义读取范围，src.read(1, window=window)只读窗口内的数据，内存占用从全国降到块级。ProcessPoolExecutor把每个窗口分给独立进程，绕开 Python 的 GIL。参数说明：block_size默认 2048，太小会增加任务调度开销，太大内存优势不明显，一般 1024 到 4096 之间。max_workers建议设为 CPU 核数的一半到全部，看 I/O 瓶颈。注意：并行写文件时要避免多个进程同时写同一个文件，这里只做统计所以安全。

另一个技巧是先把多年数据转成numpy的memmap，这样读的时候不占内存，写的时候直接落盘。但memmap对文件格式有要求，不是所有 GeoTIFF 都支持，需要先转成.npy或.dat。我一般只在需要反复随机访问时才这么做，否则窗口读取已经够用。

最后说一个我自己的习惯：每次跑完分析，先把中间结果存成parquet或csv，不要只存图。图是给人看的，表是给下一轮分析用的。土地利用变化这种工作，往往要反复调整分类和年份组合，中间结果存好了，后面改参数就是几分钟的事。希望帮到你。

本文还有配套的精品资源，点击获取