I used rio-cogeo to convert a geotiff with 1 band of Type=Float32 geotiff. To my surprise, the output generated was Type=Float64 geotiff. Any idea why that happened and how to keep the data as Type=Float32?

Here's the actual rio-cogeo command I used:

rio cogeo create --cog-profile lzw Southern_California_Topobathy_DEM_1m.tif Southern_California_Topobathy_DEM_1m_COG.tif

And here's the gdalinfo for the INPUT file:

Driver: GTiff/GeoTIFF
Files: Southern_California_Topobathy_DEM_1m.tif
       Southern_California_Topobathy_DEM_1m.tif.ovr
       Southern_California_Topobathy_DEM_1m.tif.aux.xml
Size is 316120, 225217
Coordinate System is:
PROJCS["NAD_1983_NSRS2007_UTM_Zone_11N",
    GEOGCS["GCS_NAD_1983_NSRS2007",
        DATUM["NAD_1983_NSRS2007",
            SPHEROID["GRS_1980",6378137,298.257222101]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-117],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (179523.999999998224666,3824832.000000000000000)
Pixel Size = (1.000000000000000,-1.000000000000000)
Metadata:
  AREA_OR_POINT=Area
  DataType=Generic
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  179524.000, 3824832.000) (120d29'26.74"W, 34d30'55.16"N)
Lower Left  (  179524.000, 3599615.000) (120d24'36.73"W, 32d29'15.49"N)
Upper Right (  495644.000, 3824832.000) (117d 2'50.95"W, 34d33'54.85"N)
Lower Right (  495644.000, 3599615.000) (117d 2'47.00"W, 32d32' 1.94"N)
Center      (  337584.000, 3712223.500) (118d44'57.15"W, 33d32'14.24"N)
Band 1 Block=128x128 Type=Float32, ColorInterp=Gray
  Min=-3857.312 Max=3067.828
  Minimum=-3857.312, Maximum=3067.828, Mean=-353.493, StdDev=804.656
  NoData Value=-3.40282306073709653e+38
  Overviews: 79030x56305, 39515x28153, 19758x14077, 9879x7039, 4940x3520, 2470x1760, 1235x880, 618x440, 309x220
  Metadata:
    RepresentationType=ATHEMATIC
    STATISTICS_COVARIANCES=647471.7255957157
    STATISTICS_MAXIMUM=3067.8276367188
    STATISTICS_MEAN=-353.49345874696
    STATISTICS_MINIMUM=-3857.3120117188
    STATISTICS_SKIPFACTORX=1
    STATISTICS_SKIPFACTORY=1
    STATISTICS_STDDEV=804.65627791978

And the gdalinfo for the OUTPUT file:

 Driver: GTiff/GeoTIFF
Files: tmp2cd_5blz.tif
Size is 316120, 225217
Coordinate System is:
PROJCS["NAD_1983_NSRS2007_UTM_Zone_11N",
    GEOGCS["GCS_NAD_1983_NSRS2007",
        DATUM["NAD_1983_NSRS2007",
            SPHEROID["GRS_1980",6378137,298.257222101]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",-117],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (179523.999999998224666,3824832.000000000000000)
Pixel Size = (1.000000000000000,-1.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  179524.000, 3824832.000) (120d29'26.74"W, 34d30'55.16"N)
Lower Left  (  179524.000, 3599615.000) (120d24'36.73"W, 32d29'15.49"N)
Upper Right (  495644.000, 3824832.000) (117d 2'50.95"W, 34d33'54.85"N)
Lower Right (  495644.000, 3599615.000) (117d 2'47.00"W, 32d32' 1.94"N)
Center      (  337584.000, 3712223.500) (118d44'57.15"W, 33d32'14.24"N)
Band 1 Block=512x512 Type=Float64, ColorInterp=Gray
  NoData Value=-3.40282306073709653e+38
  Overviews: 158060x112609, 79030x56305, 39515x28153, 19758x14077, 9879x7039, 4940x3520, 2470x1760, 1235x880, 618x440

closes #54 #32 #20

This PR refactor the cli and cog_translate to make the internal mask optional.

We updated the documentations to add notes about overview levels, nodata/alpha handling and internal tile size.

Breaking change

This PR adds no main breaking changes:

• internal mask creation is now optional (--add-mask).
• internal nodata or alpha channel can be forwarding to the output dataset

In all sentinel-2 directories, there is a PVI (preview) file in jp2 format that throws a blocksize exceeds raster width error when you try to convert using rio cogeo create [file].

gdalinfo confirms that the blocksize is 192x192 while the band dimension is just 171x171. So the issue seems to be more of a corner case of really small (tiled) jp2s rather than a bug with rio cogeo.

If in SAFE format, the approximate location of this file is: $(PRODUCT_NAME).SAFE/GRANULE/$(IMAGE_GRANULE)/QI_DATA/TILE_ID..._PVI.jp2

Here's a command to pull a random preview file out of the sentinel s3: aws s3 cp s3://sentinel-s2-l1c/tiles/7/J/FL/2019/2/12/0/preview.jp2 . --request-payer requester

As you can see, rio cogeo create preview.jp2 preview.tif fails.

It looks like this generic error comes when instantiating the DatasetWriterBase class. Is there a way to keep current protections but allow for cog creation on even really small preview files?

Thanks in advance for any suggestions or advice :)

Not sure whether this is the right place to report this but it seems relevant as it can be reproduced using rio cogeo create and rio cogeo validate.

Use the files in validate.zip to reproduce:

rio cogeo create temp.tif cog.tif && rio cogeo validate cog.tif with GDAL 3.0.4 creates a COG and validation outputs this:

Reading input: /temp.tif
Adding overviews...
Updating dataset tags...
Writing output to: /cog.tif
/cog.tif is a valid cloud optimized GeoTIFF

whereas running the same command with GDAL 3.1.0 (also tested 3.1.1 and 3.1.2 with the same result) makes rio cogeo validate return an error:

Reading input: /temp.tif
Adding overviews...
Updating dataset tags...
Writing output to: /cog.tif
The following errors were found:
- The offset of the main IFD should be 8 for ClassicTIFF or 16 for BigTIFF. It is 192 instead
/cog.tif is NOT a valid cloud optimized GeoTIFF

Again, I guess the error rather has its origin in either GDAL itself or rasterio, however I think it is worth at least documenting it here as it affects COG creation using rio-cogeo.

I'll investigate further on the GDAL and rasterio issue trackers but do you have a clue where the issue could be?

I am trying to use cog_translate to convert some netCDF4 datasets to COG, and it appears these files do not have proper georeferencing.

I have modified cogeo.py to accept vrt_params, so that I can manually pass the src_crs and src_transform via the WarpedVRTReaderBase mixin of WarpedVRT, as you can see here:

https://github.com/ryanjdillon/rio-cogeo/blob/cog_translate_vrt_params/rio_cogeo/cogeo.py#L87-L99

Then I'm passing the following vrt_params:

vrt_params = dict(src_transform=Affine(gdal_geotransform), src_crs='EPSG:4326')

This get's me past the CRS is None error, but it doesn't seem to acknowledge the Affine transform I pass to src_transform, resulting in the following error:

CPLE_AppDefinedError: The transformation is already "north up" or a transformation between pixel/line and georeferenced coordinates cannot be computed for /home/ryan/data/webstep/kvt/NVE_Sorlandet/NVE_Sorlandet_wslev0_xyt_synth.nc4. There is no affine transformation and no GCPs. Specify transformation option SRC_METHOD=NO_GEOTRANSFORM to bypass this check.

Any suggestions on how I might get this to work?

Because of #127, #132 and #126 we have to cut a new major release.

https://github.com/cogeotiff/rio-cogeo/blob/25f9a72f0d18b68170863d14c8690b8a550648e2/CHANGES.txt#L1-L7

changes are now in https://github.com/cogeotiff/rio-cogeo/tree/v2

Hello :hand:

I propose the addition of a "strict" flag to threat warnings as errors for the purposes of validation.

There are situations in which this might be desirable (eg. trying to determine whether a GeoTIFF is suited for streaming, overviews lacking are warnings, but affect performance significantly).

After building a COG from a vrt:

rio cogeo create my.vrt my.tif -p jpeg

the COG doesn't validate:

> rio cogeo validate /vsis3/bucket/my.tif
The following errors were found:
- The offset of the IFD for overview of index 1 is 884184, whereas it should be greater than the one of index 0, which is at byte 3167886
- The offset of the IFD for overview of index 2 is 299436, whereas it should be greater than the one of index 1, which is at byte 884184
- The offset of the IFD for overview of index 3 is 149834, whereas it should be greater than the one of index 2, which is at byte 299436
- The offset of the IFD for overview of index 4 is 107736, whereas it should be greater than the one of index 3, which is at byte 149834
- The offset of the IFD for overview of index 5 is 94524, whereas it should be greater than the one of index 4, which is at byte 107736
- The offset of the first block of the smallest overview should be after its IFD
/vsis3/bucket/my.tif is NOT a valid cloud optimized GeoTIFF

The VRT points to a group of 256px tiffs stored in S3.

I would expect to get the same set of mask_flags for Case 4 as for the rest of the cases. Am I missing something or it is a bug?

$ rio info original.tif | jq -c .mask_flags,.nodata
[["nodata"],["nodata"],["nodata"]]
0

:green_circle: Case 1

--add-mask

$ rio cogeo create --add-mask original.tif cogeo1.tif
$ rio info cogeo1.tif | jq -c .mask_flags
[["per_dataset"],["per_dataset"],["per_dataset"]]

:green_circle: Case 2

--add-mask --web-optimized

$ rio cogeo create --add-mask --web-optimized original.tif cogeo2.tif
$ rio info cogeo2.tif | jq -c .mask_flags                           
[["per_dataset"],["per_dataset"],["per_dataset"]]

:green_circle: Case 3

--add-mask --use-cog-driver

$ rio cogeo create --add-mask --use-cog-driver original.tif cogeo3.tif
$ rio info cogeo3.tif | jq -c .mask_flags  
[["per_dataset"],["per_dataset"],["per_dataset"]]

:red_circle: Case 4

--add-mask --use-cog-driver --web-optimized

$ rio cogeo create --add-mask --use-cog-driver --web-optimized original.tif cogeo4.tif
$ rio info cogeo4.tif | jq -c .mask_flags 
[["per_dataset","alpha"],["per_dataset","alpha"],["per_dataset","alpha"],["all_valid"]]

I'm opening this to document my findings, I'm not sure it's even a problem with rio-cogeo per-se, it seems to go down all the way to rasterio (1.1.0) and perhaps GDAL (2.4.3).

Input (5 bands + alpha 16bit GeoTIFF): odm_orthophoto.zip

Band 1 Block=512x512 Type=UInt16, ColorInterp=Red
Band 2 Block=512x512 Type=UInt16, ColorInterp=Green
Band 3 Block=512x512 Type=UInt16, ColorInterp=Blue
Band 4 Block=512x512 Type=UInt16, ColorInterp=Gray
Band 5 Block=512x512 Type=UInt16, ColorInterp=Gray
Band 6 Block=512x512 Type=UInt16, ColorInterp=Alpha

Command:

rio cogeo create odm_orthophoto.tif out.tif

Output:

Band 1 Block=512x512 Type=UInt16, ColorInterp=Gray  Overviews: 1000x845, 500x423
Band 2 Block=512x512 Type=UInt16, ColorInterp=Undefined  Overviews: 1000x845, 500x423
Band 3 Block=512x512 Type=UInt16, ColorInterp=Undefined  Overviews: 1000x845, 500x423
Band 4 Block=512x512 Type=UInt16, ColorInterp=Undefined  Overviews: 1000x845, 500x423
Band 5 Block=512x512 Type=UInt16, ColorInterp=Undefined  Overviews: 1000x845, 500x423
Band 6 Block=512x512 Type=UInt16, ColorInterp=Undefined  Overviews: 1000x845, 500x423

Raster values seem to have been downscaled (?) or overflowed (?) in the process:

Alpha band is out of whack too, but I won't worry about that for now.

I've narrowed this down to https://github.com/cogeotiff/rio-cogeo/blob/master/rio_cogeo/cogeo.py#L210. Once the WarpedVRT object is created, read operations on the virtual raster return erroneous values (they seem scaled down, eg. 24933 --> 97, which seems linear from 65535 to 255, but other values are still in the 65535 range, so it makes little sense). The data types are correct, although we've lost ColorInterp information (perhaps another problem).

This might be a problem with rasterio, or even GDAL? If I find a solution I will open a PR.

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}

I've a couple of rasters of pixel dimensions like 300k x 150k. Most of the pixels have NODATA values. The data looks like below

The file size with LZW compression is 7 GB. When I use following GDAL commands to convert to COG, conversion happens in around an hour:

$ gdal_translate source.tif tmp.tif -co TILED=YES -co COMPRESS=DEFLATE -co BIGTIFF=IF_SAFER -co COPY_SRC_OVERVIEWS=YES
$ gdaladdo tmp.tif 2 4 8 16 32
$ gdal_translate tmp.tif cog.tif -co TILED=YES -co BIGTIFF=IF_SAFER -co COPY_SRC_OVERVIEWS=YES

But if I run rio cogeo create source.tif cog.tif, conversion gets stuck at -

Reading input: <file path>
  [------------------------------------]    0%  2d 19:51:43

I tried multiple times and killed the process after couple of hours.

I believe the method vrt_dst.read() of class WarpedVRT is trying to unpack the raster (which is originally in LZW compression). Uncompressed file would come around 200 - 300 GB. This could be causing the conversion to take days instead of an hour.