| Help with WV_CWT [message #81684] |
Mon, 08 October 2012 14:19 |
Balthasar Indermuehle
Messages: 22
Registered: August 2012
|
Junior Member |
|
|
Hi all,
Exelis support has been dwelling on this for over 24 hours without answer thus far... does anyone here have a clue how this works?
I'm in need of a frequency power vs time analysis for a signal I'm examining. For that purpose I wanted to use the Morlet wavelet transform, using the WV_CWT function described in exelis's advanced signal processing tutorial here: http://www.exelisvis.com/portals/0/tutorials/idl/Advanced_Si gnal_Processing.pdf
Unfortunately, the tutorial gives little to no information on what it is the wavelet transform returns. It looks like the x axis of the data array it returns is the time, at least the number of elements coincides with the number of elements in the time series in my data, but the y axis which should be frequency I have no idea how to convert into actual frequency. This simple test should plot a line at 20 Hz and another one at 30 Hz:
x = DINDGEN(5000)*10d/5000 ; 5000 points, 10 seconds
y = sin(20 * 2.0D * !DPI * x) + $ ; generate a 20 Hz signal
sin(30 * 2.0D * !DPI * x) ; overlay a 20 Hz signal
wave = WV_CWT(y, 'morlet', 6, /PAD, SCALE=scales)
wavePower = ABS(wave)^2
; And plot the wavelet.
cgLOADCT, 13, RGB_TABLE=palette
cgIMAGE, wavePower, PALETTE=palette, POSITION=[0.1, 0.1, 0.9, 0.9], XRANGE=[0,90], YRANGE=[0,250 ], AXKEYWORDS={TICKLEN: -0.025}, MARGIN=0.5, /AXES
It plots two very weak lines. Scaling? Not a clue.
You can run an FFT on the y and that correctly extracts the 20 and 30 Hz components:
samplingInterval = x[1] - x[0]
samplingFreq = 1.0D / samplingInterval
Nyquist = 1.0D / (2.0D * samplingInterval)
freq = DINDGEN(N_ELEMENTS(x)) * samplingFreq / (N_ELEMENTS(x)-1)
dft = FFT(y)
; Now fold the transform
dft = dft[0:N_ELEMENTS(dft)/2+1] ; cut in half
dft /= N_ELEMENTS(dft) ; scale
dft[1:N_ELEMENTS(dft)-1] *= 2.0D ; compensate for the discarded half
magSpec = ABS(dft) ; magnitude spectrum
maxy = max(magSpec,min=miny)
CGPLOT, freq, magSpec, TITLE='FFT Magnitude', YRANGE=[miny,maxy], XRANGE=[0,nyquist], XTITLE='Freq [Hz]', /YLOG
Can anyone help me to make sense of the output of WV_CWT?
Cheers
- Balt
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
