0

I downloaded the H05 – PR-OBS-5 Accumulated precipitation at ground data from H-SAF. The Product User Manual states the dataset is in rectangular stereographic projection centred on 42°N, 10°E. I need to project it to the Web Mercator projection (EPSG:3857) to use it with my other datasets, preferably with GDAL.

This is the original dataset as downloaded in grb format and displayer in qGis.

enter image description here

As I couldn't figure the EPSG code of the aforementioned coordinate system I used simply:

'gdalwarp -t_srs ' + outSrs + ' -dstnodata ' + outNoData + ' -r near -multi -ot ' + outType + ' -of ' + outFormat + ' ' + inFile + ' ' + outFile

I format the code as string here as I'm calling gdalwarp from python with subprocess.call.

This produces multiple "ERROR 1: tolerance condition error" and an "ERROR 1: Reprojection failed, err = -20, further errors will be suppressed on the transform object.", but returns the following:

enter image description here

If I combine this image with my other data it doesn't align properly as it should. I am not familiar with grb format but it seems that GDAL cannot properly interpret the projection information from the dataset.

Does anyone knows how to properly project this dataset to the Web Mercator?

EDIT:

gdalinfo output:

ERROR 1: tolerance condition error

ERROR 1: tolerance condition error

ERROR 1: tolerance condition error

ERROR 1: tolerance condition error

Driver: GRIB/GRIdded Binary (.grb)

Files: C:/xxxx/xxxx/xxxx/xxxx/h0524/5km/daily/2018/08/07/h05_20180807_0000_24_rom.grb

Size is 1900, 900

Coordinate System is:

PROJCS["unnamed",

GEOGCS["Coordinate System imported from GRIB file",

DATUM["unknown",

SPHEROID["Spheroid imported from GRIB file",6378140,298.2529810614917]],

PRIMEM["Greenwich",0],

UNIT["degree",0.0174532925199433]],

PROJECTION["Geostationary_Satellite"],

PARAMETER["central_meridian",0],

PARAMETER["satellite_height",35785831],

PARAMETER["false_easting",0],

PARAMETER["false_northing",0]]

Origin = (-5571679.196145262569189,5574935.774084444157779)

Pixel Size = (5861.840290526315584,-12374.996168888888860)

Corner Coordinates:

Upper Left (-5571679.196, 5574935.774)

Lower Left (-5571679.196,-5562560.778)

Upper Right ( 5565817.356, 5574935.774)

Lower Right ( 5565817.356,-5562560.778)

Center ( -2930.920, 6187.498) ( 0d 1'34.78"W, 0d 3'21.45"N)

Band 1 Block=1900x1 Type=Float64, ColorInterp=Undefined

Description = 0.000[-] undefined ()

NoData Value=9999

Metadata:

GRIB_COMMENT=Estimated precipitation [kg/(m^2)]

GRIB_ELEMENT=ESTP

GRIB_FORECAST_SECONDS=0 sec

GRIB_PDS_PDTN=30

GRIB_PDS_TEMPLATE_NUMBERS=1 0 8 54 1 0 0 0 0 0 255 255 255 255 255

GRIB_REF_TIME= 1533600000 sec UTC

GRIB_SHORT_NAME=0 undefined

GRIB_UNIT=[kg/(m^2)]

GRIB_VALID_TIME= 1533600000 sec UTC

Band 2 Block=1900x1 Type=Float64, ColorInterp=Undefined

Description = 0.000[-] undefined ()

NoData Value=9999

Metadata:

GRIB_COMMENT=(prodType 3, cat 0, subcat 192) [-] [-]

GRIB_ELEMENT=unknown

GRIB_FORECAST_SECONDS=0 sec

GRIB_PDS_PDTN=30

GRIB_PDS_TEMPLATE_NUMBERS=0 192 192 0 0

GRIB_REF_TIME= 1533600000 sec UTC

GRIB_SHORT_NAME=0 undefined

GRIB_UNIT=[-]

GRIB_VALID_TIME= 1533600000 sec UTC

  • Try forcing the input CRS: '+proj=stere +lat_0=42 +lon_0=10 +k_0=1.0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs' – mkennedy Aug 7 '18 at 18:56
  • This is slightly better but still not what I am looking for. I am editing the question by adding gdalinfo of the original grb file if it will be of any help. – Alešinar Aug 8 '18 at 6:46

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.