#!/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 1076.25 750.2277901166019 29.99810041 55.9998536 -gcp 14328.0 730.676568573015 35.99817902 55.99994419 -gcp 15006.0 16485.907412548582 35.99834568 51.99987799 -gcp 427.5 16506.946973903385 29.99827426 51.99979183 -gcp 3280.5 819.0982787340367 30.99704013 55.98069036 -gcp 5485.5 867.2798445308161 31.99527314 55.96703334 -gcp 7702.5 880.787340366463 32.99828788 55.96258957 -gcp 9909.0 861.3064963908939 33.99655422 55.96691688 -gcp 12118.5 814.7815102720712 34.99672992 55.97843223 -gcp 14495.25 4662.980288728484 35.99735347 54.96599565 -gcp 14664.75 8607.852443087173 35.99730438 53.95267993 -gcp 14829.75 12537.490352581899 35.99563829 52.96690283 -gcp 12570.0 16578.17198778456 34.99472723 51.97839674 -gcp 10138.5 16629.378123264854 33.9942291 51.96684997 -gcp 7716.0 16650.365630205444 32.99798416 51.96260723 -gcp 5280.0 16630.84730705164 31.99616151 51.96816089 -gcp 2854.5 16587.41143808995 30.99811534 51.97949531 -gcp 591.0 12561.894086618546 29.99868475 52.96588855 -gcp 756.75 8632.96182676291 30.00037137 53.95136956 -gcp 916.5 4685.349111604664 29.99915 54.96527085 -gcp 7711.5 8769.327734591892 32.99909662 53.91424079 -gcp 10027.5 8751.960716268739 33.99690934 53.91780305 -gcp 7706.625 4825.544350360911 32.99853644 54.92630922 -gcp 5387.25 8757.143670183232 31.99772316 53.91817611 -gcp 7710.375 12704.029081066074 32.99710951 52.92815898 vlasenko/north/n-36.jpg n-36.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1076.25 750.2277901166019\n14328.0 730.676568573015\n15006.0 16485.907412548582\n427.5 16506.946973903385" | \
  gdaltransform -tps -output_xy n-36.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  n-36.vrt n-36.cut.vrt

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