基于Python实现nc批量转tif格式

 更新时间:2022年08月12日 14:27:24   作者:NCR1008  
做项目有时会运用到netCDF格式的气象数据,而ArcGIS中需要用栅格影像进行处理,对于较多的文件,ArcGIS一个个手动转换过于繁琐,因此我们采用Python进行转换,下面就是Python实现nc批量转tif格式的示例代码,希望对你有所帮助

由于做项目需要运用到netCDF格式的气象数据,而ArcGIS中需要用栅格影像进行处理,对于较多的文件,ArcGIS一个个手动转换过于繁琐,因此我们采用Python进行转换,当然也可以采用matlab进行转换。

首先需要安装下面几个库:

import os
import netCDF4 as nc
import numpy as np
from osgeo import gdal, osr, ogr
import glob

我们可以在下面网址中寻找对应python安装版本的安装包,下载后,在pycharm控制台中直接安装即可。例如pip install netCDF4-1.5.8-cp39-cp39-

win_amd64.whl

https://www.lfd.uci.edu/~gohlke/pythonlibs/

安装之后即可进行转换:

def nc2tif(data, Output_folder):
    tmp_data = nc.Dataset(data)  # 利用.Dataset()方法读取nc数据
    print('tmp_data', tmp_data)
 
    Lat_data = tmp_data.variables['lat'][:]
    Lon_data = tmp_data.variables['lon'][:]
    # print(Lat_data)
    # print(Lon_data)
 
    tmp_arr = np.asarray(tmp_data.variables['temp'])
 
    # 影像的左上角&右下角坐标
    Lonmin, Latmax, Lonmax, Latmin = [Lon_data.min(), Lat_data.max(), Lon_data.max(), Lat_data.min()]
    # print(Lonmin, Latmax, Lonmax, Latmin)
 
    # 分辨率计算
    Num_lat = len(Lat_data)  # 5146
    Num_lon = len(Lon_data)  # 7849
    Lat_res = (Latmax - Latmin) / (float(Num_lat) - 1)
    Lon_res = (Lonmax - Lonmin) / (float(Num_lon) - 1)
    # print(Num_lat, Num_lon)
    # print(Lat_res, Lon_res)
 
    for i in range(len(tmp_arr[:])):
        # i=0,1,2,3,4,5,6,7,8,9,...
        # 创建tif文件
        driver = gdal.GetDriverByName('GTiff')
        out_tif_name = Output_folder + '\\' + data.split('\\')[-1].split('.')[0] + '_' + str(i + 1) + '.tif'
        out_tif = driver.Create(out_tif_name, Num_lon, Num_lat, 1, gdal.GDT_Int16)
 
        # 设置影像的显示范围
        # Lat_re前需要添加负号
        geotransform = (Lonmin, Lon_res, 0.0, Latmax, 0.0, -Lat_res)
        out_tif.SetGeoTransform(geotransform)
 
        # 定义投影
        prj = osr.SpatialReference()
        prj.ImportFromEPSG(4326)  # WGS84
        out_tif.SetProjection(prj.ExportToWkt())
 
        # 数据导出
        out_tif.GetRasterBand(1).WriteArray(tmp_arr[i])  # 将数据写入内存,此时没有写入到硬盘
        out_tif.FlushCache()  # 将数据写入到硬盘
        out_tif = None  # 关闭tif文件
 
def main():
    Input_folder = r"E:\competition\航天宏图\2-meter air temperature_CMFD\Data_forcing_01yr_010deg\\"
    # Input_folder = r"E:\competition\航天宏图\2-meter air temperature_CMFD\Data_forcing_01mo_010deg\\"
    Output_folder = r"E:\competition\航天宏图\2-meter air temperature_CMFD\tif\\"
 
    # 读取所有数据
    data_list = glob.glob(os.path.join(Input_folder + '*.nc'))
    print(data_list)
 
    for i in range(len(data_list)):
        data = data_list[i]
        nc2tif(data, Output_folder)
        print(data + '转tif成功')

值得注意的是,tmp_arr = np.asarray(tmp_data.variables['temp'])中的temp可以根据需要转换的波段进行选择,我们可以在读取数据之后print一下,找到对应的波段进行替换即可。如下图中我们应该选择的就是temp。

完成上述步骤即可得到所需的tif图像:

在上述代码中,经过处理的影像是倒置的,可能是处理过程中仿射矩阵读写错误导致的。因此我们可以在写入影像的时候,进行影像的垂直镜像操作即可:WriteArray(ndvi_arr_float[i][::-1]) 

def NC_to_tiffs(data, Output_folder):
    nc_data_obj = nc.Dataset(data)
    Lon = nc_data_obj.variables['lon'][:]
    Lat = nc_data_obj.variables['lat'][:]
    ndvi_arr = np.asarray(nc_data_obj.variables['temp'])  
    ndvi_arr_float = ndvi_arr.astype(float) / 10000  之间
    # 影像的左上角和右下角坐标
    LonMin, LatMax, LonMax, LatMin = [Lon.min(), Lat.max(), Lon.max(), Lat.min()]
    # 分辨率计算
    N_Lat = len(Lat)
    N_Lon = len(Lon)
    Lon_Res = (LonMax - LonMin) / (float(N_Lon) - 1)
    Lat_Res = (LatMax - LatMin) / (float(N_Lat) - 1)
    for i in range(len(ndvi_arr[:])):
        driver = gdal.GetDriverByName('GTiff')
        out_tif_name = Output_folder + '\\' + data.split('\\')[-1].split('.')[0] + '_' + str(i + 1) + '.tif'
        out_tif = driver.Create(out_tif_name, N_Lon, N_Lat, 1, gdal.GDT_Float32)
       
        geotransform = (LonMin, Lon_Res, 0, LatMax, 0, -Lat_Res)
        out_tif.SetGeoTransform(geotransform)
      
        srs = osr.SpatialReference()
        srs.ImportFromEPSG(4326) 
        out_tif.SetProjection(srs.ExportToWkt()) 
        # 数据写出
        out_tif.GetRasterBand(1).WriteArray(ndvi_arr_float[i][::-1])  # 将数据写入内存,此时没有写入硬盘 此处[::-1]用于图像的垂直镜像对称,避免图像颠倒
        out_tif.FlushCache()  # 将数据写入硬盘
        out_tif = None  # 注意必须关闭tif文件

这样便可以得到正确输出的影像:

当然,我们除了在写入时做垂直镜像操作之外,还可以利用python对图像进行几何变换来实现翻转。具体代码如下:

图像水平翻转:

#  图像水平翻转
    im_data_hor = np.flip(im_data, axis=2)
    hor_path = train_image_path + "\\" + str(tran_num) + imageList[i][-4:]
    writeTiff(im_data_hor, im_geotrans, im_proj, hor_path)

标签水平翻转: 

 #  标签水平翻转
    Hor = cv2.flip(label, 1)
    hor_path = train_label_path + "\\" + str(tran_num) + labelList[i][-4:]
    cv2.imwrite(hor_path, Hor)
    tran_num += 1

图像垂直翻转:

 #  图像垂直翻转
    im_data_vec = np.flip(im_data, axis=1)
    vec_path = train_image_path + "\\" + str(tran_num) + imageList[i][-4:]
    writeTiff(im_data_vec, im_geotrans, im_proj, vec_path)

标签垂直翻转:

 #  标签垂直翻转
    Vec = cv2.flip(label, 0)
    vec_path = train_label_path + "\\" + str(tran_num) + labelList[i][-4:]
    cv2.imwrite(vec_path, Vec)
    tran_num += 1

图像镜像翻转:

 #  图像对角镜像
    im_data_dia = np.flip(im_data_vec, axis=2)
    dia_path = train_image_path + "\\" + str(tran_num) + imageList[i][-4:]
    writeTiff(im_data_dia, im_geotrans, im_proj, dia_path)

标签镜像翻转:

 #  标签对角镜像
    Dia = cv2.flip(label, -1)
    dia_path = train_label_path + "\\" + str(tran_num) + labelList[i][-4:]
    cv2.imwrite(dia_path, Dia)
    tran_num += 1

若是输出路径的文件夹没有建立好,则会报如下错误。当然,为了减少工作量,也可以定义一个函数,如果路径不存在则自动创建,就可以解决这个问题。

到此这篇关于基于Python实现nc批量转tif格式的文章就介绍到这了,更多相关Python nc转tif内容请搜索脚本之家以前的文章或继续浏览下面的相关文章希望大家以后多多支持脚本之家!

相关文章

  • Django Auth应用实现用户身份认证

    Django Auth应用实现用户身份认证

    Django Auth 应用一般用在用户的登录注册上,用于判断当前的用户是否合法。本文将介绍Auth的另一个功能,即认证用户身份,感兴趣的同学可以关注一下
    2021-12-12
  • python 安装impala包步骤

    python 安装impala包步骤

    这篇文章主要介绍了python 安装impala包步骤,具有很好的参考价值,希望对大家有所帮助。一起跟随小编过来看看吧
    2020-03-03
  • python pyenv多版本管理工具的使用

    python pyenv多版本管理工具的使用

    这篇文章主要介绍了python pyenv多版本管理工具的使用,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的朋友们下面随着小编来一起学习学习吧
    2019-12-12
  • Python Django 数据库的相关操作详解

    Python Django 数据库的相关操作详解

    下面小编就为大家带来一篇django数据库的相关操作,小编觉得挺不错的,现在就分享给大家,也给大家做个参考。一起跟随小编过来看看吧
    2021-11-11
  • python开发中range()函数用法实例分析

    python开发中range()函数用法实例分析

    这篇文章主要介绍了python开发中range()函数用法,以实例形式较为详细的分析了Python中range()函数遍历列表的相关技巧,需要的朋友可以参考下
    2015-11-11
  • python 中defaultdict()对字典进行初始化的用法介绍

    python 中defaultdict()对字典进行初始化的用法介绍

    这篇文章主要介绍了python 中defaultdict()对字典进行初始化,一般情况下,在使用字典时,先定义一个空字典(如dict_a = {}),然后往字典中添加元素只需要 dict_a[key] = value即可,本文通过实例代码介绍具体用法,需要的朋友可以参考下
    2022-07-07
  • Python实践之使用Pandas进行数据分析

    Python实践之使用Pandas进行数据分析

    在数据分析领域,Python的Pandas库是一个非常强大的工具。这篇文章将为大家详细介绍如何使用Pandas进行数据分析,希望对大家有所帮助
    2023-04-04
  • Python实现删除列表中满足一定条件的元素示例

    Python实现删除列表中满足一定条件的元素示例

    这篇文章主要介绍了Python实现删除列表中满足一定条件的元素,结合具体实例形式对比分析了Python针对列表元素的遍历、复制、删除等相关操作技巧,需要的朋友可以参考下
    2017-06-06
  • JS设计模式之责任链模式实例详解

    JS设计模式之责任链模式实例详解

    这篇文章主要介绍了JS设计模式之责任链模式,结合实例形式详细分析了责任链模式的概念、原理、功能、使用场景及相关操作技巧,需要的朋友可以参考下
    2018-02-02
  • python调用opencv实现猫脸检测功能

    python调用opencv实现猫脸检测功能

    这篇文章主要介绍了python调用opencv实现猫脸检测功能,具有一定的参考价值,感兴趣的小伙伴们可以参考一下
    2019-01-01

最新评论