二维改进分形海面模型及海谱分析

第34卷第2期 2013年2月

通 信 学 报 Vol.34 No. 2

Journal on Communications

February 2013

doi:10.3969/j.issn.1000-436x.2013.02.021

二维改进分形海面模型及海谱分析

陈瑜，胡云安，林涛

（海军航空工程学院 控制工程系，山东 烟台 264001）

摘 要：针对分形海面模型的功率谱在空间波数小于基波波数时不能满足正幂率的问题，提出了一种统计模型和归一化带限Weierstrass分形模型相结合的二维海面模型，确定了方向海谱的闭式解，得到了表面功率谱和方向分布函数，并且和有关文献的数据进行了对比，对分形维数、频率幅度尺度因子、风速和风向等参量对海面轮廓的影响做了进一步分析。从仿真结果可知：表面功率谱和PM谱、JONSWAP谱、文圣常谱和南海的实测数据都拟合的很好，方向分布规律和Mitsuyasu方向函数相似，海面形状具有大尺度涌浪和小尺度张力波的特性，并且统计模型和分形模型的性质在改进模型中依然适用，验证了模型的有效性。 关键词：分形模型；海谱；Longuet-Higgin模型；PM谱

中图分类号：O451 文献标识码：A 文章编号：1000-436X(2013)02-0177-09

Two-dimensional improved fractal model of the sea

surface and sea spectrum evaluation

CHEN Yu, HU Yun-an, LIN Tao

(Department of Control Engineering, Naval Aeronautical and Astronautical University, Yantai 264001,China )

Abstract: Aiming at the problem that the ominidirectional sea-wave spectrum of fractal sea-wave model could not satisfy the positive power law when spatial wave numbers was smaller than the fundamental wave number, 2-D sea-wave model was proposed by combining statistical sea model and normalization band-limited Weierstrass- Mandelbrot (WM) fractal sea model. With the proposed method, the closed form expression of the directional sea-wave spectrum was calculated and the ominidirectional sea-wave spectrum and the angular distribution function derived，the comparison with the results in the literature was made, and influences of fractal dimension, frequency and amplitude factor, wind speed and wind di-rection on sea-surface were further analyzed. Simulation results show that omnidirectional sea-wave spectrum is in good agreement with the common spectrum of Pierson-Moskowitz (PM) and JONSWAP and Wen Sheng Chang and the real data of the South China Sea, angular distribution function is very similar to the Mitsuyasu’s distribution function, sea sur-face have both characters of large-scale turbulent waves and micro-scale capillary wave, statistical model and fractal model’s characters are still suitable for improved fractal model, and the correctness of the improved model was validated. Key words: fractal model; sea-wave spectrum; Longuet-Higgin model; PM spectrum

1 引言

海面电磁散射的研究在海洋遥感、海上目标跟踪与识别等领域具有极其重要的意义[1~3]。在模拟海面电磁散射过程中提出一种合适的海面模型非常关键，一般通过基于海谱的统计模型来模拟粗糙海面，但基于海谱的统计方法不能描述海浪产生和繁衍的物理机理。最近由于分形几何在各个方面的

收稿日期：2011-11-22；修回日期：2012-03-15

逐渐应用[4~6]，人们发现二维分形海面既具有统计性又具有随机性，能够更为贴切地描述真实海面，因此利用分形模型对海浪进行研究逐渐成为国内外研究的热点。Jaggard和Sun[7]首次提出了利用一维Weierstrass带限分形函数描述海面的思想，文献[8]在Jaggard的基础上扩展到二维平面，并用此模型去描述了单一小平面的粗糙程度，该小平面在高分辨雷达中得到应用。Berizzi[9,10]在以上模型的

178 通 信 学 报 第34卷

基础上把一维扩展到二维,考虑了海浪的机动性，并引进了平台的运动速度，这种特性很适合被应用到空中和卫星雷达。实际海浪的功率谱包含正、负幂率两部分，而经典分形模型只包含负幂率部分，为了更好的描述实际海面，文献[11]应用PM谱和Weierstrass带限分形函数结合模拟一维海表面，文献[12,13]建立了基于未充分发展海谱的分形海面模型，但以上两者都没给出改进模型的海谱形式，模型的正确性无法验证。王运华等人[14,15]对Weierstrass带限分形函数一维和二维形式进行了改进，给出了方向海谱的形式，并与经典文献进行了比较，两者吻合较好，但该模型无法给出改进部分具体的物理含义。Berizzi[16]在Weierstrass带限分形函数的自相关函数中加入高斯分布的相关函数以实现对分形海谱的改进，得到的海谱与经典文献吻合较好，但从海谱得不到改进分形模型的表达式。文中根据以上问题对分形模型进行了改进。

海面上任一点的波都是由本地产生的毛细波和从其他方向传来涌浪叠加而成，而传统分形模型不考虑涌浪的情况，因而产生海谱中不包括正幂率部分，根据这一问题，基于Longuet-Higgin提出了一种改进的归一化带限Weierstrass分形海面模型，能很好的体现海浪的涌浪特性，为了进一步验证模型的正确性，推导了改进模型的方向海谱，并和有关文献[17~21]数据进行了比较，对不同参量下改进模型、经典分形模型和Longuet-Higgins模型的二维海面轮廓也进行了详细讨论。

kn= k0bn。角频率

ωn与波数kn满足色散关系为

n τs为海面表面张力，ρ为海水

密度。vx和vy为遥感平台延在x轴和y轴的速度分量，当传感器观察的时间很短时，βn(t)和Фn可看为常量，两者相互独立，Фn是[ π, π]上均匀分布的随机相位。令海浪的运动方向角βn的均值为 0，

n n 0，Ψn的概率密度函数的形式如式(2)所

示[26]

pn( ) gncos

2

2en

rect

2π

1/2

gn

2

2en 1

1 en n 1

,en T

π 1 2en n 1

2

(2)

其中，rect( )为单位矩形函数，当x 0.5时，rect(x)值为1，否则为0，Г( )为Gamma函数，nT

是波长

n接近过渡波长 T 2 时海浪的阶，

T

nT int ln kT/k0 /lnb ,kT 2 / T。

经典分形粗糙海面的功率谱为负幂率谱，正

好对应于实际海谱中的一段，并不能反映真实海面的功率谱。在以往的海浪模拟中都涉及到海浪谱的选择，而这些谱大多是根据大量实验数据，在半经验、半理论的基础上运用数学分析方法得到的。目前在国际上应用比较广泛的有Pierson- Moscowitz谱(PM谱)[17]和JONSWAP谱[18]，而在我国海域一般应用“海港水文规范”所依据的文圣常谱[19]。

PM谱、JONSWAP谱和文圣常谱中都有一个峰值，其对应的波数为km，风速越大，km越小，这表明海表面受2种波共同影响,当k＜km时，表面长重力波起主要作用，功率谱满足正幂率谱，而当k>km时，短重力波及张力波起主要作用，功率谱满足负幂率谱，这就启发笔者应用下面改进的二维分形模型来模拟实际海面，它表示为

2 改进二维分形海面模型

经典的带限Weierstrass分形模型可表示为[22]

f x,y,t b

n 0Nf 1

D 3 n

sink0bn x vxt cos n t y vyt sin n t nt n

(1)

其中，σ为海面高度起伏均方根。b为尺度因子，

b>1，b越大，海面起伏越大，海面越光滑[23]。Nf>400表示模型中含有正弦分量的个数。η为归一化因子，

x,y,t aijcos it ki[xcos j ysin j] ij

i 1j 1

Nf 1

NM

。粗糙度D决定着海

xvtcos t xn nt n b D 3 nsin k0bn

y vyt sin t n 0 n

(3) 其中，第一部分为Longuet-Higgins模型，该模型是基于海谱建立的，能很好地反映长重力波的情况。在经典分形模型中加入Longuet-Higgins模型解决

面的粗糙程度，D值越大，海面轮廓在细节方面变化得越剧烈，D的取值为2<D<3[24]。k0为基波波数，控制着空间海浪的基频，与波数kn之间的关系为

第2期 陈瑜等：二维改进分形海面模型及海谱分析 179

了分形模型不能反映涌浪的情况，在求解时，Longuet-Higgins模型中方向海谱波数的取值要小于谱峰值km，aij、ωi、 ki、θj和εij分别为组成波的振幅、角频率、波数、方向角和初相位，(x,y)为波面上某点的坐标，εij在[0,2π]范围内随机选取，波数ki与角频率ωi在深水域的关系为ki=ωi2g，aij通常由

aij 确定，其中，S(ωi,θj)为方向

求数之间的关系W(k,Ф)kdkdФ=S(k)G(k,Ф)dkdФ[27]，

得功率谱S(k)和方向分布函数G(k,Φ)。

根据步骤1求相关函数，其定义为

x1,x2,y1,y2,t1,t2

E E E E x1,y1,t1 x2,y2,t2

(5)

海谱，可写为下列形式[17]

其中，Eε( )、EФ( )、Eθ( )和Eβ( ) 分别为变量εij、 Фn、θj 和βn的平均期望，当vx=vy=0并且相位εij和Фn相互独立时，把式(3)代入式(5)可得

S(ω,θ)= S(ω)G(ω,θ) (4)

其中，S(ω)为功率谱，G(ω,θ)为方向分布函数。方向分布函数一般采用实际观测获得的方向分布函数。在本文中海浪的运动方向角βn采用Berizzi提出的概率密度函数形式[25]。在模拟海浪的过程中功率谱的频率分隔采用能量等分法，避免波浪以2π/Δω的周期重复出现，具体方法参见文献[26]。

式(3)具有2种极端情况：1) 张力波向涌浪转换的临界状态。假设海浪在有限区域内，以保证海浪没有衰减，这时基波波长 0会很长，达到几百米，根据波长和波数的关系式k=2 / ，基波波数k0和谱峰值km会相应的趋于0，这时式(3)中的aij 0，改进模型只包含分形部分，而没有Longuet-Higgins模型部分(涌浪部分)。2) 只存在涌浪的情况，即海面上没有短重力波和小尺度毛细波的状态。根据b和D的取值和适用状态，这时b ，D 2，相应式(3)第二部分中的b(D 3)n 0，这时改进模型只包含Longuet-Higgins模型部分，而没有分形部分。从以上分析可以看出，改进模型在极端情况下也同样是正确和有效的。

NMa2 ij

x, y, t E j cos[ki xcos j ysin j i t]

i 1j 12

2 2Nf 12(D 3)n

E n bcos[k(xcosysin)t] nnnn

2n 0 (6)

x1 x2, y y1 y2, t t1 t2，并且假定时其中，x

t 0，求相变海表面在短时间内为瞬时状态，这时

关函数的过程中，在Longuet-Higgins模型未知方向分布函数的情况下，假定aij只与功率谱有关，这时

aij ai ，为了得到闭合形式的相关函数，定义复空间相关函数为

ai

))] r , cos( j E jexp[iki(r i 1j 12

NM2

2

2

2Nf 1n 0

b

2D 3 n

)) cos( n E nexp ikn(r

(7)

式(7)为式(6)复空间相关函数的极坐标形式，在较短的时间内，可以认为 n和时间无关，并且满足

E( n) E( j) 0，其中，0为主风向与x轴的夹角，同时 n n 0， j j 0，并且Ψn和Ψj的概率密度函数如式(2)所示，应用关系式

3 方向海谱的求解

对改进分形模型的评价必须分析其方向海谱的形式，在求解的过程中取vx=vy=0，求解海谱的步骤如下。

步骤1 求解改进模型的自相关函数，采用极坐

, , t)，并且只考虑海面在瞬态的情况，标的形式 (r

) (r ,0)。 , , t 0，相关函数变为 (r这时令

exp(ixsinu)

m

J

m

(x)exp(imu)

(8)

Sj(m)

1

pj( )exp( im )d 2

) ( 1)mJ m(knr )Jm(knr

Sj(m) Sj( m)取式(7)的实部，可得

步骤2 对相关函数进行Fourier变换得到方向

海谱W(k,Ф)，形式如下所示

1 2 cos

)e jkr dr W(k, ) r r (,d 002

步骤3 根据方向海谱、功率谱和方向分布函

a2J kr , r i0i /2

i=1

M

N

2 a 1

2

i

i 1j 1

m 1

N

m

) 2m( Sj 2m cos J2m kir0

180

Nf 1n 0

通 信 学 报 第34卷

/2 2 2b2 D 3 nJ0 knr

k, Wn

Nf 1n 0

1

il( 0) kn k b Sn l exp kn l

(16) Jm 和Jl为不同阶的第一类Bessel函数，当

应用关系式 0时有式(9)变为 r

l 0

2

22 D 3 n

2 b

222 D 3 n

S 2l cos 1 J kr 2l( )

l

2l

n

n

1

k aim( 0) Sj m exp kn kj 1 n m

2

n

M

a2/2 2 0, i

i 1

NNf 1n 0

b

2

22D 3 n

/2 (10)

Wn k,

1 k, (17) Wn k, W n

2

其中，(·)*代表共轭，可得

又由于 S i S d ，

ai 2

i 1

W k,

1 1

2m( 0) ki k a 2Sj 2m cos ki m 12 i 1j 1

2

i

22

Nf 1

2 D 3 n

N

M

并且ω只取小于峰值角频率ωm的值，其中，ωm2/g=km，这样可以求得

b

n 0 (11)

k0bn k

1 1

Sll 22cos2( ) n0 k0bn l 12π

(18)

S d / S d (12) 根据步骤3方向海谱和功率谱、方向分布函数

0

m

之间的关系，把分形海谱写为如下形式

根据步骤2对相关函数进行Fourier变换得到

方向海谱W(k,Ф)。

在求解方向谱的过程中令i=n，式(9)可表示为

W k, W1 k 1 W2 k, /W1 k (19)

这样可以求得

其中，

, r

Nf 1n 0

n

(13) , r

S k kW1 k (20)

G k, 1 W2 k, /W1 k (21)

2π

并且

M

n

2

a2 , r n

j 1

m

Sj m cos Jm kn r m( /2 0)

2π

G k, d 1 W2 k, /W1 k d 2 (22)

2π0

上式满足了方向分布函数具有 G k, d 2 的性质。

N

b

22 D 3 n

ai2

W1 k ki k 设

(14) i 12ki

l

Nf 1n 0

Sn l cos Jl knr l( /2 0)

其中，kn’的取值与ki的取值相同，方向海谱可写成下式的形式，即

b

2 D 3 n

2 2

2k0b

n

k0bn k (23)

则根据式(18)、式(19)和式(23)可得

W k, Wn k, (15)

W2 k, n 0

Nf 1

为了得到方向海谱的闭式解，仍然采用复空间

相关函数代替空间相关函数进行求解，并应用关系式 rJl ar Jl br dr a b /a，求得的复空

0

1

2m( 0) ki k 2 a Sj 2m cos ki m 1i 1j 1

2

i

NM

2

22

Nf 1n 0

b

2 D 3 n

k, 为 间的方向海谱Wn

1

S2lcos2l( ) n0 k0b l 1

k0bn k

(24)

第2期 陈瑜等：二维改进分形海面模型及海谱分析 181

根据式(20)和式(23)可求得功率谱的形式为

ai2

,0≤k k0 2

(25) S k 22

2D 3 2 D 3 1 ,k≥k0k0 k 2lnb根据式(21)、式(23)和式(24)可得方向分布函数为

G k,

NM4 ai2

S2mcos2m()kk ji0

i 1j 1ki

1 2N

aim 1 ki k i 1ki

0≤k<k0 (26)

G k,

Nf 14 b2 D 3 n n

S2lcos2l()kbk n00 n 0k0b 1 Nf 12D 3n

bn l 1 k0b k kbn 00

k>k0 (27)

式(25)求得的功率海谱在波数0≤k<k0范围内

可根据不同的海域的涌浪情况选择不同的功率谱S( )，使改进模型的使用范围更广，并且方向海谱能表示成功率谱和方向分布函数相乘的形式，这可以证明在求解式(6)时，在未知Longuet-Higgins模型方向函数的情况下，假定aij只与功率谱有关而与方向函数无关是成立的，求解过程中海面起伏均方根及其基波波数的取值参考文献[14,28]。

4 数值计算与分析

4.1 本文功率海谱与各海域标准海谱的比较

图1为不同风速下，涌浪部分的海谱为PM谱

时，本文方法的功率谱、方法1[14,15]的功率谱、方法2[16]的功率谱与PM谱的比较。图1(a)中海面19.5m高处的风速为15m/s，尺度因子b=1.010 2，分形维数D=2.047，迭代次数Nf=400，图1(b)中风速为10m/s，其他参数不变。

从图1(a)可以看出，在和其他方法相同的仿真条件下，文中提出的方法与其他方法相比[14~16]和PM谱吻合的更好。从图1(b)可以看出，当风速改变时，其他方法的功率谱和PM谱有一定的偏差，而本文方法和PM谱大部分波数范围内吻合的较好，只在基波附近在幅值上稍有偏差。为了进一步验证本文方法的优越性和有效性，图2为本文功率海谱、文圣常谱和PM谱与不同风区JONSWAP谱的比较。

图2(a)是JONSWAP谱在风速为10m/s，峰升高因子 为3.3，风区为50km时其他谱与其的拟合情况。图2(b)是风区为20km时其他谱与其的拟合情况。其中，改进模型涌浪部分的海谱选为JONSWAP谱。

在各谱与JONSWAP谱拟合最好的情况时，从图2(a)和图2(b)可以看出，PM谱在谱峰附近与其相差很大。文圣常谱在谱峰附近与其相一致。本文海谱在高频范围内与其吻合不是很好，在其余频段吻合的很好。在图2中PM谱、文圣常谱和本文功率海谱与JONSWAP谱之间的均方差如表1所示。

从以上的分析比较可以看出，在不同风区下，PM谱与JONSWAP谱之间的误差最大，

两者不能很好地

(a) 风速为15m/s (b) 风速为10m/s

图1 不同风速下，不同改进模型的功率谱和PM谱的比较