将多个 GeoTIFF 图像的栅格时间序列转换为 NetCDF
Convert raster time series of multiple GeoTIFF images to NetCDF
我有一个栅格时间序列存储在多个 GeoTIFF
文件 (*.tif
) 中,我想将其转换为单个 NetCDF
文件。数据为uint16
。
我可能会使用 gdal_translate
将每个图像转换为 netcdf,使用:
gdal_translate -of netcdf -co FORMAT=NC4 20150520_0164.tif foo.nc
然后使用 NCO
编写一些脚本以从文件名中提取日期然后连接,但我想知道我是否可以在 Python 中使用 xarray
更有效地执行此操作并且它是新的rasterio
后端。
我可以轻松阅读文件:
import glob
import xarray as xr
f = glob.glob('*.tif')
da = xr.open_rasterio(f[0])
da
哪个returns
<xarray.DataArray (band: 1, y: 5490, x: 5490)>
[30140100 values with dtype=uint16]
Coordinates:
* band (band) int64 1
* y (y) float64 5e+05 5e+05 5e+05 5e+05 5e+05 4.999e+05 4.999e+05 ...
* x (x) float64 8e+05 8e+05 8e+05 8e+05 8.001e+05 8.001e+05 ...
Attributes:
crs: +init=epsg:32620
我可以将其中之一写入 NetCDF:
ds.to_netcdf('foo.nc')
但理想情况下,我可以使用 xr.open_mfdataset
之类的东西,写入时间值(从文件名中提取),然后将整个聚合写入 netCDF
。并让 dask
处理核心内存问题。 :-)
可以用 xarray
和 dask
来完成这样的事情吗?
Xarray 应该能够为您完成连接步骤。我在下面稍微调整了你的例子。将文件名解析为有用的内容将取决于您。
import glob
import pandas as pd
import xarray as xr
def time_index_from_filenames(filenames):
'''helper function to create a pandas DatetimeIndex
Filename example: 20150520_0164.tif'''
return pd.DatetimeIndex([pd.Timestamp(f[:8]) for f in filenames])
filenames = glob.glob('*.tif')
time = xr.Variable('time', time_index_from_filenames(filenames))
chunks = {'x': 5490, 'y': 5490, 'band': 1}
da = xr.concat([xr.open_rasterio(f, chunks=chunks) for f in filenames], dim=time)
我有一个栅格时间序列存储在多个 GeoTIFF
文件 (*.tif
) 中,我想将其转换为单个 NetCDF
文件。数据为uint16
。
我可能会使用 gdal_translate
将每个图像转换为 netcdf,使用:
gdal_translate -of netcdf -co FORMAT=NC4 20150520_0164.tif foo.nc
然后使用 NCO
编写一些脚本以从文件名中提取日期然后连接,但我想知道我是否可以在 Python 中使用 xarray
更有效地执行此操作并且它是新的rasterio
后端。
我可以轻松阅读文件:
import glob
import xarray as xr
f = glob.glob('*.tif')
da = xr.open_rasterio(f[0])
da
哪个returns
<xarray.DataArray (band: 1, y: 5490, x: 5490)>
[30140100 values with dtype=uint16]
Coordinates:
* band (band) int64 1
* y (y) float64 5e+05 5e+05 5e+05 5e+05 5e+05 4.999e+05 4.999e+05 ...
* x (x) float64 8e+05 8e+05 8e+05 8e+05 8.001e+05 8.001e+05 ...
Attributes:
crs: +init=epsg:32620
我可以将其中之一写入 NetCDF:
ds.to_netcdf('foo.nc')
但理想情况下,我可以使用 xr.open_mfdataset
之类的东西,写入时间值(从文件名中提取),然后将整个聚合写入 netCDF
。并让 dask
处理核心内存问题。 :-)
可以用 xarray
和 dask
来完成这样的事情吗?
Xarray 应该能够为您完成连接步骤。我在下面稍微调整了你的例子。将文件名解析为有用的内容将取决于您。
import glob
import pandas as pd
import xarray as xr
def time_index_from_filenames(filenames):
'''helper function to create a pandas DatetimeIndex
Filename example: 20150520_0164.tif'''
return pd.DatetimeIndex([pd.Timestamp(f[:8]) for f in filenames])
filenames = glob.glob('*.tif')
time = xr.Variable('time', time_index_from_filenames(filenames))
chunks = {'x': 5490, 'y': 5490, 'band': 1}
da = xr.concat([xr.open_rasterio(f, chunks=chunks) for f in filenames], dim=time)