《太阳影子定位-2015年全国数学建模比赛a题全国二等奖论文.doc》由会员分享，可在线阅读，更多相关《太阳影子定位-2015年全国数学建模比赛a题全国二等奖论文.doc（35页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、太阳影子定位摘要本文研究的问题是分析直杆在太阳的照射下，影子的角度和长度的变化，再结合相关地理知识和数学几何模型，推算出具体的所在地点和具体日期。该模型可以用于太阳影子定位技术中，根据物体在阳光照射下影子的变化进行定位。对于问题一，我们首先根据地球与太阳的位置关系列出太阳赤纬角，太阳高度角，太阳时角的计算式，其中需对较粗略的太阳赤纬角计算式进行修正，得出精准的计算式。再建立数学几何模型，根据太阳高度角，影长与杆长形成的角边关系，列出影长的计算式。最后建立一个太阳日照影长模型，该模型以太阳高度角计算式，太阳赤纬角计算式，太阳时角计算式为子函数，以太阳赤纬角，太阳日角，太阳时角，时间初值为中间变量
2、，以当地经纬度，从1月1日到测量日的天数，时间，杆长，年份为自变量的复合函数数学模型。然后采用由内到外计算法对此复合函数进行求解，计算出从九点到十五点的影长和太阳高度角的变化，得出直杆的太阳影子长度的变化曲线。对于问题二，我们首先分析因为时间日期已给出，所以根据太阳赤纬角计算式可知太阳赤纬角为已知量，接着我们将影长的计算式进行等式移项变换，得到一个拟合杆长及经纬度的非线性最小二乘模型，该模型将问题一中太阳日照影长模型作为参数拟合对象，以杆长和影长与太阳高度角正切值之积的差值最小误差平方和为目标函数，以太阳高度角计算式，太阳时角计算式为约束条件，以测量时间，天数，影长为已知量。将该模型在1sto
3、pt软件中运行，采用麦夸尔特算法和通用全局最优化法对该模型进行迭代计算，对实验结果统计分析后得出该直杆相应的北纬为19.29392848度，东经为108.7225248度（海南岛的西海岸）。对于问题三，除了需要拟合杆长和经纬度以外，还需拟合日期，同样参照影长等式移项变换公式，得到一个拟合杆长、经纬度及日期的非线性最小二乘模型。同样采用问题二的计算方法得到多组结果，其中附件二最优解地点为新疆维吾尔自治区喀什地区巴楚县(40.0025N,79.6587E)，附件三最优解地点为湖北省十堰市郧西县（32.9638N,110.277E）。对于问题四，我们首先将视频每隔2分钟截取图片，共取得20张截图，根
4、据图像获取直杆顶点，直杆固定点，影子顶点的坐标，以及观测影长，由于需要引入角度才能由观测影长求得实际影长，因此对问题一中的太阳日照影长模型进行改进，初始旋转角，旋转角增量，进行参数拟合，与问题二类似。在观测影长已知的条件下得到一个拟合经纬度及旋转角度的非线性最小二乘模型。经过角度分析，得出角度范围是 ，再根据最优解值分析，筛除不符合要求的项，最后剩下第四组数据即地点为湖南省永州市宁远县（25.86216N,111.9039E）。如果拍摄日期未知，我们可以在问题四的基础上，增加一个拍摄日期的拟合变量，即可解决问题。关键词：太阳日照影长模型 由内到外计算法 复合函数数学模型 非线性最小二乘模型 麦
5、夸尔特算法 通用全局最优化法 1stopt软件一、 问题的重述在视频数据分析里有两个很重要的方面，一是确定该视频的拍摄地点，二是拍摄日期，太阳影子定位技术原理就在于分析物体在阳光照射下影子的变化，从而确定该视频拍摄的时间与地点。1.建立关于影子长度变化的数学模型，并画出2015年10月22日北京时间9:00-15:00之间天安门广场高度为3米的直杆的太阳影子长度的变化曲线。2.对于固定直杆的太阳影子顶点横纵坐标数据，使用相应的数学模型确定直杆所处的地点位置。3.对于固定直杆的太阳影子顶点横纵坐标数据，使用相应的数学模型确定直杆所处的地点位置和日期。4附件4的视频描述了一根直杆在太阳照射下的影子
6、变化情况，并且已估计出此为高度2米的直杆。需建立相应的数学模型确定视频拍摄地点位置，并通过应用该模型给出若干个可能的拍摄地点。若未给出拍摄日期，你能否根据视频内容确定拍摄地点与日期？二、 模型的假设1 假设太阳射向地球的光线为平行光线。2 假设地球为均匀球体，且球面平整。3 假设忽略光传播过程中所需要的时间。4 假设忽略大气折射对光线传播路径的影响。5 假设在此过程中，忽略地球公转对影子长度、角度计算产生的影响。6 假设摄像机拍摄角度与地面平行。三、 符号及说明符号名称影长杆长太阳高度角地理纬度太阳时角太阳赤纬角从1月1日到测量日的天数测量的时间 年份四、 模型的准备从图1中，我们可以形象地看
7、出，赤道面与日地中心连线的夹角每天都在变化，这个角度就是所说的太阳赤纬角。正如我们所知，在地球自转过程中也在围绕太阳公转，而极轴与黄道面的夹角始终保持不变。正是因为以上原因，才造成每天正午时刻，太阳高度角有所差异。图1 太阳赤纬角模型太阳高度角是指太阳光的入射方向和地平面之间的夹角，它是由地理纬度，赤纬角与时角综合决定的。图2 太阳高度角和时角模型而日面中心的时角，是从观测点天球子午圈沿天赤道量至太阳所在时圈的角距离。是确定当天时刻与太阳高度角关系的纽带。五、 问题的分析及建模求解5.1 问题一的分析与建模求解5.1.1 问题一的分析针对问题一，根据提供的时间及地理位置，结合相关影响因素，变化
8、规律，得出直杆的太阳影子长度的变化曲线。为了求出长度变化曲线，需要分析有哪些变量构成影长这个变量。变化曲线的描述需要查阅资料构建模型，通过对经纬度，太阳时角，太阳赤纬角等求解修正即可计算出各个时间的影长，并绘制出太阳影子长度的变化曲线。5.1.2 模型的分析本题要求出影子长度的变化曲线，即求出影长随其他变量的变化关系，影子长度与多种因素相关联，其中太阳高度角与影长，杆长构成一种角边关系，利用这种角边关系可以求得影长。为了利用题目中所给出的经纬度，时间等条件，需引入太阳高度角的关系表达式。太阳高度角是关于纬度，太阳赤纬角，太阳时角的复合函数，进而需求出太阳赤纬角，太阳时角等物理量，得出影长即影子
9、长度的变化曲线。根据题目要求及假设，我们先建立了一个地球围绕太阳旋转的模型，并且翻阅各种文献和找出了一系列能够影响影长的因素。我们首先利用相似三角形的方法，列出了计算影长的公式：（1）其次，我们考虑到如果需要求出影长需要求出高度角和杆长，现在杆长已知3m，所以我们只需要求出太阳的高度角即可。太阳高度角是太阳相对于地平线的夹角，因为这是从太阳是盘面的几何中心到理想地平线的夹角。太阳高度我们使用了一下的算式来求解太阳的高度角的近似值1：（2）其中，为太阳高度角；为太阳时角；为当时的太阳赤纬角；为当地的纬度（天安门的维度为39度54分26秒）在求解计算高度角的过程中，需要求解太阳赤纬角、地理纬度和太
10、阳时角。现在题目中已知的是地理纬度为北纬39度54分26秒，另外太阳赤纬角也称为太阳赤纬，就是所谓的太阳直射纬度，它的计算公式为：但我们需修正该公式，要考虑多种因素。因为在周年运动中的任何时刻，太阳赤纬角的具体值均是严格已知的，所以也可以用以下的表达式来表述，即2：（3）上式中的称为日角，即2：（4）这里的又由两部分构成，即2：（5）上式当中，表示积日， 表示时间初值。所谓积日，就是所给日期在当年内的顺序号，例如，1月1日的积日为1，闰年12月31日的积日为366，而平年则为365。将式（5）代入进式（4）得式（6），即：（6）其中为2：（7）在天文学中太阳时角这个名词，意为一个天体的太阳时角
11、被定义为该天体的赤经（RA）与当地的恒星时（LST）的差值。它的计算公式为1： （8）其中，为太阳时角；为真太阳时；为当地的经度，我们考虑到问题所给的经纬度并非准确的北京天安门前的经纬度，并且有四度的偏差，所以我们小队时差来减小由于经纬度不准确的误差。5.1.3 模型的建立最后根据以上四个式子建立起一个数学模型，即：该函数是基于多个子函数的复合函数，其中以影长为因变量，以当地经度，当地纬度，从1月1日到测量日的天数，时间，杆长，年份为自变量，得出太阳日照影长模型：5.1.4 模型的求解将四个算式（1）(2)（3）（4）联立起来采用由内到外计算法进行求解。步骤一，根据已知变量求算出内侧函数。步骤
12、二，再根据已得，算出外侧函数的内侧函数。步骤三，将外侧函数带入复合函数，求得我们最后的结果。表1 影长、高度角、赤纬角随时间的变化表时间影长高度角太阳赤纬角9时7.28522.383-10.86310时5.04930.716-10.86311时4.02336.712-10.86312时3.64139.486-10.86313时3.77238.497-10.86314时4.45633.949-10.86315时5.98326.629-10.863我们首先用Microsoft visual c+6.0编写了一个C语言代码，用来计算并确定从九点到十五点的影长和太阳高度角的变化，计算结果如表1所示。图
13、2 太阳影子长度的变化曲线然后将数据放入MATLAB中，画出题目所需的影长随时间变化的曲线，如图2所示。5.2 问题二的建模及求解5.2.1 问题二的分析针对问题二，根据提供的日期，时间，影子顶点横纵坐标数据，求出若干个可能地点的经纬度，依据太阳时角、太阳高度角及太阳方位角之间的关系，根据参考问题一的求解过程和数学模型，确定已知量和变量，从而建立数学模型进行求解。 5.2.2 模型的分析本题要求出若干个可能地点的经纬度，根据影长，杆长和太阳高度角这三个量的关系式，我们采用基于多约束下的目标优化数学模型，用21组数据预测出经度，纬度，杆长，使得将量放到函数时能形成一个满足的函数形状，同时采用非线
14、性最小二乘法的思想，通过杆长和影长与太阳高度角正切值之积差值的最小误差平方和，寻找数据的最佳函数匹配。利用最小二乘法求得未知的数据，使得这些求得的计算数据与实际数据之间的误差平方和最小。建立以此平方和最小为目标函数，太阳高度角计算式，太阳时角计算式为约束条件的数学优化模型，在1stopt环境下对其进行分析影响各种因素的因子，运用麦夸尔特算法和通用全局优化法进行迭代，得出结果。根据题目要求及假设，我们先建立了一个地球围绕太阳旋转的模型，并且翻阅各种文献，找出了一系列会对影长产生影响的因素。我们首先利用相似三角形的方法，列出了计算影长的公式：接着建立优化模型，且采用非线性最小二乘法。在根据此公式设
15、定最小化误差的平方和，以此作为优化模型的目标函数，即：（1）其中为杆长和影长与太阳高度角正切值之积的差值最小，为太阳高度角,为直杆长度，为影长，后三个量为未知量。 然后在列出太阳高度角的计算公式作为约束条件之一，即1： （2）其中为太阳赤纬角，因为日期时间为已知量，所以可以求得太阳赤纬角，为已知量；为太阳高度角，为太阳时角，为未知纬度。为了表达太阳时角，列出了该量的计算表达式，也是作为约束条件表达式之一，即1：（3）其中为真太阳时，数据已给出，为已知量；为未知经度。5.2.3 模型的建立最后根据以上三个式子建立起一个以杆长和影长与太阳高度角正切值之积的差值最小为目标函数，太阳高度角计算式，太阳
16、时角计算式为约束条件的数学优化模型，已知测量时间，天数 ，影长 ，得到一个拟合杆长及经纬度的非线性最小二乘模型：5.2.4 模型的求解该模型为非线性模型，采用了非线性最小二乘法。将式子代入1stopt软件内进行求解。（1）麦夸尔特算法麦夸尔特算法类似于阻尼高斯-牛顿法，在最小二乘法的意义下，可以给出确定的非线性模型中的参数，并且增加流动性，从而解决有关非线性最小二乘法的问题。（2）通用全局最优化法同时运用通用全局最优化法进行迭代，跳出局部最优解，足够平均足够多地随机投放初值，通过一个明确的目的走出局部最优解。所以我们采用麦夸尔特算法和通用全局最优化法进行求解。3步骤一，设已知数据矩阵首先，给出
17、m个参数的初始值，由计算N组数据的残差平方和Q。步骤二，根据最小二乘法的原则，令，满足线性方程组： 步骤三，解方程组，反复迭代，投放不同的初值，直到达到要求的精度为止，当时，迭代结束。5.2.5 结果的分析因为该模型为非线性模型，所以受初值的影响即与初值的选择有关，每次最优解不同，应从多组解选择最佳对策。所以我们将程序运行20次，针对不同的初值进行检验，提高结果的可靠性，计算结果见表2：表2 第二问运行结果统计分析表运行次数 12.0336080204346619.29392875314108.72252439306122.0336080056719919.2939285728313108.7
18、225246325432.0336079946937119.2939284387685108.72252481062642.0336080135537419.2939286692245108.72252450466552.0336079946937119.2939284387685108.72252481062662.0336079871999519.2939283472871108.72252493218372.0336079890342719.2939283696747108.72252490242982.0336079853738819.2939283249479108.72252496
19、18192.0336080307917919.2939288796215108.72252422505102.0336080297674419.2939288670944108.722524241669112.0336080078170419.2939285990609108.722524597739122.0336079780441619.2939282354854108.722525080704132.0336080009588919.2939285153294108.722524708987142.0336080082750119.2939286046437108.72252459031
20、1152.0336079622024219.2939280419598108.722525337693162.0336079906726519.2939283897123108.722524875847172.0336080051400619.2939285663597108.722524641165182.0336079858916419.2939283313124108.722524953405192.0336079783460419.2939282391408108.722525075811202.0336079973425419.2939284711314108.72252476765
21、5均值2.03360799819.29392848108.7225248方差 根据表2可得出，该地点坐标为北纬19.29392848度，东经108.7225248度，均方差（RMSE）为，残差平方和（SSE）为，相关系数（R）为0.99，卡方系数（DC）为根据表2，得到结果在地图上的位置如图3所示。图3 问题二经纬度结果图示5.3 问题三的建模及求解5.3.1 问题三的分析针对问题三，提供了时间，影子顶点横纵坐标的数据，与问题二进行对比，日期变成未知量，同样求出若干个可能地点的经纬度，依据太阳时角、太阳高度角及太阳方位角之间的关系，考虑到问题二，此时仅需在此基础上将条件稍加改变，即能得出结果。
22、5.3.2 模型的分析本题同样要求出若干个可能地点的经纬度，相较于问题二，多了一个未知量，所以优化数学模型要进行相应的改变。考虑问题二已建立的模型，此时仅需将约束条件稍加改变，增加未知变量的计算式，同样采用的是非线性最小二乘法，利用最小二乘法求得未知的数据，使得这些求得的计算数据与实际数据之间的误差平方和最小。建立以此平方和最小为目标函数，太阳赤纬角计算式，太阳时角计算式，太阳高度角计算式为约束条件的数学优化模型。基于对问题三的分析，我们在问题二所建立模型的基础之上，进行数值化求解。首先由变量间的关系分析，得到关系表达式：其中，影长e、经度C、纬度、日期序号n以及杆长H为未知变量。sT是已知变
23、量。在问题二的求解中，根据相似三角形的方法，给出了影长的计算公式：由此我们采用非线性最小二乘法，建立优化模型。由于直杆的长度是不变的，所以对于第i组数据的变量H，应该与第（i+1）组数据相差尽可能小。据此，我们建立目标函数，即：其中 为杆长和影长与太阳高度角正切值之积的差值， 为第i组数据所表示的太阳高度角，是第i组数据所表示的影长。5.3.3 模型的建立最后根据上面得到的式子，以影长e、经度C、纬度、日期序号n为决策变量，已知测量时间影长得到一个拟合杆长、经纬度及日期的非线性最小二乘模型：5.3.4 模型的求解与实验结果的分析考虑到麦夸尔特算法专门应用于解决“包含多个变量的最小二乘法的非线性
24、拟合”问题，且其阶梯式拟合方式有利于求解。但是当曲线变化较小的时候，单靠麦夸尔特算法会加剧运算的迭代次数。所以，加入通用全局最优化算法可以减少迭代次数，提高效率。我们利用Levenberg-Marquardt法与UniversalGlobalOptimization（通用全局最优化算法）分别用附件的数据在1stopt软件中进行求解，并对得到的数据结果进行统计。经过分析，我们了解到所求解的值离散程度较高，所以考虑对多组计算结果进行筛选择优。对结果进行排布后，发现取值呈现概率化差异分布，某些解值完全符合我们的要求，并且出现的概率较高，而其他的解值则呈现小概率，不符常理的特点。所以，我们统计多组数据
25、，排除不可能的解值，保留解的众数，来作为我们的全局最优解，得到结果见表3：表3 问题三针对附件2的结果结果数据组别杆长时间北纬东经地理位置11.998642015.07.2739.86179.6099新疆维吾尔自治区喀什地区巴楚县22.009162015.01.2740.044479.6751新疆维吾尔自治区喀什地区巴楚县32.026572014.11.1540.363179.7719新疆维吾尔自治区阿克苏地区柯坪县42.000152014.01.2739.88279.6223新疆维吾尔自治区喀什地区巴楚县51.998962013.05.1739.861979.6145新疆维吾尔自治区喀什地区
26、巴楚县E(X)2.0067-40.002579.6587新疆维吾尔自治区喀什地区巴楚县通过对以上精选而出的数据进行分析，我们可以看出具体日期的波动较大，而其他数据拟合程度较好，根据经纬坐标进行查询，得到所求地点的经纬如图4图4 精选地点分布示意图得到结果见表4：表4 问题三针对附件3的结果结果数据组别杆长时间北纬东经地理位置13.059552014.04.0533.0262110.223湖北省十堰市郧西县23.035832014.03.2032.8585110.246湖北省十堰市郧西县33.023212013.11.1832.9026110.27湖北省十堰市郧西县43.006972013.09
27、.2433.0428110.31湖北省十堰市郧西县53.030952013.11.1832.9448110.262湖北省十堰市郧西县62.989852012.10.2732.9262110.326湖北省十堰市郧西县73.014152012.10.2733.0457110.299湖北省十堰市郧西县E（X）3.022929-32.96381110.2766湖北省十堰市郧西县对以上数据，我们标出其在地图上的位置如图5。图5 对附件3数据的处理结果5.4 问题四的建模及求解5.4.1 问题四的分析问题四中给出了一根直杆在太阳下的影子变化的视频，并且给出了直杆的高度为2米。如果确定影子长度的变化关系，就
28、可以利用麦夸尔特-通用全局最优模型进行求算。但是由于视频拍摄的角度以及图像透视效果的影响，直杆的影子长度不能由所给视频直接得出，我们考虑“如何运用合适的转换方式，将视频中的影长提取出来”，这才是解决问题的关键所在。5.4.2 模型的分析对于问题四，我们首先需要建立一个图像处理的模型，对数据进行处理。我们发现，影子围绕直杆转动时，其影子所在直线与我们所观察图像的水平方向有一个初始角度，设其为；当以相等时间间隔获取截图时，影子绕直杆底端定点所转过的角度应该相等，设其为。而实际的影长应该为图像观测长度与余弦的比值。通过对太阳高度角公式，我们可以得到图形观测中的影长与杆长的关系。其中，建立图像处理模型
29、并解决视觉透视问题是解决问题的核心。我们可以利用此模型找出直杆的顶点的坐标，直杆固定点的坐标，影子的顶点坐标。接着根据三个点的坐标的相对关系，求出杆长和影长，便可以通过上面的思路，对第四问求出最优解答。首先我们进行对图像的处理。对图像进行处理，其根本目的是更加方便地抓取上面所提及的三个点的坐标位置。我们利用MATLAB将视频中的直杆与影子的图形关系每隔3分钟截取图片，共取得14张截图，接着对其进行灰度二值化处理。由于影子的颜色深度要高于地面，所以可以将直杆影子的位置明确地展示出来。接着以直杆固定点为坐标原点，以地面水平向右方向为x轴正向，以地面上平行于场边缘线且经过远点的直线为y轴，以直杆竖直
30、向上方向为z轴正向，建立三维直角坐标系。标定图片上给出的顶点坐标，并以此计算出各个时刻的杆长和相对应的影长。针对于拍摄角度和透视原因所造成的长度测量问题，我们分析了影子角度以及观测影长的关系：其影子与水平方向的初始夹角与第n张截图影子绕直杆底端定点所转过的角度的差值的余弦值，等于图像观测影长与实际影长的比值，即：其中，且n为正整数。结合太阳高度角公式，建立出影长角度最优化模型。对于初始数据，有：对于第一张截图，影子转过，此时影子所在直线与图像水平方向的夹角为，我们将之带入模型中，可以得到：5.4.3 模型的建立分析以上关系表达式可以看出，能够将该表达式推广应用至每一时刻，已知观测影长 ，得到拟
31、合经纬度及旋转角度的非线性最小二乘模型：其中, H为杆长，为第n张截图的图像观测影长，为太阳高度角，为每隔相邻时间的影长图像角度差。通过这种模型的建立可以最大限度地减小透视效果和拍摄角度的影响，方便我们对直杆和影长的关系进一步分析，进而求解其经度纬度。5.4.4 模型的求解与计算结果分析通过上面麦夸尔特-通用全局最优化模型的建立，处理数据。由于在建立的上述模型基础之上，需要考虑影子与水平方向的初始夹角与第n张截图影子绕直杆底端定点所转过的角度的差值的关系。也就是说，当n=20，取第20张视频截图的时候，由于图像中影子方向与直杆近乎垂直。我们据此判断，需要增加的约束条件为：剔除不合理的数据，得到
32、结果见表5：表5 问题四结果结果数据组别北纬东经地理位置138.37867111.0708山西省吕梁市兴县225.82818111.341广西壮族自治区桂林市全州县338.37553111.0178山西省吕梁市兴县425.86216111.9039湖南省永州市宁远县532.0924111.0168湖北省十堰市房县625.80764110.9909广西壮族自治区桂林市全州县738.37482111.0045山西省吕梁市兴县832.10841111.2884湖北省襄阳市谷城县938.42172111.7905山西省忻州市静乐县1025.80579110.9575广西壮族自治区桂林市全州县1132.
33、10325111.2053湖北省襄阳市谷城县1225.81126111.0532广西壮族自治区桂林市全州县经过角度分析，得出角度范围是 ，再根据最优解值分析，筛除不符合要求的项，最后剩下第四组数据即地点为湖南省永州市宁远县（25.86216N,111.9039E）图6 筛选后的地点示意图如果拍摄日期未知，我们可以在问题四的基础上，增加一个拍摄日期的拟合变量，即可解决问题。六、 模型改进与评价6.1 模型的改进针对上面给出的模型，从形式上讲，他们与一般书籍中给出的并无不同。我们之所以又重新研究它，是因为以往的公式存在以下通病2：（1）我们所需要的数值，会因为所在地的地理经度及具体时刻与表值有异而
34、不同。（2）对平年和闰年不加以区分，一方面，这对闰年就不好处理，另一方面，闰年的影响有累积效应，会逐步增长具体的讲，一般需要如下三项订正：（1）经度订正：即使我们查阅的是当年天文历，也需要这项订正。在我国各个地理经度范围内，各地的订正值是：（2）年度订正：根据我们查阅资料，如果不是仅仅使用当年的天文年历值，就要对具体的时间做出修正。假设我们选用的是闰年后某年的表值，第二年再用时，就要加上-0.2（-0.2422）日的订正。这个订正到了第三年为-0.51（-0.4844）日，第四年是-0.7（-0.7266）日，但这一年为闰年，多了一日，实际应该订正应为-0.7+1=0.3（0.2734）日，以
35、此类推。-0.2日-0.3日-0.4日6.2 模型的评价6.2.1 模型的优点（1）本论文中的模型都是递进的关系，均是从简单到复杂递推出来一步步建立的。问题二是基于问题一的基础上建立的，问题三基于问题二，可读性强，文章整体上逻辑清楚，连贯。便于阅读。（2）对于问题一的模型，我们结合地理知识与简单的数学公式建立符合函数模型，易理解，利用MATLAB和C语言结合解答（3）本论文问题二三四的模型求解都使用了1stopt软件，该软件对于非线性曲线拟合，综合优化分析计算更加专业，解题速度更快。6.2.2 模型的缺点（1）本论文有些模型关于影响因素的考虑还不是很全面，会影响结果精度，产生误差。（2）参量过
36、多，会影响结果精度，产生误差。七、 参考文献1 刘琦,王德华, 建筑日照设计, 北京：水利水电出版社, 2008-09。2 王炳忠, 日期、时间和当地经纬度计算太阳天顶角和方位角，3 麦夸尔特法算法, 2012-10-01。附录一：表6 问题三针对附件2解值表DATAhDayNE最优解值11.9986415.07.2739.86179.60990.0007622.0091615.01.2740.044479.67510.0091432.0265714.11.1540.363179.77190.5445542.0001514.01.2739.88279.62230.847251.9989613.
37、05.1739.861979.61450.00056E(X)2.0067-40.002579.65870.28044VAR0.000114-0.0372260.003740.124226附录二：表7 问题三针对附件3的解值表DATAhDayNE最优解值13.0595514.04.0533.0262110.2230.964523.0358314.03.2032.8585110.2460.2922633.0232113.11.1832.9026110.270.0186443.0069713.09.2433.0428110.310.1327553.0309513.11.1832.9448110.26
38、20.006862.9898512.10.2732.9262110.3260.5787473.0141512.10.2733.0457110.2990.85005E（X）3.022929-32.96381110.27660.406248VAR0.000428-0.0047790.0011550.133707附录三：表8 第四问运算结果表angleangle1NE1.0000030.0637744.6574110.99431.0000030.17158838.37867111.07081.0000030.0702925.807110.97961.0000030.04465950.94183111
39、.01361.0000050.17061544.65541110.95731.0000030.14249425.82818111.3411.0000030.08051338.37553111.01781.0000070.02835750.99485111.90281.0000030.03520725.86216111.90391.0000030.14780332.0924111.01681.0000030.03050925.80764110.99091.0000030.1045250.94788111.11851.0000030.10499638.37482111.00451.0000030.
40、05187232.10841111.28841.0000030.14820450.94872111.13271.0000040.11106338.42172111.79051.0000060.0203425.80579110.95751.0000050.0217132.10325111.20531.0000030.14719825.81126111.05321.0000030.06676424.8886495.77293其中，angle为定义中的，angle1表示，N表示纬度，E表示经度附录四（第一小问C语言代码）：#include #include #define Pi 3.1415926#
41、define n 295#define Y 2015/影长ldouble _l(double h,double A)double l;l=h/tan(A);return l;/赤纬角S,输入日期和年份double _S(double a,double y)double S,angle,N0,t;N0=79.6764+0.2422*(y-1985)-(int)(y-1985)/4);t=a-N0;angle=2*Pi*t/365.2422;S=0.3723+23.2567*sin(angle)+0.1149*sin(2*angle)-0.1712*sin(3*angle)-0.758*cos(a
42、ngle)+0.3656*cos(2*angle)+0.0201*cos(3*angle);return S*Pi/180;/高度角Adouble _A(double N,double S,double W)double A;A=asin(sin(N)*sin(S)+cos(N)*cos(S)*cos(W);return A;/时角Wdouble _W(double t) double W;t=t-(120-116.3914)/15;W=15*(t-12);return W*Pi/180;void main()double l;double h,N,t,A,S,W;while(1)printf(请依次输入杆长h，地理纬度N，时间t：n);scanf(%lf%lf%lf,&h,&N,&t);/h=3;N=39.54;/scanf(%lf,&t);N=N*Pi/180;W=_W(t);S=_S(n,Y);A=_A(N,S,W);l=_l(h,A);A=A*180/Pi;S=S*180/Pi;printf(时间为：2015年10月22号%.0f时n,t);printf(影长l=%.3fn,l);printf(高度角A=%.3fn,A);printf(赤纬角S=%.3fn,S);printf(n);附录五（第一小问MATLAB代码）：clear allclcn=295;