Re-assess/update PySCHISM

[janahaddad]

No description provided.

[mansurjisan]

TypeError in fgrid.py: 'sjoin() got an unexpected keyword argument 'op''

Description

I encountered a TypeError in fgrid.py when trying to use the modify_by_region method. The error message indicates that sjoin() got an unexpected keyword argument 'op'. This issue seems to be caused by a change in the geopandas API, where the op parameter has been replaced by predicate.

Steps to Reproduce

1. Run the following code:

from pyschism.mesh.hgrid import Hgrid
from pyschism.mesh.fgrid import DragCoefficient

if __name__ == '__main__':

    hgrid=Hgrid.open('hgrid.gr3', crs='epsg:4326')

    depth1=-1.0
    depth2=-3.0
    bfric_river=0.0025
    bfric_land=0.025

    fgrid=DragCoefficient.linear_with_depth(hgrid, depth1, depth2, bfric_river, bfric_land)

    regions=['region1.reg', 'region2.reg']
    values=[0.001, 0.0]
    flags=[1, 0] 

    for reg, value, flag in zip(regions, values, flags):
        fgrid.modify_by_region(hgrid, f'./{reg}', value, depth1, flag)

    fgrid.write('drag.gr3', overwrite=True)

2. The error occurs on this line in fgrid.py:

gdf_in = gpd.sjoin(gdf2, gdf1, op="within")

Proposed Fix

The issue can be fixed by replacing op with predicate:

gdf_in = gpd.sjoin(gdf2, gdf1, predicate="within")

Additional Information

PySCHISM version: 0.1.15
Geopandas version: 1.0.1
Python version: 3.9.19

[mansurjisan]

I've developed a new plotting feature for pyschism that allows quick visualization of input files used in the SCHISM model. The goal here is to improve our workflow by enabling us to visually inspect our input data before running the model, potentially saving time and helping to catch issues early.

Here's what I've done:

1. Objective:

   • Develop a quick and easy way to visualize SCHISM model input files before running the model.

2. Implementation:
   (a) Created a fork of pyschism in my GitHub repo: https://github.com/mansurjisan/pyschism
   (b) Added a new 'plotting' directory inside the pyschism directory.
   (c) Created plot_gr3.py file in the plotting directory, which contains all plot-related functions, including:
   - read_gr3_files: for reading SCHISM grid files
   - Functions for downloading necessary shapefiles for coastlines and boundaries
   - plot_gr3: the main plotting function

3. Key Features:

   • Visualizes SCHISM grid elements and depth information
   • Options to show or hide land features
   • Customizable plot extent
   • Handles potential issues with input data (NaN, Inf values)
   • Uses Cartopy for map projections and coastline data

4. How to Use:
   The main function for plotting is plot_gr3. Here's an example of how to use it:

   from pyschism import read_gr3_file, plot_gr3
   
   filename = 'hgrid.gr3'
   nodes, elements = read_gr3_file(filename)
   plot_gr3(nodes, elements, "output_plot.png", min_lon=-88, max_lon=-77, min_lat=20, max_lat=30, show_land=False)

   Parameters:

   • nodes, elements: Data read from the gr3 file
   • output_file: Name of the output plot file
   • min_lon, max_lon, min_lat, max_lat: Plot extent (optional)
   • show_land: Boolean to show/hide land features (default: True)

5. New Files in Directory Structure:
   I've added new files to the pyschism directory. Here's a highlight of the relevant new additions in the directory tree:

   pyschism/
   └── plotting/
   ├── init.py
   ├── plot_gr3.py # Main plotting function for grid files
   └── plot_albedo_gr3.py # New file for plotting albedo data (in development)

   test/
   ├── plot_schism_grid.py # Example script for grid plotting
   └── plot_schism_albedo.py # Example script for albedo plotting (in development)

6. Example Script:
   I've created a sample script (plot_schism_grid.py) demonstrating how to use this new feature:

   from pyschism import read_gr3_file, plot_gr3
   
   if __name__ == "__main__":
       filename = 'hgrid.gr3'
       output_file = 'schism_grid_plot.png'
       print(f"Reading file: {filename}")
       nodes, elements = read_gr3_file(filename)
       if nodes is not None and elements is not None:
           plot_gr3(nodes, elements, "output_plot.png", min_lon=-88, max_lon=-77, min_lat=20, max_lat=30, show_land=False)
       else:
           print("Failed to read the .gr3 file. Please check the file path and format.")

Next Steps:

• Continue development of other input files visualization capability
• Implement visualization for various forcing files
• Enhance flexibility to handle different input data types
• Improve error handling

Please feel free to check out the fork, try the new feature, and let me know if you have any questions, suggestions, or if you encounter any issues.