#!/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 0.0 0.0 25.115175247192 57.419769226047 -gcp 12687.0 0.0 25.226411819458 57.428873062789 -gcp 12687.0 9627.0 25.241775512695 57.380647818129 -gcp 0.0 9627.0 25.130882263184 57.370745238986 LVVA_F7404_US1_GV1988_0002.jpg LVVA_F7404_US1_GV1988_0002.vrt

# Add cutline to VRT
gdalwarp -r $GDAL_SAMPLING_WARP -tps -dstalpha -cutline_srs EPSG:4326 \
  -cutline "POLYGON(($(echo -e \"2938.5 8110.5\n2770.5 8037.0\n2664.0 7833.0\n2488.5 7720.5\n2416.5 7695.0\n2314.5 7584.0\n2163.0 7546.5\n1779.0 7294.5\n1699.5 6996.0\n1509.0 6915.0\n1423.5 6768.0\n1288.5 6624.0\n1303.5 6522.0\n1281.0 6445.5\n1315.5 6283.5\n1252.5 6192.0\n1237.5 6075.0\n1218.0 5986.5\n1276.5 5848.5\n1258.5 5743.5\n1362.0 5515.5\n1353.0 5383.5\n1401.0 5313.0\n1555.5 5323.5\n1647.0 5367.0\n2815.5 5616.0\n2989.5 5314.5\n3163.5 5077.5\n3766.5 4164.0\n3949.5 3901.5\n4080.0 3739.5\n4312.5 3573.0\n4423.5 3397.5\n4536.0 3309.0\n5355.0 3523.5\n5571.0 3525.0\n5707.5 3432.0\n5907.0 3163.5\n6175.5 2974.5\n6369.0 2713.5\n7800.0 2034.0\n8046.0 1903.5\n8424.0 1836.0\n8569.5 2227.5\n8656.5 2308.5\n8767.5 2511.0\n8884.5 2625.0\n9015.0 2674.5\n9187.5 2713.5\n9366.0 2700.0\n9813.0 2580.0\n10036.5 2557.5\n10200.0 2590.5\n10362.0 2640.0\n10539.0 2610.0\n10632.0 2649.0\n10684.5 2778.0\n10731.0 2862.0\n10746.0 2937.0\n10720.5 3003.0\n10614.0 3046.5\n10444.5 3036.0\n10302.0 3090.0\n10317.0 3511.5\n11079.0 4006.5\n11775.0 4152.0\n11941.5 4113.0\n12034.5 4038.0\n12184.5 4090.5\n12306.0 4116.0\n12355.5 4320.0\n12168.0 4528.5\n11965.5 4696.5\n11866.5 4887.0\n11659.5 5296.5\n11491.5 5406.0\n11320.5 5464.5\n11190.0 5407.5\n11107.5 5293.5\n11002.5 5284.5\n10831.5 5359.5\n10683.0 5458.5\n10231.5 5791.5\n9886.5 5994.0\n9216.0 5550.0\n9129.0 5775.0\n9057.0 5850.0\n8952.0 6060.0\n8910.0 6049.5\n8901.0 5952.0\n8901.0 5806.5\n8457.0 6487.5\n8329.5 6592.5\n8107.5 6918.0\n8337.0 7107.0\n8005.5 7759.5\n6880.5 7498.5\n6915.75 7409.25\n6977.25 7353.0\n7037.25 7262.25\n7019.25 7146.0\n6997.5 7039.5\n7071.0 6897.75\n7002.75 6888.0\n6928.5 6858.0\n6885.0 6801.0\n6886.5 6754.5\n6879.0 6704.625\n6856.125 6675.375\n6833.25 6708.375\n6822.0 6783.0\n6793.5 6859.5\n6661.5 6982.5\n6516.0 7095.0\n6499.5 7258.5\n6417.0 7306.5\n6189.0 7233.0\n6025.5 7102.5\n5914.5 7209.0\n5832.0 7239.0\n5779.5 7341.0\n5482.5 7416.0\n5082.0 7806.0\n4845.0 7951.5\n4780.5 7996.5\n4692.0 7993.5\n4633.5 8130.0\n4450.5 8125.5\n4297.5 8076.0\n4224.0 7876.5\n4153.5 7869.0\n4089.0 7830.0\n4011.0 7467.0\n3891.0 7459.5\n3736.5 7473.0\n3745.5 7540.5\n3807.0 7671.0\n3789.0 7800.0\n3805.5 7914.0\n3763.5 7992.0\n3681.0 8058.0\n3489.0 8083.5\n3330.0 8040.0\n3303.0 8112.0\n3244.5 8148.0\" | \
  gdaltransform -tps -output_xy LVVA_F7404_US1_GV1988_0002.vrt | \
  awk 'NR==1{first=$0} {printf "%s %s,", $1,$2} END{print " " first}'))))" \
  LVVA_F7404_US1_GV1988_0002.vrt LVVA_F7404_US1_GV1988_0002.cut.vrt

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