利用太阳影子的定位技术.doc

2014高教社杯全国大学生数学建模竞赛
承 诺 书
我们仔细阅读了《全国大学生数学建模竞赛章程》和《全国大学生数学建模竞赛参赛规则》（以下简称为“竞赛章程和参赛规则”，可从全国大学生数学建模竞赛网站下载）。
我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。
我们知道，抄袭别人的成果是违反竞赛章程和参赛规则的，如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。
我们郑重承诺，严格遵守竞赛章程和参赛规则，以保证竞赛的公正、公平性。如有违反竞赛章程和参赛规则的行为，我们将受到严肃处理。
我们授权全国大学生数学建模竞赛组委会，可将我们的论文以任何形式进行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。
我们参赛选择的题号是（从A/B/C/D中选择一项填写）： A
我们的报名参赛队号为（8位数字组成的编号）： 201510033013
所属学校（请填写完整的全名）： 空军勤务学院
参赛队员 (打印并签名) ：1.
2.
3.
指导教师或指导教师组负责人 (打印并签名)：
（论文纸质版与电子版中的以上信息必须一致，只是电子版中无需签名。以上内容请仔细核对，提交后将不再允许做任何修改。如填写错误，论文可能被取消评奖资格。）
日期： 2015 年 9 月 13日

赛区评阅编号（由赛区组委会评阅前进行编号）：
2014高教社杯全国大学生数学建模竞赛
编 号 专 用 页
赛区评阅编号（由赛区组委会评阅前进行编号）：
赛区评阅记录（可供赛区评阅时使用）：
评
阅
人

评
分

备
注

全国统一编号（由赛区组委会送交全国前编号）：
全国评阅编号（由全国组委会评阅前进行编号）：
20 / 22
.
太阳影子定位
摘要
本文针对如何通过分析视频中物体的太阳影子变化，太阳影子定位技术，确定视频的拍摄地点和拍摄日期的方法为背景，提出了四个问题，并且利用二元非线性模型，三角函数和最小二乘法，同时借助matlab，成功建立了通过太阳影子位置确定地理经纬度的模型。通过代入附件1、2、3的数据，并利用Photoshop获取附件4视频中截图的数据，求解出了各个问题中的测量点地理经纬度。
针对问题一，在已知经纬度，时间的条件下，我们可以得到太阳的高度角。对于地球上的某个地点，太阳高度角是指某地太阳光线与通过该地与地心相连的地表切线的夹角（如图1）。杆的高度乘以太阳高度角的余切，可以得到影子长度变化的数学模型。最后代入数据求得，利用MATLAB可以获得2015年10月22日北京时间9:00-15:00之间天安门广场3米高的直杆的太阳影子长度的变化曲线。（如图3）。
针对问题二，在测量时刻的日期，北京时间和影子的直角坐标已知的情况下，我们可以直接求得太阳赤纬角，方位角，影长。我们对附件1的数据处理之后，通过对影长和方位角用最小二乘法拟合，得到拟合函数。通过拟合函数求得正午时刻的影长，即函数图像最低点时的方位角和影长（如图4、5）。由正午时刻的太阳高度角计算的特殊性将公式化简为，代入数据求得纬度。再根据当地时刻在正午12:00时的北京时间，求得经度。
针对问题三，通过对第二问的求得的地理经纬度，将附件二某一组数据反带求得杆长L，利用影长和杆长求得高度角。再次利用第二问的方法，拟合影长和高度角（如图6、7），联立三个方程就可以解出纬度，赤纬角，以及时角t。利用时角t和求得当地时间，和北京时间的差求出经度；利用赤纬角求出测量时刻的日期。
针对问题四，在已知日期的情况下，对附件4视频中随机取17张截图（尽量时间间隔均匀），利用Photoshop的标尺工具测量截图中的杆长度和影子长度（如图8），且已知杆2米，由等比例对照得到不同时刻影长的数据（如表1），对数据进行拟合得出正午时刻的影长（如图9）。将拍摄日期，影长代入问题三的模型当中，对模型求解得到纬度。由当地时间和北京时间的差值求得经度。在未知日期的情况下，由视频截图可知此时为早晨，太阳在东侧，则摄影机位置大概是在杆的北侧，故由Photoshop测量可得太阳方位角大致是264.6°。由方程组可以求得纬度和赤纬角，再通过赤纬角求得日期，通过当地时间求得经度。
关键词：影子定位 二元非线性方程 拟合函数 等比例对照
一．问题重述
如何确定视频的拍摄地点和拍摄日期是视频数据分析的重要方面，太阳影子定位技术就是通过分析视频中物体的太阳影子变化，确定视频拍摄的地点和日期的一种方法。
1.建立影子长度变化的数学模型，分析影子长度关于各个参数的变化规律，并应用你们建立的模型画出2015年10月22日北京时间9:00-15:00之长度的变化曲线。间天安门广场（北纬39度54分26秒,东经116度23分29秒）3米高的直杆的太阳影子
2.根据某固定直杆在水平地面上的太阳影子顶点坐标数据，建立数学模型确定直杆所处的地点。将你们的模型应用于附件1的影子顶点坐标数据，给出若干个可能的地点。
3. 根据某固定直杆在水平地面上的太阳影子顶点坐标数据，建立数学模型确定直杆所处的地点和日期。将你们的模型分别应用于附件2和附件3的影子顶点坐标数据，给出若干个可能的地点与日期。
4．附件4为一根直杆在太阳下的影子变化的视频，并且已通过某种方式估计出直杆的高度为2米。请建立确定视频拍摄地点的数学模型，并应用你们的模型给出若干个可能的拍摄地点。
如果拍摄日期未知，你能否根据视频确定出拍摄地点与日期？
二．模型假设
1. 忽略海拔，地面平整度的影响。
2. 忽略天气因素和大气对太阳光的折射造成的影响。
3. 所研究的范围不能太大，认为所研究地面为一平面。
4. 确定杆为直杆，铅直直立，且与当地地面垂直。
5. 因为视频中杆与影子距摄像机距离大致相同，所以忽略透视效果造成的测量误差。
三．符号说明
符号
释义
δ
太阳的赤纬角，即太阳直射角纬度
h
太阳高度角，太阳光的入射方向和地平面之间的夹角
φ
地理纬度

地方时(时角)以表示
L
确定直杆长度

杆的影长
M
为两个地方的经度差
H
正午的太阳高度角
四．问题分析
太阳影子定位技术是根据不同地点，不同时刻，太阳在南回归线和北回归线之间运动到不同位置时从而来确定时间，利用物体的实际长度和影长之比算出太阳的高度差，从而确定其纬度，在此过程中近视的认为所研究的范围内为一平面，在问题一中，知道时间和经纬度，从而可以根据太阳的高度角得出影子长度；而问题二，为问题一的反求，利用MATLAB的数据拟合，因为地球在绕着太阳公转的同时还在自传，即为太阳的赤纬角，所以可认为影子的顶点的位置是以定杆位置为原点的函数，利用MATLAB的数据拟合可在误差允许的范围内拟合出这个函数的方程，从而确定极小致为正午十二点时的影长，利用太阳高度角的计算公式，得出纬度值，而经度则是根据和北京时差得出；问题三则是联立太阳高度角和太阳方位角公式得出赤纬角，从而得出时间，地点则与第二问的求法相似；问题四，已知拍摄日期时，在视频中利用Photoshop随机截图若干，测量数据，利用问题二三的方法计算时间地点；未知拍摄日期时，先利用截图测量太阳方位角和高度角联立方程组求解出拍摄日期，然后照前述的方法求出地理位置。
五．模型建立与求解
5.1. 问题一
针对第一问，在已知经纬度，时间的条件下，我们可以得到太阳的高度角，即对于地球上的某个地点太阳高度角是指某地太阳光线与通过该地与地心相连的地表切线的夹角。然后杆的高度除以太阳高度角的正切，得到影子长度随时间变化的函数方程。
日升日落，同一地点一天内太阳高度角是不断变化的。时角是以正午12点为0度开始算，每一小时为15度。即14点和10点分别为30度和-30度。日出日落时角度都为0。太阳高度角随着地方时和太阳的赤纬的变化而变化。太阳赤纬（与太阳直射点纬度相等）以δ表示，观测地地理纬度用φ表示（太阳赤纬与地理纬度都是北纬为正，南纬为负），地方时(时角)以表示，有太阳高度角的计算公式：
sin h=sin φ sin δ+cos φ cos δ cos

h
（图1）

（图2）
上述式子中都涉及太阳赤纬，太阳赤纬的算法如下:
由于太阳赤纬角在周年运动中任何时刻的具体值都是严格已知的，所以它(δ)也可以用表达式表述，即:
0.3723+23.2567sinθ+0.1149sin2θ-0.1712sin3θ-0.758cosθ+0.3656cos2θ+0.0201cos3θ
式中θ称日角，即 θ=2πT/365.2422
这里T又由两部分组成，即 T=N-N0
式中N为积日，所谓积日，就是日期在年内的顺序号，例如，1月1日其积日为1，平年12月31日的积日为365，闰年则为366，等等。
N0=79.6764+0.2422×(年份-1985)-INT〔(年份-1985)/4〕
(式中INT表示取整数部分，例如INT(3.25)=3)
根据已知条件为10月22日，则 t=295
则根据公式 θ=2t/365.2422 得到θ=5.0747
0.3723+23.2567sinθ+0.1149sin2θ-0.1712sin3θ-0.758cosθ+0.3656cos2θ+0.0201cos3θ=-22.0906
根据纬度计算出 φ=40.3333
sin h=sin φ sin δ+cos φ cos δ cos t=0.0482+0.8697 cos t
所以 sin h= 0.0482+0.8697 cos t
所以
= =
由杆长的关系
所以可得影长的变化模型为：
==
=
借助MATLAB绘制出2015年10月22日北京时间9:00-15:00之长度的变化曲线。间天安门广场（北纬39度54分26秒,东经116度23分29秒）3米高的直杆的太阳影子顶点变化曲线图：

（图3）
5.2．问题二
5.2.1 模型的建立
我们以杆与地面接触的点作为原点，穿过该点的纬线在该点处的水平切线作为X轴，与X轴垂直的地面切线作为Y轴建立直角坐标系。有附件1中的数据我们可以得出在某一时刻下的影子长度和此时的方位角。由影长，杆长可以求得太阳的高度角。我们由第一问可知，知道太阳赤纬角和地理纬度，加上当地时间就可以求得太阳高度角。现在我们知道了影长，赤纬角（太阳直射点纬度），求杆长，地理纬度。所以我们将数据进行拟合估计出正午时影子长度以及正午时此刻的方位角，通过方程组联立求解杆长，地理纬度和当地时间，结合当地时间和北京时间的差值求得经度。
由公式：
0.3723+23.2567sinθ+0.1149sin2θ-0.1712sin3θ-0.758cosθ+0.3656cos2θ+0.0201cos3θ
θ=2T/365.2422(T为日期)
可知，在2015年4月14日时，22.6571，太阳赤纬角sin δ=-0.6178。
此时由太阳高度角的计算公式：
sin h=sin φ sin δ+cos φ cos δ cos
日升日落时，同一地点一天内太阳高度角是不断变化的。日出日落角度都是0度，正午太阳高度角最大，正午时时角为0，以上计算公式可以简化为：

其中，H表示正午太阳高度角。
由两点和与差的三角函数公式，可得

因此，对于北半球而言，；对于南半球而言，
其中，为当地纬度，为太阳直射点的纬度，既赤纬角。为纬度差
如果与在同一半球，则纬度差为，如果与在异半球，则纬度差为，
因为杆长未知，所以高度角可以表示为：

故上式可化为

根据太阳方位角的计算公式：

其中A为正午时刻的方位角，将sinH代入可得：

故只要用拟合出来的正午方位角，已知赤纬角就可以确定当地纬度
将当地纬度反代到

计算出杆长L，由此得出高度角H，
其中太阳高度角可以表示为 ,当地时间可表示为+B,其中t为标准北京时间，A是当地与北京时间的时角差。
最后，得出的方程为=sin φ sin δ+cos φ cos δ cos（+B）;将杆长L与地理纬度代入，计算出时差角B
时角与时间的转化公式如下：
,
地球是自西向东自传的，东边比西边先看到太阳，东边的时间比西边的早。所以规定将全球划分24个时区，每个时区相差1个小时，也就是说1度相差4分钟，故可建立以下模型：

其中为待求地地方时，t为已知地地方时，M为两地之间的经度差。
若>t则待求地在已知地的东向，既，已知地经度+M。若<t 则待求地在已知地的西向，既已知地经度—M，

其中为待求地经度，为已知地经度。
5.2.2 模型的求解
现在将附件一的数据代入到该模型当中：
根据将附件一的数据拟合，得到影长拟合函数：

当函数图像最低点时，x =13.6，y =0.5476既正午的影长为0.5476;

（图4-1）

（图4-2）

（图4-3）
方位角拟合函数：

（图5-1）

（图5-2）

（图5-3）
当函数图像最低点时，x=13.6，y=1.8836,既tanA=1.8836,A=
代入模型：

得到
将代入

计算出杆长L=1.7919m,代入下面模型，取附件一中北京时间15:00时的影长为1.3533
=sin 5.3 sin 22.66+cos 5.3 cos 22.66 cos（+B）
计算得到B=

则M=
因为北京为东经所以待测地经度为东经
综上所示，直杆的地理经纬度为（）通过卫星地图确定在海南省。
5.3 问题三
5.3.1 模型的建立
基于第二问我们知道除了需要确定直杆的地理位置外，还需要确定此时所处的日期，也就是说直杆所处的赤纬角是未知的，杆长，地理纬度，当地时间以及赤纬角，结合当地时间和北京时间的差值求得经度，根据赤纬角来确定当地日期，
此时由太阳高度角的计算公式：
sin h=sin φ sin δ+cos φ cos δ cos
其中太阳高度角可以表示为 ,当地时间可表示为+B,其中t为标准北京时间，B是当地与北京时间的时角差。
最后，得出的方程为=sin φ sin δ+cos φ cos δ cos（+B）;
由第二问可以得出杆长L，所以由太阳方位角的计算公式：

sin h=sin sin δ+cos cos δ cos

其中高度角h可以由第二问的杆长以及此刻所具有的影长，可以由此刻的北京时间确定，联立三个方程就可以解出纬度，赤纬角，以及时角差B。
由公式：
0.3723+23.2567sinθ+0.1149sin2θ-0.1712sin3θ-0.758cosθ+0.3656cos2θ+0.0201cos3θ
θ=2T/365.2422(T为日期)
就可以计算出当地日期
时角与时间的转化公式如下：
,

经度计算方法与第二问相同，即：
地球是自西向东自传的，东边比西边先看到太阳，东边的时间比西边的早。所以规定将全球划分24个时区，每个时区相差1个小时，也就是说1度相差4分钟，故可建立以下模型：

其中为待求地地方时，t为已知地地方时，M为两地之间的经度差。
若>t则待求地在已知地的东向，既，已知地经度+M。若<t 则待求地在已知地的西向，既已知地经度—M，

其中为待求地经度，为已知地经度
5.3.2 模型的求解
首先将附件二的数据进行拟合：
拟合后的影长拟合函数为：

（图6-1）

（图6-2）

（图6-3）
当函数图像在最低点时，x=14.43 ,y=0.7920
拟合后的方位角拟合函数为;

（图7-1）

（图7-2）

（图7-3）
当函数图像在最低点时，x= 14.43 ,y= -14.8716，既tanA=-14.8716
取13:02时的影子长度为1..794将数据代入方程组：

sin h=sin sin δ+cos cos δ cos

解方程可得 查表得到大约是5月22日
经度计算得
综上所述地点为（）日期为5月22日
将模型应用到附件三当中同理可得，地点为（）日期为7月13日
5.4 问题四
5.4.1 已知拍摄日期的情况
针对第四问中的附件四影子与直杆的比例关系，利用Photoshop，可以测量出视频截图中影子的长度。根据等比例对照可得若干组影长的数据。

（图8）
可得随机截得的17张视频截图中的数据如下表：
时间
图上杆长
图上影长
实际影长
时间
图上杆长
图上影长
实际影长
8-54-50
5.16
6.06
2.348837
9-11-40
5.16
5.54
2.147287
8-58-04
5.16
5.98
2.317829
9-14-29
5.16
5.4
2.093023
8-59-31
5.16
5.97
2.313953
9-19-48
5.16
5.24
2.031008
9-00-29
5.16
5.93
2.29845
9-21-59
5.16
5.16
2
9-01-57
5.16
5.91
2.290698
9-23-18
5.16
5.12
1.984496
9-03-17
5.16
5.85
2.267442
9-29-58
5.16
4.85
1.879845
9-05-09
5.16
5.78
2.24031
9-32-32
5.16
4.81
1.864341
9-08-48
5.16
5.61
2.174419
9-33-58
5.16
4.76
1.844961

9-34-26
5.16
4.74
1.837209
（表1）
通过与时间的拟合，估计出正午的影长，根据第二问中，正午太阳高度角最大，正午时时角为0，以上计算公式可以简化为：

通过日期计算出赤纬角，H为正午时的高度角，
如此就可以计算出当地纬度
选择某一时刻的影长来确定高度角h，将赤纬角，纬度以及高度角代入下面公式
=sin φ sin δ+cos φ cos δ cos（+B）
计算出时角查B
时角与时间的转化公式如下：
,

待求地经度计算模型
地球是自西向东自传的，东边比西边先看到太阳，东边的时间比西边的早。所以规定将全球划分24个时区，每个时区相差1个小时，也就是说1度相差4分钟，故可建立以下模型：

其中为待求地地方时，t为已知地地方时，M为两地之间的经度差。
若>t则待求地在已知地的东向，既，已知地经度+M。若<t 则待求地在已知地的西向，既已知地经度—M，

其中为待求地经度，为已知地经度。
综上得到直杆所在地的经纬度。
经过拟合的影长函数是

（图9-1）

（图9-2）

（图9-3）
得到拟合函数最低点，与第二问同理将数据代入到方程

求出纬度为北纬40.28°，经度为东经111.43°
5.4.2 拍摄日期未知的情况
由视频中影长可以得出太阳方位角为，高度角计算与5.4.1同理

将视频中拟合计算得到的高度角，太阳方位角代入联立解方程得直杆的纬度和赤纬角，
由公式：0.3723+23.2567sinθ+0.1149sin2θ-0.1712sin3θ-0.758cosθ+0.3656cos2θ+0.0201cos3θ
θ=2T/365.2422(T为日期)
反求得日期，确定经度与5.4.1同理。
故在日期未知的情况下只要确定太阳方位角，与高度角方程联立就可以求出日期。
六．模型评价
模型的建立是在理想化程度下，忽略海拔，地面等因素的影响利用数据的拟合得出相关数据，模型的优点在于能够快速的确定出其大致方位，数据获取与处理较为方便简单；而缺点则是由于数据不够多，由数据拟合带来的地理误差和时间误差可能较大，同时没有对因大气折射率而造成的误差作后期修正。
七．参考文献
[1] 汪晓银，周保平，数学建模与数学实验（第二版），北京：科学出版式，2012.8
[2] 同济大学数学系，高等数学（第二版）上册，上海：同济大学出版社，2009.10
[3] 薛定宇，陈阳泉，高等数学问题的MATLAB求解，北京：清华大学出版社,2008