太阳影子定位-全国大学生数学建模竞赛论文

1 太阳影子定位 摘 要 本文利用直杆影子与太阳高度角的关系 将直杆影子长度变化的模型转化为地球自 转和公转模型 使用非线性最小二乘拟合的方法 求出直杆的地理经度 纬度和日期等 参数 解决了用直杆的太阳影子来确定方位与日期的问题 对于问题一 直杆影子长度变化的模型就是直杆影子长度与其经度 纬度和日期的 函数关系 由于地球的自转和公转导致经度 纬度和日期的变化 因此将影子长度变化 的模型转化为地球的自转和公转模型 随着地球的公转 太阳在地球上的直射点发生改 变 即太阳赤纬 赤纬角 随着日期的变化而变化 于是建立太阳赤纬 赤纬角 与日 期的公转模型 地球自转时 经度在不断变化 时角也随之改变 因此可建立太阳高度 角与时角之间的自转模型 利用 MATLAB 软件对模型中的经度 纬度 日期和杆长等参数 赋值 就可得不同时间的影子长度和影子长度变化的模型 对于问题二 根据不同时刻的太阳顶点影子坐标 确定该地区的经度和纬度 经分 析可知 问题二的模型就是问题一的反向推导 反向推导的模型是超定方程组 需要利 用非线性最小二乘法拟合 求解参数直杆高度 经度和纬度 确定直杆拍摄的地点可能 为 地理纬度 N19 3355 地理经度 E116 2058 对于问题三 第三问中的测量日期未知 与问题二相比 多了一个未知参数 沿用 问题二的方法 利用非线性最小二乘法拟合 求解未知参数直杆高度 经度 纬度和日 期 在一个回归年 太阳直射点在南北回归线之间运动 同一地点对应两个日期的影子 长度相等 利用 MATLAB 软件分别对附件 2 3 求解 得到 附件 1 E115 7201 N20 7628 2015 8 20 和 2015 10 27 附件 2 E116 228 N32 847 2015 8 13 和 2015 11 1 对于问题四 首先利用米尺直接测量 将视频信息转化为数据信息 即不同时间的 影子长度已知 拍摄日期已知时 确定拍摄地点的模型与问题二中的一致 拍摄日期未 知时 确定拍摄地点的模型与问题三中的一致 利用MATLAB软件求解拍摄的时间和地点 拍摄日期已知 E118 03 N40 04 2015 7 12 拍摄日期未知 E118 03 N40 04 2015 7 12 和 2015 12 5 该模型利用非线性最小二乘拟合函数 对于小残量的问题 收敛速度很快 在查询资料之后 使用计算机对问题四中的视频进行 减小了人为测量的误差 增强了数据的准确性 关键关键词词 太阳赤纬 太阳高度角 时角 超定方程组 非线性最小二乘拟合 2 一 问题的提出 如何确定视频的拍摄地点和拍摄日期是视频数据分析的重要方面 太阳影子定位技 术就是通过分析视频中物体的太阳影子变化 确定视频拍摄的地点和日期的一种方法 1 建立影子长度变化的数学模型 分析影子长度关于各个参数的变化规律 并应用 你们建立的模型画出 2015 年 10 月 22 日北京时间 9 00 15 00 之间天安门广场 北纬 39 度 54 分 26 秒 东经 116 度 23 分 29 秒 3 米高的直杆的太阳影子长度的变化曲线 2 根据某固定直杆在水平地面上的太阳影子顶点坐标数据 建立数学模型确定直杆 所处的地点 将你们的模型应用于附件 1 的影子顶点坐标数据 给出若干个可能的地点 3 根据某固定直杆在水平地面上的太阳影子顶点坐标数据 建立数学模型确定直杆 所处的地点和日期 将你们的模型分别应用于附件 2 和附件 3 的影子顶点坐标数据 给 出若干个可能的地点与日期 4 附件 4 为一根直杆在太阳下的影子变化的视频 并且已通过某种方式估计出直杆 的高度为 2 米 请建立确定视频拍摄地点的数学模型 并应用你们的模型给出若干个可 能的拍摄地点 如果拍摄日期未知 你能否根据视频确定出拍摄地点与日期 二 问题分析 2 1 问题一的分析 要建立直杆影子长度变化的数学模型 首先要考虑影子长度与哪些参数有关 经分 析可知 影子长度与太阳高度角有关 太阳高度角与地区的经度 纬度 太阳赤纬 赤 纬角 和一天内不同时间有关 因此 要建立影子长度模型的关键就是分析其随地区经 度 纬度 太阳赤纬和一天内不同时间的变化规律 而这些参数与地球的公转和自转有 关 于是 可将影子长度模型转化为地球公转和自转模型 地球围绕太阳公转时 太阳对 地球的直射点发生改变 即太阳赤纬 2 随地球公转而变化 建立地球公转模型 就可计算 任意一天的太阳赤纬 地球自转时 随着太阳直射点在纬度上移动 同一地点在不同时 间的太阳高度角不同 一天内太阳高度角随着地方时的变化而变化 建立地球自转模型 可计算太阳高度角的变化 最后根据太阳高度角 1 的定义 经过简单的三角变换 就可计 算出影子的长度 得出影子长度变化的模型 2 2 问题二的分析 根据直杆在水平地面上的太阳影子顶点的坐标数据和测量日期 建立模型确定直杆 所处的位置 直杆所处的位置 即该地区的经度和纬度 直杆在水平地面上的太阳影子 顶点的坐标数据已知 原点在直杆的底端 根据这两个条件就可计算出不同时间的影子 长度 问题一中已知直杆所在地的经度 纬度 日期和杆长求解不同时间的影子长度 问 题二已知日期和不同时刻的影长 求解直杆所在地的经度和纬度 如果把问题一视为正 向求解 那么问题二的模型是问题一的反向模型 即问题二的模型就是将问题一的模型 3 反向推导 最后利用非线性最小二乘法对附件1中的数据进行曲线拟合 给出参数的值 就可求出直杆所在的经度和纬度 由问题二中 给出若干个可能的地点 可知可能有多个地点的影长随时间的变化与 附件1数据相同 于是可根据地理知识理论的找出有几个地点的影长随时间的变化是相 同的 结合理论给出的和实际算出的 可求出可能的地点 2 3 问题三的分析 根据直杆在水平地面上的太阳影子顶点坐标 建立模型确定直杆所处的地点和日期 问题三与问题二相比 多了一个参数 即测量日期 问题二中利用非线性最小二乘拟合 计算出参数的具体值 问题三可继续沿用问题二的模型 拟合出函数中经度 纬度等参 数的数据 再根据拟合出的参数推导出直杆的测量日期 2 4 问题四的分析 根据直杆在太阳影子下的视频和直杆的长度 建立确立视频拍摄地点的数学模型 要计算经度和纬度 最关键的问题要将视频转化为影子长度的信息 将视频转化为数据信 息时 用计算机转化 数据更准确 但由于其步骤比较困难而且比较繁琐 所以用最简单 最直接的办法 用米尺直接测量每隔一定时间影子的长度 这种方法虽然存在一定的误 差 但简单有效 将视频转化为数据信息之后 就能得到不同时间的影子长度 拍摄视频 的日期已知时 可利用问题二中的非线性最小二乘法的模型 拟合出曲线的经度参数和 纬度参数 拍摄视频的日期未知时 将模型转化为问题三中模型求解 三 符号说明 符 号 含 义 太阳高度角 太阳赤纬 赤纬角 地理纬度 地方时 时角 北京的经度 h 其它地区的经度 i l i时刻的直杆影子长度 H 直杆杆长 T 北京时间 t 其它地区的时间 N 一年的日数 从 1 月 1 日算起 L 用米尺测量视频中直杆的长度 4 四 模型的假设 1 假设太阳 6 点从东方升起 18 点从西方降落 2 假设正午 12 点的时角为 0 3 假设地球是规则的球形 4 假设地球围绕太阳做匀速圆周运动 5 假设数据的测量都是在同一水平面上 6 假设附件中给出的数据准确无误 7 假设附件中的数据测量日期都是 2015 年 五 模型的建立与求解 5 问题一模型的建立与求解 5 1 1 问题一模型的建立 由分析可知 要建立影子长度变化的模型 就是建立影子长度与直杆所在经度 纬 度 太阳赤纬和一天内不同时间等参数的函数关系 而这些参数与地球的公转和自转有 关 所以将影长变化模型转化为地球的公转和自转模型 太阳赤纬 赤纬角 是地球赤道平面与太阳和地球中心连线的夹角 是由于地球环绕 太阳运行造成的现象 它随时间而变 因为地轴方向不同 所以赤纬角随地球在运行轨 道上的不同点 具有不同的数值 如图所示赤纬角为图中的h 图 1 太阳赤纬 赤纬角 根据太阳赤纬 赤纬角 的定义可知 在地球公转时 太阳赤纬随之改变 为太 阳赤纬 赤纬角 由于地球年并不总是 365 天 所以每年夏至 冬至和春分 秋分的 日期都不固定 赤纬角的方程式为 sin N 的单位为弧度 N是日数 从每年的 1 月 1 日开始 太阳高度角是指太阳光入射方向和地平面之间的夹角 即太阳光线与通过该地与地 心相连的地表切线的夹角 太阳高度角如图中所示 其中h表示太阳高度角 5 图 2 太阳高度角 太阳方位角 3 是指太阳直射光线在地平面上的投影线与地平面正南向所夹的角 通 常以南点 S 为 0 向西为正值 向东为负值 地方时 3 时角 HA 是天文学的名词 一个天体的时角被定义为该天体的赤经 RA 与当地的恒星时 LST 的差值 HA LST RA 因此 一个天体的时角表示该天体是 否通过了当地的子午圈 中天 其数值则表示了该天体与当地子午圈的角距离 并借 用时间的单位 以小时来计量 1HA 15 度 例如 一个天体的时角是 2 5HA 就表 示他已经在 2 5 个小时之前通过当地的子午圈 并且在当地子午圈的西方 37 5 的距离 上 负数则表示在多少小时之后将通过当地的子午圈 当然 当时角为 0 时的意思就是这 个天体就在当地的子午圈上 地球自转时 同一地点一天内 太阳高度角也是不断变化的 太阳高度角随地方时 和太阳的赤纬的变化而变化 太阳 6 点从东方升起 18 点从西方落下 正午 12 点的地方 时 时角 为 0 每相差一小时 相差 15 上午时角为负值 下午时角为正值 地球 自转 24 小时 旋转 360 所以可推出每旋转 1 地球自转的时间为 4 分钟 T为北京 时间 t为其它时区的时间 为北京的经度 h为其它地区的经度 将其它地区的时间 转化为北京时间 其表达式为 经度的两侧侧或两地区在 所求地区在北京的西 经度的同侧侧或两地区在 所求地区在北京的东 hT hT t 为地方时 时角 将时间转化为时角的公式为 t 由于太阳高度角随地方时 时角 和太阳的赤纬的变化而变化 为太阳高度角 为地理纬度 太阳高度角与地方时 时角 和太阳赤纬的关系式为 coscoscossinsinsin 根据太阳高度角的定义 太阳高度角是指太阳光入射方向和地平面之间的夹角 l为 直杆影子长度 H为直杆杆长 其关系式为 6 tan H l 综上可得任意地方时 时角 的影子长度的计算方法 sin N 精度的两侧侧或两地区在 所求地区在北京的西 经度的同侧侧或两地区在 所求地区在北京的东 hT hT t t coscoscossinsinsin tan H l coscoscossin sintan arcsin H l 以上就是计算任意地方时 时角 的影子长度的模型 5 1 2 问题一模型的求解 问题一中要求 2015 年 10 月 22 日北京时间 9 00 15 00 之间天安门广场 北纬 39 度 54 分 26 秒 东经 116 度 23 分 29 秒 3 米高的直杆的太阳影子的长度变化曲线 10 月 22 日对应的日数为 295 即 N 根据太阳赤纬 赤纬角 的公式可得 问题一中所给的时间是北京时间 不需要时区的转化 只要将北京时间转化为地方 时 时角 不同的时间对应的角度不同 时间为 9 00 15 00 地方时 时角 的范 围是 45 45 最后将时角代入太阳高度角的公式 再进行三角变换就可得到影子 的长度范围是 4 25m 7 8m 利用 MATLAB 软件 5 求解不同地方时 时角 下影子长度的变化曲线如图所示 3 2 10123 4 5 5 5 5 6 6 5 7 7 5 8 时间量 影子的长度 图 3 天安门广场直杆的太阳影子长度随时间的变化曲线图 7 规定中午 12 时在图中为 0 其余时间为 t 12 北京天安门在 2015 年 10 月 22 日 9 00 15 00 直杆的太阳影子长度先变短 中午 12 点时影子最短 然后又变长 在与 12 点对称的时刻 影子的长度是对称的 直杆的太阳影子长度在 9 00 时为 7 8m 在 12 00 时为 4 25m 在 15 00 时为 7 8m 5 2 问题二模型的建立与求解 5 2 1 问题二模型的建立 不同时刻的太阳顶点影子坐标已知 需要求出直杆所在的位置 即该地区的经度和 纬度 经分析可知 问题二就是问题一模型的反向推导 已知不同时刻太阳影子顶点坐 标和原点的位置 根据两点之间的距离 就可计算出不同时刻的影子长度 i l 为i时刻影 子的长度 其表达式为 iii yxl H为直杆杆长 为太阳高度角 根据太阳高度角的定义 其表达式 tan H l T为北京时间 t为其它地区的时间 为北京的经度 h为其它地区的经度 将北 京时间转化为其它地区时间 其表达式为 hTt 当计算出的区时为负数时 应加上 24 00 日期减一天 当计算出的区时为大于或等于 24 00 时 应减去 24 00 日期加一天 将其它地区的时间转化为地方时 时角 为地方时 时角 其表达式为 t 由于测量日期已知 根据地球公转模型 为太阳赤纬 赤纬角 太阳赤纬 赤 纬角 的表达式为 sin N 为地理纬度 根据自转模型 太阳高度角的表达式为 coscoscossinsinsin 经过对以上模型的分析可知 可确定参数的个数 coscoscossinsinsin tan t hTt H l 所以在该模型中直杆的长度H 地方角 时角 为 经度h和地理纬度 为参数 根据 附件一中的 21 组数据可知 方程的个数大于参数的个数 该方程组属于超定方程组 7 8 超定方程组一般是不存在解的矛盾方程组 因此引用非线性最小二乘法拟合方程 求解 超定方程组 非线性最小二乘拟合方法 4 的理论模型为 假设已知函数 xfy 的一组测量数据 ii yx mi 要求一个关于参数 Njaj 是非线性的函数 Ni aaax 10 对一组整数 m www 使得目标函数 m i NiiiN aaaxywaaaS 1 2 1010 达到最小 则称之为非线性最小二乘法问题 属于无约束条件的最优化问题 在该问题中 时间 i t 与 i x 相对应 因此 i t 为自变量 参数直杆的长度H 地方角 时 角 为 经度h和地理纬度 与 N aaa 相对应 以上参数就是利用非线性最小二 乘法所估计的参数的值 影子长度 i l 与函数 相对应 影子长度 i l 为因变量 所以问题二 中的非线性最小二乘法模型为 已知函数 tfy 的一组测量数据 2121 iyt ii 建立关于参数直杆的长度 H 地方时 时角 经度h和地理纬度 的非线性函数 hHtll iii hTt 180 1512 t 180 15 1244 T h coscoscossinsinsin coscoscossinarcsin sin coscoscossin sintan arcsin H l 即 coscoscossin sintan arcsin H l 使得目标函数 i iiii HtllHl min 建立无约束条件的最优化模型 综上可得测量影子地点的模型为 iii yxl tan H l 9 hTt t sin N coscoscossinsinsin coscoscossin sintan arcsin H l i iiii HtllHl min 5 2 问题二模型的求解 将模型应用于附件 1 中影子顶点坐标数据 将数据通过对已知函数进行拟合利用 MATLAB 软件拟合得到图形 14 614 81515 215 415 615 8 1 1 1 2 1 3 1 4 1 5 1 6 1 7 1 8 1 9 2 时间点 影 子 的 长 度 图 4 附件 1 影子的长度与时间的拟合图形 综合考虑我们利用 MATLAB 6 非线性最小二乘拟合法的工具箱 Isqucurvefit 函数进 行拟合参数 由于利用非线性最小二乘法拟合参数 要给出适当的初值 初始值准确的 与否决定求解工作成败的关键 首先从地理方面的专业知识人工划定参数的范围 在这 个范围内设定初始值 参数范围如下杆高 0 5 5m 直杆地理纬度 S90 N90 直杆的地理经度 E0 E180 参数大致的初始值为杆高 2m 直杆地理纬度 20 直杆的地理经度 E120 有了初始值之后 利用 MATLAB 软件求解得到拟合参数 但经过多次改变初始值 会 发现直杆的地理经度是发散的 然后我们将收敛的参数杆高和地理纬度当成已知量代入 参数模型 最终得到直杆的地理经度收敛在 E116 E116 72 经过多次比较取直杆 的地理经度为 E116 2058 在问题的分析中 我们想到在地球可能存在几个点 使得影子长度符合附件 1 的数 据 一开始分析在同一经度上 纬度关于太阳赤纬对称的两个点可能随着时间的变化影 10 长也相等 但在后面的系统分析中认为只有太阳直射在赤道上这一种情况才会有这样的 两点 因为如果太阳不直射在赤道上 那么两条不同纬线上的地点白昼长不相等 进而不 可能出现影长一直相等的两点 故抛弃了这种想法 使用非线性最小二乘法 利用 MATLAB 软件求解的直杆的高度为 2 0276m 直杆 的地理纬度 N19 3355 直杆的地理经度 E116 2058 5 3 问题三模型的建立与求解 5 3 1 问题三模型的建立 问题三中已知不同时间的太阳影子顶点坐标 与问题二相比 测量日期未知 即太 阳赤纬 赤纬角 为未知参数 沿用第二问的模型 利用非线性最小二乘法求解超定方 程 问题三中最小二乘法之前的反推模型与第二问中的大体一致 唯一的不同就是在模 型中增加了太阳赤纬 赤纬角 这一未知参数 利用最小二乘法解决超定方程时 其 拟合的函数由 Htll iii 转化为 Htll iii 用非线性最小二乘拟合函数 就可求出直杆的长度H 地方时 时角 经度h 地理纬度 和太阳赤纬 赤纬角 的参数值 设拟合函数为 coscoscossin sintan arcsin H l i iiii HtllHl min 利用 MATLAB 软件求解 求出直杆的长度H 太阳高度角 经度h 地理纬度 和太阳赤纬 赤纬角 的参数值 为太阳赤纬 赤纬角 赤纬角的方程式为 sin N 将太阳赤纬 赤纬角 的公式反推 就能得到直杆的测量日期的公式 其表达式为 arcsin N 根据上述模型可求解出一个日期 然后分析太阳赤纬 赤纬角 的变化情况 也就 是地球围绕太阳公转时 太阳光在地球上的直射点的轨迹变化 太阳直射点在南回归线 和北回归线之间运动 其运动轨迹如图所示 图 5 太阳直射点的运动轨迹 11 在上图中 上半部分表示北半球 下半部分表示南半球 上半部分的虚线表示北回 归线 下半部分的虚线表示南回归线 由图分析可知 春分和秋分前后太阳直射点在赤 道上 在夏至日前后 直射点在北回归线上 在冬至日前后 直射点在南回归线上 一个 回归年太阳直射点在赤道上 3 次 在南北回归线上各一次 在赤道与南北回归线之间的 任意一个位置 太阳直射点经过两次 如图中与赤道平行的直线与直射点的运动轨迹的 两个交点可分析得出 在赤道与南北回归线之间的每一点影子长度相同时对应两个不同 的日期 而且这两个不同的日期在同一个纬度上 根据太阳赤纬 赤纬角 的公式可计算出一个地点所对应的其中一个日期 然后根 据太阳直射点的位置分类求出另一个日期 1 太阳直射点在赤道和北回归线之间 根据太阳赤纬 赤纬角 的公式计算出的日期在夏至日之前 所求日数 夏至日对应的日数 夏至日对应的日数 已知日期对应的日数 根据太阳赤纬 赤纬角 的公式计算出的日期在夏至日之后 所求日数 夏至日对应的日数 已知日期对应的日数 夏至日对应的日数 N为与日数N相对应的日数 其模型为 已知日期在夏至之后 已知日期在夏至之前 NNN NNN N N 为夏至日所对应的日数 2 太阳直射点在赤道和南回归线之间 根据太阳赤纬 赤纬角 的公式计算出的日期在冬至日之前 所求日数 冬至日对应的日数 冬至日对应的日数 已知日期对应的日数 如果计算的日数大于 365 天 则所求日数 计算的日数 365 根据太阳赤纬 赤纬角 的公式计算出的日期在冬至日之后 所求日数 冬至日对应的日数 已知日期对应的日数 冬至日对应的日数 N为与日数N相对应的日数 其表达式为 已知日期在冬至之后 已知日期在冬至之前 NNN NNN N N 为夏至日所对应的日数 每一个地点对应两个日期 根据非线性最小二乘法求出在不同太阳直射点的条件下 的地点个数 然后分别求出每一地点所对应的日期 综上可得测量影子日期的计算方法 sin N arcsin N 12 已知日期在夏至之后 已知日期在夏至之前 NNN NNN N 已知日期在冬至之后 已知日期在冬至之前 NNN NNN N 5 3 2 问题三模型的求解 对于附件 2 与第二问类似 估计出最小二乘法参数范围 直杆的高度 0 5 5 米 直 杆地理纬 S90 N90 直杆的地理经度 E0 E180 太阳赤纬角 S23 47 N23 47 参数初始值直杆的高度 直杆的地理纬度 直杆的地理经度 太阳赤纬角分别为 2m 20 115 和 10 将模型应用于附件 2 和附件 3 中数据的影子顶点坐标数据 将数据通过对已知函数 进行拟合利用 MATLAB 软件拟合得到附件 2 的图形 12 81313 213 413 613 8 0 85 0 9 0 95 1 1 05 1 1 1 15 1 2 1 25 时间点 影子 的长 度 图 6 附件 2 影子的长度与时间的拟合图形 使用 MATLAB 软件经过多次拟合 直杆的高度 直杆的地理纬度收敛 而直杆地理 经度不收敛 通过附件 3 数据 将纬度值代入参数方程 多次求解 最终可以得到直 杆地理经度为 E115 7201 直杆的高度 2 0008m 直杆的地理纬度 N20 7628 直杆的地理经度 E115 7201 太阳赤纬 N7 4176 将拟合的参数用表格表示为 表 1 问题三附件 2 中求得的参数表 直杆高度 直杆地理经度 直杆地理纬度 太阳赤纬 参数范围 0 5 5 E0 E180 S90 N90 S23 45 N23 45 参数初始值 2 E120 N20 N10 拟合参数 2 0008 E115 7201 N20 7628 N7 4176 由赤纬角与日期的时间的换算公式可以推导 N 265 33 取整为 N 265 大约为 2015 13 年 8 月 20 日 由上面模型的建立地理的知识我们知道 太阳在直杆 N20 7628 和 E115 7201 的 地方一年会直射两次 则推导另外的一个日期为 2015 年 10 月 27 日 对于附件 3 的计算与第二问类似 估计出最小二乘法参数范围 直杆的高度 0 5 5m 直杆地理纬度 S90 N90 直杆的地理经度 E0 E180 太阳赤纬角 S23 47 N23 47 参数初始值直杆的高度 直杆的地理纬度 直杆的地理经度 太阳赤纬角分别为 3m N30 E115 和 S15 将数据通过对已知函数进行拟合 利用 MATLAB 软件拟合得到附件 3 的图形 1313 213 413 613 81414 2 3 5 3 6 3 7 3 8 3 9 4 4 1 时间点 影子的长度 图 7 附件 3 影子的长度与时间的拟合图形 利用最小二乘法拟合函数解得直杆的高度为 3 0356m 直杆的地理纬度 N32 847 直杆的地理经度 S15 962 使用 MATLAB 软件经过多次拟合 直杆的高度 直杆的地理纬度收敛 而直杆地理 经度不收敛 通过附件 3 数据 将纬度值代入参数方程 多次求解 最终可以得到直 杆地理经度为 E116 228 拟合参数 直杆的高度 3 0356m 直杆的地理纬度 N32 847 直杆的地理经度 E116 228 太阳赤纬角 S15 962 将拟合的参数用表格表示为 表 2 问题三附件 3 中求得的参数表 直杆高度 直杆地理经度 直杆地理纬度 太阳赤纬 参数范围 0 5 5 E0 E180 S90 N90 S23 45 N23 45 参数初始值 3 E115 N30 S15 拟合参数 3 0356 E116 228 N32 847 S15 962 由赤纬角与日期的时间的换算公式可以推导 N 307 09 取整为 N 307 大约为 2015 14 年 11 月 1 日 由上面模型的建立地理的知识可推导太阳在直杆 N20 7628 和 E115 7201 度的地 方一年会直射两次 则另外的一个日期为 2015 年 8 月 13 日 5 4 问题四模型的建立与求解 5 4 1 问题四模型的建立 视频中直杆在太阳下 8 00 10 00 之间的影子长度不断减小 如果在北半球 视频的 拍摄是自北向南拍摄 如果在南半球 视频的拍摄时自南向北拍摄 首先对视频进行处 理 将视频信息转化数据信息 利用米尺直接测量的方法 用米尺先测量出视频中直杆 的长度L 实际直杆的长度为 2 米 利用直杆测量值与直杆实际值的比值 可计算测量 数据与实际数据的比例尺 为比例尺 其表达式为为 H L 测量直杆的长度为 11 1cm 实际杆长为 2m 所以其比例尺为 H L 然后每间隔一段时间 测量一次直杆影子的长度 视频中的时间是从上午 8 54 9 34 中间间隔大约有 40 分钟 所以每 5 分钟测量一次直杆影子的长度 就可以得到一 组在不同时间的影子长度的数据 利用比例尺 就可将测量的影子长度转化为 实际情况下的影子长度 实际影子长度 测量影子长度 比例尺 测量影子长度与实际影 子长度的表格如图所示 表 3 不同时间的影子长度表 时间 测量影子的长度 CM 实际影子的长度 M 8 54 12 95 2 333 8 59 12 51 2 254 9 04 12 11 2 182 9 09 11 71 2 110 9 14 11 35 2 045 9 19 10 98 1 978 9 24 10 62 1 945 9 29 10 26 1 849 9 34 9 92 1 787 直杆的长度已知 测量的时间已知 建立模型确定拍摄地点 因此该模型中只有两 个未知参数 经度和纬度 通过对视频的转化 可得出 9 组数据 根据问题二中的模型 建立超定方程组 利用非线性最小二乘法拟合参数 估计经度参数和纬度参数的值 15 即 coscoscossin sintan arcsin H l 使得目标函数 i iiii tlll min 5 4 2 问题四模型的求解 将模型应用于附件4中数据的影子顶点坐标数据 将数据通过对已知函数进行拟合 利用 MATLAB 软件拟合得到附件 4 的图形 8 999 19 29 39 49 59 6 1 7 1 8 1 9 2 2 1 2 2 2 3 2 4 2 5 时间点 影 子 的 长 度 图 8 附件 4 影子的长度与时间的拟合图形 利用最小二乘法拟合函数 求解直杆的地理纬度 N40 04 直杆的地理经度 N118 03 与第二问类似 估计出最小二乘法参数范围 直杆地理纬度 S90 N90 直杆的 地理经度 E0 E180 参数初始值分别为 35 和 110 经过多次拟合 直杆地理纬 度收敛 数值为 N40 04 而地理经度不收敛 通过人工测量的影子数据 将纬度值代 入参数方程多次求解 最终可以得到直杆地理经度为 E118 03 如果拍摄日期未知 则就增加了一个未知量太阳赤纬 赤纬角 有三个未知量 纬度 经度 太阳赤纬 根据问题二的模型建立超定方程组 利用非线性最小二乘法拟 合参数 估计经度参数 纬度参数 太阳赤纬 赤纬角 参数的值 即 coscoscossin sintan arcsin H l 使得目标函数 i iiii tlll min 16 利用 MATLAB 解得直杆地理纬度数值为 N40 04 直杆地理经度为 E118 03 日期为 2015 年 7 月 12 日 由问题二模型的建立中地理的知识可推导太阳在直杆地理纬度数值为 N40 04 经 度为 E118 03 的地方一年会直射两次 则另外的一个日期为 2015 年 12 月 5 日 六 模型的评价与改进 优点 1 在问题一中将影子长度的变化模型转化为地球自转和公转运动的模型 其模型 的表示精确 而且可将与影子长度变化有关的因素联系在一起 2 在问题二中运用非线性最小二乘法解决超定方程 对于小残量问题 其收敛速 度很快 3 模型依靠地理知识理论支撑 能很好的解决所给的问题有很强的说服力 缺点 1 非线性最小二乘法 对于残量很大的问题或者非线性程度很高的问题 是不收敛 的 2 在问题二到问题四中使用的模型比较单一 3 在问题四中对视频处理时 直接采用米尺测量的方法 误差较大 模型的改进 8 在问题四中 对视频的处理 采用直接用米尺测量的方法 具有一定的误差 经上 网查询我们发现一个较为方便的计算机求解过程 成为基于视频中太阳影子轨迹的经纬 度估计方法 加参考文献 这种方法属于图像处理和地理信息系统技术领域 涉及一 种基于视频中太阳影子轨迹的经纬度估计方法 先获取自然图像序列或者视频帧 并且 对每一帧图像检测出影子的轨迹点 然后确定多个灭点 并拟合出地平线 拟合互相垂 直的灭点 计算出仿射纠正和投影纠正矩阵 进而还原出经过度量纠正的世界坐标 再 拟合出经过度量纠正世界坐标中的影子点的轨迹 并且参考日晷设计 利用相似关系估 计出纬度 利用二次曲线的极值点计算时差 从而有效地恢复出所拍摄图像的经纬度信 息这种复杂度较低 精度较高 是一种利用自然图像序列或视频帧实现经纬度估计的方 法 能够根据未经过校准的影子的位置 得出所拍摄图像的经纬度信息 参考文献 1 陈青云 吴毅明 计算机在建筑环境分析中的应用 北京 北京大学出版社 1986 年 10 月 2 陈平章 星空观测原理与方法 重庆 重庆出版社 1986 年 10 月 3 A M 库普林 地球与地图知识 广东 测绘出版社 1985 年 12 月 4 韩中庚 数学建模方法及其应用 北京 高等教育出版社 2009 年 6 月 17 5 何正风 MATLAB 在数学方面的应用 北京 清华大学出版社 2012 年 1 月 6 马莉 MATLAB 语言实用教程 北京 清华大学出版社 2010 年 1 月 7 杨德平 MATLAB 基础教程 北京 机械工业出版社 2013 年 1 月 8 操 晓 春 曲 颜 龄 孙 济 洲 视 频 中 太 阳 影 子 轨 迹 的 经 纬 度 估 计 方 法 年 9 月 14 日 18 附附 录录 附录 1 问题 1 的 MATLAB 程序 clc clear t 45 1 45 t t pi 180 a 15 144 pi 180 b 39 904 pi 180 d 3 y sin b sin a cos b cos a cos t h asin y l d tan h plot t 180 pi 15 l grid on xlabel 时间量 ylabel 影子的长度 附录 2 问题 2 的 MATLAB 程序 function d fun2 b t d b 1 sqrt 1 b 2 sin 10 8795 pi 180 sqrt 1 b 2 2 cos 10 8795 pi 180 cos 15 p i t 4 b 3 116 391 12 180 2 b 2 sin 10 8795 pi 180 sqrt 1 b 2 2 cos 10 8795 pi 180 cos 15 pi t 4 b 3 116 391 12 180 clc clear t xlsread C Users Administrator Desktop 试题 A 附件 1 3 xls 附件 1 E4 E24 d xlsread C Users Administrator Desktop 试题 A 附件 1 3 xls 附件 1 D4 D24 x resn resi exitflag lsqcurvefit fun2 2 0 5 115 t d x resn exitflag plot t d p t d grid on xlabel 时间量 ylabel 影子的长度 附录 3 问题 3 的 MATLAB 程序 f