Find which points are exposed depending on the tide height and a DEM. This code is part of a larger project and is provided as an example only.

tide_access.py

from pathlib import Path
import pandas as pd
import geopandas as gpd
from shapely.geometry import Point


def main():
    poi = get_poi_csv()
    # print(poi)

    points_centroid = get_centre_of_points(poi)
    # print(points_centroid)

    # create dummy data for tide_height variable, set at Mean Sea Level:
    tide_height = 0

    # load dummy dataset for dem_data_extent
    dem_data_extent_shapefile_path = Path('DATA/NIdem_data_extent.shp')
    dem_data_extent = gpd.read_file(dem_data_extent_shapefile_path, crs={'init': 'epsg:3577'})

    # load dummy data for extent of submerged part of intertidal zone,
    # given a tide height of 0m above Mean Sea Level:
    submerged_extent_shapefile_path = Path('DATA/submerged_extent.shp')
    submerged_extent = gpd.read_file(submerged_extent_shapefile_path, crs={'init': 'epsg:3577'})

    print("\nTESTING OUT OF SCOPE DATA EXCLUSION...")
    poi = buffer_poi(poi, 1)
    analysis_poi = create_analysis_poi(poi, dem_data_extent)

    print("\nTESTING SUBMERGED poi ANALYSIS...")
    submerged_poi = find_submerged_poi(analysis_poi, submerged_extent)

    print("\nTESTING EXPOSED poi ANALYSIS...")
    analysis_poi_50m = buffer_poi(analysis_poi, 50)
    exposed_poi = find_exposed_poi(analysis_poi_50m, submerged_extent)


def get_poi_csv():
    # read poi csv file into a GeoDataFrame. Reproject into GDA94 Australian Albers coordinate system.
    try:
        poi_path = input("Please specify the path of your Points of Interest file. \n"
                         "This should be a .csv file, with three columns with the following headers: \n"
                         "name, lat, lon. Please ensure the coordinates are in the EPSG:4326 coordinate system. \n"
                         "Lat and lon values should be in decimal degrees.\n"
                         "Enter your filepath here:")

        poi_path = Path(poi_path)

    except IOError as error_details:
        print('IOError, Goodbye!\n', error_details)
    except ValueError as error_details:
        print('ValueError, Goodbye!\n', error_details)
    except TypeError as error_details:
        print('TypeError, Goodbye!\n', error_details)

    try:
        # read the csv into a GeoPandas GeoDataFrame object.
        csv_df = pd.read_csv(poi_path)
        geometry = [Point(x, y) for x, y in zip(csv_df.lon, csv_df.lat)]
        # csv_gdf = gpd.GeoDataFrame(csv_df)
        csv_gdf = gpd.GeoDataFrame(csv_df, crs={'init': 'epsg:4326'})

        csv_gdf = csv_gdf.set_geometry(geometry)
        csv_gdf_reproject = csv_gdf.to_crs(crs={'init': 'epsg:3577'})
        csv_gdf2 = gpd.GeoDataFrame(csv_gdf_reproject)
        csv_gdf2.crs = {'init': 'epsg:3577'}

        return csv_gdf2
    except IOError as error_details:
        print('IOError, Goodbye!\n', error_details)
    # informs user to check how they have named the columns in their poi csv file:
    except AttributeError as error_details:
        print('Possible column naming error in csv file, Goodbye!\n', error_details)


def get_centre_of_points(gdf):
    # gets the center of the points of interest, to help choose an appropriate tide measurement station.
    allpoints = gdf.unary_union
    centroid_x = allpoints.centroid.x
    centroid_y = allpoints.centroid.y
    return (centroid_x, centroid_y)


def buffer_poi(poi, m):
    # creates a very small buffer around pois to facilitate overlay analysis
    poi['geometry'] = poi.geometry.buffer(m)
    return poi


def create_analysis_poi(poi, extent):
    # checks the poi for points outside of the intertidal zone as specified by the NIdem dataset.
    count_poi = poi.count()['geometry']

    analysis_poi = gpd.overlay(poi, extent, how='intersection')
    count_analysis_poi = analysis_poi.count()['geometry']

    if count_poi != count_analysis_poi:
        print("Of your {} total points of interest, only {} fall within the intertidal \n"
              "zone as defined by the NIdem dataset. Points that fall outside the \n"
              "intertidal zone will be excluded from analysis.".format(count_poi, count_analysis_poi))

    else:
        print("All {} points of interest fall within the intertidal zone as defined by the NIdem dataset.\n"
              "All points included in analysis.".format(count_poi))

    return analysis_poi


def find_submerged_poi(poi, extent):
    # finds which poi will be submerged at time of analysis.
    count_poi = poi.count()['geometry']

    submerged_poi = gpd.overlay(poi, extent, how='intersection')
    count_submerged_poi = submerged_poi.count()['geometry']

    print("Of the {} points in the analysis, {} fall within the currently submerged part \n"
          "of the intertidal zone.".format(count_poi, count_submerged_poi))

    return submerged_poi


def find_exposed_poi(poi, extent):
    # finds which points will be exposed and accessible at time of analysis.
    count_poi = poi.count()['geometry']

    exposed_poi = gpd.overlay(poi, extent, how='difference')
    count_exposed_poi = exposed_poi.count()['geometry']

    print("Of the {} points in the analysis, {} fall within the currently exposed part \n"
          "of the intertidal zone.".format(count_poi, count_exposed_poi))

    return exposed_poi


main()