# -*- coding: utf-8 -*-
'''
GeophPy.geoposset
-----------------
Geographics Positioning Set management.
:copyright: Copyright 2014 Lionel Darras, Philippe Marty and contributors, see AUTHORS.
:license: GNU GPL v3.
'''
from __future__ import absolute_import
import shapefile # pyshp for reading and writing shapefiles
import geophpy.geopositioning.kml as kml
import matplotlib
import numpy as np
import matplotlib.pyplot as plt
import glob # for managing severals files thanks to "*." extension
import utm # to convert utm and wgs84 coordinates
import os
# definitions
POINT_PARAMS = 'bo'
MULTIPOINT_PARAMS = 'bo--'
LINE_PARAMS = 'r'
SURFACE_PARAMS = 'g'
refsystem_list = ['UTM', 'WGS84'] # list of reference systems availables
type_list = ['point', 'shapefile'] # list of element types availables
utm_minletter = 'E'
utm_minnumber = 1
utm_maxletter = 'X'
utm_maxnumber = 60
#---------------------------------------------------------------------------#
# Geographics Positioning Set Object #
#---------------------------------------------------------------------------#
[docs]class GeoPosSet(object):
'''
Builds a geographic positioning set to display or use to georeferencing data set.
'''
refsystem = None # Geographic positioning Reference System (None, 'UTM', 'WGS84', ...)
type = None # Type of elements contained in the array ('point', 'shapefile', ...)
list = [] # List of elements (points, shapefiles, ...)
def __init__(self):
pass
@classmethod
[docs] def from_file(cls, refsystem, type, filenameslist):
'''Builds a geographic positioning set from a list of files.
Parameters :
:refsystem: geographic positioning reference system, None, 'UTM', 'WGS84', ...
:type: type of elements contained in the files, 'point', 'shapefile', ...
:filenameslist: list of files to open,
['shapefile1', 'shapefile2', ...] or,
['shapefile*'] to open all files with a filename beginning by "shapefile" (without extensions).
Returns :
:success: true if geoposset built, fals if not.
:geoposset: geographics positioning set built with file(s).
Example :
success, geoposset = GeoPosSet.from_file('UTM', 'shapefile', "shapefile1")
'''
geoposset = cls()
geoposset.type = type
geoposset.refsystem = refsystem
geoposset.list = []
success = True # by default
if (type == 'shapefile') :
for filenames in filenameslist : # for each file in the list
try:
filenamewithoutextension, extension = os.path.splitext(filenames)
filenamesextended = glob.glob(filenamewithoutextension+".shp") # extension if the filenames if '*' is contained
for filename in filenamesextended :
geoposset.list.append(shapefile.Reader(filename)) # adds shapefile in the list
except IOError:
print("File Read Error (",filename,")")
success = False
break # exits the loop
else :
success = False
return success, geoposset
[docs] def points_getlist(self):
'''
get list of points
Returns: list of points
'''
point_num = 0 # initialisation of the point number
points_list = [] # initialisation of the points list
if (self.type == 'shapefile'):
firstpoint = True
for sf in self.list : # for each shape file in the list
shapes = sf.shapes()
for shape in shapes : # for each shape in the file
# if point type :
if ((shape.shapeType == shapefile.POINT) or (shape.shapeType == shapefile.POINTZ) or (shape.shapeType == shapefile.POINTM)) :
points_list.append([point_num, shape.points[0][0], shape.points[0][1]])
point_num = point_num + 1
else :
points_list = []
return points_list
[docs] def plot(self, picturefilename=None, dpi=None, transparent=False, i_xmin=None, i_xmax=None, i_ymin=None, i_ymax=None):
'''
display map of points, lines, and surfaces described in geographic positioning list
Parameters :
:picturefilename: name of the file to save the picture. If None, no picture saved in a file
:dpi:
:transparent: if True, picture display points not plotted as transparents
:i_xmin: x minimal value to display, None by default
:i_xmax: x maximal value to display, None by default
:i_ymin: y minimal value to display, None by default
:i_ymax: y maximal value to display, None by default
Returns :
:success: True if no error
:fig: figure object
'''
fig = plt.figure() # creation of the empty figure
point_num = 0 # initialisation of the point number
points_list = [["Num", "Lon", "Lat"]] # initialisation of the points list
firstpoint = True
success = True
if (self.type == 'shapefile'):
for sf in self.list : # for each shape file in the list
shapes = sf.shapes()
for shape in shapes : # for each shape in the file
# if point type :
if ((shape.shapeType == shapefile.POINT) or (shape.shapeType == shapefile.POINTZ) or (shape.shapeType == shapefile.POINTM)) :
plot_params = MULTIPOINT_PARAMS
points_list.append([point_num, shape.points[0][0], shape.points[0][1]])
plt.text(shape.points[0][0], shape.points[0][1],str(point_num))
point_num = point_num + 1
# if multipoint type :
elif ((shape.shapeType == shapefile.MULTIPOINT) or (shape.shapeType == shapefile.MULTIPOINTZ) or (shape.shapeType == shapefile.MULTIPOINTM)) :
plot_params = MULTIPOINT_PARAMS
# if line type :
elif ((shape.shapeType == shapefile.POLYLINE) or (shape.shapeType == shapefile.POLYLINEZ) or (shape.shapeType == shapefile.POLYLINEM)) :
plot_params = LINE_PARAMS
# if surface type :
else :
plot_params = SURFACE_PARAMS
points = np.array(shape.points)
if (firstpoint == True) : # if first point to display, initialises the x and y limits values
xmin = min(points.T[0])
xmax = max(points.T[0])
ymin = min(points.T[1])
ymax = max(points.T[1])
firstpoint = False
else : # else compares with the x and y limits established
if (min(points.T[0]) < xmin) :
xmin = min(points.T[0])
if (max(points.T[0]) > xmax) :
xmax = max(points.T[0])
if (min(points.T[1]) < ymin) :
ymin = min(points.T[1])
if (max(points.T[1]) > ymax) :
ymax = max(points.T[1])
# displays the point, the line, or the surface
plt.plot(points.T[0], points.T[1], plot_params)
else:
success = False
if (success == True):
# for each x or y limit not configured in input, initialisation of the value with the limits of points displayed
nolimits = True
if (i_xmin == None) :
i_xmin = xmin
else :
nolimits = False
if (i_xmax == None) :
i_xmax = xmax
else :
nolimits = False
if (i_ymin == None) :
i_ymin = ymin
else :
nolimits = False
if (i_ymax == None) :
i_ymax = ymax
else :
nolimits = False
if (nolimits == False) : # sets the axis limits
xmin, xmax, ymin, ymax = plt.axis([i_xmin, i_xmax, i_ymin, i_ymax])
else : # gets the axis limits
xmin, xmax, ymin, ymax = plt.axis()
# to have the same scale in X and Y axis
dy = ymax - ymin
dx = xmax - xmin
if (dy > dx):
xmax = xmin + dy
elif (dx > dy):
ymax = ymin + dx
plt.axis([xmin, xmax, ymin, ymax])
ax = plt.gca()
ax.ticklabel_format(useOffset=False)
if (picturefilename != None) : # if picturefilename in input, saves in a picture file
plt.savefig(picturefilename, dpi=dpi, transparent=transparent)
return success, fig
[docs] def to_kml(self, kmlfilename):
"""
Builds a kml file with points, lines, and surfaces described in a list of geographic positioning set
Parameters :
:kmlfilename: Name of the kml file to save
Returns :
:success: True if operation done
"""
success = True # by default
if (self.type == 'shapefile'):
kml.shapefile_to_kmlfile(self.list, kmlfilename)
else:
success = False
return success
#===========================================================================#
#---------------------------------------------------------------------------#
# General Geographic Positioning Set Environments functions #
#---------------------------------------------------------------------------#
[docs]def refsys_getlist():
"""
Returns list of reference systems availables, 'UTM', 'WGS84', ...
"""
return refsystem_list
[docs]def type_getlist():
"""
Returns list of geographic positioning type, 'point', 'shapefile', ...
"""
return type_list
[docs]def utm_to_wgs84(easting, northing, zonenumber, zoneletter):
"""
Converts UTM to WGS84 coordinates
Parameters :
:easting: east UTM coordinate
:northing: north UTM coordinate
:zonenumber: zone number
:zone letter: zone letter
Returns :
:latitude: WGS84 latitude coordinate
:longitude: WGS84 longitude coordinate
"""
return utm.to_latlon(easting, northing, zonenumber, zoneletter)
[docs]def wgs84_to_utm(latitude, longitude):
"""
Converts WGS84 to UTM coordinates
Parameters :
:latitude: WGS84 latitude coordinate
:longitude: WGS84 longitude coordinate
Returns :
:easting: east UTM coordinate
:northing: north UTM coordinate
:zonenumber: zone number
:zone letter: zone letter
"""
return utm.from_latlon(latitude, longitude)
[docs]def utm_getzonelimits():
"""
Gets the min and max numbers and letters of the UTM coordinates system
Returns:
:min_number: minimal number of the UTM zone, 1.
:min_letter: minimal letter of the UTM zone, E.
:max_numbre: maximal number of the UTM zone, 60.
:max letter: maximal letter of the UTM zone, X.
"""
return utm_minnumber, utm_minletter, utm_maxnumber, utm_maxletter