tomopy.prep.stripe

Module for pre-processing tasks.

Functions:

remove_stripe_fw(tomo[, level, wname, ...])

Remove horizontal stripes from sinogram using the Fourier-Wavelet (FW) based method [B4].

remove_stripe_ti(tomo[, nblock, alpha, ...])

Remove horizontal stripes from sinogram using Titarenko's approach [B13].

remove_stripe_sf(tomo[, size, ncore, nchunk])

Normalize raw projection data using a smoothing filter approach.

remove_stripe_based_sorting(tomo[, size, ...])

Remove full and partial stripe artifacts from sinogram using Nghia Vo's approach [B24] (algorithm 3).

remove_stripe_based_filtering(tomo[, sigma, ...])

Remove stripe artifacts from sinogram using Nghia Vo's approach [B24] (algorithm 2).

remove_stripe_based_fitting(tomo[, order, ...])

Remove stripe artifacts from sinogram using Nghia Vo's approach [B24] (algorithm 1).

remove_large_stripe(tomo[, snr, size, ...])

Remove large stripe artifacts from sinogram using Nghia Vo's approach [B24] (algorithm 5).

remove_dead_stripe(tomo[, snr, size, norm, ...])

Remove unresponsive and fluctuating stripe artifacts from sinogram using Nghia Vo's approach [B24] (algorithm 6).

remove_all_stripe(tomo[, snr, la_size, ...])

Remove all types of stripe artifacts from sinogram using Nghia Vo's approach [B24] (combination of algorithm 3,4,5, and 6).

remove_stripe_based_interpolation(tomo[, ...])

Remove most types of stripe artifacts from sinograms based on interpolation.

tomopy.prep.stripe.remove_all_stripe(tomo, snr=3, la_size=61, sm_size=21, dim=1, ncore=None, nchunk=None)[source]

Remove all types of stripe artifacts from sinogram using Nghia Vo’s approach [B24] (combination of algorithm 3,4,5, and 6).

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • snr (float) – Ratio used to locate large stripes. Greater is less sensitive.

  • la_size (int) – Window size of the median filter to remove large stripes.

  • sm_size (int) – Window size of the median filter to remove small-to-medium stripes.

  • dim ({1, 2}, optional) – Dimension of the window.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_dead_stripe(tomo, snr=3, size=51, norm=True, ncore=None, nchunk=None)[source]

Remove unresponsive and fluctuating stripe artifacts from sinogram using Nghia Vo’s approach [B24] (algorithm 6).

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • snr (float) – Ratio used to detect locations of large stripes. Greater is less sensitive.

  • size (int) – Window size of the median filter.

  • norm (bool, optional) – Remove residual stripes if True.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_large_stripe(tomo, snr=3, size=51, drop_ratio=0.1, norm=True, ncore=None, nchunk=None)[source]

Remove large stripe artifacts from sinogram using Nghia Vo’s approach [B24] (algorithm 5).

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • snr (float) – Ratio used to locate of large stripes. Greater is less sensitive.

  • size (int) – Window size of the median filter.

  • drop_ratio (float, optional) – Ratio of pixels to be dropped, which is used to reduce the false detection of stripes.

  • norm (bool, optional) – Apply normalization if True.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_based_filtering(tomo, sigma=3, size=None, dim=1, ncore=None, nchunk=None)[source]

Remove stripe artifacts from sinogram using Nghia Vo’s approach [B24] (algorithm 2).

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • sigma (float) – Sigma of the Gaussian window which is used to separate the low-pass and high-pass components of the intensity profiles of each column. Recommended values: 3->10.

  • size (int) – Window size of the median filter.

  • dim ({1, 2}, optional) – Dimension of the window.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_based_fitting(tomo, order=3, sigma=(5, 20), ncore=None, nchunk=None)[source]

Remove stripe artifacts from sinogram using Nghia Vo’s approach [B24] (algorithm 1). Suitable for removing low-pass stripes.

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • order (int) – Polynomial fit order. Recommended values: 1-> 5

  • sigma (tuple of 2 floats) – Sigmas of a 2D Gaussian window in x and y direction. Recommended values (3, 20) -> (10, 60).

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_based_interpolation(tomo, snr=3, size=31, drop_ratio=0.1, norm=True, ncore=None, nchunk=None)[source]

Remove most types of stripe artifacts from sinograms based on interpolation. Derived from algorithm 4, 5, and 6 in [B24].

New in version 1.9.

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • snr (float) – Ratio used to segment between useful information and noise.

  • size (int) – Window size of the median filter used to detect stripes.

  • drop_ratio (float, optional) – Ratio of pixels to be dropped, which is used to to reduce the possibility of the false detection of stripes.

  • norm (bool, optional) – Apply normalization if True.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_based_sorting(tomo, size=None, dim=1, ncore=None, nchunk=None)[source]

Remove full and partial stripe artifacts from sinogram using Nghia Vo’s approach [B24] (algorithm 3). Suitable for removing partial stripes.

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • size (int) – Window size of the median filter.

  • dim ({1, 2}, optional) – Dimension of the window.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_fw(tomo, level=None, wname='db5', sigma=2, pad=True, ncore=None, nchunk=None)[source]

Remove horizontal stripes from sinogram using the Fourier-Wavelet (FW) based method [B4].

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • level (int, optional) – Number of discrete wavelet transform levels.

  • wname (str, optional) – Type of the wavelet filter. ‘haar’, ‘db5’, sym5’, etc.

  • sigma (float, optional) – Damping parameter in Fourier space.

  • pad (bool, optional) – If True, extend the size of the sinogram by padding with zeros.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_sf(tomo, size=5, ncore=None, nchunk=None)[source]

Normalize raw projection data using a smoothing filter approach.

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • size (int, optional) – Size of the smoothing filter.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.remove_stripe_ti(tomo, nblock=0, alpha=1.5, ncore=None, nchunk=None)[source]

Remove horizontal stripes from sinogram using Titarenko’s approach [B13].

Parameters
  • tomo (ndarray) – 3D tomographic data.

  • nblock (int, optional) – Number of blocks.

  • alpha (int, optional) – Damping factor.

  • ncore (int, optional) – Number of cores that will be assigned to jobs.

  • nchunk (int, optional) – Chunk size for each core.

Returns

ndarray – Corrected 3D tomographic data.

tomopy.prep.stripe.stripes_detect3d(tomo, size=10, radius=3, ncore=None)[source]

Detect stripes in a 3D array. Usually applied to normalized projection data. The algorithm is presented in the paper [B10].

The method works with full and partial stripes of constant ot varying intensity.

New in version 1.14.

Parameters
  • tomo (ndarray) – 3D tomographic data of float32 data type, preferably in the [0, 1] range, although reasonable deviations accepted (e.g. the result of the normalization and the negative log taken of the raw data). The projection data should be given with [angle, detY(depth), detX(horizontal)] axis orientation. With this orientation, the stripes are features along the angle axis.

  • size (int, optional) – The pixel size of the 1D median filter orthogonal to stripes orientation to minimise false detections. Increase it if you have longer or full stripes in the data.

  • radius (int, optional) – The pixel size of the 3D stencil to calculate the mean ratio between the angular and detX orientations of the detX gradient. The larger values can affect the width of the detected stripe, use 1,2,3 values.

  • ncore (int, optional) – Number of cores that will be assigned to jobs. All cores will be used if unspecified.

Returns

ndarray – Weights in the range of [0, 1] of float32 data type where stripe’s edges are highlighted with the smaller (e.g. < 0.5) values. The weights can be manually thresholded or passed to stripes_mask3d function for further processing and a binary mask generation.

Raises

ValueError – If the tomo is not three dimensional. If the size is invalid.

tomopy.prep.stripe.stripes_mask3d(weights, threshold=0.6, min_stripe_length=20, min_stripe_depth=10, min_stripe_width=5, sensitivity_perc=85.0, ncore=None)[source]

Takes the result of the stripes_detect3d module as an input and generates a binary 3D mask with ones where stripes present.

The method tries to eliminate non-stripe features in data by checking the consistency of weights in three directions. The algorithm is presented in the paper [B10].

New in version 1.14.

Parameters
  • weights (ndarray) – 3D weights array, a result of stripes_detect3d module given in [angles, detY(depth), detX] axis orientation.

  • threshold (float, optional) – Threshold for the given weights. This parameter defines what weights will be considered as potential candidates for stripes. The lower value (< 0.5) will result in only the most prominent stripes in the data. Increase the threshold cautiously because increasing the threshold increases the probability of false detections. The good range to try is between 0.5 and 0.7.

  • min_stripe_length (int, optional) – Minimum length of a stripe in pixels with respect to the “angles” axis. If there are full stripes in the data, then this could be >50% of the size of the the “angles” axis.

  • min_stripe_depth (int, optional) – Minimum depth of a stripe in pixels with respect to the “detY” axis. The stripes do not extend very deep normally in the data. By setting this parameter to the approximate depth of the stripe more false alarms can be removed.

  • min_stripe_width (int, optional) – Minimum width of a stripe in pixels with respect to the “detX” axis. The stripes that close to each other can be merged together with this parameter.

  • sensitivity_perc (float, optional) – The value in the range [0, 100] that controls the strictness of the minimum length, depth and width parameters of a stripe. 0 is less-strict. 100 is more-strict.

  • ncore (int, optional) – Number of cores that will be assigned to jobs. All cores will be used if unspecified.

Returns

ndarray – A binary mask of bool data type with stripes highlighted as True values.

Raises

ValueError – If the input array is not three dimensional. If a min_stripe_length parameter is negative, zero or longer than its corresponding dimension (“angle”) If a min_stripe_depth parameter is negative or longer than its corresponding dimension (“detY”) If a min_stripe_width parameter is negative, zero or longer than its corresponding dimension (“detX”) If a sensitivity_perc parameter doesn’t lie in the (0,100] range