test_windows.py 54 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098
  1. import math
  2. import warnings
  3. import numpy as np
  4. from numpy import array
  5. import pytest
  6. from pytest import raises as assert_raises
  7. from scipy.fft import fft
  8. from scipy.signal import windows, get_window, resample
  9. from scipy.signal.windows._windows import _WIN_FUNC_DATA, _WIN_FUNCS
  10. from scipy._lib._array_api import (
  11. xp_assert_close, xp_assert_equal, array_namespace, is_torch, is_jax, is_cupy,
  12. assert_array_almost_equal, SCIPY_DEVICE, is_numpy, make_xp_test_case,
  13. make_xp_pytest_param, _xp_copy_to_numpy
  14. )
  15. skip_xp_backends = pytest.mark.skip_xp_backends
  16. xfail_xp_backends = pytest.mark.xfail_xp_backends
  17. lazy_xp_modules = [windows]
  18. window_funcs = [
  19. ('boxcar', ()),
  20. ('triang', ()),
  21. ('parzen', ()),
  22. ('bohman', ()),
  23. ('blackman', ()),
  24. ('nuttall', ()),
  25. ('blackmanharris', ()),
  26. ('flattop', ()),
  27. ('bartlett', ()),
  28. ('barthann', ()),
  29. ('hamming', ()),
  30. ('kaiser', (1,)),
  31. ('dpss', (2,)),
  32. ('gaussian', (0.5,)),
  33. ('general_gaussian', (1.5, 2)),
  34. ('chebwin', (1,)),
  35. ('cosine', ()),
  36. ('hann', ()),
  37. ('exponential', ()),
  38. ('taylor', ()),
  39. ('tukey', (0.5,)),
  40. ('lanczos', ()),
  41. ]
  42. @make_xp_test_case(windows.barthann)
  43. class TestBartHann:
  44. def test_basic(self, xp):
  45. xp_assert_close(windows.barthann(6, sym=True, xp=xp),
  46. xp.asarray([0, 0.35857354213752, 0.8794264578624801,
  47. 0.8794264578624801, 0.3585735421375199, 0], dtype=xp.float64),
  48. rtol=1e-15, atol=1e-15)
  49. xp_assert_close(windows.barthann(7, xp=xp),
  50. xp.asarray([0, 0.27, 0.73, 1.0, 0.73, 0.27, 0],
  51. dtype=xp.float64),
  52. rtol=1e-15, atol=1e-15)
  53. xp_assert_close(windows.barthann(6, False, xp=xp),
  54. xp.asarray([0, 0.27, 0.73, 1.0, 0.73, 0.27], dtype=xp.float64),
  55. rtol=1e-15, atol=1e-15)
  56. @make_xp_test_case(windows.bartlett)
  57. class TestBartlett:
  58. def test_basic(self, xp):
  59. xp_assert_close(windows.bartlett(6, xp=xp),
  60. xp.asarray([0, 0.4, 0.8, 0.8, 0.4, 0], dtype=xp.float64))
  61. xp_assert_close(windows.bartlett(7, xp=xp),
  62. xp.asarray([0, 1/3, 2/3, 1.0, 2/3, 1/3, 0], dtype=xp.float64))
  63. xp_assert_close(windows.bartlett(6, False, xp=xp),
  64. xp.asarray([0, 1/3, 2/3, 1.0, 2/3, 1/3], dtype=xp.float64))
  65. @make_xp_test_case(windows.blackman)
  66. class TestBlackman:
  67. def test_basic(self, xp):
  68. xp_assert_close(windows.blackman(6, sym=False, xp=xp),
  69. xp.asarray([0, 0.13, 0.63, 1.0, 0.63, 0.13], dtype=xp.float64),
  70. atol=1e-14)
  71. xp_assert_close(windows.blackman(7, sym=False, xp=xp),
  72. xp.asarray([0, 0.09045342435412804, 0.4591829575459636,
  73. 0.9203636180999081, 0.9203636180999081,
  74. 0.4591829575459636, 0.09045342435412804],
  75. dtype=xp.float64),
  76. atol=1e-8)
  77. xp_assert_close(windows.blackman(6, xp=xp),
  78. xp.asarray([0, 0.2007701432625305, 0.8492298567374694,
  79. 0.8492298567374694, 0.2007701432625305, 0],
  80. dtype=xp.float64),
  81. atol=1e-14)
  82. xp_assert_close(windows.blackman(7, True, xp=xp),
  83. xp.asarray([0, 0.13, 0.63, 1.0, 0.63, 0.13, 0],
  84. dtype=xp.float64), atol=1e-14)
  85. @make_xp_test_case(windows.blackmanharris)
  86. class TestBlackmanHarris:
  87. def test_basic(self, xp):
  88. xp_assert_close(windows.blackmanharris(6, False, xp=xp),
  89. xp.asarray([6.0e-05, 0.055645, 0.520575,
  90. 1.0, 0.520575, 0.055645], dtype=xp.float64))
  91. xp_assert_close(windows.blackmanharris(7, sym=False, xp=xp),
  92. xp.asarray([6.0e-05, 0.03339172347815117, 0.332833504298565,
  93. 0.8893697722232837, 0.8893697722232838,
  94. 0.3328335042985652, 0.03339172347815122],
  95. dtype=xp.float64))
  96. xp_assert_close(windows.blackmanharris(6, xp=xp),
  97. xp.asarray([6.0e-05, 0.1030114893456638, 0.7938335106543362,
  98. 0.7938335106543364, 0.1030114893456638, 6.0e-05],
  99. dtype=xp.float64))
  100. xp_assert_close(windows.blackmanharris(7, sym=True, xp=xp),
  101. xp.asarray([6.0e-05, 0.055645, 0.520575, 1.0, 0.520575,
  102. 0.055645, 6.0e-05], dtype=xp.float64))
  103. @make_xp_test_case(windows.taylor)
  104. class TestTaylor:
  105. def test_normalized(self, xp):
  106. """Tests windows of small length that are normalized to 1. See the
  107. documentation for the Taylor window for more information on
  108. normalization.
  109. """
  110. xp_assert_close(windows.taylor(1, 2, 15, xp=xp),
  111. xp.asarray([1.0], dtype=xp.float64))
  112. xp_assert_close(
  113. windows.taylor(5, 2, 15, xp=xp),
  114. xp.asarray([0.75803341, 0.90757699, 1.0, 0.90757699, 0.75803341],
  115. dtype=xp.float64)
  116. )
  117. xp_assert_close(
  118. windows.taylor(6, 2, 15, xp=xp),
  119. xp.asarray([
  120. 0.7504082, 0.86624416, 0.98208011, 0.98208011, 0.86624416,
  121. 0.7504082
  122. ], dtype=xp.float64)
  123. )
  124. def test_non_normalized(self, xp):
  125. """Test windows of small length that are not normalized to 1. See
  126. the documentation for the Taylor window for more information on
  127. normalization.
  128. """
  129. xp_assert_close(
  130. windows.taylor(5, 2, 15, norm=False, xp=xp),
  131. xp.asarray([
  132. 0.87508054, 1.04771499, 1.15440894, 1.04771499, 0.87508054
  133. ], dtype=xp.float64)
  134. )
  135. xp_assert_close(
  136. windows.taylor(6, 2, 15, norm=False, xp=xp),
  137. xp.asarray([
  138. 0.86627793, 1.0, 1.13372207, 1.13372207, 1.0, 0.86627793
  139. ], dtype=xp.float64)
  140. )
  141. def test_correctness(self, xp):
  142. """This test ensures the correctness of the implemented Taylor
  143. Windowing function. A Taylor Window of 1024 points is created, its FFT
  144. is taken, and the Peak Sidelobe Level (PSLL) and 3dB and 18dB bandwidth
  145. are found and checked.
  146. A publication from Sandia National Laboratories was used as reference
  147. for the correctness values [1]_.
  148. References
  149. -----
  150. .. [1] Armin Doerry, "Catalog of Window Taper Functions for
  151. Sidelobe Control", 2017.
  152. https://www.researchgate.net/profile/Armin_Doerry/publication/316281181_Catalog_of_Window_Taper_Functions_for_Sidelobe_Control/links/58f92cb2a6fdccb121c9d54d/Catalog-of-Window-Taper-Functions-for-Sidelobe-Control.pdf
  153. """
  154. M_win = 1024
  155. N_fft = 131072
  156. # Set norm=False for correctness as the values obtained from the
  157. # scientific publication do not normalize the values. Normalizing
  158. # changes the sidelobe level from the desired value.
  159. w = windows.taylor(M_win, nbar=4, sll=35, norm=False, sym=False, xp=xp)
  160. f_np = fft(_xp_copy_to_numpy(w), N_fft)
  161. spec = 20 * np.log10(np.abs(f_np / np.max(f_np)))
  162. first_zero = np.argmax(np.diff(spec) > 0)
  163. PSLL = np.max(spec[first_zero:-first_zero])
  164. BW_3dB = 2*np.argmax(spec <= -3.0102999566398121) / N_fft * M_win
  165. BW_18dB = 2*np.argmax(spec <= -18.061799739838872) / N_fft * M_win
  166. assert math.isclose(PSLL, -35.1672, abs_tol=1)
  167. assert math.isclose(BW_3dB, 1.1822, abs_tol=0.1)
  168. assert math.isclose(BW_18dB, 2.6112, abs_tol=0.1)
  169. @make_xp_test_case(windows.bohman)
  170. class TestBohman:
  171. def test_basic(self, xp):
  172. xp_assert_close(windows.bohman(6, xp=xp),
  173. xp.asarray([0, 0.1791238937062839, 0.8343114522576858,
  174. 0.8343114522576858, 0.1791238937062838, 0],
  175. dtype=xp.float64))
  176. xp_assert_close(windows.bohman(7, sym=True, xp=xp),
  177. xp.asarray([0, 0.1089977810442293, 0.6089977810442293, 1.0,
  178. 0.6089977810442295, 0.1089977810442293, 0],
  179. dtype=xp.float64))
  180. xp_assert_close(windows.bohman(6, False, xp=xp),
  181. xp.asarray([0, 0.1089977810442293, 0.6089977810442293, 1.0,
  182. 0.6089977810442295, 0.1089977810442293],
  183. dtype=xp.float64))
  184. @make_xp_test_case(windows.boxcar)
  185. class TestBoxcar:
  186. def test_basic(self, xp):
  187. xp_assert_close(windows.boxcar(6, xp=xp),
  188. xp.asarray([1.0, 1, 1, 1, 1, 1], dtype=xp.float64))
  189. xp_assert_close(windows.boxcar(7, xp=xp),
  190. xp.asarray([1.0, 1, 1, 1, 1, 1, 1], dtype=xp.float64))
  191. xp_assert_close(windows.boxcar(6, False, xp=xp),
  192. xp.asarray([1.0, 1, 1, 1, 1, 1], dtype=xp.float64))
  193. cheb_odd_true = [0.200938, 0.107729, 0.134941, 0.165348,
  194. 0.198891, 0.235450, 0.274846, 0.316836,
  195. 0.361119, 0.407338, 0.455079, 0.503883,
  196. 0.553248, 0.602637, 0.651489, 0.699227,
  197. 0.745266, 0.789028, 0.829947, 0.867485,
  198. 0.901138, 0.930448, 0.955010, 0.974482,
  199. 0.988591, 0.997138, 1.000000, 0.997138,
  200. 0.988591, 0.974482, 0.955010, 0.930448,
  201. 0.901138, 0.867485, 0.829947, 0.789028,
  202. 0.745266, 0.699227, 0.651489, 0.602637,
  203. 0.553248, 0.503883, 0.455079, 0.407338,
  204. 0.361119, 0.316836, 0.274846, 0.235450,
  205. 0.198891, 0.165348, 0.134941, 0.107729,
  206. 0.200938]
  207. cheb_even_true = [0.203894, 0.107279, 0.133904,
  208. 0.163608, 0.196338, 0.231986,
  209. 0.270385, 0.311313, 0.354493,
  210. 0.399594, 0.446233, 0.493983,
  211. 0.542378, 0.590916, 0.639071,
  212. 0.686302, 0.732055, 0.775783,
  213. 0.816944, 0.855021, 0.889525,
  214. 0.920006, 0.946060, 0.967339,
  215. 0.983557, 0.994494, 1.000000,
  216. 1.000000, 0.994494, 0.983557,
  217. 0.967339, 0.946060, 0.920006,
  218. 0.889525, 0.855021, 0.816944,
  219. 0.775783, 0.732055, 0.686302,
  220. 0.639071, 0.590916, 0.542378,
  221. 0.493983, 0.446233, 0.399594,
  222. 0.354493, 0.311313, 0.270385,
  223. 0.231986, 0.196338, 0.163608,
  224. 0.133904, 0.107279, 0.203894]
  225. @make_xp_test_case(windows.chebwin)
  226. class TestChebWin:
  227. def test_basic(self, xp):
  228. with warnings.catch_warnings():
  229. warnings.filterwarnings(
  230. "ignore", "This window is not suitable", UserWarning)
  231. xp_assert_close(windows.chebwin(6, 100, xp=xp),
  232. xp.asarray([0.1046401879356917, 0.5075781475823447,
  233. 1.0, 1.0,
  234. 0.5075781475823447, 0.1046401879356917],
  235. dtype=xp.float64),
  236. atol=1e-8
  237. )
  238. xp_assert_close(windows.chebwin(7, 100, xp=xp),
  239. xp.asarray([0.05650405062850233, 0.316608530648474,
  240. 0.7601208123539079, 1.0, 0.7601208123539079,
  241. 0.316608530648474, 0.05650405062850233],
  242. dtype=xp.float64))
  243. xp_assert_close(windows.chebwin(6, 10, xp=xp),
  244. xp.asarray([1.0, 0.6071201674458373, 0.6808391469897297,
  245. 0.6808391469897297, 0.6071201674458373, 1.0],
  246. dtype=xp.float64))
  247. xp_assert_close(windows.chebwin(7, 10, xp=xp),
  248. xp.asarray([1.0, 0.5190521247588651, 0.5864059018130382,
  249. 0.6101519801307441, 0.5864059018130382,
  250. 0.5190521247588651, 1.0], dtype=xp.float64))
  251. xp_assert_close(windows.chebwin(6, 10, False, xp=xp),
  252. xp.asarray([1.0, 0.5190521247588651, 0.5864059018130382,
  253. 0.6101519801307441, 0.5864059018130382,
  254. 0.5190521247588651], dtype=xp.float64))
  255. def test_cheb_odd_high_attenuation(self, xp):
  256. with warnings.catch_warnings():
  257. warnings.filterwarnings(
  258. "ignore", "This window is not suitable", UserWarning)
  259. cheb_odd = windows.chebwin(53, at=-40, xp=xp)
  260. assert_array_almost_equal(cheb_odd, xp.asarray(cheb_odd_true), decimal=4)
  261. def test_cheb_even_high_attenuation(self, xp):
  262. with warnings.catch_warnings():
  263. warnings.filterwarnings(
  264. "ignore", "This window is not suitable", UserWarning)
  265. cheb_even = windows.chebwin(54, at=40, xp=xp)
  266. assert_array_almost_equal(cheb_even, xp.asarray(cheb_even_true), decimal=4)
  267. def test_cheb_odd_low_attenuation(self, xp):
  268. cheb_odd_low_at_true = xp.asarray([1.000000, 0.519052, 0.586405,
  269. 0.610151, 0.586405, 0.519052,
  270. 1.000000], dtype=xp.float64)
  271. with warnings.catch_warnings():
  272. warnings.filterwarnings(
  273. "ignore", "This window is not suitable", UserWarning)
  274. cheb_odd = windows.chebwin(7, at=10, xp=xp)
  275. assert_array_almost_equal(cheb_odd, cheb_odd_low_at_true, decimal=4)
  276. def test_cheb_even_low_attenuation(self, xp):
  277. cheb_even_low_at_true = xp.asarray([1.000000, 0.451924, 0.51027,
  278. 0.541338, 0.541338, 0.51027,
  279. 0.451924, 1.000000], dtype=xp.float64)
  280. with warnings.catch_warnings():
  281. warnings.filterwarnings(
  282. "ignore", "This window is not suitable", UserWarning)
  283. cheb_even = windows.chebwin(8, at=-10, xp=xp)
  284. assert_array_almost_equal(cheb_even, cheb_even_low_at_true, decimal=4)
  285. exponential_data = {
  286. (4, None, 0.2, False):
  287. array([4.53999297624848542e-05,
  288. 6.73794699908546700e-03, 1.00000000000000000e+00,
  289. 6.73794699908546700e-03]),
  290. (4, None, 0.2, True): array([0.00055308437014783, 0.0820849986238988,
  291. 0.0820849986238988, 0.00055308437014783]),
  292. (4, None, 1.0, False): array([0.1353352832366127, 0.36787944117144233, 1.,
  293. 0.36787944117144233]),
  294. (4, None, 1.0, True): array([0.22313016014842982, 0.60653065971263342,
  295. 0.60653065971263342, 0.22313016014842982]),
  296. (4, 2, 0.2, False):
  297. array([4.53999297624848542e-05, 6.73794699908546700e-03,
  298. 1.00000000000000000e+00, 6.73794699908546700e-03]),
  299. (4, 2, 0.2, True): None,
  300. (4, 2, 1.0, False): array([0.1353352832366127, 0.36787944117144233, 1.,
  301. 0.36787944117144233]),
  302. (4, 2, 1.0, True): None,
  303. (5, None, 0.2, True):
  304. array([4.53999297624848542e-05,
  305. 6.73794699908546700e-03, 1.00000000000000000e+00,
  306. 6.73794699908546700e-03, 4.53999297624848542e-05]),
  307. (5, None, 1.0, True): array([0.1353352832366127, 0.36787944117144233, 1.,
  308. 0.36787944117144233, 0.1353352832366127]),
  309. (5, 2, 0.2, True): None,
  310. (5, 2, 1.0, True): None
  311. }
  312. @make_xp_test_case(windows.exponential)
  313. def test_exponential(xp):
  314. for k, v in exponential_data.items():
  315. if v is None:
  316. assert_raises(ValueError, windows.exponential, *k, xp=xp)
  317. else:
  318. win = windows.exponential(*k, xp=xp)
  319. xp_assert_close(win, xp.asarray(v), rtol=1e-14)
  320. @make_xp_test_case(windows.flattop)
  321. class TestFlatTop:
  322. def test_basic(self, xp):
  323. xp_assert_close(windows.flattop(6, sym=False, xp=xp),
  324. xp.asarray([-0.000421051, -0.051263156, 0.19821053, 1.0,
  325. 0.19821053, -0.051263156], dtype=xp.float64))
  326. xp_assert_close(windows.flattop(7, sym=False, xp=xp),
  327. xp.asarray([-0.000421051, -0.03684078115492348,
  328. 0.01070371671615342, 0.7808739149387698,
  329. 0.7808739149387698, 0.01070371671615342,
  330. -0.03684078115492348], dtype=xp.float64))
  331. xp_assert_close(windows.flattop(6, xp=xp),
  332. xp.asarray([-0.000421051, -0.0677142520762119,
  333. 0.6068721525762117, 0.6068721525762117,
  334. -0.0677142520762119, -0.000421051],
  335. dtype=xp.float64))
  336. xp_assert_close(windows.flattop(7, True, xp=xp),
  337. xp.asarray([-0.000421051, -0.051263156, 0.19821053, 1.0,
  338. 0.19821053, -0.051263156, -0.000421051],
  339. dtype=xp.float64))
  340. @make_xp_test_case(windows.gaussian)
  341. class TestGaussian:
  342. def test_basic(self, xp):
  343. xp_assert_close(windows.gaussian(6, 1.0, xp=xp),
  344. xp.asarray([0.04393693362340742, 0.3246524673583497,
  345. 0.8824969025845955, 0.8824969025845955,
  346. 0.3246524673583497, 0.04393693362340742],
  347. dtype=xp.float64))
  348. xp_assert_close(windows.gaussian(7, 1.2, xp=xp),
  349. xp.asarray([0.04393693362340742, 0.2493522087772962,
  350. 0.7066482778577162, 1.0, 0.7066482778577162,
  351. 0.2493522087772962, 0.04393693362340742],
  352. dtype=xp.float64))
  353. xp_assert_close(windows.gaussian(7, 3, xp=xp),
  354. xp.asarray([0.6065306597126334, 0.8007374029168081,
  355. 0.9459594689067654, 1.0, 0.9459594689067654,
  356. 0.8007374029168081, 0.6065306597126334],
  357. dtype=xp.float64))
  358. xp_assert_close(windows.gaussian(6, 3, False, xp=xp),
  359. xp.asarray([0.6065306597126334, 0.8007374029168081,
  360. 0.9459594689067654, 1.0, 0.9459594689067654,
  361. 0.8007374029168081], dtype=xp.float64))
  362. @make_xp_test_case(windows.general_cosine)
  363. class TestGeneralCosine:
  364. def test_basic(self, xp):
  365. a = xp.asarray([0.5, 0.3, 0.2])
  366. xp_assert_close(windows.general_cosine(5, a),
  367. xp.asarray([0.4, 0.3, 1, 0.3, 0.4], dtype=xp.float64))
  368. a = xp.asarray([0.5, 0.3, 0.2])
  369. xp_assert_close(windows.general_cosine(4, a, sym=False),
  370. xp.asarray([0.4, 0.3, 1, 0.3], dtype=xp.float64))
  371. @make_xp_test_case(windows.general_hamming)
  372. class TestGeneralHamming:
  373. def test_basic(self, xp):
  374. xp_assert_close(windows.general_hamming(5, 0.7, xp=xp),
  375. xp.asarray([0.4, 0.7, 1.0, 0.7, 0.4], dtype=xp.float64))
  376. xp_assert_close(windows.general_hamming(5, 0.75, sym=False, xp=xp),
  377. xp.asarray([0.5, 0.6727457514, 0.9522542486,
  378. 0.9522542486, 0.6727457514], dtype=xp.float64))
  379. xp_assert_close(windows.general_hamming(6, 0.75, sym=True, xp=xp),
  380. xp.asarray([0.5, 0.6727457514, 0.9522542486,
  381. 0.9522542486, 0.6727457514, 0.5], dtype=xp.float64))
  382. @make_xp_test_case(windows.hamming)
  383. class TestHamming:
  384. def test_basic(self, xp):
  385. xp_assert_close(windows.hamming(6, False, xp=xp),
  386. xp.asarray([0.08, 0.31, 0.77, 1.0, 0.77, 0.31],
  387. dtype=xp.float64))
  388. xp_assert_close(windows.hamming(7, sym=False, xp=xp),
  389. xp.asarray([0.08, 0.2531946911449826, 0.6423596296199047,
  390. 0.9544456792351128, 0.9544456792351128,
  391. 0.6423596296199047, 0.2531946911449826],
  392. dtype=xp.float64))
  393. xp_assert_close(windows.hamming(6, xp=xp),
  394. xp.asarray([0.08, 0.3978521825875242, 0.9121478174124757,
  395. 0.9121478174124757, 0.3978521825875242, 0.08],
  396. dtype=xp.float64))
  397. xp_assert_close(windows.hamming(7, sym=True, xp=xp),
  398. xp.asarray([0.08, 0.31, 0.77, 1.0, 0.77, 0.31, 0.08],
  399. dtype=xp.float64))
  400. @make_xp_test_case(windows.hann)
  401. class TestHann:
  402. def test_basic(self, xp):
  403. xp_assert_close(windows.hann(6, sym=False, xp=xp),
  404. xp.asarray([0, 0.25, 0.75, 1.0, 0.75, 0.25], dtype=xp.float64),
  405. rtol=1e-15, atol=1e-15)
  406. xp_assert_close(windows.hann(7, sym=False, xp=xp),
  407. xp.asarray([0, 0.1882550990706332, 0.6112604669781572,
  408. 0.9504844339512095, 0.9504844339512095,
  409. 0.6112604669781572, 0.1882550990706332],
  410. dtype=xp.float64),
  411. rtol=1e-15, atol=1e-15)
  412. xp_assert_close(windows.hann(6, True, xp=xp),
  413. xp.asarray([0, 0.3454915028125263, 0.9045084971874737,
  414. 0.9045084971874737, 0.3454915028125263, 0],
  415. dtype=xp.float64),
  416. rtol=1e-15, atol=1e-15)
  417. xp_assert_close(windows.hann(7, xp=xp),
  418. xp.asarray([0, 0.25, 0.75, 1.0, 0.75, 0.25, 0],
  419. dtype=xp.float64),
  420. rtol=1e-15, atol=1e-15)
  421. @make_xp_test_case(windows.kaiser)
  422. class TestKaiser:
  423. def test_basic(self, xp):
  424. xp_assert_close(windows.kaiser(6, 0.5, xp=xp),
  425. xp.asarray([0.9403061933191572, 0.9782962393705389,
  426. 0.9975765035372042, 0.9975765035372042,
  427. 0.9782962393705389, 0.9403061933191572],
  428. dtype=xp.float64))
  429. xp_assert_close(windows.kaiser(7, 0.5, xp=xp),
  430. xp.asarray([0.9403061933191572, 0.9732402256999829,
  431. 0.9932754654413773, 1.0, 0.9932754654413773,
  432. 0.9732402256999829, 0.9403061933191572],
  433. dtype=xp.float64))
  434. xp_assert_close(windows.kaiser(6, 2.7, xp=xp),
  435. xp.asarray([0.2603047507678832, 0.6648106293528054,
  436. 0.9582099802511439, 0.9582099802511439,
  437. 0.6648106293528054, 0.2603047507678832],
  438. dtype=xp.float64))
  439. xp_assert_close(windows.kaiser(7, 2.7, xp=xp),
  440. xp.asarray([0.2603047507678832, 0.5985765418119844,
  441. 0.8868495172060835, 1.0, 0.8868495172060835,
  442. 0.5985765418119844, 0.2603047507678832],
  443. dtype=xp.float64))
  444. xp_assert_close(windows.kaiser(6, 2.7, False, xp=xp),
  445. xp.asarray([0.2603047507678832, 0.5985765418119844,
  446. 0.8868495172060835, 1.0, 0.8868495172060835,
  447. 0.5985765418119844], dtype=xp.float64))
  448. @make_xp_test_case(windows.kaiser_bessel_derived)
  449. class TestKaiserBesselDerived:
  450. def test_basic(self, xp):
  451. # cover case `M < 1`
  452. w = windows.kaiser_bessel_derived(0.5, beta=4.0, xp=xp)
  453. xp_assert_equal(w, xp.asarray([]))
  454. M = 100
  455. w = windows.kaiser_bessel_derived(M, beta=4.0, xp=xp)
  456. w2 = windows.get_window(('kaiser bessel derived', 4.0),
  457. M, fftbins=False, xp=xp)
  458. xp_assert_close(w, w2)
  459. # Test for Princen-Bradley condition
  460. actual = w[:M // 2] ** 2 + w[-M // 2:] ** 2
  461. xp_assert_close(actual, xp.ones(actual.shape, dtype=actual.dtype))
  462. # Test actual values from other implementations
  463. # M = 2: sqrt(2) / 2
  464. # M = 4: 0.518562710536, 0.855039598640
  465. # M = 6: 0.436168993154, 0.707106781187, 0.899864772847
  466. # Ref:https://github.com/scipy/scipy/pull/4747#issuecomment-172849418
  467. actual = windows.kaiser_bessel_derived(2, beta=np.pi / 2, xp=xp)[:1]
  468. desired = xp.ones_like(actual) * math.sqrt(2) / 2.0
  469. xp_assert_close(actual, desired)
  470. xp_assert_close(windows.kaiser_bessel_derived(4, beta=np.pi / 2, xp=xp)[:2],
  471. xp.asarray([0.518562710536, 0.855039598640], dtype=xp.float64))
  472. xp_assert_close(windows.kaiser_bessel_derived(6, beta=np.pi / 2, xp=xp)[:3],
  473. xp.asarray([0.436168993154, 0.707106781187, 0.899864772847],
  474. dtype=xp.float64))
  475. def test_exceptions(self, xp):
  476. M = 100
  477. # Assert ValueError for odd window length
  478. msg = ("Kaiser-Bessel Derived windows are only defined for even "
  479. "number of points")
  480. with assert_raises(ValueError, match=msg):
  481. windows.kaiser_bessel_derived(M + 1, beta=4., xp=xp)
  482. # Assert ValueError for non-symmetric setting
  483. msg = ("Kaiser-Bessel Derived windows are only defined for "
  484. "symmetric shapes")
  485. with assert_raises(ValueError, match=msg):
  486. windows.kaiser_bessel_derived(M + 1, beta=4., sym=False, xp=xp)
  487. class TestNuttall:
  488. def test_basic(self, xp):
  489. xp_assert_close(windows.nuttall(6, sym=False, xp=xp),
  490. xp.asarray([0.0003628, 0.0613345, 0.5292298, 1.0, 0.5292298,
  491. 0.0613345], dtype=xp.float64))
  492. xp_assert_close(windows.nuttall(7, sym=False, xp=xp),
  493. xp.asarray([0.0003628, 0.03777576895352025,
  494. 0.3427276199688195,
  495. 0.8918518610776603, 0.8918518610776603,
  496. 0.3427276199688196, 0.0377757689535203],
  497. dtype=xp.float64))
  498. xp_assert_close(windows.nuttall(6, xp=xp),
  499. xp.asarray([0.0003628, 0.1105152530498718,
  500. 0.7982580969501282, 0.7982580969501283,
  501. 0.1105152530498719, 0.0003628], dtype=xp.float64))
  502. xp_assert_close(windows.nuttall(7, True, xp=xp),
  503. xp.asarray([0.0003628, 0.0613345, 0.5292298, 1.0,
  504. 0.5292298, 0.0613345, 0.0003628], dtype=xp.float64))
  505. @make_xp_test_case(windows.parzen)
  506. class TestParzen:
  507. def test_basic(self, xp):
  508. xp_assert_close(windows.parzen(6, xp=xp),
  509. xp.asarray([0.009259259259259254, 0.25, 0.8611111111111112,
  510. 0.8611111111111112, 0.25, 0.009259259259259254],
  511. dtype=xp.float64))
  512. xp_assert_close(windows.parzen(7, sym=True, xp=xp),
  513. xp.asarray([0.00583090379008747, 0.1574344023323616,
  514. 0.6501457725947521, 1.0, 0.6501457725947521,
  515. 0.1574344023323616, 0.00583090379008747],
  516. dtype=xp.float64))
  517. xp_assert_close(windows.parzen(6, False, xp=xp),
  518. xp.asarray([0.00583090379008747, 0.1574344023323616,
  519. 0.6501457725947521, 1.0, 0.6501457725947521,
  520. 0.1574344023323616], dtype=xp.float64))
  521. @make_xp_test_case(windows.triang)
  522. class TestTriang:
  523. def test_basic(self, xp):
  524. xp_assert_close(windows.triang(6, True, xp=xp),
  525. xp.asarray([1/6, 1/2, 5/6, 5/6, 1/2, 1/6], dtype=xp.float64))
  526. xp_assert_close(windows.triang(7, xp=xp),
  527. xp.asarray([1/4, 1/2, 3/4, 1, 3/4, 1/2, 1/4], dtype=xp.float64))
  528. xp_assert_close(windows.triang(6, sym=False, xp=xp),
  529. xp.asarray([1/4, 1/2, 3/4, 1, 3/4, 1/2], dtype=xp.float64))
  530. tukey_data = {
  531. (4, 0.5, True): array([0.0, 1.0, 1.0, 0.0]),
  532. (4, 0.9, True): array([0.0, 0.84312081893436686,
  533. 0.84312081893436686, 0.0]),
  534. (4, 1.0, True): array([0.0, 0.75, 0.75, 0.0]),
  535. (4, 0.5, False): array([0.0, 1.0, 1.0, 1.0]),
  536. (4, 0.9, False): array([0.0, 0.58682408883346526,
  537. 1.0, 0.58682408883346526]),
  538. (4, 1.0, False): array([0.0, 0.5, 1.0, 0.5]),
  539. (5, 0.0, True): array([1.0, 1.0, 1.0, 1.0, 1.0]),
  540. (5, 0.8, True): array([0.0, 0.69134171618254492,
  541. 1.0, 0.69134171618254492, 0.0]),
  542. (5, 1.0, True): array([0.0, 0.5, 1.0, 0.5, 0.0]),
  543. (6, 0): [1.0, 1, 1, 1, 1, 1],
  544. (7, 0): [1.0, 1, 1, 1, 1, 1, 1],
  545. (6, .25): [0.0, 1, 1, 1, 1, 0],
  546. (7, .25): [0.0, 1, 1, 1, 1, 1, 0],
  547. (6,): [0, 0.9045084971874737, 1.0, 1.0, 0.9045084971874735, 0],
  548. (7,): [0, 0.75, 1.0, 1.0, 1.0, 0.75, 0],
  549. (6, .75): [0, 0.5522642316338269, 1.0, 1.0, 0.5522642316338267, 0],
  550. (7, .75): [0, 0.4131759111665348, 0.9698463103929542, 1.0,
  551. 0.9698463103929542, 0.4131759111665347, 0],
  552. (6, 1): [0, 0.3454915028125263, 0.9045084971874737, 0.9045084971874737,
  553. 0.3454915028125263, 0],
  554. (7, 1): [0, 0.25, 0.75, 1.0, 0.75, 0.25, 0],
  555. }
  556. @make_xp_test_case(windows.tukey)
  557. class TestTukey:
  558. def test_basic(self, xp):
  559. # Test against hardcoded data
  560. for k, v in tukey_data.items():
  561. if v is None:
  562. assert_raises(ValueError, windows.tukey, *k, xp=xp)
  563. else:
  564. if is_torch(xp) and k in [(6,), (6, .75), (7, .75), (6,1)]:
  565. atol_rtol = {'rtol': 3e-8, 'atol': 1e-8}
  566. else:
  567. atol_rtol = {'rtol': 1e-15, 'atol': 1e-15 }
  568. win = windows.tukey(*k, xp=xp)
  569. xp_assert_close(win, xp.asarray(v),
  570. check_dtype=False, **atol_rtol)
  571. def test_extremes(self, xp):
  572. # Test extremes of alpha correspond to boxcar and hann
  573. tuk0 = windows.tukey(100, 0, xp=xp)
  574. box0 = windows.boxcar(100, xp=xp)
  575. xp_assert_close(tuk0, box0)
  576. tuk1 = windows.tukey(100, 1, xp=xp)
  577. han1 = windows.hann(100, xp=xp)
  578. xp_assert_close(tuk1, han1)
  579. dpss_data = {
  580. # All values from MATLAB:
  581. # * taper[1] of (3, 1.4, 3) sign-flipped
  582. # * taper[3] of (5, 1.5, 5) sign-flipped
  583. (4, 0.1, 2): ([[0.497943898, 0.502047681, 0.502047681, 0.497943898], [0.670487993, 0.224601537, -0.224601537, -0.670487993]], [0.197961815, 0.002035474]), # noqa: E501
  584. (3, 1.4, 3): ([[0.410233151, 0.814504464, 0.410233151], [0.707106781, 0.0, -0.707106781], [0.575941629, -0.580157287, 0.575941629]], [0.999998093, 0.998067480, 0.801934426]), # noqa: E501
  585. (5, 1.5, 5): ([[0.1745071052, 0.4956749177, 0.669109327, 0.495674917, 0.174507105], [0.4399493348, 0.553574369, 0.0, -0.553574369, -0.439949334], [0.631452756, 0.073280238, -0.437943884, 0.073280238, 0.631452756], [0.553574369, -0.439949334, 0.0, 0.439949334, -0.553574369], [0.266110290, -0.498935248, 0.600414741, -0.498935248, 0.266110290147157]], [0.999728571, 0.983706916, 0.768457889, 0.234159338, 0.013947282907567]), # noqa: E501
  586. (100, 2, 4): ([[0.0030914414, 0.0041266922, 0.005315076, 0.006665149, 0.008184854, 0.0098814158, 0.011761239, 0.013829809, 0.016091597, 0.018549973, 0.02120712, 0.02406396, 0.027120092, 0.030373728, 0.033821651, 0.037459181, 0.041280145, 0.045276872, 0.049440192, 0.053759447, 0.058222524, 0.062815894, 0.067524661, 0.072332638, 0.077222418, 0.082175473, 0.087172252, 0.092192299, 0.097214376, 0.1022166, 0.10717657, 0.11207154, 0.11687856, 0.12157463, 0.12613686, 0.13054266, 0.13476986, 0.13879691, 0.14260302, 0.14616832, 0.14947401, 0.1525025, 0.15523755, 0.15766438, 0.15976981, 0.16154233, 0.16297223, 0.16405162, 0.16477455, 0.16513702, 0.16513702, 0.16477455, 0.16405162, 0.16297223, 0.16154233, 0.15976981, 0.15766438, 0.15523755, 0.1525025, 0.14947401, 0.14616832, 0.14260302, 0.13879691, 0.13476986, 0.13054266, 0.12613686, 0.12157463, 0.11687856, 0.11207154, 0.10717657, 0.1022166, 0.097214376, 0.092192299, 0.087172252, 0.082175473, 0.077222418, 0.072332638, 0.067524661, 0.062815894, 0.058222524, 0.053759447, 0.049440192, 0.045276872, 0.041280145, 0.037459181, 0.033821651, 0.030373728, 0.027120092, 0.02406396, 0.02120712, 0.018549973, 0.016091597, 0.013829809, 0.011761239, 0.0098814158, 0.008184854, 0.006665149, 0.005315076, 0.0041266922, 0.0030914414], [0.018064449, 0.022040342, 0.026325013, 0.030905288, 0.035764398, 0.040881982, 0.046234148, 0.051793558, 0.057529559, 0.063408356, 0.069393216, 0.075444716, 0.081521022, 0.087578202, 0.093570567, 0.099451049, 0.10517159, 0.11068356, 0.11593818, 0.12088699, 0.12548227, 0.12967752, 0.1334279, 0.13669069, 0.13942569, 0.1415957, 0.14316686, 0.14410905, 0.14439626, 0.14400686, 0.14292389, 0.1411353, 0.13863416, 0.13541876, 0.13149274, 0.12686516, 0.12155045, 0.1155684, 0.10894403, 0.10170748, 0.093893752, 0.08554251, 0.076697768, 0.067407559, 0.057723559, 0.04770068, 0.037396627, 0.026871428, 0.016186944, 0.0054063557, -0.0054063557, -0.016186944, -0.026871428, -0.037396627, -0.04770068, -0.057723559, -0.067407559, -0.076697768, -0.08554251, -0.093893752, -0.10170748, -0.10894403, -0.1155684, -0.12155045, -0.12686516, -0.13149274, -0.13541876, -0.13863416, -0.1411353, -0.14292389, -0.14400686, -0.14439626, -0.14410905, -0.14316686, -0.1415957, -0.13942569, -0.13669069, -0.1334279, -0.12967752, -0.12548227, -0.12088699, -0.11593818, -0.11068356, -0.10517159, -0.099451049, -0.093570567, -0.087578202, -0.081521022, -0.075444716, -0.069393216, -0.063408356, -0.057529559, -0.051793558, -0.046234148, -0.040881982, -0.035764398, -0.030905288, -0.026325013, -0.022040342, -0.018064449], [0.064817553, 0.072567801, 0.080292992, 0.087918235, 0.095367076, 0.10256232, 0.10942687, 0.1158846, 0.12186124, 0.12728523, 0.13208858, 0.13620771, 0.13958427, 0.14216587, 0.14390678, 0.14476863, 0.1447209, 0.14374148, 0.14181704, 0.13894336, 0.13512554, 0.13037812, 0.1247251, 0.11819984, 0.11084487, 0.10271159, 0.093859853, 0.084357497, 0.074279719, 0.063708406, 0.052731374, 0.041441525, 0.029935953, 0.018314987, 0.0066811877, -0.0048616765, -0.016209689, -0.027259848, -0.037911124, -0.048065512, -0.05762905, -0.066512804, -0.0746338, -0.081915903, -0.088290621, -0.09369783, -0.098086416, -0.10141482, -0.10365146, -0.10477512, -0.10477512, -0.10365146, -0.10141482, -0.098086416, -0.09369783, -0.088290621, -0.081915903, -0.0746338, -0.066512804, -0.05762905, -0.048065512, -0.037911124, -0.027259848, -0.016209689, -0.0048616765, 0.0066811877, 0.018314987, 0.029935953, 0.041441525, 0.052731374, 0.063708406, 0.074279719, 0.084357497, 0.093859853, 0.10271159, 0.11084487, 0.11819984, 0.1247251, 0.13037812, 0.13512554, 0.13894336, 0.14181704, 0.14374148, 0.1447209, 0.14476863, 0.14390678, 0.14216587, 0.13958427, 0.13620771, 0.13208858, 0.12728523, 0.12186124, 0.1158846, 0.10942687, 0.10256232, 0.095367076, 0.087918235, 0.080292992, 0.072567801, 0.064817553], [0.14985551, 0.15512305, 0.15931467, 0.16236806, 0.16423291, 0.16487165, 0.16426009, 0.1623879, 0.1592589, 0.15489114, 0.14931693, 0.14258255, 0.13474785, 0.1258857, 0.11608124, 0.10543095, 0.094041635, 0.082029213, 0.069517411, 0.056636348, 0.043521028, 0.030309756, 0.017142511, 0.0041592774, -0.0085016282, -0.020705223, -0.032321494, -0.043226982, -0.053306291, -0.062453515, -0.070573544, -0.077583253, -0.083412547, -0.088005244, -0.091319802, -0.093329861, -0.094024602, -0.093408915, -0.091503383, -0.08834406, -0.08398207, -0.078483012, -0.071926192, -0.064403681, -0.056019215, -0.046886954, -0.037130106, -0.026879442, -0.016271713, -0.005448, 0.005448, 0.016271713, 0.026879442, 0.037130106, 0.046886954, 0.056019215, 0.064403681, 0.071926192, 0.078483012, 0.08398207, 0.08834406, 0.091503383, 0.093408915, 0.094024602, 0.093329861, 0.091319802, 0.088005244, 0.083412547, 0.077583253, 0.070573544, 0.062453515, 0.053306291, 0.043226982, 0.032321494, 0.020705223, 0.0085016282, -0.0041592774, -0.017142511, -0.030309756, -0.043521028, -0.056636348, -0.069517411, -0.082029213, -0.094041635, -0.10543095, -0.11608124, -0.1258857, -0.13474785, -0.14258255, -0.14931693, -0.15489114, -0.1592589, -0.1623879, -0.16426009, -0.16487165, -0.16423291, -0.16236806, -0.15931467, -0.15512305, -0.14985551]], [0.999943140, 0.997571533, 0.959465463, 0.721862496]), # noqa: E501
  587. }
  588. @make_xp_test_case(windows.dpss)
  589. class TestDPSS:
  590. def test_basic(self, xp):
  591. # Test against hardcoded data
  592. for k, v in dpss_data.items():
  593. win, ratios = windows.dpss(*k, return_ratios=True, xp=xp)
  594. xp_assert_close(win, v[0], atol=1e-7, err_msg=k)
  595. xp_assert_close(ratios, v[1], rtol=1e-5, atol=1e-7, err_msg=k)
  596. def test_unity(self, xp):
  597. # Test unity value handling (gh-2221)
  598. for M in range(1, 21):
  599. # corrected w/approximation (default)
  600. win = windows.dpss(M, M / 2.1, xp=xp)
  601. expected = M % 2 # one for odd, none for even
  602. xp_assert_equal(np.isclose(win, 1.).sum(), expected,
  603. err_msg=f'{win}')
  604. # corrected w/subsample delay (slower)
  605. win_sub = windows.dpss(M, M / 2.1, norm='subsample', xp=xp)
  606. if M > 2:
  607. # @M=2 the subsample doesn't do anything
  608. xp_assert_equal(np.isclose(win_sub, 1.).sum(), expected,
  609. err_msg=f'{win_sub}')
  610. xp_assert_close(win, win_sub, rtol=0.03) # within 3%
  611. # not the same, l2-norm
  612. win_2 = windows.dpss(M, M / 2.1, norm=2, xp=xp)
  613. expected = 1 if M == 1 else 0
  614. xp_assert_equal(np.isclose(win_2, 1.).sum(), expected,
  615. err_msg=f'{win_2}')
  616. def test_extremes(self, xp):
  617. # Test extremes of alpha
  618. lam = windows.dpss(31, 6, 4, return_ratios=True, xp=xp)[1]
  619. xp_assert_close(lam, xp.ones_like(lam))
  620. lam = windows.dpss(31, 7, 4, return_ratios=True, xp=xp)[1]
  621. xp_assert_close(lam, xp.ones_like(lam))
  622. lam = windows.dpss(31, 8, 4, return_ratios=True, xp=xp)[1]
  623. xp_assert_close(lam, xp.ones_like(lam))
  624. def test_degenerate(self, xp):
  625. # Test failures
  626. assert_raises(ValueError, windows.dpss, 4, 1.5, -1) # Bad Kmax
  627. assert_raises(ValueError, windows.dpss, 4, 1.5, -5)
  628. assert_raises(TypeError, windows.dpss, 4, 1.5, 1.1)
  629. assert_raises(ValueError, windows.dpss, 3, 1.5, 3) # NW must be < N/2.
  630. assert_raises(ValueError, windows.dpss, 3, -1, 3) # NW must be pos
  631. assert_raises(ValueError, windows.dpss, 3, 0, 3)
  632. assert_raises(ValueError, windows.dpss, -1, 1, 3) # negative M
  633. @skip_xp_backends(np_only=True)
  634. def test_degenerate_signle_samples(self, xp):
  635. # Single samples
  636. w = windows.dpss(1, 1.)
  637. xp_assert_equal(w, [1.])
  638. w, ratio = windows.dpss(1, 1., return_ratios=True)
  639. xp_assert_equal(w, [1.])
  640. assert ratio == 1.
  641. w, ratio = windows.dpss(1, 1., Kmax=4, return_ratios=True)
  642. xp_assert_equal(w, [1.])
  643. assert isinstance(ratio, np.ndarray)
  644. xp_assert_equal(ratio, [1.])
  645. assert_raises(ValueError, windows.dpss, 4, 1.5, -1, xp=xp) # Bad Kmax
  646. assert_raises(ValueError, windows.dpss, 4, 1.5, -5, xp=xp)
  647. assert_raises(TypeError, windows.dpss, 4, 1.5, 1.1, xp=xp)
  648. assert_raises(ValueError, windows.dpss, 3, 1.5, 3, xp=xp) # NW must be < N/2.
  649. assert_raises(ValueError, windows.dpss, 3, -1, 3, xp=xp) # NW must be pos
  650. assert_raises(ValueError, windows.dpss, 3, 0, 3, xp=xp)
  651. assert_raises(ValueError, windows.dpss, -1, 1, 3, xp=xp) # negative M
  652. @make_xp_test_case(windows.lanczos)
  653. class TestLanczos:
  654. def test_basic(self, xp):
  655. # Analytical results:
  656. # sinc(x) = sinc(-x)
  657. # sinc(pi) = 0, sinc(0) = 1
  658. # Hand computation on WolframAlpha:
  659. # sinc(2 pi / 3) = 0.413496672
  660. # sinc(pi / 3) = 0.826993343
  661. # sinc(3 pi / 5) = 0.504551152
  662. # sinc(pi / 5) = 0.935489284
  663. xp_assert_close(windows.lanczos(6, sym=False, xp=xp),
  664. xp.asarray([0., 0.413496672,
  665. 0.826993343, 1., 0.826993343,
  666. 0.413496672], dtype=xp.float64),
  667. atol=1e-9)
  668. xp_assert_close(windows.lanczos(6, xp=xp),
  669. xp.asarray([0., 0.504551152,
  670. 0.935489284, 0.935489284,
  671. 0.504551152, 0.], dtype=xp.float64),
  672. atol=1e-9)
  673. xp_assert_close(windows.lanczos(7, sym=True, xp=xp),
  674. xp.asarray([0., 0.413496672,
  675. 0.826993343, 1., 0.826993343,
  676. 0.413496672, 0.], dtype=xp.float64),
  677. atol=1e-9)
  678. def test_array_size(self, xp):
  679. for n in [0, 10, 11]:
  680. assert windows.lanczos(n, sym=False, xp=xp).shape[0] == n
  681. assert windows.lanczos(n, sym=True, xp=xp).shape[0] == n
  682. @make_xp_test_case(windows.get_window)
  683. class TestGetWindow:
  684. """Unit test for `scipy.signal.get_windows`. """
  685. def test_WIN_FUNC_DATA_integrity(self):
  686. """Verify that the `_windows._WIN_FUNC_DATA` dict is consistent.
  687. The keys of _WIN_FUNC_DATA are made of tuples of strings of allowed window
  688. names. Its values are 2-tuples made up of the window function and a
  689. entry characterizing the existence of window parameters as ``True``,
  690. ``False`` or ``'OPTIONAL'``.
  691. It is verified that the correct window name (i.e., corresponding to the
  692. function in the value tuple) is included in the key tuple. It is also checked
  693. that the second entry in the value tuple is either ``True``, ``False`` or
  694. ``'OPTIONAL'``.
  695. """
  696. for nn_, v_ in _WIN_FUNC_DATA.items():
  697. func_name = v_[0].__name__
  698. msg = f"Function name in {nn_} does not contain name of actual function!"
  699. assert func_name in nn_, msg
  700. assert v_[1] in (True, False, 'OPTIONAL')
  701. @make_xp_test_case(windows.boxcar)
  702. def test_boxcar(self, xp):
  703. w = windows.get_window('boxcar', 12, xp=xp)
  704. xp_assert_equal(w, xp.ones_like(w))
  705. # window is a tuple of len 1
  706. w = windows.get_window(('boxcar',), 16, xp=xp)
  707. xp_assert_equal(w, xp.ones_like(w))
  708. @make_xp_test_case(windows.chebwin)
  709. def test_cheb_odd(self, xp):
  710. with warnings.catch_warnings():
  711. warnings.filterwarnings(
  712. "ignore", "This window is not suitable", UserWarning)
  713. w = windows.get_window(('chebwin', -40), 53, fftbins=False, xp=xp)
  714. assert_array_almost_equal(
  715. w, xp.asarray(cheb_odd_true, dtype=xp.float64), decimal=4
  716. )
  717. @make_xp_test_case(windows.chebwin)
  718. def test_cheb_even(self, xp):
  719. with warnings.catch_warnings():
  720. warnings.filterwarnings(
  721. "ignore", "This window is not suitable", UserWarning)
  722. w = windows.get_window(('chebwin', 40), 54, fftbins=False, xp=xp)
  723. assert_array_almost_equal(w, xp.asarray(cheb_even_true), decimal=4)
  724. @make_xp_test_case(windows.dpss)
  725. def test_dpss(self, xp):
  726. win1 = windows.get_window(('dpss', 3), 64, fftbins=False, xp=xp)
  727. win2 = windows.dpss(64, 3, xp=xp)
  728. xp_assert_equal(win1, win2)
  729. @make_xp_test_case(windows.kaiser)
  730. def test_kaiser_float(self, xp):
  731. win1 = windows.get_window(7.2, 64, xp=xp)
  732. win2 = windows.kaiser(64, 7.2, False, xp=xp)
  733. if is_jax(xp):
  734. # On JAX with jit enabled, there is a very small discrepancy
  735. # in the results.
  736. xp_assert_close(win1, win2, rtol=xp.finfo(win1.dtype).eps)
  737. else:
  738. xp_assert_equal(win1, win2)
  739. @pytest.mark.parametrize('Nx', [-1, 3.0, np.float64(3)])
  740. @make_xp_test_case(windows.hann)
  741. def test_invalid_parameter_NX(self, Nx, xp):
  742. with pytest.raises(ValueError, match="^Parameter Nx=.*"):
  743. windows.get_window('hann', Nx, xp=xp)
  744. # noinspection PyTypeChecker
  745. def test_invalid_inputs(self, xp):
  746. """Raise all exceptions (except those concerning parameter `Nx`). """
  747. with pytest.raises(ValueError, match="^Parameter fftbins=.*"):
  748. windows.get_window('hann', 5, fftbins=1, xp=xp)
  749. with pytest.raises(ValueError, match="^Parameter window=.*"):
  750. windows.get_window(['hann',], 5, xp=xp)
  751. with pytest.raises(ValueError, match="^First tuple entry of parameter win.*"):
  752. windows.get_window((42,), 5, xp=xp)
  753. with pytest.raises(ValueError, match="^Invalid window name 'INVALID'.*"):
  754. windows.get_window('INVALID', 5, xp=xp)
  755. with pytest.raises(ValueError, match="^'hann' does not allow parameters.*"):
  756. windows.get_window(('hann', 1), 5, xp=xp)
  757. with pytest.raises(ValueError, match="^'kaiser' must have parameters.*"):
  758. windows.get_window('kaiser', 5, xp=xp)
  759. with pytest.raises(ValueError, match="^Window dpss must have one.*"):
  760. windows.get_window(('dpss', 1, 2), 5, xp=xp)
  761. with pytest.raises(ValueError, match="^'general_cosine' does not accept.*"):
  762. xp_ = xp or np # ensure parameter xp_ is not None
  763. windows.get_window(('general cosine', [1, 2]), 5, xp=xp_)
  764. @make_xp_test_case(windows.bartlett)
  765. def test_symmetric_periodic(self, xp):
  766. """Ensure that suffixes `_periodic` and `_symmetric` work for window names. """
  767. w_sym = windows.bartlett(5, sym=True, xp=xp)
  768. xp_assert_close(get_window('bartlett', 5, fftbins=False, xp=xp), w_sym)
  769. xp_assert_close(get_window('bartlett_symmetric', 5, xp=xp), w_sym)
  770. # overwrite parameter `fftbins`:
  771. xp_assert_close(get_window('bartlett_symmetric', 5, fftbins=True, xp=xp), w_sym)
  772. w_per = windows.bartlett(5, sym=False, xp=xp)
  773. xp_assert_close(get_window('bartlett', 5, xp=xp), w_per)
  774. xp_assert_close(get_window('bartlett', 5, fftbins=True, xp=xp), w_per)
  775. xp_assert_close(get_window('bartlett_periodic', 5, xp=xp), w_per)
  776. # overwrite parameter `fftbins`:
  777. xp_assert_close(get_window('bartlett_periodic', 5, fftbins=False, xp=xp),
  778. w_per)
  779. @make_xp_test_case(windows.kaiser)
  780. def test_array_as_window(self, xp):
  781. # github issue 3603
  782. osfactor = 128
  783. sig = xp.arange(128)
  784. win = windows.get_window(('kaiser', 8.0), osfactor // 2, xp=xp)
  785. mesg = "^window must" if is_cupy(xp) else "^window.shape="
  786. with assert_raises(ValueError, match=mesg):
  787. resample(sig, sig.shape[0] * osfactor, window=win)
  788. @make_xp_test_case(windows.general_cosine)
  789. def test_general_cosine(self, xp):
  790. xp_assert_close(get_window(('general_cosine', xp.asarray([0.5, 0.3, 0.2])), 4),
  791. xp.asarray([0.4, 0.3, 1, 0.3], dtype=xp.float64))
  792. xp_assert_close(get_window(('general_cosine', xp.asarray([0.5, 0.3, 0.2])), 4,
  793. fftbins=False),
  794. xp.asarray([0.4, 0.55, 0.55, 0.4], dtype=xp.float64))
  795. with pytest.raises(ValueError):
  796. get_window(('general_cosine', [0.5, 0.3, 0.2]), 4, xp=xp)
  797. @make_xp_test_case(windows.general_hamming)
  798. def test_general_hamming(self, xp):
  799. xp_assert_close(get_window(('general_hamming', 0.7), 5, xp=xp),
  800. xp.asarray([0.4, 0.6072949, 0.9427051, 0.9427051, 0.6072949],
  801. dtype=xp.float64))
  802. xp_assert_close(get_window(('general_hamming', 0.7), 5, fftbins=False, xp=xp),
  803. xp.asarray([0.4, 0.7, 1.0, 0.7, 0.4], dtype=xp.float64))
  804. @make_xp_test_case(windows.lanczos)
  805. def test_lanczos(self, xp):
  806. xp_assert_close(get_window('lanczos', 6, xp=xp),
  807. xp.asarray([0., 0.413496672, 0.826993343, 1., 0.826993343,
  808. 0.413496672], dtype=xp.float64), atol=1e-9)
  809. xp_assert_close(get_window('lanczos', 6, fftbins=False, xp=xp),
  810. xp.asarray([0., 0.504551152, 0.935489284, 0.935489284,
  811. 0.504551152, 0.], dtype=xp.float64), atol=1e-9)
  812. xp_assert_close(get_window('lanczos', 6, xp=xp),
  813. get_window('sinc', 6, xp=xp))
  814. def test_xp_default(self, xp):
  815. # no explicit xp= argument, default to numpy
  816. win = get_window('lanczos', 6)
  817. assert isinstance(win, np.ndarray)
  818. win = get_window('lanczos', 6, xp=xp)
  819. if not is_numpy(xp):
  820. assert not isinstance(win, np.ndarray)
  821. @skip_xp_backends("dask.array", reason="https://github.com/dask/dask/issues/2620")
  822. @pytest.mark.parametrize(
  823. "window,window_name,params",
  824. [
  825. make_xp_pytest_param(getattr(windows, window_name), window_name, params)
  826. for window_name, params in window_funcs
  827. ]
  828. )
  829. def test_windowfunc_basics(window, window_name, params, xp):
  830. window = getattr(windows, window_name)
  831. if is_jax(xp) and window_name in ['taylor', 'chebwin']:
  832. pytest.skip(reason=f'{window_name = }: item assignment')
  833. if window_name in ['dpss']:
  834. if is_cupy(xp):
  835. pytest.skip(reason='dpss window is not implemented for cupy')
  836. if is_torch(xp) and SCIPY_DEVICE != 'cpu':
  837. pytest.skip(reason='needs eight_tridiagonal which is CPU only')
  838. with warnings.catch_warnings():
  839. warnings.filterwarnings(
  840. "ignore", "This window is not suitable", UserWarning)
  841. # Check symmetry for odd and even lengths
  842. w1 = window(8, *params, sym=True, xp=xp)
  843. w2 = window(7, *params, sym=False, xp=xp)
  844. xp_assert_close(w1[:-1], w2)
  845. w1 = window(9, *params, sym=True, xp=xp)
  846. w2 = window(8, *params, sym=False, xp=xp)
  847. xp_assert_close(w1[:-1], w2)
  848. # Check that functions run and output lengths are correct
  849. assert window(6, *params, sym=True, xp=xp).shape[0] == 6
  850. assert window(6, *params, sym=False, xp=xp).shape[0] == 6
  851. assert window(7, *params, sym=True, xp=xp).shape[0] == 7
  852. assert window(7, *params, sym=False, xp=xp).shape[0] == 7
  853. # Check invalid lengths
  854. assert_raises(ValueError, window, 5.5, *params, xp=xp)
  855. assert_raises(ValueError, window, -7, *params, xp=xp)
  856. # Check degenerate cases
  857. xp_assert_equal(window(0, *params, sym=True, xp=xp),
  858. xp.asarray([], dtype=xp.float64))
  859. xp_assert_equal(window(0, *params, sym=False, xp=xp),
  860. xp.asarray([], dtype=xp.float64))
  861. xp_assert_equal(window(1, *params, sym=True, xp=xp),
  862. xp.asarray([1.], dtype=xp.float64))
  863. xp_assert_equal(window(1, *params, sym=False, xp=xp),
  864. xp.asarray([1.], dtype=xp.float64))
  865. # Check dtype
  866. assert window(0, *params, sym=True, xp=xp).dtype == xp.float64
  867. assert window(0, *params, sym=False, xp=xp).dtype == xp.float64
  868. assert window(1, *params, sym=True, xp=xp).dtype == xp.float64
  869. assert window(1, *params, sym=False, xp=xp).dtype == xp.float64
  870. assert window(6, *params, sym=True, xp=xp).dtype == xp.float64
  871. assert window(6, *params, sym=False, xp=xp).dtype == xp.float64
  872. # Check normalization
  873. assert xp.all(window(10, *params, sym=True, xp=xp) < 1.01)
  874. assert xp.all(window(10, *params, sym=False, xp=xp) < 1.01)
  875. assert xp.all(window(9, *params, sym=True, xp=xp) < 1.01)
  876. assert xp.all(window(9, *params, sym=False, xp=xp) < 1.01)
  877. # Check that DFT-even spectrum is purely real for odd and even
  878. res = fft(window(10, *params, sym=False, xp=xp))
  879. res = xp.imag(res)
  880. xp_assert_close(res, xp.zeros_like(res), atol=1e-14)
  881. res = fft(window(11, *params, sym=False, xp=xp))
  882. res = xp.imag(res)
  883. xp_assert_close(res, xp.zeros_like(res), atol=1e-14)
  884. @make_xp_test_case(get_window)
  885. def test_needs_params(xp):
  886. for winstr in ['kaiser', 'ksr', 'kaiser_bessel_derived', 'kbd',
  887. 'gaussian', 'gauss', 'gss',
  888. 'general gaussian', 'general_gaussian',
  889. 'general gauss', 'general_gauss', 'ggs',
  890. 'dss', 'dpss', 'general cosine', 'general_cosine',
  891. 'chebwin', 'cheb', 'general hamming', 'general_hamming',
  892. ]:
  893. assert_raises(ValueError, get_window, winstr, 7, xp=xp)
  894. _winstr = ['barthann',
  895. 'bartlett',
  896. 'blackman',
  897. 'blackmanharris',
  898. 'bohman',
  899. 'boxcar',
  900. 'cosine',
  901. 'flattop',
  902. 'hamming',
  903. 'nuttall',
  904. 'parzen',
  905. 'taylor',
  906. 'exponential',
  907. 'poisson',
  908. 'tukey',
  909. 'tuk',
  910. 'triangle',
  911. 'lanczos',
  912. 'sinc',
  913. ]
  914. @pytest.mark.parametrize(
  915. 'window,winstr',
  916. [
  917. make_xp_pytest_param(_WIN_FUNCS[winstr][0], winstr)
  918. for winstr in _winstr
  919. ]
  920. )
  921. @make_xp_test_case(get_window)
  922. def test_not_needs_params(xp, window, winstr):
  923. if is_jax(xp) and winstr in ['taylor']:
  924. pytest.skip(reason=f'{winstr}: item assignment')
  925. win = get_window(winstr, 7, xp=xp)
  926. assert win.shape[0] == 7
  927. @make_xp_test_case(windows.lanczos)
  928. def test_symmetric(xp):
  929. for win in [windows.lanczos]:
  930. # Even sampling points
  931. w = win(4096, xp=xp)
  932. flip = array_namespace(w).flip
  933. error = xp.max(xp.abs(w - flip(w)))
  934. xp_assert_equal(error, xp.asarray(0.0), check_dtype=False, check_0d=False)
  935. # Odd sampling points
  936. w = win(4097, xp=xp)
  937. error = xp.max(xp.abs(w - flip(w)))
  938. xp_assert_equal(error, xp.asarray(0.0), check_dtype=False, check_0d=False)