#!/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 759.1967845659163 209.45072182121044 53.99835636 60.00026304 -gcp 8659.76961414791 242.455719044975 59.99852535 60.00033392 -gcp 9080.553214383215 10740.849863317104 59.99868117 56.00025424 -gcp 248.45498392282957 10714.684064408662 53.99853004 56.00018634 -gcp 4705.5 318.77155052748475 56.99843871 60.00029944 -gcp 4663.6875 10834.54632842865 56.99860369 56.00022121 -gcp 2076.0 265.8770739064857 54.99895958 59.9798015 -gcp 3392.25 297.43815987933635 55.99794672 59.96933496 -gcp 6027.0 315.0 57.99831542 59.96789302 -gcp 7341.0 288.0 58.99726427 59.98067695 -gcp 8763.75 2865.476809954751 59.99771661 58.95909338 -gcp 8869.5 5493.971342383107 59.99806604 57.94392078 -gcp 8974.875 8117.552601809955 59.99827205 56.95841479 -gcp 7612.5 10794.103695324284 58.99950513 55.97950589 -gcp 6138.852281042883 10825.497248978667 57.99818219 55.96663343 -gcp 3194.25 10812.84539969834 55.99956739 55.96803475 -gcp 1722.75 10778.725301659124 55.0006183 55.97873015 -gcp 370.5 8091.505279034691 53.99466528 56.95758506 -gcp 496.125 5464.920060331825 53.99303102 57.94353148 -gcp 626.25 2834.1059577677224 53.99446723 58.95896291 -gcp 4694.25 2945.718137254902 56.99442507 58.92149548 -gcp 6082.125 5572.252262443439 57.99710719 57.91031275 -gcp 4682.625 5578.978318250377 56.99425109 57.90584772 -gcp 3289.875 5557.966440422322 55.99642339 57.91194278 -gcp 4672.6875 8209.230203619909 56.99525364 56.91925389 vlasenko/north/o-40.jpg o-40.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "759.1967845659163 209.45072182121044\n8659.76961414791 242.455719044975\n9080.553214383215 10740.849863317104\n248.45498392282957 10714.684064408662" | \
  gdaltransform -tps -output_xy o-40.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  o-40.vrt o-40.cut.vrt

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