#!/bin/bash

GDAL_TILE_PROCESSES=16
GDAL_TILE_ZOOMS=8-14
GDAL_SAMPLING_WARP=cubic
GDAL_SAMPLING_TILE=antialias

# Create virtual dataset with coordinates
gdal_translate -of VRT -a_srs EPSG:4326 -outsize 300% 300% -gcp 6766.419425629241 5085.110020545179 25.27246132 57.30012201 -gcp 38153.9406041717 4816.575791120615 25.25283426 57.30771275 -gcp 37578.814905308434 32195.222845800763 25.26566729 57.31730141 -gcp 650.777214014614 32007.38071992376 25.28906956 57.30762151 GX5131/c_GX_5131_SK_frame_15.jpg c_GX_5131_SK_frame_15.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "0.0 0.0\n41115.0 0.0\n41115.0 32961.0\n0.0 32961.0" | \
  gdaltransform -tps -output_xy c_GX_5131_SK_frame_15.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  c_GX_5131_SK_frame_15.vrt c_GX_5131_SK_frame_15.cut.vrt

# Generate tiles
gdal2tiles.py -r $GDAL_SAMPLING_TILE --xyz -z $GDAL_TILE_ZOOMS -x --processes=$GDAL_TILE_PROCESSES c_GX_5131_SK_frame_15.cut.vrt c_GX_5131_SK_frame_15.xyz