#!/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 1000.5 751.8156579677957 -0.00197134 55.99937174 -gcp 14241.0 727.8598001110494 5.99799925 55.99946845 -gcp 14931.0 16452.329816768462 5.99818235 51.99942549 -gcp 349.5 16489.871460299833 -0.00179093 51.99933351 -gcp 3212.25 822.1016102165464 0.99802225 55.99938775 -gcp 5417.25 868.6903109383675 1.99801644 55.99940381 -gcp 7619.782544564338 882.8861752341735 2.99801123 55.99941992 -gcp 9826.5 857.4340644086619 3.99800663 55.99943606 -gcp 12034.5 805.527484730705 4.99800263 55.99945224 -gcp 14415.0 4654.023598001111 5.99804931 54.9994575 -gcp 14585.25 8585.461549139367 5.99809636 53.99944669 -gcp 14757.0 12521.754927817878 5.99814064 52.99943602 -gcp 12498.75 16543.925180455302 4.99818543 51.99941007 -gcp 10068.0 16592.503609106054 3.99818906 51.99939468 -gcp 7638.0 16624.77234869517 2.99819324 51.99937933 -gcp 5211.75 16609.754025541366 1.99819797 51.99936401 -gcp 2776.5 16563.990699611328 0.99820325 51.99934874 -gcp 512.25 12551.297265408106 -0.00183203 52.99934281 -gcp 670.875 8616.528872848417 -0.00187565 53.99935228 -gcp 833.625 4678.625763464742 -0.00192201 54.99936193 -gcp 7627.125 8749.744968073292 2.99810776 53.99939931 -gcp 7622.25 4809.810174902832 2.998061 54.99940954 -gcp 7634.625 12687.487437534703 2.99815177 52.99938924 -gcp 9944.625 8724.092032204331 3.99810339 53.99941507 -gcp 5313.375 8732.345780122154 1.99811272 53.99938359 vlasenko/north/n-31.jpg n-31.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1000.5 751.8156579677957\n14241.0 727.8598001110494\n14931.0 16452.329816768462\n349.5 16489.871460299833" | \
  gdaltransform -tps -output_xy n-31.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  n-31.vrt n-31.cut.vrt

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