#!/usr/bin/env python
"""
Converts a 16 bit flood tif to 8 bit.
"""

from osgeo import gdal
gdal.TermProgress = gdal.TermProgress_nocb

import numpy as np

import sys

_, input_tif, output_tif = sys.argv

indataset = gdal.Open(input_tif, gdal.GA_ReadOnly)

out_driver = gdal.GetDriverByName("GTiff")
outdataset = out_driver.Create(output_tif, indataset.RasterXSize, indataset.RasterYSize, indataset.RasterCount, gdal.GDT_Byte)

inband = indataset.GetRasterBand(1)
outband = outdataset.GetRasterBand(1)

gt = indataset.GetGeoTransform()
if gt is not None and gt != (0.0, 1.0, 0.0, 0.0, 0.0, 1.0):
    outdataset.SetGeoTransform(gt)

prj = indataset.GetProjectionRef()
if prj is not None and len(prj) > 0:
    outdataset.SetProjection(prj)

# we should only have one band
assert(indataset.RasterCount == 1)

inband = indataset.GetRasterBand(1)
outband = outdataset.GetRasterBand(1)

arr = inband.ReadAsArray()

# move 150 and 199 down to make room for the other values
arr[arr == 150] = 31 # cloud shadow
arr[arr == 199] = 18 # normal water

# now shift all values from 200-300 down by 100
arr[arr >= 200] = arr[arr >= 200] - 100

outband.WriteArray(arr)

# explicitly close to be safe; seems to be a gdal issue?
del outdataset, outband