image = multiscaleVesselEnhancement2dImg (inImg,inMVEParams)

Detailed Description

Multiscale vessel enhancement using Frangi's approach.

The MultiscaleVesselEnhancement2dImg algorithm was proposed by Frangi et al. in 1998 [1]. It is an iterative approach used to enhance tube-like features in an image where each scale focuses on a specific tube width range.

The algorithm steps are applied for each pixel, for this reason, we avoid the notation of the coordinates $(\textbf{x})$ in the following. For a given scale $s$ , the algorithm computes at each pixel the Hessian matrix $H$ . Its eigen values $\lambda_i, i \in \left\{1, 2\right\}$ are then extracted and sorted in order to identify :

$\vert \lambda_1 \vert \leq \vert \lambda_2 \vert$

A blobness coefficient $R_B$ and a structureness coefficient $S$ are then calculated from $\lambda_1$ and $\lambda_2$ :

$\begin{eqnarray*} R_B & = & \frac{\lambda_1}{\lambda_2} \\ S & = & \sqrt{\lambda_1^2 + \lambda_2^2} \end{eqnarray*}$

The vesselness for a scale $s$ , noted $\nu_s$ , is finally computed as a combination of these coefficients. For bright features on dark background, this measure can be expressed as:

$\nu_s = \begin{cases} 0& \text{if } \lambda_2 > 0 \\ \exp \left(-\frac{R_B^2}{2 \beta^2} \right) \left( 1 - \exp \left(-\frac{S^2}{2c^2} \right) \right) &\text{otherwise} \end{cases}$

Where $\beta$ and $c$ are respectively the blobness and structureness sensitivity thresholds. $\beta$ is an input parameter and $c$ is calculated as:

$c = \frac{1}{2} \max_{\textbf{x} \in \Omega} \left( \Vert H(\textbf{x}) \Vert \right)$

Where $\Omega$ is the image domain and $H(\textbf{x})$ is the Hessian matrix of the image at pixel $\textbf{x}$ .

For dark features on bright background, the condition $\lambda_2 > 0$ becomes $\lambda_2 < 0$ .

The final result is an image containing the maximum value of $\nu$ along all the scales on each pixel:

$\nu(\textbf{x}) = \max_{s \in Scales}\left( \nu_s(\textbf{x}) \right)$

Where $Scales$ is the collection of scales given as input parameter.

The algorithm parametrization is done using the MVEParams data item, which contains:

the mandatory BlackWhite field, set to true if the user wants to enhance bright vessels on dark background ; it equals false for dark vessels on bright background,
the mandatory Scales field, a collection of Real64 values containing the scale list,
BlobSensitivityThreshold, an optionnal parameters corresponding to the $\beta$ threshold (by default, $\beta = 0.5$ ),
AreaSensitivityThresold, an optionnal parameters, not used in 2d case.

Moreover, the algorithm allows the user to provide a class image, with type UInt8. If this image is set, the algorithm will fill it by the scale index corresponding to the maximum $\nu_s$ :

$OutOptClassImg(\textbf{x}) = argmax_{s \in Scales}\left( \nu_s(\textbf{x}) \right)$

Here is an example of vessel enhancement on a 2d gray level image with the scales $\left \{ 1, 3, 5\right \}$ and its corresponding class image:

The input image has been downloaded from wikipedia: https://commons.wikimedia.org/wiki/File:Mra-mip.jpg (Attribution: By SBarnes (Own work) [CC BY-SA 3.0 (https://creativecommons.org/licenses/by-sa/3.0) or GFDL (http://www.gnu.org/copyleft/fdl.html)], via Wikimedia Commons)

References

[1] "Multiscale vessel enhancement filtering", Frangi, A. F.; Niessen, W. J.; Vincken, K. L. & Viergever, M. A., Medical Image Computing and Computer-Assisted Intervention - MICCAI '98, Springer Verlag, pages 130-137, 1998

Example of Python code :

Example imports

import PyIPSDK

import PyIPSDK.IPSDKIPLFiltering as filter

Code Example

    # open the input images
    inImg = PyIPSDK.loadTiffImageFile(inputImgPath)
    
    # Create the class image
    outVesselnessImg = PyIPSDK.createImage(PyIPSDK.eImageBufferType.eIBT_Real32, inImg.getSizeX(), inImg.getSizeY())
    outClassImg = PyIPSDK.createImage(PyIPSDK.eImageBufferType.eIBT_UInt8, inImg.getSizeX(), inImg.getSizeY())
    
    # Multiscale vessel enhancement computation
    mveParams = PyIPSDK.createMVEParameter([1, 3, 5], True)
    filter.multiscaleVesselEnhancement2dImg(inImg, mveParams, outVesselnessImg, outClassImg)

Example of C++ code :

Example informations

Header file

#include <IPSDKIPL/IPSDKIPLFiltering/Processor/MultiscaleVesselEnhancement2dImg/MultiscaleVesselEnhancement2dImg.h>

Code Example

        // Definition of the algorithm parameters
        MVEParamsPtr pMVEParams = createMVEParameter(vScales, false);
        // Calculation of the vesselness image only
        // ========================================
        // Sample with a generated output image
        // ------------------------------------
        ImagePtr pAutoOutImg = multiscaleVesselEnhancement2dImg(pInImg, pMVEParams);
        // Sample with a provided output image
        // -----------------------------------
        // create output image
        ImageGeometryPtr pOutputImageGeometry = geometry2d(eImageBufferType::eIBT_Real32, sizeX, sizeY);
        boost::shared_ptr<MemoryImage> pOutImg(boost::make_shared<MemoryImage>());
        pOutImg->init(*pOutputImageGeometry);
        // compute vesselness
        multiscaleVesselEnhancement2dImg(pInImg, pMVEParams, pOutImg);
        // Calculation of the vesselness and class images
        // ==============================================
        // create output image
        ImageGeometryPtr pOutputClassImageGeometry = geometry2d(eImageBufferType::eIBT_UInt8, sizeX, sizeY);
        boost::shared_ptr<MemoryImage> pVesselnessImg(boost::make_shared<MemoryImage>());
        boost::shared_ptr<MemoryImage> pClassImg(boost::make_shared<MemoryImage>());
        pVesselnessImg->init(*pOutputImageGeometry);
        pClassImg->init(*pOutputClassImageGeometry);
        // compute vesselness and the class image
        multiscaleVesselEnhancement2dImg(pInImg, pMVEParams, pVesselnessImg, pClassImg);