B
    q\y                 @   s   d dl Z d dlZd dlmZ d dlmZ d dlmZ yd dl	Z	dZ
d dlmZ W n ek
rh   dZ
Y nX dd	d
gZdd ZdddZdddZd ddZG dd dZedddd!dd	Zedddd"dd
ZdS )#    N)
isiterable)deprecated_renamed_argument)AstropyUserWarningT)QuantityF	SigmaClip
sigma_clipsigma_clipped_statsc                s`   t  } t fddt| jD 7  tt| j}t|  |}|d|j|d  }|S )z
    Bottleneck can only take integer axis, not tuple, so this function
    takes all the axes to be operated on and combines them into the
    first dimension of the array so that we can then use axis=0
    c             3   s   | ]}| kr|V  qd S )N ).0i)axisr	   ;lib/python3.7/site-packages/astropy/stats/sigma_clipping.py	<genexpr>"   s    z)_move_tuple_axes_first.<locals>.<genexpr>)N)lentuplerangendimnpZmoveaxisreshapeshape)arrayr   ZnaxisZdestinationZ	array_newr	   )r   r   _move_tuple_axes_first   s     r   c             C   sJ   t |trt| |d} d}t | tr8| tj| |dS tj| |dS dS )z3Bottleneck nanmean function that handle tuple axis.)r   r   N)
isinstancer   r   r   __array_wrap__
bottlenecknanmean)r   r   r	   r	   r   _nanmean1   s    

r   c             C   sJ   t |trt| |d} d}t | tr8| tj| |dS tj| |dS dS )z5Bottleneck nanmedian function that handle tuple axis.)r   r   N)r   r   r   r   r   r   	nanmedian)r   r   r	   r	   r   
_nanmedian>   s    

r   c             C   sN   t |trt| |d} d}t | tr:| tj| ||dS tj| ||dS dS )z2Bottleneck nanstd function that handle tuple axis.)r   r   )r   ddofN)r   r   r   r   r   r   nanstd)r   r   r    r	   r	   r   _nanstdK   s    


r"   c               @   sn   e Zd ZdZeddddd
dZdd Zdd Zdd Zdd Z	dddZ
d ddZd!ddZd"ddZdS )#r   aJ  
    Class to perform sigma clipping.

    The data will be iterated over, each time rejecting values that are
    less or more than a specified number of standard deviations from a
    center value.

    Clipped (rejected) pixels are those where::

        data < cenfunc(data [,axis=int]) - (sigma_lower * stdfunc(data [,axis=int]))
        data > cenfunc(data [,axis=int]) + (sigma_upper * stdfunc(data [,axis=int]))

    Invalid data values (i.e. NaN or inf) are automatically clipped.

    For a functional interface to sigma clipping, see
    :func:`sigma_clip`.

    .. note::
        `scipy.stats.sigmaclip
        <https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.sigmaclip.html>`_
        provides a subset of the functionality in this class.  Also, its
        input data cannot be a masked array and it does not handle data
        that contains invalid values (i.e.  NaN or inf).  Also note that
        it uses the mean as the centering function.

        If your data is a `~numpy.ndarray` with no invalid values and
        you want to use the mean as the centering function with
        ``axis=None`` and iterate to convergence, then
        `scipy.stats.sigmaclip` is ~25-30% faster than the equivalent
        settings here (``s = SigmaClip(cenfunc='mean', maxiters=None);
        s(data, axis=None)``).

    Parameters
    ----------
    sigma : float, optional
        The number of standard deviations to use for both the lower and
        upper clipping limit.  These limits are overridden by
        ``sigma_lower`` and ``sigma_upper``, if input.  The default is
        3.

    sigma_lower : float or `None`, optional
        The number of standard deviations to use as the lower bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    sigma_upper : float or `None`, optional
        The number of standard deviations to use as the upper bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    maxiters : int or `None`, optional
        The maximum number of sigma-clipping iterations to perform or
        `None` to clip until convergence is achieved (i.e., iterate
        until the last iteration clips nothing).  If convergence is
        achieved prior to ``maxiters`` iterations, the clipping
        iterations will stop.  The default is 5.

    cenfunc : {'median', 'mean'} or callable, optional
        The statistic or callable function/object used to compute the
        center value for the clipping.  If set to ``'median'`` or
        ``'mean'`` then having the optional `bottleneck`_ package
        installed will result in the best performance.  If using a
        callable function/object and the ``axis`` keyword is used, then
        it must be callable that can ignore NaNs (e.g. `numpy.nanmean`)
        and has an ``axis`` keyword to return an array with axis
        dimension(s) removed.  The default is ``'median'``.

        .. _bottleneck:  https://github.com/kwgoodman/bottleneck

    stdfunc : {'std'} or callable, optional
        The statistic or callable function/object used to compute the
        standard deviation about the center value.  If set to ``'std'``
        then having the optional `bottleneck`_ package installed will
        result in the best performance.  If using a callable
        function/object and the ``axis`` keyword is used, then it must
        be callable that can ignore NaNs (e.g. `numpy.nanstd`) and has
        an ``axis`` keyword to return an array with axis dimension(s)
        removed.  The default is ``'std'``.

    See Also
    --------
    sigma_clip, sigma_clipped_stats

    Examples
    --------
    This example uses a data array of random variates from a Gaussian
    distribution.  We clip all points that are more than 2 sample
    standard deviations from the median.  The result is a masked array,
    where the mask is `True` for clipped data::

        >>> from astropy.stats import SigmaClip
        >>> from numpy.random import randn
        >>> randvar = randn(10000)
        >>> sigclip = SigmaClip(sigma=2, maxiters=5)
        >>> filtered_data = sigclip(randvar)

    This example clips all points that are more than 3 sigma relative to
    the sample *mean*, clips until convergence, returns an unmasked
    `~numpy.ndarray`, and modifies the data in-place::

        >>> from astropy.stats import SigmaClip
        >>> from numpy.random import randn
        >>> from numpy import mean
        >>> randvar = randn(10000)
        >>> sigclip = SigmaClip(sigma=3, maxiters=None, cenfunc='mean')
        >>> filtered_data = sigclip(randvar, masked=False, copy=False)

    This example sigma clips along one axis::

        >>> from astropy.stats import SigmaClip
        >>> from numpy.random import normal
        >>> from numpy import arange, diag, ones
        >>> data = arange(5) + normal(0., 0.05, (5, 5)) + diag(ones(5))
        >>> sigclip = SigmaClip(sigma=2.3)
        >>> filtered_data = sigclip(data, axis=0)

    Note that along the other axis, no points would be clipped, as the
    standard deviation is higher.
    itersmaxitersz3.1      @N   medianstdc             C   sB   || _ |p|| _|p|| _|p"tj| _| || _| || _	d S )N)
sigmasigma_lowersigma_upperr   infr$   _parse_cenfunccenfunc_parse_stdfuncstdfunc)selfr)   r*   r+   r$   r.   r0   r	   r	   r   __init__   s    

zSigmaClip.__init__c             C   s    d | j| j| j| j| j| jS )Nz^SigmaClip(sigma={0}, sigma_lower={1}, sigma_upper={2}, maxiters={3}, cenfunc={4}, stdfunc={5}))formatr)   r*   r+   r$   r.   r0   )r1   r	   r	   r   __repr__   s    zSigmaClip.__repr__c          	   C   sR   d| j j d g}ddddddg}x$|D ]}|d	|t| | q(W d
|S )N<>r)   r*   r+   r$   r.   r0   z    {0}: {1}
)	__class____name__appendr3   getattrjoin)r1   linesZattrsattrr	   r	   r   __str__   s    

zSigmaClip.__str__c             C   sP   t |trL|dkr$trt}qLtj}n(|dkr>tr6t}qLtj}ntd	||S )Nr'   meanz{} is an invalid cenfunc.)
r   strHAS_BOTTLENECKr   r   r   r   r   
ValueErrorr3   )r1   r.   r	   r	   r   r-      s    
zSigmaClip._parse_cenfuncc             C   s4   t |tr0|dkr td|tr*t}ntj}|S )Nr(   z{} is an invalid stdfunc.)r   rA   rC   r3   rB   r"   r   r!   )r1   r0   r	   r	   r   r/      s    
zSigmaClip._parse_stdfuncc          	   C   sj   t  X t jdtd | j||d| _| j||d}| j|| j  | _|  j|| j	 7  _W d Q R X d S )Nignore)category)r   )
warningscatch_warningssimplefilterRuntimeWarningr.   
_max_valuer0   r*   
_min_valuer+   )r1   datar   r(   r	   r	   r   _compute_bounds
  s    
zSigmaClip._compute_boundsTFc       
   	   C   s  |  }t|tjjr$|j|j  }t|}t| rN|| }t	
dt d}d}xT|dkr|| jk r|d7 }|j}	| j|dd ||| jk|| jk@  }|	|j }qXW || _|rtjj||d}tjdd( | jt|| jk || jkO  _W dQ R X |r|| j| jfS |S dS )	z
        Sigma clip the data when ``axis`` is None.

        In this simple case, we remove clipped elements from the
        flattened array during each iteration.
        zTInput data contains invalid values (NaNs or infs), which were automatically clipped.   r   N)r   )copyrD   )invalid)Zravelr   r   maMaskedArrayrL   maskisfiniteanyrF   warnr   r$   sizerM   rK   rJ   _niterationsmasked_invaliderrstate
logical_or)
r1   rL   maskedreturn_boundsrO   filtered_dataZ	good_masknchanged	iterationrW   r	   r	   r   _sigmaclip_noaxis  s4    	

zSigmaClip._sigmaclip_noaxisc          	      s  | tt }t|r6tj|< tdt t	tj
jrbtj
 ttjt sp f tfdd D  t fddtjD }d}d}	x|dkrZ|	| jk rZ|	d7 }	tt}
| j d t| js| j|| _| j|| _tjdd	  tj| jk | jkB < W d
Q R X |
tt }qW |	| _|r|r|tj
nNtjdd	< tj
j|dd}tj
jt|| jk || jk|ddW d
Q R X |rވ| j| jfS S d
S )z
        Sigma clip the data when ``axis`` is specified.

        In this case, we replace clipped values with NaNs as placeholder
        values.
        zTInput data contains invalid values (NaNs or infs), which were automatically clipped.c             3   s$   | ]}|d k r j | n|V  qdS )r   N)r   )r
   n)r^   r	   r   r   b  s    z0SigmaClip._sigmaclip_withaxis.<locals>.<genexpr>c             3   s"   | ]\}}| krd n|V  qdS )rN   Nr	   )r
   ZdimrW   )r   r	   r   r   f  s   rN   r   )r   rD   )rP   NF)rO   )Zastypefloatr   rT   rU   nanrF   rV   r   r   rQ   rR   rY   Zfilledr   r   	enumerater   r$   Zcount_nonzeroZisnanrM   ZisscalarrK   r   rJ   rZ   rX   Zmasked_wherer[   )r1   rL   r   r\   r]   rO   Zbad_maskZmshaper_   r`   Zn_nanoutr	   )r   r^   r   _sigmaclip_withaxisF  sL    



$
zSigmaClip._sigmaclip_withaxisc             C   sf   t |}|jdkr|S t|t jjr4|j r4|S |dkrN| j||||dS | j	|||||dS dS )a	  
        Perform sigma clipping on the provided data.

        Parameters
        ----------
        data : array-like or `~numpy.ma.MaskedArray`
            The data to be sigma clipped.

        axis : `None` or int or tuple of int, optional
            The axis or axes along which to sigma clip the data.  If `None`,
            then the flattened data will be used.  ``axis`` is passed
            to the ``cenfunc`` and ``stdfunc``.  The default is `None`.

        masked : bool, optional
            If `True`, then a `~numpy.ma.MaskedArray` is returned, where
            the mask is `True` for clipped values.  If `False`, then a
            `~numpy.ndarray` and the minimum and maximum clipping
            thresholds are returned.  The default is `True`.

        return_bounds : bool, optional
            If `True`, then the minimum and maximum clipping bounds are
            also returned.

        copy : bool, optional
            If `True`, then the ``data`` array will be copied.  If
            `False` and ``masked=True``, then the returned masked array
            data will contain the same array as the input ``data`` (if
            ``data`` is a `~numpy.ndarray` or `~numpy.ma.MaskedArray`).
            The default is `True`.

        Returns
        -------
        result : flexible
            If ``masked=True``, then a `~numpy.ma.MaskedArray` is
            returned, where the mask is `True` for clipped values.  If
            ``masked=False``, then a `~numpy.ndarray` is returned.

            If ``return_bounds=True``, then in addition to the (masked)
            array above, the minimum and maximum clipping bounds are
            returned.

            If ``masked=False`` and ``axis=None``, then the output array
            is a flattened 1D `~numpy.ndarray` where the clipped values
            have been removed.  If ``return_bounds=True`` then the
            returned minimum and maximum thresholds are scalars.

            If ``masked=False`` and ``axis`` is specified, then the
            output `~numpy.ndarray` will have the same shape as the
            input ``data`` and contain ``np.nan`` where values were
            clipped.  If ``return_bounds=True`` then the returned
            minimum and maximum clipping thresholds will be be
            `~numpy.ndarray`\s.
        r   N)r\   r]   rO   )r   r\   r]   rO   )
r   Z
asanyarrayrW   r   rQ   rR   rS   allra   rg   )r1   rL   r   r\   r]   rO   r	   r	   r   __call__  s    8


zSigmaClip.__call__)r%   NNr&   r'   r(   )N)TFT)NTFT)NTFT)r9   
__module____qualname____doc__r   r2   r4   r?   r-   r/   rM   ra   rg   ri   r	   r	   r	   r   r   Y   s   w
 	

 
1 
G r#   r$   z3.1   r&   r'   r(   c             C   s&   t ||||||d}|| |||	|
dS )a8  
    Perform sigma-clipping on the provided data.

    The data will be iterated over, each time rejecting values that are
    less or more than a specified number of standard deviations from a
    center value.

    Clipped (rejected) pixels are those where::

        data < cenfunc(data [,axis=int]) - (sigma_lower * stdfunc(data [,axis=int]))
        data > cenfunc(data [,axis=int]) + (sigma_upper * stdfunc(data [,axis=int]))

    Invalid data values (i.e. NaN or inf) are automatically clipped.

    For an object-oriented interface to sigma clipping, see
    :class:`SigmaClip`.

    .. note::
        `scipy.stats.sigmaclip
        <https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.sigmaclip.html>`_
        provides a subset of the functionality in this class.  Also, its
        input data cannot be a masked array and it does not handle data
        that contains invalid values (i.e.  NaN or inf).  Also note that
        it uses the mean as the centering function.

        If your data is a `~numpy.ndarray` with no invalid values and
        you want to use the mean as the centering function with
        ``axis=None`` and iterate to convergence, then
        `scipy.stats.sigmaclip` is ~25-30% faster than the equivalent
        settings here (``sigma_clip(data, cenfunc='mean', maxiters=None,
        axis=None)``).

    Parameters
    ----------
    data : array-like or `~numpy.ma.MaskedArray`
        The data to be sigma clipped.

    sigma : float, optional
        The number of standard deviations to use for both the lower and
        upper clipping limit.  These limits are overridden by
        ``sigma_lower`` and ``sigma_upper``, if input.  The default is
        3.

    sigma_lower : float or `None`, optional
        The number of standard deviations to use as the lower bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    sigma_upper : float or `None`, optional
        The number of standard deviations to use as the upper bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    maxiters : int or `None`, optional
        The maximum number of sigma-clipping iterations to perform or
        `None` to clip until convergence is achieved (i.e., iterate
        until the last iteration clips nothing).  If convergence is
        achieved prior to ``maxiters`` iterations, the clipping
        iterations will stop.  The default is 5.

    cenfunc : {'median', 'mean'} or callable, optional
        The statistic or callable function/object used to compute the
        center value for the clipping.  If set to ``'median'`` or
        ``'mean'`` then having the optional `bottleneck`_ package
        installed will result in the best performance.  If using a
        callable function/object and the ``axis`` keyword is used, then
        it must be callable that can ignore NaNs (e.g. `numpy.nanmean`)
        and has an ``axis`` keyword to return an array with axis
        dimension(s) removed.  The default is ``'median'``.

        .. _bottleneck:  https://github.com/kwgoodman/bottleneck

    stdfunc : {'std'} or callable, optional
        The statistic or callable function/object used to compute the
        standard deviation about the center value.  If set to ``'std'``
        then having the optional `bottleneck`_ package installed will
        result in the best performance.  If using a callable
        function/object and the ``axis`` keyword is used, then it must
        be callable that can ignore NaNs (e.g. `numpy.nanstd`) and has
        an ``axis`` keyword to return an array with axis dimension(s)
        removed.  The default is ``'std'``.

    axis : `None` or int or tuple of int, optional
        The axis or axes along which to sigma clip the data.  If `None`,
        then the flattened data will be used.  ``axis`` is passed to the
        ``cenfunc`` and ``stdfunc``.  The default is `None`.

    masked : bool, optional
        If `True`, then a `~numpy.ma.MaskedArray` is returned, where the
        mask is `True` for clipped values.  If `False`, then a
        `~numpy.ndarray` and the minimum and maximum clipping thresholds
        are returned.  The default is `True`.

    return_bounds : bool, optional
        If `True`, then the minimum and maximum clipping bounds are also
        returned.

    copy : bool, optional
        If `True`, then the ``data`` array will be copied.  If `False`
        and ``masked=True``, then the returned masked array data will
        contain the same array as the input ``data`` (if ``data`` is a
        `~numpy.ndarray` or `~numpy.ma.MaskedArray`).  The default is
        `True`.

    Returns
    -------
    result : flexible
        If ``masked=True``, then a `~numpy.ma.MaskedArray` is returned,
        where the mask is `True` for clipped values.  If
        ``masked=False``, then a `~numpy.ndarray` is returned.

        If ``return_bounds=True``, then in addition to the (masked)
        array above, the minimum and maximum clipping bounds are
        returned.

        If ``masked=False`` and ``axis=None``, then the output array is
        a flattened 1D `~numpy.ndarray` where the clipped values have
        been removed.  If ``return_bounds=True`` then the returned
        minimum and maximum thresholds are scalars.

        If ``masked=False`` and ``axis`` is specified, then the output
        `~numpy.ndarray` will have the same shape as the input ``data``
        and contain ``np.nan`` where values were clipped.  If
        ``return_bounds=True`` then the returned minimum and maximum
        clipping thresholds will be be `~numpy.ndarray`\s.

    See Also
    --------
    SigmaClip, sigma_clipped_stats

    Examples
    --------
    This example uses a data array of random variates from a Gaussian
    distribution.  We clip all points that are more than 2 sample
    standard deviations from the median.  The result is a masked array,
    where the mask is `True` for clipped data::

        >>> from astropy.stats import sigma_clip
        >>> from numpy.random import randn
        >>> randvar = randn(10000)
        >>> filtered_data = sigma_clip(randvar, sigma=2, maxiters=5)

    This example clips all points that are more than 3 sigma relative to
    the sample *mean*, clips until convergence, returns an unmasked
    `~numpy.ndarray`, and does not copy the data::

        >>> from astropy.stats import sigma_clip
        >>> from numpy.random import randn
        >>> from numpy import mean
        >>> randvar = randn(10000)
        >>> filtered_data = sigma_clip(randvar, sigma=3, maxiters=None,
        ...                            cenfunc=mean, masked=False, copy=False)

    This example sigma clips along one axis::

        >>> from astropy.stats import sigma_clip
        >>> from numpy.random import normal
        >>> from numpy import arange, diag, ones
        >>> data = arange(5) + normal(0., 0.05, (5, 5)) + diag(ones(5))
        >>> filtered_data = sigma_clip(data, sigma=2.3, axis=0)

    Note that along the other axis, no points would be clipped, as the
    standard deviation is higher.
    )r)   r*   r+   r$   r.   r0   )r   r\   r]   rO   )r   )rL   r)   r*   r+   r$   r.   r0   r   r\   r]   rO   sigclipr	   r	   r   r     s     *
      @c             C   s   |dk	rt j| |} |dk	r,t j| |} t||||||d}|| |
dddd}tr~t||
d}t||
d}t||	|
d}n,t j	||
d}t j
||
d}t j||	|
d}|||fS )aN  
    Calculate sigma-clipped statistics on the provided data.

    Parameters
    ----------
    data : array-like or `~numpy.ma.MaskedArray`
        Data array or object that can be converted to an array.

    mask : `numpy.ndarray` (bool), optional
        A boolean mask with the same shape as ``data``, where a `True`
        value indicates the corresponding element of ``data`` is masked.
        Masked pixels are excluded when computing the statistics.

    mask_value : float, optional
        A data value (e.g., ``0.0``) that is ignored when computing the
        statistics.  ``mask_value`` will be masked in addition to any
        input ``mask``.

    sigma : float, optional
        The number of standard deviations to use for both the lower and
        upper clipping limit.  These limits are overridden by
        ``sigma_lower`` and ``sigma_upper``, if input.  The default is
        3.

    sigma_lower : float or `None`, optional
        The number of standard deviations to use as the lower bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    sigma_upper : float or `None`, optional
        The number of standard deviations to use as the upper bound for
        the clipping limit.  If `None` then the value of ``sigma`` is
        used.  The default is `None`.

    maxiters : int or `None`, optional
        The maximum number of sigma-clipping iterations to perform or
        `None` to clip until convergence is achieved (i.e., iterate
        until the last iteration clips nothing).  If convergence is
        achieved prior to ``maxiters`` iterations, the clipping
        iterations will stop.  The default is 5.

    cenfunc : {'median', 'mean'} or callable, optional
        The statistic or callable function/object used to compute the
        center value for the clipping.  If set to ``'median'`` or
        ``'mean'`` then having the optional `bottleneck`_ package
        installed will result in the best performance.  If using a
        callable function/object and the ``axis`` keyword is used, then
        it must be callable that can ignore NaNs (e.g. `numpy.nanmean`)
        and has an ``axis`` keyword to return an array with axis
        dimension(s) removed.  The default is ``'median'``.

        .. _bottleneck:  https://github.com/kwgoodman/bottleneck

    stdfunc : {'std'} or callable, optional
        The statistic or callable function/object used to compute the
        standard deviation about the center value.  If set to ``'std'``
        then having the optional `bottleneck`_ package installed will
        result in the best performance.  If using a callable
        function/object and the ``axis`` keyword is used, then it must
        be callable that can ignore NaNs (e.g. `numpy.nanstd`) and has
        an ``axis`` keyword to return an array with axis dimension(s)
        removed.  The default is ``'std'``.

    std_ddof : int, optional
        The delta degrees of freedom for the standard deviation
        calculation.  The divisor used in the calculation is ``N -
        std_ddof``, where ``N`` represents the number of elements.  The
        default is 0.

    axis : `None` or int or tuple of int, optional
        The axis or axes along which to sigma clip the data.  If `None`,
        then the flattened data will be used.  ``axis`` is passed
        to the ``cenfunc`` and ``stdfunc``.  The default is `None`.

    Returns
    -------
    mean, median, stddev : float
        The mean, median, and standard deviation of the sigma-clipped
        data.

    See Also
    --------
    SigmaClip, sigma_clip
    N)r)   r*   r+   r$   r.   r0   F)r   r\   r]   rO   )r   )r    r   )r   rQ   rR   Zmasked_valuesr   rB   r   r   r"   r   r   r!   )rL   rS   Z
mask_valuer)   r*   r+   r$   r.   r0   Zstd_ddofr   rn   Zdata_clippedr@   r'   r(   r	   r	   r   r     s"    Z

)N)N)Nr   )
rm   NNr&   r'   r(   NTFT)
NNro   NNr&   r'   r(   r   N)rF   Znumpyr   Zastropy.utilsr   Zastropy.utils.decoratorsr   Zastropy.utils.exceptionsr   r   rB   Zastropy.unitsr   ImportError__all__r   r   r   r"   r   r   r   r	   r	   r	   r   <module>   s:   




   
   /
   