python gdal.Warp tif裁剪、拼接、重采样

参考https://www.jianshu.com/p/821c741ff169

拼接

def rasterwarp(base_file,inpath,outpath):
    if not os.path.exists(outpath):
        os.mkdir(outpath)
    out_file = os.path.join(outpath,'mcd12q1.tif')
    inputrasfile1 = gdal.Open(base_file, gdal.GA_ReadOnly)  # 第一幅影像
    inputProj1 = inputrasfile1.GetProjection()
    options = gdal.WarpOptions(srcSRS=inputProj1, dstSRS=inputProj1, format='GTiff', resampleAlg=gdalconst.GRA_Max)
    rasterlist = []
    filelist = pathlib.Path(inpath).glob('*.tif')
    for rasterfile in filelist:
        rasterfile = str(rasterfile)
        inputrasfile2 = gdal.Open(rasterfile, gdal.GA_ReadOnly)  # 第二幅影像
        rasterlist.append(inputrasfile2)

    gdal.Warp(out_file, rasterlist, options=options)

裁剪

gdal.Warp('D:\MCD12Q1\globe\\arcpy\con_result\\final2_mcd12q1.tif', 'D:\MCD12Q1\globe\\arcpy\con_result\\re_mcd12q1.tif',outputBounds=[-180,-60,-160,60],dstSRS=Projection,resampleAlg=gdalconst.GRA_NearestNeighbour)

重采样

def warp_image(img_file: str, base_file: str, out_file: str) -> int:
    src = gdal.Open(img_file, gdalconst.GA_ReadOnly)
    src_proj = src.GetProjection()
    # src_geotrans = src.GetGeoTransform()
    nband = src.RasterCount
    dataType = src.GetRasterBand(1).DataType

    match_ds = gdal.Open(base_file, gdalconst.GA_ReadOnly)
    match_proj = match_ds.GetProjection()
    match_geotrans = match_ds.GetGeoTransform()
    wide = match_ds.RasterXSize
    high = match_ds.RasterYSize

    dst = gdal.GetDriverByName('GTiff').Create(out_file, wide, high, nband, dataType)
    dst.SetGeoTransform(match_geotrans)
    dst.SetProjection(match_proj)

    
    gdal.ReprojectImage(src, dst, src_proj, match_proj, gdalconst.GRA_NearestNeighbour)

    del dst  
    return 0