基于穷举法的太阳影子定位研究

祝长华，黄文莉，邱洪玲，郑 治

（韶关学院 数学与统计学院，广东 韶关 512005）

基于穷举法的太阳影子定位研究

祝长华，黄文莉，邱洪玲，郑 治

（韶关学院 数学与统计学院，广东 韶关 512005）

主要研究通过物体影子长度随时间的变化确定物体所处的地点和日期。通过Matlab编程运用穷举法穷举最优经纬度等信息，进而查找相应的地点和日期。首先，建立影子长度关于杆长、经度、纬度、日期和时间的函数模型，并分析各参数变化引起影子长度改变的规律。其次，建立由影子顶点坐标数据和时间确定物体的地点及日期的模型。最后，建立由物体影子随时间变化的视频确定拍摄地点和日期的模型。

太阳影子定位技术；太阳高度角；赤纬角；穷举法

1 问题提出

当今社会相机的应用非常普遍。在大数据时代，图片、视频格式数据含有海量信息。通过视频确定地理位置信息的能力和技术具有很大的应用空间，比如GPS定位系统、安全监控预测及导航系统等。由图片或视频信息确定拍摄地点的经纬度是目前计算机视觉领域的研究热点之一［1-3］，具有很深远的意义。太阳下物体影子在生活中随处可见［4-6］，特别是在户外运动时，如何通过物体影子确定当地地理位置是一个非常重要的技术。太阳影子定位技术就是通过分析视频中物体影子的变化，从而得到拍摄地的地点和日期的一种方法。现提出如下问题，希望通过太阳影子定位技术得以解决。

问题1建立影子长度变化的数学模型，分析影子长度关于各个参数的变化规律。并利用模型画出2015年10月22日北京时间9:00-15:00之间天安门广场（北纬39°54'26'',东经116°23'29''）3m高直杆太阳影子长度的变化曲线。

问题2根据某固定直杆在水平地面上的太阳影子顶点坐标数据、时间及拍摄日期，建立数学模型确定直杆所处的地点。并将模型应用于附件1［7］（2015年全国大学生数学建模竞赛a题的附件1，下同）中影子顶点坐标数据，给出若干个可能的地点。

问题3根据某固定直杆在水平地面上的太阳影子顶点坐标数据及时间，建立数学模型确定直杆所处的地点和日期。将模型应用于附件2［7］（2015年全国大学生数学建模竞赛a题的附件2，下同）、附件3［7］（2015年全国大学生数学建模竞赛a题的附件3，下同）的影子顶点坐标数据，给出若干个可能的地点与日期。

问题4根据一直杆在太阳下影子变化的视频［7］（2015年全国大学生数学建模竞赛a题的附件4，下同），建立确定视频拍摄地点的数学模型，并给出若干个可能的拍摄地点。如果拍摄日期未知，能否根据视频确定出拍摄地点与日期？

2 模型假设与符号系统

模型假设：

（1）假设地轴是固定的，地球近似为圆形，自转为匀速率；

（2）假设一年按365天计算，不考虑平年与闰年之分；

（3）假设不考虑太阳直射的影响。

符号系统如表1所示。

表1 符号系统

3 问题1模型的建立与求解

3.1 问题分析与参数选取

首先要确定影响影子长度的参数。在此，通过分析影响影子长度的各个因素，最终选定以杆长、地理位置（经纬度）、时间三个要素为参数，其中时间主要针对日期（季节）、当日时间两个要素进行细分，从而确定各个参数与影子长度的直接关系和间接关系，如图1所示。

3.2 名词解析

（1）赤纬角［8］：地球赤道平面与太阳、地球中心连线的夹角，赤纬角以年为周期，范围为±23°26'。

（2）太阳高度角：太阳光的入射方向与地平面间的夹角，范围为［0°，90°］。

（3）太阳方位角：太阳光线在地平面上的投影线与地平面正南方向所夹的角，正南方向为0°，顺时针方向为正角，表示太阳位于下午的范围，逆时针方向为负值，表示太阳位于上午的范围，取值范围为±180°。

图1 参数选取示意图

（4）时角［8］：是按照地球自转一周（24小时）相当于太阳在天球图上绕天轴一周即360°原理，不同的时角可表示在一天里不同时间的太阳位置，以其所在时圈的角度表示，即一小时相当于15°。并规定以太阳在观测点正南向，即当地时间12点的时角为0度，上午的时角为负值，下午的时角为正值。

3.3 模型建立

设H为直杆长度，L为影子长度，φ为观测点纬度，h为太阳高度角，δ为赤纬角，时角为ω，h'为正午太阳高度角，n为年积日，t地为地方时间，t北为北京时间，s地为观测地经度，t'为观测地正午时的北京时间，f为太阳方位角。太阳高度角投影示意图如图2所示，影子长度和杆长关系式为

图2 太阳高度角投影示意图

正午太阳高度角的计算公式为

假设1年按365天计算，以1月1日的年积日为1，则12月31日的年积日为365，设n为对应日期的年积日，表2为春分、秋分、夏至、冬至年积日表。

表2 春分、秋分、夏至、冬至年积日表

赤纬角的计算公式为［8］

时角的计算公式为

当地时间与北京时间的转换为s地∈［-180°，180°］，取负值表示在西经。将式（5）代入式（4）得

太阳高度角与方位角公式［9］为

综上，联立公式即得到所建模型为

其中φ∈［-90°，90°］，n∈［1，365］，H∈（0，+∞），t北∈［0，24］，s地∈［-180°，180°］。

3.4 影子长度关于各参数变化规律分析

当其他自变量不变，当日时间从日出到日落变动时，影子由长变短，再由短变长；在同一时间的不同日期，从北半球春季开始，在北回归线以南赤道以北，影子由长变短，再由短变长，夏至后再由长变短，再由短变长直至冬至；在北回归线以北，影子长度一直减小直至夏至，然后再一直增大直至冬至。南半球以此类推。在日期时间确定时，太阳直射某个纬度时，则直射点的纬度影子最短，该纬度以北或以南影子逐渐加长。某地达到正午时分，则该地经度以东以西影子都由短变长［10］。

3.5 模型的求解

将模型应用于2015年10月22日北京天安门广场时间9:00-15:00，3m高直杆影子长度，用Matlab进行样条插值得到影子随时间变化的曲线图。

对影长分析知：处于正午12:00时直杆的影子最短，为3.85m，以12:00为中心点，步长为0.5h，直杆的影子随着步长的递增而递增，在9:00-15:00这段时间内的两端点处影子达到最长，为6.97m，变化曲线如图3所示。

图3 直杆影子长度变化曲线

从图3的影长变化曲线可以看出，正午12:00的影长最短，以12:00为中轴，两边分别呈递增状态，是一条以12:00为对称轴的抛物线。

4 问题2模型的建立与求解

4.1 问题分析

本问题是基于问题1所建立的模型进行求解，确定各个参数（经度、纬度、杆长、日期、具体的时间），并将模型应用于附件1数据。已给出日期为4月18日，对应的年积日为第108天。考虑运用穷举法对所有的经纬度进行穷举，将附件1中数据依次代入问题1的公式，并将符合条件的某杆长（大于0）作为固定杆长。再按照公式计算每一经纬度及附件中每一时间点所对应的理论影长，与附件中计算得到的实际影长做相应的误差绝对值和。找到误差绝对值和最小的经纬度，查出所处的地理位置。

4.2 模型建立

分析附件1，根据每个时刻的影子顶点坐标计算影子长度，观察其变化规律，如图4所示.

图4 影子长度变化曲线

从图4可知，在北京时间14:42-15:42这个时间段内，影子的长度呈递增趋势，且增幅增大。由一天中正午影子长度最短可知，当地已过正午时分（12时）。所以不妨假设附件内影子最短时为当地正午时间，当北京时间为14：42（14.7时）时，计算经度的波动范围为

所以，可以拟定经度W在120°±40.5°间取值，纬度N在±90°间取值，在每一经度范围内对每一纬度穷举。

由式（5）将北京时间转换为地方时间，同时根据式（9），已知影长逆推计算不同纬度的杆长。将附件1中数据依次代入上面方程，当求出的杆长H大于0时，记该杆长为H0，并取杆长为该值（经计算第11个点符合），继续根据上述公式求解不同经纬度和时间所对应的影长L'，再根据附件1计算出的实际影子长L、利用穷举法求解出的不同经纬度的理论影子长L'，计算两者差的绝对值之和，求出差的绝对值之和最小的一组或几组，其对应的经纬度就是与实际值比较吻合的经纬度，在地图上找出其所处的地点。

4.3 模型的求解

通过Matlab软件，根据上述穷举法的思想，计算出一系列参数，如表3所示。

表3 附件1对应的最优参数

图5是根据Matlab的运行结果得到附件1中21个点的最优参数下做出的实际影子长度与理论影子长度的对比图。

图5 最优参数下，附件1理论影长与实际影长的对比图

从图5中可以看出，最优参数下（108.6°E，19.3°N），实际影长与理论影长的误差很小，几乎可以忽略不计。通过谷歌地图、百度地图、高德地图等查找最优参数（108.6°E，19.3°N）对应的地理位置，可以求解出该直杆可能的地点为海南省东方市附近。

5 问题3模型的建立与求解

5.1 问题分析

问题3与问题2的区别是多了一个未知参数——日期，所以在问题2的基础上多考虑一个参数即可，因此本问题继续沿用问题2的解题思路，穷举经纬度、日期的年积日数。将附件中的数据依次代入公式找到符合条件的某杆长作为定值，再将每一确定的经纬度、日期及给出的时间数据代入公式求得理论杆子影长。再将它与数据中求出的实际影子长度做差取绝对值，得出最小的绝对值之和，找到对应的经纬度、日期，从而确认直杆可能的地点与日期。

5.2 模型建立

分析附件2、附件3，根据每个时刻的影子顶点坐标计算影子长度，观察其变化规律发现：附件2的影子随时间越来越短，附件3的影子随时间越来越长。

对于附件2：在北京时间12.41-13:41这个时间段内，影子的长度呈递减趋势。由一天中正午影子长度最短可知，当地还未过正午时分（12时），所以不妨假设附件2内影子最短时为当地正午时间。北京时间为13:41（13.68时），若以此时刻作为接近当地正午的时间，计算出与北京最小的经度差β；北京时间为12:41（12.68时），根据图5可知影子长度最长，若假设此时为日出的开始时间，则可以计算出与北京最大的经度差γ，即

所以可拟定经度W在［120°-80.16°，120°-25.2°］，即W应在［39.4°，84.8°］间穷举取值，纬度N在± 90°间取值，在每一经度范围内对每一纬度穷举。

对于附件3：在北京时间13:09-14:09这个时间段内，影子的长度呈递增趋势。由一天中正午影子长度最短可知，当地已过正午时分（12时），所以不妨假设附件3内影子最短时为当地正午时间。北京时间为13:09（13.15时），若以此时刻作为接近当地正午的时间，计算出与北京最小的经度差β；北京时间为14:09（14.15时），根据图5可知影子长度最长，若假设此时为日落时间，则可以计算出与北京最大的经度差γ，即

所以可拟定经度W在［120°-17.5°，120°+57.5°］，即W应在［102.5°，177.5°］间穷举取值，纬度N在± 90°间取值，在每一经度范围内对每一纬度穷举。

由式（5）将北京时间转换为地方时间，再根据式（9），已知影长逆推计算不同纬度的杆长。将附件2、3中数据以此代入方程，求出的杆长H大于0时，记该杆长为H0，并取杆长为该值（经计算第11个点符合条件，所求杆长作为杆长的定值），继续根据上述公式及附件中数据求解不同经纬度和时间所对应的影子长L'。将由附件中计算出的实际影子长L与利用穷举法求解出的不同经纬度、日期所对应的理论影长L'，做差积的绝对值，并计算所有绝对值之和，求出差的绝对值之和最小的一组或几组，其对应的经纬度及日期就是与实际值比较吻合的经纬度，在地图上找出其所处的地点。

5.3 模型的求解

5.3.1 最优参数的求解

对模型的求解主要通过Matlab软件，根据上述穷举法的思想计算出下列参数，如表4、5所示。

表4 附件2对应的最优参数解

表5 附件3对应的最优参数解

通过画出理论影长和实际影长的对比图，可得最优参数下实际影长与理论影长的误差很小，几乎可以忽略不计。通过谷歌地图、百度地图、高德地图等查找附件2求到的经纬度（79.6°E，39.2°N）对应的地理位置，可以求解出该直杆可能的地点为新疆维吾尔自治区喀什地区巴楚县附近，日期大约为7月24日，如表6所示。附件3求到的经纬度（110°E，34.5°N）对应的地理位置，可以求解出该直杆可能的地点为陕西省渭南市附近，日期大约为10月22日，如表7所示。

5.3.2 在误差范围内其他可能的地点

在误差允许的范围内，还得到其他解。表6和表7分别给出了附件2和附件3两批数据的其他可能解，其中理想度由低到高是指实际影长误差绝对和由大到小。

表6 附件2其他可能的地点

表7 附件3其他可能的地点

6 问题4模型的建立与求解

6.1 问题分析

问题4主要是根据附件4一根直杆在太阳下影子变化的视频，建立数学模型找出可能的拍摄地点。通过视频可以获取的信息是视频的拍摄日期为7月13日的8：54-9：34，仔细观看视频，视频中出现了夏季的特征：行人穿着短袖和裙子，由此可推断出在此时太阳直射点在北半球。因为要确定视频的拍摄地点，所以通过Matlab对视频处理，对约40min的视频以2 min为步长测量出每一帧画面中的影子长度，得到21组数据，然后对这些数据进行误差修订，最终得到较精确的一组数据。本问题的第1小问是已知视频的拍摄时间，故可套用问题2的模型求解，用穷举法求解出最符合条件的几个经纬度，确认拍摄地点；第2小问没有给出具体的拍摄时间，让求出可能的拍摄时间和地点，故可套用问题3的模型求解出最符合条件的几个经纬度和拍摄时间。

6.2 数据的获取

通过Matlab读取视频，把视频转换成图片序列，对40min的视频以每间隔2 min得到一张图片，再对图片进行降噪处理，并测量出影子的长度。由于对静止图片进行测量时，所得影子长度并非真正的影长，有一个误差角度（不会太大），因此需要对测量值进行修正，得出影子的实际长度。表8是以14°角进行误差处理后得到的数据。

表8 在附件4中截取不同时间对应的影长

6.3 模型建立

获取到数据之后，由于问题4主要包括2个小问，第1小问是已知视屏拍摄时间求可能的拍摄地点，这里可通过问题2建立的模型得到；第2小问是单通过视频求解出可能的拍摄地点和拍摄时间，这里可通过问题3建立的模型得到。

6.4 模型的求解

6.4.1 已知时间求拍摄地：模型2的运用

通过Matlab软件，用表8中的数据替换附件1的数据，用相同的代码求解出的参数结果如表9所示。

表9 附件4第1问对应的最优参数

根据得到的经纬度结果查找谷歌地图、百度地图等，可得拍摄点的位置在湖北省襄阳市保康县南石线附近。

6.4.2 未知时间求拍摄地：模型3的运用

通过Matlab软件，用表8中的数据替换附件2或3的数据，用同样的代码求解出的参数结果如表10所示。

表10 附件4第2问对应的最优参数

根据得到的经纬度结果查找谷歌地图、百度地图等，可得拍摄点的位置在湖北省襄阳市保康县252省道附近，日期是7月6日。

7 模型的评价与推广

7.1 模型的评价

优点：1）问题1充分考虑影响影长变化的各个因素，通过空间几何知识，以中间变量赤纬角、太阳高度角、方位角、时角，建立影长关于自变量为杆长、经纬度、时间、日期之间的关系，得到的模型较精准，具有很好的推广作用。2）问题2至问题4的模型通过流程图解释算法通俗易懂，可以自己设定精度，按照需要可以将影长所对的地理位置和时间极为精确地输出。3）方法简单，穷举可能的经纬度即可输出想要的位置信息。4）耐心设置误差范围，便不会忽略考虑一切位置。

缺点：1）设定精度较高时，用电脑运行程序极度费时。2）输出结果包含一切地理位置，甚至不可能的地理信息（如南北极、海洋地区、山脉地区），需要人工排除。

7.2 模型的推广

问题1所建立的模型涉及的范围很广，如对球面天文学的研究、地区时差的研究、住房选址等都具有推广和借鉴意义。问题2至问题4是用Matlab编程穷举不同的经纬度得出的结果。在很多其他领域都可以应用。

［1］HAYS J,EFROS A A.IM2GPS:estimating geographic information from a single image［C］//Proceedings of the Conference on Computer Vision and Pattern Recognition.Washington DC:IEEE Computer Society,2008:1-8.

［2］JACOBS N,SATKIN S,ROMAN N,et al.Geolocating static cameras［C］//Proceedings of the 11th International Conference on Computer Vision.Washington DC:IEEE Computer Society,2007:1-6.

［3］JUNEJOI N,FOROOSH H.Estimating geo-temporal location of stationary cameras using shadow trajectories［C］//Proceedings of the Conference European on Computer Vision.Berlin:Springer,2008:318-331

［4］ANTONE M,BOSSE M.Calibration of outdoor cameras from cast shadows［C］//Proceedings of the Conference on Systems,Man and Cybernetics.Washington DC:IEEE Computer Society,2004:3040-3045.

［5］BOUGUETJ Y,PERONA P.3D photography on your desk［C］//Proceedings of the Sixth International Conference on Computer Vision.Washington DC:IEEE Computer Society,1998:43-50.

［6］CAO X,SHAH M.Camera calibration and light source estimation from images with shadows［C］//Proceedings of the Conference on Computer Vision and Pattern Recognition.Washington DC:IEEE Computer Society,2005:918-923.

［7］全国大学生数学竞赛组委会.2015年高教社杯全国大学生数学建模竞赛赛题［EB/OL］.［2015-09-12］.http://www.mcm. edu.cn/html_cn/block/8579f5fce999cdc896f78bca5d4f8237.html.

［8］吉林省建筑设计院.建筑日照设计［M］.北京:中国建筑工业出版社,1979:17-80.

［9］丹容A.球面天文学和天体力学引论［M］.北京:科学出版社,1980:5-20

［10］武琳.基于太阳阴影轨迹的经纬度估计技术研究［D］.天津:天津大学,2010.

【责任编辑：王桂珍foshanwgzh@163.com】

The study on the orientation of the sun shadow based on exhaustive method

ZHU Chang-hua,HUANG Wen-li,QIU Hong-ling,ZHENG Zhi
（Institute of Mathematics and Statistics,Shaoguan University,Shaoguan 512005,China）

In this paper,it mainly studies through the object shadow length changes over time to determine the location and data.Through MATLAB programming with exhaustive method,exhaustively the optimal information such as latitude and longitude is enumerated,and then the appropriate location and data are found.Firstly,it establishes the shadow length function model of bar length,the longitude,latitude,date and time,and analyses the shadow length changing laws as the various parameters changing.Secondly,by the shadow vertex coordinate data and time it sets up the model to determine the location of the object and the date.Finally,by video object shadow changes over time it establishes the model to determine the location and date.

the sun shadow positioning technology;the Angle of the sun;the declination Angle;exhaustive method

P128.1

A

1008-0171（2016）06-0022-09

2016-07-01

韶关学院科研项目（S201501018）；韶关学院教改项目（SYJY20151669）；韶关市科技创新基金资助项目（20157212）

祝长华（1983-），女，湖南汝城人，韶关学院讲师。