On mercredi 29 avril 2020 08:40:07 CEST Guillaume Lostis wrote:
> Like Ciaran I'm not very familiar with GDAL's python bindings, but I've
> tried running a little snippet which I believe does the same thing as what
> was done with `rasterio`: it writes a JP2 to a `vsimem` in-memory file, and
> then writes the bytes to a file on disk.
> from osgeo import gdal
> mem, out = "/vsimem/B01.jp2", "new_B01.jp2"
> gdal.Translate(mem, "B01.jp2")
> # Taken from
> https://lists.osgeo.org/pipermail/gdal-dev/2016-August/045030.html f =
> gdal.VSIFOpenL(mem, "rb")
> gdal.VSIFSeekL(f, 0, 2) # seek to end
> size = gdal.VSIFTellL(f)
> gdal.VSIFSeekL(f, 0, 0) # seek to beginning
> data = gdal.VSIFReadL(1, size, f)
> with open(out, "wb") as ff:
> When I look at the `gdalinfo` of `new_B01.jp2`, it has a CRS and a
> transform, so I am not able to reproduce the behavior we have with
> (Side note: Interestingly enough, the pixel values of `new_B01.jp2` have
> slightly changed with respect to the original `B01.jp2` file, so there is
> some data loss somewhere in the process. But maybe that is expected of the
> JP2 format and could be avoided by passing extra arguments to
By default, the JP2OpenJPEG driver uses a lossy compression (my personal take on JPEG2000 is that it has very limited interested when used in its lossless profile. Better use more conventional algorithms that are way faster).
See https://gdal.org/drivers/raster/jp2openjpeg.html#lossless-compression for options to pass to get a lossless JPEG2000 as output.
So all in all that should be something like
format = "JP2OpenJPEG", # in case you have several JP2K drivers available
options = ["REVERSIBLE=YES", "QUALITY=100"])
Spatialys - Geospatial professional services