Image Based Calculations
Table of Contents
see also: Index
Overview
The “Image-based calculations” feature allows to evaluate free-text expressions like
1.5 * (a - b) / c**2 + gt(a, 7.3) * log(d)
where a, b, c and d are variables denoting image buffers, gt() is a thresholding function. This feature can be used interactively
(File menu → Automate) and in batch-mode (documentation pending).
A family of thresholding functions can also be used in the free-text
expression: in the example above, gt() will first create a clone of
image buffer a and then set all voxels greater than the threshold
(here: 7.3) to one, all voxels below and equal to the threshold to
zero. A reference of all functions we have implemented (so far) is
given below.
You can use most of the mathematical standard functions e.g. log,
sqrt, sin, etc. A reference of standard functions we have
implemented (so far) is given below.
Simple arithmetic operations involving image buffers are voxel-based
(similar to the conventions used in the IDL programming language),
e.g. the image buffer resulting from the operation a+b is derived
from adding corresponding voxels of the image buffers referred to by
the variables a and b, resp. Multiplications with a constant
factor will multiply all voxels of an image buffer with that constant,
adding a constant value will add that value to all voxels. In case of
divisions, if a zero-voxel is encountered in the denominator, we will
assign a value of zero to the corresponding voxel in the resulting
image. Exponentiation can be done with any real number, undefined
values e.g. (-1.0)**0.5 will be set to zero.
All image-based operations (so far) work on “native” voxel volume of images. If the reslicing parameters of the images are different, Vinci will refuse computation giving the error message: Images are computed in image reset position, so the involved images MUST be co-registered in image reset position! Your images have different reslicing parameters. In this case, we prohibit the computation now as it could lead to unexpected results (maybe even unperceived). If you really have to use a resliced image (e.g. co-registered by MMM) you can first call “reset link channel” on one image, then save the other image using “current offset and rotation” and re-load the saved image.
Also, we will not automatically resample “incompatible” images: only image volumes with identical geometries (dimension and pixel size) are valid in one calculation.
Image-based calculations are among Vinci’s most versatile and powerful
features. Part of the flexibility is due to our VinciPy package,
written in the Python programming language: it allows to
remotely-control Vinci through XML-based VinciScript, also at the
heart of the high-level interfaces to C/C++, PERL and IDL. The full
source code is provided in the external/Vinci_Python directory. The
implementation of this feature is rather compact as most of the work
is left to the Python interpreter: we overload operators in an image
class, interested readers should have a look at
VinciPy/Vinci_ImageCalc.py. Clearly, adding further functions is a
straight-forward task, please let us know if you have any suggestions.
Requirements
Image-Based Calculation requires a Python interpreter (versions 2.4-2.6 will do). Python should be available on all modern Unix systems (including MacOS X) and no further installation is required. Users of MS Windows systems can get a free standard Windows installer from here: http://python.org/download (installation is just a “double-click”, the distribution’s footprint is small).
Reference
“Logical” Thresholding
gt(a, t)
Create a clone of imagebuffer a and set all voxels above threshold t to one, and all voxels below and equal to the treshold to zero.
ge(a, t)
Similar to gt() but tests for “greater than or equal”.
lt(a, t)
Create a clone of imagebuffer a and set the all voxels below threshold t to one, and all voxels above and equal to the treshold to zero.
le(a, t)
Similar to lt() but tests for “lesser than or equal”.
eq(a, t)
Create a clone of imagebuffer a and set all voxels equal (i.e. with relative difference less than 10**(-6)) to threshold t to one, and all other voxels to zero.
neq(a, t)
Create a clone of imagebuffer a and set all voxels equal (i.e. with relative difference less than 10**(-6)) to threshold t to zero, and all other voxels to one.
General Thresholding
th_u(a, t, r)
Create a clone of imagebuffer a and set the “upper part”, i.e. all voxels above threshold t to r, keep all other voxels.
th_ue(a, t, r)
Similar to th_u[pper] () but tests for “greater than or equal”.
th_l(a, t, r)
Create a clone of imagebuffer a and set the “lower part”, i.e. all voxels below threshold t to r, keep all other voxels.
th_le(a, t, r)
Similar to th_l[ower] () but tests for “lesser than or equal”.
th_eq(a, t, r)
Create a clone of imagebuffer a and set all voxels equal to threshold t to r, keep all other voxels.
th_neq(a, t, r)
Create a clone of imagebuffer a and set all voxels not equal to threshold t to r, keep all other voxels.
th_u0(a, t)
:= th_u(a, t, 0)
th_ue0(a, t)
:= th_ue(a, t, 0)
th_u1(a, t)
:= th_u(a, t, 1)
th_ue1(a, t)
:= th_ue(a, t, 1)
th_l0(a, t)
:= th_l(a, t, 0)
th_le0(a, t)
:= th_le(a, t, 0)
th_l1(a, t)
:= th_l(a, t, 1)
th_le1(a, t)
:= th_le(a, t, 1)
Mathematical Standard Functions
In the following definitions, a denotes an image buffer.
pow(a, t)
computes a**t with any real number exponent t. pow is defined for
all voxel values greater zero, for voxels values equal zero with
exponents greater equal zero, and for voxel values less than zero with
integer exponents; offending result voxels will be set to zero.
a**t
:= pow(a, t)
sqrt(a)
defined for voxel values greater or equal zero; offending result voxels will be set to zero.
log(a)
defined for voxel values greater zero; offending result voxels will be set to zero.
exp(a)
result voxels with under- and overflow will be set to zero.
sin(a)
defined for voxel values in the range [-263,+263]; offending result voxels will be set to zero. This is a technical limit set by the programming language C++, 263 is about 83.66 π.
cos(a)
defined for voxel values in the range [-263,+263]; offending result voxels will be set to zero. This is a technical limit set by the programming language C++, 263 is about 83.66 π.
tan(a)
defined for voxel values in the range [-263,+263]; offending result voxels will be set to zero. This is a technical limit set by the programming language C++, 263 is about 83.66 π.
asin(a)
defined for voxel values in the range [-1,+1]; offending result voxels will be set to zero.
acos(a)
defined for voxel values in the range [-1,+1]; offending result voxels will be set to zero.
atan(a)
results voxels will be in the range (-π/2,+π/2).
We could easily add other functions. Please, let us now if you need them.
Multi-Frame Files
sum(a)
Python’s built-in sum function works now as the index function and the length function are defined for images to give a specific frame and the total number of frames. BEWARE: frame numbers for the index function always start with zero, they are transported to the first “native” frame index according to the image type: e.g. “1” for ECAT7 images and “0” for HRRT images.
msum(a, iStart, iEnd)
Sums all frames from frame <iStart> up to and including frame <iEnd> of the multi-frame file associated with image buffer a (in case of ECAT7 data) or the multi-frame series (in case of HRRT data). “Native” frame numbers are used: the first frame of an ECAT7 multi-frame file is “1” (usually), the first frame of a HRRT multi-frame series is “0”.
msum(a, iStart)
Sums all frames of the multi-frame file, or the multi-frame series starting with frame <iStart>.
msum(a)
Sums all frames of the multi-frame file, or the multi-frame series.
mavg(a, iStart, iEnd)
Similar to msum(a, iStart, iEnd) but scales the result with
1/#frames. Please note: this function is currently only available for
ECAT7 multi-frame files.
mavg(a, iStart)
Similar to msum(a, iStart) but scales the result with
1/#frames. Please note: this function is currently only available for
ECAT7 multi-frame files.
mavg(a)
Similar to msum(a) but scales the result with 1/#frames. Please
note: this function is currently only available for ECAT7 multi-frame
files.