n = self.get_max_ir_len()

        if fs is not None:
            x = np.arange(0, n) / fs
        else:
            x = np.arange(0, n)

        impulse_response = np.zeros(n)
        impulse_response[0] = 1

        M = len(self.b)
        N = len(self.a)

        # Apply the filter operation
        for i in range(n):
            y = 0
            for k in range(M):
                if i - k >= 0:
                    y += self.b[k] * x[i - k]

            for k in range(N):
                if i - k >= 0:
                    y += self.a[k] * impulse_response[i - k]
            impulse_response[i] = y
        return x, impulse_response