From c251b8e02c111f6c32b3b6c031964873ad37335c Mon Sep 17 00:00:00 2001 From: arigdon-usgs <arigdon@usgs.gov> Date: Wed, 25 Jul 2018 11:17:35 -0600 Subject: [PATCH] Clean Timeseries fix to make sure start and end times are valid sample times --- geomagio/edge/EdgeFactory.py | 26 +++++++++++++-------- test/edge_test/EdgeFactory_test.py | 37 +++++++++++++++++++++++++++--- 2 files changed, 50 insertions(+), 13 deletions(-) diff --git a/geomagio/edge/EdgeFactory.py b/geomagio/edge/EdgeFactory.py index a26993d47..ef6d210c4 100644 --- a/geomagio/edge/EdgeFactory.py +++ b/geomagio/edge/EdgeFactory.py @@ -220,19 +220,25 @@ class EdgeFactory(TimeseriesFactory): for trace in timeseries: trace_starttime = obspy.core.UTCDateTime(trace.stats.starttime) trace_endtime = obspy.core.UTCDateTime(trace.stats.endtime) + trace_delta = trace.stats.delta + + if trace_starttime > starttime: + cnt = int((trace_starttime - starttime) / trace_delta) + if cnt > 0: + trace.data = numpy.concatenate([ + numpy.full(cnt, numpy.nan, dtype=numpy.float64), + trace.data]) + trace_starttime = trace_starttime - trace_delta * cnt + trace.stats.starttime = trace_starttime - if trace.stats.starttime > starttime: - cnt = int((trace_starttime - starttime) / trace.stats.delta) - trace.data = numpy.concatenate([ - numpy.full(cnt, numpy.nan, dtype=numpy.float64), - trace.data]) - trace.stats.starttime = starttime if trace_endtime < endtime: cnt = int((endtime - trace_endtime) / trace.stats.delta) - trace.data = numpy.concatenate([ - trace.data, - numpy.full(cnt, numpy.nan, dtype=numpy.float64)]) - trace.stats.endttime = endtime + if cnt > 0: + trace.data = numpy.concatenate([ + trace.data, + numpy.full(cnt, numpy.nan, dtype=numpy.float64)]) + trace_endtime = trace_endtime + trace_delta * cnt + trace.stats.endttime = trace_endtime def _convert_timeseries_to_decimal(self, stream): """convert geomag edge timeseries data stored as ints, to decimal by diff --git a/test/edge_test/EdgeFactory_test.py b/test/edge_test/EdgeFactory_test.py index e43055b50..4b583878f 100644 --- a/test/edge_test/EdgeFactory_test.py +++ b/test/edge_test/EdgeFactory_test.py @@ -1,10 +1,9 @@ """Tests for EdgeFactory.py""" -from obspy.core.utcdatetime import UTCDateTime -from obspy.core.stream import Stream -from obspy.core.trace import Trace +from obspy.core import Stats, Stream, Trace, UTCDateTime from geomagio.edge import EdgeFactory from nose.tools import assert_equals +import numpy as np def test__get_edge_network(): @@ -93,3 +92,35 @@ def dont_get_timeseries(): 'BOU', 'Expect timeseries to have stats') assert_equals(timeseries.select(channel='H')[0].stats.channel, 'H', 'Expect timeseries stats channel to be equal to H') + + +def test_clean_timeseries(): + """edge_test.EdgeFactory_test.test_clean_timeseries() + """ + edge_factory = EdgeFactory() + trace1 = _create_trace([1, 1, 1, 1, 1], 'H', UTCDateTime("2018-01-01")) + trace2 = _create_trace([2, 2],'E',UTCDateTime("2018-01-01")) + timeseries = Stream(traces=[trace1, trace2]) + edge_factory._clean_timeseries(timeseries,trace1.stats.starttime,trace1.stats.endtime) + assert_equals(len(trace1.data),len(trace2.data)) + assert_equals(trace1.stats.starttime,trace2.stats.starttime) + assert_equals(trace1.stats.endtime,trace2.stats.endtime) + # change starttime by less than 1 delta + starttime = trace1.stats.starttime + endtime = trace1.stats.endtime + edge_factory._clean_timeseries(timeseries,starttime-30,endtime+30) + assert_equals(trace1.stats.starttime,starttime) + # Change starttime by more than 1 delta + edge_factory._clean_timeseries(timeseries,starttime-90,endtime+90) + assert_equals(trace1.stats.starttime,starttime-60) + assert_equals(np.isnan(trace1.data[0]),np.isnan(np.NaN)) + + +def _create_trace(data,channel,starttime,delta=60.): + stats = Stats() + stats.channel = channel + stats.delta = delta + stats.starttime = starttime + stats.npts = len(data) + data = np.array(data,dtype=np.float64) + return Trace(data,stats) -- GitLab