def sort_groups(self, sort_using='mean', sort_method='no'):
        """
        Sorts and rearranges the submatrices according to the
        sorting method given.
        """
        if sort_method == 'no':
            return

        # compute the row average:
        if sort_using == 'region_length':
            matrix_avgs = np.array([x['end'] - x['start']
                                   for x in self.regions])
        else:
            matrix_avgs = np.__getattribute__(sort_using)(
                self.matrix, axis=1)

        # order per group
        _sorted_regions = []
        _sorted_matrix = []
        for idx in range(len(self.group_labels)):
            start = self.group_boundaries[idx]
            end = self.group_boundaries[idx + 1]
            order = matrix_avgs[start:end].argsort()
            if sort_method == 'descend':
                order = order[::-1]
            _sorted_matrix.append(self.matrix[start:end, :][order, :])
            # sort the regions
            _reg = self.regions[start:end]
            for idx in order:
                _sorted_regions.append(_reg[idx])

        self.matrix = np.vstack(_sorted_matrix)
        self.regions = _sorted_regions
        self.set_sorting_method(sort_method, sort_using)

    def sortMatrix(self, sort_using="mean", sort_method="no"):
        # sort the matrix using the average of values per row
        self.matrixAvgsDict = OrderedDict()
        self.lengthDict = OrderedDict()
        for label in self.matrixDict.keys():
            if sort_method != "no":
                if sort_using == "region_length":
                    matrixAvgs = np.array([x["end"] - x["start"] for x in self.regionsDict[label]])
                    b = self.parameters["upstream"] / self.parameters["bin size"]
                    # for plotting I add the upstream
                    # distance
                    self.lengthDict[label] = b + (matrixAvgs / self.parameters["bin size"])
                else:
                    matrixAvgs = np.__getattribute__(sort_using)(self.matrixDict[label], axis=1)
                SS = matrixAvgs.argsort()

                if sort_method == "descend":
                    SS = SS[::-1]
                self.matrixDict[label] = self.matrixDict[label][SS, :]
                self.regionsDict[label] = self.regionsDict[label][SS]
                self.matrixAvgsDict[label] = matrixAvgs[SS]
            try:
                self.lengthDict[label] = self.lengthDict[label][SS]
            except:
                self.lengthDict[label] = None

def add_sniff_events(h5, overwrite=True, plot=False, centering_fn='mean', *args, **kwargs):
    """

    :param h5:
    :param overwrite:
    :param plot:
    :param centering_fn: (optional) Str, ('mean', 'median'. This function name will be used to compute the offset of the
    sniff stream. Other usual value would be 'median'.
    :param args:
    :param kwargs:
    :return:
    """

    assert isinstance(h5, tb.File)
    events_group = h5.root.events
    try:
        sniff_stream_tbarray = h5.root.streams.sniff
    except tb.NoSuchNodeError:
        print('Warning, no sniff stream found.')
        return

    try:
        sniff_events_group = h5.create_group('/events', 'sniff')

    except tb.NodeError as e:
        if overwrite:
            h5.remove_node('/events/sniff', recursive=True)
            h5.flush()
            sniff_events_group = h5.create_group('/events', 'sniff')
        else:
            raise SniffExistExemption
    logging.info('Creating sniff events.')
    sniff = sniff_stream_tbarray[:, 0]
    sniff -= np.__getattribute__(centering_fn)(sniff)  # subtract the mean or median from the sniff signal.
    fs_stream = sniff_stream_tbarray._v_attrs['sample_rate_Hz']
    fs_events = events_group._v_attrs['sample_rate_Hz']
    sniff_events_array = find_inhalations_simple(sniff, fs_stream)
    sniff_events_array *= (fs_events/fs_stream)
    sniff_events_array = sniff_events_array.astype(np.int)
    sniff_events = h5.create_carray(sniff_events_group, 'inh_events', obj=sniff_events_array)
    sniff_stats = basic_sniff_stats(sniff_events_array, sniff, h5, fs_stream, fs_events)

    sniff_log = np.zeros(sniff.size, dtype=np.bool)
    for ev in sniff_events:
        ev2 = ev / 16
        sniff_log[ev2[0]:ev2[1]] = True

    m = np.max(sniff[:10000])
    plt.plot(sniff[:10000])
    plt.plot(sniff_log[:10000] * m)
    plt.grid()
    # plt.plot(plt.xlim(), [sniff.mean()]*2)
    plt.show()



    return

 def myAverage(valuesArray, avgType="mean"):
     """
     computes the mean, median, etc but only for those values
     that are not Nan
     """
     valuesArray = np.ma.masked_invalid(valuesArray)
     avg = np.__getattribute__(avgType)(valuesArray)
     if isinstance(avg, np.ma.core.MaskedConstant):
         return np.nan
     else:
         return avg

 def __getattr__(self,attr):
     """ Allow numpy overloading, called if something else is broken
         A bit hackish, but avoids poluting the namespace for the Chebfun """
     if attr in NP_OVERLOAD:
         ufunc = np.__getattribute__(attr)
         def func():
             return self._new_func( lambda x: ufunc( self._eval(x) ) )
         func.__name__ = attr
         func.__doc__  = "wraps numpy ufunc: {}".format(attr)
         return func
     return self.__getattribute__(attr)

def register_numpy_types():
    """Register the AsIs adapter for following types from numpy:
      - numpy.int8
      - numpy.int16
      - numpy.int32
      - numpy.int64

      - numpy.float16
      - numpy.float32
      - numpy.float64
      - numpy.float128
    """
    for typ in ['int8', 'int16', 'int32', 'int64',
                'float16', 'float32', 'float64', 'float128']:
        register_adapter(np.__getattribute__(typ), AsIs)

def compute_statistic(statistic, array):
    """Compute the statistic used in a Dakota response function.

    Parameters
    ----------
    statistic : str
      A string with the name of the statistic to compute ('mean',
      'median', etc.).
    array : array_like
      An array data structure, such as a numpy array.

    Returns
    -------
    float
      The value of the computed statistic.

    """
    import numpy as np
    return np.__getattribute__(statistic)(array)

def print_values_for_attributes_for_keys(dataset, attributes, keys, index=0):
    if type(attributes) is not list:
        attributes = [attributes]
    if type(index) is not list:
        index = [index]
    for key in keys:
        trajec = dataset.trajecs[key]
        print_str = key
        for a, attr in enumerate(attributes):
            attribute = trajec.__getattribute__(attr)
            if type(attribute) in [str, float, int, long, np.float64]:
                print_str += ' -- ' + attr + ': ' + str(attribute)
            else:
                if type(index[a]) is int:
                    print_str += ' -- ' + attr + ': ' + str(attribute[index[a]])
                elif type(index[a]) is str:
                    print_str += ' -- ' + attr + ': ' + str(np.__getattribute__(index[a])(attribute))
                else:
                    print_str += ' -- ' + attr + ': ' + str(attribute)
        print print_str

def test_add_sub_overload(ccd_data, operand, expect_failure, with_uncertainty,
                          operation, affects_uncertainty):
    if with_uncertainty:
        ccd_data.uncertainty = StdDevUncertainty(np.ones_like(ccd_data))
    method = ccd_data.__getattribute__(operation)
    np_method = np.__getattribute__(operation)
    if expect_failure:
        with pytest.raises(expect_failure):
            result = method(operand)
        return
    else:
        result = method(operand)
    assert result is not ccd_data
    assert isinstance(result, CCDData)
    assert (result.uncertainty is None or
            isinstance(result.uncertainty, StdDevUncertainty))
    try:
        op_value = operand.value
    except AttributeError:
        op_value = operand

    np.testing.assert_array_equal(result.data,
                                  np_method(ccd_data.data, op_value))
    if with_uncertainty:
        if affects_uncertainty:
            np.testing.assert_array_equal(result.uncertainty.array,
                                          np_method(ccd_data.uncertainty.array,
                                                    op_value))
        else:
            np.testing.assert_array_equal(result.uncertainty.array,
                                          ccd_data.uncertainty.array)
    else:
        assert result.uncertainty is None

    if isinstance(operand, u.Quantity):
        assert (result.unit == ccd_data.unit and result.unit == operand.unit)
    else:
        assert result.unit == ccd_data.unit

def test_mult_div_overload(ccd_data, operand, with_uncertainty,
                           operation, affects_uncertainty):
    if with_uncertainty:
        ccd_data.uncertainty = StdDevUncertainty(np.ones_like(ccd_data))
    method = ccd_data.__getattribute__(operation)
    np_method = np.__getattribute__(operation)
    result = method(operand)
    assert result is not ccd_data
    assert isinstance(result, CCDData)
    assert (result.uncertainty is None or
            isinstance(result.uncertainty, StdDevUncertainty))
    try:
        op_value = operand.value
    except AttributeError:
        op_value = operand

    np.testing.assert_array_equal(result.data,
                                  np_method(ccd_data.data, op_value))
    if with_uncertainty:
        if affects_uncertainty:
            np.testing.assert_array_equal(result.uncertainty.array,
                                          np_method(ccd_data.uncertainty.array,
                                                    op_value))
        else:
            np.testing.assert_array_equal(result.uncertainty.array,
                                          ccd_data.uncertainty.array)
    else:
        assert result.uncertainty is None

    if isinstance(operand, u.Quantity):
        # Need the "1 *" below to force arguments to be Quantity to work around
        # astropy/astropy#2377
        expected_unit = np_method(1 * ccd_data.unit, 1 * operand.unit).unit
        assert result.unit == expected_unit
    else:
        assert result.unit == ccd_data.unit

 def matrixAvg(matrix, avgType="mean"):
     matrix = np.ma.masked_invalid(matrix)
     return np.__getattribute__(avgType)(matrix, axis=0)

    def save_tabulated_values(self, file_handle, reference_point_label='TSS', start_label='TSS', end_label='TES', averagetype='mean'):
        """
        Saves the values averaged by col using the avg_type
        given

        Args:
            file_handle: file name to save the file
            reference_point_label: Name of the reference point label
            start_label: Name of the star label
            end_label: Name of the end label
            averagetype: average type (e.g. mean, median, std)

        """
        #  get X labels
        w = self.parameters['bin size']
        b = self.parameters['upstream']
        a = self.parameters['downstream']
        c = self.parameters.get('unscaled 5 prime', 0)
        d = self.parameters.get('unscaled 3 prime', 0)
        m = self.parameters['body']

        if b < 1e5:
            quotient = 1000
            symbol = 'Kb'
        else:
            quotient = 1e6
            symbol = 'Mb'

        if m == 0:
            xticks = [(k / w) for k in [w, b, b + a]]
            xtickslabel = ['{0:.1f}{1}'.format(-(float(b) / quotient), symbol), reference_point_label,
                           '{0:.1f}{1}'.format(float(a) / quotient, symbol)]

        else:
            xticks_values = [w]
            xtickslabel = []

            # only if upstream region is set, add a x tick
            if b > 0:
                xticks_values.append(b)
                xtickslabel.append('{0:.1f}{1}'.format(-(float(b) / quotient), symbol))

            xtickslabel.append(start_label)

            if c > 0:
                xticks_values.append(b + c)
                xtickslabel.append("")

            if d > 0:
                xticks_values.append(b + c + m)
                xtickslabel.append("")

            xticks_values.append(b + c + m + d)
            xtickslabel.append(end_label)

            if a > 0:
                xticks_values.append(b + c + m + d + a)
                xtickslabel.append('{0:.1f}{1}'.format(float(a) / quotient, symbol))

            xticks = [(k / w) for k in xticks_values]
        x_axis = np.arange(xticks[-1]) + 1
        labs = []
        for x_value in x_axis:
            if x_value in xticks:
                labs.append(xtickslabel[xticks.index(x_value)])
            else:
                labs.append("")

        with open(file_handle, 'w') as fh:
            # write labels
            fh.write("bin labels\t\t{}\n".format("\t".join(labs)))
            fh.write('bins\t\t{}\n'.format("\t".join([str(x) for x in x_axis])))

            for sample_idx in range(self.matrix.get_num_samples()):
                for group_idx in range(self.matrix.get_num_groups()):
                    sub_matrix = self.matrix.get_matrix(group_idx, sample_idx)
                    values = [str(x) for x in np.__getattribute__(averagetype)(sub_matrix['matrix'], axis=0)]
                    fh.write("{}\t{}\t{}\n".format(sub_matrix['sample'], sub_matrix['group'], "\t".join(values)))

    def plot_heatmap(self):
        matrix_flatten = None
        if self.y_min is None:
            matrix_flatten = self.hm.matrix.flatten()
            # try to avoid outliers by using np.percentile
            self.y_min = np.percentile(matrix_flatten, 1.0)
            if np.isnan(self.y_min):
                self.y_min = None

        if self.y_max is None:
            if matrix_flatten is None:
                matrix_flatten = self.hm.matrix.flatten()
            # try to avoid outliers by using np.percentile
            self.y_max = np.percentile(matrix_flatten, 98.0)
            if np.isnan(self.y_max):
                self.y_max = None

        ax_list = []
        # turn off y ticks

        for plot in range(self.numplots):
            labels = []
            col = plot % self.plots_per_row
            row = int(plot / self.plots_per_row)

            # split the ax to make room for the colorbar
            sub_grid = gridspec.GridSpecFromSubplotSpec(1, 2, subplot_spec=self.grids[row, col],
                                                        width_ratios=[0.92, 0.08], wspace=0.05)

            ax = self.fig.add_subplot(sub_grid[0])
            cax = self.fig.add_subplot(sub_grid[1])

            ax.tick_params(
                axis='y',
                which='both',
                left='off',
                right='off',
                labelleft='on')

            if self.per_group:
                title = self.hm.matrix.group_labels[plot]
            else:
                title = self.hm.matrix.sample_labels[plot]

            ax.set_title(title)
            mat = []  # when drawing a heatmap (in contrast to drawing lines)
            for data_idx in range(self.numlines):
                if self.per_group:
                    row, col = plot, data_idx
                else:
                    row, col = data_idx, plot

                sub_matrix = self.hm.matrix.get_matrix(row, col)

                if self.per_group:
                    label = sub_matrix['sample']
                else:
                    label = sub_matrix['group']
                labels.append(label)

                mat.append(np.__getattribute__(self.averagetype)(sub_matrix['matrix'], axis=0))

            img = ax.imshow(np.vstack(mat), interpolation='nearest',
                            cmap='RdYlBu_r', aspect='auto', vmin=self.y_min, vmax=self.y_max)
            self.fig.colorbar(img, cax=cax)

            ax.axes.set_xticks(self.xticks)
            ax.axes.set_xticklabels(self.xtickslabel)
            # align the first and last label
            # such that they don't fall off
            # the heatmap sides
            ticks = ax.xaxis.get_major_ticks()
            ticks[0].label1.set_horizontalalignment('left')
            ticks[-1].label1.set_horizontalalignment('right')

            # add labels as y ticks labels
            ymin, ymax = ax.axes.get_ylim()
            pos, distance = np.linspace(ymin, ymax, len(labels), retstep=True, endpoint=False)
            d_half = float(distance) / 2
            yticks = [x + d_half for x in pos]

            ax.axes.set_yticks(yticks)
            ax.axes.set_yticklabels(labels[::-1], rotation='vertical')

            ax_list.append(ax)

        plt.subplots_adjust(wspace=0.05, hspace=0.3)
        plt.tight_layout()
        plt.savefig(self.out_file_name, dpi=200, format=self.image_format)
        plt.close()

def plot_single(ax, ma, average_type, color, label, plot_type='simple'):
    """
    Adds a line to the plot in the given ax using the specified method

    Parameters
    ----------
    ax : matplotlib axis
        matplotlib axis
    ma : numpy array
        numpy array The data on this matrix is summarized according
        to the `average_type` argument.
    average_type : str
        string values are sum mean median min max std
    color : str
        a valid color: either a html color name, hex
        (e.g #002233), RGB + alpha tuple or list or RGB tuple or list
    label : str
        label
    plot_type: str
        type of plot. Either 'se' for standard error, 'std' for
        standard deviation, 'overlapped_lines' to plot each line of the matrix,
        fill to plot the area between the x axis and the value or None, just to
        plot the average line.

    Returns
    -------
    ax
        matplotlib axis

    Examples
    --------

    >>> import matplotlib.pyplot as plt
    >>> fig = plt.figure()
    >>> ax = fig.add_subplot(111)
    >>> matrix = np.array([[1,2,3],
    ...                    [4,5,6],
    ...                    [7,8,9]])
    >>> ax = plot_single(ax, matrix -2, 'mean', color=[0.6, 0.8, 0.9], label='fill light blue', plot_type='fill')
    >>> ax = plot_single(ax, matrix, 'mean', color='blue', label='red')
    >>> ax = plot_single(ax, matrix + 5, 'mean', color='red', label='red', plot_type='std')
    >>> ax =plot_single(ax, matrix + 10, 'mean', color='#cccccc', label='gray se', plot_type='se')
    >>> ax = plot_single(ax, matrix + 20, 'mean', color=(0.9, 0.5, 0.9), label='violet', plot_type='std')
    >>> ax = plot_single(ax, matrix + 30, 'mean', color=(0.9, 0.5, 0.9, 0.5), label='violet with alpha', plot_type='std')
    >>> leg = ax.legend()
    >>> plt.savefig("/tmp/test.pdf")
    >>> fig = plt.figure()


    """
    summary = np.__getattribute__(average_type)(ma, axis=0)
    # only plot the average profiles without error regions
    x = np.arange(len(summary))
    ax.plot(x, summary, color=color, label=label, alpha=0.9)
    if plot_type == 'fill':
        pass
        ax.fill_between(x, summary, facecolor=color, alpha=0.6, edgecolor='none')

    if plot_type in ['se', 'std']:
        if plot_type == 'se':  # standard error
            std = np.std(ma, axis=0) / np.sqrt(ma.shape[0])
        else:
            std = np.std(ma, axis=0)

        alpha = 0.2
        # an alpha channel has to be added to the color to fill the area
        # between the mean (or median etc.) and the std or se
        f_color = pltcolors.colorConverter.to_rgba(color, alpha)

        ax.fill_between(x, summary, summary + std, facecolor=f_color, edgecolor='none')
        ax.fill_between(x, summary, summary - std, facecolor=f_color, edgecolor='none')

    ax.set_xlim(0, max(x))

    return ax