#!/usr/bin/python
# -*- coding:Utf-8 -*-

#----------------------------------------------------------------------------------------------------
#	Visualization Project : 2D Spatial Data
#
# Author : Arthur COSTE
# Creation Date : October 14th, 2012
# Update Date :
# Version : 0
# Status : ON DEV
#
# Purpose : this code will process a set of input points to convert them into a VTK handleable points 
#
#----------------------------------------------------------------------------------------------------

# loading libraries
from array import array
import argparse
import vtk
import time

#parsing command line arguments
parser = argparse.ArgumentParser(epilog="Read input coordinates file and convert it to VTK points")
parser.add_argument("-v","--verbose", help="output verbosity", action="store_true")
parser.add_argument("-i", type=str,help="Input set of points")
parser.add_argument("-z", type=str,help="Input set of height")
parser.add_argument("-o", type=str,help="Output modified header path")
parser.add_argument("-out", type=str,help="Output PNM image")
parser.add_argument("-vtk", type=str,help="Input VTK file")
parser.add_argument("-voro", type=int,help="voronoi")
args = parser.parse_args()

print "-----------------------------------------------------------------------------------------------------------"

if args.i: 
	source_file = args.i
	if args.verbose: print "input file ok"
else: print "ERROR : Input file not specified"
if args.z: 
	source_file_2 = args.z
	if args.verbose: print "input height file ok"
else: print "ERROR : Input height file not specified"

if args.o:
	new_file = args.o
	if args.verbose: print "output file ok"
else: print "ERROR : Output file not specified"

print "-----------------------------------------------------------------------------------------------------------"
#process data
if args.i and args.z and args.o and not args.vtk: 

	points = vtk.vtkPoints()
	vertices = vtk.vtkCellArray()
	height = []
	print "processing"
	source_file = open(args.i, 'r')
	source_file_2 = open(args.z, 'r')
	new_file = open(args.o,'a')
	new_file.write("# vtk DataFile Version 2.0")
	new_file.write("\n")
	new_file.write("ASCII")
	new_file.write("\n")
	new_file.write("DATASET POLYDATA")
	new_file.write("\n")
	new_file.write("POINTS 618 float")
	new_file.write("\n")
	new_file.write("LOOKUP_TABLE default")
	new_file.write("\n")

	for i, line in enumerate (source_file_2):
		if args.verbose:print "point: ", i
		if args.verbose:print line
		height.append(float(line))

	for i, line in enumerate (source_file):
		if args.verbose:print "point: ", i
		if args.verbose:print line

		#p_x = str(line[2:19])
		p_x = str(line[0:9])
		if args.verbose:print p_x
		#p_y = str(line.split('      ')[1])
		p_y = str(line.split('    ')[1])
		if args.verbose:print p_y
		p=[round(float(p_x)),round(float(p_y)),round(height[i])]
		if args.verbose:print p

		new_file.write(str(p_x))
		new_file.write("\t")
		new_file.write(str(p_y))
		new_file.write("\t")
		new_file.write(str(height[i]))
		new_file.write("\n")
		id = points.InsertNextPoint(p)
		vertices.InsertNextCell(1)
		vertices.InsertCellPoint(id)

		# Create a polydata object
		point = vtk.vtkPolyData()
	 
		# Set the points and vertices we created as the geometry and topology of the polydata
		point.SetPoints(points)
		point.SetVerts(vertices)

	profile = vtk.vtkPolyData()
	profile.SetPoints(points)

	delny = vtk.vtkDelaunay2D()
	delny.SetInput(profile)
	delny.SetAlpha(30)
	delny.BoundingTriangulationOff()

	new_file.close()
	source_file.close()
	source_file_2.close()
	 
else: print "ERROR : MISSING ARGUMENTS"

# Define Look Up Table
lutLand = vtk.vtkLookupTable()
lutLand.SetNumberOfColors(201)
lutLand.Build()

# Visualize
mapper = vtk.vtkPolyDataMapper()

# Shrink the result to help see it better.
shrink = vtk.vtkShrinkFilter()
shrink.SetInputConnection(delny.GetOutputPort())
shrink.SetShrinkFactor(1)

map1 = vtk.vtkDataSetMapper()
map1.SetInputConnection(shrink.GetOutputPort())
#map1.SetScalarRange(-10, 10)
lut = vtk.vtkLookupTable()
lut.SetRange(0, 2000)
lut.SetNumberOfColors(15)
lut.SetValueRange(0, 1.0)
lut.SetRampToLinear()
lut.Build()
map1.SetLookupTable(lut)
#map1.SetScalarModeToUsePointData()


triangulation = vtk.vtkActor()
triangulation.SetMapper(map1)
triangulation.GetProperty().SetColor(0.2, 0.7, 0.9)

writer = vtk.vtkUnstructuredGridWriter()
writer.SetInput(delny.GetOutput());
writer.SetFileName(args.out)
writer.Update()

mapperT = vtk.vtkDataSetMapper()
mapperT.SetInputConnection(delny.GetOutputPort())
mapperT.SetScalarRange(-100, 100)
mapperT.SetScalarModeToUseCellData()
mapperT.SetLookupTable(lutLand)
mapperT.GetUseLookupTableScalarRange()
actorT = vtk.vtkActor()
actorT.SetMapper(mapperT)
#ren.AddActor(actorT)
 
if vtk.VTK_MAJOR_VERSION <= 5 and not args.vtk:
    mapper.SetInput(point)
else:
    mapper.SetInputData(point)
 
mapper.SetScalarRange(0, 250)
mapper.SetLookupTable(lutLand)
mapper.SetScalarModeToUseCellData()

actor = vtk.vtkActor()
actor.SetMapper(mapper)
actor.GetProperty().SetPointSize(1)


renderer = vtk.vtkRenderer()
renderWindow = vtk.vtkRenderWindow()
renderWindow.AddRenderer(renderer)
renderWindowInteractor = vtk.vtkRenderWindowInteractor()
renderWindowInteractor.SetRenderWindow(renderWindow)

if (args.voro == 1):
	renderer.AddActor(actorT)
else:
	renderer.AddActor(actor)
	#renderer.AddActor(triangulation)

renderWindow.Render()
renderWindowInteractor.Start()