罗兰谈MATLAB的艺术

将想法转化为MATLAB

那么闰秒是什么呢?

今天的文章由MathWorks开发团队的成员Peter Perkins提供给您。
在研究过MATLAB的一些时间和日期函数后,The MathWorks的人有时会问我一些关于日历和计时的问题。通常情况是这样的:1月初(或7月),我在午餐时间盯着手机看了一会儿,然后脱口而出:“ 公报C 发表!明年6月(或12月)不会有新的闰秒了!”,在尴尬的沉默之后,他说:“不会。 再一次 听过这句话的人会问:“闰秒到底是什么东西?”我不是专家,但我的标准回答是这样的:“所以,地球原来是在减速。”接下来是故事的其余部分。
在此过程中,我将展示对处理时间序列数据非常有用的时间表和日期时间的一些用法。

所以,事实证明地球正在减速……但不是最近

在日常生活中,我们认为一天的长度是太阳回到天顶所需要的时间,换句话说,地球绕地轴旋转一次。我们认为这是24小时,但事实上它是不同的(由于 开普勒第二定律 )由于地球轨道的偏心率(开普勒第一定律),地球的季节变化大约为+30秒/-21秒。“一天”的长度在一年的过程中变化如此之大是非常不方便的,因为如果我们认为每天是24*60*60=86400秒,那么6月的一秒长度就需要与12月的一秒长度不同,而精确测量流逝的时间是不可能的。这些季节对太阳日长度的影响可以在一年中逐渐消除,从而产生“平均太阳日”,其长度在一年中不会变化。基于平均太阳日的时间系统被称为“世界时”(UT)。
UT的问题在于,除了由偏心轨道引起的“几何效应”之外,地球的转速随着时间的推移并不是恒定的。事实上,由于潮汐与月球的相互作用,地球自转的速度长期减缓,据估计,在过去的6亿年里,一天的时间变长了几个小时。还有一个较短的部分,较小,但方向相反,据信是由于近千年来大陆冰原的融化。从历史上讲,自转变慢意味着太阳日的长度平均每世纪增加1.4-1.7毫秒。还有季节周期,以及由于地质事件和其他原因而不可预测的波动。在任何给定的几十年的时间里,甚至是一个 的意思是 太阳日可以增加,也可以减少,或者两者都有。
出于科学目的,在19世纪,秒被定义为 平均太阳日的1/86400 但是如果一天的长度随着时间的变化而变化,那么一秒实际上等于多少呢?这个定义是基于18世纪和19世纪收集的天文数据,有效地使用了1820年左右的平均太阳日。但即使到了19世纪晚期,人们也认识到当时的平均太阳日略长于86400秒。到了20世纪40年代和50年代,石英钟和原子钟使得地球自转速度的变化更加明显,并导致了 更精确的定义 一国际单位制(SI)秒的长度。毫不奇怪,地球的自转没有注意到这一点,平均太阳日的长度在20世纪继续上下波动。近几十年来,它大多比86400s长1-3毫秒。这种差异被称为“多余的一天长度(LOD)” 国际地球自转与参考系统服务
文件=“https://datacenter.iers.org/data/latestVersion/EOP_14_C04_IAU2000A_one_file_1962-now.txt”
IERSdata = readtable(文件,“NumHeaderLines”14);
IERSdata.Properties。VariableNames([4 8]) = [“MJD”“ExcessLOD”];
IERSdata。日期= datetime(IERSdata.MJD,“ConvertFrom”“mjd”);% modified Julian date -> datetime
IERSdata。ExcessLOD= seconds(IERSdata.ExcessLOD);
IERSdata = table2schedule (IERSdata(:,[ExcessLOD“日期”)));
情节(IERSdata.Date IERSdata.ExcessLOD,“b -”);
ylabel (“超额LOD”);

好吧,当然,但是什么总之是闰秒?

一些额外的ms似乎不重要,但是多余的LOD会累积起来。一个将“一天”定义为86400国际单位秒的时钟相对于地球自转有效地运行得很快,相对于太阳每天获得时间。例如,太阳在6月21日出现在天顶的时钟上的时间将连续每年晚零点几秒。不过,从20世纪60年代以来的多余LOD的总和来看,这种变化还不到一分钟,对于日常用途来说甚至都不明显,所以谁在乎呢?
IERSdata。累积elod = [0;cumsum (IERSdata.ExcessLOD (1: end-1)];
stackedplot (IERSdata [“ExcessLOD”“CumulativeELOD”]);
但在几个世纪的过程中,累积的过量LOD可能会累积到中午的太阳明显不再出现在其天顶附近。换句话说,基于每天86400国际单位秒的计时将偏离物理体验,即UT。这被认为是一个问题,尽管一个明显的事实是,在一个时区不同“端点”的城市之间,甚至在同一城市的标准时间和日光节约时间之间,中午太阳位置的变化比上个世纪的变化要大得多。1972年之前, 各种特别调整 与UT保持时间同步。到1972年, 时间系统的当前版本 被称为协调世界时(UTC),它使用原子钟定义的86400 SI秒,但偶尔会插入额外的“闰秒”调整,以保持UTC与UT ( UT1 更准确地说)。这个想法是,因为86400s相对于地球的自转来说有点太短了,所以UTC时钟相对于头顶的太阳运行得有点太快了,所以你只需要偶尔让它多走一秒,让它慢下来。也有可能地球自转在短期内会加速,使86400s也稍微加快一点 在这种情况下,闰秒可能最终是需要的 删除 .迄今为止,这种情况从未发生过;下文将详细介绍。
地球的自转是如此不可预测以至于 国际IERS仅提前6个月宣布闰秒 在1月和7月。这种特殊的时机使得在每个MATLAB版本中获取最新信息变得很困难。的第二个输出 leapseconds 功能)。

等一下!

UTC的当前定义是从1972年1月1日午夜开始采用的,因此我将只选择1971年之后的数据。
IERSdata1972 = IERSdata(时间范围(“1 - 1月- 1972”正)“ExcessLOD”);
20世纪60年代对民用计时所做的特别调整意味着1972年1月1日UTC与UT1相差约45毫秒。但是在1972年之后,如果没有闰秒调整,累积的过量LOD将导致UTC多年来偏离UT1。到2017年,这一差距将扩大到26秒左右。
DiffUT1_1972 = seconds(0.0454859);
IERSdata1972。UnadjDiff= DiffUT1_1972 + [0; cumsum(IERSdata1972.ExcessLOD(1:end-1))];
IERSdata1972 (timerange (”29日- 12月- 2016“4 - 1月- 2017”):)
ans = 6×2时间表
日期 ExcessLOD UnadjDiff
1 29日- 12月- 2016 0.0008055秒 26.402秒
2 2016年- 12月30日 0.0008525秒 26.402秒
3. 2016年- 12月31日 0.0009173秒 26.403秒
4 01 - 1月- 2017 0.001016秒 26.404秒
5 02 - 1月- 2017 0.0011845秒 26.405秒
6 03 - 1月- 2017 0.0013554秒 26.406秒
闰秒的设计就是为了缓解这种日益扩大的差异。的 leapseconds 函数枚举自1972年以来发生的每一个“闰秒事件”。发回的时间表 leapseconds 包含时间戳和事件的描述(到目前为止,所有的闰秒都是“正”插入,而不是“负”删除),以及到该事件为止并包括该事件在内的累计闰秒数。
lsEvents = leapseconds()
lsEvents = 27×2时间表
日期 类型 CumulativeAdjustment
1 1972年- 6月30日 + 1秒
2 1972年- 12月31日 + 2秒
3. 1973年- 12月31日 + 3秒
4 1974年- 12月31日 + 4秒
5 1975年- 12月31日 + 5秒
6 1976年- 12月31日 + 6秒
7 1977年- 12月31日 + 7秒
8 1978年- 12月31日 + 8秒
9 1979年- 12月31日 + 9秒
10 1981年- 6月30日 + 10秒
11 1982年- 6月30日 + 11秒
12 1983年- 6月30日 + 12秒
13 1985年- 6月30日 + 13秒
14 1987年- 12月31日 + 14秒
使用闰秒时间戳将索引到包含LOD数据的时间表中,我将添加如果没有插入闰秒,UTC和UT1之间在每个时刻的未调整差值。
lsEvents。UnadjDiff = IERSdata1972.UnadjDiff(l7ents . date)
lsEvents = 27×3的时间表
日期 类型 CumulativeAdjustment UnadjDiff
1 1972年- 6月30日 + 1秒 0.63494秒
2 1972年- 12月31日 + 2秒 1.1864秒
3. 1973年- 12月31日 + 3秒 2.2976秒
4 1974年- 12月31日 + 4秒 3.289秒
5 1975年- 12月31日 + 5秒 4.2718秒
6 1976年- 12月31日 + 6秒 5.3334秒
7 1977年- 12月31日 + 7秒 6.3469秒
8 1978年- 12月31日 + 8秒 7.398秒
9 1979年- 12月31日 + 9秒 8.3526秒
10 1981年- 6月30日 + 10秒 9.6281秒
11 1982年- 6月30日 + 11秒 10.389秒
12 1983年- 6月30日 + 12秒 11.249秒
13 1985年- 6月30日 + 13秒 12.451秒
14 1987年- 12月31日 + 14秒 13.634秒
查看每个闰秒事件中未调整的差异,很明显,每当额外的LOD累积大约一秒时,就会插入闰秒。将闰秒绘制为叠加在LOD数据上的点,可以直观地确认这一点。
情节(IERSdata1972.Date IERSdata1972.UnadjDiff);
ylabel (“累计超额LOD”);
持有
情节(lsEvents.Date IERSdata1972.UnadjDiff (lsEvents.Date),“r”。);
持有

闰秒为什么有效

每个闰秒都是发生在特定时间的瞬时事件。表示它们的另一种方法是将其作为附加到LOD数据的变量,该变量表示任意给定日期的闰秒累计和。闰秒发生在一天结束的时候,所以首先给每个事件的日期后一天赋一个1秒的值……
IERSdata1972.Adjustment(l7its . date +caldays(1)) = seconds(1);
IERSdata1972 (timerange (”29日- 12月- 2016“4 - 1月- 2017”):)
ans = 6×3的时间表
日期 ExcessLOD UnadjDiff 调整
1 29日- 12月- 2016 0.0008055秒 26.402秒 0秒
2 2016年- 12月30日 0.0008525秒 26.402秒 0秒
3. 2016年- 12月31日 0.0009173秒 26.403秒 0秒
4 01 - 1月- 2017 0.001016秒 26.404秒 1秒
5 02 - 1月- 2017 0.0011845秒 26.405秒 0秒
6 03 - 1月- 2017 0.0013554秒 26.406秒 0秒
...然后用累计和来填写 调整 所有其他时间的数据。
IERSdata1972。调整= cumsum(IERSdata1972.Adjustment);
IERSdata1972 (timerange (”29日- 12月- 2016“4 - 1月- 2017”):)
ans = 6×3的时间表
日期 ExcessLOD UnadjDiff 调整
1 29日- 12月- 2016 0.0008055秒 26.402秒 26秒
2 2016年- 12月30日 0.0008525秒 26.402秒 26秒
3. 2016年- 12月31日 0.0009173秒 26.403秒 26秒
4 01 - 1月- 2017 0.001016秒 26.404秒 27秒
5 02 - 1月- 2017 0.0011845秒 26.405秒 27秒
6 03 - 1月- 2017 0.0013554秒 26.406秒 27秒
这个新变量在数据中的每一个时间都被定义,因此它可以用来计算UTC和UT1之间的实际差值,给定所有的闰秒调整,通过从上面计算的未调整的差值中减去它。按照惯例,它被计算为UT1-UTC,相反的符号是我们所称的未调整的差值,并表示为DUT1。
IERSdata1972。DUT1= -(IERSdata1972.UnadjDiff - IERSdata1972.Adjustment);
IERSdata1972 (timerange (”29日- 12月- 2016“4 - 1月- 2017”):)
ans = 6×4时间表
日期 ExcessLOD UnadjDiff 调整 DUT1
1 29日- 12月- 2016 0.0008055秒 26.402秒 26秒 -0.40159秒
2 2016年- 12月30日 0.0008525秒 26.402秒 26秒 -0.40239秒
3. 2016年- 12月31日 0.0009173秒 26.403秒 26秒 -0.40324秒
4 01 - 1月- 2017 0.001016秒 26.404秒 27秒 0.59584秒
5 02 - 1月- 2017 0.0011845秒 26.405秒 27秒 0.59482秒
6 03 - 1月- 2017 0.0013554秒 26.406秒 27秒 0.59364秒
闰秒调整的目的是使UTC与UT1大致同步。策划 调整 随着时间的推移,分段常数阶跃函数说明了调整是如何工作的:阶跃函数大致近似累积的多余LOD,因此减去它保持DUT1较小。沿着同一时间轴绘制DUT1,可以确认调整使差值小于1秒。
stackedplot (IERSdata1972 {[“UnadjDiff”“调整”“DUT1”});

上升时间和下降时间

只有在需要时才会插入闰秒。从上面的一些数字可以看出,自1972年以来,闰秒之间的平均时间一直在增加。计算每十年发生的事件数量可以证明这一点。
十年=离散化(lsEvents。目前为止,“十年”“分类”);
柱状图(十年);
要了解原因,我们将回到多余的LOD数据。
自20世纪60年代以来,出现了几个过剩LOD在短期内下降的时期。平滑原始LOD数据可以更容易地看到局部行为。
IERSdata。smooththedelod = smoothdata(seconds(IERSdata.ExcessLOD),“黄土”“SmoothingFactor”、。4);
情节(IERSdata.Date IERSdata.ExcessLOD,“b -”);
持有
情节(IERSdata.Date IERSdata.SmoothedELOD,“g -”“线宽”2);
持有
ylabel (“超额LOD”);
从经过平滑处理的数据中,识别短期趋势改变方向的波峰和波谷,并得到每个波峰和波谷的平滑过剩LOD的日期和值。
peaks = find(islocalmax(IERSdata.SmoothedELOD));
trours = find(islocalmin(IERSdata.SmoothedELOD));
日期= IERSdata.Date([峰值;低谷]);
Value = IERSdata.SmoothedELOD([峰值;低谷]);
接下来,我将创建一个包含高峰/低谷事件的时间表,并在图上标记每个高峰/低谷事件。时间表的每一行都包含事件日期、事件类型和该事件的多余LOD值。
类型= categorical([0(大小(峰值));(大小(波谷))],[0 1]、[“峰”“槽”]);
extremaEvents = sortrows(时间表(日期,类型,值))
extremaEvents = 9×2时间表
日期 类型 价值
1 24 - 7 - 1972 0.0030
2 1975年- 7月31日 0.0028
3. 25 - 6 - 1977 0.0029
4 22日- 1月- 1987 0.0013
5 10 - 11月- 1993 0.0023
6 11月20 - - 2003 0.0003
7 10 - 2月- 2008 0.0009
8 25 - 7 - 2010 0.0007
9 28日- 12月- 2015 0.0012
持有
isPeak =(极端事件。类型= =“高峰”);
情节(extremaEvents.Date (isPeak) extremaEvents.Value (isPeak),“y ^”“MarkerFaceColor”“y”);
情节(extremaEvents.Date (~ isPeak) extremaEvents.Value (~ isPeak),“青年志愿”“MarkerFaceColor”“y”);
持有
上面是单个日子的多余LOD的图表,并有一个粗略跟踪年平均值的平滑版本。它说的是,一年中平均的白昼长度已经变得大约等于,甚至 略小于 这意味着UTC-UT1的差异实际上是持平或缩小的,而不是增加的。前面显示的累积过量LOD图表明它已经趋于平稳。这就解释了为什么自2016年12月以来没有增加闰秒。缩小的UTC-UT1增加了A 闰秒最终可能是需要的。还没有,请耐心等待,但如果这种趋势继续下去,可以想象,过多的LOD可能会在几年内累积到接近-1的值,并出现负闰秒 则根据UTC的定义是必要的。这是以前从未发生过的,并且可能会引起很多软件的大量混淆(MATLAB不会混淆)。消极的闰秒可能就是这样做的动机 最后放弃闰秒 .然而,许多专家认为,这种趋势不会持续下去,这只是一个周期的一部分。
这幅图里有很多内容。有多余的LOD(本质上是太阳日长),它的斜率(日长变化的速度有多快),以及它的积分(UTC和UT1之间的差)。

与时间赛跑

从峰值到低谷,多余的LOD正在减少,这意味着地球的自转正在加速。我将这些递减周期的开始和结束存储在一个表中。最后一个峰值没有对应的低谷;我将使用数据中的最后一个时间作为最后一个间隔的结束时间。
StartDate = IERSdata.Date(峰值);
结束日期= [IERSdata.Date(槽);IERSdata.Date(结束)];
减去事件=表(开始日期,结束日期)
decreasingEvents = 5×2表
StartDate可以 EndDate
1 24 - 7 - 1972 1975年- 7月31日
2 25 - 6 - 1977 22日- 1月- 1987
3. 10 - 11月- 1993 11月20 - - 2003
4 10 - 2月- 2008 25 - 7 - 2010
5 28日- 12月- 2015 06 - 4月- 2021
持有
startEnd =[递减事件。StartDate可以decreasingEvents。EndDate];
h = fill(startEnd(:,[1 2 2 1]),[-.]002 -。002.005.005],“红色”“FaceAlpha”2,“线型”“没有”);
持有
最后,我将计算每个间隔期间过量LOD的平均减少量(以每年“每天过量LOD的秒数”为单位),并将该信息添加到表中的间隔事件中。
dTime =递减事件。EndDate- decreasingEvents.StartDate;
dExcess = IERSdata.SmoothedELOD(reducingevents . enddate) - IERSdata.SmoothedELOD(reducingevents . startdate);
decreasingEvents。annualavgreduce = seconds(dExcess ./ years(dTime))
decreasingEvents = 5×3表
StartDate可以 EndDate AnnualAvgDecrease
1 24 - 7 - 1972 1975年- 7月31日 -8.6974 e-05秒
2 25 - 6 - 1977 22日- 1月- 1987 -0.00016705秒
3. 10 - 11月- 1993 11月20 - - 2003 -0.0001993秒
4 10 - 2月- 2008 25 - 7 - 2010 -5.6298 e-05秒
5 28日- 12月- 2015 06 - 4月- 2021 -0.00030701秒
换句话说,一整年的平均太阳日在过去几年中以每年0.3毫秒的速度减少。正如我们上面看到的,它目前非常接近86400s。如果最近的下降持续下去,从现在开始的一年后,平均太阳日比86400s少0.3ms,这意味着负的过剩LOD将以每年约十分之一秒的速度积累。但按照这个速度,需要近十年的时间才能出现负闰秒。

没有时间浪费了

地球会停止加速吗?从长远来看,是的。你不可能战胜潮汐相互作用。
近年来,相对而言,在很短的时间内,原始过剩LOD显著为负。这些只是短暂的波动,但在那些时期,地球自转了1.5毫秒 而不是86400 SI秒。
情节(IERSdata.Date IERSdata.ExcessLOD,“b -”);
ylabel (“超额LOD”);
持有
线(IERSdata。日期([1 end]),[0 0],“颜色”“k”“线型””:“
持有
ylabel (“超额LOD”);
我将确定任何一天的过量LOD小于-0.5ms的年份,找出每一年的最小过量LOD,并创建一个包含这些最小值的表。
negLOD = IERSdata(IERSdata。ExcessLOD< milliseconds(-0.5),“ExcessLOD”);
groupsummary (negLOD“日期”“年”“最小值”
ans = 11×3表
year_Date GroupCount min_ExcessLOD
1 2001 8 -0.0007064秒
2 2002 13 -0.0007436秒
3. 2003 31 -0.0009769秒
4 2004 29 -0.0010672秒
5 2005 25 -0.0010809秒
6 2007 4 -0.0006192秒
7 2010 12 -0.000784秒
8 2018 6 -0.0006457秒
9 2019 19 -0.0009571秒
10 2020 88 -0.0014663秒
11 2021 11 -0.0007026秒
那么明年会发生什么呢?之所以只提前六个月宣布闰秒,是因为很难准确预测地球转速的变化。但是,尽管 日前发表的 并不能预测未来很长一段时间内的过量LOD,但它们可以 发布DUT1的预测 有一年出去约会。
T =可读的(“https://datacenter.iers.org/data/latestVersion/finals.all.iau2000.txt”);
T = T (:,[3 10 12]);
t.Properties.VariableNames = ["日期" "DUT1" "LOD"];
t. date =日期时间(t。目前为止,“ConvertFrom”“mjd”);
t.DUT1 =秒(t.DUT1);
t.LOD =毫秒(t.LOD);
T = T (T。日期> =“1 - 1月- 2020”:)
t = 908×3表
日期 DUT1 LOD
1 01 - 1月- 2020 -0.17716秒 0.0004379秒
2 02 - 1月- 2020 -0.17763秒 0.0004915秒
3. 03 - 1月- 2020 -0.17811秒 0.0004743秒
4 04 -简- 2020 -0.17859秒 0.0004946秒
5 05 - 1月- 2020 -0.17908秒 0.0004482秒
6 06 - 1月- 2020 -0.17943秒 0.0002298秒
7 07 - 1月- 2020 -0.17955秒 4.14 e-05秒
8 08 - 1月- 2020年 -0.17953秒 -8.77 e-05秒
9 09 - 1月- 2020 -0.17938秒 -0.0002198秒
10 10 - 1月- 2020 -0.17912秒 -0.0002715秒
11 11 - 1月- 2020 -0.17888秒 -0.0001832秒
12 12 - 1月- 2020 -0.17879秒 2.29 e-05秒
13 13 - 1月- 2020 -0.17894秒 0.0002682秒
14 14 - 1月- 2020 -0.17933秒 0.0005166秒
预测=(日期>=“1 - 4月- 2021);
情节(t.Date(~)预测,t.DUT1预测(~),“b”...
t.Date(预测),t.DUT1(预测),“r”);
ylabel (“DUT1”
传奇([“实际”“预测”],“位置”“西北”
请记住,DUT1是UT1-UTC,当DUT1小于大约-0.5s时,会添加正闰秒。从过去一年的斜率可以看出,太阳日长度先是大于,然后小于,最后大致等于86400s。明年的斜率预测太阳日长度将小于,然后等于,然后大于86400s。如果正确,那就意味着……一段时间内没有闰秒,不管是积极的还是消极的!在接下来的几年里继续关注这一点将是很有趣的。如果1970年以来的模式继续下去,就会出现负闰秒 但问题是什么时候:过剩的LOD会很快作为最近十年周期的一部分再次出现,还是会像2015年以来那样继续下降?

阅读时间

我对时间系统和闰秒的了解很大一部分来自于 这些非常全面的页面 作者:Steve Allen UCO /利克天文台 .其他关于UTC和闰秒历史的非常可读的描述是由 莱文 麦卡锡 奎因 , 纳尔逊等人

是时候总结了

闰秒出现问题的主要方式是计算跨越一个或多个闰秒的两个时间瞬间之间的差值。在MATLAB中,你可以选择在你的 datetime 计算方法 UTCLeapSeconds 创建日期时间:
datetime(2016, 12日31,36岁,0,0,“时区”“UTCLeapSeconds”) - datetime(2016,12,31,12,0,0,“时区”“UTCLeapSeconds”
ans =持续时间24:00:01
但我们发现,大多数人看到明显随机的额外秒会感到恼火,所以没有明确要求闰秒计算,你会看到
datetime(2016, 12日31,36岁,0,0,“时区”“UTC”) - datetime(2016,12,31,12,0,0,“时区”“UTC”
ans =持续时间24:00:00
其他时区也一样,比如 美国/ New_York .(我们还发现,许多人觉得一年有两次多出一个小时或少了一个小时很烦人,所以你也必须选择使用某个时区的日光节约时间。)
你的工作在计算中需要闰秒吗?你用过吗 UTCLeapSeconds ?或者你只是觉得他们很烦人?应该放弃闰秒吗?与UTC不同的是,GPS时间系统不插入闰秒,因此它是一个“统一”时间系统,但每次在UTC中插入一个闰秒,就会与UTC偏离1秒。你的工作使用GPS时间吗?让我们知道你的MATLAB时间经验(痛苦?) 在这里
|

评论

如欲留言,请点击在这里登录您的MathWorks帐户或创建一个新帐户。