#!/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 1105.8309382924767 643.7840088839533 17.9977707 59.99971721 -gcp 12969.445054945054 712.91143808995 23.99781135 59.99981764 -gcp 13573.836010143703 16501.88895058301 23.99804261 55.9997597 -gcp 303.0 16437.34397556913 17.99800626 55.9996635 -gcp 7030.125 805.2176568573016 20.99778799 59.99976766 -gcp 6931.5 16615.2615907829 20.99802173 55.99971182 -gcp 3086.25 726.8820099944476 19.00006449 59.98058872 -gcp 5055.0 786.1482509716825 19.99537736 59.9678195 -gcp 9012.0 803.3982509716825 21.99486759 59.96939855 -gcp 10997.25 774.154497501388 22.99890337 59.98033161 -gcp 13125.0 4656.167823431427 24.00014064 58.95851524 -gcp 13270.5 8604.169072737368 23.99744744 57.94447229 -gcp 13420.5 12552.939061632427 23.99697928 56.95834489 -gcp 11357.25 16574.65019433648 22.99431351 55.97975693 -gcp 9144.75 16614.52679067185 21.99395401 55.96763656 -gcp 4720.5 16588.92476401999 19.99541645 55.96858919 -gcp 2521.5 16518.41143808995 19.00150091 55.98287461 -gcp 497.25 12485.451138256525 17.99487342 56.95862153 -gcp 696.75 8538.137423653527 17.9941802 57.94407481 -gcp 900.375 4588.568295391449 17.99541479 58.95823577 -gcp 4888.5 8691.985285952249 19.99580053 57.91084383 -gcp 6978.0 8717.747640199888 20.99240427 57.90716431 -gcp 9075.0 8710.7702665186 21.99276748 57.91177735 -gcp 7004.625 4760.7307745696835 20.99311661 58.92235105 -gcp 6953.25 12666.12847029428 20.99257516 56.92232528 vlasenko/north/o-34.jpg o-34.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1105.8309382924767 643.7840088839533\n12969.445054945054 712.91143808995\n13573.836010143703 16501.88895058301\n303.0 16437.34397556913" | \
  gdaltransform -tps -output_xy o-34.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  o-34.vrt o-34.cut.vrt

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