馬玉雪，令狐林玉，楊東海，孫國偉，買阿麗

(運(yùn)城學(xué)院 應(yīng)用數(shù)學(xué)系，山西 運(yùn)城 044000)

基于二次函數(shù)擬合的太陽影子定位技術(shù)研究

馬玉雪，令狐林玉，楊東海，孫國偉，買阿麗

(運(yùn)城學(xué)院 應(yīng)用數(shù)學(xué)系，山西 運(yùn)城 044000)

利用影子長度和太陽高度角的關(guān)系研究太陽影子定位技術(shù)，得到太陽高度角為45度時(shí)影子和物體一樣長；并通過對太陽影子長度變化的數(shù)據(jù)進(jìn)行分析，得出影子長度與時(shí)間的二次函數(shù)擬合方程，從而根據(jù)影長最短的時(shí)間點(diǎn)確定拍攝地點(diǎn)的經(jīng)度；然后根據(jù)太陽方位角與時(shí)角，經(jīng)緯度，高度角之間的關(guān)系，求出高度角，從而得到影子和緯度的關(guān)系；最終確定出視頻的拍攝地點(diǎn)和拍攝日期。

太陽高度角；太陽方位角；二次函數(shù)擬合；影子定位

0. 引言

太陽影子定位技術(shù)就是通過分析視頻中物體的太陽影子變化，通過測量影子方位及長度與桿的實(shí)際長度、經(jīng)緯度、時(shí)間等參數(shù)，確定視頻的拍攝地點(diǎn)和拍攝日期的一種方法。

而物體影子的變化主要源于地球自轉(zhuǎn)和繞太陽的公轉(zhuǎn)，直桿位置及長度的測量是利用CCD相機(jī)線性成像的交比不變性原理和滅點(diǎn)產(chǎn)生原理，所以通過分析影響影子變化的參數(shù)，就可確定視頻拍攝的地點(diǎn)和日期。

對于太陽影子定位技術(shù)考慮以下問題：

第一，建立影子長度變化的數(shù)學(xué)模型，分析影子長度關(guān)于各個(gè)參數(shù)的變化規(guī)律。

第二，當(dāng)日期已知，直桿高度未知時(shí)，根據(jù)水平地面上的太陽影子頂點(diǎn)坐標(biāo)數(shù)據(jù)，建立數(shù)學(xué)模型確定直桿所處的地點(diǎn)。

第三，當(dāng)日期與直桿高度都未知時(shí)，根據(jù)水平地面上的太陽影子頂點(diǎn)坐標(biāo)數(shù)據(jù)，建立數(shù)學(xué)模型確定直桿所處的地點(diǎn)和日期。

第四，根據(jù)一根2米高的直桿在太陽下的影子變化的視頻，建立確定視頻拍攝地點(diǎn)的數(shù)學(xué)模型，同時(shí)給出可能的拍攝地點(diǎn)。當(dāng)拍攝日期未知，根據(jù)數(shù)學(xué)模型確定可能的拍攝地點(diǎn)與日期。

1. 直桿影子長度變化曲線

1.1 問題分析

影響影子長度的主要因素是桿子高度和高度角。首先確立影響太陽高度角[1]變化的因素并確定高度角的計(jì)算公式，然后利用影子長度和太陽高度角的三角函數(shù)關(guān)系，得到影子長度變化規(guī)律。以天安門廣場北緯39度54分26秒,東經(jīng)116度23分29秒位置(文中記為P點(diǎn))，3米直桿在2015年10月22日北京時(shí)間9∶00-15∶00之間的太陽影子長度的變化曲線為例，說明這些量之間的變換關(guān)系。

1.2 確定影子長度與太陽高度角的關(guān)系

首先，我們確定太陽高度角[2]θ，收集資料可知太陽高度角θ的變化與太陽赤緯[3]δ、緯度φ和太陽時(shí)角[4]w三個(gè)因素有關(guān)：

(1)太陽赤緯δ與日期T有關(guān)，其關(guān)系式如下：

sinδ=0.39795cos[0.98563(T-173)]

(1)

其中T=1表示1月1日，T=365表示12月31日。

(2)緯度，其數(shù)值在0至90度之間。

(3)太陽時(shí)角w是指日面中心的時(shí)角，可以由北京時(shí)間得到，設(shè)北京時(shí)間為t,以24小時(shí)計(jì)，經(jīng)度為φ，可以得到關(guān)系式[3]：

w=15°×[t+(120°-φ)÷15°-12]

(2)

根據(jù)圖1，確定出了這三個(gè)因素關(guān)系式為[2]：

sinθ=sinφsinδ+cosφcosδcosw

(3)

最后建立函數(shù)模型，確立影子長度與太陽高度角θ的函數(shù)關(guān)系，如圖2所示，設(shè)直桿高度為h米，影子長度為l米，可以得到函數(shù)模型為：

(4)

圖1 太陽高度角計(jì)算示意圖

圖2 影子長度與太陽高度角的關(guān)系圖

1.3 畫出變化曲線圖

應(yīng)用建立的模型，確立了P點(diǎn)3米高的直桿在2015年10月22日北京時(shí)間9:00-15:00之間的太陽影子長度隨時(shí)間變化的二次函數(shù)關(guān)系，其中由計(jì)算可得φ=116.39,φ=39.9072,h=3,T=265，從而用Matlab畫出了北京時(shí)間9:00-15:00之間該地的太陽影子長度的變化曲線，如圖3所示。

分析可知，如果為東經(jīng)120度時(shí)，北京時(shí)間12點(diǎn)時(shí)影子長度最短，又因?yàn)槊客?分鐘，經(jīng)度偏西1度，而天安門廣場經(jīng)度為東經(jīng)116度23分29秒，比東經(jīng)120度少3.61度，所以該地的時(shí)間要比標(biāo)準(zhǔn)北京時(shí)間要晚14.44分鐘，所以天安門廣場影子最短的時(shí)間為12點(diǎn)14分26秒。

一天中影子變化的規(guī)律是早晚影子最長，中午最短，早上到中午影子慢慢變短，中午到晚上影子慢慢變長。從上面二次函數(shù)模型中可以得到：影子與太陽高度角的關(guān)系是，太陽高度角越小，影子越長。太陽高度角為45度時(shí)影子和物體一樣長。

2. 日期已知，直桿未知求地點(diǎn)

2.1 問題分析

已知日期時(shí)，獲得的某地桿子不同時(shí)間的頂點(diǎn)坐標(biāo)，確定拍攝地點(diǎn)所在的經(jīng)度和緯度，方法具體如下：

第一步，由已知數(shù)據(jù)得出影子長度，然后對這些數(shù)據(jù)進(jìn)行二次函數(shù)擬合，得出影子長度與時(shí)間的關(guān)系，從而根據(jù)影子最短距離的時(shí)間點(diǎn)確定可能點(diǎn)的經(jīng)度；

第二步，由獲得的經(jīng)度，根據(jù)太陽方位角[4]與時(shí)角，經(jīng)緯度，高度角之間的關(guān)系，求出高度角，然后根據(jù)影子和緯度的關(guān)系確定出可能的緯度，從而確定拍攝地點(diǎn)所在的經(jīng)度和緯度。

2.2 給定日期和固定直桿確定位置

2.2.1 確定經(jīng)度

對于該種情況，我們收集了2015年4月18日的一段時(shí)間的影子長度數(shù)據(jù)，得出影子長度與時(shí)間的關(guān)系如下表1所示。

表1 2015年4月18日某地影子長度與時(shí)間的關(guān)系

建立二次函數(shù)模型，利用Matlab對影子長度和時(shí)間t的數(shù)據(jù)進(jìn)行擬合分析，得到系數(shù)a=0.1473,b=3.7519,c=24.1325的拋物線，從而確定出影子長度與時(shí)間的關(guān)系為(如圖4)：

l=0.1473t2-3.7519t+24.1325

(5)

圖3 3米高直桿的太陽影子長度與時(shí)間的關(guān)系

圖4 表1的影子長度與時(shí)間的關(guān)系

由方程(5)可知，該地的影子最短的時(shí)間為拋物線最低點(diǎn)12.736，所以時(shí)間為12點(diǎn)44分9秒。比東經(jīng)120度的影子最短的時(shí)間晚44.16分鐘，又因?yàn)槊客?分鐘，經(jīng)度偏西1度，所以該地比東經(jīng)120度偏西11.04度，為東經(jīng)108.96度。

2.2.2 確定緯度

首先，由文獻(xiàn)[5]獲得的太陽方位角α與高度角θ、太陽赤緯δ、太陽時(shí)角w之間的關(guān)系確定出太陽方位角α的計(jì)算公式[4]：

(6)

圖5 4月18日的正午時(shí)刻影長與緯度的關(guān)系

由拋物線得到正午時(shí)刻的最短影長為0.241米，然后由緯度和影長的關(guān)系可以得出該地緯度為18.3度，由于2015年4月18日太陽直射點(diǎn)在北緯10度，又因表1的數(shù)據(jù)為下午的影子軌跡，可以得到直桿是在太陽直射點(diǎn)以北，所以應(yīng)該為北緯18.3度。

由以上分析可得到直桿所處的可能地點(diǎn)為東經(jīng)108.96度，北緯18.3度。

3. 日期與桿長都未知求地點(diǎn)

3.1 問題分析

對數(shù)據(jù)進(jìn)行二次函數(shù)擬合，根據(jù)影子最短距離的時(shí)間點(diǎn)比較分析確定兩組數(shù)據(jù)直桿所在地的經(jīng)度；對于日期不確定，遍歷一年365天，確定每一個(gè)T之后，得出日期，然后可得到正午時(shí)刻同一經(jīng)度下，緯度與影長的關(guān)系，從而確定緯度。

3.2 建立二次函數(shù)模型

利用Matlab對文獻(xiàn)[6]附件2和附件3中的兩組數(shù)據(jù)進(jìn)行分析處理，得出影子長度與時(shí)間關(guān)系的二次函數(shù)方程分別為：

l1=0.0981t2-2.8809t+21.731

l2=0.2964t2-7.5007t+47.5619

如圖6，圖7所示。

圖6 第一組影子長度與時(shí)間的關(guān)系

圖7 第二組影子長度與時(shí)間的關(guān)系

當(dāng)日期未知時(shí)，遍歷一年365天，即T∈[1,365]且T∈Z*,確定每一個(gè)T之后，用反演法確定出日期，第一組數(shù)據(jù)日期為7月20日，第二組數(shù)據(jù)日期為1月20日，可得積日T1=201,T2=20，然后可得到正午時(shí)刻同一經(jīng)度下，影長與緯度的關(guān)系分別為圖8，圖9所示。

圖8 7月20日正午時(shí)刻影長與緯度的關(guān)系

圖9 1月20日正午時(shí)刻影長與緯度的關(guān)系

所以可得該地點(diǎn)為東徑79.75度，北緯39.22度，日期為7月20日，另一個(gè)地點(diǎn)為東徑10.25度，北緯29.39度，日期為1月20日。

4. 日期未知，桿長已知求地點(diǎn)

4.1 問題分析

先將視頻導(dǎo)入Matlab，然后每2分鐘取一幀，得到21張有序的圖片，然后對每張圖片進(jìn)行灰度處理，得到每張圖片中影子頂點(diǎn)在矩陣中的位置，即像素點(diǎn)的坐標(biāo)，然后通過像素點(diǎn)坐標(biāo)及像素點(diǎn)間的距離和實(shí)際長度的比較，得到影子軌跡的實(shí)際坐標(biāo)，計(jì)算出經(jīng)、緯度，得出可能的拍攝地點(diǎn)。然后在日期未知的情況下，利用二次函數(shù)擬合模型求得可能的拍攝地點(diǎn)和日期。

5.2 根據(jù)桿長與時(shí)間求地點(diǎn)

4.2.1 確定經(jīng)度

先將視頻導(dǎo)入Matlab，并結(jié)合小孔成像原理，得到影子頂點(diǎn)軌跡的實(shí)際坐標(biāo)和實(shí)際影子長度與時(shí)間的關(guān)系如下表2所示。利用Matlab分析影子長度l和時(shí)間t的數(shù)據(jù)進(jìn)行二次函數(shù)擬合，得到系數(shù)：

a=0.0877,b=-2.4423,c=17.1244

從而得出影子長度與時(shí)間的關(guān)系拋物線方程為：

l=0.0877t2-2.4423t+17.1244

表2 2015年7月13日某地影子長度與時(shí)間的關(guān)系

用Matlab對數(shù)據(jù)做出影子長度與時(shí)間的關(guān)系符合拋物線，如圖10所示。

圖10 表2的影子長度與時(shí)間的關(guān)系

由方程可知，該地的影子最短的時(shí)間為拋物線最低點(diǎn)，時(shí)間為13點(diǎn)55分27秒。比東經(jīng)120度的影子最短的時(shí)間晚1小時(shí)55分27秒，又因?yàn)槊客?分鐘，經(jīng)度偏西1度，所以該地比東經(jīng)120度偏西28度52分，為東經(jīng)91.13度。

4.2.2 確定緯度

由上述二次函數(shù)確定出了最低點(diǎn)為13.9242，0.1209，可以得到在時(shí)間為13點(diǎn)55分27秒，經(jīng)度為東經(jīng)91度8分時(shí)，影子最短為0.1209米，而由題可知桿長h=2米，日期為2015年7月13日，可得積日T=194，然后由第一問影子長度與太陽高度角的關(guān)系，可得到正午時(shí)刻同一經(jīng)度下，緯度與影長的關(guān)系，如圖11所示。

圖11 表2的影長與緯度的關(guān)系

由拋物線得到正午時(shí)刻的最短影長為0.1209米，然后由緯度和影長的關(guān)系可以得出緯度為北緯26度和北緯20度，由于2015年7月13日太陽直射點(diǎn)在北緯23度，由于視頻中的影子往北面移動，可以得到直桿是在太陽直射點(diǎn)以北，所以應(yīng)該為北緯26度。

由模型可以得出該地點(diǎn)的經(jīng)度為東經(jīng)91.13度，北緯26度。

5.結(jié)束語

本文利用二次函數(shù)擬合方法研究太陽影子的定位技術(shù)，給出了確定視頻拍攝的地點(diǎn)和日期的一種方法。當(dāng)視頻中沒有直桿時(shí)，如何確定拍攝地點(diǎn)和拍攝日期是后續(xù)需要研究的一個(gè)方向。

[1] 屈名,王征兵,王德麾.基于交比不變性的太陽定位算法的研究[J].硅谷,2013(19).

[2] 費(fèi)云霞,王春順.對太陽高度角的了解及其計(jì)算方法[J].中小企業(yè)管理與科技,2008(1).

[3] 羅新宇.太陽赤緯的近似計(jì)算[J].鐵路航測,1993(2).

[4] 百度百科，太陽時(shí)角[EB/OL].[2015年9月12日]http://baike.baidu.com/link?url=h1Munt6Ctv0KtBIaY7-6HTmTf2w1uDOSQJ0bxzbBQnklF2ElnEQAQ8vC_0Z7raZVAUDhm5UyhPXHBdzhehmI-_

[5] 賀曉雷,于賀軍,李建英,等.太陽方位角的公式求解及其應(yīng)用[J].太陽能學(xué)報(bào),2008(1).

[6] 全國大學(xué)生數(shù)學(xué)建模網(wǎng), 2015年高教社杯全國大學(xué)生數(shù)學(xué)建模競賽賽題[EB/OL].2015年9月12日，http://www.mcm.edu.cn/html_cn/block/c61dfec317d7a5bd9b2b8efed81c8af3.html

【責(zé)任編輯 荊 瑤】

Solar Shadow Positioning Technology Based on Quadratic Function Fitting

MA Yu-xue, LINGHU Lin-yu, YANG Dong-hai, SUN Guo-wei, MAI A-li

(DepartmentofAppliedMathematics,YunchengUniversity,Yuncheng044000,China)

In this paper, solar shadow positioning technology was studied. Using the relationship between shadow length and sun elevation angle, the shadow length was obtained as long as the object when solar elevation angle is 45 degrees. By analysising the solar shadow length varying data by fitting equation, a quadratic function of the length of the shadow was obtained. According to the shortest time, this is shadow shortest time, to determine longitude. Then according to the relationship between the sun angle and azimuth, longitude latitude, elevation angle, the elevation angle is calculated, and the relationship between shadow and latitude are determined. Finally, filming locations and recording date of the video were determined.

Solar elevation angle; Solar azimuth; Quadratic function fitting; Solar shadow positioning technology

2015-08-16

山西省高等學(xué)校科技項(xiàng)目(20121113)；運(yùn)城學(xué)院科研基金項(xiàng)目(XK-2014030,XK-2014035)

馬玉雪(1994-)，女，山西潞城人，運(yùn)城學(xué)院應(yīng)用數(shù)學(xué)系1301班學(xué)生。

TP391.7

A

1008-8008(2015)06-0024-05