#!/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 0.0 0.0 24.48869705 57.07615487 -gcp 10800.0 0.0 24.49616432 56.99108267 -gcp 10800.0 10800.0 24.33621883 56.98591556 -gcp 0.0 10800.0 24.32780743 57.07295913 GX220/GX_220_SK_frame_560.jpg GX_220_SK_frame_560.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "252.0 5883.0\n376.5 5754.0\n376.5 5056.5\n246.0 4923.0\n205.5 1038.0\n1143.0 1021.5\n1140.0 525.0\n3483.0 561.0\n3534.0 150.0\n4992.0 126.0\n5136.0 258.0\n5940.0 246.0\n8676.0 246.0\n8700.0 486.0\n10764.0 462.0\n10755.0 4929.0\n10662.0 5020.5\n10663.5 5712.0\n10764.0 5808.0\n10771.5 9618.0\n10426.5 9648.0\n10260.0 9664.5\n10146.0 9733.5\n10035.0 9846.0\n9972.0 9984.0\n9963.0 10144.5\n9967.5 10618.5\n6031.5 10635.0\n5932.5 10515.0\n5224.5 10519.5\n5118.0 10647.0\n1305.0 10635.0\n1284.0 10143.0\n1270.5 10014.0\n1219.5 9904.5\n1149.0 9786.0\n1035.0 9696.0\n913.5 9637.5\n763.5 9616.5\n261.0 9639.0" | \
  gdaltransform -tps -output_xy GX_220_SK_frame_560.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  GX_220_SK_frame_560.vrt GX_220_SK_frame_560.cut.vrt

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