Source code for irfpy.util.time_panel

''' Utility for time series panels, such as E-t diagram

.. codeauthor:: Yoshifumi Futaana

The function :func:`reform` will create a new matrix considering the data gap.
'''
from irfpy.util import exception

import numpy as np

[docs]def reform(time_list, value_matrix, sampling_time): ''' Time series data is reformed considering the sampling time. :param time_list: ``np.array`` of the time when the data is sampled. (N,) shape. Float or ``datetime.datetime`` object. :param value_matrix: ``np.array`` of the data with (N, S) shape. S is the number of the spectrum data. :param sampling_time: Intrinsic sampling time. ``Float`` with a unit of day if ``time_list`` is float, but ``datetime.timedelta`` if time_list is ``datetime.datetime``. Usually, the observed data is given by the time and spectrum. For example, >>> time_list = np.array([735000, 735001, 735002, 735004, 735005, 735006]) >>> img = np.array([[1, 1, 2,], [1, 2, 2,], ... [2, 2, 1,], [0, 2, 3,], ... [1, 3, 2,], [2, 1, 0,]]) A simple way to show this data is using ``meshcolor`` function. .. code-block:: python # Note: the second argument is "4", not "3". plt.meshcolor(t, np.arange(4), img.T) Then, you may get:: - 2 1 - v v 3| -| 2 2 1 1 3 2 2| -| 1 2 2 2 2 3 1| -| 1 1 2 2 0 1 0+------------- 0|1|2|3|4|5|6 +7.35e5 Here is two problems: 1. No data is shown for 735006. This is because the "end" time of the data is not specified. Most likely you want to plot the last (2, 1, 0) data between 735006 and "735007" 2. The data at 735003 is (likely) missing, because of some reasons. However, the program does not know it, so that the data (2, 2, 1) at 735002 is plotted until 735004. Both problems can be solved if you know the "sampling time". In this sample case, 1 day could be the most reasonable sampling time. >>> time_list2, value_matrix2 = reform(time_list, img, 1) will give you a reformed time list and value matrix (np.ma.array). >>> print(time_list2) [735000 735001 735002 735003 735004 735005 735006 735007] >>> print(value_matrix2) [[1.0 1.0 2.0] [1.0 2.0 2.0] [2.0 2.0 1.0] [-- -- --] [0.0 2.0 3.0] [1.0 3.0 2.0] [2.0 1.0 0.0]] ''' n_time = time_list.shape[0] n_value = value_matrix.shape[0] if n_time != n_value: raise exception.PyanaError( 'Element number error. time_list has {} element but value has {} element. They should be the same'.format( n_time, n_value)) ### First, the end time is added. This can be done for free. time_list2 = np.hstack([time_list, [time_list[-1] + sampling_time]]) value_matrix2 = value_matrix.copy() ### Calculate the delta. delta_time_list2 = time_list2[1:] - time_list2[:-1] ### Check if the delta > sampling_time, which include data gap datagap_index = np.where(delta_time_list2 > sampling_time)[0] for idx in datagap_index[::-1]: # Iterate reversely. Index saved then. time_before = time_list2[:idx + 1] time_after = time_list2[idx + 1:] time_list2 = np.hstack([time_before, [time_before[-1] + sampling_time], time_after]) value_before = value_matrix2[:idx + 1, :] value_after = value_matrix2[idx + 1:, :] insert_values = np.zeros_like(value_matrix[0, :]) + np.nan value_matrix2 = np.vstack([value_before, insert_values, value_after]) value_matrix2 = np.ma.masked_invalid(value_matrix2) return time_list2, value_matrix2
import unittest import doctest
[docs]def doctests(): return unittest.TestSuite(( doctest.DocTestSuite(), ))
if __name__ == '__main__': unittest.main(defaultTest='doctests')