#!/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 1277.1231511254018 582.6946141032761 35.99817902 55.99994419 -gcp 14550.654019292606 762.4730705163797 41.99827759 56.00003046 -gcp 14985.65498392283 16548.38228761799 41.99843522 51.99996004 -gcp 379.3634965572496 16354.88777880364 35.99834568 51.99987799 -gcp 3489.0 683.3953359244863 36.99918271 55.98117519 -gcp 5700.0 761.4044975013882 37.99868013 55.96843113 -gcp 7909.5 804.2956690727374 38.9967408 55.96500315 -gcp 10121.625 818.6932953914493 39.9958983 55.96914685 -gcp 14652.0 4717.060244308717 41.99495446 54.9620718 -gcp 14759.625 8658.117642976125 41.99465469 53.95169678 -gcp 14870.25 12603.254441976678 41.99558368 52.96406713 -gcp 12548.25 16603.0719044975 40.99646623 51.97894568 -gcp 10112.25 16614.18739589117 39.99616206 51.96846228 -gcp 7681.5 16599.186563020543 38.99844858 51.96430001 -gcp 5251.5 16551.94280955025 38.00095324 51.96791975 -gcp 2817.0 16480.669072737368 37.00110866 51.97732374 -gcp 598.5 12410.396585230428 35.99604527 52.96466088 -gcp 827.625 8477.27977512493 35.99820831 53.95014531 -gcp 1051.875 4526.743475846752 35.99789281 54.96383529 -gcp 7849.5 4763.499444752915 38.99594551 54.92622 -gcp 7789.5 8710.866879511383 38.99511342 53.91464336 -gcp 10111.5 8734.1568573015 39.99370268 53.91663448 -gcp 5474.25 8668.331343697946 37.99937023 53.91754793 -gcp 7735.5 12649.765269294836 38.99678741 52.92901851 vlasenko/north/n-37.jpg n-37.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1277.1231511254018 582.6946141032761\n14550.654019292606 762.4730705163797\n14985.65498392283 16548.38228761799\n379.3634965572496 16354.88777880364" | \
  gdaltransform -tps -output_xy n-37.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  n-37.vrt n-37.cut.vrt

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