#!/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 1971.0 563.4857023875625 11.9940079 67.99943712 -gcp 19764.75 721.0248473070517 23.9970798 67.99993477 -gcp 21145.5 16448.413242642975 23.99750414 63.99987608 -gcp 338.25 16278.964880621876 11.99743929 63.99966641 -gcp 10873.125 1065.8334952803998 17.99702557 67.99982737 -gcp 10735.5 16833.45578845086 17.99745779 63.99977191 -gcp 4936.5 828.1432880844645 14.00452414 67.93413107 -gcp 7897.5 994.761689291101 15.99969136 67.89575969 -gcp 13836.0 1047.3016591251885 19.98930011 67.89590527 -gcp 16807.5 934.7948717948719 21.99170592 67.93399943 -gcp 20113.5 4651.920814479638 23.99956221 66.94101467 -gcp 20808.75 12511.167043740572 24.00275223 64.94327776 -gcp 17683.5 16672.892911010556 21.99104436 63.94052925 -gcp 14205.0 16815.876319758674 19.98620998 63.90038353 -gcp 7260.0 16756.710030165912 15.99640167 63.90095764 -gcp 3790.5 16573.917797888385 13.99782497 63.93711137 -gcp 744.0 12346.660633484164 11.99500848 64.94239376 -gcp 1151.25 8423.203996983408 11.99389109 65.92080317 -gcp 1559.25 4499.005656108598 11.99311732 66.93845693 -gcp 10806.0 8957.760180995476 17.99668979 65.80621573 -gcp 10842.0 5008.595776772247 17.99922013 66.82551643 -gcp 14022.0 8928.295625942685 19.98861232 65.82063517 -gcp 4365.0 8703.918174962293 14.00197093 65.85655082 -gcp 10773.0 12897.958898944195 17.99609353 64.82845493 -gcp 7581.0 8873.40874811463 15.99925294 65.82058611 vlasenko/north/q-33_34.jpg q-33_34.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e "1971.0 563.4857023875625\n19764.75 721.0248473070517\n21145.5 16448.413242642975\n338.25 16278.964880621876" | \
  gdaltransform -tps -output_xy q-33_34.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}')))" \
  q-33_34.vrt q-33_34.cut.vrt

# Generate tiles
gdal2tiles.py -r $GDAL_SAMPLING_TILE --xyz -z $GDAL_TILE_ZOOMS -x --processes=$GDAL_TILE_PROCESSES q-33_34.cut.vrt q-33_34.xyz