#!/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 3206.9542916324554 3765.057775586555 23.81251216 56.77050605 -gcp 18969.665106106127 6210.381350409725 23.82552624 56.83369119 -gcp 18935.951442721103 19760.76437017932 23.92384529 56.83616205 -gcp 1302.677902355515 19290.644329337912 23.92745554 56.76566132 GX10673/GX_10673_SK_frame_011.jpg GX_10673_SK_frame_011.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "864.0 3465.0\n20544.0 3360.4935732647828\n20256.0 20400.123393316193\n672.0 20304.691002570697" | \
  gdaltransform -tps -output_xy GX_10673_SK_frame_011.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  GX_10673_SK_frame_011.vrt GX_10673_SK_frame_011.cut.vrt

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