lanczos_wgts.py 1.04 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
""" Draw a set of Lanczos Weights based on
https://scitools.org.uk/iris/docs/v1.2/examples/graphics/SOI_filtering.html
"""

__author__ = "Carl Schreck"
__email__ = "cjschrec@ncsu.edu"
__copyright__ = "Copyright 2020, North Carolina State University"
__license__ = "BSD-3.0"

import numpy as np
import matplotlib.pyplot as plt

def low_pass_weights(window, cutoff):
    """Calculate weights for a low pass Lanczos filter.

    Args:

    window: int
        The length of the filter window.

    cutoff: float
        The cutoff frequency in inverse time steps.

    """
    order = ((window - 1) // 2 ) + 1
    nwts = 2 * order + 1
    w = np.zeros([nwts])
    n = nwts // 2
    w[n] = 2 * cutoff
    k = np.arange(1., n)
    sigma = np.sin(np.pi * k / n) * n / (np.pi * k)
    firstfactor = np.sin(2. * np.pi * cutoff * k) / (np.pi * k)
    w[n-1:0:-1] = firstfactor * sigma
    w[n+1:-1] = firstfactor * sigma
    return w[1:-1]

y = low_pass_weights(window=41, cutoff=0.10)
x = np.arange(-20,21)
plt.plot(x, y)
plt.savefig('figures/lanczos_wgts.png')
plt.show()