test_mlab.py 87 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798991001011021031041051061071081091101111121131141151161171181191201211221231241251261271281291301311321331341351361371381391401411421431441451461471481491501511521531541551561571581591601611621631641651661671681691701711721731741751761771781791801811821831841851861871881891901911921931941951961971981992002012022032042052062072082092102112122132142152162172182192202212222232242252262272282292302312322332342352362372382392402412422432442452462472482492502512522532542552562572582592602612622632642652662672682692702712722732742752762772782792802812822832842852862872882892902912922932942952962972982993003013023033043053063073083093103113123133143153163173183193203213223233243253263273283293303313323333343353363373383393403413423433443453463473483493503513523533543553563573583593603613623633643653663673683693703713723733743753763773783793803813823833843853863873883893903913923933943953963973983994004014024034044054064074084094104114124134144154164174184194204214224234244254264274284294304314324334344354364374384394404414424434444454464474484494504514524534544554564574584594604614624634644654664674684694704714724734744754764774784794804814824834844854864874884894904914924934944954964974984995005015025035045055065075085095105115125135145155165175185195205215225235245255265275285295305315325335345355365375385395405415425435445455465475485495505515525535545555565575585595605615625635645655665675685695705715725735745755765775785795805815825835845855865875885895905915925935945955965975985996006016026036046056066076086096106116126136146156166176186196206216226236246256266276286296306316326336346356366376386396406416426436446456466476486496506516526536546556566576586596606616626636646656666676686696706716726736746756766776786796806816826836846856866876886896906916926936946956966976986997007017027037047057067077087097107117127137147157167177187197207217227237247257267277287297307317327337347357367377387397407417427437447457467477487497507517527537547557567577587597607617627637647657667677687697707717727737747757767777787797807817827837847857867877887897907917927937947957967977987998008018028038048058068078088098108118128138148158168178188198208218228238248258268278288298308318328338348358368378388398408418428438448458468478488498508518528538548558568578588598608618628638648658668678688698708718728738748758768778788798808818828838848858868878888898908918928938948958968978988999009019029039049059069079089099109119129139149159169179189199209219229239249259269279289299309319329339349359369379389399409419429439449459469479489499509519529539549559569579589599609619629639649659669679689699709719729739749759769779789799809819829839849859869879889899909919929939949959969979989991000100110021003100410051006100710081009101010111012101310141015101610171018101910201021102210231024102510261027102810291030103110321033103410351036103710381039104010411042104310441045104610471048104910501051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120112111221123112411251126112711281129113011311132113311341135113611371138113911401141114211431144114511461147114811491150115111521153115411551156115711581159116011611162116311641165116611671168116911701171117211731174117511761177117811791180118111821183118411851186118711881189119011911192119311941195119611971198119912001201120212031204120512061207120812091210121112121213121412151216121712181219122012211222122312241225122612271228122912301231123212331234123512361237123812391240124112421243124412451246124712481249125012511252125312541255125612571258125912601261126212631264126512661267126812691270127112721273127412751276127712781279128012811282128312841285128612871288128912901291129212931294129512961297129812991300130113021303130413051306130713081309131013111312131313141315131613171318131913201321132213231324132513261327132813291330133113321333133413351336133713381339134013411342134313441345134613471348134913501351135213531354135513561357135813591360136113621363136413651366136713681369137013711372137313741375137613771378137913801381138213831384138513861387138813891390139113921393139413951396139713981399140014011402140314041405140614071408140914101411141214131414141514161417141814191420142114221423142414251426142714281429143014311432143314341435143614371438143914401441144214431444144514461447144814491450145114521453145414551456145714581459146014611462146314641465146614671468146914701471147214731474147514761477147814791480148114821483148414851486148714881489149014911492149314941495149614971498149915001501150215031504150515061507150815091510151115121513151415151516151715181519152015211522152315241525152615271528152915301531153215331534153515361537153815391540154115421543154415451546154715481549155015511552155315541555155615571558155915601561156215631564156515661567156815691570157115721573157415751576157715781579158015811582158315841585158615871588158915901591159215931594159515961597159815991600160116021603160416051606160716081609161016111612161316141615161616171618161916201621162216231624162516261627162816291630163116321633163416351636163716381639164016411642164316441645164616471648164916501651165216531654165516561657165816591660166116621663166416651666166716681669167016711672167316741675167616771678167916801681168216831684168516861687168816891690169116921693169416951696169716981699170017011702170317041705170617071708170917101711171217131714171517161717171817191720172117221723172417251726172717281729173017311732173317341735173617371738173917401741174217431744174517461747174817491750175117521753175417551756175717581759176017611762176317641765176617671768176917701771177217731774177517761777177817791780178117821783178417851786178717881789179017911792179317941795179617971798179918001801180218031804180518061807180818091810181118121813181418151816181718181819182018211822182318241825182618271828182918301831183218331834183518361837183818391840184118421843184418451846184718481849185018511852185318541855185618571858185918601861186218631864186518661867186818691870187118721873187418751876187718781879188018811882188318841885188618871888188918901891189218931894189518961897189818991900190119021903190419051906190719081909191019111912191319141915191619171918191919201921192219231924192519261927192819291930193119321933193419351936193719381939194019411942194319441945194619471948194919501951195219531954195519561957195819591960196119621963196419651966196719681969197019711972197319741975197619771978197919801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220132014201520162017201820192020202120222023202420252026202720282029203020312032203320342035203620372038203920402041204220432044204520462047204820492050205120522053205420552056205720582059206020612062206320642065206620672068206920702071207220732074207520762077207820792080208120822083208420852086208720882089209020912092209320942095209620972098209921002101210221032104210521062107210821092110211121122113211421152116211721182119212021212122212321242125212621272128212921302131213221332134213521362137213821392140214121422143214421452146214721482149215021512152215321542155215621572158215921602161216221632164216521662167216821692170217121722173217421752176217721782179218021812182218321842185218621872188218921902191219221932194219521962197219821992200220122022203220422052206220722082209221022112212221322142215221622172218221922202221222222232224222522262227222822292230223122322233223422352236223722382239224022412242224322442245
  1. import tempfile
  2. from numpy.testing import (assert_allclose, assert_almost_equal,
  3. assert_array_equal, assert_array_almost_equal_nulp)
  4. import numpy as np
  5. import datetime as datetime
  6. import pytest
  7. import matplotlib.mlab as mlab
  8. from matplotlib.cbook.deprecation import MatplotlibDeprecationWarning
  9. def _stride_repeat(*args, **kwargs):
  10. with pytest.warns(MatplotlibDeprecationWarning):
  11. return mlab.stride_repeat(*args, **kwargs)
  12. class TestStride:
  13. def get_base(self, x):
  14. y = x
  15. while y.base is not None:
  16. y = y.base
  17. return y
  18. def calc_window_target(self, x, NFFT, noverlap=0, axis=0):
  19. '''This is an adaptation of the original window extraction
  20. algorithm. This is here to test to make sure the new implementation
  21. has the same result'''
  22. step = NFFT - noverlap
  23. ind = np.arange(0, len(x) - NFFT + 1, step)
  24. n = len(ind)
  25. result = np.zeros((NFFT, n))
  26. # do the ffts of the slices
  27. for i in range(n):
  28. result[:, i] = x[ind[i]:ind[i]+NFFT]
  29. if axis == 1:
  30. result = result.T
  31. return result
  32. @pytest.mark.parametrize('shape', [(), (10, 1)], ids=['0D', '2D'])
  33. def test_stride_windows_invalid_input_shape(self, shape):
  34. x = np.arange(np.prod(shape)).reshape(shape)
  35. with pytest.raises(ValueError):
  36. mlab.stride_windows(x, 5)
  37. @pytest.mark.parametrize('n, noverlap',
  38. [(0, None), (11, None), (2, 2), (2, 3)],
  39. ids=['n less than 1', 'n greater than input',
  40. 'noverlap greater than n',
  41. 'noverlap equal to n'])
  42. def test_stride_windows_invalid_params(self, n, noverlap):
  43. x = np.arange(10)
  44. with pytest.raises(ValueError):
  45. mlab.stride_windows(x, n, noverlap)
  46. @pytest.mark.parametrize('shape', [(), (10, 1)], ids=['0D', '2D'])
  47. def test_stride_repeat_invalid_input_shape(self, shape):
  48. x = np.arange(np.prod(shape)).reshape(shape)
  49. with pytest.raises(ValueError):
  50. _stride_repeat(x, 5)
  51. @pytest.mark.parametrize('axis', [-1, 2],
  52. ids=['axis less than 0',
  53. 'axis greater than input shape'])
  54. def test_stride_repeat_invalid_axis(self, axis):
  55. x = np.array(0)
  56. with pytest.raises(ValueError):
  57. _stride_repeat(x, 5, axis=axis)
  58. def test_stride_repeat_n_lt_1_ValueError(self):
  59. x = np.arange(10)
  60. with pytest.raises(ValueError):
  61. _stride_repeat(x, 0)
  62. @pytest.mark.parametrize('axis', [0, 1], ids=['axis0', 'axis1'])
  63. @pytest.mark.parametrize('n', [1, 5], ids=['n1', 'n5'])
  64. def test_stride_repeat(self, n, axis):
  65. x = np.arange(10)
  66. y = _stride_repeat(x, n, axis=axis)
  67. expected_shape = [10, 10]
  68. expected_shape[axis] = n
  69. yr = np.repeat(np.expand_dims(x, axis), n, axis=axis)
  70. assert yr.shape == y.shape
  71. assert_array_equal(yr, y)
  72. assert tuple(expected_shape) == y.shape
  73. assert self.get_base(y) is x
  74. @pytest.mark.parametrize('axis', [0, 1], ids=['axis0', 'axis1'])
  75. @pytest.mark.parametrize('n, noverlap',
  76. [(1, 0), (5, 0), (15, 2), (13, -3)],
  77. ids=['n1-noverlap0', 'n5-noverlap0',
  78. 'n15-noverlap2', 'n13-noverlapn3'])
  79. def test_stride_windows(self, n, noverlap, axis):
  80. x = np.arange(100)
  81. y = mlab.stride_windows(x, n, noverlap=noverlap, axis=axis)
  82. expected_shape = [0, 0]
  83. expected_shape[axis] = n
  84. expected_shape[1 - axis] = 100 // (n - noverlap)
  85. yt = self.calc_window_target(x, n, noverlap=noverlap, axis=axis)
  86. assert yt.shape == y.shape
  87. assert_array_equal(yt, y)
  88. assert tuple(expected_shape) == y.shape
  89. assert self.get_base(y) is x
  90. @pytest.mark.parametrize('axis', [0, 1], ids=['axis0', 'axis1'])
  91. def test_stride_windows_n32_noverlap0_unflatten(self, axis):
  92. n = 32
  93. x = np.arange(n)[np.newaxis]
  94. x1 = np.tile(x, (21, 1))
  95. x2 = x1.flatten()
  96. y = mlab.stride_windows(x2, n, axis=axis)
  97. if axis == 0:
  98. x1 = x1.T
  99. assert y.shape == x1.shape
  100. assert_array_equal(y, x1)
  101. def test_stride_ensure_integer_type(self):
  102. N = 100
  103. x = np.full(N + 20, np.nan)
  104. y = x[10:-10]
  105. y[:] = 0.3
  106. # previous to #3845 lead to corrupt access
  107. y_strided = mlab.stride_windows(y, n=33, noverlap=0.6)
  108. assert_array_equal(y_strided, 0.3)
  109. # previous to #3845 lead to corrupt access
  110. y_strided = mlab.stride_windows(y, n=33.3, noverlap=0)
  111. assert_array_equal(y_strided, 0.3)
  112. # even previous to #3845 could not find any problematic
  113. # configuration however, let's be sure it's not accidentally
  114. # introduced
  115. y_strided = _stride_repeat(y, n=33.815)
  116. assert_array_equal(y_strided, 0.3)
  117. @pytest.fixture
  118. def tempcsv():
  119. with tempfile.TemporaryFile(suffix='csv', mode="w+", newline='') as fd:
  120. yield fd
  121. def test_csv2rec_names_with_comments(tempcsv):
  122. tempcsv.write('# comment\n1,2,3\n4,5,6\n')
  123. tempcsv.seek(0)
  124. array = mlab._csv2rec(tempcsv, names='a,b,c')
  125. assert len(array) == 2
  126. assert len(array.dtype) == 3
  127. @pytest.mark.parametrize('input, kwargs', [
  128. ('01/11/14\n'
  129. '03/05/76 12:00:01 AM\n'
  130. '07/09/83 5:17:34 PM\n'
  131. '06/20/2054 2:31:45 PM\n'
  132. '10/31/00 11:50:23 AM\n',
  133. {}),
  134. ('11/01/14\n'
  135. '05/03/76 12:00:01 AM\n'
  136. '09/07/83 5:17:34 PM\n'
  137. '20/06/2054 2:31:45 PM\n'
  138. '31/10/00 11:50:23 AM\n',
  139. {'dayfirst': True}),
  140. ('14/01/11\n'
  141. '76/03/05 12:00:01 AM\n'
  142. '83/07/09 5:17:34 PM\n'
  143. '2054/06/20 2:31:45 PM\n'
  144. '00/10/31 11:50:23 AM\n',
  145. {'yearfirst': True}),
  146. ], ids=['usdate', 'dayfirst', 'yearfirst'])
  147. def test_csv2rec_dates(tempcsv, input, kwargs):
  148. tempcsv.write(input)
  149. expected = [datetime.datetime(2014, 1, 11, 0, 0),
  150. datetime.datetime(1976, 3, 5, 0, 0, 1),
  151. datetime.datetime(1983, 7, 9, 17, 17, 34),
  152. datetime.datetime(2054, 6, 20, 14, 31, 45),
  153. datetime.datetime(2000, 10, 31, 11, 50, 23)]
  154. tempcsv.seek(0)
  155. array = mlab._csv2rec(tempcsv, names='a', **kwargs)
  156. assert_array_equal(array['a'].tolist(), expected)
  157. def _apply_window(*args, **kwargs):
  158. with pytest.warns(MatplotlibDeprecationWarning):
  159. return mlab.apply_window(*args, **kwargs)
  160. class TestWindow:
  161. def setup(self):
  162. np.random.seed(0)
  163. n = 1000
  164. self.sig_rand = np.random.standard_normal(n) + 100.
  165. self.sig_ones = np.ones(n)
  166. def check_window_apply_repeat(self, x, window, NFFT, noverlap):
  167. '''This is an adaptation of the original window application
  168. algorithm. This is here to test to make sure the new implementation
  169. has the same result'''
  170. step = NFFT - noverlap
  171. ind = np.arange(0, len(x) - NFFT + 1, step)
  172. n = len(ind)
  173. result = np.zeros((NFFT, n))
  174. if np.iterable(window):
  175. windowVals = window
  176. else:
  177. windowVals = window(np.ones(NFFT, x.dtype))
  178. # do the ffts of the slices
  179. for i in range(n):
  180. result[:, i] = windowVals * x[ind[i]:ind[i]+NFFT]
  181. return result
  182. def test_window_none_rand(self):
  183. res = mlab.window_none(self.sig_ones)
  184. assert_array_equal(res, self.sig_ones)
  185. def test_window_none_ones(self):
  186. res = mlab.window_none(self.sig_rand)
  187. assert_array_equal(res, self.sig_rand)
  188. def test_window_hanning_rand(self):
  189. targ = np.hanning(len(self.sig_rand)) * self.sig_rand
  190. res = mlab.window_hanning(self.sig_rand)
  191. assert_allclose(targ, res, atol=1e-06)
  192. def test_window_hanning_ones(self):
  193. targ = np.hanning(len(self.sig_ones))
  194. res = mlab.window_hanning(self.sig_ones)
  195. assert_allclose(targ, res, atol=1e-06)
  196. def test_apply_window_1D_axis1_ValueError(self):
  197. x = self.sig_rand
  198. window = mlab.window_hanning
  199. with pytest.raises(ValueError):
  200. _apply_window(x, window, axis=1, return_window=False)
  201. def test_apply_window_1D_els_wrongsize_ValueError(self):
  202. x = self.sig_rand
  203. window = mlab.window_hanning(np.ones(x.shape[0]-1))
  204. with pytest.raises(ValueError):
  205. _apply_window(x, window)
  206. def test_apply_window_0D_ValueError(self):
  207. x = np.array(0)
  208. window = mlab.window_hanning
  209. with pytest.raises(ValueError):
  210. _apply_window(x, window, axis=1, return_window=False)
  211. def test_apply_window_3D_ValueError(self):
  212. x = self.sig_rand[np.newaxis][np.newaxis]
  213. window = mlab.window_hanning
  214. with pytest.raises(ValueError):
  215. _apply_window(x, window, axis=1, return_window=False)
  216. def test_apply_window_hanning_1D(self):
  217. x = self.sig_rand
  218. window = mlab.window_hanning
  219. window1 = mlab.window_hanning(np.ones(x.shape[0]))
  220. y, window2 = _apply_window(x, window, return_window=True)
  221. yt = window(x)
  222. assert yt.shape == y.shape
  223. assert x.shape == y.shape
  224. assert_allclose(yt, y, atol=1e-06)
  225. assert_array_equal(window1, window2)
  226. def test_apply_window_hanning_1D_axis0(self):
  227. x = self.sig_rand
  228. window = mlab.window_hanning
  229. y = _apply_window(x, window, axis=0, return_window=False)
  230. yt = window(x)
  231. assert yt.shape == y.shape
  232. assert x.shape == y.shape
  233. assert_allclose(yt, y, atol=1e-06)
  234. def test_apply_window_hanning_els_1D_axis0(self):
  235. x = self.sig_rand
  236. window = mlab.window_hanning(np.ones(x.shape[0]))
  237. window1 = mlab.window_hanning
  238. y = _apply_window(x, window, axis=0, return_window=False)
  239. yt = window1(x)
  240. assert yt.shape == y.shape
  241. assert x.shape == y.shape
  242. assert_allclose(yt, y, atol=1e-06)
  243. def test_apply_window_hanning_2D_axis0(self):
  244. x = np.random.standard_normal([1000, 10]) + 100.
  245. window = mlab.window_hanning
  246. y = _apply_window(x, window, axis=0, return_window=False)
  247. yt = np.zeros_like(x)
  248. for i in range(x.shape[1]):
  249. yt[:, i] = window(x[:, i])
  250. assert yt.shape == y.shape
  251. assert x.shape == y.shape
  252. assert_allclose(yt, y, atol=1e-06)
  253. def test_apply_window_hanning_els1_2D_axis0(self):
  254. x = np.random.standard_normal([1000, 10]) + 100.
  255. window = mlab.window_hanning(np.ones(x.shape[0]))
  256. window1 = mlab.window_hanning
  257. y = _apply_window(x, window, axis=0, return_window=False)
  258. yt = np.zeros_like(x)
  259. for i in range(x.shape[1]):
  260. yt[:, i] = window1(x[:, i])
  261. assert yt.shape == y.shape
  262. assert x.shape == y.shape
  263. assert_allclose(yt, y, atol=1e-06)
  264. def test_apply_window_hanning_els2_2D_axis0(self):
  265. x = np.random.standard_normal([1000, 10]) + 100.
  266. window = mlab.window_hanning
  267. window1 = mlab.window_hanning(np.ones(x.shape[0]))
  268. y, window2 = _apply_window(x, window, axis=0, return_window=True)
  269. yt = np.zeros_like(x)
  270. for i in range(x.shape[1]):
  271. yt[:, i] = window1*x[:, i]
  272. assert yt.shape == y.shape
  273. assert x.shape == y.shape
  274. assert_allclose(yt, y, atol=1e-06)
  275. assert_array_equal(window1, window2)
  276. def test_apply_window_hanning_els3_2D_axis0(self):
  277. x = np.random.standard_normal([1000, 10]) + 100.
  278. window = mlab.window_hanning
  279. window1 = mlab.window_hanning(np.ones(x.shape[0]))
  280. y, window2 = _apply_window(x, window, axis=0, return_window=True)
  281. yt = _apply_window(x, window1, axis=0, return_window=False)
  282. assert yt.shape == y.shape
  283. assert x.shape == y.shape
  284. assert_allclose(yt, y, atol=1e-06)
  285. assert_array_equal(window1, window2)
  286. def test_apply_window_hanning_2D_axis1(self):
  287. x = np.random.standard_normal([10, 1000]) + 100.
  288. window = mlab.window_hanning
  289. y = _apply_window(x, window, axis=1, return_window=False)
  290. yt = np.zeros_like(x)
  291. for i in range(x.shape[0]):
  292. yt[i, :] = window(x[i, :])
  293. assert yt.shape == y.shape
  294. assert x.shape == y.shape
  295. assert_allclose(yt, y, atol=1e-06)
  296. def test_apply_window_hanning_2D__els1_axis1(self):
  297. x = np.random.standard_normal([10, 1000]) + 100.
  298. window = mlab.window_hanning(np.ones(x.shape[1]))
  299. window1 = mlab.window_hanning
  300. y = _apply_window(x, window, axis=1, return_window=False)
  301. yt = np.zeros_like(x)
  302. for i in range(x.shape[0]):
  303. yt[i, :] = window1(x[i, :])
  304. assert yt.shape == y.shape
  305. assert x.shape == y.shape
  306. assert_allclose(yt, y, atol=1e-06)
  307. def test_apply_window_hanning_2D_els2_axis1(self):
  308. x = np.random.standard_normal([10, 1000]) + 100.
  309. window = mlab.window_hanning
  310. window1 = mlab.window_hanning(np.ones(x.shape[1]))
  311. y, window2 = _apply_window(x, window, axis=1, return_window=True)
  312. yt = np.zeros_like(x)
  313. for i in range(x.shape[0]):
  314. yt[i, :] = window1 * x[i, :]
  315. assert yt.shape == y.shape
  316. assert x.shape == y.shape
  317. assert_allclose(yt, y, atol=1e-06)
  318. assert_array_equal(window1, window2)
  319. def test_apply_window_hanning_2D_els3_axis1(self):
  320. x = np.random.standard_normal([10, 1000]) + 100.
  321. window = mlab.window_hanning
  322. window1 = mlab.window_hanning(np.ones(x.shape[1]))
  323. y = _apply_window(x, window, axis=1, return_window=False)
  324. yt = _apply_window(x, window1, axis=1, return_window=False)
  325. assert yt.shape == y.shape
  326. assert x.shape == y.shape
  327. assert_allclose(yt, y, atol=1e-06)
  328. def test_apply_window_stride_windows_hanning_2D_n13_noverlapn3_axis0(self):
  329. x = self.sig_rand
  330. window = mlab.window_hanning
  331. yi = mlab.stride_windows(x, n=13, noverlap=2, axis=0)
  332. y = _apply_window(yi, window, axis=0, return_window=False)
  333. yt = self.check_window_apply_repeat(x, window, 13, 2)
  334. assert yt.shape == y.shape
  335. assert x.shape != y.shape
  336. assert_allclose(yt, y, atol=1e-06)
  337. def test_apply_window_hanning_2D_stack_axis1(self):
  338. ydata = np.arange(32)
  339. ydata1 = ydata+5
  340. ydata2 = ydata+3.3
  341. ycontrol1 = _apply_window(ydata1, mlab.window_hanning)
  342. ycontrol2 = mlab.window_hanning(ydata2)
  343. ydata = np.vstack([ydata1, ydata2])
  344. ycontrol = np.vstack([ycontrol1, ycontrol2])
  345. ydata = np.tile(ydata, (20, 1))
  346. ycontrol = np.tile(ycontrol, (20, 1))
  347. result = _apply_window(ydata, mlab.window_hanning, axis=1,
  348. return_window=False)
  349. assert_allclose(ycontrol, result, atol=1e-08)
  350. def test_apply_window_hanning_2D_stack_windows_axis1(self):
  351. ydata = np.arange(32)
  352. ydata1 = ydata+5
  353. ydata2 = ydata+3.3
  354. ycontrol1 = _apply_window(ydata1, mlab.window_hanning)
  355. ycontrol2 = mlab.window_hanning(ydata2)
  356. ydata = np.vstack([ydata1, ydata2])
  357. ycontrol = np.vstack([ycontrol1, ycontrol2])
  358. ydata = np.tile(ydata, (20, 1))
  359. ycontrol = np.tile(ycontrol, (20, 1))
  360. result = _apply_window(ydata, mlab.window_hanning, axis=1,
  361. return_window=False)
  362. assert_allclose(ycontrol, result, atol=1e-08)
  363. def test_apply_window_hanning_2D_stack_windows_axis1_unflatten(self):
  364. n = 32
  365. ydata = np.arange(n)
  366. ydata1 = ydata+5
  367. ydata2 = ydata+3.3
  368. ycontrol1 = _apply_window(ydata1, mlab.window_hanning)
  369. ycontrol2 = mlab.window_hanning(ydata2)
  370. ydata = np.vstack([ydata1, ydata2])
  371. ycontrol = np.vstack([ycontrol1, ycontrol2])
  372. ydata = np.tile(ydata, (20, 1))
  373. ycontrol = np.tile(ycontrol, (20, 1))
  374. ydata = ydata.flatten()
  375. ydata1 = mlab.stride_windows(ydata, 32, noverlap=0, axis=0)
  376. result = _apply_window(ydata1, mlab.window_hanning, axis=0,
  377. return_window=False)
  378. assert_allclose(ycontrol.T, result, atol=1e-08)
  379. class TestDetrend:
  380. def setup(self):
  381. np.random.seed(0)
  382. n = 1000
  383. x = np.linspace(0., 100, n)
  384. self.sig_zeros = np.zeros(n)
  385. self.sig_off = self.sig_zeros + 100.
  386. self.sig_slope = np.linspace(-10., 90., n)
  387. self.sig_slope_mean = x - x.mean()
  388. sig_rand = np.random.standard_normal(n)
  389. sig_sin = np.sin(x*2*np.pi/(n/100))
  390. sig_rand -= sig_rand.mean()
  391. sig_sin -= sig_sin.mean()
  392. self.sig_base = sig_rand + sig_sin
  393. self.atol = 1e-08
  394. def test_detrend_none_0D_zeros(self):
  395. input = 0.
  396. targ = input
  397. mlab.detrend_none(input)
  398. assert input == targ
  399. def test_detrend_none_0D_zeros_axis1(self):
  400. input = 0.
  401. targ = input
  402. mlab.detrend_none(input, axis=1)
  403. assert input == targ
  404. def test_detrend_str_none_0D_zeros(self):
  405. input = 0.
  406. targ = input
  407. mlab.detrend(input, key='none')
  408. assert input == targ
  409. def test_detrend_detrend_none_0D_zeros(self):
  410. input = 0.
  411. targ = input
  412. mlab.detrend(input, key=mlab.detrend_none)
  413. assert input == targ
  414. def test_detrend_none_0D_off(self):
  415. input = 5.5
  416. targ = input
  417. mlab.detrend_none(input)
  418. assert input == targ
  419. def test_detrend_none_1D_off(self):
  420. input = self.sig_off
  421. targ = input
  422. res = mlab.detrend_none(input)
  423. assert_array_equal(res, targ)
  424. def test_detrend_none_1D_slope(self):
  425. input = self.sig_slope
  426. targ = input
  427. res = mlab.detrend_none(input)
  428. assert_array_equal(res, targ)
  429. def test_detrend_none_1D_base(self):
  430. input = self.sig_base
  431. targ = input
  432. res = mlab.detrend_none(input)
  433. assert_array_equal(res, targ)
  434. def test_detrend_none_1D_base_slope_off_list(self):
  435. input = self.sig_base + self.sig_slope + self.sig_off
  436. targ = input.tolist()
  437. res = mlab.detrend_none(input.tolist())
  438. assert res == targ
  439. def test_detrend_none_2D(self):
  440. arri = [self.sig_base,
  441. self.sig_base + self.sig_off,
  442. self.sig_base + self.sig_slope,
  443. self.sig_base + self.sig_off + self.sig_slope]
  444. input = np.vstack(arri)
  445. targ = input
  446. res = mlab.detrend_none(input)
  447. assert_array_equal(res, targ)
  448. def test_detrend_none_2D_T(self):
  449. arri = [self.sig_base,
  450. self.sig_base + self.sig_off,
  451. self.sig_base + self.sig_slope,
  452. self.sig_base + self.sig_off + self.sig_slope]
  453. input = np.vstack(arri)
  454. targ = input
  455. res = mlab.detrend_none(input.T)
  456. assert_array_equal(res.T, targ)
  457. def test_detrend_mean_0D_zeros(self):
  458. input = 0.
  459. targ = 0.
  460. res = mlab.detrend_mean(input)
  461. assert_almost_equal(res, targ)
  462. def test_detrend_str_mean_0D_zeros(self):
  463. input = 0.
  464. targ = 0.
  465. res = mlab.detrend(input, key='mean')
  466. assert_almost_equal(res, targ)
  467. def test_detrend_detrend_mean_0D_zeros(self):
  468. input = 0.
  469. targ = 0.
  470. res = mlab.detrend(input, key=mlab.detrend_mean)
  471. assert_almost_equal(res, targ)
  472. def test_detrend_mean_0D_off(self):
  473. input = 5.5
  474. targ = 0.
  475. res = mlab.detrend_mean(input)
  476. assert_almost_equal(res, targ)
  477. def test_detrend_str_mean_0D_off(self):
  478. input = 5.5
  479. targ = 0.
  480. res = mlab.detrend(input, key='mean')
  481. assert_almost_equal(res, targ)
  482. def test_detrend_detrend_mean_0D_off(self):
  483. input = 5.5
  484. targ = 0.
  485. res = mlab.detrend(input, key=mlab.detrend_mean)
  486. assert_almost_equal(res, targ)
  487. def test_detrend_mean_1D_zeros(self):
  488. input = self.sig_zeros
  489. targ = self.sig_zeros
  490. res = mlab.detrend_mean(input)
  491. assert_allclose(res, targ, atol=self.atol)
  492. def test_detrend_mean_1D_base(self):
  493. input = self.sig_base
  494. targ = self.sig_base
  495. res = mlab.detrend_mean(input)
  496. assert_allclose(res, targ, atol=self.atol)
  497. def test_detrend_mean_1D_base_off(self):
  498. input = self.sig_base + self.sig_off
  499. targ = self.sig_base
  500. res = mlab.detrend_mean(input)
  501. assert_allclose(res, targ, atol=self.atol)
  502. def test_detrend_mean_1D_base_slope(self):
  503. input = self.sig_base + self.sig_slope
  504. targ = self.sig_base + self.sig_slope_mean
  505. res = mlab.detrend_mean(input)
  506. assert_allclose(res, targ, atol=self.atol)
  507. def test_detrend_mean_1D_base_slope_off(self):
  508. input = self.sig_base + self.sig_slope + self.sig_off
  509. targ = self.sig_base + self.sig_slope_mean
  510. res = mlab.detrend_mean(input)
  511. assert_allclose(res, targ, atol=1e-08)
  512. def test_detrend_mean_1D_base_slope_off_axis0(self):
  513. input = self.sig_base + self.sig_slope + self.sig_off
  514. targ = self.sig_base + self.sig_slope_mean
  515. res = mlab.detrend_mean(input, axis=0)
  516. assert_allclose(res, targ, atol=1e-08)
  517. def test_detrend_mean_1D_base_slope_off_list(self):
  518. input = self.sig_base + self.sig_slope + self.sig_off
  519. targ = self.sig_base + self.sig_slope_mean
  520. res = mlab.detrend_mean(input.tolist())
  521. assert_allclose(res, targ, atol=1e-08)
  522. def test_detrend_mean_1D_base_slope_off_list_axis0(self):
  523. input = self.sig_base + self.sig_slope + self.sig_off
  524. targ = self.sig_base + self.sig_slope_mean
  525. res = mlab.detrend_mean(input.tolist(), axis=0)
  526. assert_allclose(res, targ, atol=1e-08)
  527. def test_demean_0D_off(self):
  528. input = 5.5
  529. targ = 0.
  530. with pytest.warns(MatplotlibDeprecationWarning):
  531. res = mlab.demean(input, axis=None)
  532. assert_almost_equal(res, targ)
  533. def test_demean_1D_base_slope_off(self):
  534. input = self.sig_base + self.sig_slope + self.sig_off
  535. targ = self.sig_base + self.sig_slope_mean
  536. with pytest.warns(MatplotlibDeprecationWarning):
  537. res = mlab.demean(input)
  538. assert_allclose(res, targ, atol=1e-08)
  539. def test_demean_1D_base_slope_off_axis0(self):
  540. input = self.sig_base + self.sig_slope + self.sig_off
  541. targ = self.sig_base + self.sig_slope_mean
  542. with pytest.warns(MatplotlibDeprecationWarning):
  543. res = mlab.demean(input, axis=0)
  544. assert_allclose(res, targ, atol=1e-08)
  545. def test_demean_1D_base_slope_off_list(self):
  546. input = self.sig_base + self.sig_slope + self.sig_off
  547. targ = self.sig_base + self.sig_slope_mean
  548. with pytest.warns(MatplotlibDeprecationWarning):
  549. res = mlab.demean(input.tolist())
  550. assert_allclose(res, targ, atol=1e-08)
  551. def test_detrend_mean_2D_default(self):
  552. arri = [self.sig_off,
  553. self.sig_base + self.sig_off]
  554. arrt = [self.sig_zeros,
  555. self.sig_base]
  556. input = np.vstack(arri)
  557. targ = np.vstack(arrt)
  558. res = mlab.detrend_mean(input)
  559. assert_allclose(res, targ, atol=1e-08)
  560. def test_detrend_mean_2D_none(self):
  561. arri = [self.sig_off,
  562. self.sig_base + self.sig_off]
  563. arrt = [self.sig_zeros,
  564. self.sig_base]
  565. input = np.vstack(arri)
  566. targ = np.vstack(arrt)
  567. res = mlab.detrend_mean(input, axis=None)
  568. assert_allclose(res, targ,
  569. atol=1e-08)
  570. def test_detrend_mean_2D_none_T(self):
  571. arri = [self.sig_off,
  572. self.sig_base + self.sig_off]
  573. arrt = [self.sig_zeros,
  574. self.sig_base]
  575. input = np.vstack(arri).T
  576. targ = np.vstack(arrt)
  577. res = mlab.detrend_mean(input, axis=None)
  578. assert_allclose(res.T, targ,
  579. atol=1e-08)
  580. def test_detrend_mean_2D_axis0(self):
  581. arri = [self.sig_base,
  582. self.sig_base + self.sig_off,
  583. self.sig_base + self.sig_slope,
  584. self.sig_base + self.sig_off + self.sig_slope]
  585. arrt = [self.sig_base,
  586. self.sig_base,
  587. self.sig_base + self.sig_slope_mean,
  588. self.sig_base + self.sig_slope_mean]
  589. input = np.vstack(arri).T
  590. targ = np.vstack(arrt).T
  591. res = mlab.detrend_mean(input, axis=0)
  592. assert_allclose(res, targ,
  593. atol=1e-08)
  594. def test_detrend_mean_2D_axis1(self):
  595. arri = [self.sig_base,
  596. self.sig_base + self.sig_off,
  597. self.sig_base + self.sig_slope,
  598. self.sig_base + self.sig_off + self.sig_slope]
  599. arrt = [self.sig_base,
  600. self.sig_base,
  601. self.sig_base + self.sig_slope_mean,
  602. self.sig_base + self.sig_slope_mean]
  603. input = np.vstack(arri)
  604. targ = np.vstack(arrt)
  605. res = mlab.detrend_mean(input, axis=1)
  606. assert_allclose(res, targ,
  607. atol=1e-08)
  608. def test_detrend_mean_2D_axism1(self):
  609. arri = [self.sig_base,
  610. self.sig_base + self.sig_off,
  611. self.sig_base + self.sig_slope,
  612. self.sig_base + self.sig_off + self.sig_slope]
  613. arrt = [self.sig_base,
  614. self.sig_base,
  615. self.sig_base + self.sig_slope_mean,
  616. self.sig_base + self.sig_slope_mean]
  617. input = np.vstack(arri)
  618. targ = np.vstack(arrt)
  619. res = mlab.detrend_mean(input, axis=-1)
  620. assert_allclose(res, targ,
  621. atol=1e-08)
  622. def test_detrend_2D_default(self):
  623. arri = [self.sig_off,
  624. self.sig_base + self.sig_off]
  625. arrt = [self.sig_zeros,
  626. self.sig_base]
  627. input = np.vstack(arri)
  628. targ = np.vstack(arrt)
  629. res = mlab.detrend(input)
  630. assert_allclose(res, targ, atol=1e-08)
  631. def test_detrend_2D_none(self):
  632. arri = [self.sig_off,
  633. self.sig_base + self.sig_off]
  634. arrt = [self.sig_zeros,
  635. self.sig_base]
  636. input = np.vstack(arri)
  637. targ = np.vstack(arrt)
  638. res = mlab.detrend(input, axis=None)
  639. assert_allclose(res, targ, atol=1e-08)
  640. def test_detrend_str_mean_2D_axis0(self):
  641. arri = [self.sig_base,
  642. self.sig_base + self.sig_off,
  643. self.sig_base + self.sig_slope,
  644. self.sig_base + self.sig_off + self.sig_slope]
  645. arrt = [self.sig_base,
  646. self.sig_base,
  647. self.sig_base + self.sig_slope_mean,
  648. self.sig_base + self.sig_slope_mean]
  649. input = np.vstack(arri).T
  650. targ = np.vstack(arrt).T
  651. res = mlab.detrend(input, key='mean', axis=0)
  652. assert_allclose(res, targ,
  653. atol=1e-08)
  654. def test_detrend_str_constant_2D_none_T(self):
  655. arri = [self.sig_off,
  656. self.sig_base + self.sig_off]
  657. arrt = [self.sig_zeros,
  658. self.sig_base]
  659. input = np.vstack(arri).T
  660. targ = np.vstack(arrt)
  661. res = mlab.detrend(input, key='constant', axis=None)
  662. assert_allclose(res.T, targ,
  663. atol=1e-08)
  664. def test_detrend_str_default_2D_axis1(self):
  665. arri = [self.sig_base,
  666. self.sig_base + self.sig_off,
  667. self.sig_base + self.sig_slope,
  668. self.sig_base + self.sig_off + self.sig_slope]
  669. arrt = [self.sig_base,
  670. self.sig_base,
  671. self.sig_base + self.sig_slope_mean,
  672. self.sig_base + self.sig_slope_mean]
  673. input = np.vstack(arri)
  674. targ = np.vstack(arrt)
  675. res = mlab.detrend(input, key='default', axis=1)
  676. assert_allclose(res, targ,
  677. atol=1e-08)
  678. def test_detrend_detrend_mean_2D_axis0(self):
  679. arri = [self.sig_base,
  680. self.sig_base + self.sig_off,
  681. self.sig_base + self.sig_slope,
  682. self.sig_base + self.sig_off + self.sig_slope]
  683. arrt = [self.sig_base,
  684. self.sig_base,
  685. self.sig_base + self.sig_slope_mean,
  686. self.sig_base + self.sig_slope_mean]
  687. input = np.vstack(arri).T
  688. targ = np.vstack(arrt).T
  689. res = mlab.detrend(input, key=mlab.detrend_mean, axis=0)
  690. assert_allclose(res, targ,
  691. atol=1e-08)
  692. def test_demean_2D_default(self):
  693. arri = [self.sig_base,
  694. self.sig_base + self.sig_off,
  695. self.sig_base + self.sig_slope,
  696. self.sig_base + self.sig_off + self.sig_slope]
  697. arrt = [self.sig_base,
  698. self.sig_base,
  699. self.sig_base + self.sig_slope_mean,
  700. self.sig_base + self.sig_slope_mean]
  701. input = np.vstack(arri).T
  702. targ = np.vstack(arrt).T
  703. with pytest.warns(MatplotlibDeprecationWarning):
  704. res = mlab.demean(input)
  705. assert_allclose(res, targ,
  706. atol=1e-08)
  707. def test_demean_2D_none(self):
  708. arri = [self.sig_off,
  709. self.sig_base + self.sig_off]
  710. arrt = [self.sig_zeros,
  711. self.sig_base]
  712. input = np.vstack(arri)
  713. targ = np.vstack(arrt)
  714. with pytest.warns(MatplotlibDeprecationWarning):
  715. res = mlab.demean(input, axis=None)
  716. assert_allclose(res, targ,
  717. atol=1e-08)
  718. def test_demean_2D_axis0(self):
  719. arri = [self.sig_base,
  720. self.sig_base + self.sig_off,
  721. self.sig_base + self.sig_slope,
  722. self.sig_base + self.sig_off + self.sig_slope]
  723. arrt = [self.sig_base,
  724. self.sig_base,
  725. self.sig_base + self.sig_slope_mean,
  726. self.sig_base + self.sig_slope_mean]
  727. input = np.vstack(arri).T
  728. targ = np.vstack(arrt).T
  729. with pytest.warns(MatplotlibDeprecationWarning):
  730. res = mlab.demean(input, axis=0)
  731. assert_allclose(res, targ,
  732. atol=1e-08)
  733. def test_demean_2D_axis1(self):
  734. arri = [self.sig_base,
  735. self.sig_base + self.sig_off,
  736. self.sig_base + self.sig_slope,
  737. self.sig_base + self.sig_off + self.sig_slope]
  738. arrt = [self.sig_base,
  739. self.sig_base,
  740. self.sig_base + self.sig_slope_mean,
  741. self.sig_base + self.sig_slope_mean]
  742. input = np.vstack(arri)
  743. targ = np.vstack(arrt)
  744. with pytest.warns(MatplotlibDeprecationWarning):
  745. res = mlab.demean(input, axis=1)
  746. assert_allclose(res, targ,
  747. atol=1e-08)
  748. def test_demean_2D_axism1(self):
  749. arri = [self.sig_base,
  750. self.sig_base + self.sig_off,
  751. self.sig_base + self.sig_slope,
  752. self.sig_base + self.sig_off + self.sig_slope]
  753. arrt = [self.sig_base,
  754. self.sig_base,
  755. self.sig_base + self.sig_slope_mean,
  756. self.sig_base + self.sig_slope_mean]
  757. input = np.vstack(arri)
  758. targ = np.vstack(arrt)
  759. with pytest.warns(MatplotlibDeprecationWarning):
  760. res = mlab.demean(input, axis=-1)
  761. assert_allclose(res, targ,
  762. atol=1e-08)
  763. def test_detrend_bad_key_str_ValueError(self):
  764. input = self.sig_slope[np.newaxis]
  765. with pytest.raises(ValueError):
  766. mlab.detrend(input, key='spam')
  767. def test_detrend_bad_key_var_ValueError(self):
  768. input = self.sig_slope[np.newaxis]
  769. with pytest.raises(ValueError):
  770. mlab.detrend(input, key=5)
  771. def test_detrend_mean_0D_d0_ValueError(self):
  772. input = 5.5
  773. with pytest.raises(ValueError):
  774. mlab.detrend_mean(input, axis=0)
  775. def test_detrend_0D_d0_ValueError(self):
  776. input = 5.5
  777. with pytest.raises(ValueError):
  778. mlab.detrend(input, axis=0)
  779. def test_detrend_mean_1D_d1_ValueError(self):
  780. input = self.sig_slope
  781. with pytest.raises(ValueError):
  782. mlab.detrend_mean(input, axis=1)
  783. def test_detrend_1D_d1_ValueError(self):
  784. input = self.sig_slope
  785. with pytest.raises(ValueError):
  786. mlab.detrend(input, axis=1)
  787. def test_demean_1D_d1_ValueError(self):
  788. input = self.sig_slope
  789. with pytest.raises(ValueError), \
  790. pytest.warns(MatplotlibDeprecationWarning):
  791. mlab.demean(input, axis=1)
  792. def test_detrend_mean_2D_d2_ValueError(self):
  793. input = self.sig_slope[np.newaxis]
  794. with pytest.raises(ValueError):
  795. mlab.detrend_mean(input, axis=2)
  796. def test_detrend_2D_d2_ValueError(self):
  797. input = self.sig_slope[np.newaxis]
  798. with pytest.raises(ValueError):
  799. mlab.detrend(input, axis=2)
  800. def test_demean_2D_d2_ValueError(self):
  801. input = self.sig_slope[np.newaxis]
  802. with pytest.raises(ValueError), \
  803. pytest.warns(MatplotlibDeprecationWarning):
  804. mlab.demean(input, axis=2)
  805. def test_detrend_linear_0D_zeros(self):
  806. input = 0.
  807. targ = 0.
  808. res = mlab.detrend_linear(input)
  809. assert_almost_equal(res, targ)
  810. def test_detrend_linear_0D_off(self):
  811. input = 5.5
  812. targ = 0.
  813. res = mlab.detrend_linear(input)
  814. assert_almost_equal(res, targ)
  815. def test_detrend_str_linear_0D_off(self):
  816. input = 5.5
  817. targ = 0.
  818. res = mlab.detrend(input, key='linear')
  819. assert_almost_equal(res, targ)
  820. def test_detrend_detrend_linear_0D_off(self):
  821. input = 5.5
  822. targ = 0.
  823. res = mlab.detrend(input, key=mlab.detrend_linear)
  824. assert_almost_equal(res, targ)
  825. def test_detrend_linear_1d_off(self):
  826. input = self.sig_off
  827. targ = self.sig_zeros
  828. res = mlab.detrend_linear(input)
  829. assert_allclose(res, targ, atol=self.atol)
  830. def test_detrend_linear_1d_slope(self):
  831. input = self.sig_slope
  832. targ = self.sig_zeros
  833. res = mlab.detrend_linear(input)
  834. assert_allclose(res, targ, atol=self.atol)
  835. def test_detrend_linear_1d_slope_off(self):
  836. input = self.sig_slope + self.sig_off
  837. targ = self.sig_zeros
  838. res = mlab.detrend_linear(input)
  839. assert_allclose(res, targ, atol=self.atol)
  840. def test_detrend_str_linear_1d_slope_off(self):
  841. input = self.sig_slope + self.sig_off
  842. targ = self.sig_zeros
  843. res = mlab.detrend(input, key='linear')
  844. assert_allclose(res, targ, atol=self.atol)
  845. def test_detrend_detrend_linear_1d_slope_off(self):
  846. input = self.sig_slope + self.sig_off
  847. targ = self.sig_zeros
  848. res = mlab.detrend(input, key=mlab.detrend_linear)
  849. assert_allclose(res, targ, atol=self.atol)
  850. def test_detrend_linear_1d_slope_off_list(self):
  851. input = self.sig_slope + self.sig_off
  852. targ = self.sig_zeros
  853. res = mlab.detrend_linear(input.tolist())
  854. assert_allclose(res, targ, atol=self.atol)
  855. def test_detrend_linear_2D_ValueError(self):
  856. input = self.sig_slope[np.newaxis]
  857. with pytest.raises(ValueError):
  858. mlab.detrend_linear(input)
  859. def test_detrend_str_linear_2d_slope_off_axis0(self):
  860. arri = [self.sig_off,
  861. self.sig_slope,
  862. self.sig_slope + self.sig_off]
  863. arrt = [self.sig_zeros,
  864. self.sig_zeros,
  865. self.sig_zeros]
  866. input = np.vstack(arri).T
  867. targ = np.vstack(arrt).T
  868. res = mlab.detrend(input, key='linear', axis=0)
  869. assert_allclose(res, targ, atol=self.atol)
  870. def test_detrend_detrend_linear_1d_slope_off_axis1(self):
  871. arri = [self.sig_off,
  872. self.sig_slope,
  873. self.sig_slope + self.sig_off]
  874. arrt = [self.sig_zeros,
  875. self.sig_zeros,
  876. self.sig_zeros]
  877. input = np.vstack(arri).T
  878. targ = np.vstack(arrt).T
  879. res = mlab.detrend(input, key=mlab.detrend_linear, axis=0)
  880. assert_allclose(res, targ, atol=self.atol)
  881. def test_detrend_str_linear_2d_slope_off_axis0_notranspose(self):
  882. arri = [self.sig_off,
  883. self.sig_slope,
  884. self.sig_slope + self.sig_off]
  885. arrt = [self.sig_zeros,
  886. self.sig_zeros,
  887. self.sig_zeros]
  888. input = np.vstack(arri)
  889. targ = np.vstack(arrt)
  890. res = mlab.detrend(input, key='linear', axis=1)
  891. assert_allclose(res, targ, atol=self.atol)
  892. def test_detrend_detrend_linear_1d_slope_off_axis1_notranspose(self):
  893. arri = [self.sig_off,
  894. self.sig_slope,
  895. self.sig_slope + self.sig_off]
  896. arrt = [self.sig_zeros,
  897. self.sig_zeros,
  898. self.sig_zeros]
  899. input = np.vstack(arri)
  900. targ = np.vstack(arrt)
  901. res = mlab.detrend(input, key=mlab.detrend_linear, axis=1)
  902. assert_allclose(res, targ, atol=self.atol)
  903. @pytest.mark.parametrize('iscomplex', [False, True],
  904. ids=['real', 'complex'], scope='class')
  905. @pytest.mark.parametrize('sides', ['onesided', 'twosided', 'default'],
  906. scope='class')
  907. @pytest.mark.parametrize(
  908. 'fstims,len_x,NFFT_density,nover_density,pad_to_density,pad_to_spectrum',
  909. [
  910. ([], None, -1, -1, -1, -1),
  911. ([4], None, -1, -1, -1, -1),
  912. ([4, 5, 10], None, -1, -1, -1, -1),
  913. ([], None, None, -1, -1, None),
  914. ([], None, -1, -1, None, None),
  915. ([], None, None, -1, None, None),
  916. ([], 1024, 512, -1, -1, 128),
  917. ([], 256, -1, -1, 33, 257),
  918. ([], 255, 33, -1, -1, None),
  919. ([], 256, 128, -1, 256, 256),
  920. ([], None, -1, 32, -1, -1),
  921. ],
  922. ids=[
  923. 'nosig',
  924. 'Fs4',
  925. 'FsAll',
  926. 'nosig_noNFFT',
  927. 'nosig_nopad_to',
  928. 'nosig_noNFFT_no_pad_to',
  929. 'nosig_trim',
  930. 'nosig_odd',
  931. 'nosig_oddlen',
  932. 'nosig_stretch',
  933. 'nosig_overlap',
  934. ],
  935. scope='class')
  936. class TestSpectral:
  937. @pytest.fixture(scope='class', autouse=True)
  938. def stim(self, request, fstims, iscomplex, sides, len_x, NFFT_density,
  939. nover_density, pad_to_density, pad_to_spectrum):
  940. Fs = 100.
  941. x = np.arange(0, 10, 1 / Fs)
  942. if len_x is not None:
  943. x = x[:len_x]
  944. # get the stimulus frequencies, defaulting to None
  945. fstims = [Fs / fstim for fstim in fstims]
  946. # get the constants, default to calculated values
  947. if NFFT_density is None:
  948. NFFT_density_real = 256
  949. elif NFFT_density < 0:
  950. NFFT_density_real = NFFT_density = 100
  951. else:
  952. NFFT_density_real = NFFT_density
  953. if nover_density is None:
  954. nover_density_real = 0
  955. elif nover_density < 0:
  956. nover_density_real = nover_density = NFFT_density_real // 2
  957. else:
  958. nover_density_real = nover_density
  959. if pad_to_density is None:
  960. pad_to_density_real = NFFT_density_real
  961. elif pad_to_density < 0:
  962. pad_to_density = int(2**np.ceil(np.log2(NFFT_density_real)))
  963. pad_to_density_real = pad_to_density
  964. else:
  965. pad_to_density_real = pad_to_density
  966. if pad_to_spectrum is None:
  967. pad_to_spectrum_real = len(x)
  968. elif pad_to_spectrum < 0:
  969. pad_to_spectrum_real = pad_to_spectrum = len(x)
  970. else:
  971. pad_to_spectrum_real = pad_to_spectrum
  972. if pad_to_spectrum is None:
  973. NFFT_spectrum_real = NFFT_spectrum = pad_to_spectrum_real
  974. else:
  975. NFFT_spectrum_real = NFFT_spectrum = len(x)
  976. nover_spectrum = 0
  977. NFFT_specgram = NFFT_density
  978. nover_specgram = nover_density
  979. pad_to_specgram = pad_to_density
  980. NFFT_specgram_real = NFFT_density_real
  981. nover_specgram_real = nover_density_real
  982. if sides == 'onesided' or (sides == 'default' and not iscomplex):
  983. # frequencies for specgram, psd, and csd
  984. # need to handle even and odd differently
  985. if pad_to_density_real % 2:
  986. freqs_density = np.linspace(0, Fs / 2,
  987. num=pad_to_density_real,
  988. endpoint=False)[::2]
  989. else:
  990. freqs_density = np.linspace(0, Fs / 2,
  991. num=pad_to_density_real // 2 + 1)
  992. # frequencies for complex, magnitude, angle, and phase spectrums
  993. # need to handle even and odd differently
  994. if pad_to_spectrum_real % 2:
  995. freqs_spectrum = np.linspace(0, Fs / 2,
  996. num=pad_to_spectrum_real,
  997. endpoint=False)[::2]
  998. else:
  999. freqs_spectrum = np.linspace(0, Fs / 2,
  1000. num=pad_to_spectrum_real // 2 + 1)
  1001. else:
  1002. # frequencies for specgram, psd, and csd
  1003. # need to handle even and odd differentl
  1004. if pad_to_density_real % 2:
  1005. freqs_density = np.linspace(-Fs / 2, Fs / 2,
  1006. num=2 * pad_to_density_real,
  1007. endpoint=False)[1::2]
  1008. else:
  1009. freqs_density = np.linspace(-Fs / 2, Fs / 2,
  1010. num=pad_to_density_real,
  1011. endpoint=False)
  1012. # frequencies for complex, magnitude, angle, and phase spectrums
  1013. # need to handle even and odd differently
  1014. if pad_to_spectrum_real % 2:
  1015. freqs_spectrum = np.linspace(-Fs / 2, Fs / 2,
  1016. num=2 * pad_to_spectrum_real,
  1017. endpoint=False)[1::2]
  1018. else:
  1019. freqs_spectrum = np.linspace(-Fs / 2, Fs / 2,
  1020. num=pad_to_spectrum_real,
  1021. endpoint=False)
  1022. freqs_specgram = freqs_density
  1023. # time points for specgram
  1024. t_start = NFFT_specgram_real // 2
  1025. t_stop = len(x) - NFFT_specgram_real // 2 + 1
  1026. t_step = NFFT_specgram_real - nover_specgram_real
  1027. t_specgram = x[t_start:t_stop:t_step]
  1028. if NFFT_specgram_real % 2:
  1029. t_specgram += 1 / Fs / 2
  1030. if len(t_specgram) == 0:
  1031. t_specgram = np.array([NFFT_specgram_real / (2 * Fs)])
  1032. t_spectrum = np.array([NFFT_spectrum_real / (2 * Fs)])
  1033. t_density = t_specgram
  1034. y = np.zeros_like(x)
  1035. for i, fstim in enumerate(fstims):
  1036. y += np.sin(fstim * x * np.pi * 2) * 10**i
  1037. if iscomplex:
  1038. y = y.astype('complex')
  1039. # Interestingly, the instance on which this fixture is called is not
  1040. # the same as the one on which a test is run. So we need to modify the
  1041. # class itself when using a class-scoped fixture.
  1042. cls = request.cls
  1043. cls.Fs = Fs
  1044. cls.sides = sides
  1045. cls.fstims = fstims
  1046. cls.NFFT_density = NFFT_density
  1047. cls.nover_density = nover_density
  1048. cls.pad_to_density = pad_to_density
  1049. cls.NFFT_spectrum = NFFT_spectrum
  1050. cls.nover_spectrum = nover_spectrum
  1051. cls.pad_to_spectrum = pad_to_spectrum
  1052. cls.NFFT_specgram = NFFT_specgram
  1053. cls.nover_specgram = nover_specgram
  1054. cls.pad_to_specgram = pad_to_specgram
  1055. cls.t_specgram = t_specgram
  1056. cls.t_density = t_density
  1057. cls.t_spectrum = t_spectrum
  1058. cls.y = y
  1059. cls.freqs_density = freqs_density
  1060. cls.freqs_spectrum = freqs_spectrum
  1061. cls.freqs_specgram = freqs_specgram
  1062. cls.NFFT_density_real = NFFT_density_real
  1063. def check_freqs(self, vals, targfreqs, resfreqs, fstims):
  1064. assert resfreqs.argmin() == 0
  1065. assert resfreqs.argmax() == len(resfreqs)-1
  1066. assert_allclose(resfreqs, targfreqs, atol=1e-06)
  1067. for fstim in fstims:
  1068. i = np.abs(resfreqs - fstim).argmin()
  1069. assert vals[i] > vals[i+2]
  1070. assert vals[i] > vals[i-2]
  1071. def check_maxfreq(self, spec, fsp, fstims):
  1072. # skip the test if there are no frequencies
  1073. if len(fstims) == 0:
  1074. return
  1075. # if twosided, do the test for each side
  1076. if fsp.min() < 0:
  1077. fspa = np.abs(fsp)
  1078. zeroind = fspa.argmin()
  1079. self.check_maxfreq(spec[:zeroind], fspa[:zeroind], fstims)
  1080. self.check_maxfreq(spec[zeroind:], fspa[zeroind:], fstims)
  1081. return
  1082. fstimst = fstims[:]
  1083. spect = spec.copy()
  1084. # go through each peak and make sure it is correctly the maximum peak
  1085. while fstimst:
  1086. maxind = spect.argmax()
  1087. maxfreq = fsp[maxind]
  1088. assert_almost_equal(maxfreq, fstimst[-1])
  1089. del fstimst[-1]
  1090. spect[maxind-5:maxind+5] = 0
  1091. def test_spectral_helper_raises_complex_same_data(self):
  1092. # test that mode 'complex' cannot be used if x is not y
  1093. with pytest.raises(ValueError):
  1094. mlab._spectral_helper(x=self.y, y=self.y+1, mode='complex')
  1095. def test_spectral_helper_raises_magnitude_same_data(self):
  1096. # test that mode 'magnitude' cannot be used if x is not y
  1097. with pytest.raises(ValueError):
  1098. mlab._spectral_helper(x=self.y, y=self.y+1, mode='magnitude')
  1099. def test_spectral_helper_raises_angle_same_data(self):
  1100. # test that mode 'angle' cannot be used if x is not y
  1101. with pytest.raises(ValueError):
  1102. mlab._spectral_helper(x=self.y, y=self.y+1, mode='angle')
  1103. def test_spectral_helper_raises_phase_same_data(self):
  1104. # test that mode 'phase' cannot be used if x is not y
  1105. with pytest.raises(ValueError):
  1106. mlab._spectral_helper(x=self.y, y=self.y+1, mode='phase')
  1107. def test_spectral_helper_raises_unknown_mode(self):
  1108. # test that unknown value for mode cannot be used
  1109. with pytest.raises(ValueError):
  1110. mlab._spectral_helper(x=self.y, mode='spam')
  1111. def test_spectral_helper_raises_unknown_sides(self):
  1112. # test that unknown value for sides cannot be used
  1113. with pytest.raises(ValueError):
  1114. mlab._spectral_helper(x=self.y, y=self.y, sides='eggs')
  1115. def test_spectral_helper_raises_noverlap_gt_NFFT(self):
  1116. # test that noverlap cannot be larger than NFFT
  1117. with pytest.raises(ValueError):
  1118. mlab._spectral_helper(x=self.y, y=self.y, NFFT=10, noverlap=20)
  1119. def test_spectral_helper_raises_noverlap_eq_NFFT(self):
  1120. # test that noverlap cannot be equal to NFFT
  1121. with pytest.raises(ValueError):
  1122. mlab._spectral_helper(x=self.y, NFFT=10, noverlap=10)
  1123. def test_spectral_helper_raises_winlen_ne_NFFT(self):
  1124. # test that the window length cannot be different from NFFT
  1125. with pytest.raises(ValueError):
  1126. mlab._spectral_helper(x=self.y, y=self.y, NFFT=10,
  1127. window=np.ones(9))
  1128. def test_single_spectrum_helper_raises_mode_default(self):
  1129. # test that mode 'default' cannot be used with _single_spectrum_helper
  1130. with pytest.raises(ValueError):
  1131. mlab._single_spectrum_helper(x=self.y, mode='default')
  1132. def test_single_spectrum_helper_raises_mode_psd(self):
  1133. # test that mode 'psd' cannot be used with _single_spectrum_helper
  1134. with pytest.raises(ValueError):
  1135. mlab._single_spectrum_helper(x=self.y, mode='psd')
  1136. def test_spectral_helper_psd(self):
  1137. freqs = self.freqs_density
  1138. spec, fsp, t = mlab._spectral_helper(x=self.y, y=self.y,
  1139. NFFT=self.NFFT_density,
  1140. Fs=self.Fs,
  1141. noverlap=self.nover_density,
  1142. pad_to=self.pad_to_density,
  1143. sides=self.sides,
  1144. mode='psd')
  1145. assert_allclose(fsp, freqs, atol=1e-06)
  1146. assert_allclose(t, self.t_density, atol=1e-06)
  1147. assert spec.shape[0] == freqs.shape[0]
  1148. assert spec.shape[1] == self.t_specgram.shape[0]
  1149. def test_spectral_helper_magnitude_specgram(self):
  1150. freqs = self.freqs_specgram
  1151. spec, fsp, t = mlab._spectral_helper(x=self.y, y=self.y,
  1152. NFFT=self.NFFT_specgram,
  1153. Fs=self.Fs,
  1154. noverlap=self.nover_specgram,
  1155. pad_to=self.pad_to_specgram,
  1156. sides=self.sides,
  1157. mode='magnitude')
  1158. assert_allclose(fsp, freqs, atol=1e-06)
  1159. assert_allclose(t, self.t_specgram, atol=1e-06)
  1160. assert spec.shape[0] == freqs.shape[0]
  1161. assert spec.shape[1] == self.t_specgram.shape[0]
  1162. def test_spectral_helper_magnitude_magnitude_spectrum(self):
  1163. freqs = self.freqs_spectrum
  1164. spec, fsp, t = mlab._spectral_helper(x=self.y, y=self.y,
  1165. NFFT=self.NFFT_spectrum,
  1166. Fs=self.Fs,
  1167. noverlap=self.nover_spectrum,
  1168. pad_to=self.pad_to_spectrum,
  1169. sides=self.sides,
  1170. mode='magnitude')
  1171. assert_allclose(fsp, freqs, atol=1e-06)
  1172. assert_allclose(t, self.t_spectrum, atol=1e-06)
  1173. assert spec.shape[0] == freqs.shape[0]
  1174. assert spec.shape[1] == 1
  1175. def test_csd(self):
  1176. freqs = self.freqs_density
  1177. spec, fsp = mlab.csd(x=self.y, y=self.y+1,
  1178. NFFT=self.NFFT_density,
  1179. Fs=self.Fs,
  1180. noverlap=self.nover_density,
  1181. pad_to=self.pad_to_density,
  1182. sides=self.sides)
  1183. assert_allclose(fsp, freqs, atol=1e-06)
  1184. assert spec.shape == freqs.shape
  1185. def test_csd_padding(self):
  1186. """Test zero padding of csd()."""
  1187. if self.NFFT_density is None: # for derived classes
  1188. return
  1189. sargs = dict(x=self.y, y=self.y+1, Fs=self.Fs, window=mlab.window_none,
  1190. sides=self.sides)
  1191. spec0, _ = mlab.csd(NFFT=self.NFFT_density, **sargs)
  1192. spec1, _ = mlab.csd(NFFT=self.NFFT_density*2, **sargs)
  1193. assert_almost_equal(np.sum(np.conjugate(spec0)*spec0).real,
  1194. np.sum(np.conjugate(spec1/2)*spec1/2).real)
  1195. def test_psd(self):
  1196. freqs = self.freqs_density
  1197. spec, fsp = mlab.psd(x=self.y,
  1198. NFFT=self.NFFT_density,
  1199. Fs=self.Fs,
  1200. noverlap=self.nover_density,
  1201. pad_to=self.pad_to_density,
  1202. sides=self.sides)
  1203. assert spec.shape == freqs.shape
  1204. self.check_freqs(spec, freqs, fsp, self.fstims)
  1205. def test_psd_detrend_mean_func_offset(self):
  1206. if self.NFFT_density is None:
  1207. return
  1208. ydata = np.zeros(self.NFFT_density)
  1209. ydata1 = ydata+5
  1210. ydata2 = ydata+3.3
  1211. ydata = np.vstack([ydata1, ydata2])
  1212. ydata = np.tile(ydata, (20, 1))
  1213. ydatab = ydata.T.flatten()
  1214. ydata = ydata.flatten()
  1215. ycontrol = np.zeros_like(ydata)
  1216. spec_g, fsp_g = mlab.psd(x=ydata,
  1217. NFFT=self.NFFT_density,
  1218. Fs=self.Fs,
  1219. noverlap=0,
  1220. sides=self.sides,
  1221. detrend=mlab.detrend_mean)
  1222. spec_b, fsp_b = mlab.psd(x=ydatab,
  1223. NFFT=self.NFFT_density,
  1224. Fs=self.Fs,
  1225. noverlap=0,
  1226. sides=self.sides,
  1227. detrend=mlab.detrend_mean)
  1228. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1229. NFFT=self.NFFT_density,
  1230. Fs=self.Fs,
  1231. noverlap=0,
  1232. sides=self.sides)
  1233. assert_array_equal(fsp_g, fsp_c)
  1234. assert_array_equal(fsp_b, fsp_c)
  1235. assert_allclose(spec_g, spec_c, atol=1e-08)
  1236. # these should not be almost equal
  1237. with pytest.raises(AssertionError):
  1238. assert_allclose(spec_b, spec_c, atol=1e-08)
  1239. def test_psd_detrend_mean_str_offset(self):
  1240. if self.NFFT_density is None:
  1241. return
  1242. ydata = np.zeros(self.NFFT_density)
  1243. ydata1 = ydata+5
  1244. ydata2 = ydata+3.3
  1245. ydata = np.vstack([ydata1, ydata2])
  1246. ydata = np.tile(ydata, (20, 1))
  1247. ydatab = ydata.T.flatten()
  1248. ydata = ydata.flatten()
  1249. ycontrol = np.zeros_like(ydata)
  1250. spec_g, fsp_g = mlab.psd(x=ydata,
  1251. NFFT=self.NFFT_density,
  1252. Fs=self.Fs,
  1253. noverlap=0,
  1254. sides=self.sides,
  1255. detrend='mean')
  1256. spec_b, fsp_b = mlab.psd(x=ydatab,
  1257. NFFT=self.NFFT_density,
  1258. Fs=self.Fs,
  1259. noverlap=0,
  1260. sides=self.sides,
  1261. detrend='mean')
  1262. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1263. NFFT=self.NFFT_density,
  1264. Fs=self.Fs,
  1265. noverlap=0,
  1266. sides=self.sides)
  1267. assert_array_equal(fsp_g, fsp_c)
  1268. assert_array_equal(fsp_b, fsp_c)
  1269. assert_allclose(spec_g, spec_c, atol=1e-08)
  1270. # these should not be almost equal
  1271. with pytest.raises(AssertionError):
  1272. assert_allclose(spec_b, spec_c, atol=1e-08)
  1273. def test_psd_detrend_linear_func_trend(self):
  1274. if self.NFFT_density is None:
  1275. return
  1276. ydata = np.arange(self.NFFT_density)
  1277. ydata1 = ydata+5
  1278. ydata2 = ydata+3.3
  1279. ydata = np.vstack([ydata1, ydata2])
  1280. ydata = np.tile(ydata, (20, 1))
  1281. ydatab = ydata.T.flatten()
  1282. ydata = ydata.flatten()
  1283. ycontrol = np.zeros_like(ydata)
  1284. spec_g, fsp_g = mlab.psd(x=ydata,
  1285. NFFT=self.NFFT_density,
  1286. Fs=self.Fs,
  1287. noverlap=0,
  1288. sides=self.sides,
  1289. detrend=mlab.detrend_linear)
  1290. spec_b, fsp_b = mlab.psd(x=ydatab,
  1291. NFFT=self.NFFT_density,
  1292. Fs=self.Fs,
  1293. noverlap=0,
  1294. sides=self.sides,
  1295. detrend=mlab.detrend_linear)
  1296. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1297. NFFT=self.NFFT_density,
  1298. Fs=self.Fs,
  1299. noverlap=0,
  1300. sides=self.sides)
  1301. assert_array_equal(fsp_g, fsp_c)
  1302. assert_array_equal(fsp_b, fsp_c)
  1303. assert_allclose(spec_g, spec_c, atol=1e-08)
  1304. # these should not be almost equal
  1305. with pytest.raises(AssertionError):
  1306. assert_allclose(spec_b, spec_c, atol=1e-08)
  1307. def test_psd_detrend_linear_str_trend(self):
  1308. if self.NFFT_density is None:
  1309. return
  1310. ydata = np.arange(self.NFFT_density)
  1311. ydata1 = ydata+5
  1312. ydata2 = ydata+3.3
  1313. ydata = np.vstack([ydata1, ydata2])
  1314. ydata = np.tile(ydata, (20, 1))
  1315. ydatab = ydata.T.flatten()
  1316. ydata = ydata.flatten()
  1317. ycontrol = np.zeros_like(ydata)
  1318. spec_g, fsp_g = mlab.psd(x=ydata,
  1319. NFFT=self.NFFT_density,
  1320. Fs=self.Fs,
  1321. noverlap=0,
  1322. sides=self.sides,
  1323. detrend='linear')
  1324. spec_b, fsp_b = mlab.psd(x=ydatab,
  1325. NFFT=self.NFFT_density,
  1326. Fs=self.Fs,
  1327. noverlap=0,
  1328. sides=self.sides,
  1329. detrend='linear')
  1330. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1331. NFFT=self.NFFT_density,
  1332. Fs=self.Fs,
  1333. noverlap=0,
  1334. sides=self.sides)
  1335. assert_array_equal(fsp_g, fsp_c)
  1336. assert_array_equal(fsp_b, fsp_c)
  1337. assert_allclose(spec_g, spec_c, atol=1e-08)
  1338. # these should not be almost equal
  1339. with pytest.raises(AssertionError):
  1340. assert_allclose(spec_b, spec_c, atol=1e-08)
  1341. def test_psd_window_hanning(self):
  1342. if self.NFFT_density is None:
  1343. return
  1344. ydata = np.arange(self.NFFT_density)
  1345. ydata1 = ydata+5
  1346. ydata2 = ydata+3.3
  1347. ycontrol1, windowVals = _apply_window(ydata1,
  1348. mlab.window_hanning,
  1349. return_window=True)
  1350. ycontrol2 = mlab.window_hanning(ydata2)
  1351. ydata = np.vstack([ydata1, ydata2])
  1352. ycontrol = np.vstack([ycontrol1, ycontrol2])
  1353. ydata = np.tile(ydata, (20, 1))
  1354. ycontrol = np.tile(ycontrol, (20, 1))
  1355. ydatab = ydata.T.flatten()
  1356. ydataf = ydata.flatten()
  1357. ycontrol = ycontrol.flatten()
  1358. spec_g, fsp_g = mlab.psd(x=ydataf,
  1359. NFFT=self.NFFT_density,
  1360. Fs=self.Fs,
  1361. noverlap=0,
  1362. sides=self.sides,
  1363. window=mlab.window_hanning)
  1364. spec_b, fsp_b = mlab.psd(x=ydatab,
  1365. NFFT=self.NFFT_density,
  1366. Fs=self.Fs,
  1367. noverlap=0,
  1368. sides=self.sides,
  1369. window=mlab.window_hanning)
  1370. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1371. NFFT=self.NFFT_density,
  1372. Fs=self.Fs,
  1373. noverlap=0,
  1374. sides=self.sides,
  1375. window=mlab.window_none)
  1376. spec_c *= len(ycontrol1)/(np.abs(windowVals)**2).sum()
  1377. assert_array_equal(fsp_g, fsp_c)
  1378. assert_array_equal(fsp_b, fsp_c)
  1379. assert_allclose(spec_g, spec_c, atol=1e-08)
  1380. # these should not be almost equal
  1381. with pytest.raises(AssertionError):
  1382. assert_allclose(spec_b, spec_c, atol=1e-08)
  1383. def test_psd_window_hanning_detrend_linear(self):
  1384. if self.NFFT_density is None:
  1385. return
  1386. ydata = np.arange(self.NFFT_density)
  1387. ycontrol = np.zeros(self.NFFT_density)
  1388. ydata1 = ydata+5
  1389. ydata2 = ydata+3.3
  1390. ycontrol1 = ycontrol
  1391. ycontrol2 = ycontrol
  1392. ycontrol1, windowVals = _apply_window(ycontrol1,
  1393. mlab.window_hanning,
  1394. return_window=True)
  1395. ycontrol2 = mlab.window_hanning(ycontrol2)
  1396. ydata = np.vstack([ydata1, ydata2])
  1397. ycontrol = np.vstack([ycontrol1, ycontrol2])
  1398. ydata = np.tile(ydata, (20, 1))
  1399. ycontrol = np.tile(ycontrol, (20, 1))
  1400. ydatab = ydata.T.flatten()
  1401. ydataf = ydata.flatten()
  1402. ycontrol = ycontrol.flatten()
  1403. spec_g, fsp_g = mlab.psd(x=ydataf,
  1404. NFFT=self.NFFT_density,
  1405. Fs=self.Fs,
  1406. noverlap=0,
  1407. sides=self.sides,
  1408. detrend=mlab.detrend_linear,
  1409. window=mlab.window_hanning)
  1410. spec_b, fsp_b = mlab.psd(x=ydatab,
  1411. NFFT=self.NFFT_density,
  1412. Fs=self.Fs,
  1413. noverlap=0,
  1414. sides=self.sides,
  1415. detrend=mlab.detrend_linear,
  1416. window=mlab.window_hanning)
  1417. spec_c, fsp_c = mlab.psd(x=ycontrol,
  1418. NFFT=self.NFFT_density,
  1419. Fs=self.Fs,
  1420. noverlap=0,
  1421. sides=self.sides,
  1422. window=mlab.window_none)
  1423. spec_c *= len(ycontrol1)/(np.abs(windowVals)**2).sum()
  1424. assert_array_equal(fsp_g, fsp_c)
  1425. assert_array_equal(fsp_b, fsp_c)
  1426. assert_allclose(spec_g, spec_c, atol=1e-08)
  1427. # these should not be almost equal
  1428. with pytest.raises(AssertionError):
  1429. assert_allclose(spec_b, spec_c, atol=1e-08)
  1430. def test_psd_windowarray(self):
  1431. freqs = self.freqs_density
  1432. spec, fsp = mlab.psd(x=self.y,
  1433. NFFT=self.NFFT_density,
  1434. Fs=self.Fs,
  1435. noverlap=self.nover_density,
  1436. pad_to=self.pad_to_density,
  1437. sides=self.sides,
  1438. window=np.ones(self.NFFT_density_real))
  1439. assert_allclose(fsp, freqs, atol=1e-06)
  1440. assert spec.shape == freqs.shape
  1441. def test_psd_windowarray_scale_by_freq(self):
  1442. win = mlab.window_hanning(np.ones(self.NFFT_density_real))
  1443. spec, fsp = mlab.psd(x=self.y,
  1444. NFFT=self.NFFT_density,
  1445. Fs=self.Fs,
  1446. noverlap=self.nover_density,
  1447. pad_to=self.pad_to_density,
  1448. sides=self.sides,
  1449. window=mlab.window_hanning)
  1450. spec_s, fsp_s = mlab.psd(x=self.y,
  1451. NFFT=self.NFFT_density,
  1452. Fs=self.Fs,
  1453. noverlap=self.nover_density,
  1454. pad_to=self.pad_to_density,
  1455. sides=self.sides,
  1456. window=mlab.window_hanning,
  1457. scale_by_freq=True)
  1458. spec_n, fsp_n = mlab.psd(x=self.y,
  1459. NFFT=self.NFFT_density,
  1460. Fs=self.Fs,
  1461. noverlap=self.nover_density,
  1462. pad_to=self.pad_to_density,
  1463. sides=self.sides,
  1464. window=mlab.window_hanning,
  1465. scale_by_freq=False)
  1466. assert_array_equal(fsp, fsp_s)
  1467. assert_array_equal(fsp, fsp_n)
  1468. assert_array_equal(spec, spec_s)
  1469. assert_allclose(spec_s*(win**2).sum(),
  1470. spec_n/self.Fs*win.sum()**2,
  1471. atol=1e-08)
  1472. def test_complex_spectrum(self):
  1473. freqs = self.freqs_spectrum
  1474. spec, fsp = mlab.complex_spectrum(x=self.y,
  1475. Fs=self.Fs,
  1476. sides=self.sides,
  1477. pad_to=self.pad_to_spectrum)
  1478. assert_allclose(fsp, freqs, atol=1e-06)
  1479. assert spec.shape == freqs.shape
  1480. def test_magnitude_spectrum(self):
  1481. freqs = self.freqs_spectrum
  1482. spec, fsp = mlab.magnitude_spectrum(x=self.y,
  1483. Fs=self.Fs,
  1484. sides=self.sides,
  1485. pad_to=self.pad_to_spectrum)
  1486. assert spec.shape == freqs.shape
  1487. self.check_maxfreq(spec, fsp, self.fstims)
  1488. self.check_freqs(spec, freqs, fsp, self.fstims)
  1489. def test_angle_spectrum(self):
  1490. freqs = self.freqs_spectrum
  1491. spec, fsp = mlab.angle_spectrum(x=self.y,
  1492. Fs=self.Fs,
  1493. sides=self.sides,
  1494. pad_to=self.pad_to_spectrum)
  1495. assert_allclose(fsp, freqs, atol=1e-06)
  1496. assert spec.shape == freqs.shape
  1497. def test_phase_spectrum(self):
  1498. freqs = self.freqs_spectrum
  1499. spec, fsp = mlab.phase_spectrum(x=self.y,
  1500. Fs=self.Fs,
  1501. sides=self.sides,
  1502. pad_to=self.pad_to_spectrum)
  1503. assert_allclose(fsp, freqs, atol=1e-06)
  1504. assert spec.shape == freqs.shape
  1505. def test_specgram_auto(self):
  1506. freqs = self.freqs_specgram
  1507. spec, fsp, t = mlab.specgram(x=self.y,
  1508. NFFT=self.NFFT_specgram,
  1509. Fs=self.Fs,
  1510. noverlap=self.nover_specgram,
  1511. pad_to=self.pad_to_specgram,
  1512. sides=self.sides)
  1513. specm = np.mean(spec, axis=1)
  1514. assert_allclose(fsp, freqs, atol=1e-06)
  1515. assert_allclose(t, self.t_specgram, atol=1e-06)
  1516. assert spec.shape[0] == freqs.shape[0]
  1517. assert spec.shape[1] == self.t_specgram.shape[0]
  1518. # since we are using a single freq, all time slices
  1519. # should be about the same
  1520. if np.abs(spec.max()) != 0:
  1521. assert_allclose(np.diff(spec, axis=1).max()/np.abs(spec.max()), 0,
  1522. atol=1e-02)
  1523. self.check_freqs(specm, freqs, fsp, self.fstims)
  1524. def test_specgram_default(self):
  1525. freqs = self.freqs_specgram
  1526. spec, fsp, t = mlab.specgram(x=self.y,
  1527. NFFT=self.NFFT_specgram,
  1528. Fs=self.Fs,
  1529. noverlap=self.nover_specgram,
  1530. pad_to=self.pad_to_specgram,
  1531. sides=self.sides,
  1532. mode='default')
  1533. specm = np.mean(spec, axis=1)
  1534. assert_allclose(fsp, freqs, atol=1e-06)
  1535. assert_allclose(t, self.t_specgram, atol=1e-06)
  1536. assert spec.shape[0] == freqs.shape[0]
  1537. assert spec.shape[1] == self.t_specgram.shape[0]
  1538. # since we are using a single freq, all time slices
  1539. # should be about the same
  1540. if np.abs(spec.max()) != 0:
  1541. assert_allclose(np.diff(spec, axis=1).max()/np.abs(spec.max()), 0,
  1542. atol=1e-02)
  1543. self.check_freqs(specm, freqs, fsp, self.fstims)
  1544. def test_specgram_psd(self):
  1545. freqs = self.freqs_specgram
  1546. spec, fsp, t = mlab.specgram(x=self.y,
  1547. NFFT=self.NFFT_specgram,
  1548. Fs=self.Fs,
  1549. noverlap=self.nover_specgram,
  1550. pad_to=self.pad_to_specgram,
  1551. sides=self.sides,
  1552. mode='psd')
  1553. specm = np.mean(spec, axis=1)
  1554. assert_allclose(fsp, freqs, atol=1e-06)
  1555. assert_allclose(t, self.t_specgram, atol=1e-06)
  1556. assert spec.shape[0] == freqs.shape[0]
  1557. assert spec.shape[1] == self.t_specgram.shape[0]
  1558. # since we are using a single freq, all time slices
  1559. # should be about the same
  1560. if np.abs(spec.max()) != 0:
  1561. assert_allclose(np.diff(spec, axis=1).max()/np.abs(spec.max()), 0,
  1562. atol=1e-02)
  1563. self.check_freqs(specm, freqs, fsp, self.fstims)
  1564. def test_specgram_complex(self):
  1565. freqs = self.freqs_specgram
  1566. spec, fsp, t = mlab.specgram(x=self.y,
  1567. NFFT=self.NFFT_specgram,
  1568. Fs=self.Fs,
  1569. noverlap=self.nover_specgram,
  1570. pad_to=self.pad_to_specgram,
  1571. sides=self.sides,
  1572. mode='complex')
  1573. specm = np.mean(np.abs(spec), axis=1)
  1574. assert_allclose(fsp, freqs, atol=1e-06)
  1575. assert_allclose(t, self.t_specgram, atol=1e-06)
  1576. assert spec.shape[0] == freqs.shape[0]
  1577. assert spec.shape[1] == self.t_specgram.shape[0]
  1578. self.check_freqs(specm, freqs, fsp, self.fstims)
  1579. def test_specgram_magnitude(self):
  1580. freqs = self.freqs_specgram
  1581. spec, fsp, t = mlab.specgram(x=self.y,
  1582. NFFT=self.NFFT_specgram,
  1583. Fs=self.Fs,
  1584. noverlap=self.nover_specgram,
  1585. pad_to=self.pad_to_specgram,
  1586. sides=self.sides,
  1587. mode='magnitude')
  1588. specm = np.mean(spec, axis=1)
  1589. assert_allclose(fsp, freqs, atol=1e-06)
  1590. assert_allclose(t, self.t_specgram, atol=1e-06)
  1591. assert spec.shape[0] == freqs.shape[0]
  1592. assert spec.shape[1] == self.t_specgram.shape[0]
  1593. # since we are using a single freq, all time slices
  1594. # should be about the same
  1595. if np.abs(spec.max()) != 0:
  1596. assert_allclose(np.diff(spec, axis=1).max()/np.abs(spec.max()), 0,
  1597. atol=1e-02)
  1598. self.check_freqs(specm, freqs, fsp, self.fstims)
  1599. def test_specgram_angle(self):
  1600. freqs = self.freqs_specgram
  1601. spec, fsp, t = mlab.specgram(x=self.y,
  1602. NFFT=self.NFFT_specgram,
  1603. Fs=self.Fs,
  1604. noverlap=self.nover_specgram,
  1605. pad_to=self.pad_to_specgram,
  1606. sides=self.sides,
  1607. mode='angle')
  1608. assert_allclose(fsp, freqs, atol=1e-06)
  1609. assert_allclose(t, self.t_specgram, atol=1e-06)
  1610. assert spec.shape[0] == freqs.shape[0]
  1611. assert spec.shape[1] == self.t_specgram.shape[0]
  1612. def test_specgram_phase(self):
  1613. freqs = self.freqs_specgram
  1614. spec, fsp, t = mlab.specgram(x=self.y,
  1615. NFFT=self.NFFT_specgram,
  1616. Fs=self.Fs,
  1617. noverlap=self.nover_specgram,
  1618. pad_to=self.pad_to_specgram,
  1619. sides=self.sides,
  1620. mode='phase')
  1621. assert_allclose(fsp, freqs, atol=1e-06)
  1622. assert_allclose(t, self.t_specgram, atol=1e-06)
  1623. assert spec.shape[0] == freqs.shape[0]
  1624. assert spec.shape[1] == self.t_specgram.shape[0]
  1625. def test_specgram_warn_only1seg(self):
  1626. """Warning should be raised if len(x) <= NFFT."""
  1627. with pytest.warns(UserWarning, match="Only one segment is calculated"):
  1628. mlab.specgram(x=self.y, NFFT=len(self.y), Fs=self.Fs)
  1629. def test_psd_csd_equal(self):
  1630. Pxx, freqsxx = mlab.psd(x=self.y,
  1631. NFFT=self.NFFT_density,
  1632. Fs=self.Fs,
  1633. noverlap=self.nover_density,
  1634. pad_to=self.pad_to_density,
  1635. sides=self.sides)
  1636. Pxy, freqsxy = mlab.csd(x=self.y, y=self.y,
  1637. NFFT=self.NFFT_density,
  1638. Fs=self.Fs,
  1639. noverlap=self.nover_density,
  1640. pad_to=self.pad_to_density,
  1641. sides=self.sides)
  1642. assert_array_almost_equal_nulp(Pxx, Pxy)
  1643. assert_array_equal(freqsxx, freqsxy)
  1644. def test_specgram_auto_default_equal(self):
  1645. '''test that mlab.specgram without mode and with mode 'default' and
  1646. 'psd' are all the same'''
  1647. speca, freqspeca, ta = mlab.specgram(x=self.y,
  1648. NFFT=self.NFFT_specgram,
  1649. Fs=self.Fs,
  1650. noverlap=self.nover_specgram,
  1651. pad_to=self.pad_to_specgram,
  1652. sides=self.sides)
  1653. specb, freqspecb, tb = mlab.specgram(x=self.y,
  1654. NFFT=self.NFFT_specgram,
  1655. Fs=self.Fs,
  1656. noverlap=self.nover_specgram,
  1657. pad_to=self.pad_to_specgram,
  1658. sides=self.sides,
  1659. mode='default')
  1660. assert_array_equal(speca, specb)
  1661. assert_array_equal(freqspeca, freqspecb)
  1662. assert_array_equal(ta, tb)
  1663. def test_specgram_auto_psd_equal(self):
  1664. '''test that mlab.specgram without mode and with mode 'default' and
  1665. 'psd' are all the same'''
  1666. speca, freqspeca, ta = mlab.specgram(x=self.y,
  1667. NFFT=self.NFFT_specgram,
  1668. Fs=self.Fs,
  1669. noverlap=self.nover_specgram,
  1670. pad_to=self.pad_to_specgram,
  1671. sides=self.sides)
  1672. specc, freqspecc, tc = mlab.specgram(x=self.y,
  1673. NFFT=self.NFFT_specgram,
  1674. Fs=self.Fs,
  1675. noverlap=self.nover_specgram,
  1676. pad_to=self.pad_to_specgram,
  1677. sides=self.sides,
  1678. mode='psd')
  1679. assert_array_equal(speca, specc)
  1680. assert_array_equal(freqspeca, freqspecc)
  1681. assert_array_equal(ta, tc)
  1682. def test_specgram_complex_mag_equivalent(self):
  1683. specc, freqspecc, tc = mlab.specgram(x=self.y,
  1684. NFFT=self.NFFT_specgram,
  1685. Fs=self.Fs,
  1686. noverlap=self.nover_specgram,
  1687. pad_to=self.pad_to_specgram,
  1688. sides=self.sides,
  1689. mode='complex')
  1690. specm, freqspecm, tm = mlab.specgram(x=self.y,
  1691. NFFT=self.NFFT_specgram,
  1692. Fs=self.Fs,
  1693. noverlap=self.nover_specgram,
  1694. pad_to=self.pad_to_specgram,
  1695. sides=self.sides,
  1696. mode='magnitude')
  1697. assert_array_equal(freqspecc, freqspecm)
  1698. assert_array_equal(tc, tm)
  1699. assert_allclose(np.abs(specc), specm, atol=1e-06)
  1700. def test_specgram_complex_angle_equivalent(self):
  1701. specc, freqspecc, tc = mlab.specgram(x=self.y,
  1702. NFFT=self.NFFT_specgram,
  1703. Fs=self.Fs,
  1704. noverlap=self.nover_specgram,
  1705. pad_to=self.pad_to_specgram,
  1706. sides=self.sides,
  1707. mode='complex')
  1708. speca, freqspeca, ta = mlab.specgram(x=self.y,
  1709. NFFT=self.NFFT_specgram,
  1710. Fs=self.Fs,
  1711. noverlap=self.nover_specgram,
  1712. pad_to=self.pad_to_specgram,
  1713. sides=self.sides,
  1714. mode='angle')
  1715. assert_array_equal(freqspecc, freqspeca)
  1716. assert_array_equal(tc, ta)
  1717. assert_allclose(np.angle(specc), speca, atol=1e-06)
  1718. def test_specgram_complex_phase_equivalent(self):
  1719. specc, freqspecc, tc = mlab.specgram(x=self.y,
  1720. NFFT=self.NFFT_specgram,
  1721. Fs=self.Fs,
  1722. noverlap=self.nover_specgram,
  1723. pad_to=self.pad_to_specgram,
  1724. sides=self.sides,
  1725. mode='complex')
  1726. specp, freqspecp, tp = mlab.specgram(x=self.y,
  1727. NFFT=self.NFFT_specgram,
  1728. Fs=self.Fs,
  1729. noverlap=self.nover_specgram,
  1730. pad_to=self.pad_to_specgram,
  1731. sides=self.sides,
  1732. mode='phase')
  1733. assert_array_equal(freqspecc, freqspecp)
  1734. assert_array_equal(tc, tp)
  1735. assert_allclose(np.unwrap(np.angle(specc), axis=0), specp,
  1736. atol=1e-06)
  1737. def test_specgram_angle_phase_equivalent(self):
  1738. speca, freqspeca, ta = mlab.specgram(x=self.y,
  1739. NFFT=self.NFFT_specgram,
  1740. Fs=self.Fs,
  1741. noverlap=self.nover_specgram,
  1742. pad_to=self.pad_to_specgram,
  1743. sides=self.sides,
  1744. mode='angle')
  1745. specp, freqspecp, tp = mlab.specgram(x=self.y,
  1746. NFFT=self.NFFT_specgram,
  1747. Fs=self.Fs,
  1748. noverlap=self.nover_specgram,
  1749. pad_to=self.pad_to_specgram,
  1750. sides=self.sides,
  1751. mode='phase')
  1752. assert_array_equal(freqspeca, freqspecp)
  1753. assert_array_equal(ta, tp)
  1754. assert_allclose(np.unwrap(speca, axis=0), specp,
  1755. atol=1e-06)
  1756. def test_psd_windowarray_equal(self):
  1757. win = mlab.window_hanning(np.ones(self.NFFT_density_real))
  1758. speca, fspa = mlab.psd(x=self.y,
  1759. NFFT=self.NFFT_density,
  1760. Fs=self.Fs,
  1761. noverlap=self.nover_density,
  1762. pad_to=self.pad_to_density,
  1763. sides=self.sides,
  1764. window=win)
  1765. specb, fspb = mlab.psd(x=self.y,
  1766. NFFT=self.NFFT_density,
  1767. Fs=self.Fs,
  1768. noverlap=self.nover_density,
  1769. pad_to=self.pad_to_density,
  1770. sides=self.sides)
  1771. assert_array_equal(fspa, fspb)
  1772. assert_allclose(speca, specb, atol=1e-08)
  1773. # extra test for cohere...
  1774. def test_cohere():
  1775. N = 1024
  1776. np.random.seed(19680801)
  1777. x = np.random.randn(N)
  1778. # phase offset
  1779. y = np.roll(x, 20)
  1780. # high-freq roll-off
  1781. y = np.convolve(y, np.ones(20) / 20., mode='same')
  1782. cohsq, f = mlab.cohere(x, y, NFFT=256, Fs=2, noverlap=128)
  1783. assert_allclose(np.mean(cohsq), 0.837, atol=1.e-3)
  1784. assert np.isreal(np.mean(cohsq))
  1785. #*****************************************************************
  1786. # These Tests where taken from SCIPY with some minor modifications
  1787. # this can be retrieved from:
  1788. # https://github.com/scipy/scipy/blob/master/scipy/stats/tests/test_kdeoth.py
  1789. #*****************************************************************
  1790. class TestGaussianKDE:
  1791. def test_kde_integer_input(self):
  1792. """Regression test for #1181."""
  1793. x1 = np.arange(5)
  1794. kde = mlab.GaussianKDE(x1)
  1795. y_expected = [0.13480721, 0.18222869, 0.19514935, 0.18222869,
  1796. 0.13480721]
  1797. np.testing.assert_array_almost_equal(kde(x1), y_expected, decimal=6)
  1798. def test_gaussian_kde_covariance_caching(self):
  1799. x1 = np.array([-7, -5, 1, 4, 5], dtype=float)
  1800. xs = np.linspace(-10, 10, num=5)
  1801. # These expected values are from scipy 0.10, before some changes to
  1802. # gaussian_kde. They were not compared with any external reference.
  1803. y_expected = [0.02463386, 0.04689208, 0.05395444, 0.05337754,
  1804. 0.01664475]
  1805. # set it to the default bandwidth.
  1806. kde2 = mlab.GaussianKDE(x1, 'scott')
  1807. y2 = kde2(xs)
  1808. np.testing.assert_array_almost_equal(y_expected, y2, decimal=7)
  1809. def test_kde_bandwidth_method(self):
  1810. np.random.seed(8765678)
  1811. n_basesample = 50
  1812. xn = np.random.randn(n_basesample)
  1813. # Default
  1814. gkde = mlab.GaussianKDE(xn)
  1815. # Supply a callable
  1816. gkde2 = mlab.GaussianKDE(xn, 'scott')
  1817. # Supply a scalar
  1818. gkde3 = mlab.GaussianKDE(xn, bw_method=gkde.factor)
  1819. xs = np.linspace(-7, 7, 51)
  1820. kdepdf = gkde.evaluate(xs)
  1821. kdepdf2 = gkde2.evaluate(xs)
  1822. assert kdepdf.all() == kdepdf2.all()
  1823. kdepdf3 = gkde3.evaluate(xs)
  1824. assert kdepdf.all() == kdepdf3.all()
  1825. class TestGaussianKDECustom:
  1826. def test_no_data(self):
  1827. """Pass no data into the GaussianKDE class."""
  1828. with pytest.raises(ValueError):
  1829. mlab.GaussianKDE([])
  1830. def test_single_dataset_element(self):
  1831. """Pass a single dataset element into the GaussianKDE class."""
  1832. with pytest.raises(ValueError):
  1833. mlab.GaussianKDE([42])
  1834. def test_silverman_multidim_dataset(self):
  1835. """Use a multi-dimensional array as the dataset and test silverman's
  1836. output"""
  1837. x1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
  1838. with pytest.raises(np.linalg.LinAlgError):
  1839. mlab.GaussianKDE(x1, "silverman")
  1840. def test_silverman_singledim_dataset(self):
  1841. """Use a single dimension list as the dataset and test silverman's
  1842. output."""
  1843. x1 = np.array([-7, -5, 1, 4, 5])
  1844. mygauss = mlab.GaussianKDE(x1, "silverman")
  1845. y_expected = 0.76770389927475502
  1846. assert_almost_equal(mygauss.covariance_factor(), y_expected, 7)
  1847. def test_scott_multidim_dataset(self):
  1848. """Use a multi-dimensional array as the dataset and test scott's output
  1849. """
  1850. x1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
  1851. with pytest.raises(np.linalg.LinAlgError):
  1852. mlab.GaussianKDE(x1, "scott")
  1853. def test_scott_singledim_dataset(self):
  1854. """Use a single-dimensional array as the dataset and test scott's
  1855. output"""
  1856. x1 = np.array([-7, -5, 1, 4, 5])
  1857. mygauss = mlab.GaussianKDE(x1, "scott")
  1858. y_expected = 0.72477966367769553
  1859. assert_almost_equal(mygauss.covariance_factor(), y_expected, 7)
  1860. def test_scalar_empty_dataset(self):
  1861. """Use an empty array as the dataset and test the scalar's cov factor
  1862. """
  1863. with pytest.raises(ValueError):
  1864. mlab.GaussianKDE([], bw_method=5)
  1865. def test_scalar_covariance_dataset(self):
  1866. """Use a dataset and test a scalar's cov factor
  1867. """
  1868. np.random.seed(8765678)
  1869. n_basesample = 50
  1870. multidim_data = [np.random.randn(n_basesample) for i in range(5)]
  1871. kde = mlab.GaussianKDE(multidim_data, bw_method=0.5)
  1872. assert kde.covariance_factor() == 0.5
  1873. def test_callable_covariance_dataset(self):
  1874. """Use a multi-dimensional array as the dataset and test the callable's
  1875. cov factor"""
  1876. np.random.seed(8765678)
  1877. n_basesample = 50
  1878. multidim_data = [np.random.randn(n_basesample) for i in range(5)]
  1879. def callable_fun(x):
  1880. return 0.55
  1881. kde = mlab.GaussianKDE(multidim_data, bw_method=callable_fun)
  1882. assert kde.covariance_factor() == 0.55
  1883. def test_callable_singledim_dataset(self):
  1884. """Use a single-dimensional array as the dataset and test the
  1885. callable's cov factor"""
  1886. np.random.seed(8765678)
  1887. n_basesample = 50
  1888. multidim_data = np.random.randn(n_basesample)
  1889. kde = mlab.GaussianKDE(multidim_data, bw_method='silverman')
  1890. y_expected = 0.48438841363348911
  1891. assert_almost_equal(kde.covariance_factor(), y_expected, 7)
  1892. def test_wrong_bw_method(self):
  1893. """Test the error message that should be called when bw is invalid."""
  1894. np.random.seed(8765678)
  1895. n_basesample = 50
  1896. data = np.random.randn(n_basesample)
  1897. with pytest.raises(ValueError):
  1898. mlab.GaussianKDE(data, bw_method="invalid")
  1899. class TestGaussianKDEEvaluate:
  1900. def test_evaluate_diff_dim(self):
  1901. """
  1902. Test the evaluate method when the dim's of dataset and points have
  1903. different dimensions.
  1904. """
  1905. x1 = np.arange(3, 10, 2)
  1906. kde = mlab.GaussianKDE(x1)
  1907. x2 = np.arange(3, 12, 2)
  1908. y_expected = [
  1909. 0.08797252, 0.11774109, 0.11774109, 0.08797252, 0.0370153
  1910. ]
  1911. y = kde.evaluate(x2)
  1912. np.testing.assert_array_almost_equal(y, y_expected, 7)
  1913. def test_evaluate_inv_dim(self):
  1914. """
  1915. Invert the dimensions; i.e., for a dataset of dimension 1 [3, 2, 4],
  1916. the points should have a dimension of 3 [[3], [2], [4]].
  1917. """
  1918. np.random.seed(8765678)
  1919. n_basesample = 50
  1920. multidim_data = np.random.randn(n_basesample)
  1921. kde = mlab.GaussianKDE(multidim_data)
  1922. x2 = [[1], [2], [3]]
  1923. with pytest.raises(ValueError):
  1924. kde.evaluate(x2)
  1925. def test_evaluate_dim_and_num(self):
  1926. """Tests if evaluated against a one by one array"""
  1927. x1 = np.arange(3, 10, 2)
  1928. x2 = np.array([3])
  1929. kde = mlab.GaussianKDE(x1)
  1930. y_expected = [0.08797252]
  1931. y = kde.evaluate(x2)
  1932. np.testing.assert_array_almost_equal(y, y_expected, 7)
  1933. def test_evaluate_point_dim_not_one(self):
  1934. x1 = np.arange(3, 10, 2)
  1935. x2 = [np.arange(3, 10, 2), np.arange(3, 10, 2)]
  1936. kde = mlab.GaussianKDE(x1)
  1937. with pytest.raises(ValueError):
  1938. kde.evaluate(x2)
  1939. def test_evaluate_equal_dim_and_num_lt(self):
  1940. x1 = np.arange(3, 10, 2)
  1941. x2 = np.arange(3, 8, 2)
  1942. kde = mlab.GaussianKDE(x1)
  1943. y_expected = [0.08797252, 0.11774109, 0.11774109]
  1944. y = kde.evaluate(x2)
  1945. np.testing.assert_array_almost_equal(y, y_expected, 7)
  1946. def test_psd_onesided_norm():
  1947. u = np.array([0, 1, 2, 3, 1, 2, 1])
  1948. dt = 1.0
  1949. Su = np.abs(np.fft.fft(u) * dt)**2 / (dt * u.size)
  1950. P, f = mlab.psd(u, NFFT=u.size, Fs=1/dt, window=mlab.window_none,
  1951. detrend=mlab.detrend_none, noverlap=0, pad_to=None,
  1952. scale_by_freq=None,
  1953. sides='onesided')
  1954. Su_1side = np.append([Su[0]], Su[1:4] + Su[4:][::-1])
  1955. assert_allclose(P, Su_1side, atol=1e-06)
  1956. def test_psd_oversampling():
  1957. """Test the case len(x) < NFFT for psd()."""
  1958. u = np.array([0, 1, 2, 3, 1, 2, 1])
  1959. dt = 1.0
  1960. Su = np.abs(np.fft.fft(u) * dt)**2 / (dt * u.size)
  1961. P, f = mlab.psd(u, NFFT=u.size*2, Fs=1/dt, window=mlab.window_none,
  1962. detrend=mlab.detrend_none, noverlap=0, pad_to=None,
  1963. scale_by_freq=None,
  1964. sides='onesided')
  1965. Su_1side = np.append([Su[0]], Su[1:4] + Su[4:][::-1])
  1966. assert_almost_equal(np.sum(P), np.sum(Su_1side)) # same energy