#!/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 1793.6385246750715 5127.48659110771 25.15238285 57.26702455 -gcp 19305.7853262072 4969.005777267554 25.1530695 57.33653993 -gcp 16696.327226350768 19131.370758708297 25.2584374 57.32487665 -gcp 3696.9196206570314 18281.393257723787 25.24918914 57.27353317 GX8464/GX_8464_SK_frame_279.jpg GX_8464_SK_frame_279.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "720.0 2783.893036211699\n768.0 19632.722005571028\n20736.0 19344.294150417823\n21120.0 2400.080222841226" | \
  gdaltransform -tps -output_xy GX_8464_SK_frame_279.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  GX_8464_SK_frame_279.vrt GX_8464_SK_frame_279.cut.vrt

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