Merge pull request #1 from JuliaDSP/gfa/types

API and type hierarchy change

README.md

@@ -5,7 +5,7 @@
 [![Build Status](https://travis-ci.org/JuliaDSP/Wavelets.jl.svg?branch=master)](https://travis-ci.org/JuliaDSP/Wavelets.jl)
 [![Coverage Status](https://coveralls.io/repos/JuliaDSP/Wavelets.jl/badge.png?branch=master)](https://coveralls.io/r/JuliaDSP/Wavelets.jl?branch=master)
 
-A [Julia](https://github.com/JuliaLang/julia) package for fast wavelet transforms (1-d, 2-d, by filtering or lifting).
+A [Julia](https://github.com/JuliaLang/julia) package for fast wavelet transforms (1-D, 2-D, by filtering or lifting).
 
 * 1st generation wavelets using filter banks (periodic and orthogonal). Filters are included for the following types: Haar, Daubechies, Coiflet, Symmlet, Battle-Lemarie, Beylkin, Vaidyanathan.
 
@@ -15,7 +15,7 @@ A [Julia](https://github.com/JuliaLang/julia) package for fast wavelet transform
 
 * Wavelet utilities e.g. indexing and size calculation, scaling and wavelet functions computation, test functions, up and down sampling, filter mirrors, coefficient counting, inplace circshifts, and more.
 
-* Plotting/visualization utilities for 1-d and 2-d signals.
+* Plotting/visualization utilities for 1-D and 2-D signals.
 
 Roughly 20x speedup and 50x less memory usage than [this](https://github.com/tomaskrehlik/Wavelets) implementation of `dwt`. Loosely inspired by [this](https://github.com/tomaskrehlik/Wavelets) and [this](http://statweb.stanford.edu/~wavelab). 
 
@@ -32,36 +32,34 @@ julia> Pkg.add("Wavelets")
 julia> using Wavelets
 ```
 
-A rough idea of the API:
+A few usage examples:
 
 ```julia
 # the simplest way to transform a signal x is
-xt = fwt(x, wavelet("db2"))
+xt = dwt(x, wavelet("db2"))
 
 # the transform type can be more explicitly specified
 # set up wavelet type (filter, Periodic, Orthogonal, 4 vanishing moments)
-wt = wavelet("Coiflet", 4)  # or
-wt = wavelet("Coiflet", 4, transform="filter")  # or
-wt = waveletfilter("coif4", boundary="periodic", t="orthogonal")  # or even
-wt = wavelet("coif4", transform="filter", boundary="P", t="O")
-# which equivalent to 
-wt = POfilter("coif4")
+wt = wavelet("Coiflet", 4, transform="filter", boundary="per")  # or
+wt = waveletfilter("coif4")
+# which is equivalent to 
+wt = OrthoFilter("coif4")
 # the object wt determines the transform type 
 # wt now contains instructions for a periodic biorthogonal CDF 9/7 lifting scheme
 wt = waveletls("cdf9/7")
 # xt is a 5 level transform of vector x
-xt = fwt(x, 5, wt)
+xt = dwt(x, wt, 5)
 # inverse tranform
-x2 = iwt(xt, 5, wt)
+xti = idwt(xt, wt, 5)
 # a full transform
-xt = fwt(x, wt)
+xt = dwt(x, wt)
 
 # scaling filters is easy
 wt = scale(wt, 1/sqrt(2))
 # signals can be transformed inplace with a low-level command
 # requiring very little memory allocation (especially for L=1 for filters)
-dwt!(x, L, wt, true)      # inplace (lifting)
-dwt!(xt, x, L, wt, true)  # write to xt (filter)
+dwt!(x, wt, L, true)      # inplace (lifting)
+dwt!(xt, x, wt, L, true)  # write to xt (filter)
 
 # denoising with default parameters (VisuShrink hard thresholding)
 x0 = testfunction(n, "HeaviSine")
@@ -70,23 +68,44 @@ y = denoise(x)
 
 # plotting utilities 1-d (see images and code in /example)
 x = testfunction(n, "Bumps")
-y = fwt(x, GPLS("cdf9/7"))
+y = dwt(x, waveletls("cdf9/7"))
 d,l = wplotdots(y, 0.1, n)
 A = wplotim(y)
 # plotting utilities 2-d
 img = imread("lena.png")
 x = permutedims(img.data, [ndims(img.data):-1:1])
 L = 2
-xts = wplotim(x, L, POfilter("db3"))
+xts = wplotim(x, L, waveletfilter("db3"))
 ```
 
 ![Bumps](/example/transform1d_bumps.png)
 
 ![Lena](/example/transform2d_lena.jpg)
 
+
+API
+---------
+
+#### Wavelet transforms and types
+```julia
+# Type construction,
+# also accept (class::String, n::Union(Integer,String); ...)
+wavelet(name::String; transform::String="filter", boundary::String="per")
+waveletfilter(name::String; boundary::String="per")
+waveletls(name::String; boundary::String="per")
+# DWT (discrete wavelet transform)
+dwt(x::AbstractArray, wt::DiscreteWavelet, L::Integer=nscales(size(x,1)))
+idwt(x::AbstractArray, wt::DiscreteWavelet, L::Integer=nscales(size(x,1)))
+dwt!(y::AbstractArray, x::AbstractArray, filter::OrthoFilter, L::Integer, fw::Bool)
+dwt!(y::AbstractArray, scheme::GLS, L::Integer, fw::Bool)
+# DWTC (column-wise discrete wavelet transform)
+dwtc(x::AbstractArray, wt::DiscreteWavelet, L::Integer=nscales(size(x,1)))
+idwtc(x::AbstractArray, wt::DiscreteWavelet, L::Integer=nscales(size(x,1)))
+```
+
 #### Wavelet class information
 
-| Long name | Short | Type | Numbers |
+| Class | Short | Type | Numbers |
 |:------- |:------ |:----- |:----- |
 | `Haar` | `haar` | Ortho |   |
 | `Coiflet` | `coif` | Ortho | 2:2:8 |
@@ -101,22 +120,22 @@ xts = wplotim(x, L, POfilter("db3"))
 Benchmarks
 ---------
 
-Timing of `fwt` (using `db2` filter of length 4) and `fft`. The wavelet transforms are faster and use less memory than `fft` in all cases. `fwt` by lifting is currently an order of magnitude faster than by filtering.
+Timing of `dwt` (using `db2` filter of length 4) and `fft`. The wavelet transforms are faster and use less memory than `fft` in all cases. `dwt` by lifting is currently an order of magnitude faster than by filtering.
 
 ```julia
 # 10 iterations
-fwt by filtering (N=1048576), 18 levels
+dwt by filtering (N=1048576), 18 levels
 elapsed time: 1.262672396 seconds (125866088 bytes allocated, 2.35% gc time)
-fwt by lifting (N=1048576), 18 levels
+dwt by lifting (N=1048576), 18 levels
 elapsed time: 0.153162937 seconds (104927608 bytes allocated, 13.18% gc time)
 fft (N=1048576), (FFTW)
 elapsed time: 2.836224168 seconds (587236088 bytes allocated, 4.83% gc time)
 ```
 
-For 2D transforms (the `fwt` using a CDF 9/7 lifting scheme):
+For 2D transforms (the `dwt` using a CDF 9/7 lifting scheme):
 ```julia
 # 10 iterations
-fwt by lifting (N=1024x1024), 10 levels
+dwt by lifting (N=1024x1024), 10 levels
 elapsed time: 0.952415385 seconds (109075368 bytes allocated, 2.14% gc time)
 fft (N=1024 x 1024), (FFTW)
 elapsed time: 2.945895417 seconds (587236728 bytes allocated, 4.69% gc time)
@@ -127,7 +146,7 @@ By using the low-level function `dwt!` and pre-allocating temporary arrays, sign
 # 1000 iterations
 dwt! (N=32768), 13 levels
 elapsed time: 5.081105993 seconds (701512 bytes allocated)
-fwt (N=32768), 13 levels
+dwt (N=32768), 13 levels
 elapsed time: 5.151621451 seconds (395317512 bytes allocated, 1.62% gc time)
 ```
 
@@ -159,7 +178,6 @@ To-do list
 ---------
 
 * don't use sub arrays in 2d
-* error: y[ind...] = fwt(xc, L, wt)
 * Boundary orthogonal wavelets
 * Define more lifting schemes
 * Redundant transforms and wavelet packets

benchmark/bm_denoise.jl

@@ -2,8 +2,8 @@
 #include("setup_1d.jl")
 using Wavelets
 
-wf = POfilter("db2")
-wl = GPLS("db2")
+wf = waveletfilter("db2")
+wl = waveletls("db2")
 N = 8*1024;
 x0 = rand(N);
 L = nscales(N)  # int(log2(N)-2)

benchmark/bm_dwt.jl

@@ -1,7 +1,7 @@
 
 using Wavelets
 
-wf = POfilter("db6")
+wf = waveletfilter("db6")
 N = 1024*32;
 x0 = rand(N);
 L = int(log2(N)-2)
@@ -24,7 +24,7 @@ y=x0*0.0
 
 println("dwt! (N=",N,"), ", L, " levels")
 #f(x0, L, wf) = for i = 1:tn; fwt(x0, L, wf); end
-f(x0, L, wf) = for i = 1:tn; dwt!(y, x0, L, wf, fw, dcfilter, scfilter, si, snew); end
+f(x0, L, wf) = for i = 1:tn; dwt!(y, x0, wf, L, fw, dcfilter, scfilter, si, snew); end
 f(x0, L, wf);
 @time f(x0, L, wf);
 

benchmark/bm_fwt.jl

@@ -2,7 +2,7 @@
 include("setup_1d.jl")
 
 println("fwt by filtering (N=",N,"), ", L, " levels")
-f(x0, L, wf) = for i = 1:tn; fwt(x0, L, wf); end
+f(x0, L, wf) = for i = 1:tn; dwt(x0, wf, L); end
 f(x0, L, wf);
 @time f(x0, L, wf);
 

benchmark/bm_fwt2.jl

@@ -2,7 +2,7 @@
 include("setup_2d.jl")
 
 println("fwt (N=",N,"x",N,"), ", L, " levels")
-f(x0, L, wf) = for i = 1:tn; fwt(x0, L, wf); end
+f(x0, L, wf) = for i = 1:tn; dwt(x0, wf, L); end
 f(x0, L, wf);
 @time f(x0, L, wf);
 

benchmark/bm_iwt.jl

@@ -1,14 +1,14 @@
 
 using Wavelets
 
-wf = POfilter("db6")
+wf = waveletfilter("db6")
 N = 1024*32;
 x0 = rand(N);
 L = int(log2(N)-2)
 tn = 1000
 
 println("iwt (N=",N,"), ", L, " levels")
-f(x0, L, wf) = for i = 1:tn; iwt(x0, L, wf); end
+f(x0, L, wf) = for i = 1:tn; idwt(x0, wf, L); end
 f(x0, L, wf);
 @time f(x0, L, wf);
 

benchmark/bm_iwt2.jl

@@ -1,15 +1,15 @@
 
 using Wavelets
 
-wf = POfilter("db4")
+wf = waveletfilter("db4")
 N = 1024;
 x0 = rand(N,N);
 L = int(log2(N)-2)
 tn = 10
 #y=x0*0.0
 
 println("iwt (N=",N,"x",N,"), ", L, " levels")
-f(x0, L, wf) = for i = 1:tn; iwt(x0, L, wf); end
+f(x0, L, wf) = for i = 1:tn; idwt(x0, wf, L); end
 f(x0, L, wf);
 @time f(x0, L, wf);
 

benchmark/bm_lsfwt.jl

@@ -2,7 +2,7 @@
 include("setup_1d.jl")
 
 println("fwt by lifting (N=",N,"), ", L, " levels")
-f(x0, L, wl) = for i = 1:tn; fwt(x0, L, wl); end
+f(x0, L, wl) = for i = 1:tn; dwt(x0, wl, L); end
 f(x0, L, wl);
 @time f(x0, L, wl);
 

benchmark/bm_lsfwt2.jl

@@ -2,7 +2,7 @@
 include("setup_2d.jl")
 
 println("fwt by lifting (N=",N,"x",N,"), ", L, " levels")
-f(x0, L, wl) = for i = 1:tn; fwt(x0, L, wl); end
+f(x0, L, wl) = for i = 1:tn; dwt(x0, wl, L); end
 f(x0, L, wl);
 @time f(x0, L, wl);
 

benchmark/setup_1d.jl

@@ -1,8 +1,8 @@
 
 using Wavelets
 
-wf = POfilter("db2")
-wl = GPLS("db2")
+wf = waveletfilter("db2")
+wl = waveletls("db2")
 N = 1024*1024;
 x0 = rand(N);
 L = nscales(N)  # int(log2(N)-2)

benchmark/setup_2d.jl

@@ -1,8 +1,8 @@
 
 using Wavelets
 
-wf = POfilter("db4")
-wl = GPLS("cdf9/7")
+wf = waveletfilter("db4")
+wl = waveletls("cdf9/7")
 N = 1024;
 x0 = rand(N,N);
 L = nscales(N)  # int(log2(N)-2)

example/transform1d.jl

@@ -5,7 +5,7 @@ using PyPlot
 J = 11
 n = 2^J
 x = testfunction(n,"Bumps")
-y = fwt(x, GPLS("cdf9/7"))
+y = dwt(x, waveletls("cdf9/7"))
 d,l = wplotdots(y, 0.1, n)
 A = wplotim(y)
 

example/transform2d.jl

@@ -5,7 +5,7 @@ using Images  # for imread and imwrite
 img = imread("lena.tiff")
 x = permutedims(img.data, [ndims(img.data):-1:1])
 L = 2
-xts = wplotim(x, L, POfilter("db3"))
+xts = wplotim(x, L, waveletfilter("db3"))
 
 imwrite(xts, "lena_2d.png")
 # convert with ImageMagick

src/plot.jl

@@ -1,6 +1,5 @@
 module Plot
-using ..Util, ..Transforms
-import ..WaveletTypes: WaveletType
+using ..Util, ..WaveletTypes, ..Transforms
 export wplotdots, wplotim
 
 # PLOTTING UTILITIES
@@ -52,7 +51,7 @@ end
 
 # return an array of scaled detail coefficients and unscaled scaling coefficients
 # ready to be plotted as an image
-function wplotim(x::AbstractArray, L::Integer, wt::Union(WaveletType,Nothing)=nothing; wabs::Bool=true, power::Real=0.7, pnorm::Real=1)
+function wplotim(x::AbstractArray, L::Integer, wt::Union(DiscreteWavelet,Nothing)=nothing; wabs::Bool=true, power::Real=0.7, pnorm::Real=1)
     dim = ndims(x)
     (dim == 2 || dim == 3) || error("dimension ",dim," not supported")
     n = size(x,1)
@@ -65,9 +64,9 @@ function wplotim(x::AbstractArray, L::Integer, wt::Union(WaveletType,Nothing)=no
     # do wavelet transform
     if wt != nothing
         if size(x,3)>1
-            x = fwtc(x, L, wt)
+            x = dwtc(x, wt, L)
         else
-            x = fwt(x, L, wt)
+            x = dwt(x, wt, L)
         end
     end