Python nc文件批量转tif再转投影
import warnings import nc import netCDF4 warnings.filterwarnings('ignore') warnings.filterwarnings('ignore', category=DeprecationWarning) import netCDF4 import pandas as pd import numpy as np from osgeo import gdal import matplotlib.pyplot as plt import math import h5py from pyhdf.SD import * from osgeo import osr import numpy as np ### def array2raster(newRasterfn, rasterOrigin, xsize, ysize, array): """ newRasterfn: 输出tif路径 rasterOrigin: 原始栅格数据Extent xsize: x方向像元大小 ysize: y方向像元大小 array: 计算后的栅格数据 """ cols = array.shape[1] # 矩阵列数 rows = array.shape[0] # 矩阵行数 originX = rasterOrigin[0] # 起始像元经度 originY = rasterOrigin[1] # 起始像元纬度 driver = gdal.GetDriverByName('GTiff') outRaster = driver.Create(newRasterfn, cols, rows, 1, gdal.GDT_Float32) # 括号中两个0表示起始像元的行列号从(0,0)开始 outRaster.SetGeoTransform((originX, xsize, 0, originY, 0, ysize)) # 获取数据集第一个波段,是从1开始,不是从0开始 outband = outRaster.GetRasterBand(1) outband.WriteArray(array) outRasterSRS = osr.SpatialReference() # 代码4326表示WGS84坐标 outRasterSRS.ImportFromEPSG(4326) outRaster.SetProjection(outRasterSRS.ExportToWkt()) outband.FlushCache() ### from os import listdir import re inputf = 'E:/FSC_AVHRR/train_nc/' output_TIF = 'E:/FSC_AVHRR/train_tif_4326/' output_TIF_RE = 'E:/FSC_AVHRR/train_tif_6933/' files = listdir(inputf) pattern = '.nc$' timer = 0 for i in files: temp_input = inputf + i temp_out_name = re.split(pattern, i)[0] + '.tif' temp_out_T_name = re.split(pattern, i)[0] + '_T.tif' temp_output = output_TIF + temp_out_name temp_output_T = output_TIF_RE + temp_out_T_name nc_data = netCDF4.Dataset(temp_input) FSC = nc_data.variables['scfg'][:] # plt.imshow(C) FSC_data = np.flipud(FSC[0, :, :]) xsize = 0.05 ysize = -0.05 # 对于全球尺度rasterOrigin 为[-180(xmin),-90(ymin)] array2raster(temp_output, [-180, 90], xsize, ysize, FSC_data) option4 = gdal.WarpOptions(format='GTiff', srcSRS='EPSG:4326', dstSRS='EPSG:6933') FSC_RE = gdal.Warp(temp_output_T, temp_output, options=option4) timer = timer + 1 print(timer)