2020年全国各省、地级市、县的10米分辨率的土地利用数据的制作方法与分享

#引入模块包
import os
from osgeo import gdal, gdalconst
import rasterio as rio
import rasterio.mask
import rasterio
from tqdm import tqdm

def Land_Cover(Mask_PATH,Data_path,OutPut_path):
    tifPaths_folder_SHENG = os.listdir(Mask_PATH)
    for mask_path in tqdm(tifPaths_folder_SHENG):
        try:
            #省目录
            Land_Cover_SHENG_PATH = os.path.join(mask_path, Mask_PATH)
            #获得文件中的名字
            SHENG_PATH=Land_Cover_SHENG_PATH+"\\\\"+mask_path+"\\\\"+mask_path+".shp"
            #裁剪函数
            City = gpd.read_file(SHENG_PATH)
            Land_Cover = gpd.read_file(Data_path)
            City=City.to_crs(Land_Cover .crs)
            Sheng_Land_Cover= gpd.clip(Parking_shp,City)
            mkdir(OutPut_path+mask_path+"_2017年土地利用数据_10米分辨率")
            #保存文件SHP
            save_Land_Cover_path = OutPut_path + mask_path + "\\\\" + mask_path +"土地利用数据"+ ".tif"
            Sheng_Land_Cover.to_file(save_SHP_path)
        except:
            pass

#引入模块包
import os
from osgeo import gdal, gdalconst
import rasterio as rio
import rasterio.mask
import rasterio
from tqdm import tqdm

def mosaic_tif(tifPaths_folder)#tifPaths_folder是输入文件的文件夹路径
    # 循环目录
    for path in tqdm(tifPaths_folder):
        Land_Cover_SMALL_PATH = os.path.join(tifPath, path)
        for (pathname, dirs, files) in os.walk(Land_Cover__SMALL_PATH):
            son_Paths_file=files
        # 如果影像数量大于一
        if len(son_Paths_file) >= 2:
            Land_Cover_SMALL_PATH2=Land_Cover_SMALL_PATH+"\\\\"
            # 循环子目录，进行镶嵌
            #循环同一个文件下的tif文件
            inputFiles = []
            for path_small in son_Paths_file:
                # 每一个栅格的路径
                son_Paths_PATH = os.path.join(Land_Cover_SMALL_PATH2, path_small)
                #读取影像
                inputrasfile = gdal.Open(son_Paths_PATH, gdal.GA_ReadOnly)  # 读取影像
                inputProj = inputrasfile.GetProjection()  # 获取坐标系
                inputFiles.append(inputrasfile)  # 推入列表
                options = gdal.WarpOptions(srcSRS=inputProj,  # 输入坐标系
                                           dstSRS=inputProj,  # 输出坐标系
                                           format='GTiff',  # 图像格式
                                           resampleAlg=gdalconst.GRIORA_NearestNeighbour,  # 重采样算法，这里是双线性内插
                                           dstNodata=Nodata,  # 缺省值
                                           cutlineLayer=outline,  # 输出范围，这里可以是一个外轮廓shp数据
                                           outputType=gdalconst.GDT_Int16)  # 数据类型，这里是有符号32位整型
            #输出文件名
            outputfilePath ="输出文件夹"+"\\\\"+path+"_Land_Cover_"+"_10m分辨率_Esri数据集"+".tif"
            #写栅格
            gdal.Warp(outputfilePath, inputFiles, options=options)  # 图像镶嵌