diff --git a/ehfheatwaves/__init__.py b/ehfheatwaves/__init__.py index f8f10b3..828ae0f 100644 --- a/ehfheatwaves/__init__.py +++ b/ehfheatwaves/__init__.py @@ -4,5 +4,5 @@ from ehfheatwaves.ehfheatwaves import * -__version__ = '1.1.5' +__version__ = '1.2.0' diff --git a/ehfheatwaves/ehfheatwaves.py b/ehfheatwaves/ehfheatwaves.py old mode 100644 new mode 100755 index ac7290b..b6af65f --- a/ehfheatwaves/ehfheatwaves.py +++ b/ehfheatwaves/ehfheatwaves.py @@ -106,8 +106,14 @@ def identify_hw(ehfs:np.ndarray)->tuple: """ # Agregate consecutive days with EHF>0 # First day contains duration - events = (ehfs>0.).astype(int) - events[events.mask==True] = 0 + if np.isnan(ehfs).any(): # then ehfs is a view to tmax or tmin + # This is handled differently to EHFs because tmax or tmin could theoretically hold -ve + # values. Therefore we identify values of the exceedence index that are not nan values as + # heatwaves. + events = np.logical_not(np.isnan(ehfs)).astype(int) + else: # ehfs is a view to actual EHF values + events = (ehfs>0.0).astype(int) + events[events.mask==True] = 0 for i in range(events.shape[0] - 2, -1, -1): events[i,events[i,...]>0] = events[i+1,events[i,...]>0]+1 @@ -147,8 +153,14 @@ def identify_semi_hw(ehfs:np.ndarray)->tuple: """ # Agregate consecutive days with EHF>0 # First day contains duration - events = (ehfs>0.).astype(int) - events[events.mask==True] = 0 + if np.isnan(ehfs).any(): # then ehfs is a view to tmax or tmin + # This is handled differently to EHFs because tmax or tmin could theoretically hold -ve + # values. Therefore we identify values of the exceedence index that are not nan values as + # heatwaves. + events = np.logical_not(np.isnan(ehfs)).astype(int) + else: # ehfs is a view to actual EHF values + events = (ehfs>0.0).astype(int) + events[events.mask==True] = 0 for i in range(events.shape[0] - 2, -1, -1): events[i,events[i,...]>0] = events[i+1,events[i,...]>0]+1 @@ -473,12 +485,14 @@ def main(): idoy = i-timedata.daysinyear*int((i+1)/timedata.daysinyear) txexceed[i,...] = tmax[i,...]>txpct[idoy,...] txexceed[txexceed>0] = tmax[txexceed>0] + txexceed[txexceed==0] = np.nan # need nans to be identified correctly if options.keeptmin: tnexceed = np.ma.ones(tmin.shape)*np.nan for i in range(0, tmin.shape[0]): idoy = i-timedata.daysinyear*int((i+1)/timedata.daysinyear) tnexceed[i,...] = tmin[i,...]>tnpct[idoy,...] tnexceed[tnexceed>0] = tmin[tnexceed>0] + tnexceed[tnexceed==0] = np.nan # need nans to be identified correctly # Calculate daily output if options.dailyout or options.dailyout: event, ends = identify_hw(EHF) diff --git a/tests/run_tests.py b/tests/run_tests.py index e51575f..9f56a0c 100755 --- a/tests/run_tests.py +++ b/tests/run_tests.py @@ -182,13 +182,24 @@ class TestIdentifyHW(unittest.TestCase): """Tests for the identify_hw function.""" # Some example EHF index data - ehfdata = np.ma.array([-1.3, -0.3, 3.4, 8.5, -0.4, # Not a heatwave - 5.6, 10.2, 20.4, -1.4, # a three day heatwave starting on day 6 - 7.8, 15.5, 16.9, 17.9, 30.2, -3.3]) # a 5 day heatwave starting on day 10 + ehfdata = np.ma.array([ + -1.3, -0.3, 3.4, 8.5, -0.4, # Not a heatwave + 5.6, 10.2, 20.4, -1.4, # a three day heatwave starting on day 6 + 7.8, 15.5, 16.9, 17.9, 30.2, -3.3, # a 5 day heatwave starting on day 10 + ]) ehfdata[ehfdata<0] = 0 known_events = np.ma.array([0,0,0,0,0,1,1,1,0,1,1,1,1,1,0]) known_ends = np.ma.array([0,0,0,0,0,3,0,0,0,5,0,0,0,0,0]) + # Example data for a heatwave with negative temperatures. + tmindata = np.ma.array([ + np.nan,np.nan,np.nan,np.nan, # does not exceed threshold so is masked out. + -2.5,-1.9,-0.2,0.5,2.1, # Heatwaves with negative values. + np.nan,np.nan, # does not exceed threshold. + ]) + tmin_known_events = np.ma.array([0,0,0,0,1,1,1,1,1,0,0]) + tmin_known_ends = np.ma.array([0,0,0,0,5,0,0,0,0,0,0]) + def testReturnTupple(self): """Should return a tupple containging the event indicator and the durations (numpy.ndarray).""" result = ehfheatwaves.identify_hw(self.ehfdata) @@ -209,6 +220,12 @@ def testShape(self): self.assertEqual(events.shape, input_shape) self.assertEqual(ends.shape, input_shape) + def testNegativeTemperatures(self): + """Tmin or tmax exceedances with negative values should be identified correctly.""" + events, ends = ehfheatwaves.identify_hw(self.tmindata) + self.assertTrue((events==self.tmin_known_events).all()) + self.assertTrue((ends==self.tmin_known_ends).all()) + class TestIdentifySemiHW(unittest.TestCase): """Tests for the identify_semi_hw function."""