Feature/rmax psurgev29 (#154)

* carry over subsetted elements in adcirc subset parsing

* adding in the Psurgev2.9 RMW errors as a Gaussian distribution same as the other variables and renaming legacy method for RMW to RadiustoMaximumWindsPersistent

* reformatting and changing perturbationtype of isotachs to None

* correcting index of forecast times for new RMW/Rmax errors

* name for RMW persistent has to have persisent removed to correspond to right column of the data

* transfer variable instead of variable_names to write so that can distinguish between the two RadiusOfMaximumWinds

* corrections to allow the paralle processing to work as well as delineate between the two RadiusOfMaximumWind options

* adding new references for best_track_ensemble test with gaussian RMW errors

* making sure test_existing_advisory test works with persistent RMW option as before

* adding test for the online OFCL advisory using new RMW forecasts in StormEvents

* reformat

* adding dependency on stormevents >= v2.3 due to RMW regression forecast

* first commit for new implementation that changes isotach radii based on changed Vmax, Rmax, and track using the GAHM profile

* working method for full GAHM profile, needs checking for non-convergence

* using bi-section method guess for GAHM root-finding algo to avoid non-convergence

* cleaning up with function for GAHM root finding algo

* finishing cleanup by moving the isotach_rad finding algo call into the GAHM profile function

* update to test reference original track for Rmax bias-correction and fix to not exceed the maximum isotach

* ensuring order of Rmax is last in the perturb loop

* changes to make GAHM functions neater

* improve accuracy of function return outputs

* ensuring correct units for the translation_speed and preserving the initial speed

* removing surplus isotach lines if they all go to zero (for 50-kt and 64-kt only)

* updating reference files for the online advisory test with the GAHM isotach perturbation method

* updating besttrack perturbation test reference files for the GAHM based isotach perturbation upgrade

* updating existing advisory test reference files for the GAHM based isotach perturbation upgrade

* changing translation speed impact on Vmax to use LC12 0.66 factor

* oitnb format

* cleaning of radial Vmax computation by moving to a function of the MaximumSustainedWindSpeed class and including a lower_bound to prevent low and negative values

* cleaning up radial Vmax computation in the isotach class and adding translation_speed output from the radial Vmax function to be used in Vr calculation

* moving central pressure computation and isotach radius computation outside of loop to make it independent of loop order

* updates to the GAHM isotach perturbator to make sure original track data can be matched, especially relevant for the Rmax over isotach_rad ratio > 1

* adding a radius dependent factor for translation_speed based on LC12 figure 2

* adding check for small or zero Vr_test which may not go negative

* new advisory reference test files with latest changes

* new best_track reference test files with latest changes

* ensuring no changes to the original file by only changing pressure and isotach if perturbations are made

* changes for original test files given preservation

* differentiate between a zero perturbation and the original file by using nan for the original

* changes to json for original file (now NaN instead of 0)

* correction to call to the ts_factor, and adjustment of B where it goes nan but should be a solution

* corrections to allow isotachs that are originally zero to be filled in for the perturbed storm in the case it can fit with the GAHM profile

* Correct pyproject.toml

* Update test refs

* Fix xarray compatiblitiy issue

* Pin scipy to at or below 1.13

* Revert back to scipy 1.12

---------

Co-authored-by: William Pringle <[email protected]>

ensembleperturbation/parsing/adcirc.py

@@ -758,6 +758,20 @@ def combine_outputs(
                 if elevation_selection is not None:
                     elevation_subset = subset
 
+                if 'element' in file_data:
+                    # keep only elements where all nodes are present
+                    elements = file_data['element'].values
+                    element_mask = numpy.isin(elements, file_data['node'].values).all(axis=1)
+                    elements = elements[element_mask]
+                    # map nodes in element table to local numbering system (start at 0)
+                    node_mapper = numpy.zeros(file_data['node'].max().values + 1, dtype=int)
+                    node_index = numpy.arange(len(file_data['node']))
+                    node_mapper[file_data['node'].values] = node_index
+                    elements = node_mapper[elements]
+                    # update element table in dataset
+                    ele_da = DataArray(data=elements, dims=['nele', 'nvertex'])
+                    file_data = file_data.drop('element').assign_coords({'element': ele_da})
+
             output_data[basename] = file_data
 
     for basename, file_data in output_data.items():