网上基本没有这些代码，需要的可以看看，仅供参考
HDF直接转成Tif的也有代码，可以看看我的其他博文
HDF转CSV
from pyhdf.SD import SDimport numpy as npimport osfilepath="E:CO2" #文件位置，存放hdf文件的文件夹filenames=os.listdir(filepath)print(filenames)for i,filename in enumerate(filenames): pathhdf=filepath+'\'+filename print(pathhdf) hdfarr = pathhdf.split("\")[-1].split("_")[2].split(".") csvname = pathhdf.split("\")[-1].split(".")[0] + "." + pathhdf.split("\")[-1].split(".")[1] print(hdfarr) y = hdfarr[0] m = hdfarr[1] ylist = [] mlist = [] for i in range(41 * 71): ylist.append(y) mlist.append(m) hdf = SD(pathhdf) print(hdf.info()) # 信息类别数 data = hdf.datasets() a = [] b = [] CO2Amount_grid = [] for i in data: print(i) # 具体类别 Lat = hdf.select(i)[:] lat = np.array(Lat) print(lat.shape) if i == 'lat_grid': a.append(lat) if i == 'lon_grid': b.append(lat) if i == 'CO2Amount_grid': CO2Amount_grid.append(lat) print(type(CO2Amount_grid)) CO2Amount_grid = np.array(list(CO2Amount_grid)).reshape(-1) # print(a) grid = [] grid_1 = [] for i in range(41): for k in range(71): grid_1.append(a[0][i]) grid_1 = np.array(grid_1).squeeze() print(grid_1.shape) grid_1 = grid_1.reshape(-1) print(grid_1.shape) grid_2 = [] for k in range(41): grid_2.append(b) grid_2 = np.array(grid_2).squeeze() print(grid_2.shape) grid_2 = grid_2.reshape(-1) print(grid_2.shape) grid.append(ylist) grid.append(mlist) grid.append(grid_1) grid.append(grid_2) grid.append(CO2Amount_grid) head = ["year", "month", "lat", "lon", "CO2Amount"] grid = np.array(grid, dtype=np.str) print(grid.shape) grid = grid.T print(grid.shape) grid = np.insert(grid, 0, np.array(head), axis=0) print(grid.shape) # print(CO2Amount_grid) np.savetxt("E:\temp\"+csvname + ".csv", grid, delimiter=",", fmt='%s')
效果如下
CSV转tif
import numpy as npfrom osgeo import gdal, osrfrom glob import globimport osimport datetimeimport pandas as pdos.environ['PROJ_LIB'] = r'D:pythonpythonDatavenvLibsite-packagesosgeodataproj'def csv2array(filen): df = pd.read_csv(filen) lat = (np.array((df.loc[:,'lat']).tolist())).reshape(-1,1) lon = (np.array((df.loc[:,'lon']).tolist())).reshape(-1,1) co2 = (np.array((df.loc[:,'CO2Amount']).tolist())).reshape(-1,1) ns = [] nl = [] for line in lon: ns.append(int(line)) ns = sorted(ns) unique_data1 = np.unique(ns) for line in lat: nl.append(int(line)) nl = sorted(nl) unique_data2 = np.unique(nl) xx,yy = [len(unique_data2),len(unique_data1)] outdata = co2.reshape(xx,yy) res_lon = abs(unique_data1[1] - unique_data1[0]) res_lat = -abs(unique_data2[1] - unique_data2[0]) latmax = np.max(unique_data2) lonmin = np.min(unique_data1) return outdata,res_lon,res_lat,latmax,lonmindef writetif(data,outname,geotransform): nl,ns = [data.shape[0],data.shape[1]] bands = 1 driver = gdal.GetDriverByName("GTiff") out_tif = driver.Create(outname, ns, nl, bands, gdal.GDT_Float32) out_tif.SetGeoTransform(geotransform) srs = osr.SpatialReference() proj_type = 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]' out_tif.SetProjection(proj_type) # 给新建图层赋予投影信息 out_tif.GetRasterBand(1).WriteArray(data) del out_tifinpath = 'E:/temp'outpath = 'E:/temp'dirs = os.listdir(inpath) # 获取指定路径下的文件file_name = []for x in dirs: # 循环读取路径下的文件并筛选输出 if os.path.splitext(x)[1] == '.csv' or os.path.splitext(x)[1] == '.CSV': # 筛选.HDF文件 file_name.append(x)for k in range(len(file_name)): filename = inpath+'/'+file_name[k] tt = file_name[k] outdata,res_lon,res_lat,latmax,lonmin = csv2array(filename) geotransform = (lonmin, res_lon, 0, latmax, 0, res_lat) outname = outpath + '/' + tt[0:len(tt)-3] + 'tif' writetif(outdata,outname,geotransform)
注意，有些库 python3版本不能直接下载，例如proj这个包，pychram可能识别不到，需要下载下来进行本地安装，我的另一篇博客解决了这个问题