#!/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 2559.75 38905.84477487493 24.837216666666666 52 -gcp 2733.75 1734.7766814897163 26.337216666666666 52 -gcp 22375.5 1854.9230127848805 26.337216666666666 51.5 -gcp 22196.25 39120.11527237354 24.837216666666666 51.5 Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.JPG Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "12540.0 1809.0083287062744\n2731.5 1734.108273181566\n2652.0 14945.995002776235\n2617.5 28176.024153248196\n2558.25 38903.19627984453\n12375.0 39018.0029150472\n22196.25 39122.247917823435\n22225.5 28015.541643531367\n22302.0 14390.719322598556\n22374.75 1852.4958356468628" | \
  gdaltransform -tps -output_xy Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.vrt Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.cut.vrt

# Generate tiles
gdal2tiles.py -r $GDAL_SAMPLING_TILE --xyz -z $GDAL_TILE_ZOOMS -x --processes=$GDAL_TILE_PROCESSES Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.cut.vrt Lietuv.un_Latv.dienvidu_dalas_top.kartes_Vetli-Lubaseva-Serniki_6828_A1_L35_1_02.xyz