#!/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 987.0 789.1020942408377 23.99781135 59.99981764 -gcp 12856.5 784.8726003490402 29.99787597 59.99991567 -gcp 13575.0 16537.92844677138 29.99810041 55.9998536 -gcp 298.5 16566.385689354276 23.99804261 55.9997597 -gcp 6924.375 910.174555802332 26.9978407 59.99986702 -gcp 6933.375 16685.712867851194 26.99806887 55.999807 -gcp 2966.25 854.3062187673515 24.99816116 59.98217464 -gcp 4950.0 898.5807884508606 26.00009528 59.97013693 -gcp 8901.75 900.827456968351 27.99602321 59.96917686 -gcp 10880.0625 846.5931079955581 28.99670955 59.98353606 -gcp 13035.0 4716.094947251527 29.99754162 58.95948072 -gcp 13212.75 8652.0622570794 29.99676417 57.94622073 -gcp 13397.25 12595.223278734036 29.99901466 56.95918937 -gcp 11361.75 16612.981815657968 28.99563069 55.98298292 -gcp 9149.4375 16671.24049139367 27.99515533 55.97015112 -gcp 4722.75 16680.483273181566 25.9962631 55.9701935 -gcp 2502.75 16638.857648528596 24.99396318 55.9812883 -gcp 473.25 12622.569336479733 23.99921664 56.95775744 -gcp 644.25 8679.865769017213 23.99869245 57.94305551 -gcp 817.125 4736.309550249862 23.99903909 58.95666955 -gcp 6927.0 8801.725152692949 26.99468103 57.91021807 -gcp 4835.4375 8789.862264019988 25.99762732 57.91387793 -gcp 6923.25 4857.393114936146 26.99465277 58.92387157 -gcp 9022.5 8784.106607440312 27.99410371 57.91402817 -gcp 6930.0 12742.178095502499 26.99436158 56.92562508 vlasenko/north/o-35.jpg o-35.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "987.0 789.1020942408377\n12856.5 784.8726003490402\n13575.0 16537.92844677138\n298.5 16566.385689354276" | \
  gdaltransform -tps -output_xy o-35.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  o-35.vrt o-35.cut.vrt

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