Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Document coefficient vector layout #37

Open
mschauer opened this issue May 30, 2018 · 4 comments
Open

Document coefficient vector layout #37

mschauer opened this issue May 30, 2018 · 4 comments

Comments

@mschauer
Copy link
Contributor

I a hope I am not missing anything but I could not find the layout of the coefficient vector. Inspecting the size of the coefficient also does not immediately suggest a certain order to me.

@gummif
Copy link
Member

gummif commented May 30, 2018

You are right, this is not well documented. There are some utiliities in https://github.com/JuliaDSP/Wavelets.jl/blob/master/src/mod/Util.jl , especially the dyadicdetailrange and dyadicscalingrange functions (assuming the signal is dyadic, a power of 2). There is also detailrange for other more general sizes of signals, but no scalingrange function. These give you the range of the vector that corresponds to each level of the coefficients. I will take a look at adding missing functions and documenting these better.

@gugushvili
Copy link

What would be very helpful is to have equivalents of the functions accessC, accessD, putC and putD from the wavethresh package in R. These allow to access or replace smooth and detail coefficients in entire resolution levels of DWT. See https://cran.r-project.org/web/packages/wavethresh/

@mschauer
Copy link
Contributor Author

Ah, so my first guess for the memory layout was actually right according to

dyadicdetailrange(j::Integer) = (2^j+1):(2^(j+1))

I just did not believe it: i started with smooth signal of length 128 and expected the coefficients xt[65:128] small and of the same order, but that was not the case for the example I looked at:

julia> t = linspace(0,1,128)
0.0:0.007874015748031496:1.0

julia> x = sin.(2pi*t)

julia> using Wavelets

julia> xt = dwt(x, wavelet(WT.db4, WT.Periodic))
128-element Array{Float64,1}:
 -5.48173e-16
  5.31657    
 -4.19097    
  4.11042    
 -0.290043   
  0.489353   
  0.33211    
 -0.574409   
  0.0127308  
  0.00217123 
  ⋮          
 -2.50345e-6 
 -2.30573e-6 
 -2.08546e-6 
 -1.84479e-6 
 -1.58608e-6 
 -1.31184e-6 
 -1.02478e-6 
 -7.27685e-7 
 -4.23475e-7 



julia> xt[65:80]
16-element Array{Float64,1}:
  0.0239865 
 -0.00864036
  0.00214933
  6.5407e-7 
  9.53168e-7
  1.24294e-6
  1.52056e-6
  1.78329e-6
  2.02859e-6
  2.25404e-6
  2.45744e-6
  2.63679e-6
  2.79036e-6
  2.91662e-6
  3.01436e-6
  3.0826e-6 

Maybe these are boundary effects?

@mschauer
Copy link
Contributor Author

mschauer commented Jun 1, 2018

Ah, I found my mistake, the time points t = linspace(0,1,128) were not at equal distance on a circle of length 1, this works.

julia> t = (1:128)./128
julia> x = sin.(2pi*t)
julia> xt = dwt(x, wavelet(WT.db4, WT.Periodic))
julia> xt[65:80]
16-element Array{Float64,1}:
 -1.10697e-7
  1.86888e-7
  4.82674e-7
  7.73811e-7
  1.0575e-6 
  1.331e-6  
  1.59168e-6
  1.83703e-6
  2.06469e-6
  2.27247e-6
  2.45836e-6
  2.62058e-6
  2.75756e-6
  2.86798e-6
  2.95079e-6
  3.00517e-6

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants