#!/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 1213.0228680799362 1124.2354603447734 25.48350692 57.38933863 -gcp 9597.358364175907 639.8714272164166 25.49692869 57.31723006 -gcp 8872.025056352213 9858.320819896744 25.353688 57.31874196 -gcp 2227.3634790277433 9734.21571635263 25.34775496 57.37641988 GX8473/GX_8473_SK_frame_203_copy.jpg GX_8473_SK_frame_203_copy.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1122.0 471.0016713091922\n1104.0 1170.0\n657.0 1172.9966573816155\n612.0 2879.9949860724237\n242.7271728515625 3311.9998168945312\n351.2728271484375 9274.908782958984\n981.0 9266.984958217272\n1275.0 9572.998328690806\n1257.0 10133.998328690806\n9141.0 10128.0\n9141.0 9629.978272980501\n9439.090576171875 9308.727172851562\n9939.272827148438 9293.454711914062\n9960.0 336.02005571030645" | \
  gdaltransform -tps -output_xy GX_8473_SK_frame_203_copy.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  GX_8473_SK_frame_203_copy.vrt GX_8473_SK_frame_203_copy.cut.vrt

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