PixelMath Reference

By Juan Conejero (PTeam)

Introduction

PixelMath is part of the standard set of process modules of the PixInsight platform. A first version of this essential tool was introduced along with build 207 of the core PixInsight application, and has experienced important additions and improvements since then.

This work includes a reference documentation for PixelMath's expression syntax, a description of PixelMath's features and capabilities, and a number of examples. The PixelMath process is available on the Process Explorer window, under the PixelMath category, and installs on the Favorite Processes folder by default. It can also be found at Process > PixelMath > PixelMath in the core application's main menu.

PixelMath features

These are the main characteristics of the PixelMath process in PixInsight:

High-performance, proprietary expression parser/interpreter. PixelMath expressions are plain text strings in common algebraic notation, following a set of syntax rules very similar to the basic C syntax, with extensions and some new operators adapted to manage 2-D images efficiently.

Per-channel expressions. Three different expressions can be defined and applied to each individual RGB channel, or a single expression can be used for all nominal channels of the target image (default setting). This is a user-selectable option.

Variables. User-defined symbols can be declared as temporary variables. This allows to store intermediate results, and is especially useful when combined with multiple statements (see next point).

Multiple statements. A PixelMath expression can include multiple statements to perform temporary calculations that can be assigned to variables. This greatly improves execution performance and facilitates writing very complex expressions.

Parallelized process. PixelMath can use all available processors on multicore and multiprocessor systems.

Optional automatic rescaling of the result to an arbitrary output range.

49 functions available (as of writing this document), including statistical, algebraic and transcendental functions, as well as pixel access functions in rectangular and polar coordinates.

Optional usage of a 64-bit floating point temporary image to accumulate intermediate results. By default, a 32-bit floating point temporary image is used.

Recursive detection and precalculation of invariant subexpressions. For example, in the expression Min(img)*x, the function call Min(img) is pre-evaluated as an invariant subexpression equal to the minimum sample value of an image img (see Section 2.6). All invariants are precalculated prior expression execution.

Operand images can have arbitrary dimensions and/or color spaces. If images with different sizes are involved in a PixelMath execution, interpolation is carried out on the fly as appropriate to generate the required pixel values. If different color spaces are involved (e.g. RGB and grayscale), the process automatically manages grayscale images as RGB images with three identical channels.

Metasymbols allow defining PixelMath expressions that are independent on a specific target image. This allows to write, test and refine a PixelMath expression once, and then apply it to any number of images without modifications.

Optional generation of a new image with the result of the PixelMath calculation. For newly generated images, a complete set of parameters can be specified: dimensions, color space, sample format and image identifier. Optionally, newly generated images can adopt the same properties of the target image, if available.

Maskable process. PixelMath can be applied through mask images without limitations.

Previewable process. PixelMath can be tried out on previews, even if generation of a new image has been selected to store the result.

Expression editor interface with pre-evaluation of expressions for syntax check and diagnostics of potential errors.

Optional execution in the global context (to generate new images), along with normal execution on a target view (to replace its contents).

1. How PixelMath Works

Target image

PixelMath expressions are always executed on a target image. The target image can be either an existing image to which the PixelMath instance is applied, or a newly created image when PixelMath is executed in the global context. In both, cases, the target image defines the range of pixel coordinates where PixelMath runs.

Current pixel

PixelMath expressions are always executed pixel by pixel. It is a common mistake to think that a PixelMath expression can be applied to an image as a whole. Whereas some PixelMath functions calculate values from entire images (extreme pixel values and standard deviation, for example), a PixelMath expression evaluates to a unique value that is always assigned to a single pixel of the target image.

The following pseudocode describes how PixelMath works with the simple expression A+B:

for c = 0 to number_of_nominal_channels_in_target - 1 for y = 0 to height_of_target - 1 for x = 0 to width_of_target - 1 a = A[c,x,y] b = B[c,x,y] target[c,x,y] = a + b end_for end_for end_for

where the expression img[c,x,y] refers to the sample value at the x,y pixel coordinates of the channel c in an image img. In the pseudocode above, target refers to the target image, as we have described it in the preceding section.

When we talk about the current pixel or current sample, we refer to the pixel or sample corresponding to the current channel and pixel coordinates on the target image (perhaps mapped from a mismatching image geometry; see the next point).

Pixel sample values are always expressed as real numbers in the normalized range [0,1] (0=black, 1=white) for all supported data types: 8, 16 and 32-bit integers, 32 and 64-bit floating point. However, intermediate and temporary values, as well as variables, can have any real numerical value during PixelMath execution.

Automatic adaptation between different image geometries and color spaces

In the example above, what happens if the A, B and target images have different dimensions in pixels? In these cases, PixelMath automatically interpolates sample values by mapping the actual range of coordinates of each operand image to that of the target image. In other words, when image dimensions differ, PixelMath maps all operand images, as necessary, to the coordinate range of the target image.

And how about images in different color spaces? This poses two possible scenarios:

The target image is a RGB color image, and an operand is a grayscale image. In this case, the grayscale operand image is treated as a RGB color image with identical red, green and blue channels. For each nominal channel of the target image, the same operand values are applied.

The target image is a grayscale image, and an operand is a RGB color image. In this special case, the red channel of the operand RGB image is applied to the nominal channel of the grayscale target image, and the green and blue channels of the RGB operand are ignored unless explicitly referenced (see section 2.3 about channel references).

Invariant subexpressions

There are many PixelMath functions and subexpressions that can be precalculated, that is, evaluated only once before applying a PixelMath expression to a target image. This happens when such subexpressions remain constant during PixelMath execution. For example, in the expression:

A + 2 * 3

The subexpression "2*3" is an invariant. Obviously, the actual expression that PixelMath would execute is:

A + 6

This has been a pretty straightforward example. Here's a more interesting one:

A + StdDev( B )/Med( A )

In the expression above, the StdDev and Med functions are being used, respectively, to calculate the standard deviation and median values of the images B and A. Since these values won't change during PixelMath execution, they can be precalculated and stored in temporary variables. Furthermore, the quotient between StdDev( B ) and Med( A ) becomes also an invariant once both functions have been pre-evaluated, so it yields a single value that remains constant during PixelMath execution. The final expression that PixelMath would apply is:

A + k

where k would be precalculated as k = StdDev( B )/Med( A ). Recursive precalculation of invariant subexpressions is a feature of PixelMath's interpreter that greatly improves its performance, especially for complex expressions. In the two examples that we have seen here, the complexities of the applied expressions are identical, even if the second original expression is much more complicated than the first one.

Automatic rescaling and truncation

The PixelMath process can optionally rescale the resulting target image to a user-defined output range. This feature can be controlled with the standard PixelMath interface by means of the Rescale result check box and the two bound fields (lower and upper bounds). By default, rescaling is enabled and the output range is the standard normalized range [0,1].

In general, if v0 and v1 are, respectively, the lower and upper bounds of the output range, then for each sample value v of the resulting target image, we have:

v' = v0 + (v - m)×(v1 - v0)/(M - m),

where m and M are, respectively, the resulting minimum and maximum sample values in the resulting image, and v' is the final rescaled sample value.

It is important to point out that the m and M extreme sample values in the expression above, are calculated for the whole resulting target image. This means that the automatic rescaling operation does not change in any way the relative values of the nominal channels in the target image (the relations between red, green and blue) that result after PixelMath expression evaluation.

When the automatic rescaling feature is disabled, PixelMath truncates the resulting values to the [0,1] range.

2. PixelMath Syntax

This section describes the main elements of PixelMath expressions and their syntax rules. A valid PixelMath expression is a plain text string such that:

Follows PixelMath syntax rules, as described in this document.

Represents a set of functional, arithmetic, logical or bitwise logical operations performed between images, pixel sample values, and constant values.

Evaluates to a single sample value that is necessarily assigned to the current sample of the target image during PixelMath execution (see Section 1 for definitions).

In this document, we use a practical, mostly informal approach to describe PixelMath syntax, instead of a formal, rigorous language specification. This is because PixelMath syntax is very close to common algebraic notation, as used by most programmable calculators for example, or in programming languages such as C, FORTRAN or BASIC. We think that most users will have at least a basic knowledge of that notation, so understanding PixelMath expressions shouldn't be difficult.

Employed metalanguage

Even if we are not doing a formal language description, we need at least a basic metalanguage to describe PixelMath syntax rules and syntactical elements. From now on, we'll use

italics-serif-text

to denote one or more metasymbols. The italicized serif text is used for elements that are not included in a PixelMath expression literally; rather, metasymbols represent elements of a PixelMath expression in formal terms. For example, in the following description of the StdDev function:

StdDev( img )

img is a metasymbol that represents a reference to an image. The following metasymbols are frequently used in this document:

expr	Represents an expression (or subexpression) of any kind.
img	Represents an image reference, optionally followed by a channel reference (see Section 2.3).
[]	Items enclosed by italicized square brackets are optional.
[c]	A channel reference with channel index c (see Section 2.3). Note that the square brackets are not part of the metalanguage in this case (they are not italicized), so they must appear literally.
x1[, x2[, ... xn]]	A sequence of function arguments (separated by commas) of length n >= 1. This construction is used to denote a function that takes a variable number of arguments.
[, opt-arg=def-value]	An optional function argument, denoted by opt-arg, whose default value (i.e. the value assumed when the argument is not specified) is denoted by def-value.

Immediate numerical expressions

Numerals can be specified in the usual way, including full support for scientific notation. Examples:

`123456`	The integer constant 123456
`1.23456`	The real number 1.23456
`0.123e-3`	Scientific notation. This is equivalent to 0.123×10^-3, or 0.000123.
`-0.12345d+05`	FORTRAN exponent separators (d, f) are also supported. This is equivalent to -12345.0

Immediate numerical expressions are always decoded and internally stored as 64-bit IEEE 754 floating point values.

Image identifiers

Image identifiers follow the usual rules in the PixInsight/PCL platform: they are case-sensitive, standard C identifiers. An image identifier can start with an alphabetic or underscore character, followed by any combination of alphabetic, decimal digit and underscore characters, without specific length limit. In the expression

Image01+Image02

both Image01 and Image02 are considered as image identifiers by the PixelMath expression parser. If no images with such identifiers are found upon expression execution, a runtime error occurs and the process fails. Since image identifiers are case-sensitive,

MyFinestImage
myfinestimage
myFinestImage

are considered as three different identifiers.

In a PixelMath expression, an image identifier is an image reference expression that evaluates to the current sample value in the referenced image. Channel references, covered in the next section, can be used to limit the value of an image reference expression to a particular channel.

Channel references

An image reference can optionally be followed by a channel reference in a PixelMath expression. A channel reference is an integer number enclosed by square brackets, following an image identifier, as in myImage[1]. Channels are indexed starting from zero, so [0] refers to the red and gray channels, for RGB color images and grayscale images, respectively. [1] and [2] refer to the green and blue channels of a color image, respectively, or to the first and second alpha channels of a grayscale image; [4] refers to the first alpha channel of a color image, or to the third alpha channel of a grayscale image, and so on.

When a channel reference is found, PixelMath will take only pixel samples from a specific channel of the referenced image. In absence of a channel reference, an image reference will cause PixelMath to take pixels from all channels of the referenced image, as appropriate.

For example, in the expression

Image01*Image02[2]

PixelMath will multiply each pixel of Image01 by the blue channel sample value of Image02 (i.e. the sample in the third channel of Image02), assuming that Image02 is a RGB color image.

Metaidentifiers

The PixelMath expression parser supports metasymbols. A metasymbol starts with the '$' character. Metasymbols are like special variables that are resolved at execution time. At the time of writing this document, PixelMath accepts only metaidentifiers, which are specified with the $target metasymbol.

In the following expression

$target*0.723

PixelMath will multiply the target image by the real constant 0.723. Note that the $target matasymbol will be replaced by the appropriate image reference when the above PixelMath expression is executed on a particular target image.

Metaidentifiers can have channel references, as in

Min( $target[1], $target[2] )

More metasymbols will be implemented in incoming versions of PixelMath.

Operators

Below is a list of all supported operators in precedence order. Operators with higher precedence are evaluated first. Two or more consecutive subexpressions with equal precedence are evaluated from left to right.

The highest precedence is 0 and corresponds to grouping operators (parentheses and channel selector). The lowest precedence is 12 and corresponds to the assignment operator.

Operator	Precedence	Syntax	Description
()	0	(expr)	Parentheses The expression expr is evaluated with the highest precedence.
[]	0	img[c]	Channel Selector A channel of an image img, whose index is c, is referenced. The channel index c must be an integer literal >= 0. If the specified channel does not exist, a runtime error occurs. Example: `Image01[2]` selects the blue channel of the `Image01` image.
~	1	~x	Pixel inversion ~x is equivalent to 1 - x, where x is expressed in the normalized real range from 0=black to 1=white.
-	1	-x	Unary Minus (Sign Change) The result is the operand x with negative sign if x is positive; positive x if x is negative.
+	1	+x	Unary Plus (No Change) The result is the operand x.
!	1	!x	Logical Negation The result is 1 if the operand is zero; 0 if the operand is nonzero.
^	2	x ^ y	Exponentiation The result if the left operand raised to the right operand.
*	3	x * y	Multiplication The left operand is multiplied by the right operand.
/	3	x / y	Division The left operand is divided by the right operand. If the right operand is zero or insignificant, the result is the maximum value in the normalized range (1=white) with the same sign as the numerator.
%	3	x % y	Modulus The result is the floating-point remainder of the x/y division, such that x = i×y + f, where i is an integer, f has the same sign as x, and the absolute value of f is less than the absolute value of y.
+	4	x + y	Addition The right operand is added to the left operand.
-	4	x - y	Subtraction The right operand is subtracted from the left operand.
--	4	x -- y	Absolute Difference The result is the absolute value of the difference x-y.
<	5	x < y	Less Than The result is 1 if x is less than y; 0 if x is greater than or equal to y.
<=	5	x <= y	Less Than Or Equal To The result is 1 if x is less than or equal to y; 0 if x is greater than y.
>	5	x > y	Greater Than The result is 1 if x is greater than y; 0 if x is less than or equal to y.
>=	5	x >= y	Greater Than Or Equal To The result is 1 if x is greater than or equal to y; 0 if x is less than y.
==	6	x == y	Equal To The result is 1 if x is equal to y; 0 if x is not equal to y.
!=	6	x != y	Not Equal To The result is 1 if x is not equal to y; 0 if x is equal to y.
&	7	x & y	Bitwise AND Both operands are transformed to 8-bit integers, their bitwise AND operation is performed, and the result is transformed back to the normalized real range [0,1].
!&	7	x !& y	Bitwise NAND Both operands are transformed to 8-bit integers, their bitwise NAND operation is performed, and the result is transformed back to the normalized real range [0,1].
&\|	8	x &\| y	Bitwise XOR Both operands are transformed to 8-bit integers, their bitwise XOR operation is performed, and the result is transformed back to the normalized real range [0,1].
!&\|	8	x !&\| y	Bitwise XNOR Both operands are transformed to 8-bit integers, their bitwise XNOR operation is performed, and the result is transformed back to the normalized real range [0,1].
\|	9	x \| y	Bitwise OR Both operands are transformed to 8-bit integers, their bitwise OR operation is performed, and the result is transformed back to the normalized real range [0,1].
!\|	9	x !\| y	Bitwise NOR Both operands are transformed to 8-bit integers, their bitwise NOR operation is performed, and the result is transformed back to the normalized real range [0,1].
&&	10	x && y	Logical AND The result is 1 if both operands are nonzero; 0 if any of the operands is zero.
\|\|	11	x \|\| y	Logical OR The result is 0 if both operands are zero; 1 if any of the operands is nonzero.
=	12	var = expr	Assignment The expression expr is evaluated, and its resulting value is assigned to the variable var. See Section 2.7 for a description of variables in PixelMath expression.

Functions

A function is an object that receives a set of zero or more arguments and returns a single value. Syntactically, a function is an expression that evaluates to the value that it returns, so functions can be used to build expressions just like images and operators.

Function arguments are specified as a list of zero or more expressions enclosed by parentheses. When a function takes more than one argument, arguments must be separated by commas.

The following table lists all functions currently supported by PixelMath.

Function syntax	Description	Invariant subexpression
Abs( x )	Absolute value of x	For invariant x
ArcCos( x )	Arc cosine of x, in radians.	For invariant x
ArcCosh( x )	Hyperbolic arc cosine of x	For invariant x
ArcSin( x )	Arc sine of x, in radians.	For invariant x
ArcSinh( x )	Hyperbolic arc sine of x	For invariant x
ArcTan( y[, x=1] )	Arc tangent of y/x, in radians. If the argument x is not specified, a value of one is assumed. This function always gives arc tangent values in the correct quadrant. The result is in the range [-π,+π].	For invariant x, y
ArcTanh( x )	Hyperbolic arc tangent of x	For invariant x
Avg( x1[, x2[, ... xn]] )	Arithmetic mean of a set x1, x2, ..., xn If a single argument is specified, as in Avg( img ), img must be a reference to an existing image. In this case, the Avg function is an invariant that is precalculated as the average sample value of the referenced image prior expression execution. If an optional channel selector is specified, as in Avg( img[c] ), the invariant average sample value is precalculated for the specified channel c exclusively.	Only with a single image reference argument, as in Avg( img ) or Avg( img[c] ).
AvgDev( img )	Average deviation of an image img.	Always
Ceil( x )	Lowest integer >= x	For invariant x
Cos( x )	Cosine of x, x in radians.	For invariant x
Cosh( x )	Hyperbolic cosine of x	For invariant x
Exp( x )	Exponential function, e^x	For invariant x
Floor( x )	Highest integer <= x	For invariant x
Frac( x )	Fractional part of x in the ]-1,+1[ interval.	For invariant x
Height( [img=$target] )	Height of an image img in pixels. If no argument is specified, the result is the height in pixels of the target image.	Always
iif( condition, expr-true, expr-false )	Conditional funcion (inline if function) The result is expr-true if the condition expression evaluates to a nonzero value. The result is expr-false if condition evaluates to zero.
IsColor( [img=$target] )	Color space of an image img The result is a nonzero value if img is a RGB color image. The result is zero if img is a grayscale image.	Always
Ln( x )	Natural logarithm of x (base e logarithm)	For invariant x
Log( x )	Base 10 logarithm of x	For invariant x
Log2( x )	Base 2 logarithm of x	For invariant x
Lum( img )	Luminance of the current pixel in an image img. The luminance is calculated in the current RGB working space (RGBWS) of the specified image. For grayscale images, the current pixel value in the nominal gray channel is returned.
Max( x1[, x2[, ... xn]] )	Maximum element of a set x1, x2, ..., xn If a single argument is specified, as in Max( img ), img must be a reference to an existing image. In this case, the Max function is an invariant that is precalculated as the maximum sample value of the referenced image prior expression execution. If an optional channel selector is specified, as in Max( img[c] ), the invariant maximum sample value is precalculated for the specified channel c exclusively.	Only with a single image reference argument, as in Max( img ) or Max( img[c] ).
MaxSample( img )	Maximum sample of the current pixel in a color image img. The argument img must be a reference to an existing RGB color image. This function is equivalent to: Max( img[0], img[1], img[2] ).
Med( x1[, x2[, ... xn]] )	Median value of a set x1, x2, ..., xn If a single argument is specified, as in Med( img ), img must be a reference to an existing image. In this case, the Med function is an invariant that is precalculated as the median sample value of the referenced image prior expression execution. If an optional channel selector is specified, as in Med( img[c] ), the invariant median sample value is precalculated for the specified channel c exclusively.	Only with a single image reference argument, as in Med( img ) or Med( img[c] ).
Min( x1[, x2[, ... xn]] )	Minimum element of a set x1, x2, ..., xn If a single argument is specified, as in Min( img ), img must be a reference to an existing image. In this case, the Min function is an invariant that is precalculated as the minimum sample value of the referenced image prior expression execution. If an optional channel selector is specified, as in Min( img[c] ), the invariant minimum sample value is precalculated for the specified channel c exclusively.	Only with a single image reference argument, as in Min( img ) or Min( img[c] ).
MinSample( img )	Minimum sample of current pixel in a color image img. The argument img must be a reference to an existing RGB color image. This function is equivalent to: Min( img[0], img[1], img[2] ).
MTF( m, x )	Midtones Transfer Function (MTF) of x for a midtones balance 0 <= m <= 1. The calculated MTF is identical to the same functional value in the standard HistogramTransform process. The MTF in PixInsight is a rational interpolation on a curve that passes through the points {0,0}, {m,0.5} and {1,1}, where the X axis represents input pixel sample values, and the Y axis represents output MTF values. The MTF is much more controllable and yields more accurate results than a standard gamma curve.
NumberOfChannels( [img=$target] )	Number of channels in an image img, including nominal and alpha channels.	Always
PAngle( [xc=Width()/2, yc=Height()/2] )	Polar angle of the current pixel in the target image, measured in radians in the range [-π,+π] with respect to a center point with the specified {xc, yc} coordinates. Positive position angles are reckoned clockwise, starting from the three o'clock direction. Negative position angles are reckoned counter-clockwise starting from the same direction. The reason why positive position angles cover the third and fourth quadrants, instead of the first and second as usual, is that vertical Y-axis coordinates grow from top to bottom. If no arguments are specified, polar angle is measured with respect to the geometric center of the target image.
Pi()	The π constant 3.141592...	Always
Pixel( img, x, y[, c] )	Sample value of an image img at the specified {x, y} pixel coordinates of the channel at index c. It is legal to specify pixel coordinates that are out of range, e.g. negative x or y coordinates, x >= Width(), or y >= Height(). In these cases, the returned value is zero. Channels are indexed in the range [0,NumberOfChannels()-1]. If no channel index is specified, sample values are taken from the current channel in the specified image img (i.e. the channel where PixelMath is being executed). Nonexistent channel indexes are not allowed. If an out-of-range channel index is specified, a runtime error occurs.
RDist( [xc=Width()/2, yc=Height()/2] )	Radial distance of the current pixel in the target image, in pixels, with respect to a center point with the specified {xc, yc} coordinates. If no arguments are specified, radial distance is measured with respect to the geometric center of the target image.
Round( x[, n=0] )	Round function: x rounded to n decimal digits. If no second argument is specified, the Round function rounds x to its nearest integer. Rounding is performed towards ∓∞, so that: Round( 1.5 ) returns 2 Round( -1.5 ) returns -2	For invariant x, n
Sign( x )	Sign of x: +1 if x is positive or zero; -1 if x is negative.	For invariant x
Sin( x )	Sine of x, x in radians.	For invariant x
Sinh( x )	Hyperbolic sine of x.	For invariant x
Sqrt( x )	Square root of x.	For invariant x
StdDev( img )	Standard deviation of an image img.	Always
Sum( x1[, x2[, ... xn]] )	Summatory of a set x1, x2, ..., xn
Tan( x )	Tangent of x, x in radians.	For invariant x
Tanh( x )	Hyperbolic tangent of x.	For invariant x
Trunc( x )	Truncated integer part of x.	For invariant x
Var( img )	Variance of an image img.	Always
Width( [img=$target] )	Width of an image img in pixels. If no argument is specified, the result is the width in pixels of the target image.	Always
XPos()	Horizontal coordinate (X-axis) of the current pixel in the target image.
X()	Normalized horizontal coordinate (X-axis) of the current pixel in the target image. The returned value is in the range [0,1], where 0 corresponds to the zero horizontal pixel coordinate and 1 corresponds to Width()-1.
YPos()	Vertical coordinate (Y-axis) of the current pixel in the target image.
Y()	Normalized vertical coordinate (Y-axis) of the current pixel in the target image. The returned value is in the range [0,1], where 0 corresponds to the zero vertical pixel coordinate and 1 corresponds to Height()-1.

Variables

An unlimited number of variables can be declared in an instance of the PixelMath process. A variable is denoted by its identifier, and can be thought of as a placeholder to store and manipulate a single real numerical value.

Variable identifiers follow the same exact syntax rules as image identifiers (see Section 2.2). For this reason, variables must be explicitly declared: any identifier that has not been declared as a variable is treated as a reference to an image by the PixelMath interpreter.

To declare one or more variables, the PixelMath interface provides a set of specific controls: the Variables field and the Variable Editor dialog.

Variables can be used as right-hand-side operands in PixelMath expressions, or as left-hand-side operands in assignment statements (see Section 2.8). For example, in the statement:

x = 1 + y

the variable x is being assigned the result of the 1+y expression, where y could be either a variable or an image reference.

We'll see a good example of variables in the next section, working along with multiple statements.

Multiple statements

A PixelMath expression can include any number of statements. Statements are syntactical constructs that cannot be evaluated, in contrast with expressions. As in most programming languages, a statement is executed for its side-effects in a PixelMath expression. Currently the only effect that statements can produce is variable assignment in PixelMath.

Statements must be separated with semicolons. The general syntax of a multiple-statement PixelMath expression is:

[statement1;[ statement2;[ ... statementn;]]] expr[;]

where the last expression expr must evaluate to a single value, which will be the value assigned to the current sample of the target image. The terminating semicolon is optional.

Statements are extremely useful to build complex PixelMath expressions. They allow to perform temporary calculations and store their results in variables. Those variables can then be used in the expression that evaluates to replace the current sample value of the target image.

Consider the following PixelMath expression:

d = Abs( RDist() - 250 );
f = (1.25 - iif( d < 1.25, d, 1.25 ))/1.25;
f + (1 - f)*$target

This expression draws an antialiased circumference of radius 250 pixels, using a white pen whose width is 1.25 pixels. The circumference is drawn centered at the geometric center of the target image. Figure 1 shows an application example.

Figure 1— A PixelMath expression to draw an antialiased circumference.

Click on the image to enlarge.

As you can easily verify by yourself (just try to obtain an equivalent single expression by doing the necessary replacements), the same task would require an extremely complex and confusing expression without the help of variables and multiple statements.

Note that the last statement, namely

f + (1 - f)*$target

in this example, is always the expression whose value is assigned to the current sample of the target image.

In the expression above, we are using two variables: d and f. However, the radius and pen width are being specified literally (250 and 1.25 pixels, respectively). We can use two additional variables to store those constant quantities, which will lead to a much more readable and manageable PixelMath expression:

radius = 250;
penWidth = 1.25;
d = Abs( RDist() - radius );
f = (penWidth - iif( d < penWidth, d, penWidth ))/penWidth;
f + (1 - f)*$target

Now if you want to change the pen width, for example, you just have to edit a single item (the 1.25 constant in the second statement). Without these additional variables, four values would have to be edited, in a much more error-prone operation. In addition, the radius and penWidth variable names make the expression much more readable. Using additional variables in this way only has a negligible impact on PixelMath's execution times.

3. PixelMath Examples

Averaging a set of images

The best way to increase the signal-to-noise ratio of an image is to average several shots of the same scene that have been previously registered (mutually aligned). For example, if we want to average five images, a PixelMath expression could be as simple as:

Image01 + Image02 + Image03 + Image04 + Image05

assuming that automatic rescaling is enabled (the Rescale result option on PixelMath interface; see Section 1). With automatic rescaling disabled, the correct expression would be:

(Image01 + Image02 + Image03 + Image04 + Image05)/5

Since rescaling is enabled by default, you usually will prefer expressions like the former one, so there's no need to care about the number of averaged images.

2a— Original image

2b—Random noise added, uniform distribution, 25%

2c—Average of 2 noisy images

2d—Average of 4 noisy images

2e—Average of 8 noisy images

2f—Average of 16 noisy images

Figure 2— A PixelMath image integration experiment.

Linear color corrections

By multiplying each channel of a RGB color image by a constant value, we apply a linear color correction. For example:

0.983 * $target[0]
0.870 * $target[1]
0.955 * $target[2]

The above three expressions correspond to the red, green and blue channels of the target image. For these expressions to work, the Use a single expression for all channels option must be disabled on the PixelMath interface. The Rescale result option can be enabled or disabled, depending on whether you want to redistribute the resulting values over the full dynamic range or not.

If you specify each of the above expressions to the red, green and blue channels respectively, in that order, channel references can be omitted, and the following expressions yield exactly the same result:

0.983 * $target
0.870 * $target
0.955 * $target

since PixelMath reads sample values from the corresponding channels by default.

As long as you apply only linear operations (addition, subtraction, multiplication and division), the linearity (or nonlinearity) of the target image won't change. This includes the automatic rescaling operation in PixelMath, which is also linear. This is very important and must be taken into account when manipulating raw linear image data.

Posterization

Posterization consists in reducing the number of existing colors or gray levels in an image. Some image structures can be easier to detect in this way because after posterization they are clearly isolated by hard brightness steps. This is very useful to define isophotes over deep-sky images, for example to draw the contours of nebulae and other nonstellar objects, to generate sky charts.

It is very easy to write a PixelMath expression to produce the desired posterization. Here is an example:

Round( $target*10 )/10

for ten posterization steps or, using variables to make things clearer:

n = 10; Round( $target*n )/n

where n is the number of desired colors or gray levels after posterization. Figure 3 shows an example.

Figure 3— Posterization with the PixelMath expression:

Round( $target*8 )/8

To the left, the same expression applied to a grayscale version of the image.

Nonlinear transformations

Nonlinear functions can be applied to redistribute existing pixel values over the available dynamic range. This allows us to change their relative mapping to the shadows, midtones and highlights parts of the histogram. The particular nonlinear transform to apply in each case depends on the effect that we want to achieve; e.g. to increase the midtones without changing the highlights, or vice-versa, etc. Here is a brief survey with some of the nonlinear functions that have proved to be useful, and also some that we like for some reason, perhaps because we find them elegant.

Gamma stretch

The most well-known nonlinear transform is the gamma stretch function, consisting on raising pixel values to a given constant γ. For example:

$target ^ 0.75

where γ = 0.75. When γ < 1, the gamma stretch function brightens the target image because it pushes the midtones up. The inverse effect (darkening) is achieved for γ > 1.

A frequently-used gamma stretch function can be conveniently expressed as the square root:

Sqrt( $target )

which is equivalent to γ = 1/2.

Figure 4— Gamma stretch functions. Top-left: original image. Top-right: γ = 0.75. Bottom-left: γ = 1.5. The applied PixelMath expressions are, respectively:

$target ^ 0.75$target ^ 1.5

Exponential functions

Another transform with similar properties is the exponential function, namely e^x. The exponential function can be very useful to reveal details that are hidden in the shadows (Figure 5). Here is a PixelMath expression to implement it:

~Exp( -$target )

Note that the argument of Exp() must be a negative number for the function to give a meaningful result, and that the result must be inverted. Compared to a similar gamma stretch function, the exponential function improves the midtones with superior preservation of the shadows and highlights.

Figure 5— Exponential function. Left: Original image. Right: After applying the PixelMath expression:

~Exp( -$target ).

Logarithmic functions

PixelMath provides three logarithmic functions: natural (base e), base 10 and base 2 logarithms. Logarithmic functions are a bit more difficult to use. One must make sure that PixelMath will never try to calculate a logarithm of zero, since these functions are undefined for that value. Here is an useful PixelMath expression involving a logarithmic function:

alpha = 0.7;
Ln( 1 + $target )/alpha

This expression is intended to work with automatic rescaling turned off. The alpha variable allows to control the contrast manipulation that the expression applies (the growth of the logarithmic function, to be more precise). The appropriate alpha value depends on the image and on the desired result; you can try in the range from 0.5 to 1.5 for most images with reasonable contrast/brightness ratios.

If you ever need a logarithm of an arbitrary base b > 0 that is not one of e, 10, or 2, you can make use of the following relation:

Base 10 logarithms can also be used (any base, in fact), but base 2 logarithms are slightly faster.

Figure 6— Logarithmic function. Left: Original image. Right: After applying the PixelMath expression:

alpha = 0.5; Ln( 1 + $target )/alpha.

Trigonometric functions

The arguments of the sine, cosine and tangent functions must be expressed in radians. Assuming that the target image has no values outside the normalized range [0,1], you can apply sine and cosine in the range [0, π/2] with the following expressions:

Sin( Pi()/2*$target )
~Cos( Pi()/2*$target )

The result of the cosine function must always be inverted, because cos( 0 ) = 1 and cos( π/2 ) = 0. The sine function brightens the image, and the cosine function has the inverse effect. The tangent function is not well behaved because it is singular for angles close to π/2 in absolute value. You must apply the tangent in a range where it won't produce too large results, that is. For example, in the range [0, π/3], and assuming again that the target image is normalized:

Tan( Pi()/3*$target )

Figure 7— Trigonometric functions. Top-left: original image. Top-right: sine function. Bottom-left: inverse of cosine function. Bottom-right: tangent function restricted to the range [0,π/3]. The applied PixelMath expressions are, respectively:

Sin( Pi()/2*$target )
~Cos( Pi()/2*$target )Tan( Pi()/3*$target )

Sigmoidal functions

Sigmoidal functions can be very useful because they allow us to control the midtones, shadows and highlights independently. This is because sigmoids are S-shaped curves.

The logistic function, namely

can be implemented in PixelMath as

1/(1 + Exp( -$target ))

and should be used with automatic rescaling enabled. More interesting transforms involve trigonometric functions; for example, this one

improves overall contrast by reducing the shadows and increasing the highlights (the typical S-curve). The α variable (0 < α ≤ 1) controls the function's strength (more aggressive for smaller α values). This is the corresponding PixelMath expression for α=1/2:

a = 0.5*Pi();
(1 + 1/Sin( a/2 ) * Sin( a*($target - 1/2) ))/2

The above expression can be used without automatic rescaling, since it is bounded to [0,1].

A similar function with the tangent:

with the α variable having the same meaning as before, has the inverse effect: the shadows are increased and the highlights are decreased, which diminishes overall contrast. This is the PixelMath implementation for α=1/2:

a = 0.5*Pi();
(1 + 1/Tan( a/2 ) * Tan( a*($target - 1/2) ))/2

Figure 8— Sigmoidal functions. Top-left: original image. Top-right: sigmoidal/sine function. Bottom-left: sigmoidal/tangent function. The applied PixelMath expressions are, respectively:

a = 0.5*Pi(); (1 + 1/Sin( a/2 ) * Sin( a*($target - 1/2) ))/2
a = 0.5*Pi(); (1 + 1/Tan( a/2 ) * Tan( a*($target - 1/2) ))/2

Hyperbolic arc sine

The hyperbolic arc sine function is also interesting because it is naturally bounded to the normalized [0,1] range. Lupton et al. (arXiv:astro-ph/0312483v1) propose the following formulation to implement an hyperbolic arc sine stretch:

where the α parameter applies a linear stretch, Q governs the improvement of bright features, and m is the minimum displayed value (we always have m=0 in PixInsight/PixelMath). This is the PixelMath expression that implements this transform in PixInsight:

a = 0.01;
Q = 5;
ArcSinh( a*Q*$target )/Q

As a general rule of thumb, useful values for α seem to be in the range from 0.1 to 0.001, and for Q between 1 and 20.

Adaptive midtones transfer functions

A midtones transfer function (MTF) can be easily controlled to apply contrast and saturation transforms that depend on current pixel values. A very simple PixelMath expression:

MTF( ~$target, $target )

can yield very good results with many poorly-contrasted images. The above expression must be applied to all channels of a RGB color image, or to the nominal channel of a grayscale image, with the Rescale result option enabled. Figure 9 shows a good example.

Figure 9— Adaptive MTF. Left: Original image. Right: After applying the PixelMath expression: MTF( ~$target, $target ).

Merging images

There are many ways and techniques to combine images with PixelMath. For example, the Max and Min functions are useful to combine two or more images selectively, or to replace certain pixels with constant values (see Section 2.6 for a reference of all supported functions). Figure 10 shows several interesting examples. In the cases where images are being composed without mixing pixel values (Figures 10a, 10b), the Crop tool (Geometry category) has been used to generate free black areas very easily.

Exercise: ¿Can you figure out a way to generate a r-g-b image with PixelMath, as on Figure 10a, without needing three separate images?

Another exercise, again about Figure 10a. ¿How would you obtain the same result if the three source images (Red, Green and Blue on the figure) were grayscale images instead of RGB color?

10a— Merging three color images with the Max function.

Click on the image to enlarge.

10b— Mirroring an image with the iif and Pixel functions.

Click on the image to enlarge.

10c— Merging two images with PixelMath. From top to bottom: 1—PixelMath setup in PixInsight (click to enlarge); 2—Manhattan source image; 3—M45 (Pleiades) source image, conveniently cropped; 4—the resulting composite image.

Figure 10c deserves a more in-depth analysis. Here is the expression used to generate the Alien Manhattan landscape:

k1 = 0.87;
k2 = 0.92;
m = Min( 1, Max( 0, ($target[2] - k1)/(k2 - k1) ) );
$target*(1 - m) + Pleiades*m

The k1 and k2 constants are the two break points that define how both images are being merged. Below the k1 value, the target image (Manhattan in this case) will be fully preserved. Above the k2 value, the Pleiades image will replace the target image completely. For values between k1 and k2, the m variable is calculated as a masking value that will control how both images are merged proportionally (a simple linear proportion is being calculated in this case, but more refined mixing strategies can be easily defined), so that the transitions between them will be smooth. Finally, the last expression calculates the mixed sample value between both images.

In the example of Figure 10c, the blue channel of the target image ($target[2] in the expression) has been used as the reference value that decides how both images are combined. This is because in the Manhattan image that we have used the sky has been represented with relatively strong blue values.

This merging technique, and similar ones, can be used with many images. The key is to find some property of one of the images that allows us to shoot through it with a PixelMath expression; then it's easy to fill the gaps with other image.

Interchanging color channels

These expressions, applied to the red, green and blue channels respectively, swap the red and blue channels of a target RGB image:

$target[2]
$target[1]
$target[0]

Conditional expressions (using the iif function, see Section 2.6) allow us to implement more elaborated channel swapping routines. For example, we applied the following expressions to the red, green and blue channels, respectively,

$target
iif( $target[2] < 0.77, $target[2], $target )
iif( $target < 0.77, $target[1], $target )

to exchange green and blue selectively in the image shown on Figure 11. Note that in the expression above we are using a technique that is conceptually similar to what we applied in the example of Figure 10c. Conditional expressions can be much more sophisticated to adapt these techniques to more difficult cases.

Figure 11— Left: Original image. Right: After applying the following PixelMath expressions (red, green and blue channels, respectively):

$target iif( $target[2] < 0.77, $target[2], $target ) iif( $target < 0.77, $target[1], $target )

Convolutions

It is relatively straightforward to write PixelMath expressions implementing spatial convolutions with small filters. Although PixelMath is definitely not the best tool for this task (the JavaScript interface provides predefined functions that run at native machine speed), this is an interesting exercise to gain experience with moderately complex PixelMath expressions.

First a low-pass filter. The box average filter can be implemented as a direct convolution with the following kernel:

1 1 1 1 1 1 1 1 1

This is the simplest low-pass filter. It replaces each pixel in the target image with the arithmetic mean of it and its eight nearest neighbors. Here is a PixelMath expression that implements this convolution:

x = XPos();
y = YPos();
xm1 = x-1;
xp1 = x+1;
ym1 = y-1;
yp1 = y+1;
p11 = Pixel( $target, xm1, ym1 );
p12 = Pixel( $target, x,   ym1 );
p13 = Pixel( $target, xp1, ym1 );
p21 = Pixel( $target, xm1, y   );
p22 = Pixel( $target, x,   y   );
p23 = Pixel( $target, xp1, y   );
p31 = Pixel( $target, xm1, yp1 );
p32 = Pixel( $target, x,   yp1 );
p33 = Pixel( $target, xp1, yp1 );
(p11+p12+p13+p21+p22+p23+p31+p32+p33)/9

The p_ij variables are used to store neighbor sample values at the i-th rows and j-th columns of the kernel. p22 is the central value, that is, the sample value below the central kernel element, which is also the current sample value during PixelMath execution. The xm1, xp1, ym1 and yp1 variables are used for optimization by suppressing redundant calculations. The final expression is just the average of the nine sample values that fall below the filter kernel.

Now a high-pass filter. Here we can choose between several options, but we'll implement a classical mild sharpen filter, whose kernel is:

0, -1, 0 -1, 6, -1 0, -1, 0

Here is the PixelMath expression for a mild sharpen filter:

x = XPos();
y = YPos();
p12 = Pixel( $target, x,   y-1 );
p21 = Pixel( $target, x-1, y   );
p22 = Pixel( $target, x,   y   );
p23 = Pixel( $target, x+1, y   );
p32 = Pixel( $target, x,   y+1 );
(6*p22 - (p12 + p21 + p23 + p32))/2

This expression is much simpler than the box average filter, mainly because there are only five nonzero kernel elements.

Of course, since we haven't implemented any correction for boundary artifacts, they will occur on processed images. The low-pass filter will generate a one-pixel dark border artifact at the four edges of the target image, and the high-pass filter a one-pixel bright border artifact.

Compared to equivalent, optimized JavaScript routines, these PixelMath expressions are quite slow. This is because PixInsight's JavaScript runtime provides access to high-performance native convolution routines. However, a JavaScript implementation without calls to runtime convolution functions is much slower than PixelMath, the difference being dramatical on multiprocessor systems.

Figure 12— PixelMath convolutions. Top-left: original image. Top-right: box average filter. Bottom-left: mild sharpen filter.

Drawing a projected torus

If you don't know what a torus is, think on a ring, or a doughnut; a torus is a similar 3-D surface. We'll write a PixelMath expression to draw the normal projection of a torus over an image. Call it a ring, for simplicity.

Again, JavaScript can be used to accomplish this task very efficiently. However, the PixelMath implementation has several niceties, as the possibility to generate graphics directly with floating point accuracy (with JavaScript we generate all graphics on 32-bit RGBA bitmaps, which provide 8 bits per color, then bitmaps are blended with images of any data type).

Without more preambles, this is the expression:

w = Min( Width(), Height() )/10;
R = Min( Width(), Height() )/3;
r = RDist();
d = Abs( r - R );
f = Sin( Pi()/2*(w - iif( d < w, d, w ))/w );
f + (1 - f)*$target

In the expression above, w and R are, respectively, the thickness and radius of the ring in pixels. The 10 and 3 constant values can be varied to yield figures with different proportions. Note that the ring is drawn centered at the geometric center of the target image. This is because we are using the RDist function (radial distance) without parameters.

The call to Sin() and the corresponding scaling to the range [0, π/2] simulate the illumination on a (circular) toroidal profile in terms of the transparency of the generated pixel values. This yields a translucent rendition. Figure 13 shows an example.

Note that this expression draws white (or gray) pixels. Drawing with arbitrary colors is straightforward; it is left as an exercise to the reader.

And a further exercise: write a PixelMath expression to draw the same figure centered at an arbitrary location.

Figure 13— A translucent ring drawn with PixelMath.