Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Multiple bugfixes for CARQ correction of OFCL track #68

Merged
merged 10 commits into from
Mar 24, 2023
+21,873 −21,822
Merged

Multiple bugfixes for CARQ correction of OFCL track #68

Show file tree
Hide file tree
Changes from 5 commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
cb53e04
Fix multiple bugs
SorooshMani-NOAA Mar 15, 2023
a7a3825
Style fixes
SorooshMani-NOAA Mar 15, 2023
c17afdb
Test for carq fix of ofcl track
SorooshMani-NOAA Mar 15, 2023
c5751f7
Update test references after carq correction bugfix
SorooshMani-NOAA Mar 15, 2023
2ebcc72
Fix invalid test vaule after track order fix #67
SorooshMani-NOAA Mar 15, 2023
4097d56
corrections for adding forecast background pressure usin CARQ backgro…
WPringle Mar 23, 2023
7db6cb0
Fix more correction bugs & update tests
SorooshMani-NOAA Mar 23, 2023
aa94d10
Style fix
SorooshMani-NOAA Mar 23, 2023
3bc58c4
Bugfix for background and central pressure correction on write
SorooshMani-NOAA Mar 24, 2023
e6b777d
Style fix
SorooshMani-NOAA Mar 24, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
133 changes: 84 additions & 49 deletions stormevents/nhc/track.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -981,55 +981,9 @@ def unfiltered_data(self, dataframe: DataFrame):
# fill missing values of MRD and MSLP in the OFCL advisory
if "OFCL" in self.advisories:
tracks = separate_tracks(dataframe)

if "OFCL" in tracks:
ofcl_tracks = tracks["OFCL"]
carq_tracks = tracks["CARQ"]

for initial_time, forecast in ofcl_tracks.items():
if initial_time in carq_tracks:
carq_forecast = carq_tracks[initial_time]
else:
carq_forecast = carq_tracks[list(carq_tracks)[0]]

relation = HollandBRelation()
holland_b = relation.holland_b(
max_sustained_wind_speed=carq_forecast[
"max_sustained_wind_speed"
],
background_pressure=carq_forecast["background_pressure"],
central_pressure=carq_forecast["central_pressure"],
)

holland_b[holland_b == numpy.inf] = numpy.nan
holland_b = numpy.nanmean(holland_b)

mrd_missing = pandas.isna(forecast["radius_of_maximum_winds"])
mslp_missing = pandas.isna(forecast["central_pressure"])
radp_missing = pandas.isna(forecast["background_pressure"])
Comment on lines -1007 to -1009
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Original values for these variables from the ATCF file were not empty, but they were 0 instead. So, the correction logic would ignore them and didn't fix the incorrect 0 value. #64


# fill OFCL maximum wind radius with the first entry from the CARQ advisory
forecast.loc[
mrd_missing, "radius_of_maximum_winds"
] = carq_forecast["radius_of_maximum_winds"].iloc[0]
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The first CARQ track could be a -24 hour "forecast", which has missing values for variables of interest. iloc[0] just gets the first entry of carq_forecast dataframe instead of getting the first hour-0 entry


# fill OFCL background pressure with the first entry from the CARQ advisory central pressure (at sea level)
forecast.loc[radp_missing, "background_pressure"] = carq_forecast[
"central_pressure"
].iloc[0]

# fill OFCL central pressure (at sea level) with the 3rd hour entry, preserving Holland B
forecast.loc[
mslp_missing, "central_pressure"
] = relation.central_pressure(
max_sustained_wind_speed=forecast.loc[
mslp_missing, "max_sustained_wind_speed"
],
background_pressure=forecast.loc[
mslp_missing, "background_pressure"
],
holland_b=holland_b,
)
Comment on lines -1022 to -1032
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Prior to the fixes in this pull, this update forecast data frame was not stored anywhere and the corrections were discarded. #66

if all(adv in tracks for adv in ["OFCL", "CARQ"]):
tracks = correct_ofcl_based_on_carq_n_hollandb(tracks)
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The logic for CARQ correction is now moved to a standalone function

dataframe = combine_tracks(tracks)
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The corrected tracks are now converted back into a dataframe which replaces the original one #66


if len(self.__advisories_to_remove) > 0:
dataframe = dataframe[
Expand Down Expand Up @@ -1175,6 +1129,7 @@ def separate_tracks(data: DataFrame) -> Dict[str, Dict[str, DataFrame]]:
tracks = {}
for advisory in pandas.unique(data["advisory"]):
advisory_data = data[data["advisory"] == advisory]
advisory_data["forecast_hours"] = advisory_data.forecast_hours.astype(int)
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The forecast_hours data is used for sorting the dataframe. Because the data was originally of string type, the sort used to happen at character level (i.e. first 1s and then 2s, etc. so that 120 would come before 24). Now this issue is fixed by using int type #67


if advisory == "BEST":
advisory_data = advisory_data.sort_values("datetime")
Expand All @@ -1196,3 +1151,83 @@ def separate_tracks(data: DataFrame) -> Dict[str, Dict[str, DataFrame]]:
] = track_data

return tracks


def combine_tracks(tracks: Dict[str, Dict[str, DataFrame]]) -> DataFrame:
"""
combine tracks separated using `separate_tracks`

:param tracks: dictionary of forecasts for each advisory (aside from best track ``BEST``, which only has one hindcast)
:return: data frame of track
"""

return pandas.concat([df for adv_trk in tracks.values() for df in adv_trk.values()])


def correct_ofcl_based_on_carq_n_hollandb(
tracks: Dict[str, Dict[str, DataFrame]]
) -> Dict[str, Dict[str, DataFrame]]:
"""
Correct official forecast using consensus track along with holland-b
relation

:param tracks: dictionary of forecasts for each advisory (aside from best track ``BEST``, which only has one hindcast)
:return: dictionary of forecasts for each advisory (aside from best track ``BEST``, which only has one hindcast) with corrected OFCL
"""

ofcl_tracks = tracks["OFCL"]
carq_tracks = tracks["CARQ"]

corr_ofcl_tracks = dict()

for initial_time, forecast in ofcl_tracks.items():
if initial_time in carq_tracks:
carq_forecast = carq_tracks[initial_time]
else:
carq_forecast = carq_tracks[list(carq_tracks)[0]]

relation = HollandBRelation()
holland_b = relation.holland_b(
max_sustained_wind_speed=carq_forecast["max_sustained_wind_speed"],
background_pressure=carq_forecast["background_pressure"],
central_pressure=carq_forecast["central_pressure"],
)

holland_b[holland_b == numpy.inf] = numpy.nan
holland_b = numpy.nanmean(holland_b)

# Get CARQ from forecast hour 0 and isotach 34kt (i.e. the first item)
carq_ref = carq_forecast.loc[carq_forecast.forecast_hours == 0].iloc[0]
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

CARQ track entry to be used is now the first entry of hour-0 forecast (i.e. 34 kt) #65


columns_of_interest = forecast[
["radius_of_maximum_winds", "central_pressure", "background_pressure"]
]
columns_of_interest[columns_of_interest == 0] = pandas.NA
missing = columns_of_interest.isna()
# Order of columns is the same as columns_of_interest
mrd_missing = missing.iloc[:, 0]
mslp_missing = missing.iloc[:, 1]
radp_missing = missing.iloc[:, 2]
Comment on lines +1204 to +1209
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The values are now checked for both NA as well as 0 to indicate the missing data. #64


# fill OFCL maximum wind radius with the first entry from the CARQ advisory
forecast.loc[mrd_missing, "radius_of_maximum_winds"] = carq_ref[
"radius_of_maximum_winds"
]

# fill OFCL background pressure with the first entry from the CARQ advisory central pressure (at sea level)
forecast.loc[radp_missing, "background_pressure"] = carq_ref["central_pressure"]

# fill OFCL central pressure (at sea level) with the 3rd hour entry, preserving Holland B
forecast.loc[mslp_missing, "central_pressure"] = relation.central_pressure(
max_sustained_wind_speed=forecast.loc[
mslp_missing, "max_sustained_wind_speed"
],
background_pressure=forecast.loc[mslp_missing, "background_pressure"],
holland_b=holland_b,
)

corr_ofcl_tracks[initial_time] = forecast
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Now the correction is stored in the corr_ofcl_tracks which replaces the values in the original forecast track. #66


tracks["OFCL"] = corr_ofcl_tracks

return tracks
Loading