From ff21ebcea3209ca1fe1d848e64cd207655ebf64a Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Dion=20H=C3=A4fner?= <dion.haefner@nbi.ku.dk>
Date: Fri, 19 Mar 2021 11:08:43 +0100
Subject: [PATCH] fix metadata; suppress warnings

---
 fowd/operators.py  | 47 +++++++++++++++++++++++++++-------------------
 fowd/output.py     | 17 ++++++++++++++---
 fowd/processing.py |  2 +-
 3 files changed, 43 insertions(+), 23 deletions(-)

diff --git a/fowd/operators.py b/fowd/operators.py
index 6cb181e..541fcdc 100644
--- a/fowd/operators.py
+++ b/fowd/operators.py
@@ -7,6 +7,7 @@
 import os
 import math
 import hashlib
+import warnings
 
 import numpy as np
 import scipy.stats
@@ -81,25 +82,33 @@ def compute_dynamic_window_size(t, elevation, min_length, max_length, num_window
         window_stationarity = []
         window_offsets = np.linspace(0, window_length, num_samples, endpoint=False, dtype='int')
 
-        for window_offset in window_offsets:
-            if window_offset >= len(t):
-                continue
-
-            current_time_idx = window_offset
-            window_stds = []
-            while True:
-                next_time = t[current_time_idx] + np.timedelta64(window_length, 's')
-                if next_time > t[-1]:
-                    break
-
-                next_time_idx = get_time_index(next_time, t)
-                window_stds.append(np.nanstd(elevation[current_time_idx:next_time_idx], ddof=1.5))
-                current_time_idx = next_time_idx
-
-            window_stds = np.asarray(window_stds)
-            window_stationarity.append(np.nanstd(window_stds[1:] / window_stds[:-1] - 1, ddof=1))
-
-        mean_window_stationarity = np.nanmean(window_stationarity)
+        # suppress nan warnings
+        with warnings.catch_warnings():
+            warnings.simplefilter('ignore', category=RuntimeWarning)
+
+            for window_offset in window_offsets:
+                if window_offset >= len(t):
+                    continue
+
+                current_time_idx = window_offset
+                window_stds = []
+                while True:
+                    next_time = t[current_time_idx] + np.timedelta64(window_length, 's')
+                    if next_time > t[-1]:
+                        break
+
+                    next_time_idx = get_time_index(next_time, t)
+                    window_stds.append(
+                        np.nanstd(elevation[current_time_idx:next_time_idx], ddof=1.5)
+                    )
+                    current_time_idx = next_time_idx
+
+                window_stds = np.asarray(window_stds)
+                window_stationarity.append(
+                    np.nanstd(window_stds[1:] / window_stds[:-1] - 1, ddof=1)
+                )
+
+            mean_window_stationarity = np.nanmean(window_stationarity)
 
         if mean_window_stationarity < best_window_stationarity:
             best_window_length = window_length
diff --git a/fowd/output.py b/fowd/output.py
index 0764385..255063b 100644
--- a/fowd/output.py
+++ b/fowd/output.py
@@ -191,7 +191,7 @@
                 dtype='int64',
                 attrs=dict(
                     long_name='Length of dynamically computed sea state window',
-                    units='seconds',
+                    units='minutes',
                 )
             )
         })
@@ -622,9 +622,20 @@ def write_records(wave_record_iterator, filename, station_name, extra_metadata=N
                 if np.issubdtype(data.dtype, np.datetime64):
                     data = (data - np.datetime64(TIME_ORIGIN)) / np.timedelta64(1, 'ms')
 
-                # convert timedelta64 to seconds
+                # convert timedelta64 to target unit
                 if np.issubdtype(data.dtype, np.timedelta64):
-                    data = data / np.timedelta64(1, 's')
+                    unit = meta.get('attrs', {}).get('units')
+
+                    if unit == 'seconds':
+                        target_unit = 's'
+                    elif unit == 'minutes':
+                        target_unit = 'm'
+                    elif unit == 'hours':
+                        target_unit = 'h'
+                    else:
+                        raise ValueError(f'Got unrecognized time unit {unit}')
+
+                    data = data / np.timedelta64(1, target_unit)
 
                 v[0, chunk_slice, ...] = data
 
diff --git a/fowd/processing.py b/fowd/processing.py
index 8f4772f..4e41ae6 100644
--- a/fowd/processing.py
+++ b/fowd/processing.py
@@ -224,7 +224,7 @@ def get_dynamic_window_size(current_window_size, current_time, time, elevation,
         num_samples=NUM_DYNAMIC_WINDOW_SAMPLES,
     )
 
-    return dynamic_sea_state_period // 60, current_time
+    return np.timedelta64(dynamic_sea_state_period, 's'), current_time
 
 
 def compute_wave_records(time, elevation, elevation_normalized, outfile, statefile,