Attempts to intelligently load the .mat file and take average of left and right edges

        :return: left and right averages
        :return: times for each column
        :return: accept/reject for each column
        :return: pixel-inch ratio
354548
def _loadedges(self) -> typing.Tuple[np.ndarray, np.ndarray, np.ndarray, float, np.ndarray]:
        
        data = sco.loadmat(self.filename)
        datakeys = [k for k in data.keys()
                    if ('right' in k) or ('left' in k) or ('edge' in k)]
        averagedata = ((data[datakeys[0]] + data[datakeys[1]]) / 2)

        try:
            times = (data['times'] - data['times'].min())[0]
        except KeyError:
            times = np.arange(len(data[datakeys[0]][0]))

        try:
            accept = data['accept']
        except KeyError:
            accept = np.zeros(len(times))

        try:
            ratio = data['ratio']
        except KeyError:
            ratio = 1

        try:
            viscosity = data['viscosity']
        except KeyError:
            viscosity = np.ones(len(times))
        return averagedata, times, accept, ratio, viscosity
252706
def align_unwrapped(sino):
    
    samples = []
    if len(sino.shape) == 2:
        # 2D
        # take 1D samples at beginning and end of array
        samples.append(sino[:, 0])
        samples.append(sino[:, 1])
        samples.append(sino[:, 2])
        samples.append(sino[:, -1])
        samples.append(sino[:, -2])

    elif len(sino.shape) == 3:
        # 3D
        # take 1D samples at beginning and end of array
        samples.append(sino[:, 0, 0])
        samples.append(sino[:, 0, -1])
        samples.append(sino[:, -1, 0])
        samples.append(sino[:, -1, -1])
        samples.append(sino[:, 0, 1])

    # find discontinuities in the samples
    steps = np.zeros((len(samples), samples[0].shape[0]))
    for i in range(len(samples)):
        t = np.unwrap(samples[i])
        steps[i] = samples[i] - t

    # if the majority believes so, add a step of PI
    remove = mode(steps, axis=0)[0][0]

    # obtain divmod min
    twopi = 2*np.pi
    minimum = divmod_neg(np.min(sino), twopi)[0]
    remove += minimum*twopi

    for i in range(len(sino)):
        sino[i] -= remove[i]
1229940
def lengths( self ):
        
        return( np.array( [ math.sqrt( sum( row**2 ) ) for row in self.matrix ] ) )
1363787
def _setup(self):
        
        dtype = [('inv', 'f8'), ('min', 'f8'), ('max', 'f8'),('N', 'i4'), ('spacing','object')]
        dtype        = numpy.dtype(dtype)
        self._info   = numpy.empty(self.Ndim, dtype=dtype)
        self.min     = self._info['min']
        self.max     = self._info['max']
        self.N       = self._info['N']
        self.inv     = self._info['inv']
        self.spacing = self._info['spacing']

        for i, dim in enumerate(self.dims):

            self.N[i] = len(self.edges[i])-1
            self.min[i] = self.edges[i][0]
            self.max[i] = self.edges[i][-1]

            # determine the type of spacing
            self.spacing[i] = None
            lin_diff = numpy.diff(self.edges[i])
            with numpy.errstate(divide='ignore', invalid='ignore'):
                log_diff = numpy.diff(numpy.log10(self.edges[i]))
            if numpy.allclose(lin_diff, lin_diff[0]):
                self.spacing[i] = 'linspace'
                self.inv[i] = self.N[i] * 1.0 / (self.max[i] - self.min[i])
            elif numpy.allclose(log_diff, log_diff[0]):
                self.spacing[i] = 'logspace'
                self.inv[i] = self.N[i] * 1.0 / numpy.log10(self.max[i] / self.min[i])

        self.shape = self.N + 2

        # store Rmax
        self.Rmax = self.max[0]
1110039
def get_axis_grid(self, ind):
        
        ng = self.dim
        num_pts = ng[ind]
        lengths = self.structure.lattice.abc
        return [i / num_pts * lengths[ind] for i in range(num_pts)]
306839
def avgMultiplePlots(data, calc_mean=True, calc_std=False,
                     calc_density=False, nsample=None):
    
    xArr, yArr = bringPlotOverSameX(data, nsample)
    out = [xArr]
    if calc_mean:
        out.append(np.nanmean(yArr, axis=0))
    if calc_std:
        out.append(np.nanstd(yArr, axis=0))
    if calc_density:
        out.append(len(yArr) - np.isnan(yArr).sum(axis=0))
    return tuple(out)
701245
def _calculate_edges(self):
        
        # Calculate the left edge of the label.
        left = self.specs.left_margin
        left += (self.specs.label_width * (self._position[1] - 1))
        if self.specs.column_gap:
            left += (self.specs.column_gap * (self._position[1] - 1))
        left *= mm

        # And the bottom.
        bottom = self.specs.sheet_height - self.specs.top_margin
        bottom -= (self.specs.label_height * self._position[0])
        if self.specs.row_gap:
            bottom -= (self.specs.row_gap * (self._position[0] - 1))
        bottom *= mm

        # Done.
        return float(left), float(bottom)
1347687
def normalize_dense(M, norm="frag", order=1, iterations=3):
    

    s = np.array(M, np.float64)
    floatorder = np.float64(order)

    if norm == "SCN":
        for _ in range(0, iterations):

            sumrows = s.sum(axis=1)
            maskrows = (sumrows != 0)[:, None] * (sumrows != 0)[None, :]
            sums_row = sumrows[:, None] * np.ones(sumrows.shape)[None, :]
            s[maskrows] = 1. * s[maskrows] / sums_row[maskrows]

            sumcols = s.sum(axis=0)
            maskcols = (sumcols != 0)[:, None] * (sumcols != 0)[None, :]
            sums_col = sumcols[None, :] * np.ones(sumcols.shape)[:, None]
            s[maskcols] = 1. * s[maskcols] / sums_col[maskcols]

    elif norm == "mirnylib":
        try:
            from mirnylib import numutils as ntls
            s = ntls.iterativeCorrection(s, iterations)[0]
        except ImportError as e:
            print(str(e))
            print("I can't find mirnylib.")
            print("Please install it from "
                  "https://bitbucket.org/mirnylab/mirnylib")
            print("I will use default norm as fallback.")
            return normalize_dense(M, order=order, iterations=iterations)

    elif norm == "frag":
        for _ in range(1, iterations):
            s_norm_x = np.linalg.norm(s, ord=floatorder, axis=0)
            s_norm_y = np.linalg.norm(s, ord=floatorder, axis=1)
            s_norm = np.tensordot(s_norm_x, s_norm_y, axes=0)
            s[s_norm != 0] = 1. * s[s_norm != 0] / s_norm[s_norm != 0]

    elif norm == "global":
        s_norm = np.linalg.norm(s, ord=floatorder)
        s /= 1. * s_norm

    elif callable(norm):
        s = norm(M)

    else:
        print("Unknown norm. Returning input as fallback")

    return (s + s.T) / 2
229958
def loadSignalFromDelimitedFile(filename, timeColumnIdx=0, dataColumnIdx=1, delimiter=',', skipHeader=0) -> Signal:
    
    data = np.genfromtxt(filename, delimiter=delimiter, skip_header=skipHeader)
    columnCount = data.shape[1]
    if columnCount < timeColumnIdx + 1:
        raise ValueError(
            filename + ' has only ' + columnCount + ' columns,  time values can\'t be at column ' + timeColumnIdx)
    if columnCount < dataColumnIdx + 1:
        raise ValueError(
            filename + ' has only ' + columnCount + ' columns,  data values can\'t be at column ' + dataColumnIdx)

    t = data[:, [timeColumnIdx]]
    samples = data[:, [dataColumnIdx]]
    # calculate fs as the interval between the time samples
    fs = int(round(1 / (np.diff(t, n=1, axis=0).mean()), 0))
    source = Signal(samples.ravel(), fs)
    return source
602042
function calculateColumnAverage(inputMatrix)
{
	var rows = inputMatrix.rows();
	
	var averages = new Array();
	for(var i = 1; i <= inputMatrix.cols(); i++)
	{
		var sum = 0;
		for(var j = 1; j <= inputMatrix.rows(); j++)
		{
			sum += inputMatrix.e(j, i);
		}
		averages[i-1] = sum/rows;
	}
	
	return $V(averages);
}
336640
def bg_compensate(img, sigma, splinepoints, scale):
    ''''''
    
    from PIL import Image
    import pylab
    from matplotlib.image import pil_to_array
    from centrosome.filter import canny
    import matplotlib

    img = Image.open(img)
    if img.mode=='I;16':
        # 16-bit image
        # deal with the endianness explicitly... I'm not sure
        # why PIL doesn't get this right.
        imgdata = np.fromstring(img.tostring(),np.uint8)
        imgdata.shape=(int(imgdata.shape[0]/2),2)
        imgdata = imgdata.astype(np.uint16)
        hi,lo = (0,1) if img.tag.prefix == 'MM' else (1,0)
        imgdata = imgdata[:,hi]*256 + imgdata[:,lo]
        img_size = list(img.size)
        img_size.reverse()
        new_img = imgdata.reshape(img_size)
        # The magic # for maximum sample value is 281
        if 281 in img.tag:
            img = new_img.astype(np.float32) / img.tag[281][0]
        elif np.max(new_img) < 4096:
            img = new_img.astype(np.float32) / 4095.
        else:
            img = new_img.astype(np.float32) / 65535.
    else:
        img = pil_to_array(img)
    
    pylab.subplot(1,3,1).imshow(img, cmap=matplotlib.cm.Greys_r)
    pylab.show()
    
    if len(img.shape)>2:
        raise ValueError('Image must be grayscale')

## Create mask that will fix problem when image has black areas outside of well
    edges = canny(img, np.ones(img.shape, bool), 2, .1, .3)
    ci = np.cumsum(edges, 0)
    cj = np.cumsum(edges, 1)
    i,j = np.mgrid[0:img.shape[0], 0:img.shape[1]]
    mask = ci > 0
    mask = mask & (cj > 0)
    mask[1:,:] &= (ci[0:-1,:] < ci[-1,j[0:-1,:]])
    mask[:,1:] &= (cj[:,0:-1] < cj[i[:,0:-1],-1])
    
    import time
    t0 = time.clock()
    bg = backgr(img, mask, MODE_AUTO, sigma, splinepoints=splinepoints, scale=scale)
    print("Executed in %f sec" % (time.clock() - t0))
    bg[~mask] = img[~mask]

    pylab.subplot(1,3,2).imshow(img - bg, cmap=matplotlib.cm.Greys_r)
    pylab.subplot(1,3,3).imshow(bg, cmap=matplotlib.cm.Greys_r)
    pylab.show()
26619
def run_norm(net, df=None, norm_type='zscore', axis='row', keep_orig=False):
  '''
  
  '''

  # df here is actually a dictionary of several dataframes, 'mat', 'mat_orig',
  # etc
  if df is None:
    df = net.dat_to_df()

  if norm_type == 'zscore':
    df = zscore_df(df, axis, keep_orig)

  if norm_type == 'qn':
    df = qn_df(df, axis, keep_orig)

  net.df_to_dat(df)
830454
def block_widths(self):
        
        if self._widths_cache is None:
            try:
                # The first column will have the correct lengths. We have an
                # invariant that requires that all blocks be the same width in a
                # column of blocks.
                self._widths_cache = np.array(
                    ray.get([obj.width().oid for obj in self._partitions_cache[0]])
                    if len(self._partitions_cache) > 0
                    else []
                )
            except RayTaskError as e:
                handle_ray_task_error(e)
        return self._widths_cache
142500
def build_disagg_matrix(bdata, bin_edges, sid, mon=Monitor):
    
    with mon('build_disagg_matrix'):
        mag_bins, dist_bins, lon_bins, lat_bins, eps_bins = bin_edges
        dim1, dim2, dim3, dim4, dim5 = shape = get_shape(bin_edges, sid)

        # find bin indexes of rupture attributes; bins are assumed closed
        # on the lower bound, and open on the upper bound, that is [ )
        # longitude values need an ad-hoc method to take into account
        # the 'international date line' issue
        # the 'minus 1' is needed because the digitize method returns the
        # index of the upper bound of the bin
        mags_idx = numpy.digitize(bdata.mags+pmf.PRECISION, mag_bins) - 1
        dists_idx = numpy.digitize(bdata.dists[:, sid], dist_bins) - 1
        lons_idx = _digitize_lons(bdata.lons[:, sid], lon_bins[sid])
        lats_idx = numpy.digitize(bdata.lats[:, sid], lat_bins[sid]) - 1

        # because of the way numpy.digitize works, values equal to the last bin
        # edge are associated to an index equal to len(bins) which is not a
        # valid index for the disaggregation matrix. Such values are assumed
        # to fall in the last bin
        mags_idx[mags_idx == dim1] = dim1 - 1
        dists_idx[dists_idx == dim2] = dim2 - 1
        lons_idx[lons_idx == dim3] = dim3 - 1
        lats_idx[lats_idx == dim4] = dim4 - 1

        out = {}
        cache = {}
        cache_hit = 0
        num_zeros = 0
        for k, allpnes in bdata.items():
            pnes = allpnes[:, sid, :]  # shape (U, N, E)
            cache_key = pnes.sum()
            if cache_key == pnes.size:  # all pnes are 1
                num_zeros += 1
                continue  # zero matrices are not transferred
            try:
                matrix = cache[cache_key]
                cache_hit += 1
            except KeyError:
                mat = numpy.ones(shape)
                for i_mag, i_dist, i_lon, i_lat, pne in zip(
                        mags_idx, dists_idx, lons_idx, lats_idx, pnes):
                    mat[i_mag, i_dist, i_lon, i_lat] *= pne
                matrix = 1. - mat
                cache[cache_key] = matrix
            out[k] = matrix

        # operations, hits, num_zeros
        if hasattr(mon, 'cache_info'):
            mon.cache_info += numpy.array([len(bdata), cache_hit, num_zeros])
        else:
            mon.cache_info = numpy.array([len(bdata), cache_hit, num_zeros])
    return out
427810
def mat_to_numpy_arr(self):
  '''  '''
  import numpy as np
  self.dat['mat'] = np.asarray(self.dat['mat'])
576696
def aggregate_periods(self, periods):
        
        try:
            fieldname = self.raster_field.name
        except TypeError:
            raise exceptions.FieldDoesNotExist('Raster field not found')
        arrays = self.arrays(fieldname)
        arr = arrays[0]
        if len(arrays) > 1:
            if getattr(arr, 'ndim', 0) > 2:
                arrays = np.vstack(arrays)
            fill = getattr(arr, 'fill_value', None)
            arr = np.ma.masked_values(arrays, fill, copy=False)
        # Try to reshape using equal sizes first and fall back to unequal
        # splits.
        try:
            means = arr.reshape((periods, -1)).mean(axis=1)
        except ValueError:
            means = np.array([a.mean() for a in np.array_split(arr, periods)])
        obj = self[0]
        setattr(obj, fieldname, means)
        return [obj]
294565
def backgroundMMD(self, catalog, method='cloud-in-cells', weights=None):
        

        # Select objects in annulus
        cut_annulus = self.roi.inAnnulus(catalog.lon,catalog.lat)
        mag_1 = catalog.mag_1[cut_annulus]
        mag_2 = catalog.mag_2[cut_annulus]

        # Units are (deg^2)
        solid_angle = ugali.utils.binning.take2D(self.solid_angle_mmd, mag_2, mag_1,
                                                 self.roi.bins_mag, self.roi.bins_mag)

        # Weight each object before binning
        # Divide by solid angle and bin size in magnitudes to get number density 
        # [objs / deg^2 / mag^2]
        if weights is None:
            number_density = (solid_angle*self.roi.delta_mag**2)**(-1)
        else:
            number_density = weights*(solid_angle*self.roi.delta_mag**2)**(-1)

        method = str(method).lower()
        if method == 'cloud-in-cells':
            # Apply cloud-in-cells algorithm
            mmd_background = ugali.utils.binning.cloudInCells(mag_2,mag_1,
                                                              [self.roi.bins_mag,self.roi.bins_mag],
                                                              weights=number_density)[0]
        elif method == 'bootstrap':
            # Not implemented
            raise ValueError("Bootstrap method not implemented")
            mag_1 + (mag_1_err * np.random.normal(0, 1., len(mag_1)))
            mag_2 + (mag_2_err * np.random.normal(0, 1., len(mag_2)))

        elif method == 'histogram':
            # Apply raw histogram
            mmd_background = np.histogram2d(mag_1,mag_2,bins=[self.roi.bins_mag,self.roi.bins_mag],
                                            weights=number_density)[0]

        elif method == 'kde':
            # Gridded kernel density estimator
            logger.warning("### KDE not implemented properly")
            mmd_background = ugali.utils.binning.kernelDensity(mag_2,mag_1,
                                                               [self.roi.bins_mag,self.roi.bins_mag],
                                                               weights=number_density)[0]
        elif method == 'uniform':
            logger.warning("### WARNING: Uniform CMD")
            hist = np.histogram2d(mag_1,mag_2,bins=[self.roi.bins_mag,self.roi.bins_mag], 
                                  weights=number_density)[0]
            mmd_background = np.mean(hist)*np.ones(hist.shape)
            observable = (self.solid_angle_mmd > self.minimum_solid_angle)
            mmd_background *= observable
            return mmd_background
        else:
            raise ValueError("Unrecognized method: %s"%method)
        ## Account for the objects that spill out of the observable space
        ## But what about the objects that spill out to red colors??
        #for index_color in range(0, len(self.roi.centers_color)):
        #    for index_mag in range(0, len(self.roi.centers_mag)):
        #        if self.solid_angle_cmd[index_mag][index_color] < self.minimum_solid_angle:
        #            cmd_background[index_mag - 1][index_color] += cmd_background[index_mag][index_color]
        #            cmd_background[index_mag][index_color] = 0.
        #            break

        mmd_area = self.solid_angle_mmd*self.roi.delta_mag**2 # [deg^2 * mag^2]

        # ADW: This accounts for leakage to faint magnitudes
        # But what about the objects that spill out to red colors??
        # Maximum obsevable magnitude index for each color (uses the fact that
        # np.argmin returns first minimum (zero) instance found.
        # NOTE: More complicated maps may have holes causing problems

        observable = (self.solid_angle_mmd > self.minimum_solid_angle)
        index_mag_1 = observable.argmin(axis=0) - 1
        index_mag_2 = np.arange(len(self.roi.centers_mag))
        # Add the cumulative leakage back into the last bin of the CMD
        leakage = (mmd_background * ~observable).sum(axis=0)
        ### mmd_background[[index_mag_1,index_mag_2]] += leakage
        # Zero out all non-observable bins
        ### mmd_background *= observable

        # Avoid dividing by zero by setting empty bins to the value of the 
        # minimum filled bin of the CMD. This choice is arbitrary and 
        # could be replaced by a static minimum, some fraction of the 
        # CMD maximum, some median clipped minimum, etc. However, should 
        # be robust against outliers with very small values.
        min_mmd_background = max(mmd_background[mmd_background > 0.].min(),
                                 1e-4*mmd_background.max())
        mmd_background[observable] = mmd_background[observable].clip(min_mmd_background)

        ### # ADW: This is a fudge factor introduced to renormalize the CMD
        ### # to the number of input stars in the annulus. While leakage
        ### # will still smooth the distribution, it shouldn't result in 
        ### fudge_factor = len(mag) / float((cmd_background*cmd_area).sum())
        ### cmd_background *= fudge_factor

        return mmd_background