|
8 | 8 |
|
9 | 9 |
|
10 | 10 | def ideal_bandpass(data, sample_period, bandpass_freqs): |
11 | | - # Derived from YAN Chao-Gan 120504 based on REST. |
| 11 | + """ |
| 12 | + Apply ideal bandpass filtering to a 1D time series data using FFT. Derived from YAN Chao-Gan 120504 based on REST. |
| 13 | +
|
| 14 | + Parameters |
| 15 | + ---------- |
| 16 | + data : NDArray |
| 17 | + 1D time series data to be filtered. |
| 18 | + sample_period : float |
| 19 | + Length of sampling period in seconds. |
| 20 | + bandpass_freqs : tuple |
| 21 | + Tuple containing the bandpass frequencies (LowCutoff, HighCutoff). |
| 22 | +
|
| 23 | + Returns |
| 24 | + ------- |
| 25 | + NDArray |
| 26 | + Filtered time series data. |
| 27 | +
|
| 28 | + """ |
12 | 29 | sample_freq = 1.0 / sample_period |
13 | 30 | sample_length = data.shape[0] |
| 31 | + nyquist_freq = sample_freq / 2.0 |
14 | 32 |
|
15 | | - data_p = np.zeros(int(2 ** np.ceil(np.log2(sample_length)))) |
| 33 | + # Length of zero-padded data for efficient FFT |
| 34 | + N = int(2 ** np.ceil(np.log2(len(data)))) |
| 35 | + data_p = np.zeros(N) |
16 | 36 | data_p[:sample_length] = data |
17 | 37 |
|
18 | 38 | LowCutoff, HighCutoff = bandpass_freqs |
19 | 39 |
|
20 | 40 | if LowCutoff is None: # No lower cutoff (low-pass filter) |
21 | 41 | low_cutoff_i = 0 |
22 | | - elif LowCutoff > sample_freq / 2.0: |
| 42 | + elif LowCutoff > nyquist_freq: |
23 | 43 | # Cutoff beyond fs/2 (all-stop filter) |
24 | | - low_cutoff_i = int(data_p.shape[0] / 2) |
| 44 | + low_cutoff_i = int(N / 2) |
25 | 45 | else: |
26 | | - low_cutoff_i = np.ceil(LowCutoff * data_p.shape[0] * sample_period).astype( |
27 | | - "int" |
28 | | - ) |
| 46 | + low_cutoff_i = np.ceil(LowCutoff * N * sample_period).astype("int") |
29 | 47 |
|
30 | | - if HighCutoff > sample_freq / 2.0 or HighCutoff is None: |
| 48 | + if HighCutoff > nyquist_freq or HighCutoff is None: |
31 | 49 | # Cutoff beyond fs/2 or unspecified (become a highpass filter) |
32 | | - high_cutoff_i = int(data_p.shape[0] / 2) |
| 50 | + high_cutoff_i = int(N / 2) |
33 | 51 | else: |
34 | | - high_cutoff_i = np.fix(HighCutoff * data_p.shape[0] * sample_period).astype( |
35 | | - "int" |
36 | | - ) |
| 52 | + high_cutoff_i = np.fix(HighCutoff * N * sample_period).astype("int") |
37 | 53 |
|
38 | 54 | freq_mask = np.zeros_like(data_p, dtype="bool") |
39 | 55 | freq_mask[low_cutoff_i : high_cutoff_i + 1] = True |
40 | | - freq_mask[data_p.shape[0] - high_cutoff_i : data_p.shape[0] + 1 - low_cutoff_i] = ( |
41 | | - True |
42 | | - ) |
| 56 | + freq_mask[N - high_cutoff_i : N + 1 - low_cutoff_i] = True |
43 | 57 |
|
44 | 58 | f_data = fft(data_p) |
45 | | - f_data[freq_mask is not True] = 0.0 |
| 59 | + f_data[~freq_mask] = 0.0 |
46 | 60 | return np.real_if_close(ifft(f_data)[:sample_length]) |
47 | 61 |
|
48 | 62 |
|
|
0 commit comments