{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Using Special Ops" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This notebook dives into some technical details of the ImageJ Ops library.\n", "\n", "It is recommended that you first read and understand the ImageJ Ops notebook before tackling this one." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Added new repo: imagej.public\n" ] }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "", "version_major": 2, "version_minor": 0 }, "method": "display_data" }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "2e223d23-b953-4e6f-8db7-3c5a0b5687ec", "version_major": 2, "version_minor": 0 }, "method": "display_data" }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "2.0.0-rc-71" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "%classpath config resolver imagej.public \"https://maven.imagej.net/content/groups/public\"\n", "%classpath add mvn net.imagej imagej 2.0.0-rc-71\n", "ij = new net.imagej.ImageJ()\n", "ij.getVersion()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This tutorial covers the so-called `SPECIAL` ops. Every special op produces one primary output from some fixed number of primary inputs. The number of such inputs for a given op is known as its `ARITY`. The framework provides interfaces for three arities:\n", "* 1 input = UNARY e.g.: eight = math.sqrt(64)\n", "* 2 inputs = BINARY e.g.: eight = math.add(3, 5)\n", "* 0 inputs = NULLARY e.g.: zero = math.zero()\n", "Additional arities (e.g., ternary) are feasible, but not implemented.\n", "\n", "The `SPECIAL` ops can be classified into three main categories: `COMPUTER`, `FUNCTION` and `INPLACE`:\n", "* `COMPUTER`: An op which will write its result to pre-allocated output, i.e. an empty (zero) image\n", "* `FUNCTION`: An op which will allocate on its own and then return it\n", "* `INPLACE`: An op which will write its output to one of the input arguments, hence mutating one of the inputs\n", "\n", "Special ops may also have any number of additional secondary inputs, which do not affect the arity of the op. More on that below.\n", "\n", "Let's try some example op calls, to illustrate the difference between these various kinds of special ops.\n", "\n", "First, we will create some test images to work with. In order to do that, we define and use a helper function, which creates an image of specified width and height using the `create.img` op, and then populates the contents of that image using the `image.equation` op, which sets the value of each sample according to a specified formula. In this case, we create four images depicting different patterns based on the location of a sample, which may be accessed via the variable `p[]` inside the formula." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Created image: Sinusoid of size [150,150]\n", "Created image: Gradient of size [150,150]\n", "Created image: Circle of size [150,150]\n", "Created image: Diamond of size [150,150]\n" ] }, { "data": { "text/html": [ "
sinusoidgradientspotcross
" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imglib2.type.numeric.real.DoubleType\n", "\n", "// Helper function to create DoubleType images from a formula\n", "imageFromFormula = { name, width, height, formula -> \n", " img = ij.op().run(\"create.img\", [width, height], new DoubleType())\n", " ij.op().run(\"image.equation\", img, formula)\n", " \n", " println(\"Created image: \" + name + \" of size [\" + width + \",\" + height + \"]\")\n", " \n", " return img\n", "}\n", "\n", "sinusoid = imageFromFormula(\"Sinusoid\", 150, 150, \"10 * (Math.cos(0.3*p[0]) + Math.sin(0.5*p[1]))\")\n", "gradient = imageFromFormula(\"Gradient\", 150, 150, \"p[0] + p[1]\")\n", "spot = imageFromFormula(\"Circle\", 150, 150, \"Math.sqrt((p[0]-75)*(p[0]-75) + (p[1]-75)*(p[1]-75))\")\n", "cross = imageFromFormula(\"Diamond\", 150, 150, \"Math.abs(p[0]-75) + Math.abs(p[1]-75)\")\n", "\n", "ij.notebook().display([[\"sinusoid\": sinusoid, \"gradient\": gradient, \"spot\": spot, \"cross\": cross]])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## BINARY OPS\n", "First we will look at binary ops. A `BINARY` op is one with two primary inputs, which produces an output. Suppose we want to add two images together. There are several versions of the `math.add` op. In the following example we will try the functional version:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "//Pointwise add FUNCTION\n", "functionOutput = ij.op().run(\"math.add\", sinusoid, gradient)\n", "ij.notebook().display(functionOutput)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The above invocation calls a `math.add` `FUNCTION`, which allocates and returns a newly created output image, without mutating any of the input data. This is great from a functional standpoint, but costs space and time. Now suppose we already have a properly dimensioned output image allocated and ready to go. In that case, it would be wasteful to call math.add as a function, only to then copy the data from the newly allocated result back into our already-existing image. Much better would be to tell `math.add` to use our existing output buffer directly. Fortunately, this approach is made possible by `COMPUTER` ops, which require you to specify a pre-allocated output buffer where the result will be stored. \n", "\n", "For example, imagine you want to do some analysis on a large set of images which would require to run e.g. a Gaussian smoothing as an intermediate operation. In this case, this intermediate step would not have to be store as it is not part of the end result. Thus, we can make use of the `COMPUTER` paradigm. We pre-allocate a single output buffer to store the results of the smoothing and continue with our analysis for the current image. For the next image, we simply overwrite the same buffer and continue. This approach is much faster compared to calculating the smoothing using a `FUNCTION` which would create a new output buffer for every intermediate step.\n", "\n", "In the following example, we will first pre-allocate an output image and then write the result of the addition to this image:\n" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "// Perallocate output image\n", "computerOutput = ij.op().run(\"create.img\", sinusoid)\n", "// Pointwise add using a COMPUTER\n", "ij.op()run(\"math.add\", computerOutput, spot, gradient)\n", "ij.notebook().display(computerOutput)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In the above call, the same operation is performed, but the output is written into the provided \"result\" buffer. `COMPUTER` ops do have an important restriction: the output buffer __MUST NOT__ be the same reference as any of the inputs. Otherwise, the computer would end up mutating the input data as a consequence of writing to the output buffer. Hence, computers cannot be used to calculate results \"in place\" -- that's what `INPLACE` ops are for.\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imagej.ops.special.inplace.Inplaces;\n", "import net.imagej.ops.Ops\n", "\n", "// Create a copy of one of the example images which we want to mutate\n", "gradientCopy = ij.op().run(\"copy.rai\", gradient)\n", "// Create an INPLACE math.add Op\n", "inplaceOp = Inplaces.binary1(ij.op(), Ops.Math.Add.class, gradientCopy, sinusoid)\n", "\"\"" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The maximum pixel value is: 308.6507787063433\n" ] }, { "data": { "text/html": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "inplaceOp.mutate(gradientCopy)\n", "max = ij.op()run(\"stats.max\", gradientCopy)\n", "println(\"The maximum pixel value is: \" + max)\n", "ij.notebook().display(gradientCopy)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The above invocation differs from earlier: the `Inplaces.binary1` method requests a `math.add` binary op which mutates its first argument inplace. A BinaryInplace1Op instance is returned, upon which we then call mutate(). In order to show that the content of the mutated image actually changes, we also calculate and print the maximum pixel value of the image using the `stats.max` op. If the above cell is executed more than once, one can see that the maximum of the image actually changes as values are added to the pixels. Note that we stored the `INPLACE` op in a variable. Such ops can be used in other ways too; see `REUSING SPECIAL OPS` below.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## UNARY and NULLARY OPS\n", "In addition to binary ops which have the `ARITY` 2, there are also `UNARY` (one primary input, one primary output) and `NULLARY` (no primary inputs, one primary output). Next, we are going to look at an example for both. First, we look at the `UNARY` `math.sqrt` op which simply calculates the square root of a number. Note, that this op is also a `COMPUTER`, meaning we first create an output container to write the result to." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "The sqrt of 64 is 8.0" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imglib2.type.numeric.real.DoubleType\n", "\n", "// Create preallocated output\n", "number = new DoubleType();\n", "// Calculate the square root of 64\n", "ij.op().run(\"math.sqrt\", number, new DoubleType(64))\n", "\n", "\"The sqrt of \" + 64 + \" is \" + number" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Second, we look at the `NULLARY` `math.zero` op. This op simply sets a number to zero. Note, that the variable `number` still contains the square root of the calculation from the cell above." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "number == 0.0? : true" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imglib2.type.numeric.real.DoubleType\n", "\n", "ij.op().run(\"math.zero\", number)\n", "\n", "\"number == 0.0? : \" + number.equals(new DoubleType(0))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Reusing Special Ops\n", "One of the most useful characteristics of special ops is the ability to look them up in a type-safe way, obtaining a reference to an op instance which can then be used however you need -- one time or many." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "net.imagej.ops.math.IIToIIOutputII$Add@3ce19b95" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imagej.ops.special.computer.Computers\n", "import net.imagej.ops.Ops\n", "\n", "// Create output for the COMPUTER\n", "computerOutput = ij.op().run(\"create.img\", sinusoid)\n", "// Get an instance of the \"math.add\" op and store it for later\n", "addOp = Computers.binary(ij.op(), Ops.Math.Add.class, computerOutput, gradient, cross)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this case, we use the `Computers.binary` utility method to require a `BINARY` `COMPUTER` `Op`, so that we can take advantage of the BinaryComputerOp API, which is type-safe -- whereas the Op interface's `run(\"op.name\", inputs..)` method is not.\n", "\n", "As usual, compile-time type safety is a two-edged sword: helpful to avoid coding errors, but sometimes rather verbose. As this notebook uses Groovy, which is a Scripting language, we can get around some of Java's verboseness. In java, the line from above creating the op would look like this:\n", "\n", "```java\n", "final BinaryComputerOp, IterableInterval, IterableInterval> addOp = Computers.binary(ops, Ops.Math.Add.class, computerOutput, gradient, cross);\n", "```\n", "However, the verbose way of doing it helps to understand what it going on underneath.\n", "\n", "In the cell below we can now just call the `run()` method to execute the op using the inputs specified above." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "addOp.run()\n", "ij.notebook().display(computerOutput)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Or we just run the op specifying different arguments by calling the `compute()` method. Note that for the `compute()` method the order of arguments will change. In this case the output to write to will be the last argument. Also note that the name of the method coincides with the type of op we are using. Here, we are using a `COMPUTER` op, hence the method is called `compute()`. In the case of an `INPLACE` op, the method will be called `mutate()` as one of the inputs is mutated inplace." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "// Call the op using different inputs\n", "addOp.compute(gradient, spot, computerOutput)\n", "ij.notebook().display(computerOutput)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "An op may additionally have secondary inputs. These are usually algorithm parameters or an `ENUM` that specifies the behaviour of the op. One example would be the sigma value or the out-of-bounds strategy used for a Gaussian filter op, which would look like this (sigma = 2, out-of-bounds = Boundary.SINGLE):\n", "\n", "```java\n", "gaussSigma2ComputerOp = Computers.unary(ij.op(), Ops.Filter.Gauss.class, input, input, 2, new OutOfBoundsMirrorFactory(Boundary.SINGLE))\n", "```\n", "Note the difference between primary and secondary inputs: the latter are not given when calling the special op's type-safe API; hence, secondary inputs remain fixed as given during the original matching (when we requested the op). Thus, the `UNARY` `COMPUTER` we saved in the variable `gaussSigma2ComputerOp` will always perform a Gaussian filtering with sigma 2. Example usage could look like this:\n", "\n", "```java\n", "gaussSigma2ComputerOp.compute(image, smoothedImage)\n", "```\n", "\n", "The variable `smoothedImage` will now contain the smoothed image using a sigma value of 2 as defined above when we requested the op.\n", "\n", "Which inputs are considered \"primary\" is a decision left to the author of each op. But here are two rules of thumb:\n", "* The primary inputs and outputs should correspond to the main relevant data structure -- e.g., ops which operate on images should declare the image argument(s) as primary.\n", "* There is usually an anticipated main use case involving variation of specific input(s) -- those should be chosen as the main inputs, so that the op can be reused across different values of those inputs. Hence, the inputs that change when executing the op several times are the primary ones, i.e. in most cases the input/output images." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Reuse an op to improve performance\n", "It is important to understand that each call into the Ops framework requires the system to search for the best-matching op amongst all those available. This is very useful for reasons of extensibility: the best op for the job is always selected, and anyone can extend the system with new ops optimized for particular scenarios. However, matching can have substantial performance implications:" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ij.op().run() needs\t2812 ms for 200 iterations\n", "addOp.compute() needs \t14 ms for 200 iterations\n" ] }, { "data": { "text/plain": [] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imglib2.type.numeric.real.DoubleType\n", "import net.imglib2.IterableInterval\n", "\n", "slowStart = System.currentTimeMillis();\n", "200.times {\n", " // Search for the op every time and then run it\n", " ij.op().run(\"math.add\", computerOutput, gradient, spot)\n", "}\n", "slowEnd = System.currentTimeMillis();\n", "println(\"ij.op().run() needs\\t\" + (slowEnd - slowStart) + \" ms for 200 iterations\")\n", "\n", "fastStart = System.currentTimeMillis();\n", "200.times {\n", " // Execute the same op instance repeatedly\n", " addOp.compute(gradient, spot, computerOutput)\n", "}\n", "fastEnd = System.currentTimeMillis();\n", "\n", "println(\"addOp.compute() needs \\t\" + (fastEnd - fastStart) + \" ms for 200 iterations\")\n", "\"\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As one can see, reusing the same matched op every time is much faster. However, it is the responsibility of the caller to ensure that subsequent arguments passed in are compatible with that same op." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Ops Chaining\n", "Some ops accept other ops as arguments, forming chains or trees of ops. For example, the map op executes the given op over all elements of an iteration. In the following example, we want to take the pixelwise square root of an image using the `math.sqrt` op. Unfortunately, this op is only applicable for singe numbers compared to the `math.add` op which can already handle `IterableInterval` for convenience. However, we can make use of the `map` op which expects another op as input and executes it on every element of an iterable; in this case every pixel. " ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The maximum pixel value:\n", "* before pointwise square root is:\t\t19.999118601072674\n", "* after pointwise square root should be:\t4.472037410518015\n", "* after pointwise square root is:\t\tNaN\n" ] }, { "data": { "text/plain": [] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imagej.ops.special.computer.Computers\n", "import net.imglib2.type.numeric.real.DoubleType\n", "import net.imagej.ops.Ops\n", "\n", "// Create output for the COMPUTER\n", "mapOutput = ij.op().run(\"create.img\", sinusoid)\n", "\n", "imgMax = ij.op().run(\"stats.max\", sinusoid)\n", "println(\"The maximum pixel value:\")\n", "println(\"* before pointwise square root is:\\t\\t\" + imgMax)\n", "println(\"* after pointwise square root should be:\\t\" + ij.op().run(\"math.sqrt\", imgMax.get()))\n", "\n", "// Create and save a sqrt op operating on DoubleType\n", "sqrtOp = Computers.unary(ij.op(), Ops.Math.Sqrt.class, DoubleType.class, DoubleType.class)\n", "\n", "// 'Map' the op on the 'sinusoid' image to calculate the square root of every pixel.\n", "ij.op().run(\"map\", mapOutput, sinusoid, sqrtOp);\n", "\n", "println(\"* after pointwise square root is:\\t\\t\" + ij.op().stats().max(mapOutput))\n", "\"\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In the above example, we additionally calculate some image statistics, namely the maximum pixel value in order to verify that the `map` op is actually doing what we'd expect. First we print the maximum of the `sinusoid` image and the square root of that value. After the map operation, the new maximum of the map output should be the same as previously calculated one, which is indeed the case :)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, we want to look at a slightly more advanced example of chaining. Here, we want to calculate the Difference of Gaussian of an image which is a classical image processing operation. If we look at the method signatures of the `dog` op, we see that there is a version with the following inputs:\n", "* RandomAccessibleInterval out? : the output to write to\n", "* RandomAccessibleInterval in : the input to calculate the dog on\n", "* UnaryComputerOp gauss1 : an op that calculates the first Gaussian filtered version of the image\n", "* UnaryComputerOp gauss2 : an op that calculates the second Gaussian filtered version of the image\n", "* UnaryFunctionOp outputCreator : an op to allocate output if the op is used as a FUNCTION\n", "* UnaryFunctionOp tmpCreator : an op to create a temporary storage needed by the algorithm in order to calculate the first Gaussian\n", "\n", "Therefore, we will create all of these ops and call the `dog` op using them:" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import net.imglib2.type.numeric.real.DoubleType\n", "import net.imagej.ops.special.computer.Computers\n", "import net.imagej.ops.Ops;\n", "import net.imglib2.RandomAccessibleInterval\n", "import java.util.Collections;\n", "import net.imagej.ops.special.function.Functions;\n", "import net.imglib2.outofbounds.OutOfBoundsMirrorFactory;\n", "import net.imglib2.outofbounds.OutOfBoundsMirrorFactory.Boundary;\n", "\n", "// For this example we will use a different image as this will illustrate the result of the `dog` more nicely\n", "input = ij.scifio().datasetIO().open(\"https://imagej.net/images/clown.png\")\n", "// Convert the input to float in order to get a corerct result\n", "input = ij.op().run(\"convert.float32\", input)\n", "\n", "Class raiClass = RandomAccessibleInterval.class\n", "// Create two Gaussian filter ops with different sigmas as secondary inputs. Furthermore, we need to specify an out-of-bounds strategy\n", "gaussSigma2Op = Computers.unary(ij.op(), Ops.Filter.Gauss.class, input, input, 2, new OutOfBoundsMirrorFactory(Boundary.SINGLE))\n", "gaussSigma8Op = Computers.unary(ij.op(), Ops.Filter.Gauss.class, input, input, 8, new OutOfBoundsMirrorFactory(Boundary.SINGLE))\n", "\n", "// Create the tmp and output creator ops\n", "outputCreatorFunc = Functions.unary(ij.op(), Ops.Create.Img.class, raiClass, input)\n", "tmpCreatorFunc = Functions.unary(ij.op(), Ops.Create.Img.class, raiClass, input)\n", "\n", "// Allocate output\n", "dogOutput = ij.op().run(\"create.img\", input)\n", "// Run the `dog` op\n", "ij.op().run(\"dog\", dogOutput, input, gaussSigma2Op, gaussSigma8Op, outputCreatorFunc, tmpCreatorFunc)\n", "\n", "ij.notebook().display(dogOutput)" ] } ], "metadata": { "kernelspec": { "display_name": "Groovy", "language": "groovy", "name": "groovy" }, "language_info": { "codemirror_mode": "groovy", "file_extension": ".groovy", "mimetype": "", "name": "Groovy", "nbconverter_exporter": "", "version": "2.4.3" }, "toc": { "base_numbering": 1, "nav_menu": { "height": "12px", "width": "214px" }, "number_sections": false, "sideBar": false, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": false, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }