#!/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 25.49759388 56.29176372 -gcp 21690.0 0.0 25.45915016 56.31081931 -gcp 21690.0 21123.0 25.49331039 56.33226602 -gcp 0.0 21123.0 25.53189397 56.31205592 GX751/GX_751_SD_frame_118.jpg GX_751_SD_frame_118.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "354.0 11643.0\n489.0 11472.0\n720.0 11331.0\n726.0 9990.0\n480.0 9819.0\n369.0 9633.0\n429.0 1785.0\n2359.5 1792.5\n2412.0 55.5\n6228.0 109.5\n10071.0 129.0\n10228.5 159.0\n10306.5 223.5\n10401.0 441.0\n10828.5 511.5\n11368.5 501.0\n11764.5 469.5\n11896.5 228.0\n11986.5 175.5\n12078.0 138.0\n15819.0 180.0\n19248.0 198.0\n21480.0 210.0\n21583.5 274.5\n21661.5 382.5\n21622.5 4890.0\n21579.0 9777.0\n21543.0 9954.0\n21436.5 10038.0\n21240.0 10141.5\n21207.0 10519.5\n21199.5 11083.5\n21210.0 11479.5\n21297.0 11572.5\n21424.5 11619.0\n21513.0 11725.5\n21540.0 11958.0\n21475.5 21120.0\n2355.0 21120.0\n2400.0 20175.0\n2352.0 19926.0\n2247.0 19683.0\n1977.0 19362.0\n1731.0 19236.0\n1554.0 19167.0\n1356.0 19140.0\n276.0 19218.0" | \
  gdaltransform -tps -output_xy GX_751_SD_frame_118.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  GX_751_SD_frame_118.vrt GX_751_SD_frame_118.cut.vrt

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