#!/usr/bin/env python
#
# kismet2shp
#
# 8/04/2006
# Dylan Beaudette
# Simple conversion script: Kismet XML report --> shapefile (network center or network range)
#
# projection stuff from http://www.larsen-b.com/Article/204.html
# shapefile field addition
# adapted from Sean Guiles
# http://zcologia.com/news/16
# notes from perrygeo : http://www.perrygeo.net/wordpress/?p=4
#
# command line parsing ideas from: http://docs.python.org/lib/module-optparse.html

# default output to latlong WGS84

#	range polygons
# ./kismet2shp -r -s --file=Kismet-Sep-03-2006-1.xml --outfile=range.shp

#	network centers
# ./kismet2shp -s --file=Kismet-Sep-03-2006-1.xml --outfile=points.shp  


#TODO : 
# compare output with gpsmap
# encode all information from Kismet
# figure out how to encode power data from Kismet-xxx.gps

'''usage: %prog kismet-xxx.xml [options]
   -r, --range: create range polygons instead of network center points
   -s --keep-srs leave points / polys in cartesian coordinates '''


import ogr, osr
from optparse import OptionParser
from elementtree import ElementTree
import sys, math

#process command line (optparse)
parser = OptionParser()
parser.add_option("-f", "--file", dest="infile", help="input Kismet-xxx.xml file", metavar="FILE")

parser.add_option("-o", "--outfile", dest="outfile", help="full name of shapefile to create", metavar="FILE")

parser.add_option("-r", "--range", action="store_true", dest="range_polys", default=False, help="compute range polygons instead of network centers [false]")

parser.add_option("-s", "--keep-srs", action="store_true", dest="keep_srs", default=False,  help="leave points / polygons in cartesian coordinates [false]")

(options, args) = parser.parse_args()


#hack
#require an input file
if not options.infile:
	print "ERROR: must supply an input file"
	sys.exit(1)

#hack
#require an output file
if not options.outfile:
	print "ERROR: must supply an output file"
	sys.exit(1)

#open input file
try:
   infile = options.infile
   #input data source
   data = open(infile,'r').read()
except:
   print "ERROR: Cannot open: " + infile
   sys.exit(1)
   


#setup output type
if options.range_polys == True:
	the_geom_type = ogr.wkbPolygon
else:
	the_geom_type = ogr.wkbPoint	


# Set up spatial reference systems
#something planimetric, units in meters: using our AEA projection
proj = osr.SpatialReference()
proj.ImportFromProj4('+proj=aea +x_0=0.0 +y_0=0 +lon_0=-96 +lat_0=40.0 +lat_1=20 +lat_2=60.0 +datum=NAD83')

#original lat/lon
latlong = osr.SpatialReference()
latlong.ImportFromProj4('+proj=latlong +datum=WGS84')

#create range output layer
filename = options.outfile
driver = ogr.GetDriverByName('ESRI Shapefile')
driver.DeleteDataSource(filename)
ds = driver.CreateDataSource(filename)


#output layer type
layer = ds.CreateLayer('wifi', geom_type=the_geom_type)

# Create a field for the SSID
fd = ogr.FieldDefn('ssid', ogr.OFTString)
fd.SetWidth(25)
fd.SetPrecision(0)
layer.CreateField(fd)

# Create a field for the approx range
fd = ogr.FieldDefn('range', ogr.OFTReal)
fd.SetWidth(10)
fd.SetPrecision(0)
layer.CreateField(fd)

# Create a field for the encryption range
fd = ogr.FieldDefn('encryption', ogr.OFTString)
fd.SetWidth(10)
fd.SetPrecision(0)
layer.CreateField(fd)

# Create a field for the packet info
fd = ogr.FieldDefn('num_packets', ogr.OFTInteger)
fd.SetWidth(10)
fd.SetPrecision(0)
layer.CreateField(fd)

#setup XML tree
detection = ElementTree.XML(data)

for node in detection:
	#need this in case there is a node with a missing SSID : i.e. an individual computer
	try:
		#extract SSID
		ssid = node.find('SSID').text
		encryption = node.find('encryption').text
		
		#extract GPS data
		gps = node.find('gps-info')
		
		#extract packet data
		packets = node.find('packets')
		num_packets = int(packets.find('LLC').text)
		
		lon_min = float(gps.find('min-lon').text)
		lon_max = float(gps.find('max-lon').text)
		
		#hack to prevent incorporation of (180,90) points ? GPS error?
		if (lon_min == 180):
			continue
		else:	
			lat_min = float(gps.find('min-lat').text)
			lat_max = float(gps.find('max-lat').text)
			
			#locate the center in latlong
			lon_center = (abs(lon_max - lon_min)/2) + lon_min
			lat_center = (abs(lat_max - lat_min)/2) + lat_min
					
			#init the feature
			f = ogr.Feature(feature_def=layer.GetLayerDefn())
						
			#work on the network center
			#make a text version of the point 
			wkt = 'POINT(%f %f)' % (lon_center, lat_center)
			p = ogr.CreateGeometryFromWkt(wkt)
			#assign latlon WG84 srs
			p.AssignSpatialReference(latlong)
			
			#now work on the min,max values
			wkt = 'POINT(%f %f)' % (lon_min, lat_max)
			p_1 = ogr.CreateGeometryFromWkt(wkt)
			wkt = 'POINT(%f %f)' % (lon_min, lat_min)
			p_2 = ogr.CreateGeometryFromWkt(wkt)
			#assign latlon WG84 srs
			p_1.AssignSpatialReference(latlong)
			p_2.AssignSpatialReference(latlong)
			
			#transform to planimetric projection:
			p.TransformTo(proj)
			p_1.TransformTo(proj)
			p_2.TransformTo(proj)
			
			#calculate range radius:
			#simple approximation via calculation of hypotenuse connecting center-side bounded right triangle- done in cartesian coordinate space
#			
#			working with Western Longitudes (negative values) and Northern Latitudes (positive)
#			
#	p_1(lon_min,lat_max)
#			o--------+--------+
#			|        |      / |
#			|        |  r /   | dy
#			|      c |  /     |
#			+--------o--------+
#			|        |   dx   |
#			|        |        |
#			|        |        |
#			+--------+--------o 
#						p_2(lon_max,lat_min)
			
			#note method for accessing x,y coords
			dx = abs(p_1.GetX() - p_2.GetX())/2
			dy = abs(p_1.GetY() - p_2.GetY())/2
			
			#cleanup
			p_1.Destroy()
			p_2.Destroy()
			
			# a^2 + b^2 = r^2
			approx_range = (dx**2 + dy**2) ** .5
			
			#setup output type
			if options.range_polys == True:
				#create buffer in XY space:
				poly_range = p.Buffer(approx_range)
				poly_range.AssignSpatialReference(proj)
				
				if options.keep_srs == False:
					#transform back to latlong
					poly_range.TransformTo(latlong)
				
				#register feature: the range buffer
				f.SetGeometryDirectly(poly_range)
				
			else:
				if options.keep_srs == False:
					#back transform to latlong
					p.TransformTo(latlong)
				#register feature: the center point
				f.SetGeometryDirectly(p)	
				
			
			#assign attribute data:
			f.SetField(0,ssid)
			f.SetField(1,approx_range)
			f.SetField(2,encryption)
			f.SetField(3,num_packets)
			
			#save feature
			layer.CreateFeature(f)
			f.Destroy()
		
	#catch a missing SSID
	except AttributeError:pass  #hidden SSID
	#print "BSSID: " +  node.find('BSSID').text,
	
# destroying data source closes the output file
ds.Destroy()

#write a prj file:
if options.keep_srs == False:
	wkt_proj = latlong.ExportToWkt()
else:
	wkt_proj = proj.ExportToWkt()

#open a filehandle
prj_file_name = filename[:-4] + '.prj'
prj_file = open(prj_file_name, 'w')
prj_file.write(wkt_proj)
prj_file.close()