UCD4IDS / ContinuousWavelets.jl

wide array of continuous wavelet transforms using Julia
https://dsweber2.github.io/ContinuousWavelets.jl/dev/
Other
19 stars 6 forks source link

Continuous Wavelet Transform Crashes for Complex Inputs #51

Open mattcbro opened 2 months ago

mattcbro commented 2 months ago

My use case for wavelets involves analyzing complex baseband signals like chirps and tones. This crashes the wavelet transform, e.g.

 yt = cwt(chirp, wave)
┌ Warning: the lowest frequency wavelet has more than 1% its max at zero, so it may not be analytic. Think carefully
│   lowAprxAnalyt = 0.050211
└ @ ContinuousWavelets ~/.julia/packages/ContinuousWavelets/pycif/src/sanityChecks.jl:7
ERROR: DimensionMismatch: arrays could not be broadcast to a common size; got a dimension with lengths 4096 and 2049
Stacktrace:
  [1] _bcs1
    @ ./broadcast.jl:555 [inlined]
  [2] _bcs
    @ ./broadcast.jl:549 [inlined]
  [3] broadcast_shape
    @ ./broadcast.jl:543 [inlined]
  [4] combine_axes
    @ ./broadcast.jl:524 [inlined]
  [5] instantiate
    @ ./broadcast.jl:306 [inlined]
  [6] materialize
    @ ./broadcast.jl:903 [inlined]
  [7] analyticTransformComplex!(wave::Array{…}, daughters::Matrix{…}, x̂::Matrix{…}, fftPlan::FFTW.cFFTWPlan{…}, averagingType::NoAve)
    @ ContinuousWavelets ~/.julia/packages/ContinuousWavelets/pycif/src/apply.jl:172
  [8] cwt(Y::Vector{…}, cWav::CWT{…}, daughters::Matrix{…}, fftPlans::Int64)
    @ ContinuousWavelets ~/.julia/packages/ContinuousWavelets/pycif/src/apply.jl:56
  [9] cwt
    @ ~/.julia/packages/ContinuousWavelets/pycif/src/apply.jl:23 [inlined]
 [10] #cwt#21
    @ ~/.julia/packages/ContinuousWavelets/pycif/src/apply.jl:240 [inlined]
 [11] cwt(Y::Vector{ComplexF64}, c::CWT{SymBoundary, Float64, Morlet, Float64, true})
    @ ContinuousWavelets ~/.julia/packages/ContinuousWavelets/pycif/src/apply.jl:236
 [12] top-level scope
    @ /data/projects/Airanaculus/GPSSignalClass/signalclass.jl:38
Some type information was truncated. Use `show(err)` to see complete types.

chirp is a length 2048 complex chirp signal.

I think I can work around this by separately taking the cwt() operation of the real and imaginary parts, but I don't see why these transforms don't work out of the box on a complex vector.

dsweber2 commented 2 months ago

Could you provide your definitions for chirp and wave? When I just take the front page example and do cwt(f * (1.0+1.0im), c), it works just fine, so I suspect its either the particular wavelet, or there may be something odd about the shape of chirp.

mattcbro commented 2 months ago

It's nothing to do with the chirp. It's the fact that the waveform is complex. As a matter of fact this problem occurs throughout a lot of the DSP related packages. See the simple example below:

using ContinuousWavelets

x = randn(1024) + 1im * randn(1024)
wave = wavelet(Morlet(π), averagingType=NoAve(), β=2)

# works
yt = cwt(real.(x), wave)

# doesn't work
yt = cwt(x, wave)
dsweber2 commented 2 months ago

Ok this narrows it down to averagingType=NoAve(), since the same example works if I set averagingType=Father(). Shouldn't take too long to fix. Thanks for the issue