commit 6f3cb48c90143e20d022fec24bfe759bda3eadbb Author: Andreas Demmelbauer Date: Wed Sep 15 15:25:00 2021 +0200 initial commit diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..7b1c5d6 --- /dev/null +++ b/.gitignore @@ -0,0 +1,10 @@ +source_images* +sources +output +playing-around +__pycache__ +*.png +*.R +*.tif +*.tiff +*.jpg \ No newline at end of file diff --git a/README.md b/README.md new file mode 100644 index 0000000..e69de29 diff --git a/Requirements.txt b/Requirements.txt new file mode 100644 index 0000000..0dc4673 --- /dev/null +++ b/Requirements.txt @@ -0,0 +1,3 @@ +lxml +wand +PyExifTool \ No newline at end of file diff --git a/canvas.svg b/canvas.svg new file mode 100644 index 0000000..878d9ea --- /dev/null +++ b/canvas.svg @@ -0,0 +1,110 @@ + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + diff --git a/export-gis-jpg.py b/export-gis-jpg.py new file mode 100644 index 0000000..c14a5c3 --- /dev/null +++ b/export-gis-jpg.py @@ -0,0 +1,263 @@ +import os +import argparse +import lxml.etree as ET +import subprocess +import flirimageextractor +import cv2 +import numpy as np +from pathlib import Path +from wand.image import Image +from osgeo import gdal +from osgeo import osr + +arg_parser = argparse.ArgumentParser(description='Export SVG composition of FLIR images as TIFF with thermo layer') + +arg_parser.add_argument('Input', + metavar='input_svg', + type=str, + help='Path to the input SVG file cotaining xlinks to FLIR images') + +arg_parser.add_argument('Output', + metavar='output_tiff', + type=str, + help='Output filename') + + +args = arg_parser.parse_args() +dirname = os.path.dirname(__file__) +INPUT_PATH = os.path.join(dirname, args.Input) +INPUT_DIR = os.path.split(INPUT_PATH)[0] + +TEMP_MAP_THERMALPNG_SVG_PATH = os.path.join(INPUT_DIR, 'map_thermalpng.svg') +TEMP_MAP_THERMALPNG_PATH = os.path.join(INPUT_DIR, 'map_thermalpng.png') +TEMP_MAP_PREVIEW_PATH = os.path.join(INPUT_DIR, 'map_preview.png') +THERMALPNG_DIR = 'thermalpngs' + +OUTPUT_PATH = os.path.join(dirname, args.Output) + + + +def make_thermalpng_tiles(): + """ + Extract thermal infomration as greyscale PNG-16 (temp * 1000 to retain some decimals) + and save the png tiles + """ + Path(os.path.join(INPUT_DIR, THERMALPNG_DIR)).mkdir(parents=True, exist_ok=True) + png_output_dir = os.path.join(INPUT_DIR, THERMALPNG_DIR) + + for root_path, directories, file in os.walk(os.path.join(dirname, INPUT_DIR)): + for file in file: + if(file.endswith(".jpg")): + print('Extracting thermal info from ' + file) + full_filepath = os.path.join(root_path, file) + flir = flirimageextractor.FlirImageExtractor() + flir.process_image(full_filepath) + thermal_img_np = flir.thermal_image_np + multiplied_image = cv2.multiply(thermal_img_np, 1000) + output_file_path = os.path.join(png_output_dir, file + '.thermal.png') + print(output_file_path) + cv2.imwrite(output_file_path, multiplied_image.astype(np.uint16)) + + + + +def make_thermalpng_svg(): + """ + replaces the image paths with the thermal pngs + and creates new SVG file + """ + # print("svg_file") + # print(dir(svg_file)) + tree = ET.parse(INPUT_PATH) + root = tree.getroot() + # print(ET.tostring(root)) + # tile_rows = root.xpath('//image', namespaces={'n': "http://www.w3.org/2000/svg"}) + # print(dir(root)) + tile_elements = root.xpath('//*[@class="thermal_image"]') + linkattrib ='{http://www.w3.org/1999/xlink}href' + for tile in tile_elements: + tile.attrib[linkattrib] = os.path.join(THERMALPNG_DIR, tile.attrib[linkattrib] + '.thermal.png') + # newxml = ET.tostring(tree, encoding="unicode") + # print(newxml) + # return newxml + + with open(TEMP_MAP_THERMALPNG_SVG_PATH, 'wb') as f: + tree.write(f, encoding='utf-8') + + return tree + + + +def make_thermalpng(): + """ + exports the SVG canvas as Gray_16 PNG + """ + command = [ + '/snap/bin/inkscape', + '--pipe', + '--export-type=png', + '--export-png-color-mode=Gray_16' + ], + input_file = open(TEMP_MAP_THERMALPNG_SVG_PATH, "rb") + output_file = open(TEMP_MAP_THERMALPNG_PATH, "wb") + completed = subprocess.run( + *command, + cwd=INPUT_DIR, # needed for reative image links + stdin=input_file, + stdout=output_file + ) + return completed + +def make_thermalpreview(): + """ + exports the preview image + """ + command = [ + '/snap/bin/inkscape', + '--pipe', + '--export-type=png', + '--export-png-color-mode=Gray_8' + ], + input_file = open(TEMP_MAP_THERMALPNG_SVG_PATH, "rb") + output_file = open(TEMP_MAP_PREVIEW_PATH, "wb") + completed = subprocess.run( + *command, + cwd=INPUT_DIR, # needed for reative image links + stdin=input_file, + stdout=output_file + ) + return completed + +# def make_thermalpreview(): +# """ +# exports the preview image +# """ +# command = [ +# '/snap/bin/inkscape', +# '--pipe', +# '--export-type=png', +# '--export-png-color-mode=Gray_8' +# ] +# input_file = open(TEMP_MAP_THERMALPNG_SVG_PATH, "rb") +# output_file = open(TEMP_MAP_PREVIEW_PATH, "wb") +# completed = subprocess.run( +# *command, +# cwd=INPUT_DIR, # needed for reative image links +# stdin=input_file, +# stdout=output_file +# ) +# return completed + + +def get_thermal_numpy_array(): + # input_file = open(TEMP_MAP_THERMALPNG_PATH, "rb") + image = cv2.imread(TEMP_MAP_THERMALPNG_PATH, cv2.IMREAD_ANYDEPTH) + image_float = image.astype(np.float32) + image_float_normalized = cv2.divide(image_float, 1000) + print(image_float_normalized[1000][905]) + # cv2.imshow("OpenCV Image Reading", image) + return image_float_normalized + +def get_used_tiles_relpaths(): + """ + outputs an array of all used tile filenames in the input SVG + (relative filepaths like they appear in the svg.) + """ + images = [] + tree = ET.parse(INPUT_PATH) + root = tree.getroot() + tile_elements = root.xpath('//*[@class="thermal_image"]') + linkattrib ='{http://www.w3.org/1999/xlink}href' + for tile in tile_elements: + images.append(tile.attrib[linkattrib]) + return images + +def deg_coordinates_to_decimal(coordStr): + coordArr = coordStr.split(', ') + calculatedCoordArray = [] + for calculation in coordArr: + calculationArr = calculation.split('/') + calculatedCoordArray.append(int(calculationArr[0]) / int(calculationArr[1])) + degrees = calculatedCoordArray[0] + minutes = calculatedCoordArray[1] + seconds = calculatedCoordArray[2] + decimal = (degrees + (minutes * 1/60) + (seconds * 1/60 * 1/60)) + # print(decimal) + return decimal + + +def read_coordinates_from_tile(filename): + full_filepath = os.path.join(INPUT_DIR, filename) + with Image(filename=full_filepath) as image: + for key, value in image.metadata.items(): + if key == 'exif:GPSLatitude': + # print('latstr', value) + lat = deg_coordinates_to_decimal(value) # lat -> Y vertical + if key == 'exif:GPSLongitude': + # print('lonstr', value) + lon = deg_coordinates_to_decimal(value) # lon -> X horizontal + if key == 'exif:GPSImgDirection': + direction = value.split('/') + print(int(direction[0])/int(direction[1])/2, ' ', (value)) + + return [lat, lon] + +def get_coordinate_boundaries(): + image_names = get_used_tiles_relpaths() + coordinates = { + 'lat': [], + 'lon': [] + } + for filename in image_names: + tile_coordinates = read_coordinates_from_tile(filename) + coordinates['lat'].append(tile_coordinates[0]) + coordinates['lon'].append(tile_coordinates[1]) + + boundaries = { + 'xmin': min(coordinates['lon']), + 'xmax': max(coordinates['lon']), + 'ymin': min(coordinates['lat']), + 'ymax': max(coordinates['lat']), + } + return boundaries + + +def make_geotiff_image(): + thermal_numpy_array = get_thermal_numpy_array() + + # coordinates of all tiles + geo_bound = get_coordinate_boundaries() + print('boundaries', geo_bound) + + np_shape = thermal_numpy_array.shape + image_size = (np_shape[0], np_shape[1]) + + # set geotransform + nx = image_size[0] + ny = image_size[1] + xres = (geo_bound['xmax'] - geo_bound['xmin']) / float(nx) + yres = (geo_bound['ymax'] - geo_bound['ymin']) / float(ny) + geotransform = (geo_bound['xmin'], xres, 0, geo_bound['ymax'], 0, -yres) + + # create the 3-band raster file + dst_ds = gdal.GetDriverByName('GTiff').Create(OUTPUT_PATH, ny, nx, 1, gdal.GDT_Float32) + dst_ds.SetGeoTransform(geotransform) # specify coords + srs = osr.SpatialReference() # establish encoding + res = srs.SetWellKnownGeogCS( "WGS84" ) # WGS84 lat/long + dst_ds.SetProjection(srs.ExportToWkt()) # export coords to file + dst_ds.GetRasterBand(1).WriteArray(thermal_numpy_array) # write thermal-band to the raster + dst_ds.FlushCache() # write to disk + + + + +# make_thermalpng_tiles() +# make_thermalpng_svg() +# make_thermalpreview() +# make_thermalpng() +make_geotiff_image() + + +# dataset = gdal.Open("working_result_example.tif", gdal.GA_ReadOnly) +# print(dir(dataset)) +# print(dataset.GetMetadata_List()) diff --git a/gis-svg-stitcher.py b/gis-svg-stitcher.py new file mode 100644 index 0000000..3853cc0 --- /dev/null +++ b/gis-svg-stitcher.py @@ -0,0 +1,238 @@ +from wand.image import Image +import PIL.Image +import io +import exiftool +import subprocess +import os +import lxml.etree as ET +import copy + +import cv2 +import flirimageextractor +from matplotlib import cm +import numpy as np +import urllib.request + + + +dirname = os.path.dirname(__file__) + +working_dir = 'source_images_full' + + +filename = os.path.join(dirname, 'canvas.svg') +tree = ET.parse(filename) + +root = tree.getroot() +d = root.nsmap + +main_layer = root.xpath('//*[@id="tiles"]', namespaces={'n': "http://www.w3.org/2000/svg"})[0] + +tile_rows = root.xpath('//*[@id="tile_rows"]', namespaces={'n': "http://www.w3.org/2000/svg"})[0] + + +def deg_coordinates_to_decimal(coordStr): + coordArr = value.split(', ') + calculatedCoordArray = [] + for calculation in coordArr: + calculationArr = calculation.split('/') + calculatedCoordArray.append(int(calculationArr[0]) / int(calculationArr[1])) + degrees = calculatedCoordArray[0] + minutes = calculatedCoordArray[1] + seconds = calculatedCoordArray[2] + return (degrees + (minutes * 1/60) + (seconds * 1/60 * 1/60)) + +# # extracting TIF Data +# for root, directories, file in os.walk(os.path.join(dirname, working_dir)): +# for file in file: +# if(file.endswith(".jpg")): +# print(os.path.join(root, file)) +# full_filepath = os.path.join(root, file) +# with exiftool.ExifTool() as et: +# cmd = ['exiftool', full_filepath, "-b", "-RawThermalImage"] +# tif_data = subprocess.check_output(cmd) +# tif_image = PIL.Image.open(io.BytesIO(tif_data)) +# tif_filepath = os.path.join(dirname, working_dir, file.split('.')[0] + '_thermal.tif') +# tif_image.save(tif_filepath) +# print(tif_filepath) + + + +# finding the boundaries of the whole canvas +latsArr = [] +lonsArr = [] +for root_path, directories, file in os.walk(os.path.join(dirname, working_dir)): + for file in file: + if(file.endswith(".jpg")): + print(os.path.join(root_path, file)) + full_filepath = os.path.join(root_path, file) + with Image(filename=full_filepath) as image: + print(image.width) + print(image.height) + for key, value in image.metadata.items(): + if key == 'exif:GPSLatitude': + lat = deg_coordinates_to_decimal(value) # lat -> Y vertical + latsArr.append(lat) + print("{}: {}".format(key, value)) + print('lat '+ str(lat)) + if key == 'exif:GPSLongitude': + lon = deg_coordinates_to_decimal(value) # lon -> X horizontal + lonsArr.append(lon) + print("{}: {}".format(key, value)) + print('lon '+ str(lon)) + + +minLat = min(latsArr) +minLon = min(lonsArr) +maxLat = max(latsArr) +maxLon = max(lonsArr) +width = maxLon - minLon +height = maxLat- minLat + + + + +# placing the images into the svg + +rotation = 125 + +y_scale = -1800000 #-400000 +x_scale = 655000 #-950000 + +# y_scale = 2600000 #-400000 +# x_scale = 1200000 #-950000 + +image_rotation_up = rotation #32 +image_rotation_down = rotation + 180 #192 + +for root_path, directories, file in os.walk(os.path.join(dirname, working_dir)): + for file in file: + if(file.endswith(".jpg")): + print(os.path.join(root_path, file)) + full_filepath = os.path.join(root_path, file) + with Image(filename=full_filepath) as image: + print(image.width) + print(image.height) + for key, value in image.metadata.items(): + # print("{}: {}".format(key, value)) + if key == 'exif:GPSLatitude': + lat = deg_coordinates_to_decimal(value) - minLat + print('lat '+ str(lat)) + if key == 'exif:GPSLongitude': + lon = deg_coordinates_to_decimal(value) - minLon + print('lon '+ str(lon)) + if key == 'exif:GPSImgDirection': + direction = value.split('/') + rotation = int(direction[0])/int(direction[1])/2 + + g_pos_el_attributes = { + # 'x': str(lat*scale), + # 'y': str(lon*scale), + 'transform': "translate({}, {})".format(format(lon*x_scale, '.20f'), format(lat*y_scale, '.20f')), + 'data-lat': format(lat, '.20f'), + 'data-lon': format(lon, '.20f'), + 'class': 'tile', + 'id': 'tile_{}'.format(file.split('.')[0]), + # 'style': 'opacity:.6', + } + g_pos_el = ET.SubElement(main_layer, 'g', attrib=g_pos_el_attributes) + + g_offset_corr_el_attributes = { + 'transform': "translate(150, 0)", + 'class': 'tile-offset-corr', + } + g_offset_corr_el = ET.SubElement(g_pos_el, 'g', attrib=g_offset_corr_el_attributes) + + + g_center_el_attributes = { + 'class': 'tile-center', + 'transform': 'translate({}, {})'.format(str(image.width/2*-1), str(image.height/2*-1)) + } + g_center_el = ET.SubElement(g_offset_corr_el, 'g', attrib=g_center_el_attributes) + + g_rot_el_attributes = { + 'class': 'tile-rotate', + 'data-image-rotation': str(image_rotation_up), + 'data-image-dimensions': str(image.width/2) + ' ' + str(image.height/2), + 'transform': 'rotate({} {} {})'.format(str(image_rotation_up), str(image.width/2), str(image.height/2)) + } + g_rot_el = ET.SubElement(g_center_el, 'g', attrib=g_rot_el_attributes) + + xlinkns ="http://www.w3.org/1999/xlink" + image_el = ET.SubElement(g_rot_el, 'image', { + "class": 'thermal_image', + "{%s}href" % xlinkns: file, + "width": str(image.width), + "height": str(image.height), + "mask" : "url(#tilefademask)", + }) + +# transform_str = "translate(-{}, -{})".format(str(min(latsArr)*scale), str(min(lonsArr)*scale)) +# print(transform_str) +# main_layer.attrib['transform'] = transform_str + + + +# sort elements +def getkey(elem): + # Used for sorting elements by @LIN. + # returns a tuple of ints from the exploded @LIN value + # '1.0' -> (1,0) + # '1.0.1' -> (1,0,1) + return float(elem.get('id').split('_')[2]) + +main_layer[:] = sorted(main_layer, key=getkey) + + +# rotate image if previous element is under the current one +last_state = 'down' +for index, el in enumerate(main_layer): + if(el.getprevious() is not None): + if (el.getprevious().attrib['data-lon'] > el.attrib['data-lon'] or (el.getprevious().attrib['data-lon'] == el.attrib['data-lon'] and last_state == 'up')): + print('up') + rot_el = el[0][0][0] + # print(rot_el.attrib['data-image-rotation']) + # print(rot_el.attrib['data-image-dimensions']) + el.attrib['data-direction'] = 'up' + + # print(el.attrib['data-lat'], el.getprevious().attrib['data-lat']) + else: + rot_el = el[0][0][0] + el.attrib['data-direction'] = 'down' + # el.attrib['style'] = 'opacity:0' + new_rotation = image_rotation_down #float(rot_el.attrib['data-image-rotation']) + 180 + rot_el.attrib['transform'] = "rotate({} {})".format(str(new_rotation), rot_el.attrib['data-image-dimensions']) + print('down') + # print(rot_el.attrib['data-image-rotation']) + # print(rot_el.attrib['data-image-dimensions']) + + # merge tiles into groups + print(index) + print("el.attrib['data-direction'] " + el.attrib['data-direction']) + print("last_state " + last_state) + if index is 1 or last_state != el.attrib['data-direction']: + current_row = ET.SubElement(tile_rows, 'g', attrib={ 'class': 'tile-row' }) + copyElem = copy.deepcopy(el) + current_row.insert(0, copyElem) + last_state = el.attrib['data-direction'] + +root.remove(main_layer) + + +with open(os.path.join(working_dir,'map.svg'), 'wb') as f: + tree.write(f, encoding='utf-8') + + + +# # get some base satellite map for reference +# apikey = "MYaMHCLtPz1fUfe0FzZqOMI35m893jIV80oeHG19Piw" +# lon_center = +# lat_center = +# zoom = +# map_width = +# request = "https://image.maps.ls.hereapi.com/mia/1.6/mapview?apiKey={}&c={},{}&sb=mk&t=1&z={}&w={}&nodot".format(apikey, lon_center, lat_center, zoom, map_width) + +# svg = ET.tostring(tree, encoding="unicode") +# print(svg) + +print('Done!') diff --git a/mask.svg b/mask.svg new file mode 100644 index 0000000..b3c8604 --- /dev/null +++ b/mask.svg @@ -0,0 +1,197 @@ + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + + + + + + + diff --git a/thermal-image-extractor.py b/thermal-image-extractor.py new file mode 100644 index 0000000..e1104e5 --- /dev/null +++ b/thermal-image-extractor.py @@ -0,0 +1,43 @@ +import cv2 +# import flirimageextractor +# flir = flirimageextractor.FlirImageExtractor() +# flir.process_image('source_images_full/20200621_125936_R.jpg',RGB=True) +# c = flir.get_metadata('source_images_full/20200621_125936_R.jpg') +# img = flir.extract_embedded_image() + +import flirimageextractor +from matplotlib import cm +import numpy as np +# flir = flirimageextractor.FlirImageExtractor() +# flir.process_image('source_images_full/20200621_125936_R.jpg') +# flir.save_images() +# flir.plot() + + + +flir = flirimageextractor.FlirImageExtractor() +flir.process_image('source_images_full/20200621_125936_R.jpg') +# img = flir.save_images() #min=0,max=64, +# print(dir(img)) +# print(img) +thermal_img_np = flir.thermal_image_np +h,w = thermal_img_np.shape +# print(flir.flir_img_bytes()) + +print(thermal_img_np) +print(type(thermal_img_np[0][0])) +print(h) +print(w) + +# cv2.imdecode(thermal_img_np) + +multiplied_image = cv2.multiply(thermal_img_np, 1000) + +cv2.imwrite('color_img.png', multiplied_image.astype(np.uint16)) + + + +# img.astype('uint8') + +# image = cv2.imdecode(thermal_img_np,4) +# print('Image Dimensions :', image.shape)