Use Bessels.jl for Bessel function calculations

Replace calls to `SpecialFunctions.besseli(0, ...)` with calls to
`Bessels.besseli0(...)`, adding Bessels.jl
(https://github.com/JuliaMath/Bessels.jl) as a dependency.

When resampling signals with a very low or high arbitrary rate, the call
to `Windows.kaiser` is what dominates.  Profiling this, most of the time
is taken up by the call to `SpecialFunctions.besseli`. This allocates,
but also the computation is just fundamentally slow for these arguments.

Replacing these calls with those to `Bessels.besseli0` gives a factor
of 17 speedup when resampling a 1000-point random trace to a rate of
0.02.

Closes JuliaDSP#506.

Project.toml

@@ -3,6 +3,7 @@ uuid = "717857b8-e6f2-59f4-9121-6e50c889abd2"
 version = "0.8.0"
 
 [deps]
+Bessels = "0e736298-9ec6-45e8-9647-e4fc86a2fe38"
 Compat = "34da2185-b29b-5c13-b0c7-acf172513d20"
 FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
 IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e"
@@ -14,6 +15,7 @@ SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 
 [compat]
+Bessels = "0.2"
 Compat = "4.1"
 DelimitedFiles = "1.6"
 FFTW = "1.8"

src/windows.jl

@@ -1,6 +1,6 @@
 module Windows
 using ..Util
-import SpecialFunctions: besseli
+using Bessels: besseli0
 using LinearAlgebra: Diagonal, SymTridiagonal, eigen!, mul!, rmul!
 using FFTW
 
@@ -501,9 +501,9 @@ $(twoD_docs("α"))
 $zerophase_docs
 """
 function kaiser(n::Integer, α::Real; padding::Integer=0, zerophase::Bool=false)
-    pf = 1.0/besseli(0,pi*α)
+    pf = 1.0/besseli0(pi*α)
     makewindow(n, padding, zerophase) do x
-        pf*besseli(0, pi*α*(sqrt(1 - (2x)^2)))
+        pf*besseli0(pi*α*(sqrt(1 - (2x)^2)))
     end
 end