GenericMappingTools / gmtserver-admin

Cache data and script for managing the GMT data server
GNU Lesser General Public License v3.0
7 stars 3 forks source link

Error in mars_relief.recipe #189

Closed Esteban82 closed 1 year ago

Esteban82 commented 1 year ago

I run the code to process Mars MOLA grid (see below). The range of grid has values that do not match those of the recipe. As a result, the higher grid resolution weren't created (up to 30s included).

Info from mars_relief.recipe: # We use a precision of 0.5 m with a 6000 m offset to fit the range of -8528 to +21226 in 16-bit ints Info from gmt grdinfo: v_min: -8528 v_max: 31060 name: z

I think the Mount Olimpus is higher with a higher resolution.

Solution

So, I think that I have to change these values so the script would work.

# DST_SCALE=0.5
# DST_OFFSET=6000

I will try to figure out the new values.

Full message when I run the code:

bash scripts/srv_downsampler_grid.sh mars_relief
Convert Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_12s_p.grd=ns+s0.5+o6000
grdconvert [ERROR]: Cannot write format ns = GMT netCDF format (16-bit integer), CF-1.7.
grdconvert [ERROR]: The packed z-range, [-29056,50120], exceeds the maximum representable size. Adjust scale and offset parameters or remove out-of-range values.
grdconvert (gmt_api.c:5690(gmtapi_export_grid)): NetCDF: Numeric conversion not representable [mars/mars_relief/mars_relief_12s_p.grd=ns+s0.5+o6000]
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
grdedit [ERROR]: Cannot find file mars/mars_relief/mars_relief_12s_p.grd
Down-filter Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_15s_g.grd=ns+s0.5+o6000
grdfilter [ERROR]: Cannot write format ns = GMT netCDF format (16-bit integer), CF-1.7.
grdfilter [ERROR]: The packed z-range, [-29056,50120], exceeds the maximum representable size. Adjust scale and offset parameters or remove out-of-range values.
grdfilter (gmt_api.c:5690(gmtapi_export_grid)): NetCDF: Numeric conversion not representable [mars/mars_relief/mars_relief_15s_g.grd=ns+s0.5+o6000]
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
grdedit [ERROR]: Cannot find file mars/mars_relief/mars_relief_15s_g.grd
Down-filter Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_15s_p.grd=ns+s0.5+o6000
grdfilter [ERROR]: Cannot write format ns = GMT netCDF format (16-bit integer), CF-1.7.
grdfilter [ERROR]: The packed z-range, [-29052,50120], exceeds the maximum representable size. Adjust scale and offset parameters or remove out-of-range values.
grdfilter (gmt_api.c:5690(gmtapi_export_grid)): NetCDF: Numeric conversion not representable [mars/mars_relief/mars_relief_15s_p.grd=ns+s0.5+o6000]
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
grdedit [ERROR]: Cannot find file mars/mars_relief/mars_relief_15s_p.grd
Down-filter Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_30s_g.grd=ns+s0.5+o6000
grdfilter [ERROR]: Cannot write format ns = GMT netCDF format (16-bit integer), CF-1.7.
grdfilter [ERROR]: The packed z-range, [-29055,38147], exceeds the maximum representable size. Adjust scale and offset parameters or remove out-of-range values.
grdfilter (gmt_api.c:5690(gmtapi_export_grid)): NetCDF: Numeric conversion not representable [mars/mars_relief/mars_relief_30s_g.grd=ns+s0.5+o6000]
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
grdedit [ERROR]: Cannot find file mars/mars_relief/mars_relief_30s_g.grd
Down-filter Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_30s_p.grd=ns+s0.5+o6000
grdfilter [ERROR]: Cannot write format ns = GMT netCDF format (16-bit integer), CF-1.7.
grdfilter [ERROR]: The packed z-range, [-29044,42693], exceeds the maximum representable size. Adjust scale and offset parameters or remove out-of-range values.
grdfilter (gmt_api.c:5690(gmtapi_export_grid)): NetCDF: Numeric conversion not representable [mars/mars_relief/mars_relief_30s_p.grd=ns+s0.5+o6000]
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
[Session gmt.exe (0)]: Error returned from GMT API: GMT_GRID_WRITE_ERROR (19)
grdedit [ERROR]: Cannot find file mars/mars_relief/mars_relief_30s_p.grd
Down-filter Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc to mars/mars_relief/mars_relief_01m_g.grd=ns+s0.5+o6000
Esteban82 commented 1 year ago

If I undertand the logic behing packing closer via scale/offsets, I think that in this case there is no case to do it.

So, I should use these values.

# DST_SCALE=1
# DST_OFFSET=0

Am I right?

PaulWessel commented 1 year ago

But 31060 cannot be right since just googling what the peak is gives 21,229. It may have changed a tiny bit though. Should still fit with a scale of 0.5 but may need to adjust the offset so that value = (max - offset ) /scale and (min - offset) / scale both fit inside ±32767 m

PaulWessel commented 1 year ago

Something is odd with my gmt convert to nc from TIFF anyway:

Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: v_min: -6575 v_max: 0 name: z

Esteban82 commented 1 year ago

Should still fit with a scale of 0.5 but may need to adjust the offset so that value = (max - offset ) /scale and (min - offset) / scale both fit inside ±32767 m

Are you sure? (31060 + 8528 ) * 0.5 = 79176.

I have to change the scale to a value greater than 0.6 to have less than 65534 values.

PaulWessel commented 1 year ago

But where is 31060 coming from if Wikipedia says Olympus Mons is just over 20 km?

Esteban82 commented 1 year ago

From the grid itself. I use gmt grdinfo and I got: v_min: -8528 v_max: 31060 name: z. I got the same result when I did it on the tiff file.

PaulWessel commented 1 year ago

Well, here is a quote:

"The full range of topography on Mars is about 19 miles (30 kilometers), one and a half times the range of elevations found on Earth," noted Dr. David Smith of NASA's Goddard Space Flight Center, Greenbelt, MD, the principal investigator for MOLA

. S as long as the range is < 32767 m then we can use a scale of 0.5 and a suitable offset and store in units of 50 cm.

Esteban82 commented 1 year ago

So, I guess that there is problem with the grid.

Or maybe I used the wrong one.

Esteban82 commented 1 year ago

This is from the 11 GB tif.

$ gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Title: Grid imported via GDAL
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Command:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Remark:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Pixel node registration used [Geographic grid]
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Grid file format: gd = Import/export through GDAL
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: x_min: -180 x_max: 179.998447904 x_inc: 0.00337412083064 name: x n_columns: 106694
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: y_min: -89.9992239522 y_max: 90 y_inc: 0.00337412083064 name: y n_rows: 53347
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: v_min: -8528 v_max: 31060 name: z
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: scale_factor: 1 add_offset: 0
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Default CPT:
Esteban82 commented 1 year ago

And this from the 3.4 GB NC grid:

$ gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Title: Grid imported via GDAL
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Command:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Remark:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Pixel node registration used [Geographic grid]
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Grid file format: nf = GMT netCDF format (32-bit float), CF-1.7
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: x_min: -180 x_max: 180 x_inc: 0.00337413537781 name: longitude n_columns: 106694
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: y_min: -90 y_max: 90 y_inc: 0.00337413537781 name: latitude n_rows: 53347
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: v_min: -8528 v_max: 31060 name: z
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: scale_factor: 1 add_offset: 0
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: format: netCDF-4 chunk_size: 129,129 shuffle: on deflation_level: 3
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Default CPT:
GEOGCS["unknown",
    DATUM["unknown",
        SPHEROID["unknown",3396190,0]],
    PRIMEM["Reference meridian",0],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AXIS["Longitude",EAST],
    AXIS["Latitude",NORTH]]
joa-quim commented 1 year ago

Does gdalinfo say the same thing?

Esteban82 commented 1 year ago

Well, I made a map of Olympus and its top is barely over 21 km.

gmt grdimage Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc -R222/232/14/21 -png Olympus -JW15 -Baf

Olympus

Esteban82 commented 1 year ago

Does gdalinfo say the same thing?

Yes, the same.

actual_range={-8528,31060}

Full Info:

$ gdalinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc 
Driver: netCDF/Network Common Data Format
Files: Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc
Size is 106694, 53347
Coordinate System is:
GEOGCRS["unknown",
    DATUM["unknown",
        ELLIPSOID["unknown",3396190,0,
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]]],
    PRIMEM["Reference meridian",0,
        ANGLEUNIT["degree",0.0174532925199433,
            ID["EPSG",9122]]],
    CS[ellipsoidal,2],
        AXIS["longitude",east,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]],
        AXIS["latitude",north,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]]
Data axis to CRS axis mapping: 1,2
Origin = (-180.000000000000000,90.000000000000000)
Pixel Size = (0.003374135377809,-0.003374135377809)
Metadata:
  grid_mapping#spatial_ref=GEOGCS["unknown",
    DATUM["unknown",
        SPHEROID["unknown",3396190,0]],
    PRIMEM["Reference meridian",0],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AXIS["Longitude",EAST],
    AXIS["Latitude",NORTH]]
  lat#actual_range={-90,90}
  lat#axis=Y
  lat#long_name=latitude
  lat#standard_name=latitude
  lat#units=degrees_north
  lon#actual_range={-180,180}
  lon#axis=X
  lon#long_name=longitude
  lon#standard_name=longitude
  lon#units=degrees_east
  NC_GLOBAL#Conventions=CF-1.7
  NC_GLOBAL#description=
  NC_GLOBAL#GMT_version=6.4.0 [64-bit] [MP]
  NC_GLOBAL#history=
  NC_GLOBAL#node_offset=1
  NC_GLOBAL#title=Grid imported via GDAL
  z#actual_range={-8528,31060}
  z#grid_mapping=grid_mapping
  z#long_name=z
  z#_FillValue=-nan
Corner Coordinates:
Upper Left  (-180.0000000,  90.0000000) (180d 0' 0.00"W, 90d 0' 0.00"N)
Lower Left  (-180.0000000, -90.0000000) (180d 0' 0.00"W, 90d 0' 0.00"S)
Upper Right ( 180.0000000,  90.0000000) (180d 0' 0.00"E, 90d 0' 0.00"N)
Lower Right ( 180.0000000, -90.0000000) (180d 0' 0.00"E, 90d 0' 0.00"S)
Center      (   0.0000000,   0.0000000) (  0d 0' 0.01"E,  0d 0' 0.01"N)
Band 1 Block=129x129 Type=Float32, ColorInterp=Undefined
  NoData Value=nan
  Metadata:
    actual_range={-8528,31060}
    grid_mapping=grid_mapping
    long_name=z
    NETCDF_VARNAME=z
    _FillValue=-nan
joa-quim commented 1 year ago

And what about?

gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc -L0

Esteban82 commented 1 year ago

I got the same.

Later I will make a map to see where the high values are.

PaulWessel commented 1 year ago

I would have swore I did run through Mars many months ago and got sensible results. I have not refreshed the original TIFF file nor nc file I created. So I was going to blame GDAL but not sure. How come I get max = 0 while @Esteban82 gets 31k?

joa-quim commented 1 year ago

Later I will make a map to see where the high values are.

grdinfo -M tells you that.

PaulWessel commented 1 year ago

gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc -M Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Title: Grid imported via GDAL Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Command: Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Remark: Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Pixel node registration used [Geographic grid] Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Grid file format: nf = GMT netCDF format (32-bit float), CF-1.7 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: x_min: -180 x_max: 180 x_inc: 0.00337413537781 name: longitude n_columns: 106694 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: y_min: -90 y_max: 90 y_inc: 0.00337413537781 name: latitude n_rows: 53347 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: v_min: -6575 at x = -7.92078279941 y = 69.9154591636 v_max: 0 at x = -179.998312932 y = 68.1541604964 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: scale_factor: 1 add_offset: 0 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: 6474 nodes (0.0%) set to NaN

Esteban82 commented 1 year ago

I found some strange values near the north pole!! They are also in the tiff file.

Mars_NorthPole_tif2

$ gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc -M
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Title: Grid imported via GDAL
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Command:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Remark:
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Pixel node registration used [Geographic grid]
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Grid file format: nf = GMT netCDF format (32-bit float), CF-1.7
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: x_min: -180 x_max: 180 x_inc: 0.00337413537781 name: longitude n_columns: 106694
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: y_min: -90 y_max: 90 y_inc: 0.00337413537781 name: latitude n_rows: 53347
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: v_min: -8528 at x = 62.1397641854 y = -32.7999700077 v_max: 31060 at x = -133.779406527 y = 88.4360882524
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: scale_factor: 1 add_offset: 0
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: 39956 nodes (0.0%) set to NaN
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: format: netCDF-4 chunk_size: 129,129 shuffle: on deflation_level: 3
Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc: Default CPT:
PaulWessel commented 1 year ago

Maybe we both had transmission issues or theirs has errors. I find it hard to believe theirs has errors given widespread usage?

PaulWessel commented 1 year ago

@joa-quim any ideas?

joa-quim commented 1 year ago

Nope, but I guess I'll have to untie it. Where is that Mars tiff file?

Does it have a checksum?

PaulWessel commented 1 year ago

Recipe says # SRC_FILE=https://planetarymaps.usgs.gov/mosaic/Mars/HRSC_MOLA_Blend/Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif

Esteban82 commented 1 year ago

Sorry. I checked the wrong file.

e98092a06bea74d5ec9a704dcd1fe318 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif

https://planetarymaps.usgs.gov/mosaic/Mars/HRSC_MOLA_Blend/Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif.md5

joa-quim commented 1 year ago

So you can check if the checksum of the download file is the same.

Esteban82 commented 1 year ago

Is it ok with this command md5sum Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif ?

If so, I got a different hash. b18c673dddd8c22a1ae1995f20c1540a *Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif

Esteban82 commented 1 year ago

When I download the file, I had some issues. Do you have any type on how to download it? Any software to use?

joa-quim commented 1 year ago

I'm having troubles downloading that file with Firefox. Connection keeps breaking. Maybe try with wget

PaulWessel commented 1 year ago

I am trying curl which is what our src_downloader_grid.sh script does. Good you are trying wget

joa-quim commented 1 year ago

I can't download it. Permanently breaking connection.

PaulWessel commented 1 year ago

Yep:

curl: (18) transfer closed with 11023654838 bytes remaining to read

PaulWessel commented 1 year ago

I am trying another source. It has an original GMT grid that is 43 Gb and the 11 Gb tiff.

PaulWessel commented 1 year ago

This went very fast and gave this output which seems sane:

gmt grdinfo Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Title: Grid imported via GDAL Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Command: Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Remark: Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Pixel node registration used [Geographic grid] Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: Grid file format: gd = Import/export through GDAL Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: x_min: -180 x_max: 179.998447904 x_inc: 0.00337412083064 name: x n_columns: 106694 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: y_min: -89.9992239522 y_max: 90 y_inc: 0.00337412083064 name: y n_rows: 53347 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif: v_min: -8528 v_max: 21226 name: z

PaulWessel commented 1 year ago

Any luck for you to download this file? If should I place it in my SOEST ftp just in case?

Esteban82 commented 1 year ago

Sorry, I forgot to ask you if was going to do it. I can do it. I have no problem.

I can't donwload it. Yes, put it your SOEST, please.

Esteban82 commented 1 year ago

Now I could enter. I will let you know if I have any problem.

PaulWessel commented 1 year ago

I was unable to copy over the juge TIF file but running the grdedit command (see _marsrecipe) I was able to copy over the NC file. So, you can get that from my public ftp. File is Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.nc. If you place that in the staging area then the recipe should detect that and skip on to the conversions.

Esteban82 commented 1 year ago

Thanks Paul. I manage to download it. I will try to process it on Monday.

Esteban82 commented 1 year ago

I download the file. It is ok.

The cheksum is the same that USGS's.

md5sum staging/Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif
e98092a06bea74d5ec9a704dcd1fe318 *staging/Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif

USGS:
e98092a06bea74d5ec9a704dcd1fe318 Mars_HRSC_MOLA_BlendDEM_Global_200mp_v2.tif
Esteban82 commented 1 year ago

I run this code and I got this image. So this issus can be closed. gmt grdimage @mars_relief -Rg -Bf -Cgeo -png Mars

Mars