code: def lfilter(b, a, x, axis=-1, zi=None):
      """
      Filter data along one-dimension with an IIR or FIR filter.
  
      Filter a data sequence, `x`, using a digital filter.  This works for many
      fundamental data types (including Object type).  The filter is a direct
      form II transposed implementation of the standard difference equation
      (see Notes).
  
      The function `sosfilt` (and filter design using ``output='sos'``) should be
      preferred over `lfilter` for most filtering tasks, as second-order sections
      have fewer numerical problems.
  
      Parameters
      ----------
      b : array_like
          The numerator coefficient vector in a 1-D sequence.
      a : array_like
          The denominator coefficient vector in a 1-D sequence.  If ``a[0]``
          is not 1, then both `a` and `b` are normalized by ``a[0]``.
      x : array_like
          An N-dimensional input array.
      axis : int, optional
          The axis of the input data array along which to apply the
          linear filter. The filter is applied to each subarray along
          this axis.  Default is -1.
      zi : array_like, optional
          Initial conditions for the filter delays.  It is a vector
          (or array of vectors for an N-dimensional input) of length
          ``len(b) + len(a) - 2``. The first ``len(b)`` numbers correspond to the
          last elements of the previous input and the last ``len(a)`` to the last
          elements of the previous output. If `zi` is None or is not given then
          initial rest is assumed.  See `lfiltic` for more information.
  
          **Note**: This argument differs from dimensions from the SciPy
          implementation! However, as long as they are chained from the same
          library, the output result will be the same. Please make sure to use
          the `zi` from CuPy calls and not from SciPy. This due to the parallel
          nature of this implementation as opposed to the serial one in SciPy.
  
      Returns
      -------
      y : array
          The output of the digital filter.
      zf : array, optional
          If `zi` is None, this is not returned, otherwise, `zf` holds the
          final filter delay values.
  
      See Also
      --------
      lfiltic : Construct initial conditions for `lfilter`.
      lfilter_zi : Compute initial state (steady state of step response) for
                   `lfilter`.
      filtfilt : A forward-backward filter, to obtain a filter with zero phase.
      savgol_filter : A Savitzky-Golay filter.
      sosfilt: Filter data using cascaded second-order sections.
      sosfiltfilt: A forward-backward filter using second-order sections.
  
      Notes
      -----
      The filter function is implemented as a direct II transposed structure.
      This means that the filter implements::
  
            a[0]*y[n] = b[0]*x[n] + b[1]*x[n-1] + ... + b[M]*x[n-M]
                                  - a[1]*y[n-1] - ... - a[N]*y[n-N]
  
      where `M` is the degree of the numerator, `N` is the degree of the
      denominator, `n` is the sample number and `L` denotes the length of the
      input.  It is implemented by computing first the FIR part and then
      computing the IIR part from it::
  
               a[0] * y = r(f(x, b), a)
               f(x, b)[n] = b[0]*x[n] + b[1]*x[n-1] + ... + b[M]*x[n-M]
               r(y, a)[n] = - a[1]*y[n-1] - ... - a[N]*y[n-N]
  
      The IIR result is computed in parallel by first dividing the input signal
      into chunks (`g_i`) of size `m`. For each chunk, the IIR recurrence
      equation is applied to each chunk (in parallel). Then the chunks are merged
      based on the last N values of the last chunk::
  
               nc = L/m
               x = [g_0, g_1, ..., g_nc]
  
               g_i = [x[i * m], ..., x[i * m + m - 1]]
               p_i = r(g_i, a)
  
               o_i = r(p_i, c(p_{i - 1}))   if i > 1,
                     r(p_i, zi)             otherwise
  
               y = [o_0, o_1, ..., o_nc]
  
      where `c` denotes a function that takes a chunk, slices the last `N` values
      and adjust them using a correction factor table computed using the
      (1, 2, ..., N)-fibonacci sequence. For more information see [1]_.
  
      The rational transfer function describing this filter in the
      z-transform domain is::
  
                               -1              -M
                   b[0] + b[1]z  + ... + b[M] z
           Y(z) = -------------------------------- X(z)
                               -1              -N
                   a[0] + a[1]z  + ... + a[N] z
  
      References
      ----------
      .. [1] Sepideh Maleki and Martin Burtscher.
             2018. Automatic Hierarchical Parallelization of Linear Recurrences.
             SIGPLAN Not. 53, 2 (February 2018), 128-138.
             `10.1145/3173162.3173168 <https://doi.org/10.1145/3173162.3173168>`_
      """
      a0 = a[0]
      a_r = -a[1:] / a0
      b = b / a0
  
      num_b = b.size - 1
      num_a = a_r.size
      x_ndim = x.ndim
      axis = internal._normalize_axis_index(axis, x_ndim)
      n = x.shape[axis]
      fir_dtype = cupy.result_type(x, b)
  
      prev_in = None
      prev_out = None
      pad_shape = list(x.shape)
      pad_shape[axis] += num_b
  
      x_full = cupy.zeros(pad_shape, dtype=fir_dtype)
      if zi is not None:
          zi = cupy.atleast_1d(zi)
          if num_b > 0:
              prev_in = axis_slice(zi, 0, num_b, axis=axis)
          if num_a > 0:
              prev_out = axis_slice(
                  zi, zi.shape[axis] - num_a, zi.shape[axis], axis=axis
              )
  
      if prev_in is not None:
          x_full = axis_assign(x_full, prev_in, 0, num_b, axis=axis)
  
      x_full = axis_assign(x_full, x, num_b, axis=axis)
      origin = -num_b // 2
      out = cupy.empty_like(x_full, dtype=fir_dtype)
      out = _filters.convolve1d(
          x_full,
          b,
          axis=axis,
          mode="constant",
          origin=origin,
          output=out,
      )
  
      if num_b > 0:
          out = axis_slice(
              out, out.shape[axis] - n, out.shape[axis], axis=axis
          )
  
      if a_r.size > 0:
          iir_dtype = cupy.result_type(fir_dtype, a)
          const_dtype = cupy.dtype(a.dtype)
          if const_dtype.kind == "u":
              const_dtype = cupy.dtype(const_dtype.char.lower())
              a = a.astype(const_dtype)
  
          out = apply_iir(
              out, a_r, axis=axis, zi=prev_out, dtype=iir_dtype
          )
  
      if zi is not None:
          zi = cupy.empty(zi.shape, dtype=out.dtype)
          if num_b > 0:
              prev_in = axis_slice(
                  x, x.shape[axis] - num_b, x.shape[axis], axis=axis
              )
              zi = axis_assign(zi, prev_in, 0, num_b, axis=axis)
          if num_a > 0:
              prev_out = axis_slice(
                  out,
                  out.shape[axis] - num_a,
                  out.shape[axis],
                  axis=axis,
              )
              zi = axis_assign(
                  zi,
                  prev_out,
                  zi.shape[axis] - num_a,
                  zi.shape[axis],
                  axis=axis,
              )
          return out, zi
      else:
          return out
  
imports: import warnings
  import cupy
  from cupy._core import internal
  from cupy.linalg import lstsq
  import cupyx.scipy.fft
  from cupyx.scipy.ndimage import _util
  from cupyx.scipy.ndimage import _filters
  from cupyx.scipy.signal import _signaltools_core
  from cupyx.scipy.signal._arraytools import const_ext, even_ext, odd_ext, axis_reverse, axis_slice, axis_assign
  from cupyx.scipy.signal._iir_utils import apply_iir, apply_iir_sos, compute_correction_factors, compute_correction_factors_sos
label_desc: Use image node_id 15 for calling a global function with example usage: lfilter(b, a, x, axis, zi) and returns: out,  zi, out