diff --git a/docs/_posts/2019-10-02-rebuilding-finn-for-open-source.md b/docs/_posts/2019-10-02-rebuilding-finn-for-open-source.md
index 95381f289e5936057c86a67db9dd397054c30aad..324576ee46ecd044f75072976ae644d420e8fd59 100644
--- a/docs/_posts/2019-10-02-rebuilding-finn-for-open-source.md
+++ b/docs/_posts/2019-10-02-rebuilding-finn-for-open-source.md
@@ -11,21 +11,21 @@ to be more modular, more usable and more open-source!
 
 Over the past few years, the team at Xilinx Research Labs Ireland has done quite a bit of research of Quantized Neural Networks
 (QNNs).
-Starting with <a href="https://arxiv.org/abs/1612.07119">Binarized Neural Networks (BNNs) on FPGAs</a> back in 2016, we've since 
+Starting with <a href="https://arxiv.org/abs/1612.07119">Binarized Neural Networks (BNNs) on FPGAs</a> back in 2016, we've since
 looked at many aspects of quantized deep learning, ranging from
 at <a href ="https://arxiv.org/abs/1807.00301">better quantization methods</a> and
-<a href="https://arxiv.org/abs/1709.06262">mixing quantization and pruning</a>, 
+<a href="https://arxiv.org/abs/1709.06262">mixing quantization and pruning</a>,
 to <a href="https://arxiv.org/pdf/1807.10577.pdf">accuracy-throughput tradeoffs</a> and
 <a href="https://arxiv.org/pdf/1807.04093.pdf">recurrent topologies</a>.
 
-Although some <a href="https://github.com/Xilinx/BNN-PYNQ">demonstrators</a> of our work has been open source for some time, 
+Although some <a href="https://github.com/Xilinx/BNN-PYNQ">demonstrators</a> of our work has been open source for some time,
 we want to take things a step further.
-We love QNNs and the high-performance, high-efficiency dataflow accelerators we can build for them on Xilinx FPGAs, and we want you and 
+We love QNNs and the high-performance, high-efficiency dataflow accelerators we can build for them on Xilinx FPGAs, and we want you and
 the FPGA/ML community to be able to do the same.
 The (co-)design process for making this happen is actually quite involved, starting from customizing a neural network in a machine
-learning framework, going through multiple design steps that involve many optimizations, HLS code generation and Vivado synthesis, and 
+learning framework, going through multiple design steps that involve many optimizations, HLS code generation and Vivado synthesis, and
 ending up with an FPGA bitstream that you can deploy as part of some application.
-Many of those steps require some manual effort, but having a modular, flexible solution stack to support you through this process is greatly 
+Many of those steps require some manual effort, but having a modular, flexible solution stack to support you through this process is greatly
 helpful.
 This is why we are rebulding our FINN solution stack from the ground-up to make it more modular, and we hope to build a community
 around it that shares our excitement around QNNs for FPGAs.
@@ -40,31 +40,31 @@ frameworks like <a href="http://llvm.org">LLVM</a>.
 This stack breaks down the complex co-design problem into parts, and each layer focuses on a different sub-problem, consuming
 the artifacts produced by the previous one.
 The diagram on the left illustrates this briefly, and over the next few months we hope to make a first few QNNs go through all
-the layers of this stack to produce cool FPGA dataflow accelerators. 
+the layers of this stack to produce cool FPGA dataflow accelerators.
 In fact, some of these components are already available today for you to explore!
 
 Let's have a look at the main parts:
 
 * <b>Brevitas</b> is a PyTorch library that lets you do quantization-aware training. It gives you a set of `torch.nn` building
-blocks to explore different forms of weight, activation and accumulator quantization schemes. You can also learn the bitwidths for 
+blocks to explore different forms of weight, activation and accumulator quantization schemes. You can also learn the bitwidths for
 different layers with backpropagation! See the <a href="https://xilinx.github.io/brevitas/">Brevitas page</a> for more information.
-* <b>Frontend</b>. Once you are happy with the accuracy of your quantized neural network in Brevitas, you'll be able to export it into a custom 
+* <b>Frontend</b>. Once you are happy with the accuracy of your quantized neural network in Brevitas, you'll be able to export it into a custom
 <a href="https://onnx.ai">ONNX</a> representation that FINN uses internally to represent QNNs. More details about this custom ONNX
 representation will be available in an upcoming blog post.
 * The <b>FINN Compiler</b> will then import this ONNX representation, and go through several steps of optimizations such as the
-  <a href="https://arxiv.org/pdf/1709.04060.pdf">streamlining transform</a> to make the QNN simpler. 
+  <a href="https://arxiv.org/pdf/1709.04060.pdf">streamlining transform</a> to make the QNN simpler.
 * The <b>FPGA dataflow backend</b> will then convert the optimized QNN into a series of streaming HLS library calls. An important
- part of the stack is the <a href="https://github.com/Xilinx/finn-hlslib">FINN HLS library</a>, which provides optimized Vivado HLS 
+ part of the stack is the <a href="https://github.com/Xilinx/finn-hlslib">FINN HLS library</a>, which provides optimized Vivado HLS
  descriptions of several common layer types (convolutions, thresholding, pooling...) found in QNNs.
  * <b>Synthesis</b>. Once the HLS calls are generated, the next steps are to call Vivado HLS and Vivado to generate a bitstream for the target
  Xilinx FPGA. We have plans to support Vivado IPI block design code generation as well for increased agility and modularity.
  * <b>PYNQ deployment</b>. Finally, you will be able to use any of the supported <a href="http://www.pynq.io/">PYNQ</a> platforms to directly call the
  generated accelerator from Python and integrate it with other functionality. Since FINN-generated dataflow accelerators expose
- streaming interfaces, we think it will be exciting to use streaming-oriented Python frameworks such as 
+ streaming interfaces, we think it will be exciting to use streaming-oriented Python frameworks such as
  <a href="https://github.com/ray-project/ray">Ray</a> to create heterogeneous, high-performance task graphs incorporating QNNs.
- 
+
  ### Getting started
- 
+
  More will be available in the coming weeks and months, but if you want to get your hands dirty there's already plenty to start with!
  If you haven't done so already, we recommend starting with <a href="https://github.com/Xilinx/BNN-PYNQ">BNN-PYNQ</a> to see what
  dataflow QNN accelerators look and feel like.
@@ -72,4 +72,3 @@ representation will be available in an upcoming blog post.
  put together a streaming pipeline with the <a href="https://github.com/Xilinx/finn-hlslib">FINN HLS library</a>.
  We have also created a <a href="https://gitter.im/xilinx-finn/community">Gitter channel</a> to make it easier to get in touch with
  the community, and hope to see many of you there! :)
- 
diff --git a/docs/index.md b/docs/index.md
index dc6c5033233e970b2a8032978558b06ab5fd92f3..ac635db3474f9deb8cca13da84adb4f79f50b222 100644
--- a/docs/index.md
+++ b/docs/index.md
@@ -21,7 +21,7 @@ Depending on what you would like to do, we have different suggestions on where t
 * **I want to try out prebuilt QNN accelerators on real hardware.** Head over to <a href="https://github.com/Xilinx/BNN-PYNQ" target="_blank">BNN-PYNQ</a> repository to try out some image
 classification accelerators, or to <a href="https://github.com/Xilinx/LSTM-PYNQ" target="_blank">LSTM-PYNQ</a>
 to try optical character recognition with LSTMs.
-* **I want to train new quantized networks for FINN.** Check out <a href="https://github.com/Xilinx/brevitas">Brevitas</a>, 
+* **I want to train new quantized networks for FINN.** Check out <a href="https://github.com/Xilinx/brevitas">Brevitas</a>,
 our PyTorch library for training quantized networks. The Brevitas-to-FINN part of the flow is coming soon!
 * **I want to understand the computations involved in quantized inference.** Check out these Jupyter notebooks on <a href="https://github.com/maltanar/qnn-inference-examples">QNN inference</a>. This repo contains simple Numpy/Python layer implementations and a few pretrained QNNs for instructive purposes.
 * **I want to understand how it all fits together.** Check out our [publications](#publications),
diff --git a/src/finn/core/execute_custom_node.py b/src/finn/core/execute_custom_node.py
new file mode 100644
index 0000000000000000000000000000000000000000..c941b5011bbae20cef7f751351dd581086d05f9d
--- /dev/null
+++ b/src/finn/core/execute_custom_node.py
@@ -0,0 +1,25 @@
+# import onnx.helper as helper
+
+import finn.core.multithreshold as multiThresh
+
+
+def execute_custom_node(node, context, graph):
+    """Call custom implementation to execute a single custom node.
+    Input/output provided via context."""
+
+    if node.op_type == "MultiThreshold":
+
+        # save inputs
+        v = context[node.input[0]]
+        thresholds = context[node.input[1]]
+
+        # calculate output
+        output = multiThresh.execute(v, thresholds)
+
+        # setting context according to output
+        context[node.output[0]] = output
+        print(output)
+
+    else:
+        # exception if op_type is not supported
+        raise Exception("This custom node is currently not supported.")
diff --git a/src/finn/core/multithreshold.py b/src/finn/core/multithreshold.py
new file mode 100755
index 0000000000000000000000000000000000000000..23b5cca5a137ad9d6be5e06436ffdef212e7828c
--- /dev/null
+++ b/src/finn/core/multithreshold.py
@@ -0,0 +1,56 @@
+import numpy as np
+
+
+def compare(x, y):
+    if x >= y:
+        return 1.0
+    else:
+        return 0.0
+
+
+def execute(v, thresholds):
+
+    # the inputs are expected to be in the shape (N,C,H,W)
+    # N : Batch size
+    # C : Number of channels
+    # H : Heigth of the input images
+    # W : Width of the input images
+    #
+    # the thresholds are expected to be in the shape (C, B)
+    # C : Number of channels (must be the same value as C in input tensor)
+    # B : Desired activation steps => i.e. for 4-bit activation, B=7 (2^(n)-1 and n=4)
+
+    # assert if channel sizes do not match
+    assert v.shape[1] == thresholds.shape[0]
+
+    # save the required shape sizes for the loops (N, C and B)
+    num_batch = v.shape[0]
+    num_channel = v.shape[1]
+
+    num_act = thresholds.shape[1]
+
+    # reshape inputs to enable channel-wise reading
+    vr = v.reshape((v.shape[0], v.shape[1], -1))
+
+    # save the new shape size of the images
+    num_img_elem = vr.shape[2]
+
+    # initiate output tensor
+    ret = np.zeros_like(vr)
+
+    # iterate over thresholds channel-wise
+    for t in range(num_channel):
+        channel_thresh = thresholds[t]
+
+        # iterate over batches
+        for b in range(num_batch):
+
+            # iterate over image elements on which the thresholds should be applied
+            for elem in range(num_img_elem):
+
+                # iterate over the different thresholds that correspond to one channel
+                for a in range(num_act):
+                    # apply successive thresholding to every element of one channel
+                    ret[b][t][elem] += compare(vr[b][t][elem], channel_thresh[a])
+
+    return ret.reshape(v.shape)
diff --git a/src/finn/core/onnx_exec.py b/src/finn/core/onnx_exec.py
index 5521bccbdc277bd87f52ad6f27b0b6207c39f70f..0e4f0237d64331a02cf6b85d5d0ea07645957abe 100644
--- a/src/finn/core/onnx_exec.py
+++ b/src/finn/core/onnx_exec.py
@@ -29,34 +29,53 @@ import copy
 import onnx.helper as helper
 import onnxruntime as rt
 
+import finn.core.execute_custom_node as ex_cu_node
+
 
 def execute_node(node, context, graph):
     """Call onnxruntime to execute a single node. Input/output provided via context."""
 
-    # onnxruntime unfortunately does not implement run_node as defined by ONNX,
-    # it can only execute entire models -- so we create a model which solely
-    # consists of our current node.
-    node_inputs = list(filter(lambda x: x.name in node.input, graph.input))
-    node_inputs += list(filter(lambda x: x.name in node.input, graph.value_info))
-    node_outputs = list(filter(lambda x: x.name in node.output, graph.output))
-    node_outputs += list(filter(lambda x: x.name in node.output, graph.value_info))
-    node_graph = helper.make_graph(
-        nodes=[node], name="single-node-exec", inputs=node_inputs, outputs=node_outputs
-    )
-    node_model = helper.make_model(node_graph)
-    input_dict = dict()
-    for inp in node.input:
-        input_dict[inp] = context[inp]
-    sess = rt.InferenceSession(node_model.SerializeToString())
-    output_list = sess.run(None, input_dict)
-    for output_ind in range(len(node.output)):
-        outp = node.output[output_ind]
-        if output_list[output_ind].shape != context[outp].shape:
-            raise Exception(
-                "Output shapes disagree after node execution: found %s vs expected %s"
-                % (str(output_list[output_ind].shape.shape), str(context[outp].shape))
-            )
-        context[outp] = output_list[output_ind]
+    # run node with custom function or by using onnxruntime
+
+    if node.domain == "finn":
+
+        ex_cu_node.execute_custom_node(node, context, graph)
+
+    else:
+
+        # onnxruntime unfortunately does not implement run_node as defined by ONNX,
+        # it can only execute entire models -- so we create a model which solely
+        # consists of our current node.
+        node_inputs = list(filter(lambda x: x.name in node.input, graph.input))
+        node_inputs += list(filter(lambda x: x.name in node.input, graph.value_info))
+        node_outputs = list(filter(lambda x: x.name in node.output, graph.output))
+        node_outputs += list(filter(lambda x: x.name in node.output, graph.value_info))
+        node_graph = helper.make_graph(
+            nodes=[node],
+            name="single-node-exec",
+            inputs=node_inputs,
+            outputs=node_outputs,
+        )
+        node_model = helper.make_model(node_graph)
+        input_dict = dict()
+        for inp in node.input:
+            input_dict[inp] = context[inp]
+
+        sess = rt.InferenceSession(node_model.SerializeToString())
+        output_list = sess.run(None, input_dict)
+
+        for output_ind in range(len(node.output)):
+            outp = node.output[output_ind]
+            if output_list[output_ind].shape != context[outp].shape:
+                raise Exception(
+                    """Output shapes disagree after node execution:
+                    found %s vs expected %s"""
+                    % (
+                        str(output_list[output_ind].shape.shape),
+                        str(context[outp].shape),
+                    )
+                )
+            context[outp] = output_list[output_ind]
 
 
 def execute_onnx(model, input_dict, return_full_exec_context=False):
diff --git a/tests/test_custom_onnx_exec.py b/tests/test_custom_onnx_exec.py
new file mode 100644
index 0000000000000000000000000000000000000000..0d07b9888d1330753446d589619149c1f8f316cf
--- /dev/null
+++ b/tests/test_custom_onnx_exec.py
@@ -0,0 +1,213 @@
+import numpy as np
+from onnx import TensorProto, helper
+
+import finn.core.execute_custom_node as ex_cu_node
+
+
+def test_execute_custom_node_multithreshold():
+    inputs = np.ndarray(
+        shape=(6, 3, 2, 2),
+        buffer=np.array(
+            [
+                4.8,
+                3.2,
+                1.2,
+                4.9,
+                7.8,
+                2.4,
+                3.1,
+                4.7,
+                6.2,
+                5.1,
+                4.9,
+                2.2,
+                6.2,
+                0.0,
+                0.8,
+                4.7,
+                0.2,
+                5.6,
+                8.9,
+                9.2,
+                9.1,
+                4.0,
+                3.3,
+                4.9,
+                2.3,
+                1.7,
+                1.3,
+                2.2,
+                4.6,
+                3.4,
+                3.7,
+                9.8,
+                4.7,
+                4.9,
+                2.8,
+                2.7,
+                8.3,
+                6.7,
+                4.2,
+                7.1,
+                2.8,
+                3.1,
+                0.8,
+                0.6,
+                4.4,
+                2.7,
+                6.3,
+                6.1,
+                1.4,
+                5.3,
+                2.3,
+                1.9,
+                4.7,
+                8.1,
+                9.3,
+                3.7,
+                2.7,
+                5.1,
+                4.2,
+                1.8,
+                4.1,
+                7.3,
+                7.1,
+                0.4,
+                0.2,
+                1.3,
+                4.3,
+                8.9,
+                1.4,
+                1.6,
+                8.3,
+                9.4,
+            ]
+        ),
+    )
+
+    threshold_values = np.ndarray(
+        shape=(3, 7),
+        buffer=np.array(
+            [
+                0.8,
+                1.4,
+                1.7,
+                3.5,
+                5.2,
+                6.8,
+                8.2,
+                0.2,
+                2.2,
+                3.5,
+                4.5,
+                6.6,
+                8.6,
+                9.2,
+                1.3,
+                4.1,
+                4.5,
+                6.5,
+                7.8,
+                8.1,
+                8.9,
+            ]
+        ),
+    )
+
+    v = helper.make_tensor_value_info("v", TensorProto.FLOAT, [6, 3, 2, 2])
+    thresholds = helper.make_tensor_value_info("thresholds", TensorProto.FLOAT, [3, 7])
+    out = helper.make_tensor_value_info("out", TensorProto.FLOAT, [6, 3, 2, 2])
+
+    node_def = helper.make_node(
+        "MultiThreshold", ["v", "thresholds"], ["out"], domain="finn"
+    )
+
+    graph_def = helper.make_graph([node_def], "test_model", [v, thresholds], [out])
+
+    execution_context = {}
+    execution_context["v"] = inputs
+    execution_context["thresholds"] = threshold_values
+
+    ex_cu_node.execute_custom_node(node_def, execution_context, graph_def)
+
+    outputs = np.ndarray(
+        shape=(6, 3, 2, 2),
+        buffer=np.array(
+            [
+                4.0,
+                3.0,
+                1.0,
+                4.0,
+                5.0,
+                2.0,
+                2.0,
+                4.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                5.0,
+                0.0,
+                1.0,
+                4.0,
+                1.0,
+                4.0,
+                6.0,
+                7.0,
+                7.0,
+                1.0,
+                1.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                3.0,
+                4.0,
+                2.0,
+                3.0,
+                7.0,
+                3.0,
+                3.0,
+                1.0,
+                1.0,
+                7.0,
+                5.0,
+                4.0,
+                6.0,
+                2.0,
+                2.0,
+                1.0,
+                1.0,
+                2.0,
+                1.0,
+                3.0,
+                3.0,
+                2.0,
+                5.0,
+                3.0,
+                3.0,
+                4.0,
+                5.0,
+                7.0,
+                3.0,
+                1.0,
+                3.0,
+                2.0,
+                1.0,
+                4.0,
+                6.0,
+                6.0,
+                0.0,
+                1.0,
+                1.0,
+                3.0,
+                6.0,
+                1.0,
+                1.0,
+                6.0,
+                7.0,
+            ]
+        ),
+    )
+
+    assert (execution_context["out"] == outputs).all()
diff --git a/tests/test_mixed_onnx_exec.py b/tests/test_mixed_onnx_exec.py
new file mode 100644
index 0000000000000000000000000000000000000000..912d62b5cb3962b881f5267c3c425231f0f8ee32
--- /dev/null
+++ b/tests/test_mixed_onnx_exec.py
@@ -0,0 +1,224 @@
+import numpy as np
+import onnx
+from onnx import TensorProto, helper
+
+import finn.transformation.infer_shapes as si
+from finn.core.modelwrapper import ModelWrapper
+import finn.core.onnx_exec as oxe
+
+
+def test_execute_mixed_model():
+
+    out0 = helper.make_tensor_value_info("out0", TensorProto.FLOAT, [6, 3, 2, 2])
+
+    graph_def = helper.make_graph(
+        nodes=[
+            helper.make_node(
+                "MultiThreshold", ["v", "thresholds"], ["out0"], domain="finn"
+            ),
+            helper.make_node("Relu", ["out0"], ["out1"]),
+        ],
+        name="test-model",
+        inputs=[
+            helper.make_tensor_value_info("v", TensorProto.FLOAT, [6, 3, 2, 2]),
+            helper.make_tensor_value_info("thresholds", TensorProto.FLOAT, [3, 7]),
+        ],
+        outputs=[helper.make_tensor_value_info("out1", TensorProto.FLOAT, [6, 3, 2, 2])],
+        value_info = [out0],
+    )
+    model_def = helper.make_model(graph_def, producer_name='onnx-example')
+
+    model = ModelWrapper(model_def)
+    model = model.transform_single(si.infer_shapes)
+
+    inputs = np.asarray(
+            [
+                4.8,
+                3.2,
+                1.2,
+                4.9,
+                7.8,
+                2.4,
+                3.1,
+                4.7,
+                6.2,
+                5.1,
+                4.9,
+                2.2,
+                6.2,
+                0.0,
+                0.8,
+                4.7,
+                0.2,
+                5.6,
+                8.9,
+                9.2,
+                9.1,
+                4.0,
+                3.3,
+                4.9,
+                2.3,
+                1.7,
+                1.3,
+                2.2,
+                4.6,
+                3.4,
+                3.7,
+                9.8,
+                4.7,
+                4.9,
+                2.8,
+                2.7,
+                8.3,
+                6.7,
+                4.2,
+                7.1,
+                2.8,
+                3.1,
+                0.8,
+                0.6,
+                4.4,
+                2.7,
+                6.3,
+                6.1,
+                1.4,
+                5.3,
+                2.3,
+                1.9,
+                4.7,
+                8.1,
+                9.3,
+                3.7,
+                2.7,
+                5.1,
+                4.2,
+                1.8,
+                4.1,
+                7.3,
+                7.1,
+                0.4,
+                0.2,
+                1.3,
+                4.3,
+                8.9,
+                1.4,
+                1.6,
+                8.3,
+                9.4,
+            ]
+        , dtype=np.float32,
+    ).reshape(6, 3, 2, 2)
+
+    threshold_values = np.asarray(
+            [
+                0.8,
+                1.4,
+                1.7,
+                3.5,
+                5.2,
+                6.8,
+                8.2,
+                0.2,
+                2.2,
+                3.5,
+                4.5,
+                6.6,
+                8.6,
+                9.2,
+                1.3,
+                4.1,
+                4.5,
+                6.5,
+                7.8,
+                8.1,
+                8.9,
+            ]
+        , dtype=np.float32,
+    ).reshape(3, 7)
+
+    input_dict = {}
+    input_dict['v']=inputs
+    input_dict['thresholds']=threshold_values
+
+    output_dict = oxe.execute_onnx(model, input_dict)
+
+    outputs = np.asarray(
+            [
+                4.0,
+                3.0,
+                1.0,
+                4.0,
+                5.0,
+                2.0,
+                2.0,
+                4.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                5.0,
+                0.0,
+                1.0,
+                4.0,
+                1.0,
+                4.0,
+                6.0,
+                7.0,
+                7.0,
+                1.0,
+                1.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                3.0,
+                4.0,
+                2.0,
+                3.0,
+                7.0,
+                3.0,
+                3.0,
+                1.0,
+                1.0,
+                7.0,
+                5.0,
+                4.0,
+                6.0,
+                2.0,
+                2.0,
+                1.0,
+                1.0,
+                2.0,
+                1.0,
+                3.0,
+                3.0,
+                2.0,
+                5.0,
+                3.0,
+                3.0,
+                4.0,
+                5.0,
+                7.0,
+                3.0,
+                1.0,
+                3.0,
+                2.0,
+                1.0,
+                4.0,
+                6.0,
+                6.0,
+                0.0,
+                1.0,
+                1.0,
+                3.0,
+                6.0,
+                1.0,
+                1.0,
+                6.0,
+                7.0,
+            ]
+        , dtype=np.float32,
+    ).reshape(6, 3, 2, 2)
+
+    assert (output_dict["out1"] == outputs).all()
+ 
diff --git a/tests/test_multi_thresholding.py b/tests/test_multi_thresholding.py
new file mode 100644
index 0000000000000000000000000000000000000000..49305e5572f35eb4a2e5f7678c73038777eb8b92
--- /dev/null
+++ b/tests/test_multi_thresholding.py
@@ -0,0 +1,199 @@
+import numpy as np
+
+import finn.core.multithreshold as multi_thresh
+
+
+def test_execute_multi_thresholding():
+
+    inputs = np.ndarray(
+        shape=(6, 3, 2, 2),
+        buffer=np.array(
+            [
+                4.8,
+                3.2,
+                1.2,
+                4.9,
+                7.8,
+                2.4,
+                3.1,
+                4.7,
+                6.2,
+                5.1,
+                4.9,
+                2.2,
+                6.2,
+                0.0,
+                0.8,
+                4.7,
+                0.2,
+                5.6,
+                8.9,
+                9.2,
+                9.1,
+                4.0,
+                3.3,
+                4.9,
+                2.3,
+                1.7,
+                1.3,
+                2.2,
+                4.6,
+                3.4,
+                3.7,
+                9.8,
+                4.7,
+                4.9,
+                2.8,
+                2.7,
+                8.3,
+                6.7,
+                4.2,
+                7.1,
+                2.8,
+                3.1,
+                0.8,
+                0.6,
+                4.4,
+                2.7,
+                6.3,
+                6.1,
+                1.4,
+                5.3,
+                2.3,
+                1.9,
+                4.7,
+                8.1,
+                9.3,
+                3.7,
+                2.7,
+                5.1,
+                4.2,
+                1.8,
+                4.1,
+                7.3,
+                7.1,
+                0.4,
+                0.2,
+                1.3,
+                4.3,
+                8.9,
+                1.4,
+                1.6,
+                8.3,
+                9.4,
+            ]
+        ),
+    )
+
+    thresholds = np.ndarray(
+        shape=(3, 7),
+        buffer=np.array(
+            [
+                0.8,
+                1.4,
+                1.7,
+                3.5,
+                5.2,
+                6.8,
+                8.2,
+                0.2,
+                2.2,
+                3.5,
+                4.5,
+                6.6,
+                8.6,
+                9.2,
+                1.3,
+                4.1,
+                4.5,
+                6.5,
+                7.8,
+                8.1,
+                8.9,
+            ]
+        ),
+    )
+
+    outputs = np.ndarray(
+        shape=(6, 3, 2, 2),
+        buffer=np.array(
+            [
+                4.0,
+                3.0,
+                1.0,
+                4.0,
+                5.0,
+                2.0,
+                2.0,
+                4.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                5.0,
+                0.0,
+                1.0,
+                4.0,
+                1.0,
+                4.0,
+                6.0,
+                7.0,
+                7.0,
+                1.0,
+                1.0,
+                3.0,
+                3.0,
+                3.0,
+                1.0,
+                3.0,
+                4.0,
+                2.0,
+                3.0,
+                7.0,
+                3.0,
+                3.0,
+                1.0,
+                1.0,
+                7.0,
+                5.0,
+                4.0,
+                6.0,
+                2.0,
+                2.0,
+                1.0,
+                1.0,
+                2.0,
+                1.0,
+                3.0,
+                3.0,
+                2.0,
+                5.0,
+                3.0,
+                3.0,
+                4.0,
+                5.0,
+                7.0,
+                3.0,
+                1.0,
+                3.0,
+                2.0,
+                1.0,
+                4.0,
+                6.0,
+                6.0,
+                0.0,
+                1.0,
+                1.0,
+                3.0,
+                6.0,
+                1.0,
+                1.0,
+                6.0,
+                7.0,
+            ]
+        ),
+    )
+
+    results = multi_thresh.execute(inputs, thresholds)
+
+    assert (results == outputs).all()