python 土地分类动态变化

数据来源：立方数据学社地表覆盖数据

from osgeo import gdal
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib as mlp

LC=gdal.Open('./数据/上海/2000上海市.tif')
data_2000=LC.ReadAsArray(0,0)

LC=gdal.Open('./数据/上海/2010上海市.tif')
data_2010=LC.ReadAsArray(0,0)

LC=gdal.Open('./数据/上海/2020上海市.tif')
data_2020=LC.ReadAsArray(0,0)

data_pd = pd.DataFrame(data_2000)
data1 = data_pd[data_pd != 65535]

data_pd1 = pd.DataFrame(data_2010)
data2 = data_pd1[data_pd1 != 65535]

data_pd2 = pd.DataFrame(data_2020)
data3 = data_pd2[data_pd2 != 0]

#创建画图，以及camera
plt.rcParams['font.sans-serif'] = ['KaiTi']  #全局字体替换_中文显示
plt.rcParams['axes.unicode_minus'] = False  #坐标轴负数的负号显示问题
from celluloid import Camera
fig= plt.figure()
camera = Camera(fig)
labels = ['2000年','2010年','2020年']

n=0
#定义映射颜色
cmap = mpl.colors.ListedColormap(['red', 'green', 'c','lightcoral','w','b','k','cyan','y','m'])
tick = ['耕地','森林','草地','灌木地','湿地','水体','苔原','人造地表','裸地','冰川']
bounds = [10,20,30,40,50,60,70,80,90,100]
norm = mpl.colors.BoundaryNorm(bounds, cmap.N)
for i in Data:  
    plt.imshow(i,cmap,vmin=10,vmax=100)
    plt.text(200,500,labels[n],fontsize=20)
    n=n+1
    camera.snap()
    
bar1 = plt.colorbar(ticks=bounds)
bar1.ax.set_yticklabels(tick)
#plt.colorbar()
animation = camera.animate(interval=1000,blit=True)
animation.save('上海市土地变化.gif')