Add polyline_to_points.py

polyline_to_points.py

@@ -0,0 +1,92 @@
+"""polyline_to_points.py
+Extract interpolated points along transects in a shapefile.
+
+D. Howe
+d.howe@wrl.unsw.edu.au
+2020-02-19
+"""
+import sys
+import argparse
+import numpy as np
+import pandas as pd
+import geopandas as gpd
+from shapely.geometry import LineString
+
+
+def extract(shp_path, spacing=1, field=None):
+
+    rows = []
+    shp = gpd.read_file(shp_path)
+
+    if field is None:
+        # Assume profile name is second field in shapefile
+        field = shp.columns[1]
+
+    for i, line in shp.iterrows():
+        g = line.geometry
+        chainages = np.arange(0, g.length, step=spacing)
+        for chainage in chainages:
+            easting, northing = g.interpolate(chainage).xy
+
+            row = {
+                'ProfileNum': line[field],
+                'Easting': easting[0],
+                'Northing': northing[0],
+                'Chainage': chainage,
+            }
+            rows.append(row)
+
+    # Create output dataframe
+    df = pd.DataFrame(rows)
+
+    # Re-order columns
+    df = df[['ProfileNum', 'Easting', 'Northing', 'Chainage']]
+
+    # Export to csv
+    csv_path = shp_path.replace('.shp', '.csv')
+    df.to_csv(csv_path, index=False)
+
+
+def main():
+    example_text = """examples:
+    # Extract points at default spacing (1m)
+    $ python polyline_to_points.py path/to/shp
+
+    # Extract points at 5m increments
+    $ python polyline_to_points.py -s 5 path/to/shp
+
+    # Extract points from field "ProfileID"
+    $ python polyline_to_points.py -f ProfileID path/to/shp
+    """
+    # Set up command line arguments
+    parser = argparse.ArgumentParser(
+        epilog=example_text,
+        formatter_class=argparse.RawDescriptionHelpFormatter)
+    parser.add_argument('shp_path',
+                        metavar='SHAPEFILE',
+                        help='path to input shapefile')
+    parser.add_argument('-s',
+                        '--spacing',
+                        metavar='SPACING',
+                        default=1,
+                        type=int,
+                        help='space between points (default=1)')
+    parser.add_argument('-f',
+                        '--field',
+                        metavar='FIELDNAME',
+                        help='profile field name in attribute table')
+
+    # Print usage if no arguments are provided
+    if len(sys.argv) == 1:
+        parser.print_help(sys.stderr)
+        sys.exit(1)
+
+    # Parse arguments
+    args = parser.parse_args()
+
+    # Combine images
+    extract(**vars(args))
+
+
+if __name__ == '__main__':
+    main()