when using jpeg profile silently creates empty output if bit depth != 8

[alon]

Vrersion:
$ pip freeze | grep rio
rasterio==1.0.22
rio-cogeo==1.0b3

Problem: When creating from a unsigned int 16 bit per pixel tif a jpeg profile / jpeg compressed cogeo, successfully creates an empty file.

What I expected: an error message indicating this can only work for 8 bits per pixel.

Details:

When doing the same thing with gdal directly I get warnings that JPEG compression is not compatible with 16 bit pixel depth (using UINT16 Tiff bands), but rio cogeo silently swallows these warnings:

gdal_translate temp.tif gdal_cog.tif -co TILED=YES -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co COPY_SRC_OVERVIEWS=YES
Input file size is 4171, 4061
0ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
..ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
.10...20.ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
..30...40...50...60...70...80...90...100 - done.
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG
ERROR 1: WriteEncodedTile/Strip() failed.
ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG

 $ rio cogeo create --cog-profile jpeg temp.tif cog_jpeg.tif
Reading input: /home/hiro/git/tilingtest/temp.tif
  [####################################]  100%
Adding overviews...
Updating dataset tags...
Writing output to: /home/hiro/git/tilingtest/cog_jpeg.tif
$ ls -l cog_jpeg.tif 
-rw-------. 1 hiro hiro 4992 Mar 31 13:02 cog_jpeg.tif
$ ls -l temp.tif
-rw-------. 1 hiro hiro 138236380 Mar 31 12:32 temp.tif
$ rio info temp.tif
{"bounds": [694696.0, 3674712.0, 736406.0, 3715322.0], "colorinterp": ["red", "green", "blue"], "count": 3, "crs": "EPSG:32636", "descriptions": [null, null, null], "driver": "GTiff", "dtype": "uint16", "height": 4061, "indexes": [1, 2, 3], "interleave": "pixel", "lnglat": [35.316978660693565, 33.37281747525934], "mask_flags": [["all_valid"], ["all_valid"], ["all_valid"]], "nodata": null, "res": [10.0, 10.0], "shape": [4061, 4171], "tiled": false, "transform": [10.0, 0.0, 694696.0, 0.0, -10.0, 3715322.0, 0.0, 0.0, 1.0], "units": [null, null, null], "width": 4171}

[vincentsarago]

Thanks @alon this is interesting. I'll have a look at it

[alon]

Thank you. If I may, I have a very slow tiling work flow I want to replace with a COG. Is it possible to serve a COG without a tiling service, i.e. have leaflet or such directly GE it,, like vsicurl w/ gdal would? That is the goal, right?

…

On March 31, 2019 7:54:27 PM GMT+03:00, Vincent Sarago ***@***.***> wrote: Thanks @alon this is interesting. I'll have a look at it -- You are receiving this because you were mentioned. Reply to this email directly or view it on GitHub: #73 (comment)

-- Sent from my Android device with K-9 Mail. Please excuse my brevity.

[vincentsarago]

@alon

Is it possible to serve a COG without a tiling service

Absolutely, I'm usually using https://github.com/cogeotiff/rio-tiler (and https://github.com/vincentsarago/lambda-tiler) to do this but you could also check https://github.com/DHI-GRAS/terracotta or https://github.com/mojodna/marblecutter which provide the full architecture.

That's said, if you have int16 dataset you'll need custom script to rescale your data to int8.

[vincentsarago]

@alon This is a rasterio problem and I'm not sure if there is a way around it.

Basically rasterio.shutil.copy https://github.com/mapbox/rasterio/blob/master/rasterio/shutil.pyx#L56 is swallowing GDAL warnings but I don't think there is a way to fix that.

cc @sgillies

[sgillies]

@vincentsarago the errors occur during the execution of GDALCreateCopy at https://github.com/mapbox/rasterio/blob/master/rasterio/shutil.pyx#L137, right? Because we've released Python's global interpreter lock, the errors cannot be turned into exceptions. It might be possible to create a GDAL error handler that would break out of GDALCreateCopy and the nogil block using setjmp and longjmp, but I've never done that in Cython before and I'm hesitant to do so 😰 Validating the data type before copying is probably the best solution.

[vincentsarago]

@sgillies thanks for the investigation.
That's said, I think GDAL doesn't throw error but just warnings (ERROR 1: JPEGSetupEncode:BitsPerSample 16 not allowed for JPEG)

Validating the data type before copying is probably the best solution

Yeah I think this is what we should do 🙏

[vincentsarago]

🤔 @sgillies I've just encountered this with other options

# This command runs find but the result is empty (data is not copied)
$ rio cogeo create float64.tif float64cog.tif --co BLOCKXSIZE=256 --co BLOCKXSIZE=256 --cog-profile deflate --co ZLEVEL=6 --co PREDICTOR=2

while with gdal 
gdal_translate  float64.tif float64cog.tif -co BLOCKXSIZE=256 -co BLOCKYSIZE=256 -co COMPRESS=DEFLATE -co ZLEVEL=6 -co PREDICTOR=2
Input file size is 301, 520
0...10...20...30...40...50...60...70...80...90...100 - done.
ERROR 1: PredictorSetup:Horizontal differencing "Predictor" not supported with 64-bit samples
ERROR 1: PredictorSetup:Horizontal differencing "Predictor" not supported with 64-bit samples
ERROR 1: WriteEncodedTile/Strip() failed.
ERROR 1: PredictorSetup:Horizontal differencing "Predictor" not supported with 64-bit samples
ERROR 1: WriteEncodedTile/Strip() failed.
ERROR 1: PredictorSetup:Horizontal differencing "Predictor" not supported with 64-bit samples

I feel there are a lot of corner case I'll need to check ...

[alon]

@alon

Is it possible to serve a COG without a tiling service

Absolutely, I'm usually using https://github.com/cogeotiff/rio-tiler (and https://github.com/vincentsarago/lambda-tiler) to do this but you could also check https://github.com/DHI-GRAS/terracotta or https://github.com/mojodna/marblecutter which provide the full architecture.

Right, I understand it is doable with a service which creates the tiles on demand, and then serves them:
COG -> TILER -> FRONTEND
And showing the artifacts:
COG: single file on disk
TILER: in memory / cached on disk tiles (png per X/Y/Z)

But I'd like to forgoe the tiler, and have
COG -> FRONTEND
where the arrow signifies HTTP traffic.

My current flow is:
Preprocessing: TIFF (not COG) -> tiles -> cloud (azure)
CLOUD -> FRONTEND (via leaflet)
The tiling stage is very slow, I'm using gdal2tiles with zoom levels 18, which is brutal for a 50k x 20k image, even though most of it is NODATA.

@vincentsarago I would happily discuss this offline, or publically but on a more appropriate place than this ticket :) I'm alon at pobox dotcom

[vincentsarago]

@alon

COG -> FRONTEND

checkout https://github.com/geotiffjs/geotiff.js

I would happily discuss this offline, or publically but on a more appropriate place than this ticket :) I'm alon at pobox dotcom

We have a #cogeotiff Chanel on thespatialcommunity slack -> https://thespatialcommunity.org