diff --git a/.gitignore b/.gitignore index 796644d..44d6ac3 100644 --- a/.gitignore +++ b/.gitignore @@ -131,8 +131,9 @@ dmypy.json # others .idea logs/* -geobhumi/__pycache__ -geobhumi/geo_app/__pycache__ -geobhumi/raster_ops/__pycache__ -geobhumi/vector_ops/__pycache__ -tests/__pycache__ \ No newline at end of file +rsgtools/__pycache__ +rsgtools/geo_app/__pycache__ +rsgtools/raster_ops/__pycache__ +rsgtools/vector_ops/__pycache__ +tests/__pycache__ +tests/results/* diff --git a/README.md b/README.md index 9a1a9cc..1ba6859 100644 --- a/README.md +++ b/README.md @@ -48,6 +48,7 @@ The `geoutils-rsg` Python package aims to simplify geospatial tasks and automate | :--------------------------------- | :--------------------------------------------------------------------------------------------- | | generate_shoreline_raster() | Shoreline Extraction | | generate_morphometric_parameters() | Morphometric Analysis for Prioritizing Sub-watershed and Management Using Geospatial Technique | +| shoreline_change_analysis() | Digital Shoreline Change-Rate Analysis (example: `notebooks/shoreline_change_rate.ipynb`) | # Installation diff --git a/consts.py b/consts.py new file mode 100644 index 0000000..08ac7de --- /dev/null +++ b/consts.py @@ -0,0 +1,7 @@ +# constants values + +DEG_TO_KM = 111.0 + + +# References +# https://education.nationalgeographic.org/resource/latitude/ diff --git a/docs/functions.md b/docs/functions.md index 2d9d880..327240e 100644 --- a/docs/functions.md +++ b/docs/functions.md @@ -14,6 +14,10 @@ offset_random_point read_raster_as_array np2gdal_dtype write_geotiff_file +timeit +Histogram +OTSU +get_shapefile_epsg_code ``` ### ref_scripts.py diff --git a/docs/notes.md b/docs/notes.md index 490e70b..a13a00c 100644 --- a/docs/notes.md +++ b/docs/notes.md @@ -1,9 +1,41 @@ # Notes: -- convert ipynb to py file in vscode: ctrl+shift+p --> type Export --> select Export to Python Script -- project layouts: https://realpython.com/python-application-layouts/ -- python project structure: https://github.com/yngvem/python-project-structure +- Convert ipynb to py file in vscode: ctrl+shift+p --> type Export --> select Export to Python Script +- Project layouts: https://realpython.com/python-application-layouts/ +- Python project structure: https://github.com/yngvem/python-project-structure - Reset GitHub repository URL on local system if you changed it: https://stackoverflow.com/questions/30443333/error-with-renamed-repo-in-github-remote-this-repository-moved-please-use-th - To check the current one: git remote -v - If it was origin, Then you change it like: git remote set-url origin https://github.com/YOUR-USERNAME/YOUR-REPO - If it was upstream, you change it like: git remote set-url upstream https://github.com/YOUR-USERNAME/YOUR-REPO +- Git Commands: + - Run this command to see what branch you're on: `git status` + - Run this command to see local branches: `git branch` + - Run this command to see remote branches: `git branch -r` + - Run this command to see all local and remote branches: `git branch -a` + - Create a new branch: + - `git checkout -b ` + - Switch to a branch in your local repo: `git checkout ` + - Run this command to get a list of all branches from the remote: `git pull` + - Run this command to switch to a branch that came from a remote repo: `git checkout --track origin/` + - Push to a branch if local branch already exists on the remote: `git push origin ` or `git push` + - Push to a branch if local branch does not exist on the remote: `git push -u origin ` or `git push -u origin HEAD` + - Rename git branch: + - To rename the current branch: `git branch -m ` + - To rename a branch while pointed to any branch: `git branch -m ` + - Merge a branch: + - `git status` + - `git checkout ` + - `git merge ` + - Delete branches: + - Run this command to delete a remote branch: `git push origin --delete ` + - Run this command to delete a local branch: `git branch -d ` or `git branch -D ` + - More commands: + - `git stash` when you want to change to a different branch, and you want to store changes that are not ready to be committed. + - Steps to follow after every edit in the project: + - `git status` + - `git add .` or `git add ` + - `git commit -m ""` + - `git stash` (if needed) + - `git pull` (if needed) + - `git push origin ` + - `git status` diff --git a/main.py b/main.py new file mode 100644 index 0000000..716c062 --- /dev/null +++ b/main.py @@ -0,0 +1,6 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 18:41:31 2024 + +@author: USER +""" diff --git a/notebooks/data/shorelines/Baseline.cpg b/notebooks/data/shorelines/Baseline.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/Baseline.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/Baseline.dbf b/notebooks/data/shorelines/Baseline.dbf new file mode 100644 index 0000000..9caae50 Binary files /dev/null and b/notebooks/data/shorelines/Baseline.dbf differ diff --git a/notebooks/data/shorelines/Baseline.prj b/notebooks/data/shorelines/Baseline.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/Baseline.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/Baseline.qpj b/notebooks/data/shorelines/Baseline.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/Baseline.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/Baseline.shp b/notebooks/data/shorelines/Baseline.shp new file mode 100644 index 0000000..92e865d Binary files /dev/null and b/notebooks/data/shorelines/Baseline.shp differ diff --git a/notebooks/data/shorelines/Baseline.shx b/notebooks/data/shorelines/Baseline.shx new file mode 100644 index 0000000..161167d Binary files /dev/null and b/notebooks/data/shorelines/Baseline.shx differ diff --git a/notebooks/data/shorelines/OffShoreline.cpg b/notebooks/data/shorelines/OffShoreline.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/OffShoreline.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/OffShoreline.dbf b/notebooks/data/shorelines/OffShoreline.dbf new file mode 100644 index 0000000..3aafc9c Binary files /dev/null and b/notebooks/data/shorelines/OffShoreline.dbf differ diff --git a/notebooks/data/shorelines/OffShoreline.prj b/notebooks/data/shorelines/OffShoreline.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/OffShoreline.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/OffShoreline.qpj b/notebooks/data/shorelines/OffShoreline.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/OffShoreline.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/OffShoreline.shp b/notebooks/data/shorelines/OffShoreline.shp new file mode 100644 index 0000000..cbd42d4 Binary files /dev/null and b/notebooks/data/shorelines/OffShoreline.shp differ diff --git a/notebooks/data/shorelines/OffShoreline.shx b/notebooks/data/shorelines/OffShoreline.shx new file mode 100644 index 0000000..80a20ee Binary files /dev/null and b/notebooks/data/shorelines/OffShoreline.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1973.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1973.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1973.dbf new file mode 100644 index 0000000..669dcb3 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1973.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.prj b/notebooks/data/shorelines/shorelines/Shoreline_1973.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1973.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1973.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1973.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.shp b/notebooks/data/shorelines/shorelines/Shoreline_1973.shp new file mode 100644 index 0000000..08ec7f1 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1973.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1973.shx b/notebooks/data/shorelines/shorelines/Shoreline_1973.shx new file mode 100644 index 0000000..dc68244 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1973.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1977.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1977.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1977.dbf new file mode 100644 index 0000000..3b435b3 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1977.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.prj b/notebooks/data/shorelines/shorelines/Shoreline_1977.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1977.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1977.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1977.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.shp b/notebooks/data/shorelines/shorelines/Shoreline_1977.shp new file mode 100644 index 0000000..c7c055d Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1977.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1977.shx b/notebooks/data/shorelines/shorelines/Shoreline_1977.shx new file mode 100644 index 0000000..37b13bf Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1977.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1980.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1980.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1980.dbf new file mode 100644 index 0000000..15639e3 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1980.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.prj b/notebooks/data/shorelines/shorelines/Shoreline_1980.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1980.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1980.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1980.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.shp b/notebooks/data/shorelines/shorelines/Shoreline_1980.shp new file mode 100644 index 0000000..165c721 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1980.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1980.shx b/notebooks/data/shorelines/shorelines/Shoreline_1980.shx new file mode 100644 index 0000000..1487df1 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1980.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1989.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1989.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1989.dbf new file mode 100644 index 0000000..3db3371 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1989.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.prj b/notebooks/data/shorelines/shorelines/Shoreline_1989.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1989.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1989.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1989.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.shp b/notebooks/data/shorelines/shorelines/Shoreline_1989.shp new file mode 100644 index 0000000..6b9d292 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1989.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1989.shx b/notebooks/data/shorelines/shorelines/Shoreline_1989.shx new file mode 100644 index 0000000..d8d2a65 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1989.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1993.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1993.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1993.dbf new file mode 100644 index 0000000..7e45e3e Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1993.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.prj b/notebooks/data/shorelines/shorelines/Shoreline_1993.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1993.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1993.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1993.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.shp b/notebooks/data/shorelines/shorelines/Shoreline_1993.shp new file mode 100644 index 0000000..fcf31c7 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1993.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1993.shx b/notebooks/data/shorelines/shorelines/Shoreline_1993.shx new file mode 100644 index 0000000..1249ddb Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1993.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.cpg b/notebooks/data/shorelines/shorelines/Shoreline_1997.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1997.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.dbf b/notebooks/data/shorelines/shorelines/Shoreline_1997.dbf new file mode 100644 index 0000000..bd2aedc Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1997.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.prj b/notebooks/data/shorelines/shorelines/Shoreline_1997.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1997.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.qpj b/notebooks/data/shorelines/shorelines/Shoreline_1997.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_1997.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.shp b/notebooks/data/shorelines/shorelines/Shoreline_1997.shp new file mode 100644 index 0000000..1ae9f51 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1997.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_1997.shx b/notebooks/data/shorelines/shorelines/Shoreline_1997.shx new file mode 100644 index 0000000..4ce0c45 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_1997.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.cpg b/notebooks/data/shorelines/shorelines/Shoreline_2003.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2003.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.dbf b/notebooks/data/shorelines/shorelines/Shoreline_2003.dbf new file mode 100644 index 0000000..f6d263b Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2003.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.prj b/notebooks/data/shorelines/shorelines/Shoreline_2003.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2003.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.qpj b/notebooks/data/shorelines/shorelines/Shoreline_2003.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2003.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.shp b/notebooks/data/shorelines/shorelines/Shoreline_2003.shp new file mode 100644 index 0000000..5a093ea Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2003.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2003.shx b/notebooks/data/shorelines/shorelines/Shoreline_2003.shx new file mode 100644 index 0000000..afbb04f Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2003.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.cpg b/notebooks/data/shorelines/shorelines/Shoreline_2008.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2008.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.dbf b/notebooks/data/shorelines/shorelines/Shoreline_2008.dbf new file mode 100644 index 0000000..7407d93 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2008.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.prj b/notebooks/data/shorelines/shorelines/Shoreline_2008.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2008.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.qpj b/notebooks/data/shorelines/shorelines/Shoreline_2008.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2008.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.shp b/notebooks/data/shorelines/shorelines/Shoreline_2008.shp new file mode 100644 index 0000000..77eb5e1 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2008.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2008.shx b/notebooks/data/shorelines/shorelines/Shoreline_2008.shx new file mode 100644 index 0000000..249489f Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2008.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.cpg b/notebooks/data/shorelines/shorelines/Shoreline_2014.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2014.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.dbf b/notebooks/data/shorelines/shorelines/Shoreline_2014.dbf new file mode 100644 index 0000000..982598d Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2014.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.prj b/notebooks/data/shorelines/shorelines/Shoreline_2014.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2014.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.qpj b/notebooks/data/shorelines/shorelines/Shoreline_2014.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2014.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.shp b/notebooks/data/shorelines/shorelines/Shoreline_2014.shp new file mode 100644 index 0000000..bdc8e7e Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2014.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2014.shx b/notebooks/data/shorelines/shorelines/Shoreline_2014.shx new file mode 100644 index 0000000..42039df Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2014.shx differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.cpg b/notebooks/data/shorelines/shorelines/Shoreline_2018.cpg new file mode 100644 index 0000000..3ad133c --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2018.cpg @@ -0,0 +1 @@ +UTF-8 \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.dbf b/notebooks/data/shorelines/shorelines/Shoreline_2018.dbf new file mode 100644 index 0000000..69068e0 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2018.dbf differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.prj b/notebooks/data/shorelines/shorelines/Shoreline_2018.prj new file mode 100644 index 0000000..b9cf562 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2018.prj @@ -0,0 +1 @@ +PROJCS["WGS_1984_UTM_Zone_45N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]] \ No newline at end of file diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.qpj b/notebooks/data/shorelines/shorelines/Shoreline_2018.qpj new file mode 100644 index 0000000..e9f9151 --- /dev/null +++ b/notebooks/data/shorelines/shorelines/Shoreline_2018.qpj @@ -0,0 +1 @@ +PROJCS["WGS 84 / UTM zone 45N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",87],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32645"]] diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.shp b/notebooks/data/shorelines/shorelines/Shoreline_2018.shp new file mode 100644 index 0000000..010a0c1 Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2018.shp differ diff --git a/notebooks/data/shorelines/shorelines/Shoreline_2018.shx b/notebooks/data/shorelines/shorelines/Shoreline_2018.shx new file mode 100644 index 0000000..190ac5c Binary files /dev/null and b/notebooks/data/shorelines/shorelines/Shoreline_2018.shx differ diff --git a/notebooks/shoreline_change_rate.ipynb b/notebooks/shoreline_change_rate.ipynb new file mode 100644 index 0000000..f552e40 --- /dev/null +++ b/notebooks/shoreline_change_rate.ipynb @@ -0,0 +1,188 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "'d:\\\\work\\\\geoutils-rsg'" + ] + }, + "execution_count": 1, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "import os\n", + "os.chdir('../.')\n", + "%pwd" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "baseline = 'notebooks/data/shorelines/Baseline.shp'\n", + "offshoreline = 'notebooks/data/shorelines/OffShoreline.shp'\n", + "shorelines_path = \"notebooks/data/shorelines/shorelines\"\n", + "transects = 'tests/results/transects.shp'\n", + "temp_dir = 'tests/results'\n", + "out_changerate_csvfile = 'tests/results/shoreline_changerate.csv'" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "# generate transects\n", + "from rsgtools.geo_app.generate_shoreline_transects import create_shoreline_transects\n", + "create_shoreline_transects(\n", + " onshore_line=baseline,\n", + " offshore_line=offshoreline,\n", + " out_transect_line=transects.replace('.shp', '_temp.shp'),\n", + " x_interval=100\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "[2024-04-13 00:21:22,116: WARNING: collection]: Normalized/laundered field name: 'index_right' to 'index_righ'\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "d:\\work\\geoutils-rsg\\rsgtools\\utils.py:481: UserWarning: Column names longer than 10 characters will be truncated when saved to ESRI Shapefile.\n", + " selected_features.to_file(out_selected_features)\n" + ] + } + ], + "source": [ + "from rsgtools.geo_app.generate_transects_seqid import yield_transects_seqid\n", + "\n", + "yield_transects_seqid(\n", + " temp_dir, \n", + " baseline, \n", + " transects.replace('.shp', '_temp.shp'),\n", + " transects,\n", + " col_name='SeqID'\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "0\n", + "1\n", + "2\n", + "3\n", + "4\n", + "5\n", + "6\n", + "7\n", + "8\n", + "9\n", + "10\n", + "11\n", + "12\n", + "13\n", + "14\n", + "15\n", + "16\n", + "17\n", + "18\n", + "19\n", + "20\n", + "21\n", + "22\n", + "23\n", + "24\n", + "25\n", + "26\n", + "27\n", + "28\n", + "29\n", + "30\n", + "31\n", + "32\n", + "33\n", + "34\n", + "35\n", + "36\n", + "37\n", + "38\n", + "39\n", + "40\n", + "41\n", + "42\n", + "43\n", + "44\n", + "45\n", + "46\n", + "47\n", + "48\n", + "49\n" + ] + } + ], + "source": [ + "from rsgtools.geo_app.shoreline_change_analysis import generate_shoreline_change_rate\n", + "\n", + "generate_shoreline_change_rate(\n", + " tmp=temp_dir, \n", + " offshoreline=offshoreline, \n", + " baseline=baseline, \n", + " shorelines_path=shorelines_path, \n", + " transects=transects, \n", + " out_csvfile=out_changerate_csvfile, \n", + " interval=100, \n", + " date_sep='-', \n", + " date_format='yyyy-mm-dd'\n", + ")" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.10.10" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/pyproject.toml b/pyproject.toml index a226b6c..b107e9d 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -18,7 +18,7 @@ classifiers = [ ] [tool.setuptools.packages.find] -where = ["geobhumi"] +where = ["rsgtools"] [project.urls] "Homepage" = "https://github.com/dghorai/geoutils-rsg" diff --git a/rsgtools/__init__.py b/rsgtools/__init__.py index cb2dc3c..3d06d53 100644 --- a/rsgtools/__init__.py +++ b/rsgtools/__init__.py @@ -33,8 +33,8 @@ def error_handling(err, sys_err: sys): # create log folder and file project_dir = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) log_dir = os.path.join(project_dir, "logs") -log_file_name = f"{datetime.datetime.now().strftime('%m_%d_%Y_%H_%M_%S')}.log" -log_file = os.path.join(log_dir, log_file_name) +# log_file_name = f"{datetime.datetime.now().strftime('%m_%d_%Y_%H_%M_%S')}.log" +log_filepath = os.path.join(log_dir, "running_logs.log") os.makedirs(log_dir, exist_ok=True) # setup logging @@ -45,9 +45,9 @@ def error_handling(err, sys_err: sys): level=logging.INFO, format=logging_str, handlers=[ - logging.FileHandler(log_file), + logging.FileHandler(log_filepath), logging.StreamHandler(sys.stdout) ] ) -logger = logging.getLogger("GeoBhumiLogger") +logger = logging.getLogger("RSGToolsLogger") diff --git a/rsgtools/config_entity.py b/rsgtools/config_entity.py index 35f8470..623ec83 100644 --- a/rsgtools/config_entity.py +++ b/rsgtools/config_entity.py @@ -7,7 +7,6 @@ class PointConfig: inpoint: Path single_point: tuple point_offset: float - point_file_coordinate_sys: str @dataclass(frozen=True) @@ -24,7 +23,6 @@ class PolylineConfig: linefield: str interval: int offset: int - coordsys: str @dataclass(frozen=True) @@ -52,7 +50,3 @@ class RasterdataConfig: lmin_list: list qcal_min: float qcal_max: float - prefix: str - lulc_raster_file: Path - lulc_out_raster_file: Path - lilc_class_code: int diff --git a/rsgtools/geo_app/calc_shoreline_changerate.py b/rsgtools/geo_app/calc_shoreline_changerate.py new file mode 100644 index 0000000..b2e85e0 --- /dev/null +++ b/rsgtools/geo_app/calc_shoreline_changerate.py @@ -0,0 +1,285 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 18:41:31 2024 + +@author: Debabrata Ghorai, Ph.D. + +SHORELINE CHANGE RATE ANALYSIS +------------------------------ +Measurement of shoreline change rate using following statistical models: +1) Shoreline Change Envelope (SCE) +2) Net Shoreline Movement (NSM) +3) End Point Rate (EPR) +4) Linear Regression Rate-of-Change(LRR) +5) Weighted Linear Regression (WLR) + +Notes: +1) Shoraline Change-Rate Computation - (Reference - DSAS) +2) Considering Leap Year Days (366) and Not a Leap Year Days (365) + +Model references: +http://woodshole.er.usgs.gov/project-pages/DSAS/version4/images/pdf/DSASv4.pdf +http://www.soest.hawaii.edu/coasts/nps/erosionHazards.php +https://woodshole.er.usgs.gov/discus/messages/2/785.html?1404737665 +http://answers.google.com/answers/threadview/id/761806.html + +""" + +# Import modules +import re +import numpy + +from rsgtools.utils import leap_year +from rsgtools import CustomException + + +class ShorelineChangerate: + """ + Required Fields and Index: + a) Transect ID (Index 0) + b) Year ID (Index 1) + c) Date (Index 2) + d) Distance (Index 3) + e) E (Index 4) + + Date Format: MM/DD/YYYY + """ + + def __init__(self): + pass + + def get_data(self, infile): + # Read transect data for the temporal shorelines + fileopen = open(infile, 'r') + datalist1 = list() + for f in fileopen: + dt = re.sub('\n', '', f) + dt = re.split(',', dt) + datalist1.append(dt) + datalist2 = datalist1[1:] + transect = list() + for ts in datalist2: + transect.append(int(ts[0])) + unique = list(set(transect)) + return datalist2, unique + + def _sce(self, y_list): + # ------------------------------ + # SCE - Shoreline Change Envelope + # SCE = Greatest distance between all shorelines + # ----------------------------------------------- + max_dist = max(y_list) + diff_list = list() + for dist in y_list: + diff_list.append(max_dist - dist) + SCE = max(diff_list) + return SCE + + def old_and_yound_time(self, yr_list, yr_dist_list): + # -------------------------------------------- + # NSM - Net Shoreline Movement + # NSM = Distance between oldest and youngest shorelines + # ------------------------------------------------------ + old_year = min(yr_list) + young_year = max(yr_list) + for yr in yr_dist_list: + if old_year == yr[0]: + old_d = yr[1] + old_t = yr[2] + if young_year == yr[0]: + young_d = yr[1] + young_t = yr[2] + return old_d, old_t, young_d, young_t + + def _nsm(self, yr_list, yr_dist_list): + # -------------------------------------------- + # NSM - Net Shoreline Movement + # NSM = Distance between oldest and youngest shorelines + # ------------------------------------------------------ + old_d, _, young_d, _ = self.old_and_yound_time(yr_list, yr_dist_list) + NSM = (young_d - old_d) + return NSM + + def _epr(self, yr_list, yr_dist_list, nsm): + # ----------------------------------------- + # EPR - End Point Rate + # EPR = Distance between two years/Time between oldest and most recent shoreline + # ------------------------------------------------------------------------------- + _, old_t, _, young_t = self.old_and_yound_time(yr_list, yr_dist_list) + EPR = (nsm/(young_t - old_t)) + return EPR + + def _lrr(self, yr_list, xy_list, x_list, y_list, xsqr_list, ysqr_list): + # ---------------------------------------- + # LRR - Linear Regression Rate-of-Change + # LRR = Plot of shoreline date (year) v/s distance from baseline + # --------------------------------------------------------------- + sum_x = sum(x_list) + sum_y = sum(y_list) + sum_xsqr = sum(xsqr_list) + sum_ysqr = sum(ysqr_list) + sum_xy = sum(xy_list) + n = int(numpy.size(yr_list)) + A = float((n*sum_xy) - (sum_x*sum_y)) + B = float((n*sum_xsqr) - (sum_x**2)) + slope = A/B + intercept = float(sum_y - (slope*sum_x))/n + if numpy.sign(intercept) == -1.0: + slp = slope + intp = intercept + else: + slp = slope + intp = intercept + D = float((n*sum_ysqr) - (sum_y**2)) + r = float(A/(numpy.sqrt(B*D))) # Correlation co-efficient + R2 = r**2 # Coefficient of determination + return slp, intp, R2 + + def _wlr(self, w_list, wx_list, wxsqr_list, wy_list, wxy_list, wgt_list): + # -------------------------------------------------- + # WLR - Weighted Linear Regression + # WLR = Plot of shoreline date (year) v/s distance from baseline with weighted values + # ------------------------------------------------------------------------------------ + I = sum(w_list) + II = sum(wx_list) + III = sum(wxsqr_list) + IV = sum(wy_list) + V = sum(wxy_list) + # If m = the slope of the line and b = its intercept, then they are the solutions + # to the following two equations: + # I*b + II*m = IV + # II*b + III*m = V + E = I*III-II**2 + b = (III*IV - II*V)/E + m = (I*V - II*IV)/E + if numpy.sign(b) == -1.0: + slp = m + intp = b + else: + slp = m + intp = b + mean_x = II/I + mean_y = IV/I + sxx_list = list() + syy_list = list() + sxy_list = list() + for j in wgt_list: + sxx_list.append(j[2]*(j[0]-mean_x)**2) + syy_list.append(j[2]*(j[1]-mean_y)**2) + sxy_list.append(j[2]*(j[0]-mean_x)*(j[1]-mean_y)) + s_xx = sum(sxx_list) + s_yy = sum(syy_list) + s_xy = sum(sxy_list) + wr = float(s_xy/(numpy.sqrt(s_xx*s_yy))) + WR2 = wr**2 + return slp, intp, WR2 + + def extract_shoreline_charngerate_params(self, datalist2, u, date_sep='/', date_format='yyyy-mm-dd'): + x_list = list() + y_list = list() + w_list = list() + xsqr_list = list() + ysqr_list = list() + xy_list = list() + yr_list = list() + yr_dist_list = list() + wx_list = list() + wxsqr_list = list() + wy_list = list() + wxy_list = list() + wgt_list = list() + # Loop over the datalist + for q in datalist2: + if u == int(q[0]): + if date_format == 'yyyy-mm-dd': + yyyy = int(re.split(date_sep, q[2])[0]) + mm = int(re.split(date_sep, q[2])[1]) + dd = int(re.split(date_sep, q[2])[2]) + else: + mm = int(re.split(date_sep, q[2])[0]) + dd = int(re.split(date_sep, q[2])[1]) + yyyy = int(re.split(date_sep, q[2])[2]) + time = None + check = leap_year(yyyy) + if check[0] == "Leap year": + time = yyyy + ((((mm-1)*(366.0/12.0)) + dd)/366.0) + elif check[0] == "Not a leap year": + time = yyyy + ((((mm-1)*(365.0/12.0)) + dd)/365.0) + else: + raise CustomException("time is out of bound") + e = float(q[4]) # Index 4 + x, y, w = time, float(q[3]), 1.0/(e**2) + xsqr, ysqr, xy = x**2, y**2, x*y + wx, wxsqr, wy, wxy = w*x, w*x**2, w*y, w*x*y + # append values + x_list.append(x) + y_list.append(y) + w_list.append(w) + xsqr_list.append(xsqr) + ysqr_list.append(ysqr) + xy_list.append(xy) + yr_list.append(yyyy) + yr_dist_list.append((yyyy, y, x)) + wx_list.append(wx) + wxsqr_list.append(wxsqr) + wy_list.append(wy) + wxy_list.append(wxy) + wgt_list.append((x, y, w)) + res1 = [x_list, y_list, w_list] + res2 = [xsqr_list, ysqr_list] + res3 = [xy_list, yr_list, yr_dist_list] + res4 = [wx_list, wxsqr_list, wy_list, wxy_list, wgt_list] + return res1, res2, res3, res4 + + def write_shoreline_changerate(self, crList, out_csvfile): + w = open(out_csvfile, 'w') + w.writelines( + 'TrID,SCE,NSM,EPR,LRR_M,LRR_B,LRR_R2,WLR_M,WLR_B,WLR_WR2\n') + for ii in crList: + w.writelines(",".join(str(jj) for jj in ii)+"\n") + w.close() + return + + def calc_shoreline_changerate_values(self, infile, out_csvfile, date_sep=None, date_format=None): + # Read transect data for the temporal shorelines + datalist2, unique = self.get_data(infile) + # Loop over the unique transects + crList = list() + for u in unique: + print(u) + # get params + res1, res2, res3, res4 = self.extract_shoreline_charngerate_params(datalist2, u, date_sep=date_sep, date_format=date_format) + x_list, y_list, w_list = res1 + xsqr_list, ysqr_list = res2 + xy_list, yr_list, yr_dist_list = res3 + wx_list, wxsqr_list, wy_list, wxy_list, wgt_list = res4 + # shoreline change-rate calculation + SCE = self._sce(y_list) + NSM = self._nsm(yr_list, yr_dist_list) + EPR = self._epr(yr_list, yr_dist_list, NSM) + m1, b1, R2 = self._lrr(yr_list, xy_list, x_list, + y_list, xsqr_list, ysqr_list) + m2, b2, WR2 = self._wlr(w_list, wx_list, wxsqr_list, + wy_list, wxy_list, wgt_list) + # append values + crt = list() + crt.append(u) + crt.append(round(SCE, 2)) + crt.append(round(NSM, 2)) + crt.append(round(EPR, 2)) + crt.append(round(m1, 2)) + crt.append(round(b1, 2)) + crt.append(round(R2, 2)) + crt.append(round(m2, 2)) + crt.append(round(b2, 2)) + crt.append(round(WR2, 2)) + crList.append(crt) + # save + self.write_shoreline_changerate(crList, out_csvfile) + return + + +def generate_shoreline_chargerate_values(in_params_csvfile, out_slcr_csvfile, date_sep=None, date_format=None): + scrt = ShorelineChangerate() + scrt.calc_shoreline_changerate_values(in_params_csvfile, out_slcr_csvfile, date_sep=date_sep, date_format=date_format) + return diff --git a/rsgtools/geo_app/calc_shoreline_dist_from_baseline.py b/rsgtools/geo_app/calc_shoreline_dist_from_baseline.py new file mode 100644 index 0000000..d5458de --- /dev/null +++ b/rsgtools/geo_app/calc_shoreline_dist_from_baseline.py @@ -0,0 +1,112 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 18:41:31 2024 + +@author: Debabrata Ghorai, Ph.D. + +Distance measurement from baseline to past shorelines. + +""" + +from osgeo import ogr +from rsgtools.utils import dist_calc + + +# historical shoreline's distance measurement from baseline +def distance_from_baseline(qset, alist): + dists = [] + for seq in qset: + tlist = [] + for f in alist: + if seq == f[0]: + tlist.append(f[1:]) + base = None + for t in tlist: + if t[1] == 'Baseline': + base = t[3] + cnt = 1 + for t2 in tlist: + if t2[1] != 'Baseline': + dist = dist_calc(base, t2[3]) + dists.append((seq, cnt, t2[0], dist, t2[2])) + cnt += 1 + return dists + + +def write_historical_shoreline_distance(dists, out_csvfile): + # write csv file + w = open(out_csvfile, 'w') + w.writelines('TrID,YearID,Date,Distance,E\n') + for ii in dists: + w.writelines(",".join(str(i) for i in ii)+"\n") + w.close() + return + + +def extract_distance_parameters(layer, fields): + qlist = [] + alist = [] + for i in range(layer.GetFeatureCount()): + feat = layer.GetFeature(i) + # get seqid field + seqid = feat.GetField(fields[0]) + # get date field + date = feat.GetField(fields[1]) + # get E field + e = feat.GetField(fields[2]) + # get name field + name = feat.GetField(fields[3]) + # get geometry + geom = feat.GetGeometryRef() + point = geom.GetPoints() + # append values + alist.append((seqid, date, name, e, point[0])) + qlist.append(seqid) + return qlist, alist + + +def measure_historical_shoreline_distance(shore_transect_xpoint_file, out_csvfile): + """ + Historical shoreline's distance measurement from baseline. + + Required Datasets: + A merged point file of the following layers: + - Baseline and Transect intersects points + - Historical shorelines and Transect intersect points + + Required Fields: + a) SeqID (Index 1) + b) DATE (Index 2) + c) E (Index 3) + d) NAME (Index 4) + """ + + # data source + ds = ogr.Open(shore_transect_xpoint_file) + layer = ds.GetLayer(0) + layerDefn = layer.GetLayerDefn() + + # get spatial reference + sr = layer.GetSpatialRef() + sr.ExportToWkt() + + # get all the fields + fields = [] + for i in range(layerDefn.GetFieldCount()): + field_name = layerDefn.GetFieldDefn(i).GetName() + fields.append(field_name) + + # calculate historical distances + try: + # see the total no. of fields should be four + if len(fields) == 4: + qlist, alist = extract_distance_parameters(layer, fields) + # calculate unique seqid + qset = list(set(qlist)) + # calculate distances + dists = distance_from_baseline(qset, alist) + # write csv file + write_historical_shoreline_distance(dists, out_csvfile) + except: + raise Exception("Check the input datasets and it's attributes") + return diff --git a/rsgtools/geo_app/calc_shoreline_uncertainty.py b/rsgtools/geo_app/calc_shoreline_uncertainty.py new file mode 100644 index 0000000..8d1de7b --- /dev/null +++ b/rsgtools/geo_app/calc_shoreline_uncertainty.py @@ -0,0 +1,4 @@ + +def shoreline_uncertainty_value(): + # page 10-11: https://www.nccr.gov.in/sites/default/files/schangenew.pdf + pass \ No newline at end of file diff --git a/rsgtools/geo_app/generate_shoreline_baseline_intersects.py b/rsgtools/geo_app/generate_shoreline_baseline_intersects.py new file mode 100644 index 0000000..9c3976d --- /dev/null +++ b/rsgtools/geo_app/generate_shoreline_baseline_intersects.py @@ -0,0 +1,94 @@ +import os +import glob +import geopandas as gpd +import pandas as pd + +from shapely.geometry import LineString +from rsgtools.utils import extend_line_shapefile, get_shapefile_epsg_code, unlink_files + + +def get_intersects(baseline, transects, save=False, intersects_points=None): + line1 = gpd.read_file(baseline) + line2 = gpd.read_file(transects) + epsg_code = get_shapefile_epsg_code(transects) + + columns_data = [] + geoms = [] + for i in line1.itertuples(): + date = i.DATE + + try: + name = i.NAME + except: + name = '' + + try: + e_val = i.E + except: + e_val = 0.01 + + r = i.geometry + for j in line2.itertuples(): + seqid = j.SeqID + xg = j.geometry + l_coords = list(xg.coords) + # extent at start + c = LineString(l_coords) + intersect = r.intersection(c) + if intersect.geom_type == 'MultiPoint': + points = [p for p in intersect.geoms] + for pnt in points: + columns_data.append((seqid, date, e_val, name)) + geoms.append(pnt) + else: + columns_data.append((seqid, date, e_val, name)) + geoms.append(intersect) + + gdf = gpd.GeoDataFrame( + columns_data, + crs='epsg:'+str(epsg_code), + geometry=geoms, + columns=['SeqID', 'DATE', 'E', 'NAME'] + ) + + if save: + gdf.to_file(intersects_points) + else: + return gdf + return + + +def create_changerate_analysis_points(temp_dir, transects_with_seqid, baseline, shorelines_path, merge_intersect_points): + # define temporary files + transects_extended = os.path.join(temp_dir, 'tests', 'results', 'transects_seqid_extd.shp') + base_and_tr_intx = os.path.join(temp_dir, 'tests', 'results', 'base_with_transects_intx.shp') + shorelines_and_tr_intx = os.path.join(temp_dir, 'tests', 'results', 'shorelines_with_transects_intx.shp') + + # extent transects lines + extend_line_shapefile(transects_with_seqid, transects_extended, max_snap_dist=25.0, keep_field='SeqID', flip=True) + # get intersection points between baseline and transaction line + get_intersects(baseline, transects_extended, save=True, intersects_points=base_and_tr_intx) + # get intersection points between shorelines and transaction line + shorelines = glob.glob(os.path.join(shorelines_path, '*.shp')) + + gdfs = [] + for shoreline in shorelines: + gdf = get_intersects(shoreline, transects_extended) + gdfs.append(gdf) + + fgdf = pd.concat(gdfs).pipe(gpd.GeoDataFrame) + fgdf.reset_index(drop=True, inplace=True) + fgdf.to_file(shorelines_and_tr_intx) + + # merge both the point files + points1 = gpd.read_file(base_and_tr_intx) + points2 = gpd.read_file(shorelines_and_tr_intx) + intersect_points = pd.concat([points1, points2]).pipe(gpd.GeoDataFrame) + intersect_points.to_file(merge_intersect_points) + + # remove temporary files + unlink_files(is_file=True, file_path=transects_extended) + unlink_files(is_file=True, file_path=base_and_tr_intx) + unlink_files(is_file=True, file_path=shorelines_and_tr_intx) + + return diff --git a/rsgtools/geo_app/generate_shoreline_transects.py b/rsgtools/geo_app/generate_shoreline_transects.py new file mode 100644 index 0000000..4e77304 --- /dev/null +++ b/rsgtools/geo_app/generate_shoreline_transects.py @@ -0,0 +1,61 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 21:29:13 2024 + +@author: Debabrata Ghorai, Ph.D. + +Generate shoreline transects. + +""" + +from osgeo import ogr +from rsgtools.ref_scripts import fixed_interval_points +from rsgtools.vector_ops.find_nearest_point import get_nearest_point + + +def create_shoreline_transects(onshore_line=None, offshore_line=None, out_transect_line=None, x_interval=None): + driver = ogr.GetDriverByName("ESRI Shapefile") + ds_base = driver.Open(onshore_line, 0) + lyr = ds_base.GetLayer(0) + sr = lyr.GetSpatialRef() + sr.ExportToWkt() + ds_xline = driver.CreateDataSource(out_transect_line) + lyr = ds_xline.CreateLayer('myLyr', sr, ogr.wkbMultiLineString) + lyrDef = lyr.GetLayerDefn() + + baseline_points = fixed_interval_points( + onshore_line, x_interval, flipline=True) + + for ix, pnt in enumerate(baseline_points): + # get nearest point + near_pnt = get_nearest_point(pnt, offshore_line) + # get the forward node (projected node) + ix_f = near_pnt[0]+(near_pnt[0]-pnt[0]) + iy_f = near_pnt[1]+(near_pnt[1]-pnt[1]) + fnd = [ix_f, iy_f] + # transect drawing + x_line = ogr.Geometry(ogr.wkbLineString) + x_line.AddPoint(pnt[0], pnt[1]) + x_line.AddPoint(fnd[0], fnd[1]) + # create feature + feature = ogr.Feature(lyrDef) + feature.SetGeometry(x_line) + feature.SetFID(ix+1) + lyr.CreateFeature(feature) + # flush + ds_xline.Destroy() + driver = None + return + + +# if __name__ == '__main__': +# point_interval = 500 +# shore_baseline = '' # baseline +# shore_referenceline = '' # reference line +# out_crosssection_line = "" # Output Transect File +# create_shoreline_transects( +# onshore_line=shore_baseline, +# offshore_line=shore_referenceline, +# out_transect_line=out_crosssection_line, +# x_interval=point_interval +# ) diff --git a/rsgtools/geo_app/generate_transects_seqid.py b/rsgtools/geo_app/generate_transects_seqid.py new file mode 100644 index 0000000..e6e232f --- /dev/null +++ b/rsgtools/geo_app/generate_transects_seqid.py @@ -0,0 +1,160 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 21:29:13 2024 + +@author: Debabrata Ghorai, Ph.D. + +Generate sequence ID of shoreline transects. + +""" + +import os +import geopandas as gpd + +from pathlib import Path +from osgeo import ogr +from rsgtools.vector_ops.lines_middle_points import generate_lines_middle_point +from rsgtools.utils import ( + select_by_location, + reading_polyline, + get_shapefile_epsg_code, + split_lines_at_points, + unlink_files +) +from rsgtools.vector_ops.create_buffer import buffer_feature + + +def create_line_seqid(line_feature, out_line_feature): + # Input shapefile layer + shapefile = ogr.Open(line_feature) + layer = shapefile.GetLayer(0) + # Get input file spatial reference + sr = layer.GetSpatialRef() + sr.ExportToWkt() + # Create new shapefile + driver = ogr.GetDriverByName('ESRI Shapefile') + if os.path.exists(out_line_feature): + driver.DeleteDataSource(out_line_feature) + + # create data source + ds = driver.CreateDataSource(out_line_feature) + # Create new shapefile layer + lyr = ds.CreateLayer('mylyr1', sr, ogr.wkbPoint) + lyrDef = lyr.GetLayerDefn() + # Add the Sequence ID Field + lyr.CreateField(ogr.FieldDefn("SeqID", ogr.OFTInteger)) + + # Loop over the input features + for i in range(layer.GetFeatureCount()): + feature = layer.GetFeature(i) + geom = feature.GetGeometryRef() + pnts = geom.GetPoints() + sgmnt_end = pnts[-1] + # Point file generation + point = ogr.Geometry(ogr.wkbPoint) + point.SetPoint(0, sgmnt_end[0], sgmnt_end[1]) + featIndex = 0 + feat = ogr.Feature(lyrDef) + feat.SetGeometry(point) + feat.SetFID(featIndex) + feat.SetField("SeqID", i+1) + lyr.CreateFeature(feat) + # Flush + ds.Destroy() + ds = None + return + + +def yield_transects_seqid(temp_dir, baseline, transects_lines, out_transects_with_seqid, col_name=None): + """Generate sequence ID of shoreline transects""" + # define temporary files + out_transects = Path(transects_lines) + out_union_line = Path(os.path.join(temp_dir, 'tests', 'results', 'union_xline.shp')) + out_mid_points = Path(os.path.join(temp_dir, 'tests', 'results', 'union_xline_midpoints.shp')) + out_union_line_selected = Path(os.path.join(temp_dir, 'tests', 'results', 'union_xline_selected.shp')) + out_seqid_xpoints = Path(os.path.join(temp_dir, 'tests', 'results', 'xpoints_seqid.shp')) + + # check exist or not + out_transects.parent.mkdir(parents=True, exist_ok=True) + out_union_line.parent.mkdir(parents=True, exist_ok=True) + out_mid_points.parent.mkdir(parents=True, exist_ok=True) + out_union_line_selected.parent.mkdir(parents=True, exist_ok=True) + out_seqid_xpoints.parent.mkdir(parents=True, exist_ok=True) + + # convert to realpath + out_transects = os.path.realpath(out_transects) + out_union_line = os.path.realpath(out_union_line) + out_mid_points = os.path.realpath(out_mid_points) + out_union_line_selected = os.path.realpath(out_union_line_selected) + out_seqid_xpoints = os.path.realpath(out_seqid_xpoints) + + # get src epsg code + src_epsg_code = get_shapefile_epsg_code(out_transects) + + # merge and split shoreline and trsancets lines + # https://geopandas.org/en/stable/docs/reference/api/geopandas.GeoDataFrame.html + _, geoms = reading_polyline(out_transects) + + points = [] + for geom in geoms: + points.append(geom[0]) + points.append(geom[-1]) + + gdf_segments = split_lines_at_points( + baseline, points, snap_dist=10, src_epsg_code=src_epsg_code) + gdf2 = gpd.read_file(out_transects) + + # union line layers + res_union = gpd.overlay( + gdf_segments, gdf2, how='union', keep_geom_type=True) + res_union.to_file(out_union_line) + + # create mid-point of lines + generate_lines_middle_point( + line_feature=out_union_line, out_mid_points=out_mid_points) + + # Create 1 meter buffer for Baseline + baseline_buffer = buffer_feature( + baseline, buffer_offset=1.0, src_epsg_code=src_epsg_code) + baseline_buffer_gdf = gpd.GeoDataFrame(geometry=gpd.GeoSeries(baseline_buffer)) + # Select points based on the Buffer Baseline Polygon + point_in_poly = select_by_location( + select_feature_from=out_mid_points, feature_overlap_with=baseline_buffer_gdf, overlap_method='within') + # Create 1 meter buffer for selected point file + spoints_buffer = buffer_feature( + point_in_poly, buffer_offset=1.0, src_epsg_code=src_epsg_code) + spoints_buffer_gdf = gpd.GeoDataFrame(geometry=gpd.GeoSeries(spoints_buffer)) + + # Select Base-Trans-Splited files based on selected points file polygon + select_by_location( + select_feature_from=out_union_line, + feature_overlap_with=spoints_buffer_gdf, + overlap_method='intersects', + save=True, + out_selected_features=out_union_line_selected + ) + + # assign sequence id to shoreline points + create_line_seqid(out_union_line_selected, out_seqid_xpoints) + + # k) Create 1 meter buffer on the (j) created point file + seqid_point_buffer = buffer_feature( + out_seqid_xpoints, buffer_offset=1.0, src_epsg_code=src_epsg_code) + seqid_point_buffer_gdf = gpd.GeoDataFrame(geometry=gpd.GeoSeries(seqid_point_buffer)) + + # l) Spatial join between transects and (k) feature class + select_by_location( + select_feature_from=out_transects, + feature_overlap_with=seqid_point_buffer_gdf, + overlap_method='intersects', + save=True, + out_selected_features=out_transects_with_seqid, + col_name=col_name + ) + + # remove temporary files + unlink_files(is_file=True, file_path=out_union_line) + unlink_files(is_file=True, file_path=out_mid_points) + unlink_files(is_file=True, file_path=out_union_line_selected) + unlink_files(is_file=True, file_path=out_seqid_xpoints) + return diff --git a/rsgtools/geo_app/shoreline_change_analysis.py b/rsgtools/geo_app/shoreline_change_analysis.py new file mode 100644 index 0000000..1dbf369 --- /dev/null +++ b/rsgtools/geo_app/shoreline_change_analysis.py @@ -0,0 +1,121 @@ +#!/usr/bin/env python +""" +**************************************************************************** +* File Name: shoreline_change_analysis.py +* Created: April 12, 2024 +* Author: Debabrata Ghorai, Ph.D. +* Purpose: Digital Shoreline Change Rate Analysis System. +* Description: Create shoreline transects, generate sequence id of the transects, create change rate analysis inputs, and + calculate shoreline change rate value with different models (SCE, NSM, EPR, LRR_M, LRR_B, LRR_R2, WLR_M, WLR_B, WLR_WR2). +* Research Paper: Ghorai, D., & Bhunia, G. S. (2022). Automatic shoreline detection and its forecast: a case study + on Dr. Abdul Kalam Island in the section of Bay of Bengal. Geocarto International, 37(8), 2273-2292. + Available at: https://www.tandfonline.com/doi/abs/10.1080/10106049.2020.1815868 +* Revesions: N/A +**************************************************************************** +""" + + +import os +from rsgtools.utils import unlink_files +from rsgtools.geo_app.generate_shoreline_transects import create_shoreline_transects +from rsgtools.geo_app.generate_transects_seqid import yield_transects_seqid +from rsgtools.geo_app.generate_shoreline_baseline_intersects import create_changerate_analysis_points +from rsgtools.geo_app.calc_shoreline_dist_from_baseline import measure_historical_shoreline_distance +from rsgtools.geo_app.calc_shoreline_changerate import generate_shoreline_chargerate_values + + +class ShorelineChangeAnalysis: + """ + baseline: + DATE: Null + NAME: Baseline + + past shorelines: + DATE: YYYY-MM-DD + E: Shoreline Uncertainty Value + + transects: + SeqID: Value of SequenceID + """ + + def __init__(self, tmp, offshoreline, baseline, shorelines_path, transects, out_csvfile, interval=100, date_sep='-', date_format='yyyy-mm-dd'): + self.temp_dir = tmp + self.baseline = baseline + self.offshoreline = offshoreline + self.shorelines_path = shorelines_path + self.x_interval = interval + self.transects_with_seqid = transects + self.out_changerate_csvfile = out_csvfile + self.transects = os.path.join(tmp, transects.replace('.shp', '_temp.shp')) + self.merge_intersect_points = os.path.join(tmp, 'tests', 'results', 'merge_points.shp') + self.shoreline_change_params = os.path.join(tmp, 'tests', 'results', 'shoreline_distances.csv') + self.date_sep = date_sep + self.date_format = date_format + + def run(self): + # generate transects + # create_shoreline_transects( + # onshore_line=self.baseline, + # offshore_line=self.offshoreline, + # out_transect_line=self.transects, + # x_interval=self.x_interval + # ) + # add sequence id + # yield_transects_seqid( + # self.temp_dir, + # self.baseline, + # self.transects, + # self.transects_with_seqid, + # col_name='SeqID' + # ) + # generate intersects points for shorelines and baseline + create_changerate_analysis_points( + self.temp_dir, + self.transects_with_seqid, + self.baseline, + self.shorelines_path, + self.merge_intersect_points + ) + # calculate shoreline distance from baseline + measure_historical_shoreline_distance( + self.merge_intersect_points, + self.shoreline_change_params + ) + # calculate shoreline change rate + generate_shoreline_chargerate_values( + self.shoreline_change_params, + self.out_changerate_csvfile, + date_sep=self.date_sep, + date_format=self.date_format + ) + # remove temporary files + unlink_files(is_file=True, file_path=self.transects) + unlink_files(is_file=True, file_path=self.merge_intersect_points) + unlink_files(is_file=True, file_path=self.shoreline_change_params) + return + + +def generate_shoreline_change_rate( + tmp=None, + offshoreline=None, + baseline=None, + shorelines_path=None, + transects=None, + out_csvfile=None, + interval=None, + date_sep=None, + date_format=None +): + sca = ShorelineChangeAnalysis( + tmp, + offshoreline, + baseline, + shorelines_path, + transects, + out_csvfile, + interval=interval, + date_sep=date_sep, + date_format=date_format + ) + sca.run() + return \ No newline at end of file diff --git a/rsgtools/ref_scripts.py b/rsgtools/ref_scripts.py index 9655f37..4285272 100644 --- a/rsgtools/ref_scripts.py +++ b/rsgtools/ref_scripts.py @@ -5,12 +5,14 @@ from osgeo import ogr, gdal from pyproj import Proj from scipy.spatial import ConvexHull +from consts import DEG_TO_KM from rsgtools.utils import ( dist_calc, flip_line, reading_polyline, line_slope, - read_raster_as_array + read_raster_as_array, + get_shapefile_epsg_code ) @@ -56,43 +58,34 @@ def intersect_point_to_line(point, line_start, line_end): ix = dist_calc(point, line_start) iy = dist_calc(point, line_end) if ix > iy: - return line_end + res = line_end else: - return line_start + res = line_start else: ix = line_start[0]+u*(line_end[0]-line_start[0]) iy = line_start[1]+u*(line_end[1]-line_start[1]) - return ix, iy + res = (ix, iy) + + return res # [3, 4] -def fixed_interval_points(infc, interval, flipline=False, save=False, outfc=None, lineslope=False): - """Generate fixed interval points along polyline. - """ +def mid_point(p1, p2, l1, l2): + iX = p1[0]+((p2[0]-p1[0])*(l1/l2)) + iY = p1[1]+((p2[1]-p1[1])*(l1/l2)) + pts = [iX, iY] + return pts - def mid_point(p1, p2, l1, l2): - iX = p1[0]+((p2[0]-p1[0])*(l1/l2)) - iY = p1[1]+((p2[1]-p1[1])*(l1/l2)) - pts = [iX, iY] - return pts - # get nodes and objectids - if flipline == True: - openfile = ogr.Open(infc) - in_layer = openfile.GetLayer(0) - objects, nodes = flip_line(in_layer) - else: - objects, nodes = reading_polyline(infc) - - # generate points +def interval_point(objects, nodes, interval): points = [] slopes = [] for n, i in zip(objects, nodes): tDist = 0 # total distance vPoint = None # previous point for pnt in i: - if not (vPoint is None): - thisDist = dist_calc(vPoint, pnt) + if not isinstance(vPoint, type(None)): + thisDist = dist_calc(list(vPoint), list(pnt)) maxAddDist = interval - tDist if (tDist+thisDist) > interval: pCnt = int((tDist+thisDist)/interval) @@ -100,11 +93,11 @@ def mid_point(p1, p2, l1, l2): maxAddDist = interval - tDist nPoint = mid_point( vPoint, pnt, maxAddDist, thisDist) - slope = line_slope(vPoint, nPoint) + slope = line_slope(list(vPoint), nPoint) points.append(nPoint+[n]) slopes.append(slope) vPoint = nPoint - thisDist = dist_calc(vPoint, pnt) + thisDist = dist_calc(list(vPoint), list(pnt)) tDist = 0 tDist += thisDist else: @@ -112,8 +105,31 @@ def mid_point(p1, p2, l1, l2): else: tDist = 0 vPoint = pnt + return points, slopes + + +def fixed_interval_points(infc, interval, flipline=False, save=False, outfc=None, lineslope=False): + """Generate fixed interval points along polyline. + """ + + # check scr projection and change interval accordingly + src_epsg_code = get_shapefile_epsg_code(infc) + if src_epsg_code == 4326: + interval = interval/DEG_TO_KM/1000.0 + + # get nodes and objectids + if flipline == True: + openfile = ogr.Open(infc) + in_layer = openfile.GetLayer(0) + objects, nodes = flip_line(in_layer) + else: + objects, nodes = reading_polyline(infc) + + # generate points + points, slopes = interval_point(objects, nodes, interval) # save if needed + res = None if save == True: openfile = ogr.Open(infc) in_layer = openfile.GetLayer(0) @@ -141,12 +157,11 @@ def mid_point(p1, p2, l1, l2): # Flush shapeData.Destroy() - - # return points - if lineslope == True: - return points, slopes + elif lineslope == True: + res = (points, slopes) else: - return points + res = points + return res # [5] @@ -328,6 +343,19 @@ def extract_pixel_value(imgfile, point_list): return pnt_values +# [13] +def find_wgs2utm_epsg_code(lon: float, lat: float): + """Based on lat and lng, return best utm epsg-code""" + utm_band = str((math.floor((lon + 180) / 6) % 60) + 1) + if len(utm_band) == 1: + utm_band = '0'+utm_band + if lat >= 0: + epsg_code = '326' + utm_band + return epsg_code + epsg_code = '327' + utm_band + return epsg_code + + # Reference: # [1] https://gis.stackexchange.com/questions/57964/get-vector-features-inside-a-specific-extent # [2] https://gis.stackexchange.com/questions/396/nearest-neighbor-between-point-layer-and-line-layer @@ -341,6 +369,7 @@ def extract_pixel_value(imgfile, point_list): # [10] https://stackoverflow.com/questions/24467972/calculate-area-of-polygon-given-x-y-coordinates # [11] http://stackoverflow.com/questions/13416764/clipping-raster-image-with-a-polygon-suggestion-to-resolve-an-error-related-to # [12] http://gis.stackexchange.com/questions/46893/how-do-i-get-the-pixel-value-of-a-gdal-raster-under-an-ogr-point-without-numpy +# [13] https://gis.stackexchange.com/questions/269518/auto-select-suitable-utm-zone-based-on-grid-intersection # https://gis.stackexchange.com/questions/392515/create-a-shapefile-from-geometry-with-ogr # https://www.gis.usu.edu/~chrisg/python/2009/lectures/ospy_slides2.pdf diff --git a/rsgtools/utils.py b/rsgtools/utils.py index 08a4c98..488ca75 100644 --- a/rsgtools/utils.py +++ b/rsgtools/utils.py @@ -12,21 +12,32 @@ import random import datetime import numpy as np +import geopandas as gpd +import pandas as pd +import glob +import collections +import operator from osgeo import ogr, osr, gdal, gdalconst +from shapely.geometry import Point, mapping, LineString +from shapely.ops import split, snap +from fiona import collection +from fiona.crs import from_epsg from ensure import ensure_annotations -from typing import Any +# from typing import Any from rsgtools import CustomException, logger +from consts import DEG_TO_KM +# Set GDAL/Geopandas Configuration +gdal.SetConfigOption('SHAPE_RESTORE_SHX', 'YES') -@ensure_annotations -def dist_calc(p1: list, p2: list) -> float: + +def dist_calc(p1, p2): """Calculate distance between two points""" return math.sqrt((p1[0]-p2[0])**2+(p1[1]-p2[1])**2) -@ensure_annotations -def line_slope(p1: list, p2: list) -> Any: +def line_slope(p1, p2): """Calculate slope of a line""" return "Inf" if p1[0] == p2[0] else ((p2[1]-p1[1])/(p2[0]-p1[0])) @@ -109,50 +120,54 @@ def reading_polyline(inshp, shptype='polyline', fieldname=None, isnode=False): # return results if fieldname != None: - if isnode == True: - return nodes, objects, fdisc + if isnode: + res = (nodes, objects, fdisc) else: - return fdisc + res = fdisc else: - return objects, nodes + res = (objects, nodes) + + return res -def writting_polyline(xylists, outfile, inlayer='', epsgno=0): +def write_line_shapefile(line_coordinates, out_line_file, ref_line_file=None, dst_epsg_code=None): """ Create new line shapefile. """ # Set up the shapefile driver driver = ogr.GetDriverByName("ESRI Shapefile") + if os.path.exists(out_line_file): + driver.DeleteDataSource(out_line_file) # create the data source - ds = driver.CreateDataSource(outfile) + ds = driver.CreateDataSource(out_line_file) - if len(inlayer) > 0: + if ref_line_file != None: + openfile = ogr.Open(ref_line_file) + ref_lyr = openfile.GetLayer(0) # get the spatial reference system - srs = inlayer.GetSpatialRef() - # create one layer - layer = ds.CreateLayer("line", srs, ogr.wkbLineString) - elif epsgno > 0: + srs = ref_lyr.GetSpatialRef() + elif dst_epsg_code != None: # create the spatial reference system srs = osr.SpatialReference() - srs.ImportFromEPSG(epsgno) - # create one layer - layer = ds.CreateLayer("line", srs, ogr.wkbLineString) + srs.ImportFromEPSG(dst_epsg_code) else: raise CustomException( - "SpatialReference not define, pass inlayer or epsgno arguments to the function") + "SpatialReference not define, pass inlayer or dst_epsg_code arguments to the function") + + # create one layer + dst_layer = ds.CreateLayer('mylyr', srs, ogr.wkbLineString) # Add an ID field idField = ogr.FieldDefn("id", ogr.OFTInteger) - layer.CreateField(idField) - + dst_layer.CreateField(idField) # Create the feature and set values - featureDefn = layer.GetLayerDefn() + featureDefn = dst_layer.GetLayerDefn() feature = ogr.Feature(featureDefn) # loop over line_list - for n, xylist in enumerate(xylists): + for n, xylist in enumerate(line_coordinates): # Creating a line geometry linegeom = ogr.Geometry(ogr.wkbLineString) @@ -161,8 +176,8 @@ def writting_polyline(xylists, outfile, inlayer='', epsgno=0): feature.SetGeometry(linegeom) feature.SetField("id", n+1) - - layer.CreateFeature(feature) + # create feature + dst_layer.CreateFeature(feature) feature = None # Save and close DataSource @@ -305,7 +320,7 @@ def offset_random_point(lon, lat, offset=10, unit='km'): if offset > 0: if unit == 'km': # 1 degree = 111 km. - a0 = offset*(1/111.0) + a0 = offset*(1/DEG_TO_KM) else: raise Exception("convert offset to km.") @@ -442,3 +457,198 @@ def get_shapefile_epsg_code(in_shp): shp_lyr = shp_ds.GetLayer() src_epsg_code = int(shp_lyr.GetSpatialRef().GetAttrValue("AUTHORITY", 1)) return src_epsg_code + + +def select_by_location(select_feature_from=None, feature_overlap_with=None, overlap_method=None, save=False, out_selected_features=None, col_name=None): + if isinstance(select_feature_from, str): + get_selected = gpd.read_file(select_feature_from) + else: + get_selected = select_feature_from + + if isinstance(feature_overlap_with, str): + selected_by = gpd.read_file(feature_overlap_with) + else: + selected_by = feature_overlap_with + + selected_features = gpd.sjoin( + get_selected, selected_by, predicate=overlap_method) + + if isinstance(col_name, str): + selected_features.rename(columns={'index_right':col_name}, inplace=True) + + res = None + if save: + selected_features.to_file(out_selected_features) + else: + res = selected_features + + return res + + +def split_lines_at_points(line_feature, points, snap_dist=None, src_epsg_code=None, save=False, out_line_feature=None): + if src_epsg_code == 4326: + tolerance = snap_dist/DEG_TO_KM/1000 + else: + tolerance = snap_dist + # create point dataframe + df = pd.DataFrame( + { + "Latitude": [y for _, y in points], + "Longitude": [x for x, _ in points], + } + ) + # create points geodataframe + pnt_gdf = gpd.GeoDataFrame(df, geometry=gpd.points_from_xy( + df.Longitude, df.Latitude), crs=f"EPSG:{src_epsg_code}") + + if isinstance(line_feature, str): + ln_gdf = gpd.read_file(line_feature) + else: + ln_gdf = line_feature + + # get geometry + point_gdf = pnt_gdf.geometry.unary_union + line_gdf = ln_gdf.geometry.unary_union + + # split lines + split_line = split(snap(line_gdf, point_gdf, tolerance), point_gdf) + + # convert geometry collection to geo-dataframe + segments = [feat for feat in split_line.geoms] + gdf_segments = gpd.GeoDataFrame( + list(range(len(segments))), geometry=segments, crs=f"EPSG:{src_epsg_code}") + gdf_segments.columns = ['objectid', 'geometry'] + + res = None + if save: + gdf_segments.to_file(out_line_feature) + else: + res = gdf_segments + + return res + + +def leap_year(y): + if (y % 4.0 == 0) and (y % 100.0 == 0): + if (y % 400.0 == 0): + result = "Leap year", y + else: + result = "Not a leap year", y + elif (y % 4.0 == 0) and (y % 100.0 != 0): + result = "Leap year", y + else: + result = "Not a leap year", y + return result + + +def unlink_files(dir_path=None, file_extent=None, is_file=False, file_path=None): + # search files + if is_file: + flist = glob.glob( + os.path.join( + os.path.dirname(file_path), + os.path.basename(file_path).split('.')[0]+'*' + ) + ) + else: + flist = glob.glob( + os.path.join( + dir_path, + '*.'+file_extent + ) + ) + # remove files + for f in flist: + os.remove(f) + return + + +def write_point_file(points, out_point_shp, dst_epsg_cide): + schema = {'geometry': 'Point', 'properties': {}} + with collection( + out_point_shp, + "w", "ESRI Shapefile", + schema, + crs=from_epsg(dst_epsg_cide) + ) as output: + for x, y in points: + output.write({'properties': {}, 'geometry': mapping(Point(x, y))}) + + +def extend_line_shapefile(line_shapefile, extended_shapefile, max_snap_dist=None, keep_field=None, flip=False): + # compute new coordinates + def new_node(nodes, dist): + x1, y1 = nodes[0] + x2, y2 = nodes[1] + dx = x2 - x1 + dy = y2 - y1 + line_ln = math.hypot(dx, dy) + new_x = x2 + dx/line_ln * dist + new_y = y2 + dy/line_ln * dist + return new_x, new_y + + def accumulate(i): + it = iter(i) + total = next(it) + yield total + for elm in it: + total = operator.add(total, elm) + yield total + return + + def extract_line_geoms(line_gdf): + nodes = [] + line_ids = [] + for i in line_gdf.itertuples(): + objid = i.Index + geom = i.geometry + coords = list(geom.coords) + nodes.append(coords) + line_ids.append(objid) + return line_ids, nodes + + # read line geometry + line2 = gpd.read_file(line_shapefile) + epsg_code = get_shapefile_epsg_code(line_shapefile) + seq_ids = line2[keep_field].tolist() + + # get line ids and line nodes + if flip: + shp_ds = ogr.Open(line_shapefile) + shp_lyr = shp_ds.GetLayer() + line_ids, line_nodes = flip_line(shp_lyr) + shp_lyr = None + shp_ds = None + else: + line_ids, line_nodes = extract_line_geoms(line2) + + if epsg_code == 4326: + max_snap_dist = max_snap_dist*(1/DEG_TO_KM)/1000 + + # count vectices + vert_counts = list(accumulate(collections.Counter(line_ids).values())) + # get the end node of the lines + end_pnts = [pnt for i, pnt in enumerate(line_nodes) if i+1 in vert_counts or i+2 in vert_counts] + # calculate new nodes + new_nodes = [new_node(geom, max_snap_dist) for geom in end_pnts] + # extended line nodes + ext_nodes = [p_nds + [n_nd] for p_nds, n_nd in zip(line_nodes, new_nodes)] + + # construct line feature + line_strings = [] + line_index = [] + for line_coords, line_id in zip(ext_nodes, seq_ids): + line = LineString(line_coords) + line_strings.append(line) + line_index.append(line_id) + + # save line feature + col_dict = {keep_field: [ix for ix in line_index]} + gdf = gpd.GeoDataFrame( + col_dict, + index=range(len(line_index)), + crs='epsg:'+str(epsg_code), + geometry=line_strings + ) + gdf.to_file(extended_shapefile) + return diff --git a/rsgtools/vector_ops/create_buffer.py b/rsgtools/vector_ops/create_buffer.py new file mode 100644 index 0000000..90a271b --- /dev/null +++ b/rsgtools/vector_ops/create_buffer.py @@ -0,0 +1,52 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 19:29:41 2024 + +@author: USER +""" + +import geopandas as gpd +import pandas as pd +import shapely.ops as sp_ops + +from pyproj import Transformer +from rsgtools.ref_scripts import find_wgs2utm_epsg_code + +def buffer_feature(input_feature, buffer_offset=None, src_epsg_code=None, save=False, out_buffer_feature=None): + if isinstance(input_feature, str): + gdf = gpd.read_file(input_feature) + else: + gdf = input_feature + + # check src projection + if src_epsg_code == 4326: + minx, miny, maxx, maxy = gdf.geometry.total_bounds + dst_epsg_code = int(find_wgs2utm_epsg_code(minx, maxy)) + transformer1 = Transformer.from_crs(f'EPSG:{src_epsg_code}', f'EPSG:{dst_epsg_code}', always_xy=True) + # set geometry + geom_transformed1 = pd.Series([sp_ops.transform(transformer1.transform, geom) for geom in gdf['geometry'].tolist()]) + gdf.set_geometry(geom_transformed1, inplace=True, crs=dst_epsg_code) + # create buffer + gdf_buffer_series = gdf.buffer(buffer_offset) + gdf_buffer = gpd.GeoDataFrame(geometry=gpd.GeoSeries(gdf_buffer_series)) + # re-project to src coordinates + transformer2 = Transformer.from_crs(f'EPSG:{dst_epsg_code}', f'EPSG:{src_epsg_code}', always_xy=True) + # set geometry + geom_transformed2 = pd.Series([sp_ops.transform(transformer2.transform, geom) for geom in gdf_buffer['geometry'].tolist()]) + gdf_buffer.set_geometry(geom_transformed2, inplace=True, crs=src_epsg_code) + + else: + gdf_buffer = gdf.buffer(buffer_offset) + + res = gdf_buffer + if save: + gdf_buffer.to_file(out_buffer_feature) + else: + res = gdf_buffer + + return res + + + +# Reference +# https://github.com/geopandas/geopandas/issues/1175 diff --git a/rsgtools/vector_ops/find_nearest_point.py b/rsgtools/vector_ops/find_nearest_point.py index 3cc025a..106f72c 100644 --- a/rsgtools/vector_ops/find_nearest_point.py +++ b/rsgtools/vector_ops/find_nearest_point.py @@ -43,12 +43,25 @@ def find_closest_point(self): return nearest_pnt -def get_nearest_point(point: str, ref_line_file: Path): +def get_nearest_point(point: list, ref_line_file: Path): pntconfig = PointConfig( - single_point=point + inpoint=None, + single_point=point, + point_offset=None ) plconfig = PolylineConfig( - inline=ref_line_file + inline=ref_line_file, + fnodeid=None, + tnodeid=None, + lineid=None, + seqid=None, + is_xscl_available=None, + outxsclline=None, + uidfield=None, + objectidfield=None, + linefield=None, + interval=None, + offset=None ) obj = FindNearestPoint(pntconfig=pntconfig, plineconfig=plconfig) res = obj.find_closest_point() diff --git a/rsgtools/vector_ops/lines_middle_points.py b/rsgtools/vector_ops/lines_middle_points.py new file mode 100644 index 0000000..ae58c4b --- /dev/null +++ b/rsgtools/vector_ops/lines_middle_points.py @@ -0,0 +1,77 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Mar 25 21:29:13 2024 + +@author: Debabrata Ghorai, Ph.D. + +Generate middle point of a line. + +""" + +import os +import glob + +from osgeo import ogr +from rsgtools.ref_scripts import interval_point +from rsgtools.utils import dist_calc + + +def generate_lines_middle_point(line_feature=None, out_mid_points=None): + src_ds = ogr.Open(line_feature) + layer = src_ds.GetLayer(0) + # Get input file spatial reference + sr = layer.GetSpatialRef() + # sr.ExportToWkt() + # Create new shapefile + driver = ogr.GetDriverByName('ESRI Shapefile') + if os.path.exists(out_mid_points): + try: + driver.DeleteDataSource(out_mid_points) + except: + pass + flist = glob.glob( + os.path.join( + os.path.dirname(out_mid_points), + os.path.basename(out_mid_points).split('.')[0]+'*' + ) + ) + for f in flist: + os.remove(f) + # create ds + ds = driver.CreateDataSource(out_mid_points) + # Create new shapefile layer + lyr = ds.CreateLayer('mylyr', sr, ogr.wkbPoint) + lyrDef = lyr.GetLayerDefn() + + # calculate mid point of a line + for i in range(layer.GetFeatureCount()): + feature = layer.GetFeature(i) + geom = feature.GetGeometryRef() + pnts = geom.GetPoints() + if not isinstance(pnts, type(None)): + # calc segment dist + lineDist = list() + for p in range(len(pnts)-1): + lineDist.append(dist_calc(pnts[p], pnts[p+1])) + # get mid-length + midLength = sum(lineDist)*0.5 + objects = [i+1] + nodes = [pnts] + points, _ = interval_point(objects, nodes, midLength) + # Point file generation + point = ogr.Geometry(ogr.wkbPoint) + point.SetPoint(0, points[0][0], points[0][1]) + featIndex = 0 + feat = ogr.Feature(lyrDef) + feat.SetGeometry(point) + feat.SetFID(featIndex) + lyr.CreateFeature(feat) + # Flush + layer = None + sr = None + ds.Destroy() + ds = None + driver = None + lyrDef = None + lyr = None + return diff --git a/setup.py b/setup.py index 00000bd..d358d71 100644 --- a/setup.py +++ b/setup.py @@ -14,8 +14,8 @@ def get_requirements(): setup( name='geoutils-rsg', version='0.0.3', - package_dir={"": "geobhumi"}, - packages=find_packages(where='geobhumi'), + package_dir={"": "rsgtools"}, + packages=find_packages(where='rsgtools'), url='https://github.com/dghorai', license='MIT License', author='Debabrata Ghorai, Ph.D.', diff --git a/tests/results/merge_lines.dbf b/tests/results/merge_lines.dbf deleted file mode 100644 index 4ff327d..0000000 Binary files a/tests/results/merge_lines.dbf and /dev/null differ diff --git a/tests/results/merge_lines.prj b/tests/results/merge_lines.prj deleted file mode 100644 index f45cbad..0000000 --- a/tests/results/merge_lines.prj +++ /dev/null @@ -1 +0,0 @@ -GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]] \ No newline at end of file diff --git a/tests/results/merge_lines.shp b/tests/results/merge_lines.shp deleted file mode 100644 index ca13aba..0000000 Binary files a/tests/results/merge_lines.shp and /dev/null differ diff --git a/tests/results/merge_lines.shx b/tests/results/merge_lines.shx deleted file mode 100644 index 5e05b07..0000000 Binary files a/tests/results/merge_lines.shx and /dev/null differ diff --git a/tests/results/out_cross_section.dbf b/tests/results/out_cross_section.dbf deleted file mode 100644 index f8eb8ae..0000000 Binary files a/tests/results/out_cross_section.dbf and /dev/null differ diff --git a/tests/results/out_cross_section.prj b/tests/results/out_cross_section.prj deleted file mode 100644 index f45cbad..0000000 --- a/tests/results/out_cross_section.prj +++ /dev/null @@ -1 +0,0 @@ -GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]] \ No newline at end of file diff --git a/tests/results/out_cross_section.shp b/tests/results/out_cross_section.shp deleted file mode 100644 index 65a2f6d..0000000 Binary files a/tests/results/out_cross_section.shp and /dev/null differ diff --git a/tests/results/out_cross_section.shx b/tests/results/out_cross_section.shx deleted file mode 100644 index 3d3dd67..0000000 Binary files a/tests/results/out_cross_section.shx and /dev/null differ diff --git a/tests/results/point_to_grid.dbf b/tests/results/point_to_grid.dbf deleted file mode 100644 index 58ff07e..0000000 Binary files a/tests/results/point_to_grid.dbf and /dev/null differ diff --git a/tests/results/point_to_grid.prj b/tests/results/point_to_grid.prj deleted file mode 100644 index e5b6987..0000000 --- a/tests/results/point_to_grid.prj +++ /dev/null @@ -1 +0,0 @@ -GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]] \ No newline at end of file diff --git a/tests/results/point_to_grid.shp b/tests/results/point_to_grid.shp deleted file mode 100644 index 410ca0f..0000000 Binary files a/tests/results/point_to_grid.shp and /dev/null differ diff --git a/tests/results/point_to_grid.shx b/tests/results/point_to_grid.shx deleted file mode 100644 index 65ccc71..0000000 Binary files a/tests/results/point_to_grid.shx and /dev/null differ diff --git a/tests/test_geobhumi.py b/tests/test_rsgtools.py similarity index 100% rename from tests/test_geobhumi.py rename to tests/test_rsgtools.py