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.

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
160 changes: 96 additions & 64 deletions stormevents/nhc/track.py
Original file line number Diff line number Diff line change
Expand Up @@ -564,7 +564,6 @@ def atcf(self, advisory: ATCF_Advisory = None) -> DataFrame:
atcf["max_sustained_wind_speed"] = (
atcf["max_sustained_wind_speed"].astype("string").str.pad(5)
)
atcf["central_pressure"] = atcf["central_pressure"].astype("string").str.pad(5)
atcf["development_level"] = atcf["development_level"].str.pad(3)
atcf["isotach_radius"] = atcf["isotach_radius"].astype("string").str.pad(4)
atcf["isotach_quadrant_code"] = atcf["isotach_quadrant_code"].str.pad(4)
Expand All @@ -582,22 +581,20 @@ def atcf(self, advisory: ATCF_Advisory = None) -> DataFrame:
)

atcf["background_pressure"].fillna(method="ffill", inplace=True)
atcf.loc[
~pandas.isna(self.data["central_pressure"])
& (self.data["background_pressure"] <= self.data["central_pressure"])
& (self.data["central_pressure"] < 1013),
"background_pressure",
] = "1013"
atcf.loc[
~pandas.isna(self.data["central_pressure"])
& (self.data["background_pressure"] <= self.data["central_pressure"])
& (self.data["central_pressure"] < 1013),
"background_pressure",
] = (
self.data["central_pressure"] + 1
atcf["background_pressure"] = atcf["background_pressure"].astype(int)
atcf["central_pressure"] = atcf["central_pressure"].astype(int)
press_cond = (
~pandas.isna(atcf["central_pressure"])
& (atcf["background_pressure"] <= atcf["central_pressure"])
& (atcf["central_pressure"] < 1013)
)
atcf.loc[press_cond, "background_pressure"] = 1013
atcf.loc[press_cond, "background_pressure"] = (
atcf.loc[press_cond, "central_pressure"] + 1
WPringle marked this conversation as resolved.
Show resolved Hide resolved
)
atcf["central_pressure"] = atcf["central_pressure"].astype("string").str.pad(5)
atcf["background_pressure"] = (
atcf["background_pressure"].astype(int).astype("string").str.pad(5)
atcf["background_pressure"].astype("string").str.pad(5)
)

atcf["radius_of_last_closed_isobar"] = (
Expand Down Expand Up @@ -981,55 +978,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 +1126,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 +1148,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 0-hr CARQ
forecast.loc[mrd_missing, "radius_of_maximum_winds"] = carq_ref[
"radius_of_maximum_winds"
]

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

# fill OFCL central pressure (at sea level), preserving Holland B from 0-hr CARQ
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
38 changes: 19 additions & 19 deletions tests/data/reference/test_vortex_track/ike2008.fort.22
Original file line number Diff line number Diff line change
Expand Up @@ -61,18 +61,18 @@ AL, 09, 2008090706, , BEST, 0, 211N,716W, 115, 947, HU, 64, NEQ, 40,
AL, 09, 2008090712, , BEST, 0, 210N,728W, 110, 947, HU, 34, NEQ, 140, 120, 100, 125, 1012, 200, 15, 140, 0, L, 0, ,265, 6, IKE, 26
AL, 09, 2008090712, , BEST, 0, 210N,728W, 110, 947, HU, 50, NEQ, 90, 75, 50, 75, 1012, 200, 15, 140, 0, L, 0, ,265, 6, IKE, 26
AL, 09, 2008090712, , BEST, 0, 210N,728W, 110, 947, HU, 64, NEQ, 50, 40, 30, 50, 1012, 200, 15, 140, 0, L, 0, ,265, 6, IKE, 26
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 34, NEQ, 140, 120, 100, 125,,,, ,,,,,270, 12,, 27
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 50, NEQ, 90, 75, 50, 75,,,, ,,,,,270, 12,, 27
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 64, NEQ, 50, 40, 30, 50,,,, ,,,,,270, 12,, 27
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 34, NEQ, 140, 120, 100, 125, 948,,, ,,,,,270, 12,, 27
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 50, NEQ, 90, 75, 50, 75, 948,,, ,,,,,270, 12,, 27
AL, 09, 2008090713, , BEST, 0, 210N,732W, 110, 947, HU, 64, NEQ, 50, 40, 30, 50, 948,,, ,,,,,270, 12,, 27
AL, 09, 2008090718, , BEST, 0, 210N,740W, 105, 946, HU, 34, NEQ, 145, 125, 100, 125, 1012, 200, 15, 130, 0, L, 0, ,270, 5, IKE, 28
AL, 09, 2008090718, , BEST, 0, 210N,740W, 105, 946, HU, 50, NEQ, 90, 90, 50, 75, 1012, 200, 15, 130, 0, L, 0, ,270, 5, IKE, 28
AL, 09, 2008090718, , BEST, 0, 210N,740W, 105, 946, HU, 64, NEQ, 50, 40, 30, 50, 1012, 200, 15, 130, 0, L, 0, ,270, 5, IKE, 28
AL, 09, 2008090800, , BEST, 0, 211N,752W, 115, 945, HU, 34, NEQ, 150, 140, 100, 125, 1009, 200, 15, 130, 0, L, 0, ,275, 6, IKE, 29
AL, 09, 2008090800, , BEST, 0, 211N,752W, 115, 945, HU, 50, NEQ, 90, 105, 50, 75, 1009, 200, 15, 130, 0, L, 0, ,275, 6, IKE, 29
AL, 09, 2008090800, , BEST, 0, 211N,752W, 115, 945, HU, 64, NEQ, 50, 40, 30, 50, 1009, 200, 15, 130, 0, L, 0, ,275, 6, IKE, 29
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 34, NEQ, 150, 140, 100, 125,,,, ,,,,,270, 7,, 30
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 50, NEQ, 90, 105, 50, 75,,,, ,,,,,270, 7,, 30
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 64, NEQ, 50, 40, 30, 50,,,, ,,,,,270, 7,, 30
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 34, NEQ, 150, 140, 100, 125, 946,,, ,,,,,270, 7,, 30
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 50, NEQ, 90, 105, 50, 75, 946,,, ,,,,,270, 7,, 30
AL, 09, 2008090802, , BEST, 0, 211N,757W, 115, 945, HU, 64, NEQ, 50, 40, 30, 50, 946,,, ,,,,,270, 7,, 30
AL, 09, 2008090806, , BEST, 0, 211N,765W, 100, 950, HU, 34, NEQ, 155, 175, 100, 125, 1010, 200, 15, 120, 0, L, 0, ,270, 6, IKE, 31
AL, 09, 2008090806, , BEST, 0, 211N,765W, 100, 950, HU, 50, NEQ, 90, 120, 50, 75, 1010, 200, 15, 120, 0, L, 0, ,270, 6, IKE, 31
AL, 09, 2008090806, , BEST, 0, 211N,765W, 100, 950, HU, 64, NEQ, 50, 40, 30, 50, 1010, 200, 15, 120, 0, L, 0, ,270, 6, IKE, 31
Expand All @@ -91,9 +91,9 @@ AL, 09, 2008090906, , BEST, 0, 220N,814W, 70, 965, HU, 64, NEQ, 30,
AL, 09, 2008090912, , BEST, 0, 224N,824W, 70, 965, HU, 34, NEQ, 155, 155, 100, 170, 1008, 275, 15, 85, 0, L, 0, ,293, 5, IKE, 36
AL, 09, 2008090912, , BEST, 0, 224N,824W, 70, 965, HU, 50, NEQ, 90, 75, 30, 90, 1008, 275, 15, 85, 0, L, 0, ,293, 5, IKE, 36
AL, 09, 2008090912, , BEST, 0, 224N,824W, 70, 965, HU, 64, NEQ, 30, 0, 0, 30, 1008, 275, 15, 85, 0, L, 0, ,293, 5, IKE, 36
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 34, NEQ, 155, 155, 100, 170,,,, ,,,,,293, 8,, 37
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 50, NEQ, 90, 75, 30, 90,,,, ,,,,,293, 8,, 37
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 64, NEQ, 30, 0, 0, 30,,,, ,,,,,293, 8,, 37
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 34, NEQ, 155, 155, 100, 170, 966,,, ,,,,,293, 8,, 37
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 50, NEQ, 90, 75, 30, 90, 966,,, ,,,,,293, 8,, 37
AL, 09, 2008090914, , BEST, 0, 226N,829W, 70, 965, HU, 64, NEQ, 30, 0, 0, 30, 966,,, ,,,,,293, 8,, 37
AL, 09, 2008090918, , BEST, 0, 227N,833W, 65, 966, HU, 34, NEQ, 155, 150, 105, 170, 1008, 275, 15, 80, 0, L, 0, ,285, 3, IKE, 38
AL, 09, 2008090918, , BEST, 0, 227N,833W, 65, 966, HU, 50, NEQ, 100, 75, 30, 90, 1008, 275, 15, 80, 0, L, 0, ,285, 3, IKE, 38
AL, 09, 2008090918, , BEST, 0, 227N,833W, 65, 966, HU, 64, NEQ, 20, 0, 0, 30, 1008, 275, 15, 80, 0, L, 0, ,285, 3, IKE, 38
Expand Down Expand Up @@ -139,20 +139,20 @@ AL, 09, 2008091300, , BEST, 0, 283N,940W, 95, 952, HU, 64, NEQ, 110,
AL, 09, 2008091306, , BEST, 0, 291N,946W, 95, 951, HU, 34, NEQ, 225, 200, 125, 125, 1007, 325, 30, 115, 0, L, 0, ,327, 5, IKE, 52
AL, 09, 2008091306, , BEST, 0, 291N,946W, 95, 951, HU, 50, NEQ, 150, 160, 80, 75, 1007, 325, 30, 115, 0, L, 0, ,327, 5, IKE, 52
AL, 09, 2008091306, , BEST, 0, 291N,946W, 95, 951, HU, 64, NEQ, 110, 90, 55, 45, 1007, 325, 30, 115, 0, L, 0, ,327, 5, IKE, 52
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 34, NEQ, 225, 200, 125, 125,,,, ,,,,,336, 7,, 53
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 50, NEQ, 150, 160, 80, 75,,,, ,,,,,336, 7,, 53
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 64, NEQ, 110, 90, 55, 45,,,, ,,,,,336, 7,, 53
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 34, NEQ, 225, 200, 125, 125, 951,,, ,,,,,336, 7,, 53
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 50, NEQ, 150, 160, 80, 75, 951,,, ,,,,,336, 7,, 53
AL, 09, 2008091307, , BEST, 0, 293N,947W, 95, 950, HU, 64, NEQ, 110, 90, 55, 45, 951,,, ,,,,,336, 7,, 53
AL, 09, 2008091312, , BEST, 0, 303N,952W, 85, 959, HU, 34, NEQ, 125, 180, 125, 60, 1007, 325, 30, 105, 0, L, 0, ,337, 7, IKE, 54
AL, 09, 2008091312, , BEST, 0, 303N,952W, 85, 959, HU, 50, NEQ, 75, 90, 60, 45, 1007, 325, 30, 105, 0, L, 0, ,337, 7, IKE, 54
AL, 09, 2008091312, , BEST, 0, 303N,952W, 85, 959, HU, 64, NEQ, 50, 45, 30, 20, 1007, 325, 30, 105, 0, L, 0, ,337, 7, IKE, 54
AL, 09, 2008091318, , BEST, 0, 317N,953W, 50, 974, TS, 34, NEQ, 75, 150, 60, 40, 1007, 325, 40, 65, 0, L, 0, ,357, 7, IKE, 55
AL, 09, 2008091318, , BEST, 0, 317N,953W, 50, 974, TS, 50, NEQ, 40, 45, 0, 0, 1007, 325, 40, 65, 0, L, 0, ,357, 7, IKE, 55
AL, 09, 2008091400, , BEST, 0, 335N,949W, 35, 980, TS, 34, NEQ, 60, 90, 0, 0, 1007, 300, 50, 45, 0, L, 0, , 11, 9, IKE, 56
AL, 09, 2008091406, , BEST, 0, 355N,937W, 35, 985, TS, 34, NEQ, 45, 90, 20, 0, 1007, 300, 50, 40, 0, L, 0, , 26, 11, IKE, 57
AL, 09, 2008091412, , BEST, 0, 376N,910W, 40, 987, EX, 34, NEQ, 0, 90, 40, 0,,,, ,,,,, 45, 16,, 58
AL, 09, 2008091418, , BEST, 0, 403N,872W, 50, 988, EX, 34, NEQ, 0, 180, 150, 0,,,, ,,,,, 47, 21,, 59
AL, 09, 2008091418, , BEST, 0, 403N,872W, 50, 988, EX, 50, NEQ, 0, 160, 0, 0,,,, ,,,,, 47, 21,, 59
AL, 09, 2008091500, , BEST, 0, 433N,815W, 50, 988, EX, 34, NEQ, 0, 180, 150, 0,,,, ,,,,, 53, 27,, 60
AL, 09, 2008091500, , BEST, 0, 433N,815W, 50, 988, EX, 50, NEQ, 0, 160, 0, 0,,,, ,,,,, 53, 27,, 60
AL, 09, 2008091506, , BEST, 0, 458N,753W, 40, 986, EX, 34, NEQ, 0, 180, 150, 0,,,, ,,,,, 58, 26,, 61
AL, 09, 2008091512, , BEST, 0, 472N,711W, 35, 986, EX, 34, NEQ, 0, 180, 150, 0,,,, ,,,,, 63, 17,, 62
AL, 09, 2008091412, , BEST, 0, 376N,910W, 40, 987, EX, 34, NEQ, 0, 90, 40, 0, 988,,, ,,,,, 45, 16,, 58
AL, 09, 2008091418, , BEST, 0, 403N,872W, 50, 988, EX, 34, NEQ, 0, 180, 150, 0, 989,,, ,,,,, 47, 21,, 59
AL, 09, 2008091418, , BEST, 0, 403N,872W, 50, 988, EX, 50, NEQ, 0, 160, 0, 0, 989,,, ,,,,, 47, 21,, 59
AL, 09, 2008091500, , BEST, 0, 433N,815W, 50, 988, EX, 34, NEQ, 0, 180, 150, 0, 989,,, ,,,,, 53, 27,, 60
AL, 09, 2008091500, , BEST, 0, 433N,815W, 50, 988, EX, 50, NEQ, 0, 160, 0, 0, 989,,, ,,,,, 53, 27,, 60
AL, 09, 2008091506, , BEST, 0, 458N,753W, 40, 986, EX, 34, NEQ, 0, 180, 150, 0, 987,,, ,,,,, 58, 26,, 61
AL, 09, 2008091512, , BEST, 0, 472N,711W, 35, 986, EX, 34, NEQ, 0, 180, 150, 0, 987,,, ,,,,, 63, 17,, 62
Loading