Python遥感开发之dat批量转tif掩膜裁剪

Python遥感开发之dat批量转tif掩膜裁剪

前言:使用numpy读取dat数据,再使用GDAL转成tif数据,再使用GDAL中的Warp方法进行掩膜裁剪。


了解dat数据

  • 遥感中的dat数据一般是由ENVI软件生成的遥感影像数据,每一条数据由xxx.dat和xxx.hdr两个文件组成,其中hdr里面存储了一些影像的遥感信息。

  • 一个遥感数据由dat和hdr两部分组成
    在这里插入图片描述

  • hdr存储的信息包括
    在这里插入图片描述

代码实现

  • 注意1:model中的tif数据是dat数据用argmap软件转的tif,目的是为了拿到遥感信息
  • 注意2:shp是掩膜裁剪的矢量图
  • 注意3:samples是列数,lines是行数,值可以在hdr文件里面拿到
import glob
import os

import numpy as np
from osgeo import gdal

def read_dat(file,samples,lines):
    data = np.fromfile(file, dtype=np.uint16)
    data.shape = [lines, samples]
    data = data.astype(np.float32)
    return data

def dat_to_tifdata(data,model,out_all_put):
    ds = gdal.Open(model)
    shape = data.shape
    driver = gdal.GetDriverByName("GTiff")
    dataset = driver.Create(out_all_put, shape[1], shape[0], 1, gdal.GDT_Float32)  # 以float类型进行存储
    dataset.SetGeoTransform(ds.GetGeoTransform())
    dataset.SetProjection(ds.GetProjection())
    dataset.GetRasterBand(1).WriteArray(data)
    return dataset

def get_filepath_dat(file_path):
    list1 = []  # 文件的完整路径
    list2 = [] #文件的名字
    if os.path.isdir(file_path):
        os.chdir(file_path)
        names = glob.glob("*.dat")  # 读取文件
        for f in names:
            list2.append(f)
            pre_data = file_path + '\\' + f  # 文件的完整路径
            list1.append(pre_data)
    return list1,list2

if __name__ == '__main__':
    model = r"C:\Users\Administrator\Desktop\model\dat1.tif"
    shp = r"D:\fly\研究方向\矢量图\江苏省\江苏省.shp"
    # dat = r"C:\Users\Administrator\Desktop\dat遥感\rad2000008.dat"
    dat = r"C:\Users\Administrator\Desktop\dat遥感"
    samples = 6164
    lines = 3540
    out = r"C:\Users\Administrator\Desktop\out"
    out_all = r"C:\Users\Administrator\Desktop\all"
    dat_path_list,dat_name_list = get_filepath_dat(dat)

    for dat_path,dat_name in zip(dat_path_list,dat_name_list):
        data= read_dat(dat_path,samples,lines)

        out_put = out + "\\" + dat_name + ".tif"
        out_all_put = out_all + "\\" + dat_name + ".tif"

        tif_data = dat_to_tifdata(data,model,out_all_put)
        gdal.Warp(out_put, tif_data, cutlineDSName=shp, cropToCutline=True)  # 按掩膜提取
        print(dat_name,"完成")

在这里插入图片描述


版权声明:本文为qq_32306361原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接和本声明。