mercantile.InvalidLatitudeError: Y can not be computed: lat=90

[raj-techprescient]

I am trying to serve GeoTIFFs using terracotta and using MySQL as the datastore. However, I am encountering an error during the process.

The GeoTIFF I am using has the following bounds

Upper Left (-180.0000051, 89.9999975) (180° 0' 0.02"W, 89° 59' 59.99"N)
Lower Left (-180.0000051, -90.0000051) (180° 0' 0.02"W, 90° 0' 0.02"S)
Upper Right ( 180.0000051, 89.9999975) (180° 0' 0.02"E, 89° 59' 59.99"N)
Lower Right ( 180.0000051, -90.0000051) (180° 0' 0.02"E, 90° 0' 0.02"S)
Center ( 0.0000000, -0.0000038) ( 0° 0' 0.01"E, 0° 0' 0.01"S)

Traceback

Traceback (most recent call last):
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\mercantile\__init__.py", line 391, in _xy
    y = 0.5 - 0.25 * math.log((1.0 + sinlat) / (1.0 - sinlat)) / math.pi
                              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
ZeroDivisionError: float division by zero

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 1536, in __call__
    return self.wsgi_app(environ, start_response)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 1514, in wsgi_app
    response = self.handle_exception(e)
               ^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask_cors\extension.py", line 194, in wrapped_function
    return cors_after_request(app.make_response(f(*args, **kwargs)))
                                                ^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 1511, in wsgi_app
    response = self.full_dispatch_request()
               ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 919, in full_dispatch_request
    rv = self.handle_user_exception(e)
         ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask_cors\extension.py", line 194, in wrapped_function
    return cors_after_request(app.make_response(f(*args, **kwargs)))
                                                ^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 917, in full_dispatch_request
    rv = self.dispatch_request()
         ^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\flask\app.py", line 902, in dispatch_request
    return self.ensure_sync(self.view_functions[rule.endpoint])(**view_args)  # type: ignore[no-any-return]
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\terracotta\server\singleband.py", line 149, in get_singleband
    return _get_singleband_image(keys, tile_xyz)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\terracotta\server\singleband.py", line 197, in _get_singleband_image
    image = singleband(parsed_keys, tile_xyz=tile_xyz, **options)
            ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Python312\Lib\contextlib.py", line 81, in inner
    return func(*args, **kwds)
           ^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\terracotta\handlers\singleband.py", line 50, in singleband
    tile_data = xyz.get_tile_data(
                ^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\terracotta\xyz.py", line 41, in get_tile_data
    if not tile_exists(wgs_bounds, tile_x, tile_y, tile_z):
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\terracotta\xyz.py", line 60, in tile_exists
    mintile = mercantile.tile(bounds[0], bounds[3], tile_z)
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\mercantile\__init__.py", line 415, in tile
    x, y = _xy(lng, lat, truncate=truncate)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\rajsh\projects\POC\poc-terra\venv\Lib\site-packages\mercantile\__init__.py", line 393, in _xy
    raise InvalidLatitudeError("Y can not be computed: lat={!r}".format(lat))

Upon inspecting the metadata table in MySQL, I noticed that the bounds are being rounded off and stored as

bound_noth -180,
bound_east -90,
bound_south 180,
bound_west 90

Note:- This works fine with SQLite DB.

[dionhaefner]

That's funny. Can you dump a schema of that table? What's the type of those bounds?

[raj-techprescient]

That's funny. Can you dump a schema of that table? What's the type of those bounds?
Here is the table schema.

Code reference.
(https://github.com/DHI/terracotta/blob/main/terracotta/drivers/relational_meta_store.py#L79)

[dionhaefner]

That's interesting. The data type is float even though I think it should be double. We specify the type like this:

SQLA_METADATA_TYPE_LOOKUP: Dict[str, Any] = {
        "real": functools.partial(sqla.types.Float, precision=8),
        "text": sqla.types.Text,
        "blob": sqla.types.LargeBinary,
    }

MySQL manual says this:

For FLOAT, the SQL standard permits an optional specification of the precision (but not the range of the exponent) in bits following the keyword FLOAT in parentheses, that is, FLOAT(p). MySQL also supports this optional precision specification, but the precision value in FLOAT(p) is used only to determine storage size. A precision from 0 to 23 results in a 4-byte single-precision FLOAT column. A precision from 24 to 53 results in an 8-byte double-precision DOUBLE column.

On the other hand, the SQLalchemy docs clarify that the precision arg corresponds to the number of decimal places, not binary precision:

Backends should attempt to ensure this precision indicates a number of digits for the generic Float datatype.

So to me it sounds like the correct solution would be to increase that precision by a lot (for example to 53, as suggested here). @charalamm any thoughts?