You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
93 lines
2.4 KiB
Python
93 lines
2.4 KiB
Python
5 years ago
|
"""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()
|