import json

# Load cat data
db = json.load(open('cat-position-data.json'))

# Print scale factors of map images to get to 1 pixel per metre
# This is used to resize each layer in road.xcf which is 1PPM.
for cat in db['cats']:
    print('{0}: {1:.2f}%'.format(cat['name'], 100.0/cat['pixelsPerMeter']))

from PIL import Image

# Get the size of each map image in pixels and hence pixels in road.tiff
image_sizes = {}
image_scales = {}
for cat in db['cats']:
    im = Image.open(cat['name'].lower() + '.jpg')
    w, h = im.size
    scale = 1.0/cat['pixelsPerMeter']
    w = int(w*scale); h = int(h * scale)
    print('{0}: {1}x{2}'.format(cat['name'], w, h))
    image_sizes[cat['name'].lower()] = (w,h)
    image_scales[cat['name'].lower()] = scale

# This is a list of the X, Y positions of the top-left (i.e. minimum) pixel co-ord of each image in road.xcf
image_locs = {
  'chip': (2159, 2034),
  'coco': (1636, 2007),
  'deebee': (1807, 2038),
  'ginger': (1294, 2174),
  'hermie': (1619, 1880),
  'kato': (1948, 2348),
  'phoebe': (1948, 2348),
  'orlando': (0, 0), # not geo-referenced
  'rosie': (1870, 2165),
  'sooty': (1481, 2051),
}

from osgeo import gdal

# Open road map and get the geo transform
map_ds = gdal.Open('road.tiff')
ox, xs, _, oy, _, ys = map_ds.GetGeoTransform()
def pixel_to_grid(x, y):
    return (ox + xs * x, oy + ys * y)

pixel_to_grid(0,0)

pixel_to_grid(2159, 2034)

from osgeo import osr

# Setup an OS grid -> WGS 84 projection
bng = osr.SpatialReference()
bng.ImportFromEPSG(27700)
wgs84 = osr.SpatialReference()
wgs84.ImportFromEPSG(4326)

bng2wgs84 = osr.CoordinateTransformation(bng, wgs84)

# Project the points for each cat into os grid and thence WGS84
for cat in db['cats']:
    name = cat['name'].lower()
    dx, dy = image_locs[name]
    
    # sckip unlocated images
    if dx == 0 and dy == 0:
        continue
    
    scale = image_scales[name]
    
    coords = []
    for pt in cat['points']:
        if not 'projected' in pt:
            continue
        x, y = tuple(float(v) for v in pt['projected'].split(','))
        x *= scale; y *= scale
        gx, gy = pixel_to_grid(x+dx, y+dy)
        
        # To WGS84
        x, y, _ = bng2wgs84.TransformPoint(gx, gy)
        coords.append((x, y))

    geometry = { 'type': 'LineString', 'coordinates': coords }
    feature = { 'type': 'Feature', 'properties': cat, 'geometry': geometry }
    json.dump(feature, open(name + '.geojson', 'w'))