#!/usr/bin/env python
# -*- coding: utf-8 -*-
import math
from shapely.geometry import shape, Point
import json


def pairs(lst):
    """
    yield iterator of two coordinates of linestring
    :param lst: list object
    :return: yield iterator of two coordinates
    """
    for i in range(1, len(lst)):
        yield lst[i - 1], lst[i]


def calc_3d_distance_2pts(x1, y1, z1, x2, y2, z2):
    """
    :input two point coordinates (x1,y1,z1),(x2,y2,2)
    :param x1: x coordinate first segment
    :param y1: y coordiante first segment
    :param z1: z height value first coordinate
    :param x2: x coordinate second segment
    :param y2: y coordinate second segment
    :param z2: z height value seconc coordinate
    :return: 3D distance between two input 3D coordinates
    """
    d = math.sqrt((x2 - x1) ** 2 + (y2 - y1) ** 2 + (z2 - z1) ** 2)
    return d


def readin_json(jsonfile):
    """
    input: geojson or json file
    """
    with open(jsonfile) as json_data:
        d = json.load(json_data)
        return d


geoj_27563_file = "../geodata/velowire_stage_16_27563.geojson"

# create python dict type from geojson file object
json_load = readin_json(geoj_27563_file)

# set start lengths
length_3d = 0.0
length_2d = 0.0

# go through each geometry in our linestring
for f in json_load['features']:
    # create shapely shape from geojson
    s = shape(f['geometry'])
    geo = f['geometry']

    # calculate 2D total length
    length_2d = s.length

    # set start elevation
    elevation_gain = 0

    # go through each coordinate pair
    for vert_start, vert_end in pairs(s.coords):
        line_start = Point(vert_start)
        line_end = Point(vert_end)

        # create
        x1 = line_start.coords[0][0]
        y1 = line_start.coords[0][1]
        z1 = line_start.coords[0][2]
        x2 = line_end.coords[0][0]
        y2 = line_end.coords[0][1]
        z2 = line_end.coords[0][2]

        distance = calc_3d_distance_2pts(x1, y1, z1, x2, y2, z2)
        length_3d += distance

        # calculate total elevation gain
        if z1 > z2:
            elevation_gain = ((z1 - z2) + elevation_gain )
            z2 = z1
        else:
            elevation_gain = elevation_gain # no height change
            z2 = z1


print "total elevation gain is: " + str(elevation_gain) + " meters"

# print coord_pair
print "3D line distance is: " + str(length_3d / 1000)
print "2D line distance is: " + str(length_2d / 1000)
print "difference 2D and 3D is: " + str(length_3d - length_2d) + " meters"
@@ -804,25 +804,24 @@
 (x2 - x1) **
- 
 2 + (y2 - y1
@@ -824,17 +824,16 @@
 - y1) **
- 
 2 + (z2 
@@ -840,17 +840,16 @@
 - z1) **
- 
 2)%0A    r
@@ -1690,16 +1690,34 @@
 # create
+ input coordinates
 %0A       
@@ -1926,16 +1926,48 @@
 %5B0%5D%5B2%5D%0A%0A
+        # calculate 3d distance%0A
         
@@ -2022,16 +2022,63 @@
  y2, z2)
+%0A%0A        # sum distances from vertex to vertex
 %0A       
@@ -2395,16 +2395,24 @@
 on_gain)
+ %5C%0A     
  + %22 met
@@ -2555,31 +2555,31 @@
 nt %22
-difference 2D and 3D is
+3D-2D length difference
 : %22 
@@ -2606,16 +2606,24 @@
 ngth_2d)
+ %5C%0A     
  + %22 met

clean code