@@ -2316,9 +2316,13 @@ def deconvolve(signal, divisor):
2316
2316
2317
2317
2318
2318
def hilbert (x , N = None , axis = - 1 ):
2319
- """
2320
- Compute the analytic signal, using the Hilbert transform.
2319
+ r"""FFT-based computation of the analytic signal.
2321
2320
2321
+ The analytic signal is calculated by filtering out the negative frequencies and
2322
+ doubling the amplitudes of the positive frequencies in the FFT domain.
2323
+ The imaginary part of the result is the hilbert transform of the real-valued input
2324
+ signal.
2325
+
2322
2326
The transformation is done along the last axis by default.
2323
2327
2324
2328
Parameters
@@ -2337,24 +2341,27 @@ def hilbert(x, N=None, axis=-1):
2337
2341
2338
2342
Notes
2339
2343
-----
2340
- The analytic signal ``x_a(t)`` of signal ``x(t)`` is:
2344
+ The analytic signal ``x_a(t)`` of a real-valued signal ``x(t)``
2345
+ can be expressed as [1]_
2341
2346
2342
- .. math:: x_a = F^{-1}(F(x) 2U) = x + i y
2347
+ .. math:: x_a = F^{-1}(F(x) 2U) = x + i y\ ,
2343
2348
2344
2349
where `F` is the Fourier transform, `U` the unit step function,
2345
- and `y` the Hilbert transform of `x`. [1 ]_
2350
+ and `y` the Hilbert transform of `x`. [2 ]_
2346
2351
2347
2352
In other words, the negative half of the frequency spectrum is zeroed
2348
- out, turning the real-valued signal into a complex signal. The Hilbert
2353
+ out, turning the real-valued signal into a complex-valued signal. The Hilbert
2349
2354
transformed signal can be obtained from ``np.imag(hilbert(x))``, and the
2350
2355
original signal from ``np.real(hilbert(x))``.
2351
2356
2352
2357
References
2353
2358
----------
2354
2359
.. [1] Wikipedia, "Analytic signal".
2355
2360
https://en.wikipedia.org/wiki/Analytic_signal
2356
- .. [2] Leon Cohen, "Time-Frequency Analysis", 1995. Chapter 2.
2357
- .. [3] Alan V. Oppenheim, Ronald W. Schafer. Discrete-Time Signal
2361
+ .. [2] Wikipedia, "Hilbert Transform".
2362
+ https://en.wikipedia.org/wiki/Hilbert_transform
2363
+ .. [3] Leon Cohen, "Time-Frequency Analysis", 1995. Chapter 2.
2364
+ .. [4] Alan V. Oppenheim, Ronald W. Schafer. Discrete-Time Signal
2358
2365
Processing, Third Edition, 2009. Chapter 12.
2359
2366
ISBN 13: 978-1292-02572-8
2360
2367
@@ -2363,41 +2370,38 @@ def hilbert(x, N=None, axis=-1):
2363
2370
In this example we use the Hilbert transform to determine the amplitude
2364
2371
envelope and instantaneous frequency of an amplitude-modulated signal.
2365
2372
2373
+ Let's create a chirp of which the frequency increases from 20 Hz to 100 Hz and
2374
+ apply an amplitude modulation:
2375
+
2366
2376
>>> import numpy as np
2367
2377
>>> import matplotlib.pyplot as plt
2368
2378
>>> from scipy.signal import hilbert, chirp
2369
-
2370
- >>> duration = 1.0
2371
- >>> fs = 400.0
2372
- >>> samples = int(fs*duration)
2373
- >>> t = np.arange(samples) / fs
2374
-
2375
- We create a chirp of which the frequency increases from 20 Hz to 100 Hz and
2376
- apply an amplitude modulation.
2377
-
2379
+ ...
2380
+ >>> duration, fs = 1, 400 # 1 s signal with sampling frequency of 400 Hz
2381
+ >>> t = np.arange(int(fs*duration)) / fs # timestamps of samples
2378
2382
>>> signal = chirp(t, 20.0, t[-1], 100.0)
2379
2383
>>> signal *= (1.0 + 0.5 * np.sin(2.0*np.pi*3.0*t) )
2380
2384
2381
- The amplitude envelope is given by magnitude of the analytic signal. The
2385
+ The amplitude envelope is given by the magnitude of the analytic signal. The
2382
2386
instantaneous frequency can be obtained by differentiating the
2383
2387
instantaneous phase in respect to time. The instantaneous phase corresponds
2384
2388
to the phase angle of the analytic signal.
2385
2389
2386
2390
>>> analytic_signal = hilbert(signal)
2387
2391
>>> amplitude_envelope = np.abs(analytic_signal)
2388
2392
>>> instantaneous_phase = np.unwrap(np.angle(analytic_signal))
2389
- >>> instantaneous_frequency = ( np.diff(instantaneous_phase) /
2390
- ... (2.0*np.pi) * fs)
2391
-
2392
- >>> fig, ( ax0, ax1) = plt.subplots(nrows=2 )
2393
- >>> ax0.plot(t, signal, label='signal' )
2394
- >>> ax0.plot(t, amplitude_envelope , label='envelope ')
2395
- >>> ax0.set_xlabel("time in seconds" )
2393
+ >>> instantaneous_frequency = np.diff(instantaneous_phase) / (2.0*np.pi) * fs
2394
+ ...
2395
+ >>> fig, (ax0, ax1) = plt.subplots(nrows=2, sharex='all', tight_layout=True)
2396
+ >>> ax0.set_title("Amplitude-modulated Chirp Signal" )
2397
+ >>> ax0.set_ylabel("Amplitude" )
2398
+ >>> ax0.plot(t, signal , label='Signal ')
2399
+ >>> ax0.plot(t, amplitude_envelope, label='Envelope' )
2396
2400
>>> ax0.legend()
2397
- >>> ax1.plot(t[1:], instantaneous_frequency )
2398
- >>> ax1.set_xlabel("time in seconds" )
2399
- >>> ax1.set_ylim(0.0, 120.0 )
2400
- >>> fig.tight_layout ()
2401
+ >>> ax1.set(xlabel="Time in seconds", ylabel="Phase in rad", ylim=(0, 120) )
2402
+ >>> ax1.plot(t[1:], instantaneous_frequency, 'C2-', label='Instantaneous Phase' )
2403
+ >>> ax1.legend( )
2404
+ >>> plt.show ()
2401
2405
2402
2406
"""
2403
2407
x = np .asarray (x )
0 commit comments