Как преобразовать файл .grib в GeoTIFF с правильной проекцией, используя Python, GDAL, ArcPy - PullRequest
Новая
       31

Как преобразовать файл .grib в GeoTIFF с правильной проекцией, используя Python, GDAL, ArcPy

0 голосов
/ 08 января 2019

Я пытаюсь преобразовать файл .grib в GeoTIFF для использования в ГИС (точнее, в ArcGIS), но у меня возникают проблемы с правильным отображением изображения. Я смог создать GeoTIFF, используя GDAL в Python, который показывает данные, но не отображается в правильном месте при переносе в ArcGIS. Полученное изображение ниже.

What the data looks like when I put it into ArcGIS

Данные, с которыми я работаю, можно загрузить с: https://gimms.gsfc.nasa.gov/SMOS/SMAP/L05/

Я пытаюсь проецировать данные в WGS84 Web Mercator (Вспомогательная сфера), EPSG: 3857

Примечание. Я попытался ввести данные через ArcMap, создав растровую мозаику, которая должна работать с данными .grib, но мне не повезло.

Обновление: я также пытался использовать инструмент Project Raster, но ArcGIS не нравится проекция по умолчанию, которая берется из файла .grib, и выдает ошибку.

Код, который я использую: 

import gdal

src_filename = r"C:\att\project\arcshare\public\disaster_response\nrt_products\smap\20150402_20150404_anom1.grib"
dst_filename = r"C:\att\project\arcshare\public\disaster_response\nrt_products\smap\smap_py_test1.tif"

#Open existing dataset
src_ds = gdal.Open(src_filename)

#Open output format driver, see gdal_translate --formats for list
format = "GTiff"
driver = gdal.GetDriverByName( format )

#Output to new format
dst_file = driver.CreateCopy( dst_filename, src_ds, 0 )

#Properly close the datasets to flush to disk
dst_ds = None
src_ds = None

Я не очень хорошо разбираюсь в использовании GDAL или GDAL в Python, поэтому любая помощь или советы будут с благодарностью.

1 Ответ

0 голосов
/ 04 августа 2019

Примерно так должно преобразовать ваши родные координаты в желаемую проекцию. Это еще не проверено. (Может по широте вместо широт). 

from cfgrib import xarray_store
from pyproj import Proj, transform

grib_data = xarray_store.open_dataset('your_grib_file.grib')
lat = grib_data.latitudes.value
lon = grib_data.longitudes.value
lon_transformed, lat_transformed = transform (Proj(init='init_projection'), 
 Proj(init='target_projection', lon, lat)
...