Add Events from External Data to Timetable
This example shows how to use events with timestamped data stored in timetables. An event consists of a timestamp (when something happened), often a description (what happened), and sometimes additional information about the event. The timestamp for an instantaneous event is a single point in time, either a datetime
or duration
value. The timestamp for an interval event consists of two points in time that define the start and end of the event.
You can use events to find or mark data of interest for plotting and analysis. This example shows how you can define events in data using information that is external to your data set. To find events that are within your timestamped data, see Find Events in Timestamped Data.
You can store events separately from the main data set. A timetable is a convenient way to store events because it can hold the event timestamps along with their descriptions and any other information. You can use the event timestamps as subscripts into the main data set to select data at, before, or after an event, or between two events. You can also store events in the same timetable as the data set they describe by using one or more additional timetable variables. This example shows how to use both representations. The representation that is more useful depends on the data analysis you plan to perform.
Define Data Set Using Variations in Length of Day
By definition a day is 86,400 seconds long, where the second has a precise definition in the International System of Units (SI). However, the length of a day actually varies due to several physical causes. It varies with the seasons by as much as 30 seconds over and 21 seconds under the SI definition because of the eccentricity of Earth's orbit and the tilt of its axis. Averaging these seasonal effects enables the definition of the mean solar day, which does not vary in length over a year.
Also, there is a very long-term slowing in the rotational speed of the Earth due to tidal interaction with the moon; a smaller, opposite, shorter-term component believed to be due to melting of continental ice sheets; very short-term cycles on the order of decades; and unpredictable fluctuations due to geological events and other causes. Because of those effects, the length of a mean solar day might increase or decrease. In recent decades, it has fluctuated up and down, but has mostly been 1–3 milliseconds longer than 86,400 seconds. That difference is known as the excess Length of Day, or excess LOD.
For this example, create a timetable that contains the excess LOD for every day from January 1, 1962, to the present. The International Earth Rotation and Reference Systems Service (IERS) collects and publishes this data. However, this data needs preprocessing before storing in a MATLAB timetable because the dates are modified Julian dates. To read the IERS data into a table, use the readtable
function. Rename the two variables of interest to MJD
and ExcessLOD
.
file = "https://datacenter.iers.org/data/latestVersion/223_EOP_C04_14.62-NOW.IAU1980223.txt"; IERSdata = readtable(file,"NumHeaderLines",14); IERSdata.Properties.VariableNames([4 8]) = ["MJD","ExcessLOD"];
To store the excess LOD values in a timetable, convert the modified Julian dates to datetime
values. Use the datetime
function with the "ConvertFrom","mjd"
name-value argument. Convert the excess LOD values to duration
values by using the seconds
function. Then convert IERSdata
from a table to a timetable using the table2timetable
function.
IERSdata.Date = datetime(IERSdata.MJD,"ConvertFrom","mjd"); IERSdata.ExcessLOD = seconds(IERSdata.ExcessLOD); IERSdata = table2timetable(IERSdata(:,["Date","ExcessLOD"]))
IERSdata=22017×1 timetable
Date ExcessLOD
___________ ____________
01-Jan-1962 0.001723 sec
02-Jan-1962 0.001669 sec
03-Jan-1962 0.001582 sec
04-Jan-1962 0.001496 sec
05-Jan-1962 0.001416 sec
06-Jan-1962 0.001382 sec
07-Jan-1962 0.001413 sec
08-Jan-1962 0.001505 sec
09-Jan-1962 0.001628 sec
10-Jan-1962 0.001738 sec
11-Jan-1962 0.001794 sec
12-Jan-1962 0.001774 sec
13-Jan-1962 0.001667 sec
14-Jan-1962 0.00151 sec
15-Jan-1962 0.001312 sec
16-Jan-1962 0.001112 sec
⋮
Plot the excess LOD as a function of time. The excess LOD is currently decreasing on average but has remained positive except during very brief periods.
plot(IERSdata.Date,IERSdata.ExcessLOD,"b-"); ylabel("Excess LOD");
Calculate Cumulative Excess Length of Day
A clock that defines a day as 86,400 SI seconds is effectively running fast with respect to the Earth's rotation, gaining time with respect to the sun each day. A few extra milliseconds per day might seem unimportant, but the excess LOD since 1962 accumulates.
Add the cumulative excess LOD to IERSdata
as another table variable. To calculate the cumulative sum of the excess LOD, use the cumsum
function. Then use the stackedplot
function to plot excess LOD and cumulative excess LOD together. This plot shows that the shift due to the accumulated excess LOD since the 1960s has been less than one minute.
IERSdata.CumulativeELOD = [0; cumsum(IERSdata.ExcessLOD(1:end-1))]; h = stackedplot(IERSdata,["ExcessLOD","CumulativeELOD"])
h = StackedLineChart with properties: SourceTable: [22017×2 timetable] DisplayVariables: {'ExcessLOD' 'CumulativeELOD'} Color: [0 0.4470 0.7410] LineStyle: '-' LineWidth: 0.5000 Marker: 'none' MarkerSize: 6 Show all properties
Define Leap Seconds as Events in Separate TImetable
At midnight on January 1, 1972, the current version of the system of time known as Coordinated Universal Time (UTC) was enacted, which uses 86,400 SI seconds per day. On that date, UTC was defined to be roughly in sync with solar time at the Greenwich Meridian, or more precisely, with the time known as UT1. However, timekeeping based on exactly 86,400 SI seconds per day would drift away from our physical experience of solar time because of accumulated excess LOD. So when needed, an extra leap second adjustment is inserted to keep UTC approximately in sync with solar time. Without these leap second adjustments after 1972, accumulated excess LOD would have caused UTC to drift away from UT1 over the years. Within decades, the difference would have grown to tens of seconds.
Show the difference since 1972. First, select the post-1971 excess LOD data by subscripting into IERSdata
with a time range starting on January 1, 1972. To create that time range, use the timerange
function.
IERSdata1972 = IERSdata(timerange("1972-01-01",Inf),"ExcessLOD")
IERSdata1972=18365×1 timetable
Date ExcessLOD
___________ ____________
01-Jan-1972 0.002539 sec
02-Jan-1972 0.002708 sec
03-Jan-1972 0.002897 sec
04-Jan-1972 0.003065 sec
05-Jan-1972 0.003177 sec
06-Jan-1972 0.00322 sec
07-Jan-1972 0.003201 sec
08-Jan-1972 0.003137 sec
09-Jan-1972 0.003033 sec
10-Jan-1972 0.002898 sec
11-Jan-1972 0.002772 sec
12-Jan-1972 0.002672 sec
13-Jan-1972 0.002621 sec
14-Jan-1972 0.002642 sec
15-Jan-1972 0.00274 sec
16-Jan-1972 0.002937 sec
⋮
Add another variable with the unadjusted difference since 1972. Start with the cumulative excess LOD on January 1, 1972, which was 0.0454859 second. For every date that follows, calculate the unadjusted difference. By the beginning of 2017, the unadjusted difference would have grown to about 26 seconds. Display IERSdata1972
using a time range that includes the start of 2017.
DiffUT1_1972 = seconds(0.0454859); IERSdata1972.UnadjustedDiff = DiffUT1_1972 + [0; cumsum(IERSdata1972.ExcessLOD(1:end-1))]; IERSdata1972(timerange("2016-12-29","2017-01-04"),:)
ans=6×2 timetable
Date ExcessLOD UnadjustedDiff
___________ _____________ ______________
29-Dec-2016 0.0008055 sec 26.402 sec
30-Dec-2016 0.0008525 sec 26.402 sec
31-Dec-2016 0.0009173 sec 26.403 sec
01-Jan-2017 0.001016 sec 26.404 sec
02-Jan-2017 0.0011845 sec 26.405 sec
03-Jan-2017 0.0013554 sec 26.406 sec
This drift is what leap seconds are designed to mitigate. You can think of each leap second that is inserted into the UTC timeline as an event. These events are not part of the LOD data. However, the leapseconds
function lists each of the leap second events that has occurred since 1972. (The leap second events listed by leapseconds
come from another data set provided by the IERS.) The timetable returned by leapseconds
contains a timestamp, a description of each event (+
for leap second insertions, -
for removals), and the cumulative number of leap seconds up to and including the event. To date, all leap second events have been insertions.
lsEvents = leapseconds()
lsEvents=27×2 timetable
Date Type CumulativeAdjustment
___________ ____ ____________________
30-Jun-1972 + 1 sec
31-Dec-1972 + 2 sec
31-Dec-1973 + 3 sec
31-Dec-1974 + 4 sec
31-Dec-1975 + 5 sec
31-Dec-1976 + 6 sec
31-Dec-1977 + 7 sec
31-Dec-1978 + 8 sec
31-Dec-1979 + 9 sec
30-Jun-1981 + 10 sec
30-Jun-1982 + 11 sec
30-Jun-1983 + 12 sec
30-Jun-1985 + 13 sec
31-Dec-1987 + 14 sec
31-Dec-1989 + 15 sec
31-Dec-1990 + 16 sec
⋮
At this point, the main LOD data is in the IERSdata1972
timetable. The leap second events are in the lsEvents
timetable.
Use Event Timestamps to Select Data
Add the unadjusted difference, or cumulative excess LOD, at each date in lsEvents
. To find the unadjusted differences on these dates, index into the UnadjustedDiff
variable of IERSdata1972
and select the unadjusted difference at the timestamps from lsEvents
. Then append those differences to each event in lsEvents
as additional information about the event.
lsEvents.UnadjustedDiff = IERSdata1972.UnadjustedDiff(lsEvents.Date)
lsEvents=27×3 timetable
Date Type CumulativeAdjustment UnadjustedDiff
___________ ____ ____________________ ______________
30-Jun-1972 + 1 sec 0.63494 sec
31-Dec-1972 + 2 sec 1.1864 sec
31-Dec-1973 + 3 sec 2.2976 sec
31-Dec-1974 + 4 sec 3.289 sec
31-Dec-1975 + 5 sec 4.2718 sec
31-Dec-1976 + 6 sec 5.3334 sec
31-Dec-1977 + 7 sec 6.3469 sec
31-Dec-1978 + 8 sec 7.398 sec
31-Dec-1979 + 9 sec 8.3526 sec
30-Jun-1981 + 10 sec 9.6281 sec
30-Jun-1982 + 11 sec 10.389 sec
30-Jun-1983 + 12 sec 11.249 sec
30-Jun-1985 + 13 sec 12.451 sec
31-Dec-1987 + 14 sec 13.634 sec
31-Dec-1989 + 15 sec 14.669 sec
31-Dec-1990 + 16 sec 15.379 sec
⋮
Plot Events Against Data
The lsEvents
timetable shows that at each event, the IERS inserted a leap second whenever the excess LOD accumulated by roughly one additional second. To confirm this observation visually, plot the instantaneous leap second events as points overlaid on the LOD data over time.
plot(IERSdata1972.Date,IERSdata1972.UnadjustedDiff); hold on plot(lsEvents.Date,IERSdata1972.UnadjustedDiff(lsEvents.Date),"r+"); hold off
Represent Events in Additional Timetable Variable
Each leap second is an instantaneous event that occurs at a specific time. Another way to represent them is as a state variable appended to the original LOD data. A state variable describes the "state" of the process being measured at each time point in the data rather than tagging a specific instant. One possible state variable for this leap second data is an indicator variable that records an adjustment on each event date, with missing values on all the other dates. Assign the leap second event types to a state variable in IERSdata1972
, using timestamps from lsEvents
to assign event types to the correct rows. This approach works because all timestamps from lsEvents
are also timestamps in IERSdata1972
. The assignment automatically fills the other elements of the new Adjustment
variable with missing values.
IERSdata1972.Adjustment(lsEvents.Date) = lsEvents.Type; IERSdata1972(timerange("2016-12-29","2017-01-04"),:)
ans=6×3 timetable
Date ExcessLOD UnadjustedDiff Adjustment
___________ _____________ ______________ ___________
29-Dec-2016 0.0008055 sec 26.402 sec <undefined>
30-Dec-2016 0.0008525 sec 26.402 sec <undefined>
31-Dec-2016 0.0009173 sec 26.403 sec +
01-Jan-2017 0.001016 sec 26.404 sec <undefined>
02-Jan-2017 0.0011845 sec 26.405 sec <undefined>
03-Jan-2017 0.0013554 sec 26.406 sec <undefined>
Another possibility, more useful in this example, is to add a state variable that indicates the cumulative sum of leap seconds at any given date in the LOD data. Leap seconds occur at the end of a day, so first assign a value of one second to the day after the date for each event.
IERSdata1972.Adjustment = []; % remove the previous state variable IERSdata1972.Adjustment(lsEvents.Date+caldays(1)) = seconds(1); IERSdata1972(timerange("2016-12-29","2017-01-04"),:)
ans=6×3 timetable
Date ExcessLOD UnadjustedDiff Adjustment
___________ _____________ ______________ __________
29-Dec-2016 0.0008055 sec 26.402 sec 0 sec
30-Dec-2016 0.0008525 sec 26.402 sec 0 sec
31-Dec-2016 0.0009173 sec 26.403 sec 0 sec
01-Jan-2017 0.001016 sec 26.404 sec 1 sec
02-Jan-2017 0.0011845 sec 26.405 sec 0 sec
03-Jan-2017 0.0013554 sec 26.406 sec 0 sec
Then use a cumulative sum to fill in the Adjustment
state variable at all other times in the data. Instead of containing values of +
and <undefined>
, the state variable Adjustment
now lists the number of leap seconds added as of the date listed in each row of IERSdata1972
.
IERSdata1972.Adjustment = cumsum(IERSdata1972.Adjustment); IERSdata1972(timerange("2016-12-29","2017-01-04"),:)
ans=6×3 timetable
Date ExcessLOD UnadjustedDiff Adjustment
___________ _____________ ______________ __________
29-Dec-2016 0.0008055 sec 26.402 sec 26 sec
30-Dec-2016 0.0008525 sec 26.402 sec 26 sec
31-Dec-2016 0.0009173 sec 26.403 sec 26 sec
01-Jan-2017 0.001016 sec 26.404 sec 27 sec
02-Jan-2017 0.0011845 sec 26.405 sec 27 sec
03-Jan-2017 0.0013554 sec 26.406 sec 27 sec
Although not used here, the retime
function is often a useful tool for converting instantaneous events to a state variable.
While the new state variable represents the instantaneous events, it is just like all the other variables in IERSdata1972
, with a value that is defined at each time. You can use it to compute the actual differences between UTC and UT1, given all the leap second adjustments, by subtracting it from the unadjusted difference. The convention is to compute the actual difference with the opposite sign as UT1 – UTC and denote it as DUT1.
IERSdata1972.DUT1 = IERSdata1972.UnadjustedDiff - IERSdata1972.Adjustment; IERSdata1972(timerange("2016-12-29","2017-01-04"),:)
ans=6×4 timetable
Date ExcessLOD UnadjustedDiff Adjustment DUT1
___________ _____________ ______________ __________ ____________
29-Dec-2016 0.0008055 sec 26.402 sec 26 sec 0.40159 sec
30-Dec-2016 0.0008525 sec 26.402 sec 26 sec 0.40239 sec
31-Dec-2016 0.0009173 sec 26.403 sec 26 sec 0.40324 sec
01-Jan-2017 0.001016 sec 26.404 sec 27 sec -0.59584 sec
02-Jan-2017 0.0011845 sec 26.405 sec 27 sec -0.59482 sec
03-Jan-2017 0.0013554 sec 26.406 sec 27 sec -0.59364 sec
Plot State Variable Against Data
The IERS makes leap second adjustments to keep UTC roughly in sync with UT1. To show how this adjustment works, plot the Adjustment
state variable as a piecewise-constant step function over time. The step function approximates the accumulated excess LOD, so subtracting it keeps DUT1 small. For visual confirmation, plot DUT1 along the same time axis by using the stackedplot
function. Plot UnadjustedDiff
and Adjustment
together along one y-axis, with a legend, and plot DUT1
along a second y-axis.
h = stackedplot(IERSdata1972,{["UnadjustedDiff","Adjustment"],"DUT1"})
h = StackedLineChart with properties: SourceTable: [18365×4 timetable] DisplayVariables: {{1×2 cell} 'DUT1'} Color: [0 0.4470 0.7410] LineStyle: '-' LineWidth: 0.5000 Marker: 'none' MarkerSize: 6 Show all properties
Add the value of DUT1 before and after each leap second to lsEvents
. To get these values, use the instantaneous form of the leap second events. Index into the DUT1
variable of IERSdata1972
using the timestamps from lsEvents
. Get the value of DUT1
on the date each leap second was added and on the date of the day after. In general, DUT1
has a small negative value the day after a leap second is added.
lsEvents.UnadjustedDiff = []; % remove the previous additional event information
lsEvents.DUT1_before = IERSdata1972.DUT1(lsEvents.Date);
lsEvents.DUT1_after = IERSdata1972.DUT1(lsEvents.Date+caldays(1))
lsEvents=27×4 timetable
Date Type CumulativeAdjustment DUT1_before DUT1_after
___________ ____ ____________________ ___________ ____________
30-Jun-1972 + 1 sec 0.63494 sec -0.36236 sec
31-Dec-1972 + 2 sec 0.18641 sec -0.8107 sec
31-Dec-1973 + 3 sec 0.29763 sec -0.69952 sec
31-Dec-1974 + 4 sec 0.28898 sec -0.70817 sec
31-Dec-1975 + 5 sec 0.2718 sec -0.72578 sec
31-Dec-1976 + 6 sec 0.33338 sec -0.66411 sec
31-Dec-1977 + 7 sec 0.34693 sec -0.64986 sec
31-Dec-1978 + 8 sec 0.39802 sec -0.59892 sec
31-Dec-1979 + 9 sec 0.35261 sec -0.64508 sec
30-Jun-1981 + 10 sec 0.6281 sec -0.37053 sec
30-Jun-1982 + 11 sec 0.38922 sec -0.60911 sec
30-Jun-1983 + 12 sec 0.24852 sec -0.75008 sec
30-Jun-1985 + 13 sec 0.45106 sec -0.54797 sec
31-Dec-1987 + 14 sec 0.63443 sec -0.36417 sec
31-Dec-1989 + 15 sec 0.66947 sec -0.32879 sec
31-Dec-1990 + 16 sec 0.37936 sec -0.61885 sec
⋮
Choose Events or State Variables
The two representations, a separate list of instantaneous events and a state variable defined at all times, are conceptually different. In some cases, there might not be a useful definition of a "state" that corresponds to the periods between instantaneous events. But when there is, the two representations are equivalent and useful in similar ways. For example, to select all the data after the 22nd leap second but before the 23rd, you can use event times and timerange
.
from22to23 = timerange(lsEvents.Date(22),lsEvents.Date(23),"openleft");
data22to23 = IERSdata1972(from22to23,:)
data22to23=2557×4 timetable
Date ExcessLOD UnadjustedDiff Adjustment DUT1
___________ _____________ ______________ __________ ____________
01-Jan-1999 0.0009738 sec 21.283 sec 22 sec -0.71702 sec
02-Jan-1999 0.000888 sec 21.284 sec 22 sec -0.71604 sec
03-Jan-1999 0.0008605 sec 21.285 sec 22 sec -0.71516 sec
04-Jan-1999 0.0008798 sec 21.286 sec 22 sec -0.71429 sec
05-Jan-1999 0.0008935 sec 21.287 sec 22 sec -0.71342 sec
06-Jan-1999 0.0009693 sec 21.287 sec 22 sec -0.71252 sec
07-Jan-1999 0.0010658 sec 21.288 sec 22 sec -0.71155 sec
08-Jan-1999 0.0010857 sec 21.29 sec 22 sec -0.71049 sec
09-Jan-1999 0.0010715 sec 21.291 sec 22 sec -0.7094 sec
10-Jan-1999 0.0010519 sec 21.292 sec 22 sec -0.70833 sec
11-Jan-1999 0.0010021 sec 21.293 sec 22 sec -0.70728 sec
12-Jan-1999 0.0008986 sec 21.294 sec 22 sec -0.70628 sec
13-Jan-1999 0.0007891 sec 21.295 sec 22 sec -0.70538 sec
14-Jan-1999 0.0007236 sec 21.295 sec 22 sec -0.70459 sec
15-Jan-1999 0.0006842 sec 21.296 sec 22 sec -0.70386 sec
16-Jan-1999 0.0006091 sec 21.297 sec 22 sec -0.70318 sec
⋮
As an alternative, you can create a logical subscript from the state variable that selects the same data.
from22to23 = (IERSdata1972.Adjustment == seconds(22)); data22to23 = IERSdata1972(from22to23,:)
data22to23=2557×4 timetable
Date ExcessLOD UnadjustedDiff Adjustment DUT1
___________ _____________ ______________ __________ ____________
01-Jan-1999 0.0009738 sec 21.283 sec 22 sec -0.71702 sec
02-Jan-1999 0.000888 sec 21.284 sec 22 sec -0.71604 sec
03-Jan-1999 0.0008605 sec 21.285 sec 22 sec -0.71516 sec
04-Jan-1999 0.0008798 sec 21.286 sec 22 sec -0.71429 sec
05-Jan-1999 0.0008935 sec 21.287 sec 22 sec -0.71342 sec
06-Jan-1999 0.0009693 sec 21.287 sec 22 sec -0.71252 sec
07-Jan-1999 0.0010658 sec 21.288 sec 22 sec -0.71155 sec
08-Jan-1999 0.0010857 sec 21.29 sec 22 sec -0.71049 sec
09-Jan-1999 0.0010715 sec 21.291 sec 22 sec -0.7094 sec
10-Jan-1999 0.0010519 sec 21.292 sec 22 sec -0.70833 sec
11-Jan-1999 0.0010021 sec 21.293 sec 22 sec -0.70728 sec
12-Jan-1999 0.0008986 sec 21.294 sec 22 sec -0.70628 sec
13-Jan-1999 0.0007891 sec 21.295 sec 22 sec -0.70538 sec
14-Jan-1999 0.0007236 sec 21.295 sec 22 sec -0.70459 sec
15-Jan-1999 0.0006842 sec 21.296 sec 22 sec -0.70386 sec
16-Jan-1999 0.0006091 sec 21.297 sec 22 sec -0.70318 sec
⋮
Both representations have their uses. For example, while each event can be plotted as a point, the state variable is the more convenient form for highlighting regions between events in a plot. To highlight the region between the 22nd and 23rd leap second, use the isbetween
function with the Adjustment
state variable from IERSdata1972
.
plot(IERSdata1972.Date,IERSdata1972.UnadjustedDiff); hold on from22to23 = isbetween(IERSdata1972.Adjustment,seconds(22),seconds(23)); plot(IERSdata1972.Date(from22to23),IERSdata1972.UnadjustedDiff(from22to23),"r-","LineWidth",4); hold off
You can switch between the two representations. In this case, to get the leap second event dates from the adjustment state variable, find the locations where the adjustment changes, and subtract one day. The eventTimes
output represents the dates on which leap seconds were added, which are instantaneous events.
eventTimes = IERSdata1972.Date(diff(IERSdata1972.Adjustment) ~= 0) - caldays(1)
eventTimes = 27×1 datetime
29-Jun-1972
30-Dec-1972
30-Dec-1973
30-Dec-1974
30-Dec-1975
30-Dec-1976
30-Dec-1977
30-Dec-1978
30-Dec-1979
29-Jun-1981
29-Jun-1982
29-Jun-1983
29-Jun-1985
30-Dec-1987
30-Dec-1989
30-Dec-1990
29-Jun-1992
29-Jun-1993
29-Jun-1994
30-Dec-1995
29-Jun-1997
30-Dec-1998
30-Dec-2005
30-Dec-2008
29-Jun-2012
29-Jun-2015
30-Dec-2016
Whether you use instantaneous events, state variables, or switch between them depends on which representation is more convenient and useful for the data analysis that you plan to perform. Both representations are useful ways to add information about events to your timestamped data in a timetable.
See Also
timetable
| datetime
| seconds
| leapseconds
| timerange
| isbetween
| stackedplot
| readtable
| table2timetable
| cumsum