You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
194 lines
6.5 KiB
Python
194 lines
6.5 KiB
Python
import wafo.spectrum.models as sm
|
|
import wafo.transform.models as wtm
|
|
import wafo.objects as wo
|
|
from wafo.spectrum import SpecData1D
|
|
import numpy as np
|
|
from numpy import NAN
|
|
from numpy.testing import assert_array_almost_equal, assert_array_equal
|
|
import unittest
|
|
|
|
|
|
def slow(f):
|
|
f.slow = True
|
|
return f
|
|
|
|
|
|
class TestSpectrumHs7(unittest.TestCase):
|
|
def setUp(self):
|
|
self.Sj = sm.Jonswap(Hm0=7.0, Tp=11)
|
|
self.S = self.Sj.tospecdata()
|
|
|
|
def test_tocovmatrix(self):
|
|
acfmat = self.S.tocov_matrix(nr=3, nt=256, dt=0.1)
|
|
vals = acfmat[:2, :]
|
|
true_vals = np.array([[3.06073383, 0.0000000, -1.67748256, 0.],
|
|
[3.05235423, -0.1674357, -1.66811444,
|
|
0.18693242]])
|
|
assert_array_almost_equal(vals, true_vals)
|
|
|
|
def test_tocovdata(self):
|
|
|
|
Nt = len(self.S.data) - 1
|
|
acf = self.S.tocovdata(nr=0, nt=Nt)
|
|
vals = acf.data[:5]
|
|
|
|
true_vals = np.array(
|
|
[3.06090339, 2.22658399, 0.45307391, -1.17495501, -2.05649042])
|
|
assert_array_almost_equal(vals, true_vals)
|
|
assert((np.abs(vals - true_vals) < 1e-6).all())
|
|
|
|
def test_to_t_pdf(self):
|
|
f = self.S.to_t_pdf(pdef='Tc', paramt=(0, 10, 51), speed=7, seed=100)
|
|
|
|
truth = [0.0, 0.014068786046738972, 0.027384724577108947, 0.039538002584522454,
|
|
0.050183061144017056, 0.05948762020247726, 0.0669017098497974,
|
|
0.07251637759775977, 0.07729759248201125, 0.08151306823047058]
|
|
assert_array_almost_equal(f.data[:10], truth, decimal=1e-3)
|
|
|
|
# estimated error bounds
|
|
truth = [0., 0.00028327, 0.00027281, 0.00042283, 0.00058736, 0.00083936,
|
|
0.00160774, 0.00186591, 0.00196073, 0.00213102]
|
|
self.assertTrue(np.all(0.5 * f.err[:10] <= truth))
|
|
|
|
@slow
|
|
def test_sim(self):
|
|
S = self.S
|
|
|
|
import scipy.stats as st
|
|
x2 = S.sim(20000, 20)
|
|
truth1 = [0, np.sqrt(S.moment(1)[0]), 0., 0.]
|
|
funs = [np.mean, np.std, st.skew, st.kurtosis]
|
|
for fun, trueval in zip(funs, truth1):
|
|
res = fun(x2[:, 1::], axis=0)
|
|
m = res.mean()
|
|
sa = res.std()
|
|
assert(np.abs(m - trueval) < sa)
|
|
|
|
@slow
|
|
def test_sim_nl(self):
|
|
S = self.S
|
|
|
|
import scipy.stats as st
|
|
x2, _x1 = S.sim_nl(ns=20000, cases=40)
|
|
truth1 = [0, np.sqrt(S.moment(1)[0][0])] + S.stats_nl(moments='sk')
|
|
truth1[-1] = truth1[-1] - 3
|
|
|
|
# truth1
|
|
# [0, 1.7495200310090633, 0.18673120577479801, 0.061988521262417606]
|
|
|
|
funs = [np.mean, np.std, st.skew, st.kurtosis]
|
|
for fun, trueval in zip(funs, truth1):
|
|
res = fun(x2.data, axis=0)
|
|
m = res.mean()
|
|
sa = res.std()
|
|
# trueval, m, sa
|
|
assert(np.abs(m - trueval) < 2 * sa)
|
|
|
|
def test_stats_nl(self):
|
|
S = self.S
|
|
me, va, sk, ku = S.stats_nl(moments='mvsk')
|
|
assert(me == 0.0)
|
|
assert_array_almost_equal(va, 3.0608203389019537)
|
|
assert_array_almost_equal(sk, 0.18673120577479801)
|
|
assert_array_almost_equal(ku, 3.0619885212624176)
|
|
|
|
def test_testgaussian(self):
|
|
Hs = self.Sj.Hm0
|
|
S0 = self.S
|
|
# ns =100; dt = .2
|
|
# x1 = S0.sim(ns, dt=dt)
|
|
S = S0.copy()
|
|
me, _va, sk, ku = S.stats_nl(moments='mvsk')
|
|
S.tr = wtm.TrHermite(mean=me, sigma=Hs / 4, skew=sk, kurt=ku, ysigma=Hs / 4)
|
|
ys = wo.mat2timeseries(S.sim(ns=2 ** 13))
|
|
g0, _gemp = ys.trdata()
|
|
t0 = g0.dist2gauss()
|
|
t1 = S0.testgaussian(ns=2 ** 13, test0=None, cases=50)
|
|
assert(sum(t1 > t0) <= 5)
|
|
|
|
|
|
class TestSpectrumHs5(unittest.TestCase):
|
|
def setUp(self):
|
|
self.Sj = sm.Jonswap(Hm0=5.0)
|
|
self.S = self.Sj.tospecdata()
|
|
|
|
def test_moment(self):
|
|
S = self.S
|
|
vals, txt = S.moment()
|
|
true_vals = [1.5614600345079888, 0.95567089481941048]
|
|
true_txt = ['m0', 'm0tt']
|
|
|
|
assert_array_almost_equal(vals, true_vals)
|
|
for tv, v in zip(true_txt, txt):
|
|
assert(tv == v)
|
|
|
|
def test_nyquist_freq(self):
|
|
S = self.S
|
|
assert_array_almost_equal(S.nyquist_freq(), 3.0)
|
|
|
|
def test_sampling_period(self):
|
|
S = self.S
|
|
assert_array_almost_equal(S.sampling_period(), 1.0471975511965976)
|
|
|
|
def test_normalize(self):
|
|
S = self.S
|
|
mom, txt = S.moment(2)
|
|
assert_array_almost_equal(mom,
|
|
[1.5614600345079888, 0.95567089481941048])
|
|
assert_array_equal(txt, ['m0', 'm0tt'])
|
|
vals, _txt = S.moment(2)
|
|
true_vals = [1.5614600345079888, 0.95567089481941048]
|
|
assert_array_almost_equal(vals, true_vals)
|
|
|
|
Sn = S.copy()
|
|
Sn.normalize()
|
|
|
|
# Now the moments should be one
|
|
new_vals, _txt = Sn.moment(2)
|
|
assert_array_almost_equal(new_vals, np.ones(2))
|
|
|
|
def test_characteristic(self):
|
|
S = self.S
|
|
ch, R, txt = S.characteristic(1)
|
|
assert_array_almost_equal(ch, 8.59007646)
|
|
assert_array_almost_equal(R, 0.03040216)
|
|
self.assert_(txt == ['Tm01'])
|
|
|
|
ch, R, txt = S.characteristic([1, 2, 3]) # fact a vector of integers
|
|
assert_array_almost_equal(ch, [8.59007646, 8.03139757, 5.62484314])
|
|
assert_array_almost_equal(R,
|
|
[[0.03040216, 0.02834263, NAN],
|
|
[0.02834263, 0.0274645, NAN],
|
|
[NAN, NAN, 0.01500249]])
|
|
assert_array_equal(txt, ['Tm01', 'Tm02', 'Tm24'])
|
|
|
|
ch, R, txt = S.characteristic('Ss') # fact a string
|
|
assert_array_almost_equal(ch, [0.04963112])
|
|
assert_array_almost_equal(R, [[2.63624782e-06]])
|
|
assert_array_equal(txt, ['Ss'])
|
|
|
|
# fact a list of strings
|
|
ch, R, txt = S.characteristic(['Hm0', 'Tm02'])
|
|
|
|
assert_array_almost_equal(ch,
|
|
[4.99833578, 8.03139757])
|
|
assert_array_almost_equal(R, [[0.05292989, 0.02511371],
|
|
[0.02511371, 0.0274645]])
|
|
assert_array_equal(txt, ['Hm0', 'Tm02'])
|
|
|
|
|
|
class TestSpectrumHs3(unittest.TestCase):
|
|
def test_bandwidth(self):
|
|
|
|
Sj = sm.Jonswap(Hm0=3, Tp=7)
|
|
w = np.linspace(0, 4, 256)
|
|
S = SpecData1D(Sj(w), w) # Make spectrum object from numerical values
|
|
vals = S.bandwidth([0, 1, 2, 3])
|
|
true_vals = [0.73062845, 0.34476034, 0.68277527, 2.90817052]
|
|
assert_array_almost_equal(vals, true_vals)
|
|
|
|
|
|
if __name__ == '__main__':
|
|
import nose
|
|
nose.run()
|