屈耀紅, 邢哲文, 袁冬莉, 張友民

(1.西北工業(yè)大學(xué) 自動(dòng)化學(xué)院, 陜西 西安　710072; 2.康科迪亞大學(xué) 機(jī)械與工業(yè)工程系, 蒙特利爾，QC H3G 1M8, 加拿大)

?

基于懸停四旋翼位置姿態(tài)信息的風(fēng)場(chǎng)估計(jì)方法研究

屈耀紅1, 邢哲文1, 袁冬莉1, 張友民2

(1.西北工業(yè)大學(xué) 自動(dòng)化學(xué)院, 陜西 西安710072; 2.康科迪亞大學(xué) 機(jī)械與工業(yè)工程系, 蒙特利爾，QC H3G 1M8, 加拿大)

摘要：基于四旋翼飛行器懸停狀態(tài)下的位置及姿態(tài)信息，提出了一種離線的風(fēng)場(chǎng)估計(jì)方法。首先根據(jù)Dryden大氣紊流模型建立了四旋翼飛行器所處的風(fēng)場(chǎng)環(huán)境，并通過(guò)分析有風(fēng)情況下旋翼升力的變化，得到旋翼升力與風(fēng)場(chǎng)信息(風(fēng)速、風(fēng)向)的函數(shù)關(guān)系式；接著利用牛頓-歐拉方法推導(dǎo)出有風(fēng)擾動(dòng)下的四旋翼動(dòng)力學(xué)方程，并進(jìn)一步設(shè)計(jì)了用于保持飛行器懸停狀態(tài)的PID控制器；最后，基于懸停狀態(tài)下四旋翼飛行器的位置姿態(tài)信息，計(jì)算得到飛行器所處的風(fēng)場(chǎng)環(huán)境信息。MATLAB仿真結(jié)果表明所提方法在有紊流干擾的情況下，能夠有效地提取出風(fēng)場(chǎng)環(huán)境里的主風(fēng)信息。

關(guān)鍵詞：控制器;MATLAB;四旋翼飛行器;旋翼升力;Dryden模型;定點(diǎn)懸停;風(fēng)場(chǎng)估計(jì)

在未知近海岸低空區(qū)域中，準(zhǔn)確的風(fēng)場(chǎng)估計(jì)對(duì)該區(qū)域里風(fēng)能的探測(cè)、開發(fā)和利用有著至關(guān)重要的作用[1]。在這樣的應(yīng)用背景下，傳統(tǒng)陸基測(cè)風(fēng)儀器設(shè)備的昂貴和平臺(tái)依賴性等弊端就體現(xiàn)了出來(lái)。相反的，以無(wú)人飛行器(UAV)為主要平臺(tái)的風(fēng)場(chǎng)估計(jì)方法因其廉價(jià)可回收、方便靈活等特點(diǎn)，成為了對(duì)傳統(tǒng)風(fēng)環(huán)境估計(jì)方法的良好替代[2]。

根據(jù)所選用的飛行器類型，無(wú)人機(jī)風(fēng)場(chǎng)估計(jì)方法可分為固定翼無(wú)人機(jī)風(fēng)場(chǎng)估計(jì)方法和旋翼無(wú)人機(jī)風(fēng)場(chǎng)估計(jì)方法。雖然以固定翼無(wú)人機(jī)[3]為平臺(tái)的風(fēng)場(chǎng)估計(jì)方法具有測(cè)量范圍大、續(xù)航時(shí)間長(zhǎng)等特點(diǎn)，但是其無(wú)法做到對(duì)風(fēng)信息的定點(diǎn)測(cè)量。并且因?yàn)楣潭ㄒ頍o(wú)人機(jī)需要保持較高的空速以保證其正常的飛行，所以低風(fēng)速和小范圍的風(fēng)信息變化因占比小而無(wú)法被測(cè)量[4]。

以四旋翼無(wú)人機(jī)為平臺(tái)的風(fēng)場(chǎng)估計(jì)方法能夠很好的彌補(bǔ)固定翼無(wú)人機(jī)測(cè)風(fēng)的不足。首先，四旋翼飛行器能夠保持定點(diǎn)懸停狀態(tài)，從而實(shí)現(xiàn)風(fēng)信息的定點(diǎn)測(cè)量，得到固定位置氣流信息(風(fēng)速和風(fēng)向)隨時(shí)間的變化關(guān)系。其次，由于定點(diǎn)懸停狀態(tài)下旋翼的氣動(dòng)特性造成飛行器橫側(cè)向擾流較小，低風(fēng)速和小范圍的風(fēng)信息變化的測(cè)量精度得以提高。利用四旋翼無(wú)人機(jī)的風(fēng)場(chǎng)估計(jì)方法有2個(gè)關(guān)鍵的問(wèn)題亟待解決。第一是如何實(shí)現(xiàn)有風(fēng)條件下四旋翼飛行器的懸停保持。四旋翼作為一個(gè)典型的欠驅(qū)動(dòng)強(qiáng)耦合非線性系統(tǒng)，能夠?qū)崿F(xiàn)六自由度的飛行運(yùn)動(dòng)。在有風(fēng)擾動(dòng)情況下，飛行器運(yùn)動(dòng)的數(shù)學(xué)模型會(huì)發(fā)生變化，簡(jiǎn)單的控制方法難以實(shí)現(xiàn)理想的懸?？刂菩Ч?。第二是如何從四旋翼飛行器的飛行數(shù)據(jù)中提取有效部分，從而計(jì)算得到風(fēng)場(chǎng)環(huán)境信息。風(fēng)擾動(dòng)與四旋翼的氣動(dòng)耦合關(guān)系不是簡(jiǎn)單的線性疊加，很難顯式的求解得到風(fēng)場(chǎng)信息的表達(dá)式。

針對(duì)上述問(wèn)題，本文提出了一種基于懸停四旋翼飛行器位置姿態(tài)信息的風(fēng)場(chǎng)估計(jì)方法。該方法通過(guò)六通道PID控制器實(shí)現(xiàn)了四旋翼在風(fēng)擾情況下的懸?？刂?，并利用已知的四旋翼氣動(dòng)參數(shù)，推導(dǎo)出風(fēng)場(chǎng)信息的近似計(jì)算方程。最后利用MATLAB/M-File編程對(duì)本文所提出風(fēng)場(chǎng)估計(jì)方法的有效性進(jìn)行了仿真測(cè)試。

1四旋翼飛行器動(dòng)力學(xué)模型

1.1旋翼空氣動(dòng)力分析

有風(fēng)情況下，旋翼的空氣動(dòng)力學(xué)分析如圖1所示。

圖1　旋翼氣動(dòng)力分析圖

圖中,v1表示由旋翼旋轉(zhuǎn)而產(chǎn)生的誘導(dǎo)速度,vw表示風(fēng)速,vt表示誘導(dǎo)速度和風(fēng)速的合速度,α表示風(fēng)與旋翼旋轉(zhuǎn)平面形成的角度,風(fēng)自上而下流過(guò)旋翼旋轉(zhuǎn)平面角度為正,FT表示旋翼的升力,沿yD軸正方向?yàn)檎?與v1方向相反。

根據(jù)旋翼滑流理論,無(wú)風(fēng)條件下旋翼升力可以表示為

(1)

式中,ρ表示空氣密度,A表示槳盤的有效面積。

在無(wú)風(fēng)條件下,旋翼的升力和扭矩均與轉(zhuǎn)速的平方成正比[5]

(2)

式中,kF、kM分別表示旋翼的升力比和扭矩比。

旋翼的誘導(dǎo)速度可以表示為

(3)

式中,R表示槳葉半徑,槳盤的有效面積近似與槳葉的旋轉(zhuǎn)面積相等。

在有風(fēng)場(chǎng)情況下,旋翼的總升力[6]可表示為

(4)

根據(jù)旋翼葉素理論,經(jīng)過(guò)簡(jiǎn)化的旋翼升力和扭矩公式可表示為

(5)

式中,CT、CM分別表示旋翼的升力系數(shù)和扭矩系數(shù),二者僅與旋翼的形狀和槳葉數(shù)量有關(guān)。

考慮圖1中v1、vw與vt的幾何關(guān)系,由余弦定理

(6)

(4)式可以表示為

(7)

1.2四旋翼飛行器建模

假設(shè)四旋翼飛行器是由剛體機(jī)身和4個(gè)理想旋翼組成。4個(gè)旋翼被分成2組,旋翼1和旋翼3的旋轉(zhuǎn)方向與zb軸正方向一致,旋翼2和旋翼4的旋轉(zhuǎn)方向與zb軸的正方向相反,用以抵消旋翼1和旋翼3因旋轉(zhuǎn)產(chǎn)生的反扭矩。

圖2　四旋翼飛行器坐標(biāo)系示意圖

通過(guò)控制旋翼2和旋翼4的轉(zhuǎn)速差,能夠?qū)崿F(xiàn)機(jī)體繞xb軸的滾轉(zhuǎn)運(yùn)動(dòng)和沿yb的平動(dòng);通過(guò)控制旋翼1和旋翼3的轉(zhuǎn)速差,能夠?qū)崿F(xiàn)機(jī)體繞yb軸的俯仰運(yùn)動(dòng)和沿xb軸的平動(dòng);通過(guò)控制旋翼1、2、3、4的轉(zhuǎn)速不同,能夠?qū)崿F(xiàn)機(jī)體繞zb軸的偏航運(yùn)動(dòng)和沿zb軸的平動(dòng)。

在導(dǎo)航坐標(biāo)系下,根據(jù)牛頓-歐拉方法,對(duì)飛行器進(jìn)行受力和運(yùn)動(dòng)分析。

(8)

根據(jù)空氣阻力的計(jì)算公式

(9)

(10)

在機(jī)體坐標(biāo)系下,建立四旋翼飛行器旋轉(zhuǎn)方程[7]

(11)

(12)

(13)

式中,l表示四旋翼的臂長(zhǎng),kyaw是偏航系數(shù)。

由機(jī)體坐標(biāo)系下的傾轉(zhuǎn)角速度到導(dǎo)航坐標(biāo)系下的姿態(tài)角速度的變換矩陣如下

(14)

綜合(8)～(14)式,得到有風(fēng)擾動(dòng)下的四旋翼飛行器的運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)模型

(15)

2紊流風(fēng)場(chǎng)環(huán)境建模

大氣擾動(dòng)中包含各種時(shí)間和空間尺度的運(yùn)動(dòng),其產(chǎn)生的機(jī)理和發(fā)展過(guò)程各不相同。采用復(fù)雜的大氣動(dòng)力學(xué)方程來(lái)研究大氣環(huán)境對(duì)飛行影響的規(guī)律是不方便的。因而廣泛使用簡(jiǎn)化大氣擾動(dòng)模型,它主要描述風(fēng)場(chǎng)基本物理參數(shù)之間的關(guān)聯(lián),而忽略一些次要的影響因素。本文中采用Dryden大氣紊流模型,將標(biāo)準(zhǔn)高斯白噪聲序列通過(guò)成形濾波器,形成有色噪聲序列,完成對(duì)大氣紊流的模擬。

根據(jù)文獻(xiàn)[8],Dryden模型的時(shí)間譜函數(shù)為

(16)

式中,σx、σy和σz表示3個(gè)方向的紊流強(qiáng)度,Lx、Ly和Lz表示3個(gè)方向的紊流尺度,v表示飛行器的飛行速度。

由于四旋翼飛行器的飛行高度有限,根據(jù)文獻(xiàn)[9],在低空條件下,紊流強(qiáng)度和紊流尺度可由(17)、(18)式計(jì)算

(17)

(18)

式中,h表示飛行高度,u20表示6.096 m高度的風(fēng)速。

將單位強(qiáng)度的白噪聲序列r(t)通過(guò)一個(gè)傳遞函數(shù)為G(s)的濾波器,產(chǎn)生輸出序列x(t),輸出頻譜函數(shù)為

(19)

式中,*表示共軛。將(16)式按照(19)式分解,對(duì)于3個(gè)紊流速度,可以得到所需的成形濾波器的傳遞函數(shù),為方便實(shí)現(xiàn),均簡(jiǎn)化為一階形式

(20)

以h0作為步長(zhǎng)(即時(shí)間間隔),將成形濾波器傳遞函數(shù)離散化,生成x(t)的離散序列

x(i+1)=Px(i)+Qr(i+1)

(21)

3基于懸停過(guò)程的風(fēng)場(chǎng)估計(jì)

3.1懸停控制器設(shè)計(jì)

本文采用六通道的PID控制器,實(shí)現(xiàn)了四旋翼飛行器的懸?？刂?。將4個(gè)旋翼的轉(zhuǎn)速ω1、ω2、ω3、ω4定義為系統(tǒng)的輸入量。為了簡(jiǎn)化控制器的設(shè)計(jì)過(guò)程,做出如下假設(shè):

1) 槳盤各處的誘導(dǎo)速度相等,忽略旋翼旋轉(zhuǎn)的邊緣效應(yīng),不考慮因旋翼工作而產(chǎn)生的擾流的相互影響。

2) 與飛行器的傾轉(zhuǎn)扭矩相比,有旋翼的旋轉(zhuǎn)力矩是小量,可以忽略。

3) 與四旋翼飛行器的質(zhì)量相比,旋翼的質(zhì)量為小量,忽略旋翼的轉(zhuǎn)動(dòng)慣量。

4) 本文只考慮小角度的四旋翼運(yùn)動(dòng),因此,導(dǎo)航坐標(biāo)系下的姿態(tài)角速度與機(jī)體坐標(biāo)系下的傾轉(zhuǎn)角速度近似相等。

控制系統(tǒng)結(jié)構(gòu)如圖3所示。

圖3　懸?？刂葡到y(tǒng)結(jié)構(gòu)示意圖

其中,xd、yd、zd、φd、θd、ψd分別表示期望控制四旋翼飛行器達(dá)到的位置和姿態(tài)。T、Mφ、Mθ、Mψ分別表示為達(dá)到期望懸停狀態(tài)而計(jì)算得到的總升力、滾轉(zhuǎn)力矩、俯仰力矩和偏航力矩。

該控制系統(tǒng)分為2個(gè)部分:內(nèi)環(huán)的姿態(tài)控制和外環(huán)的位置控制。首先利用3個(gè)通道的PID控制器將期望位置信息(xd,yd,zd)和反饋得到的當(dāng)前位置信息(x,y,z)送入逆解算單元,計(jì)算得到合升力T,期望滾轉(zhuǎn)角φd和期望俯仰角θd;然后利用期望姿態(tài)信息(φd,θd,ψd)和當(dāng)前姿態(tài)信息(φ,θ,ψ),通過(guò)另外3個(gè)通道的PID控制器,計(jì)算得到目標(biāo)滾轉(zhuǎn)力矩Mφ、俯仰力矩Mθ和偏航力矩Mψ;最后利用之前計(jì)算得到的總升力和3個(gè)方向的力矩(T,Mφ,Mθ,Mψ)通過(guò)轉(zhuǎn)速分配模塊,計(jì)算得到系統(tǒng)的輸入量(ω1,ω2,ω3,ω4)。

① 位置控制器設(shè)計(jì)

通道1、通道2、通道3的PID控制器分別用作計(jì)算控制四旋翼飛行器沿導(dǎo)航坐標(biāo)系下3個(gè)坐標(biāo)軸的加速度。數(shù)學(xué)表達(dá)式如下

(22)

式中,kpx、kix、kdx、kpy、kiy、kdy、kpz、kiz、kdz是這3個(gè)通道的PID控制器的控制參數(shù)。

(23)

(24)

② 姿態(tài)控制器設(shè)計(jì)

通道4、通道5、通道6的PID控制器分別用作計(jì)算控制四旋翼飛行器繞機(jī)體坐標(biāo)系下3個(gè)坐標(biāo)軸的角加速度。數(shù)學(xué)表達(dá)式如下

(25)

式中,kpφ、kiφ、kdφ、kpθ、kiθ、kdθ、kdθ、kiψ、kdψ是這3個(gè)通道的PID控制器的控制參數(shù)。

考慮到姿態(tài)角對(duì)整個(gè)控制系統(tǒng)非線性性質(zhì)造成的影響,在控制器設(shè)計(jì)時(shí),對(duì)滾轉(zhuǎn)角和俯仰角進(jìn)行限幅,-10°≤φ≤10°,-10°≤θ≤10°。由于在姿態(tài)控制環(huán)中引入積分項(xiàng)會(huì)引起系統(tǒng)超調(diào),造成控制系統(tǒng)超越限幅和不穩(wěn)定,所以,令kiφ=0,kiθ=0,kiψ=0,3個(gè)通道的PID控制退化為PD控制。

繞3個(gè)坐標(biāo)軸的旋轉(zhuǎn)力矩的數(shù)學(xué)表達(dá)式如下

(26)

3.2基于懸停過(guò)程的風(fēng)場(chǎng)估計(jì)方法

四旋翼飛行器的懸停過(guò)程總是要保持導(dǎo)航坐標(biāo)系下3個(gè)方向的位移X=0,偏航角ψ=0。根據(jù)(15)式,可以得到風(fēng)擾力的表達(dá)式

(27)

考慮到四旋翼飛行器的低風(fēng)擾動(dòng)飛行環(huán)境和小角度懸停控制過(guò)程,橫側(cè)向方程可做如下簡(jiǎn)化

(28)

根據(jù)(9)式,飛行器橫側(cè)向分解空速的表達(dá)式為

(29)

又由(10)式,橫側(cè)向風(fēng)場(chǎng)的估計(jì)表達(dá)式為

(30)

4仿真與分析

4.1紊流風(fēng)場(chǎng)仿真

取6.096m高度的風(fēng)速為5m/s,四旋翼飛行器的飛行高度為10m。紊流風(fēng)場(chǎng)的參數(shù)如表1所示,生成的紊流風(fēng)場(chǎng)如圖4所示。

表1　紊流風(fēng)場(chǎng)參數(shù)

圖4　紊流風(fēng)場(chǎng)圖

在xn軸方向紊流的基礎(chǔ)上疊加了一個(gè)大小4 m/s,方向與xn軸正方向相反的常值風(fēng);在yn軸方向紊流的基礎(chǔ)上疊加了一個(gè)大小3 m/s,方向與yn軸正方向相同的常值風(fēng)。

4.2懸?？刂品抡?/p>

模型參數(shù)如表2所示,PID控制器參數(shù)如表3所示?？刂菩Ч鐖D5、圖6所示。

表2　模型參數(shù)

表3　PID控制器參數(shù)

由仿真結(jié)果圖可看出,PID懸?？刂破魅〉昧死硐氲目刂菩ЧＮ灰破屏枯^小,角度偏移量小于10°。

4.3風(fēng)場(chǎng)環(huán)境估計(jì)

利用仿真得到的有風(fēng)擾動(dòng)下的四旋翼懸停飛行數(shù)據(jù),結(jié)合之前提出的風(fēng)場(chǎng)估計(jì)方法,對(duì)四旋翼的飛行環(huán)境重新估計(jì)。仿真結(jié)果如圖7所示。

圖5　四旋翼懸停橫側(cè)向位移圖　　　圖6　四旋翼懸停滾轉(zhuǎn)角/俯仰角偏移圖　　　圖7　風(fēng)場(chǎng)估計(jì)效果圖

圖8　估計(jì)風(fēng)場(chǎng)環(huán)境誤差分析圖

比較圖4和圖7,估計(jì)得到的風(fēng)場(chǎng)環(huán)境信息與仿真風(fēng)場(chǎng)環(huán)境信息具有相同的變化趨勢(shì)。誤差分析如圖8所示,xn方向與yn方向的估計(jì)誤差均小于1 m/s。造成誤差的原因主要是因?yàn)槔?28)式對(duì)旋翼升力做了近似化處理,而風(fēng)擾動(dòng)對(duì)旋翼升力的影響不可忽略。

5結(jié)論

針對(duì)四旋翼的定點(diǎn)懸停動(dòng)力學(xué)特性,本文提出了基于懸停四旋翼位置及姿態(tài)信息的風(fēng)場(chǎng)環(huán)境估計(jì)方法。仿真驗(yàn)證了在低風(fēng)速擾動(dòng)情況下,該方法的可行性。

此研究工作考慮了風(fēng)擾對(duì)四旋翼飛行器氣動(dòng)建模的影響,推導(dǎo)得到了估計(jì)風(fēng)場(chǎng)環(huán)境的近似方程,為無(wú)人機(jī)探測(cè)風(fēng)場(chǎng)的研究提供了一種新的思路。

參考文獻(xiàn)：

[1]Fichaux N, Beurskens J, Jensen P H, et al. Upwind: Design Limits and Solutions for Very Large Wind Turbines[R]. Brussels, Sixth Framework Programme, 2011

[2]Moyano Cano J. Quadrotor UAV for Wind Profile Characterization[D]. Madrid, Universidad Carlos Ⅲ de Madrid, 2013

[3]屈耀紅, 于麗麗, 周迎,等. 小型固定翼無(wú)人機(jī)多步長(zhǎng)在線風(fēng)場(chǎng)估計(jì)[J]. 西北工業(yè)大學(xué)學(xué)報(bào), 2015, 33(5): 750-756

Qu Yaohong, Yu Lili, Zhou Ying, et al. Multi-Step Wind Field Estimation Online of Small Fixed-Wing UAVs[J]. Journal of Northwestern Polytechnical University, 2015, 33(5): 750-756 (in Chinese)

[4]Velasco-Carrau J, García-Nieto S, Salcedo J V, et al. Multi-Objective Optimization for Wind Estimation and Aircraft Model Identification[J]. Journal of Guidance, Control, and Dynamics, 2015, 39: 1-18

[5]Nagaty A, Saeedi S, Thibault C, et al. Control and Navigation Framework for Quadrotor Helicopters[J]. Journal of intelligent & robotic systems, 2013, 70(1/2/3/4): 1-12

[6]Tran N K, Bulka E, Nahon M. Quadrotor Control in a Wind Field[C]∥2015 International Conference on Unmanned Aircraft Systems, 2015: 320-328

[7]He Y, Chen Y, Zhou M. Modeling and Control of a Quadrotor Helicopter Under Impact of Wind Disturbance[J]. Journal of Chinese Inertial Technology, 2013, 5: 013

[8]Beal T R. Digital Simulation of Atmospheric Turbulence for Dryden and Von Karman Models[J]. Journal of Guidance, Control and Dynamics, 1993, 16(1): 132-138

[9]肖葉倫,金長(zhǎng)江. 大氣擾動(dòng)中的飛行原理[M]. 北京:國(guó)防工業(yè)出版社,1993: 73-74

XiaoYelun, Jin Changjiang. Flight theory in Atmospheric Disturbance[M]. Beijing, National Defense Industry Press,1993: 73-74 (in Chinese)

Wind Field Estimation Based on Position and AttitudeInformation of Quadrotor in Hover

Qu Yaohong1, Xing Zhewen1, Yuan Dongli1, Zhang Youmin2

1.School of Automation, Northwestern Polytechnical University, Xi′an 710072, China2.Department of Mechanical and Industrial Engineering, Concordia University, Montreal, QC H3G 1M8, Canada

Abstract:A wind field estimation method is proposed based on the position and attitude information of a quadrotor in this paper. First, the flight environment of the quadrotor aircraft is built based on Dryden atmospheric turbulence model, and the function relation between the rotor lift and the wind field information (wind speed and the wind direction) is obtained by analyzing the change of rotor lift; then the dynamic equation of quadrotor with wind disturbance is deduced using the Newton Euler method, and further a PID controller is designed in order to keep the aircraft flying in fixed-point hover; finally, the wind field information is calculated. The simulation results in Matlab/M-file show that the proposed method can effectively extract the main wind information of flight environment.

Keywords:controllers, MATLAB, quadrotor aircraft, rotor lift, Dryden model, fixed-point hover, wind field estimation

收稿日期：2016-03-10

基金項(xiàng)目：國(guó)家自然科學(xué)基金(61473229、60974146)資助

作者簡(jiǎn)介：：屈耀紅(1971—)，西北工業(yè)大學(xué)副教授，主要從事無(wú)人機(jī)組合導(dǎo)航與風(fēng)場(chǎng)估計(jì)研究。

中圖分類號(hào)：V19

文獻(xiàn)標(biāo)志碼：A

文章編號(hào):1000-2758(2016)04-0684-07