search

geotiffjs / geotiff.js

geotiff.js is a small library to parse TIFF files for visualization or analysis. It is written in pure JavaScript, and is usable in both the browser and node.js applications.
https://geotiffjs.github.io/
MIT License
861 stars 181 forks source link

GDAL_NODATA value being rounded #259

Closed rooby closed 2 years ago

rooby commented 2 years ago

Generally the GDAL_NODATA metadata is working for me but I have some files that get an incorrect value.

For example with this file, when I run gdalinfo I get a nodata value of 9.99999998050644787e+18 but geotiff.js says 1.0E19 (gdalinfo is correct).

I'm not sure what the difference is with how gdalinfo determines the value, so I'm not sure if it's possible that the nodata metadata is incorrect but gdalinfo gets its info from somewhere else, or if geotiff.js is not reading the value correctly.

$ gdalinfo example.tif

Driver: GTiff/GeoTIFF
Files: example.tif
Size is 279, 416
Coordinate System is:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Origin = (135.033453258464164,-23.596259085710940)
Pixel Size = (0.003889256068923,-0.003789922154798)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_RESOLUTIONUNIT=1 (unitless)
  TIFFTAG_XRESOLUTION=1
  TIFFTAG_YRESOLUTION=1
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  ( 135.0334533, -23.5962591) (135d 2' 0.43"E, 23d35'46.53"S)
Lower Left  ( 135.0334533, -25.1728667) (135d 2' 0.43"E, 25d10'22.32"S)
Upper Right ( 136.1185557, -23.5962591) (136d 7' 6.80"E, 23d35'46.53"S)
Lower Right ( 136.1185557, -25.1728667) (136d 7' 6.80"E, 25d10'22.32"S)
Center      ( 135.5760045, -24.3845629) (135d34'33.62"E, 24d23' 4.43"S)
Band 1 Block=272x416 Type=Float32, ColorInterp=Gray
  NoData Value=9.99999998050644787e+18
constantinius commented 2 years ago

Hi @rooby

Could you share the tiff tags as read from geotiff.js or give us the output of tiffdump? (Or share the original file if possible)

This sounds like a floating point rounding issue to me. Generally, it is not a good idea to compare equality of float values, but always compare with a slight relative +- epsilon.

rooby commented 2 years ago

@constantinius thanks for the quick reply. Here is the output of the tiffdump, and a zipped version of the example.tif is attached.

example.tif:
Magic: 0x4d4d <big-endian> Version: 0x2a <ClassicTIFF>
Directory 0: offset 8 (0x8) next 0 (0)
ImageWidth (256) SHORT (3) 1<279>
ImageLength (257) SHORT (3) 1<416>
BitsPerSample (258) SHORT (3) 1<32>
Compression (259) SHORT (3) 1<1>
Photometric (262) SHORT (3) 1<1>
SamplesPerPixel (277) SHORT (3) 1<1>
XResolution (282) RATIONAL (5) 1<1>
YResolution (283) RATIONAL (5) 1<1>
ResolutionUnit (296) SHORT (3) 1<1>
TileWidth (322) SHORT (3) 1<272>
TileLength (323) SHORT (3) 1<416>
TileOffsets (324) LONG (4) 2<433 453041>
TileByteCounts (325) LONG (4) 2<452608 452608>
SampleFormat (339) SHORT (3) 1<3>
34264 (0x85d8) DOUBLE (12) 16<0.00388926 0 0 135.033 0 -0.00378992 0 -23.5963 0 0 0 0 0 0 0 1>
34735 (0x87af) SHORT (3) 16<1 1 2 3 1024 0 1 2 1025 0 1 1 2048 0 1 4326>
42113 (0xa481) ASCII (2) 21<9.999999980506448E18\0>

example.zip

In this case I think it is maybe an issue with the reading of the value as opposed to a calculation rounding issue. An epsilon value to make these two numbers match would is very large (19493552128).

constantinius commented 2 years ago

@rooby I can't reproduce this issue. I wrote this small script to get the nodata value from the image:

  const tiff = await fromUrl('data/example.tif');
  const image = await tiff.getImage();
  console.log(image.getGDALNoData());

Which gives me this output:

9999999980506448000

This is in the Browser (Chrome).

Which is the same as the exponential notation 9.999999980506448E18. So I now believe that you retrieve the right value from the TIFF, but when it is displayed, the number is somehow rounded.

Can you share the context in which the value is used?

rooby commented 2 years ago

@constantinius Oh sorry, I should have given more context initially. I'm in an Angular app, fetching the data from a WCS GetCoverage request into a Blob and then using fromBlob().

After some more investigation, it seems changing my request to return an array buffer and using fromArrayBuffer doesn't make any difference. When I request from the WCS directly in my browser and downloading the file gives me a file with the correct number, but requesting from the WCS in Angular and saving the blob as a file gives me a file with the incorrect number.

constantinius commented 2 years ago

Have you checked the nodata value from the WCS response with GDAL?

I'm very certain, that geotiff.js is not the culprit on this case, but a wrong nodata value is included in the response.

DanielJDufour commented 2 years ago

I've run into a similar issue as well with other Float32 GeoTIFFs. I'm also not sure if it's an issue with geotiff.js or the program that created the GeoTIFF file (not sure it's GDAL in my case).

Here's a Float32 GeoTIFF of Sea Ice Concentration (originally via https://github.com/GeoTIFF/georaster-layer-for-leaflet/issues/47#issuecomment-716449294) https://github.com/GeoTIFF/test-data/blob/main/files/nt_20201024_f18_nrt_s.tif

A couple issues come to mind:

Hope this helps a little.

constantinius commented 2 years ago

@rooby Is there an update on this? Can you confirm that the WCS result has the same nodata value set as the original image?

rooby commented 2 years ago

@constantinius
I have a workaround for my issue for now. Still doing a bit more investigation to confirm exactly what is the source of the inaccurate number, however I don't think the issue is with geotiff.js. Happy to close this now if you like, otherwise I'll close it off when I've finished the full investigation.

rooby commented 2 years ago

Thanks for the info @DanielJDufour

constantinius commented 2 years ago

I'll close it now. Please re-open if the issue remains.

zakjan commented 2 years ago

I have a similar issue with CMEMS sea_ice_fraction scaled NetCDF. It contains non-sea points with descaled value 1.27999997138977, that I have to override as nodata.

gdal_translate -a_nodata 1.27999997138977 stores it as 1.27999997138977051 string. QGIS shows these points correctly as nodata. JS parses this string as 1.2799999713897705, so strict equality in JS doesn't match the points.

My current workaround:

function maskData(data, nodata = undefined) {
  if (!nodata) {
    return data;
  }

  const maskedData = new data.constructor(Array.from(data).map(value => {
    return Math.abs(value - nodata) > Number.EPSILON * 2 ? value : NaN;
  }));

  return maskedData;
}
zakjan commented 2 years ago

It seems that the issue is with specifically Float32 GeoTIFF bands, because geotiff.js parses GDAL_NODATA into a Float64 number.

How about applying Math.fround for Float32 bands here? https://github.com/geotiffjs/geotiff.js/blob/ae4eb19864879c894ef89e6499d9a08c18b0248f/src/geotiffimage.js#L801

rooby commented 2 years ago

@zakjan your issue seems different to the one I was reporting. My issue is not one relating to inaccuracies in floating point numbers, and the problem occurs before getGDALNoData() ever gets called, since the number parsed from the file header into the fileDirectory.GDAL_NODATA is already incorrect and was either already incorrect in the source file or was incorrectly parsed from the header.

Probably best to create a separate issue for your findings.

zakjan commented 2 years ago

I think my issue is similar to @DanielJDufour's issue. As fileDirectory.GDAL_NODATA is stored as ASCII text, is it possible that it is a GDAL bug (or a "feature")?