太阳影子定位

太阳影子定位问题摘要本文主要研究的是利用视频中太阳影子的影长，方位等有效信息，建立合理精确的模型定位拍摄时间和地点。太阳影子定位技术对于刑侦、记录等方面有着重要意义。问题一，经过分析和查阅资料，得到太阳影长和太阳高度角、直杆高度的关系;s=h*tan,而太阳高度角与日期、经纬度、太阳时、时角、赤纬角、太阳方位角等参数存在数学关系，建立数学模型，在已知时间、经纬度的情况下，可以求得影长变化曲线，通过数据的导入与计算，我们得出天安门广场9:00-10:00的影长数据，并绘制为曲线图，最小值为在12:14时的3.8640米，最大值为在9:00时的7.7389米。问题二，我们先由X，Y坐标算出影长，并将影长与时间的关系进行傅立叶拟合，得出最短影长的时间，可通过与北京时间的时间差算出经度，在已知大概的经度范围内，再用最小二乘法的优化拟合进行检索，在用傅立叶拟合估算出经度后再进行最小二乘法优化大大提高了精度，也优于直接用最小二乘法计算的方法，经过数据导入与计算，算出杆所在地为海南省东方市，云南省陇川县和印度洋；误差最小的最优解为海南省东方市。问题三，在未知日期并只给出X，Y的情况下，我们先算出影长数据，再继续利用最小二乘法进行计算，将日期与赤纬角的关系表示出来并直接带入到总的公式中，则出现比模型二多检索一个变量的问题，我们仍然运用傅立叶拟合算出大致经度，减少一个变量，然后再进行检索的方法，最后分别检索出两组数据分别是4月25日与6月20日的新疆喀什地区和1月6日与11月9日的湖北省十堰市。问题四，在给出一个视频附件的情况下，我们先应用tracker软件对视频中的影长进行坐标追踪，利用tracker进行影长坐标的提取，又因为提取坐标不能直接应用于算法中，于是我们再建立模型，由小孔成像的原理得出提取坐标与实际坐标的关系，最后算出真实坐标，再在模型二与模型三的算法基础上进行数据的计算，最后计算得出地点为呼和浩特市玉泉区。关键词：最小二乘优化遍历搜索傅里叶拟合MATLAB一、问题重述视频数据分析在刑侦破案，记录生活等具有重要意义。本文主要是基于视频数据中的太阳影子的长度和方位等有效数据的挖掘，分析影子长度的地理影响因素，建立太阳影长的数学模型，确定视频拍摄的时间和地点。问题一：分析影响影长的各个参数的变化规律，如：太阳时、时角、太阳高度角、赤纬角、太阳方位角，建立合理的数学模型，得到影子的长度变化，方位变化。进而可以确定，天安门广场3米高直杆在2015年10月22日9:00-15:00的影子长度变化曲线。问题二：由附件一给出的各时刻影子顶点坐标数据，得到影长、方位等有效信息，建立数学模型求解附件一的可能位置。问题三：由附件二、三给出的各时刻影子顶点坐标数据，得到影长、方位等有效信息的基础上，建立数学模型得到附件二、三的可能位置和拍摄时间。问题四：已知一直杆在阳光下的影子变化视频和直杆长度，根据视频确定拍摄地点和日期。二、问题分析问题一分析，影子长度的变化规律集中体现在太阳高度角的变化规律上，经过查阅资料我们得到，影响参数有：直杆高度、太阳时、时角、赤纬角、太阳方位角、地理纬度等。我们将这些因素对太阳高度角的影响关系得到，就可以建立数学模型确定某地点一定时间内的影长变化曲线。问题二分析，因为直杆高度对于采用太阳高度角计算直杆位置至关重要，我们首先采用最小二乘优化拟合得到直杆的高度。由第一问的求解结果知道某一位置的太阳影子的变化与该地的位置该地所处时间有着密切联系，因此可根据第一问影子长度与时间、地点的关系，反推出纬度；已知北京的经度为120E，利用地方正午时间与北京正午时间的差异，可求出该地经度。问题三分析，在问题二的基础上，拍摄日期未知，即n未知。又由,3652842sin45.23度n将用n表示，我们选择了MATLAB中lsqucurvefit函数来对算式进行优化，用最小二乘法中的遍历搜索，筛选出可能的拍摄时间和地点。问题四分析，我们先应用tracker软件对视频的影长进行坐标提取，再由小孔成像的原理建立模型得出的提取坐标与实际坐标的关系，最后可以算出真实坐标，再在模型二与模型三的基础上进行数据的计算。三、模型假设1.假设分析在一天内地球的运动时，不考虑地球的公转；2.假设在计算日期序号时，二月都当做28天考虑；3.假设太阳光在空气介质中不发生折射；4.假设附件中给出的数据信息真实有效；四、符号说明太阳时影长时角某一时刻赤纬角影子坐标太阳高度角对应影长太阳方位角经度日期序号纬度杆长北京时间杆距平面的距离像中杆的长度像中杆中心距离stAnhsit(,)iixyilBtBtl'hr五、模型的建立与求解5.1.问题一5.1.1模型的建立设杆子的长度为h,竖直杆在太阳光照射下影子的长度为s,太阳高度角为图1由图1知tanhs即tanhs太阳高度角是太阳相对于地平线的高度角，这是以太阳视盘面的几何中心和理想地平线所夹的角度,其与时角，当时的太阳赤纬，当地的纬度有着密切的联系。太阳高度角可以使用下面的算式，经由计算得到很好的近似值：,coscoscossinsinsin下面是各个参量的定义与算法。1.太阳时st时间的计量以地球自转为依据，地球自转一周，计24太阳时，当太阳达到正南处为12:00。钟表所指的时间也称为平太阳时（简称为平时），我国采用东经120度经圈上的平太阳时作为全国的标准时间，即“北京时间”。在计算地方xyhsa真太阳时st的计算公式为(120)/15sBtt其中Bt为北京时间，为当地度经度值2.时角记为，时角以正午12点为零度开始计算，每一小时为15度，时角的值可正可负（上午为正，下午为负）；,1215度st3.赤纬角记为，赤纬角记太阳直射纬度，近似计算公,3652842sin45.23度n其中，日期序号记为n，例如：1n，表示1月1日；85n，表示3月26日。由以上参量的关系，已知日期n可以求出赤纬角,已知时间t可以求出时角，影长s可以由坐标x,y求出；则在方程式：tan((sinsincoscoscos))hsarc中只有影长s一个未知量，则问题求解。5.1.2模型的求解由题意知杆的长度：3hm该地的经度：=116.391E该地的纬度39.907N由日期值得日期序号273n,将其代入式(4)求得4.067(度)将北京时t=9:00-15:00的代入式(3)求得每一时刻对应的时角tw最终将这些参量值代入高度角公式得到太阳高度角在不同时刻值t将时刻t作为横轴，太阳高度角作为纵轴利用MATLAB数学软件对数据处理得到太阳高度角随时间的变化曲线（具体程序见附录附件1），得到关系图如下：从图中可以直观的看出，大约在12：14左右，该地的太阳高度角达到一天中的最大值，为38。又因为tansh得到太阳高度角变化曲线，即可绘制出影长变化曲线。如下图所示：图2图3由图可以直观的看出，在该地12:14时刻，因为天安门广场经度与北京时间的120度经度有一定的经度差，则取当地时间的正午会与北京的12:00有一定的时间差，此时太阳高度角达到一天中的最大值，直杆影长达到一天中的最小值，约为3.8640米。影子最长的时刻为9：00，约为7.7389米。5.2.问题二时间影长9:007.738910:005.336911:004.257412:003.862813:003.998214:004.716615:006.34485.2.1.求经度的模型问题二我们采用最小二乘法求解，但由于在计算过程中，未知量过多，误差较大，所以我们采用傅里叶先拟合出大致经度数据，再让其进行检索，极大的提高了最后算出结果的精度，以下为傅里叶拟合建模方法：设在某天某地时刻为：123,,ntttt分别取一固定杆影子的坐标：112233(,),(,),(,)(,)nnxyxyxyxy，对应影长为；22iiilxy(1,2,3)in根据地理知识知当一天中影子最短的时刻即为该地地方时间的12:00，由于所有函数都可展开成傅立叶级数,因此可采用MATLAB中的拟合工具将时刻与影长进行傅立叶拟合，（拟合方法见附录2）进行人工遍历搜索找到最短影长及其对应的北京时间0t，在中国某一地理位置的经度可根据下式求得012015(12)t（度）我们用matlab将影长和时间进行一次拟合,来研究影长与时间的关系。图4从图中可以看出，影长与时间基本上是曲线关系傅里叶拟合方程为：011?fxaacosxwbsinxwa0=25.82(-310,361.6)a1=-4.251(-304.8,296.3)b1=-24.94(-314.7,264.8)w=2.656(-15.69,21.01)将傅立叶拟合出来的系数带入，再将时间x的值进行扩范围选取，便能拟合出其他时刻的影长值，拟合结果如下：图55.2.2求纬度、杆长的模型由第一问的求解结果知道某一位置的太阳影子的变化与该地的位置和当地所处时间有着密切联系，因此可根据第一问某地某一时间段影子长度变化与时间地点的关系的思想反推出直杆所处的位置。在地球表面上的位置是由经度和纬度确定的。我们已知北京的经度为120E利用地方正午时间与北京正午时间的差异，可由此求出该地经度。最小二乘优化：在无约束最优化问题中,有一些重要的特殊情形,比如目标函数由若干个函数的平方和构成,这类函一般可以写成：其中:1[,]Tnxxx一般假设mn。把极小化这类函数的问题：21min()()miiFxfx称为最小二乘优化问题。最小二乘优化是一类比较特殊的优化问题，在处理这类问题时，MATLAB提供了一些强大的函数有lsqlin，lsqucurvefit，lsqmonlin，lsqmonneg.在此我们选择了MATLAB中lsqucurvefit函数来对算式进行优化，以确定的经度值为背景，作为最小二乘优化中便利搜索的起始范围。由题意我们构建了拟合函数tan((sinsincoscoscos)hyarc其中y为影长，h为杆的长度，为时角，为当时的太阳赤纬，为当地的纬度。由问题一知赤纬角的计算公式为,3652842sin45.23度nn为日期序号；将附件一的数据得：=10.511时角的计算公式为,1215度st其中真太阳st的计算公式为(120)/15sBtt其中Bt为北京时间，为当地纬度值。21()(),mniiFxfxxR将st的计算公式代入的计算公式得15((120)/1512)Bt(度)根据上式公式可将具体数值带入拟合函数得到tan((sinsincoscoscos(15((120)/1512)))Bhyarct已知的观测数据如下表，求拟合函数的参数h,,,（lsqucurvefit函数算法程序见附录模型二程序）求得下列几组可能的地点值（不合理的位置点已经去除）杆长纬度经度地点误差范围2.036519.2732108.6451海南省东方市3.52E-082.765524.577198.1058云南省德宏傣族景颇族自治州陇川县9.70E-062.2500-3.6110102.7380印度洋1.09E-05（默认正值表示北纬，东经）5.3.问题三我们沿用第二问的思路，采用无约束拟合中的最小二乘优化根据题意得拟合函数（观测数据与程序见附录3）将赤纬角，时角，高度角公式联立得出的总方程式：2+n2+n(120)tan((sinsin(23.45sincoscos(23.45sincos(15(12)))36536515Bhyarct（284）（284）（）)（）)继续应用模型二最小二乘法将n与赤纬角的关系式列出将模型二中的n换为X（4）将关系式直接带入总的方程在模型二的基础上将日期n变化为方