# coding: utf-8 # /*########################################################################## # # Copyright (c) 2015-2018 European Synchrotron Radiation Facility # # Permission is hereby granted, free of charge, to any person obtaining a copy # of this software and associated documentation files (the "Software"), to deal # in the Software without restriction, including without limitation the rights # to use, copy, modify, merge, publish, distribute, sublicense, and/or sell # copies of the Software, and to permit persons to whom the Software is # furnished to do so, subject to the following conditions: # # The above copyright notice and this permission notice shall be included in # all copies or substantial portions of the Software. # # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR # IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, # FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE # AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER # LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, # OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN # THE SOFTWARE. # # ###########################################################################*/ """This module provides median filter function for 1D and 2D arrays. """ __authors__ = ["H. Payno", "J. Kieffer"] __license__ = "MIT" __date__ = "02/05/2017" from cython.parallel import prange cimport cython cimport median_filter import numpy cimport numpy as cnumpy from libcpp cimport bool ctypedef unsigned long uint64 ctypedef unsigned int uint32 ctypedef unsigned short uint16 MODES = {'nearest':0, 'reflect':1, 'mirror':2, 'shrink':3} def medfilt1d(data, kernel_size=3, bool conditional=False, mode='nearest'): """Function computing the median filter of the given input. Behavior at boundaries: the algorithm is reducing the size of the window/kernel for pixels at boundaries (there is no mirroring). :param numpy.ndarray data: the array for which we want to apply the median filter. Should be 1d. :param kernel_size: the dimension of the kernel. :type kernel_size: int :param bool conditional: True if we want to apply a conditional median filtering. :returns: the array with the median value for each pixel. """ return medfilt(data, kernel_size, conditional, mode) def medfilt2d(image, kernel_size=3, bool conditional=False, mode='nearest'): """Function computing the median filter of the given input. Behavior at boundaries: the algorithm is reducing the size of the window/kernel for pixels at boundaries (there is no mirroring). :param numpy.ndarray data: the array for which we want to apply the median filter. Should be 2d. :param kernel_size: the dimension of the kernel. :type kernel_size: For 1D should be an int for 2D should be a tuple or a list of (kernel_height, kernel_width) :param bool conditional: True if we want to apply a conditional median filtering. :returns: the array with the median value for each pixel. """ return medfilt(image, kernel_size, conditional, mode) def medfilt(data, kernel_size=3, bool conditional=False, mode='nearest'): """Function computing the median filter of the given input. Behavior at boundaries: the algorithm is reducing the size of the window/kernel for pixels at boundaries (there is no mirroring). :param numpy.ndarray data: the array for which we want to apply the median filter. Should be 1d or 2d. :param kernel_size: the dimension of the kernel. :type kernel_size: For 1D should be an int for 2D should be a tuple or a list of (kernel_height, kernel_width) :param bool conditional: True if we want to apply a conditional median filtering. :param str mode: the algorithm used to determine how values at borders are determined. :returns: the array with the median value for each pixel. """ if mode not in MODES: err = 'Requested mode %s is unknowed.' % mode raise ValueError(err) reshaped = False if len(data.shape) == 1: data = data.reshape(data.shape[0], 1) reshaped = True elif len(data.shape) > 2: raise ValueError("Invalid data shape. Dimemsion of the arary should be 1 or 2") # simple median filter apply into a 2D buffer output_buffer = numpy.zeros_like(data) check(data, output_buffer) if type(kernel_size) in (tuple, list): if(len(kernel_size) == 1): ker_dim = numpy.array(1, [kernel_size[0]], dtype=numpy.int32) else: ker_dim = numpy.array(kernel_size, dtype=numpy.int32) else: ker_dim = numpy.array([kernel_size, kernel_size], dtype=numpy.int32) if data.dtype == numpy.float64: medfilterfc = _median_filter_float64 elif data.dtype == numpy.float32: medfilterfc = _median_filter_float32 elif data.dtype == numpy.int64: medfilterfc = _median_filter_int64 elif data.dtype == numpy.uint64: medfilterfc = _median_filter_uint64 elif data.dtype == numpy.int32: medfilterfc = _median_filter_int32 elif data.dtype == numpy.uint32: medfilterfc = _median_filter_uint32 elif data.dtype == numpy.int16: medfilterfc = _median_filter_int16 elif data.dtype == numpy.uint16: medfilterfc = _median_filter_uint16 else: raise ValueError("%s type is not managed by the median filter" % data.dtype) medfilterfc(input_buffer=data, output_buffer=output_buffer, kernel_size=ker_dim, conditional=conditional, mode=MODES[mode]) if reshaped: data = data.reshape(data.shape[0]) output_buffer = output_buffer.reshape(data.shape[0]) return output_buffer def check(input_buffer, output_buffer): """Simple check on the two buffers to make sure we can apply the median filter """ if (input_buffer.flags['C_CONTIGUOUS'] is False): raise ValueError(' must be a C_CONTIGUOUS numpy array.') if (output_buffer.flags['C_CONTIGUOUS'] is False): raise ValueError(' must be a C_CONTIGUOUS numpy array.') if not (len(input_buffer.shape) <= 2): raise ValueError(' dimension must mo higher than 2.') if not (len(output_buffer.shape) <= 2): raise ValueError(' dimension must mo higher than 2.') if not(input_buffer.dtype == output_buffer.dtype): raise ValueError('input buffer and output_buffer must be of the same type') if not (input_buffer.shape == output_buffer.shape): raise ValueError('input buffer and output_buffer must be of the same dimension and same dimension') ######### implementations of the include/median_filter.hpp function ############ @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def reflect(int index, int length_max): """find the correct index into [0, length_max-1] for index in reflect mode :param int index: the index to move into [0, length_max-1] in reflect mode :param int length_max: the higher bound limit """ return median_filter.reflect(index, length_max) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def mirror(int index, int length_max): """find the correct index into [0, length_max-1] for index in mirror mode :param int index: the index to move into [0, length_max-1] in mirror mode :param int length_max: the higher bound limit """ return median_filter.mirror(index, length_max) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_float32(float[:, ::1] input_buffer not None, float[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[float]( & input_buffer[0,0], & output_buffer[0,0], & kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_float64(double[:, ::1] input_buffer not None, double[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[double]( & input_buffer[0, 0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_int64(cnumpy.int64_t[:, ::1] input_buffer not None, cnumpy.int64_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[long]( & input_buffer[0,0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_uint64( cnumpy.uint64_t[:, ::1] input_buffer not None, cnumpy.uint64_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[uint64]( & input_buffer[0,0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_int32(cnumpy.int32_t[:, ::1] input_buffer not None, cnumpy.int32_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[int]( & input_buffer[0,0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_uint32(cnumpy.uint32_t[:, ::1] input_buffer not None, cnumpy.uint32_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[uint32]( & input_buffer[0,0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_int16(cnumpy.int16_t[:, ::1] input_buffer not None, cnumpy.int16_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[short]( & input_buffer[0,0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode) @cython.cdivision(True) @cython.boundscheck(False) @cython.wraparound(False) @cython.initializedcheck(False) def _median_filter_uint16( cnumpy.uint16_t[:, ::1] input_buffer not None, cnumpy.uint16_t[:, ::1] output_buffer not None, cnumpy.int32_t[::1] kernel_size not None, bool conditional, int mode): cdef: int y = 0 int image_dim = input_buffer.shape[1] - 1 int[2] buffer_shape, buffer_shape[0] = input_buffer.shape[0] buffer_shape[1] = input_buffer.shape[1] for y in prange(input_buffer.shape[0], nogil=True): median_filter.median_filter[uint16]( & input_buffer[0, 0], & output_buffer[0, 0], &kernel_size[0], buffer_shape, y, 0, image_dim, conditional, mode)