太陽影子定位模型

佘梓航，馬曉娜

（韓山師范學院數學與統計學院，廣東潮州521041）

太陽影子定位模型

佘梓航，馬曉娜

（韓山師范學院數學與統計學院，廣東潮州521041）

摘要：以2015年全國大學生數學建模競賽A題的附件1和附件4為對象.通過建立太陽影子定位的模型，分析確定影子對應物體所在的經緯度以及拍攝時間.在問題1中，利用直角三角形的性質，建立物體影子長度關于太陽赤緯、觀測點緯度、真太陽時和物體高度變化模型，在MATLAB環境下繪制影子長度與各參數的三維相關圖，得到其變化規律；在問題2中,基于問題1中建立的影子長度變化模型，結合真太陽時與觀測點經度的關系，建立影子長度關于觀測點緯度和觀測點經度關系模型，通過非線性擬合估計觀測點經緯度；在問題3中，將視頻數據進行轉化，獲取太陽影子頂點坐標數據,基于問題2中建立的影子長度變化模型，結合日期與太陽赤緯的關系，通過非線性擬合估計視頻拍攝的可能地點.

關鍵詞：控制變量法；數據可視化；非線性最小二乘；非線性擬合；真太陽時

太陽影子定位技術就是通過分析物體的太陽影子變化，確定視頻拍攝地點經緯度和日期的一種方法.如何確定拍攝地點和拍攝日期是數據分析的重要方面，本文采用2015年全國大學生數學建模競賽A題附件1和附件4的數據，嘗試解決以下問題.

問題1：建立影子長度變化的數學模型，分析影子長度關于各個參數的變化規律，并畫出2015年10 月22日北京時間9:00-15:00之間天安門廣場（北緯39°54'26"，東經116°23'29"）3 m高的直桿的太陽影子長度的變化曲線.

問題2：根據附件1中某固定直桿在水平地面上的太陽影子頂點坐標數據，建立數學模型確定直桿所處的地點.

問題3：根據附件4中一根2 m高直桿在太陽下的影子變化的視頻，建立數學模型確定若干個可能的視頻拍攝地點.

注：題目及數據附件都可到全國大學生數學建模競賽網站（www.mcm.edu.cn）下載.

1模型的假設

（1）假設物體較規則且物體直立；

（2）拍攝當天天氣狀況良好，受太陽光持續照射；

（3）物體四周無其余物體遮擋，保證影子軌跡不受外界影響；（4）不考慮空氣折射和太陽視半徑.

2模型的建立與求解

2.1問題1分析與建模求解

問題1要求建立影子長度變化的數學模型，分析影子長度關于各個參數的變化規律，并繪制太陽影子長度的變化曲線.本文利用控制變量法及數據可視化處理分析太陽影子長度與各相關參數的變化關系建立模型1，繪制2015年10月22日北京時間9:00-15:00之間天安門廣場3 m高的直桿的太陽影子長度的變化曲線.

2.1.1模型1的建立

太陽高度角的正弦值［1］計算如下：

根據物體影子長度與太陽高度角的余切關系，影子長度公式為：

根據（2）式，計算各個參數，最終得到影子長度變化的余切模型：

其中L為物體影子長度，H為物體的高度或長度，φ表示當地緯度，δ為太陽赤緯角，Ω為太陽時角.

2.1.2模型1的求解

保持太陽赤緯角δ不變，取3 m高的較規則物體，并假設其直立.結合（3）式，利用MATLAB編程，繪制該物體的影子長度在太陽赤緯角δ為0度時，隨太陽時角和觀測點緯度變化的三維圖，如圖1所示.

由圖1觀測到：影子長度在日出到日落的變化趨勢先是由長變短再變長，以太陽時角為0度即真太陽時為12點作為分界點，長度最短.利用（3）式，將題目中所給數據代入對應公式，得到天安門廣場2015年10月22日北京時間9:00-15:00之間天安門廣場（北緯39°54'26"，東經116°23'29"）3 m高的直桿的太陽影子長度，繪制太陽影子長度與真太陽時的曲線圖.

圖1　影子長度隨太陽時角和觀測點緯度變化的三維圖

在（3）式中，求影子長度需要得到太陽赤緯角δ和太陽時角Ω，其計算步驟如下.

（1）赤緯角的計算

根據有關文獻，太陽赤緯角的計算公式［2］為：

其中，θ表示日角，計算公式為：

此處t由兩部分組成，即：

式中N為積日，積日表示日期在年內的順序號.且：

INT（X）為取不大于X的最大整數，Y為年份.

（2）太陽時角的計算

規定正午太陽時角Ω為0，上午Ω為負值，下午Ω為正值.由于地球自轉一周360°需要24小時，每小時對應時角為15°.太陽時角的計算公式為：

其中tO為真太陽時，以24小時計.在某地要取得正確的當地真太陽時tO，由北京時間來換算需要修正兩種差.

由于各地經度不同而存在時間差.用Bt表示北京時間，用Et表示時差，那么平太陽時tm為：

其中時差與當地經度有關，其計算公式為：

其中λ表示當地經度，時差Et的單位為小時.

平太陽時tm假設地球繞太陽是標準的圓形，即一年中每天都是24小時平均的.但地球的公轉軌道是一個橢圓，也就是說地球的公轉并不均勻，并非每天都為24小時.因此要將平太陽時調整為真太陽時，需要修正真太陽時與平太陽時之差η.時差只與日期相關，與地區所處位置無關.即：

其中：

式中，t表示每年自1月1日后的天數，其中不包括1月1日.η單位為分鐘.

根據上述公式，且真太陽時為平太陽時與真平太陽時差相加之和［3］，得到真太陽時為：

結合（3）式，應用以上的求解步驟，利用MATLAB編程實現，繪制2015年10月22日北京時間9:00-15:00之間天安門廣場3 m高的直桿的太陽影子長度的變化曲線.以3 min為時間間隔取點，計算120個時間點各個參數的值，最終得到直桿影子長度隨時間變化的圖像如圖2所示.

由圖2曲線可知，天安門廣場3 m高的直桿的太陽影子在指定的時間段北京時間9:00-15:00內，隨著時間變化，影子長度總體呈現出先變短后變長的趨勢，在9:00時，影子長度值達到最高點為6.438 6 m；接著長度值逐漸減少，約在12:00處到達最低點；12:00過后長度值曲線回升，最后在15:00處.

2.2問題2分析與建模求解

問題2要求根據某固定直桿在水平地面上的太陽影子頂點坐標數據（附件1），建立數學模型確定直桿所處地的緯經度.由（3）式可知影子長度與參數桿長H，緯度φ，太陽赤緯角δ，真太陽時tO存在確定的非線性關系且真太陽時與觀測點經度存在確定關系，此問題可以轉化為非線性最小二乘估計問題.已知數據的測量時間為2015年4月18日，可以算出太陽赤緯.利用非線性最小二乘法［4］擬合估計出觀測點的經度和緯度，再根據經緯度在地圖上查詢直桿的具體位置.

2.2.1模型2的建立

真太陽時與經度關系為：

圖2　直桿影子長度曲線變化圖

其中Bt表示北京時間，η為真平太陽時時差，且：

太陽時角與真太陽時關系：

將（14）和（15）式代入（16）式可得太陽時角與經度的關系式：

要使擬合殘差平方之和最小，即目標函數為：

其中n=2，Li為直桿真實影長，為每一次估計的直桿擬合影長，E表示殘差平方之和計算如下：

綜合以上結果，可以得到使得殘差平方之和最小的最優模型：

2.2.2模型2的求解

Step2：利用MATLAB編程調用非線性擬合函數lsqnonlin擬合參數桿長L，緯度φ，觀測點經度λ的估計值，且設置的參數初始值，可得到多個擬合結果.

Step3：繪制擬合圖，利用殘差檢驗擬合優度.若擬合效果差，應返回Step2，重新擬合；擬合效果好，通過檢驗.仍返回Step2重新設置區間，得到多組可能的觀測點經緯度，如表1所示.

表1　直桿可能的所在地

根據本文計算的數據所得的擬合圖如圖3所示.由圖3可知，影長真實值和估計值基本吻合，擬合優度較好.

2.3問題3分析與建模求解

本題要求根據一根2 m高的直桿在水平地面上的太陽影子變化的視頻（附件4），建立數學模型確定拍攝地點經緯度.本題給的是視頻數據，則需要將視頻的數據進行轉化，獲取太陽影子頂點坐標數據，建立太陽高難度角相對于觀測點緯度和經度的變化模型，其中3個參數分別為：太陽高度角h，觀測點緯度φ，觀測點經度λ，可以利用非線性最小二乘法擬合估計，得到經度和緯度數據，從而直桿所在地得以確定.

圖3　影子長度與真實值的擬合圖

2.3.1模型3的建立

由于桿長和時間是已知的，太陽赤緯可以根據時間代入公式（4）～（7）求出；通過讀取視頻處理得到的坐標數據和已知的桿長，可以得到太陽高度角，結合（1）式和（17）式得到太陽高度關于觀測點緯度與經度的關系式為：

此問中利用太陽高度角和觀測點緯度的非線性最小二乘估計，要使擬合殘差之和最小，則有基本模型：

2.3.2模型3的求解

Step1：利用MATLAB軟件中的VideoReader函數，將視頻轉化為圖片，由于視頻的時長約為40 min，且附件1～3中均以21個數據點作為分析數據.不失一般性，每間隔2 min取圖片，從而獲取21張圖片.

Step2：在MATLAB環境下，調用rgb2gray編程實現將Step1獲取的21張真彩色圖像轉換為灰度圖像，記錄下每一時間間隔桿頂影子的像素坐標a（x1，y1），桿頂的像素坐標b（x2，y2），桿底的像素坐標c（x3，y3）.

Step3：將Step2得到的桿頂影子像素坐標a（x1，y1）轉化為以桿底像素坐標為原點的像素坐標（x1-y3，y1-y3），對附件4中的數據進行處理，結果如表2所示.

表2　桿頂點陰影的像素坐標數據　（像素）

Step4：根據桿頂點與底端的像素坐標，其縱坐標之差為d=y3-y2，結合題目中給出桿的高度為2 m，則可將桿頂影子像素坐標轉化為以米為單位的坐標，即可得桿影的頂點坐標數據，將附件4中的數據代入，以桿底為原點坐標，結果如表3所示.

表3　桿影頂點坐標數據 （m）

類比模型2，對附件4中視頻數據轉化后，結合 （21）式，通過MATLAB編程實現，給出直桿可能的地點，見表4.

表4　拍攝地點經緯度

2.3.3誤差分析

通過分析得到的影子坐標計算影子長度L，則在其基礎上添加正態隨機數Z（0，0.005，1，21），即影子長度"，再進行非線性擬合求解，對擬合的結果進行白噪聲檢驗［5］，影子長度添加正態隨機數后，R值（0.999）較為接近，且地點多為內蒙古，則說明通過視頻讀取的數據具有合理性.

3結語

對于問題1，本文給出了影子長度同各個參數之間的變化規律關系圖，為解答接下去的問題提供依據.在畫出影子長度變化圖時，能夠考慮到晚上沒有影子這個事實，畫出來的圖更具有代表性.同時，利用非線性最小二乘法給出的結果，接近數模組委會給出來的結果，正確率高.對于最后一個問題，由于視頻中給出的影子長度受到很多不確定性因素的影響（如：拍攝角度，光的折射等），因此，最后一個問題不可能給出準確的結果.即便如此，為了驗證本文所給方法的正確性，本文通過對影子坐標數據加入白噪聲，再進行非線性擬合求解，得到同沒加白噪聲相近的擬合精度，說明了從視頻中讀取數據具有合理性.

參考文獻：

［1］鄭鵬飛，林大鈞，劉小羊，等.基于影子軌跡線反求采光效果的技術研究［J］.華東理工大學學報，2010，36（3）：459-460.

［2］陳曉勇，鄭科科.對建筑日照計算中太陽赤緯角公式的探討［J］.浙江建筑，2011，28（9）：7-8.

［3］孫吉山，洪薇.北京標準時間與各地真太陽時的換算法［J］.上海針灸雜志，1994，13（1）：25-26.

［4］薛毅.最優化理論與算法［M］.北京：北京工業大學出版社，2004：201-202.

［5］王燕.應用時間序列分析［M］.3版.北京：中國人民大學出版社，2012：27-31.

（責任編輯：邵曉軍）

中圖分類號：O29

文獻標識碼：A

文章編號：1007-5348（2016）04-0006-07

［收稿日期］2016-03-12

［作者簡介］佘梓航（1989-），男，廣東潮州人，韓山師范學院數學與統計學院助教，碩士；研究方向：高頻數據估計量構建與數學建模.

Position Models of the Sun Shadow

SHE Zi-hang,MA Xiao-na
（DepartmentofMathematicsandStatistics,HanshanNormal University,Chaozhou521041,Guangdong,China）

Abstract:In this paper,a solution is presented for the ProblemA of Contemporary Undergraduate Mathematical Contest in Modeling in 2015.In question 1,by using the property of right triangle,it built the sun shadow model about the solar declination,latitude,apparent solar time and the length of object.It used MATLAB todrawthe 3D correlation diagramof the parameters and shadowlength.In question 2,based on the model of question 1 and the relationshipbetweentheapparent solar timeandthelongitudeof theobservationpoint,it built lengthof sunshadow model about the latitude and longitude by using the observation point.Therefore,it estimated the latitude and longitudeby NonlinearFittingMethod（NFM）.Inquestion3,itusedthevideodatatoobtainthevertex coordinatesof the sun shadow.Based on and the relationship between the date and the solar declination,it estimated the possible locationsofthevideobyNFM.

Key words:control variable method;data visualization;nonlinear least square method;nonlinear fitting method （NFM）;solardeclination