#!/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 8203.5 3232.6951083935523 25.837216666666666 56.5 -gcp 19150.5 3424.567954419122 26.337216666666666 56.5 -gcp 19081.5 33641.05364091162 26.337216666666666 55.75 -gcp 7978.5 33631.83185102835 25.837216666666666 55.75 Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.JPG Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "10426.5 3271.5041689827685\n1704.0 3072.0\n1279.5 33568.49916712937\n10177.5 33637.51834352418\n19073.25 33640.49916712937\n19007.25 23461.498334258744\n18942.0 18546.00666296502\n19008.0 13605.089950027761\n19150.5 3423.0" | \
  gdaltransform -tps -output_xy Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.vrt Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.cut.vrt

# Generate tiles
gdal2tiles.py -r $GDAL_SAMPLING_TILE --xyz -z $GDAL_TILE_ZOOMS -x --processes=$GDAL_TILE_PROCESSES Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.cut.vrt Liet.un_Latv.dienvidu_dalas_top.kartes_Lielzalve-Birzi-Jekabp._6828_A1_L35_3_03.xyz