Search code examples
pythonrasteriojpeg2000

How do I mount these JPEG2000 R, G, B files, with rasterio?

So my problem consists of the trying of getting 3 different .jp2 files - one corresponding to each band R, G, B - and mounting them as a solo RGB file. I'm using rasterio for the issue, and have followed the following tutorial: https://towardsdatascience.com/satellite-imagery-access-and-analysis-in-python-jupyter-notebooks-387971ece84b

My code got to this:

import os
import geopandas
import rasterio as rio
from PIL import Image

Image.MAX_IMAGE_PIXELS = 300000000

## RGB Transform

# Image paths: B04 = Red, B03 = Green, B02 = Blue
redPath = "assets/geodata/dump_r2/S2A_MSIL2A_20220512T133231_N0400_R081_T22JCP_20220512T202012.SAFE/GRANULE/L2A_T22JCP_A035971_20220512T133948/IMG_DATA/R10m/T22JCP_20220512T133231_B04_10m.jp2"
greenPath = "assets/geodata/dump_r2/S2A_MSIL2A_20220512T133231_N0400_R081_T22JCP_20220512T202012.SAFE/GRANULE/L2A_T22JCP_A035971_20220512T133948/IMG_DATA/R10m/T22JCP_20220512T133231_B03_10m.jp2"
bluePath = "assets/geodata/dump_r2/S2A_MSIL2A_20220512T133231_N0400_R081_T22JCP_20220512T202012.SAFE/GRANULE/L2A_T22JCP_A035971_20220512T133948/IMG_DATA/R10m/T22JCP_20220512T133231_B02_10m.jp2"

## Bands

"""
redConv = Image.open(redPath)
redConv.save("red.tiff", "TIFF")
greenConv = Image.open(greenPath)
greenConv.save("green.tiff", "TIFF")
blueConv = Image.open(bluePath)
blueConv.save("blue.tiff", "TIFF")
"""
red = rio.open(redPath)
green = rio.open(greenPath)
blue = rio.open(bluePath) 

# Creates RGB file.
rgb = rio.open("RGB.tiff", 'w+', driver="Gtiff", width=red.width, height=red.height, count=3, crs=red.crs, transform=red.transform, dtype=red.dtypes[0])
rgb.write(blue.read(1), 1)
rgb.write(green.read(1), 2)
rgb.write(red.read(1), 3)

print(red.shape)
print(green.shape)
print(blue.shape)

rgb.close()

The result however seems to be just gray-scale images, with no channels in reality. My think-processes were I. There's something wrong with the drivers, be it the JP2-TIFF conversion or the JP2 reading or the .tiff saving. II. I don't have enough imagery knowledge, I try opening the "rgb" file through GIMP and it only shows gray-scale images, perhaps there's extra steps I need to take, but as far as I reached, nothing. III. The rgb.read().shape is [3, 10980, 10980] (like "3 matrices"), I've read other file types (like jpeg) with like matplotlib.image, but it gives a different shape, something like (1980, 1080, 3) ("1 matrix with 3 values in each cell"). But I suppose that's just how .tiff works? Are .tiff even supposed to be visualized? enter image description here Image opened by gimp

Meta-data of the image (exiftool)

ExifTool Version Number         : 12.16
File Name                       : RGB.tiff
Directory                       : .
File Size                       : 690 MiB
File Modification Date/Time     : 2022:08:24 11:56:03-03:00
File Access Date/Time           : 2022:08:24 11:55:18-03:00
File Inode Change Date/Time     : 2022:08:24 11:56:03-03:00
File Permissions                : rw-r--r--
File Type                       : TIFF
File Type Extension             : tif
MIME Type                       : image/tiff
Exif Byte Order                 : Little-endian (Intel, II)
Image Width                     : 10980
Image Height                    : 10980
Bits Per Sample                 : 16 16 16
Compression                     : Uncompressed
Photometric Interpretation      : BlackIsZero
Strip Offsets                   : (Binary data 108116 bytes, use -b option to extract)
Samples Per Pixel               : 3
Rows Per Strip                  : 1
Strip Byte Counts               : (Binary data 65879 bytes, use -b option to extract)
Planar Configuration            : Chunky
Extra Samples                   : Unknown (0 0)
Sample Format                   : Unsigned; Unsigned; Unsigned
Pixel Scale                     : 10 10 0
Model Tie Point                 : 0 0 0 300000 6900040 0
Geo Tiff Version                : 1.1.0
GT Model Type                   : Projected
GT Raster Type                  : Pixel Is Area
GT Citation                     : WGS 84 / UTM zone 22S
Geog Citation                   : WGS 84
Geog Angular Units              : Angular Degree
Projected CS Type               : WGS84 UTM zone 22S
Proj Linear Units               : Linear Meter
Image Size                      : 10980x10980
Megapixels                      : 120.6
Solution

  • The "Photometric Interpretation" looks wrong with "BlackIsZero".

    Try changing to:

    rgb = rio.open("RGB.tiff", 'w+', ..., photometric="RGB")

    and then write the bands in RGB order - rather than BGR.

    Alternatively, there is a tool called tiffset that is installed alongside libtiff and you could try changing the "Photometric Interpretation" from MinIsBlack to RGB on your rasterio-generated image like this:

    tiffset -s 262 2 RGB.TIF