Huge GeoTIFF file after reprojecting file

Question

I have some geotiff file. Their CRS ise EPSG:32635. I want to change crs to EPSG:4326. I use below script from Rasteri doc (link). In addition I added 'compress': 'lzw' to compress new geotiff file. But my output file size is huge. For example The size of a 12 MB file is 945 MB after conversion.

How can I solve this problem?

import numpy as np
import rasterio
from rasterio.warp import calculate_default_transform, reproject, Resampling

dst_crs = 'EPSG:4326'

with rasterio.open('rasterio/tests/data/RGB.byte.tif') as src:
    transform, width, height = calculate_default_transform(
        src.crs, dst_crs, src.width, src.height, *src.bounds)
    kwargs = src.meta.copy()
    kwargs.update({
        'crs': dst_crs,
        'transform': transform,
        'width': width,
        'height': height,
        'compress': 'lzw'
    })

    with rasterio.open('/tmp/RGB.byte.wgs84.tif', 'w', **kwargs) as dst:
        for i in range(1, src.count + 1):
            reproject(
                source=rasterio.band(src, i),
                destination=rasterio.band(dst, i),
                src_transform=src.transform,
                src_crs=src.crs,
                dst_transform=transform,
                dst_crs=dst_crs,
                resampling=Resampling.nearest)

The output {i: dtype for i, dtype in zip(dataset.indexes, dataset.dtypes)} is below.

Original img output 
{1: 'uint8', 2: 'uint8', 3: 'uint8'} 
------------- 
Reproject img output 
{1: 'uint8', 2: 'uint8', 3: 'uint8'} 

gdalinfo for original image

Size is 13691, 16259 Coordinate System is: PROJCS["WGS 84 / UTM zone 35N", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0, AUTHORITY["EPSG","8901"]], UNIT["degree",0.0174532925199433, AUTHORITY["EPSG","9122"]], AUTHORITY["EPSG","4326"]], PROJECTION["Transverse_Mercator"], PARAMETER["latitude_of_origin",0], PARAMETER["central_meridian",27], PARAMETER["scale_factor",0.9996], PARAMETER["false_easting",500000], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001"]], AXIS["Easting",EAST], AXIS["Northing",NORTH], AUTHORITY["EPSG","32635"]] Origin = (543529.770710000069812,5977099.678410000167787) Pixel Size = (0.065900000000000,-0.065900000000000) Metadata: AREA_OR_POINT=Area TIFFTAG_SOFTWARE=pix4dmapper Image Structure Metadata: COMPRESSION=YCbCr JPEG INTERLEAVE=PIXEL SOURCE_COLOR_SPACE=YCbCr ERROR 1: ***/oam_data/aerial_imgs/01.tif: Can not read TIFF directory count ERROR 1: TIFFReadDirectory:Failed to read directory at offset 15655350 Corner Coordinates: Upper Left ( 543529.771, 5977099.678) ( 27d39'47.28"E, 53d56'25.62"N) Lower Left ( 543529.771, 5976028.210) ( 27d39'46.73"E, 53d55'50.95"N) Upper Right ( 544432.008, 5977099.678) ( 27d40'36.76"E, 53d56'25.34"N) Lower Right ( 544432.008, 5976028.210) ( 27d40'36.19"E, 53d55'50.67"N) Center ( 543980.889, 5976563.944) ( 27d40'11.74"E, 53d56' 8.15"N) Band 1 Block=512x512 Type=Byte, ColorInterp=Red NoData Value=-10000 Band 2 Block=512x512 Type=Byte, ColorInterp=Green NoData Value=-10000 Band 3 Block=512x512 Type=Byte, ColorInterp=Blue NoData Value=-10000

gdalinfo for reproject image

Driver: GTiff/GeoTIFF Files: ***/oam_data/aerial_imgs/01_crs_changed.tif Size is 17689, 12355 Coordinate System is: GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]] Origin = (27.662980139607967,53.940449109530391) Pixel Size = (0.000000785608906,-0.000000785608906) Metadata: AREA_OR_POINT=Area Image Structure Metadata: COMPRESSION=LZW INTERLEAVE=PIXEL Corner Coordinates: Upper Left ( 27.6629801, 53.9404491) ( 27d39'46.73"E, 53d56'25.62"N) Lower Left ( 27.6629801, 53.9307429) ( 27d39'46.73"E, 53d55'50.67"N) Upper Right ( 27.6768768, 53.9404491) ( 27d40'36.76"E, 53d56'25.62"N) Lower Right ( 27.6768768, 53.9307429) ( 27d40'36.76"E, 53d55'50.67"N) Center ( 27.6699285, 53.9355960) ( 27d40'11.74"E, 53d56' 8.15"N) Band 1 Block=17689x1 Type=Byte, ColorInterp=Red Band 2 Block=17689x1 Type=Byte, ColorInterp=Green Band 3 Block=17689x1 Type=Byte, ColorInterp=Blue

When I try with GDAL using below code the original image is 12 MB and reproject image is 655 MB.

def gdal_reproject(): import gdal gdal.Warp('***/reproject.tif', '***/original.tif', dstSRS='EPSG:4326')

And this is the output error messages

ERROR 1: TIFFFetchDirectory:***/aerial_imgs/original.tif: Can not read TIFF directory count ERROR 1: TIFFReadDirectory:Failed to read directory at offset 15655350 ERROR 1: TIFFFillTile:Read error at row 2048, col 8192, tile 648; got 0 bytes, expected 1664 ERROR 1: TIFFReadEncodedTile() failed. ERROR 1: ***/original.tif, band 1: IReadBlock failed at X offset 0, Y offset 24: TIFFReadEncodedTile() failed.

When I use "compress": "JPEG", photometric": "YCBCR" parameters instead of 'compress': 'lzw' the file size has shrunk.The size of original and reproject image almost same.

Answer 1

Your input tif is lossily compressed with JPEG YCBCR which gives much higher compression than lossless LZW. Try "compress": "JPEG", "photometric": "YCBCR" instead.

Answer 2

import rasterio
from rasterio.warp import calculate_default_transform, reproject, Resampling

dst_crs = 'EPSG:4326'

with rasterio.open('rasterio/tests/data/RGB.byte.tif') as src:
    transform, width, height = calculate_default_transform(
        src.crs, dst_crs, src.width, src.height, *src.bounds)
    kwargs = src.meta.copy()
    kwargs.update({
        'crs': dst_crs,
        'transform': transform,
        'width': width,
        'height': height,
        'compress': 'lzw',
         'dtype': 'uint8'
    })

    with rasterio.open('/tmp/RGB.byte.wgs84.tif', 'w', **kwargs) as dst:
        for i in range(1, src.count + 1):
            reproject(
                source=rasterio.band(src, i),
                destination=rasterio.band(dst, i),
                src_transform=src.transform,
                src_crs=src.crs,
                dst_transform=transform,
                dst_crs=dst_crs,
                resampling=Resampling.nearest)

Sometimes it can happen that when warping a raster dataset to other projection, end up having pixels within the new extent, and they will be assigned to no data to keep the data values. Does that solve your problem? Or do you have the same results?