2015年度全国数学建模竞赛A题全国一等奖专业论文12.pdf

《2015年度全国数学建模竞赛A题全国一等奖专业论文12.pdf》由会员分享，可在线阅读，更多相关《2015年度全国数学建模竞赛A题全国一等奖专业论文12.pdf（24页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、 1 太阳影子定位研究模型太阳影子定位研究模型 摘要摘要 本文针对太阳影子定位问题,通过非线性最小二乘拟合模型、最优化问题、 CAD 软件制图、零点函数等建立数学模型，综合分析求解固定杆的经纬度定位和 视频拍摄地点及日期等问题。 针对问题一，利用地理学中的太阳高度角公式和几何学知识，建立影子长度 变化规律的数学模型， 分析影子长度的影响参数有太阳赤纬、 地理经纬度和时角， 考虑天安门广场当地时间和北京时间的时差， 天安门广场正午时刻的北京时间为 12.24 时，利用此模型画出 2015 年 10 月 22 日北京时间 9:00-15:00 之间天安门 广场 3 米高的直杆的太阳影子长度的变化曲

2、线。 针对问题二，预处理附录 1 中的坐标数据,计算出不同时间下的影子长度， 利用 lsqcurvefit 函数计算出当地时差，再使用经度公式计算得到当地的经度； 非线性拟合影子长度关于纬度、杆长和赤纬三个参数的曲线，将问题转化为最优 化问题，利用 lsqcurvefit 函数算出纬度、杆长和赤纬，位置在北纬 19.2858 度，东经 108.7215 度（海南省） 。 针对问题三， 经度和纬度的计算方法与问题二中的一致， 利用赤纬计算公式， 由赤纬值和零点函数求得日期， 得到附件二的测量地位置为新疆省维吾尔族自治 区，经纬度坐标是北纬 39.8951，经度 79.7490，测量日期为 201

3、5 年 6 月 20 日 或 2015 年 4 月 25 日；附件三测量地位置为湖北省，经度坐标是北纬 32.8472， 经度 110.2455，测量日期为 2015 年 10 月 23 日或 2015 年 2 月 6 日。 针对问题四，通过每两分钟截取一次视频，得到了 21 张杆子影子图片。使 用 CAD 制图软件分别测量图片中杆子影子的长度， 借助比例公式求出实际影子的 长度， 得到影子长度随时间变化数据；使用问题三中所建立的数学模型对数据进 行处理，算出拍摄点的经纬度，位置约在东经111.0135度，北纬42.0732度（内 蒙古） 。 针对问题五，由问题四计算得到的杆子的影长随时间变化

4、数据，使用问题二 中建立的数学模型，得到了两个拍摄地点经纬度，由于太阳直射点在春分-秋分 作周期运动，对于一个地理坐标可以得到两个日期。从而在拍摄日期未知时，仍 然可以根据视频确定出拍摄地点与日期。得出位置经纬度为北纬 41.9462 度，东 经 112.8705 度（内蒙古） ，拍摄日期为日期为 6 月 21 日或 4 月 25 日。 关键词：太阳高度角公式 最小二乘拟合 CAD 软件制图 零点函数 2 一、问题重述一、问题重述 1.1 1.1 问题背景问题背景 影子是由于物体遮住了光线的传播， 不能穿过不透明物体而形成的较暗区域， 就是我们常说的影子，它是一种光学现象。影子长度与季节、当日

5、时间、物体高 度和地理位置等有关， 如何确定视频的拍摄地点和拍摄日期是视频数据分析的重 要方面，太阳影子定位技术就是通过分析视频中物体的太阳影子变化，确定视频 拍摄的地点和日期的一种方法。 1.2 1.2 要解决的问题要解决的问题 1.建立影子长度变化的数学模型，分析影子长度关于各个参数的变化规律， 并应用建立的模型画出 2015 年10月22 日北京时间 9:00-15:00 之间天安门广场 （北纬 39 度 54 分 26 秒,39.907222 东经 116 度 23 分 29 秒）3 米高的直杆的太 阳影子长度的变化曲线。 2.根据某固定直杆在水平地面上的太阳影子顶点坐标数据， 建立数

6、学模型确 定直杆所处的地点。将模型应用于附件 1 的影子顶点坐标数据，给出若干个可能 的地点。 3.根据某固定直杆在水平地面上的太阳影子顶点坐标数据， 建立数学模型确 定直杆所处的地点和日期。 将模型分别应用于附件 2 和附件 3 的影子顶点坐标数 据，给出若干个可能的地点与日期。 4附件 4 为一根直杆在太阳下的影子变化的视频，并且已通过某种方式估 计出直杆的高度为 2 米。请建立确定视频拍摄地点的数学模型，并应用模型给出 若干个可能的拍摄地点。 5如果拍摄日期未知，根据视频确定出拍摄地点与日期。 二、问题分析二、问题分析 2.1 2.1 问题的总体分析问题的总体分析 本题主要通过影子长度的

7、变化和时间确定拍摄地点和日期。 建立影子长度随 时间的变化曲线，分析曲线最低点，确定时差，从而求出经度；根据影子长度和 时间的多组数据，拟合影长关于物体长度、赤纬和纬度的曲线，确定物体长度、 赤纬和纬度。而赤纬和日期有关系，通过赤纬可求得日期。通过实际和视频中的 影长和物体长度的比例相等的关系，计算出实际物体影长，继而按照上述分析确 定拍摄地点的经纬度和拍摄日期。 2.2 2.2 具体问题分析具体问题分析 （1） 针对问题一， 利用天文学知识将太阳光线的照射情况转化为几何知识， 建立影子长度变化规律模型，分析影子长度的影响因素，建立影子长度和北京时 间的函数关系，用 matlab 画出 201

8、5 年 10 月 22 日北京时间 9:00-15:00 之间天 安门广场 3 米高的直杆的太阳影子长度的变化曲线。 3 （2）针对问题二，根据附件一中影子长度和北京时间的数据，对影子长度 和时间进行拟合，计算出正午时刻的北京时间，算出经度，通过拟合影子长度关 于纬度的曲线，从而计算出纬度。 （3）针对问题三，问题三比问题二少一个已知量日期，通过非线性拟合影 子长度关于纬度和赤纬的曲线，计算出经纬度和赤纬，通过日期和赤纬的关系， 求出日期。 （4）针对问题四，问题四的难点在于实际影长的确定，通过 CAD 制图软件 计算出视频中影长和物体长度的比例， 利用实际中和视频中的影长和物体长度的 比例相

9、等的关系，实际物体长度已知，求出各时间下的影长。利用分析（3）中 相同方法，确定拍摄地点。 （5）针对问题五，在未知拍摄时间下，还是可以求出实际影长和时间，通 过非线性拟合影长和赤纬和经纬度，确定拍摄地点。 三、模型假设三、模型假设 1.假设模型中均已修正大气折射等因素。 2.假设量取视频中的长度误差很小，可以忽略。 3.忽略平闰年的影响。 4.忽略天气对影子长度的影响。 四、符号说明四、符号说明 符号 意义 1 t 日角 太阳赤纬 地理纬度 地理经度 t 时角 纬度差 h 太阳高度角 z 影子长度 d 距离 t , ,a b c 时差 系数 f 零点函数 4 五、模型的建立与求解五、模型的建

10、立与求解 5.1 5.1 问题一：问题一：画出影子长度变化曲线 5.1.1 5.1.1 模型建立模型建立 利用地理天文学知识，将太阳光照射地面的情况转化为几何图形，如图 1。 图 1 太阳高度角 对于地球上的某个地点P， 太阳高度角是指太阳光的入射方向和地平面之间 的夹角，记为h，地点的纬度和赤纬的差值记为，记太阳赤纬为，地点 的地理纬度为，地理经度为，根据物理中光线的平行性和几何学中平行、垂 直性质，有 2 h 推出太阳高度角 sinhcos cossinsin 考虑到当日具体时间的变化，同理可导出准确太阳高度角公式 sinhsin sincos coscos ( ) t 太阳高度角 arc

11、sin sinsincoscoscos ( )ht PP 太阳光线 P Q 5 其中t为北京时间（单位：h） ，( ) t表示当地时角，它与当地经度、时差 有关，时角以小时来计量（1HA = 15 度） ，当地正午时（12 时）的时角为 0，可 由以下公式计算出，当地时角公式： 120 15() 180 15 ( )2 360 t t 太阳赤纬可由公式获取。太阳赤纬计算公式为： 11 1111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=ED tt tttt 其中 1 t称日角， 2422.365

12、2 1 x t 这里x由两部分组成，即0NNx 式中 N 为积日，积日即为日期在年内的顺序号。 079.67640.2422(年份-1985)(年份-1985)/ 4NINT 式中INT表示进行取整运算。 物体成像示意图见图 2，h为太阳高度角，y为物体高度，z为影子长度。 图 2 物体成像示意图 根据三角函数求得影子长度 cothzy 将求得的太阳高度角公式代入上式，得出影子长度与北京时间的zt函数 5.1.2 5.1.2 模型求解模型求解 北京时间经度为 E120， 天安门广场经纬度为北纬 39 度 54 分 26 秒,东经 116 度 23 分 29 秒为 116 度 23 分 29 秒

13、，即纬度 39.907222 度、经度 116.39133 度。 cotarcsin sinsincoscoscos ( )zyt 6 利用太阳赤纬公式 11 1111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=ED tt tttt 计算出天安门广场 2015 年 10 月 22 日的太阳赤纬-10.8627， 程序见附录一。 利用时角公式 120 15() 180 15 ( )2 360 t t 计算出表达式 15(0.2406) 180 ( )2 360 t t 代入模型的函数得出 3cota

14、rcsinsin(-10.8627 )sin39.907222 15(0.2406)180 cos(-10.8627 )cos39.907222 cos2 360 z t 利用 matlab 画出 2015 年 10 月 22 日北京时间 9:00-15:00 之间天安门广场 3 米高的直杆的太阳影子长度的变化曲线，如图 3. 图 3 天安门广场杆子的影子长度随时间变化图 5.25.2 问题二：根据影子顶点坐标确定直杆地点问题二：根据影子顶点坐标确定直杆地点 5.2.15.2.1 模型建立模型建立 （1）数据预处理 zt 7 根据给出的时间和影子顶点坐标, ii x y，将时间全部转化为小时，

15、利用原点 到某点的距离公式 22 ii dxy，算出不同时间下的影子长度。 北京时间 x 坐标(米) y 坐标(米) 北京时间（时） 杆子影长 （米） a 点 b 分 i x i y /60ab 22 ii xy （2）经度计算，先使用经度公式 120 1512 m t 计算出经度，再用时角公式： （3）纬度计算，使用赤纬计算公式 111111 0.0201cos3+2+0.3656cos0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=EDtttttt 计算出 2015 年 4 月 18 日该地点的太阳赤纬为 10.6305 度。 根据已求的经度，

16、预处理北京时间，利用时角公式 120 15() 180 15 ( )2 360 t t 可计算出当地的时角。 利用非线性曲线拟合，拟合影长关于纬度和杆长的曲线 cotarcsin sinsincoscoscos ( )zyt 则将问题转化为最优化问题： 2 21 1 min( , ) cotarcsin sinsincoscoscos ( ) iiii i F yytz 使用 lsqcurvefit 函数，求出纬度和杆长y，程序见附录五、附录六。 （4）综上，可计算出直杆地点的经纬度，实现定位。 5.2.2 5.2.2 模型求解模型求解 （1）对附件 1 数据进行预处理，计算出时间和影长，如表

17、 1. 表 1 预处理影长和时间数据表格 北京时间 x 坐标 （米） y 坐标(米) 北京时间(时) 杆子的影长 120 15() 180 15 ( )2 360 t t 8 14:42 1.0365 0.4973 14.7 1.149626 14:45 1.0699 0.5029 14.75 1.182199 14:48 1.1038 0.5085 14.8 1.215297 14:51 1.1383 0.5142 14.85 1.249051 14:54 1.1732 0.5198 14.9 1.283195 14:57 1.2087 0.5255 14.95 1.317993 15:00

18、 1.2448 0.5311 15 1.353364 15:03 1.2815 0.5368 15.05 1.389387 15:06 1.3189 0.5426 15.1 1.426153 15:09 1.3568 0.5483 15.15 1.463400 15:12 1.3955 0.5541 15.2 1.501482 15:15 1.4349 0.5598 15.25 1.540232 15:18 1.4751 0.5657 15.3 1.579853 15:21 1.516 0.5715 15.35 1.620145 15:24 1.5577 0.5774 15.4 1.66127

19、1 15:27 1.6003 0.5833 15.45 1.703291 15:30 1.6438 0.5892 15.5 1.746206 15:33 1.6882 0.5952 15.55 1.790051 15:36 1.7337 0.6013 15.6 1.835014 15:39 1.7801 0.6074 15.65 1.880875 15:42 1.8277 0.6135 15.7 1.927918 （2）利用解最优化约束问题最优解的方法，使用 lsqcurvefit 函数，解出 了纬度 ，时差 m t 和杆长 y ，分别为 m 19.2858 ,2.0336,t0.6755y

20、程序见附录五、附录六。 （3）经度计算 程序见附录三、附录四。 （4）综上，杆子位置在(19.2858, 108.7215)NE 北纬 19.2858 度，东经 108.7215 度为，即位置在：海南省。 5.2.3 5.2.3 模型检验模型检验 利用问题二的模型算出的经纬度、 赤纬和杆长， 使用问题一的影子长度变化 模型， 画出影子长度变化曲线，在曲线上标出附件一所给点的坐标来验证拟合的 优度， 影子长度变化的曲线和附件一所有点的坐标标注见图 2， 程序见附录十一。 0 120 1512108.7215 m t 9 图 2 模型检验 分析图 2 可以观察到曲线拟合效果较好。 5.35.3 问

21、题三：根据影子顶点坐标确立直杆日期与地点模型问题三：根据影子顶点坐标确立直杆日期与地点模型 5.35.31 1 模型建立模型建立 （1）数据预处理 根据给出的时间和影子顶点坐标, ii x y，将时间全部转化为以小时为单位， 利用原点到某点的距离公式 22 ii dxy，算出不同时间下的影子长度。 北京时间 x 坐标(米) y 坐标(米) 北京时间（时） 杆子影长 （米） a 点 b 分 i x i y /60ab 22 ii xy （2）经度计算 据此可算出时差(12) m tt ，从而计算出经度差15Et，经度公式 单位为：度 （3）纬度、赤纬和杆长的计算 拟合影长关于纬度、杆长和赤纬的曲

22、线 120 1512 m t 10 cotarcsin sinsincoscoscos ( )zyt 则将问题转化为最优化问题： 2 21 1 min( , , ,) cotarcsin sinsincoscoscos ( ) miiiiii i F ytytz 使用 lsqcurvefit 函数，求出纬度、赤度和杆长y。 (4)日期的确定 太阳赤纬计算公式： 用赤纬与日期的关系，由赤纬反求日期。利用函数零点求出（积日） 。 令赤纬-时间函数 111111 11 ( )ED(t ) 0.3723+23.2567sin +0.1149sin20.1712sin30.758cos +0.3656c

23、os2 +0.0201cos3 0 f ttttt tt 代入赤纬 1 ED(t )，通过求函数 1 f t的零点，求出 1 t确定日期。 （5）综上，可计算出直杆地点的经纬度和日期。 5.35.32 2 模型求解模型求解 5.3.2.1 5.3.2.1 对附件二求解对附件二求解 （1）对附件 2 数据进行预处理，计算出时间和影长，如表 3. 表 3 预处理影长和时间数据表格（附件 2） 北京时间 x 坐标(米) y 坐标(米) 北京时间（时） 杆子的影长 12:41 -1.2352 0.173 12.683 1.247256 12:44 -1.2081 0.189 12.733 1.2227

24、95 12:47 -1.1813 0.2048 12.783 1.198921 12:50 -1.1546 0.2203 12.833 1.175429 12:53 -1.1281 0.2356 12.883 1.152440 12:56 -1.1018 0.2505 12.933 1.129917 12:59 -1.0756 0.2653 12.983 1.107835 13:02 -1.0496 0.2798 13.033 1.086254 13:05 -1.0237 0.294 13.083 1.065081 13:08 -0.998 0.308 13.133 1.044446 13:1

25、1 -0.9724 0.3218 13.183 1.024264 13:14 -0.947 0.3354 13.233 1.004640 13:17 -0.9217 0.3488 13.283 0.985491 11 11111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=)ED(t tt tttt 1 t 11 13:20 -0.8965 0.3619 13.333 0.966790 13:23 -0.8714 0.3748 13.383 0.948585 13:26 -0.8464 0.3876

26、 13.433 0.930928 13:29 -0.8215 0.4001 13.483 0.913752 13:32 -0.7967 0.4124 13.533 0.897109 13:35 -0.7719 0.4246 13.583 0.880974 13:38 -0.7473 0.4366 13.633 0.865492 13:41 -0.7227 0.4484 13.683 0.850504 （2）利用解最优化约束问题最优解的方法，使用 lsqcurvefit 函数，程序 运行结果：（见附录七、八） 运行结果 x = 0.6414 2.0008 0.3545 2.6834 则有： 纬度

27、为 2 arcsin(0.6414)/236039.8951 经度为 杆长 2 2.0008y 米 赤纬为 2 arcsin(0.3545)/236020.7640 经纬位置在东经 71.960 度，北纬 30.9168 度，大概位置。 （4）将已求的赤纬值代入赤纬-日期零点函数，得出日期。 111111 11 ( )ED(t ) 0.3723+23.2567sin +0.1149sin20.1712sin30.758cos +0.3656cos2 +0.0201cos3 0 f ttttt tt 首先假设年份为 2015 年， 0NNx . 0N=79.9424 代入零点公式，运行下列程序：

28、 fzero(6.7151-(0.3723+23.2567*sin(2*pi*x/365.2422-1.3752)+0.1149* sin(2*(2*pi*x/365.2422-1.3752)-0.1712*sin(3*(2*pi*x/365.2422-1.3752) )-0.758*cos(2*pi*x/365.2422-1.3752)+0.3656*cos(2*(2*pi*x/365.2422-1.3 752)+0.0201*cos(3*(2*pi*x/365.2422-1.3752),180)% 求零点值 当赋初值为 180 时，即当天的积日为 201，日期为 6 月 20 日； 当赋初

29、值为 150 时，即积日为 145，日期为 4 月 25 日。 120 151279.7490 m t 12 （5）综上得到附件二给出的地点为东经 79.7490 度，北纬 39.8951 度，位 置为新疆，日期为 6 月 20 日或日期为 4 月 25 日。 5.3.2.2 5.3.2.2 对附件三求解对附件三求解 （1）对附件 3 数据进行预处理，计算出时间和影长，如表 4. 表 4 预处理影长和时间数据表格（附件 3） 北京时间 x 坐标(米) y 坐标(米 北京时间（时） 杆子的影长 13:09 1.1637 3.336 13.15 3.533142184 13:12 1.2212 3

30、.3299 13.2 3.546768029 13:15 1.2791 3.3242 13.25 3.561797643 13:18 1.3373 3.3188 13.3 3.578100715 13:21 1.396 3.3137 13.35 3.595750783 13:24 1.4552 3.3091 13.4 3.61493428 13:27 1.5148 3.3048 13.45 3.635425983 13:30 1.575 3.3007 13.5 3.657218272 13:33 1.6357 3.2971 13.55 3.680541115 13:36 1.697 3.293

31、7 13.6 3.705167836 13:39 1.7589 3.2907 13.65 3.731278025 13:42 1.8215 3.2881 13.7 3.758917911 13:45 1.8848 3.2859 13.75 3.788087888 13:48 1.9488 3.284 13.8 3.818701015 13:51 2.0136 3.2824 13.85 3.850809619 13:54 2.0792 3.2813 13.9 3.88458522 13:57 2.1457 3.2805 13.95 3.919911828 14:00 2.2131 3.2801

32、14 3.956875992 14:03 2.2815 3.2801 14.05 3.99553479 14:06 2.3508 3.2804 14.1 4.035750835 14:09 2.4213 3.2812 14.15 4.077863059 （2）利用解最优化约束问题最优解的方法，使用 lsqcurvefit 函数，程序运 行结果：通过运行 M 文件和主程序得到结果为： x = 0.5424 3.0356 -0.2750 0.6503 则有： 纬度为 2 arcsin(0.5424)/236032.8472 经度为 13 杆长 2 3.0356y 米 赤纬为 2 arcsin( 0

33、.2750)/236015.9620 经纬度坐标位置为湖北省十堰市。 （4）将已求的赤纬值代入赤纬-日期零点函数，得出日期。 111111 11 ( )ED(t )0.3723+23.2567sin +0.1149sin20.1712sin30.758cos +0.3656cos2 +0.0201cos3 0 f ttttt tt 首先假设年份为 2015 年 fzero(-15.9620-(0.3723+23.2567*sin(2*pi*x/365.2422-1.3752)+0.114 9*sin(2*(2*pi*x/365.2422-1.3752)-0.1712*sin(3*(2*pi*x

34、/365.2422-1.375 2)-0.758*cos(2*pi*x/365.2422-1.3752)+0.3656*cos(2*(2*pi*x/365.2422-1 .3752)+0.0201*cos(3*(2*pi*x/365.2422-1.3752),180)% 求零点值 运行结果： 当赋初值为 50 时，结果为 37，日期为 2 月 6 日， 当赋初值为 200 时，结果为 311，日期为 10 月 23 日。 （5）得到附件三的地点为东经 110.2455 度，北纬 32.8472 度，日期为 2 月 6 日或 10 月 23 日。 5.4 5.4 问题四：根据视频问题四：根据视频

35、预测若干个可能的拍摄地点预测若干个可能的拍摄地点 5.45.41 1 模型建立模型建立 （1）赤纬计算 根据问题一中的模型，赤纬与日期的关系式 11 1111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=ED tt tttt 由视频中的日期可以计算出拍摄地点的赤纬。 (2)实际影长 利用CAD制图软件测量得出从视频中截取的21张图像中的杆长和影长的长 度。 记视频中的杆长和影长分别为,y l ,实际的杆长和影长分别为, y l，其中 2y ,根据放射性质，有 120 1512110.2455 m t

36、 yy ll 14 则实际影长2 yll l yy . （3）纬度的计算 拟合影长关于纬度曲线 2cotarcsin sinsincoscoscos ( )zt 则将问题转化为最优化问题： 2 21 1 min( ,)2cotarcsin sinsincoscoscos ( ) miiii i Fttz 使用 lsqcurvefit 函数，求出纬度和时差 m t . . 5.45.42 2 模型求解模型求解 （1）赤纬计算公式 11 1111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=ED tt

37、tttt , 2015 年 7 月 13 日积日为 184，代入赤纬计算公式，计算出 赤纬 4 23.0049。 （2）利用 CAD 制图软件，计算视频中的杆长和影长，见图 3。 图 3 CAD 测量的影长和杆长 使用实际影长2 yll l yy 公式，算出各时间点的实际影长，见表 5. 表 5 利用杆长与影长的比例算出各时间点的实际影长 15 时间 图片中的影长 图片的杆长 实际杆长 实际影长 8:54:10 109.19 93.63 2 2.332372 8:56:10 108.2 93.63 2 2.311225 8:58:10 106.23 93.63 2 2.269145 9:00:

38、10 104.91 93.63 2 2.240948 9:02:10 103.93 93.63 2 2.220015 9:04:10 102.61 93.63 2 2.191819 9:06:10 101.95 93.63 2 2.177721 9:08:10 99.65 93.63 2 2.128591 9:10:10 99.32 93.63 2 2.121542 9:12:10 97.68 93.63 2 2.086511 9:14:10 96.36 93.63 2 2.058315 9:16:10 95.7 93.63 2 2.044217 9:18:10 93.73 93.63 2 2

39、.002136 9:20:10 93.07 93.63 2 1.988038 9:22:10 91.43 93.63 2 1.953007 9:24:10 90.44 93.63 2 1.931859 9:26:10 88.47 93.63 2 1.889779 9:28:10 87.48 93.63 2 1.868632 9:30:10 86.82 93.63 2 1.854534 9:32:10 83.54 93.63 2 1.784471 9:34:10 81.86 93.63 2 1.748585 (3）利用解最优化约束问题最优解的方法，使用 lsqcurvefit 函数，算得 x =

40、 0.6701 0.5991 纬度为 4 arcsin(0.6701)/236042.0723 经度为： 综上，拍摄地点的经纬度为东经111.0135度，北纬42.0732度（内蒙古） 。 5.5 5.5 问题五：在时间未知下，根据视频问题五：在时间未知下，根据视频推断若干个可能的拍摄地点推断若干个可能的拍摄地点 在视频拍摄日期未知的情况下，我们可以从视频显示屏上知道拍摄的时间， 即时角。直杆高度没变，仍是 2 米。可以使用问题四中方法截取视频中图片，测 量得到图片中直杆影子长度与直杆长度，根据公式计算出实际的直杆影子长度。 再用最小二乘拟合求出影子长度最小值时间， 即正午时间， 计算出观测点

41、的经度。 再使用最优化约束问题求出纬度。 下面由最小二乘拟合求出纬度和赤纬的最优解： 120 1512111.0135 m t 16 利用问题二的模型得出: 当赋初值为（40，30，0.6）时，得出纬度为 41.4602 度，经度为 112.8705 度，赤纬为 20.7444 度，利用赤纬与积日的关系得出日期， 11 1111 0.0201cos3+0.3656cos2+ 0.758cos0.1712sin30.1149sin2+23.2567sin+0.3723=ED tt tttt 日期为 6 月 21 日或 4 月 25 日。 六模型评价六模型评价 6.16.1 模型缺点模型缺点 1.

42、第二、三、四、五问中，本文使用了非线性最小二乘拟合求出了时差、经 度和纬度，无法计算或者估计经度误差，无法确定拟合精度是否精确。 2.第四问中，由于 matlab 版本过低，是 2009a 的版本，读取视频函数无法 使用，例如 mmread,VideoReader 函数。因此，本文采用了对附件视频进行截屏 的方式来获得每个时刻 t 变动的图像，这里 t 取的值为 2 分钟，共取到了 21 组 数据。这样截取的缺点是截取得到的图片大小不一定一样，而使用 matlab 处理 得到的帧的大小的相同，这样会测量更多次数的杆长和影长顶点坐标，使得过程 更加复杂。 3、由模型算出第五问视频拍摄时间为 6

43、月 21 日，而由第四问视频可知真正 拍摄时间为 7 月 13 日，因此模型存在一定误差。 6.26.2 模型优点模型优点 1.第一问中，在利用时差计算太阳高度角时，进行了时差修正，根据天安门 广场的经度计算出天安门广场和北京时间的时间差， 得到了天安门广场的正午时 间为 12 点 26 分，使得天安门广场直杆影子随时间变化曲线更准确。 2.第四问中，无法利用 matlab 函数编程来求得图中的影子和影子的长度， 且杆子影子是立体三维空间图， 而在图中得到的影子顶点和杆子顶点以杆子和影 子顶点交点坐标是平面坐标，运用欧氏距离计算出长度会造成测量误差。我们使 用了 CAD 制图软件，建立三维立体

44、空间来测量杆长、影长，误差更小。 七、参考文献七、参考文献 1王国安、米洪涛，太阳高度角和日出日落时刻太阳方位角一年变化范围的计 算，气象与环境科学，30（增刊），161-164,2007 年。 2刘魁敏，计算机绘图，北京：机械工业出版社，2008 年。 3赵静、但琦，数学建模与数学实验，北京：高等教育出版社，2003 年。 4郑鹏飞、林大均、刘小洋、吴志庭，基于影子轨迹线反求采光效果的技术研 究，华东理工大学学报,36，459-460，2010 年。 5肖智勇、刘宇翔，一种新的纬度测量方法，大学物理，29,51-53,2010 年。 17 八八附录附录 附录一： %计算出天安门广场的太阳赤纬

45、 n=295;%10 月 22 日在 2015 年内的顺序号号 no=79.6764+0.2422*(2015-1985)-7;%7 为（年份-1985）/4 取整后得到的值 x=n-no; t=2*pi*x/365.2422; ED=0.3723+23.2567*sin(t)+0.1149*sin(2*t)-0.1712*sin(3*t)-0.758*cos(t)+0.3656*co s(2*t)+0.0201*cos(3*t) 附录二： t=9:0.01:15; h=sin(39.907222/360*2*pi)*sin(-10.8627/360*2*pi)+cos(39.907222/3

46、60*2*pi)*cos(-1 0.8627/360*2*pi)*cos(180-15*(t-0.2406)./360*2*pi)%计算太阳高度角的正弦值 h1=asin(h)%求出太阳高度角 y=3./(tan(h1)%计算影子长度 plot(t,y)%画出影子长度变化曲线 xlabel(时间);ylabel(影子长度); title(影子长度变化曲线) 附录三：%第二问求经度 function f=ff(x,tdata) f=2*(1./(x(1)*(0.3908)+(1-(x(1)2)(1/2)*(0.9205)*cos(180-15.*(tdata-x(2)./360.*2.*pi).

47、2- 1).(1/2);%其中 x(1)为太阳高度角的正弦值 附录四： tdata=8.903 8.936 8.969 9.002 9.035 9.068 9.101 9.134 9.167 9.2 9.233 9.266 9.299 9.332 9.365 9.398 9.431 9.464 9.497 9.53 18 9.563 ydata=2.332372103 2.311225035 2.269144505 2.240948414 2.220014952 2.191818861 2.177720816 2.128591263 2.121542241 2.086510734 2.058314643 2.044216597 2.002136067 1.988038022 1.953006515 1.931859447 1.889778917 1.868631849 1.854533803 1.784470789 1.748584855 x0=39/360*2*pi,0.5; x=lsqcurvefit(ff,x0,tdata,ydata) 附录五： %第二问求纬度的 m 文件 %建立 ff.m 文件 %建立影长关于纬度和杆长的函数 function f=ff(x,tdata) f=