numpy.logical_and

Here are the examples of the python api numpy.logical_and taken from open source projects. By voting up you can indicate which examples are most useful and appropriate.

117 Examples 7

Example 1

Project: worldengine Source File: generation.py
def harmonize_ocean(ocean, elevation, ocean_level):
    """
    The goal of this function is to make the ocean floor less noisy.
    The underwater erosion should cause the ocean floor to be more uniform
    """

    shallow_sea = ocean_level * 0.85
    midpoint = shallow_sea / 2.0

    ocean_points = numpy.logical_and(elevation < shallow_sea, ocean)

    shallow_ocean = numpy.logical_and(elevation < midpoint, ocean_points)
    elevation[shallow_ocean] = midpoint - ((midpoint - elevation[shallow_ocean]) / 5.0)

    deep_ocean = numpy.logical_and(elevation > midpoint, ocean_points)
    elevation[deep_ocean] = midpoint + ((elevation[deep_ocean] - midpoint) / 5.0)

Example 2

Project: pyphi Source File: utils.py
def causally_significant_nodes(cm):
    """Returns a tuple of all nodes indices in the connectivity matrix which
    are causally significant (have inputs and outputs)."""
    inputs = cm.sum(0)
    outputs = cm.sum(1)
    nodes_with_inputs_and_outputs = np.logical_and(inputs > 0, outputs > 0)
    return tuple(np.where(nodes_with_inputs_and_outputs)[0])

Example 3

Project: TheCannon Source File: plot_apogee_lamost_cannon.py
def plot_row(name, SNR, ax, x_all, y_all, xlabel, ylabel, cut1, cut2, lim):
    #axarr[0].set_ylabel("LAMOST %s" %name)
    #ax[0].set_ylabel(ylabel)
    #ax[1].set_xlabel(xlabel)

    cond = SNR < cut1
    x = x_all[cond]
    y = y_all[cond]
    plot_one(title, ax[0], x, y, lim)

    cond = np.logical_and(SNR>cut1,SNR<cut2)
    x = x_all[cond]
    y = y_all[cond]
    plot_one(title, ax[1], x, y, lim)

    cond = SNR > cut2
    x = x_all[cond]
    y = y_all[cond]
    plot_one(title, ax[2], x, y, lim)

Example 4

Project: nupic.research Source File: sp_metrics.py
def calculateStability(activeColumnsCurrentEpoch, activeColumnsPreviousEpoch):
  activeColumnsStable = np.logical_and(activeColumnsCurrentEpoch,
                                       activeColumnsPreviousEpoch)
  stability = np.mean(np.sum(activeColumnsStable, 1))/\
              np.mean(np.sum(activeColumnsCurrentEpoch, 1))
  return stability

Example 5

Project: generativebot Source File: random.py
Function: random_points_in_rectangle
def random_points_in_rectangle(n, xx, yy, w, h):
  rnd = random((n,2))
  height_mask = logical_and(rnd[:,1]>yy-h*0.5, rnd[:,1]<yy+h*0.5)
  width_mask = logical_and(rnd[:,0]>xx-w*0.5, rnd[:,0]<xx+w*0.5)
  mask = logical_and(height_mask, width_mask)
  return rnd[mask,:]

Example 6

Project: sci-wms Source File: data_handler.py
def ugrid_lat_lon_subset_idx(lon, lat, bbox, padding=None):
    """
    Assumes the size of lat/lon are equal (UGRID variables).
    Returns a boolean mask array of the indexes, suitable for slicing
    """

    padding = padding or 0.18

    minlon = bbox[0] - padding
    minlat = bbox[1] - padding
    maxlon = bbox[2] + padding
    maxlat = bbox[3] + padding

    land = np.logical_and
    return land(land(lon >= minlon, lon <= maxlon),
                land(lat >= minlat, lat <= maxlat))

Example 7

Project: mcedit2 Source File: __init__.py
Function: box_mask
    def box_mask(self, box):
        masks = [s.box_mask(box) for s in self.selections]
        if any(m is None for m in masks):
            return None

        m = masks.pop()

        while len(masks):
            numpy.logical_and(m, masks.pop(), m)

        return m

Example 8

Project: kaggle-seeclickfix-ensemble Source File: munge.py
Function: segment_data
def segment_data(dfTrn,dfTest,segment):
    if segment == 'remote_api_created':
        dfTrn = dfTrn[dfTrn.source == segment]
        dfTest = dfTest[dfTest.source == segment]
    if segment in ['Oakland']:
        #For Oakland, clean remote_api issues off of the test set.  Leave them in training set
        #due to marginally better CV performance with remote_api issues left in training set
        dfTrn = dfTrn[dfTrn.city == segment]
        dfTest = dfTest[dfTest.city == segment]
        dfTest = dfTest[np.logical_and(dfTest.city == segment,dfTest.source != 'remote_api_created')]
    if segment in ['Chicago']:
        dfTrn = dfTrn[np.logical_and(dfTrn.city == segment,dfTrn.source != 'remote_api_created')]
        dfTest = dfTest[np.logical_and(dfTest.city == segment,dfTest.source != 'remote_api_created')]
    if segment in ['New Haven','Richmond']:
        dfTrn = dfTrn[dfTrn.city == segment]
        dfTest = dfTest[dfTest.city == segment]
    return dfTrn, dfTest

Example 9

Project: qiime Source File: shared_phylotypes.py
def _calc_shared_phylotypes_pairwise(otu_table, i, j):
    """Calculate shared otus between two samples in column i and j.

    otu_table: OTU tables as a OTUtable subclass

    i: a sample id in the OTU table

    j: a sample id in the OTU table
    """
    shared_phylos = logical_and(
        otu_table.data(i, 'sample'),
        otu_table.data(j, 'sample'))

    return shared_phylos.sum()

Example 10

Project: paegan Source File: n_cell.py
Function: get_xyind_from_bbox
    def get_xyind_from_bbox(self, var, bbox):
        grid = self.getgridobj(var)
        xbool = grid.get_xbool_from_bbox(bbox)
        ybool = grid.get_ybool_from_bbox(bbox)
        inds = np.where(np.logical_and(xbool, ybool))
        return inds, inds #xinds, yinds

Example 11

Project: ptsa Source File: cluster.py
Function: get_components
def _get_components(x_in, connectivity):
    """get connected components from a mask and a connectivity matrix"""
    cs_graph_components = sparse.cs_graph_components

    mask = np.logical_and(x_in[connectivity.row], x_in[connectivity.col])
    data = connectivity.data[mask]
    row = connectivity.row[mask]
    col = connectivity.col[mask]
    shape = connectivity.shape
    idx = np.where(x_in)[0]
    row = np.concatenate((row, idx))
    col = np.concatenate((col, idx))
    data = np.concatenate((data, np.ones(len(idx), dtype=data.dtype)))
    connectivity = sparse.coo_matrix((data, (row, col)), shape=shape)
    _, components = cs_graph_components(connectivity)
    # print "-- number of components : %d" % np.unique(components).size
    return components

Example 12

Project: differential-lattice Source File: differentialLattice.py
Function: decay
  def decay(self, age, ratio=0.01):

    num = self.num
    die = logical_and(
        self.age[:num,0]<self.itt-age,
        self.link_counts[:num,0]<1
        )
    alive = logical_not(logical_and(die, random(len(die))<ratio))

    inds = alive.nonzero()[0]

    new_num = len(inds)

    self.xy[:new_num,:] = self.xy[inds,:]
    self.num = new_num

Example 13

Project: SimpleCV Source File: demo_cv_threshold.py
def show_depth():
    global threshold
    global current_depth

    depth, timestamp = freenect.sync_get_depth()
    depth = 255 * np.logical_and(depth >= current_depth - threshold,
                                 depth <= current_depth + threshold)
    depth = depth.astype(np.uint8)
    image = cv.CreateImageHeader((depth.shape[1], depth.shape[0]),
                                 cv.IPL_DEPTH_8U,
                                 1)
    cv.SetData(image, depth.tostring(),
               depth.dtype.itemsize * depth.shape[1])
    cv.ShowImage('Depth', image)

Example 14

Project: ramp Source File: selectors.py
Function: sets
    def sets(self, x, y, n_keep):
        cnts = y.value_counts()
        assert(len(cnts) == 2)
        logging.info("Computing IG scores...")
        scores = []
        for c in x.columns:
            true_positives = sum(np.logical_and(x[c], y))
            false_positives = sum(np.logical_and(x[c], np.logical_not(y)))
            score = abs(norm.ppf(true_positives / float(cnts[1])) - norm.ppf(false_positives / float(cnts[0])))
            scores.append((score, c))
        scores.sort(reverse=True)
        return [s[1] for s in scores[:n_keep]]

Example 15

Project: landlab Source File: fault_facet_finder.py
    def show_possible_nodes(self):
        """
        Once the subsets by aspect and slope have been set, call this function
        to see both the whole elevation map, and the subset of nodes that
        will be searched.
        """
        possible_core_nodes = np.logical_and(
            self.steep_nodes, self.aspect_close_nodes)
        figure(1)
        gridshow.imshow_grid_at_node(self.grid, self.elevs)
        figure(2)
        gridshow.imshow_grid_at_node(self.grid, self.slopes)
        figure(3)
        gridshow.imshow_grid_at_node(self.grid, self.aspect)
        figure(4)
        gridshow.imshow_grid_at_node(self.grid, possible_core_nodes)
        show()

Example 16

Project: TheCannon Source File: mass_age_comparison.py
def load_comparison():
    direc = "/Users/annaho/Data/LAMOST/Mass_And_Age"
    hdulist = pyfits.open("%s/age_vs_age_catalog.fits" %direc)
    tbdata = hdulist[1].data
    hdulist.close()
    snr = tbdata.field("snr")
    mass_raw = tbdata.field("cannon_mass")
    choose = np.logical_and(snr > 30, mass_raw > 0)
    mass = np.log10(tbdata.field("cannon_mass")[choose])
    age = tbdata.field("cannon_age")[choose]
    age_err = tbdata.field("cannon_age_err")[choose]
    mass_err = tbdata.field("cannon_mass_err")[choose]
    ness_age = np.log10(np.exp(tbdata.field("lnAge")[choose]))
    ness_age_err = tbdata.field("e_logAge")[choose]
    ness_mass = np.log10(np.exp(tbdata.field("lnM")[choose]))
    ness_mass_err = tbdata.field("e_logM")[choose]
    return mass, age, ness_mass, ness_age

Example 17

Project: statsmodels Source File: utils.py
def _get_dataset_meta(dataname, package, cache):
    # get the index, you'll probably want this cached because you have
    # to download info about all the data to get info about any of the data...
    index_url = ("https://raw.github.com/vincentarelbundock/Rdatasets/master/"
                 "datasets.csv")
    data, _ = _urlopen_cached(index_url, cache)
    # Python 3
    if PY3:  # pragma: no cover
        data = data.decode('utf-8', 'strict')
    index = read_csv(StringIO(data))
    idx = np.logical_and(index.Item == dataname, index.Package == package)
    dataset_meta = index.ix[idx]
    return dataset_meta["Title"].item()

Example 18

Project: mne-python Source File: urls.py
Function: url_match
def url_match(condition, data_format, data_type):
    """Function to match MEGSIM data files."""
    inds = np.logical_and(conditions == condition, data_formats == data_format)
    inds = np.logical_and(inds, data_types == data_type)
    inds = np.logical_and(inds, data_formats == data_format)
    good_urls = list(urls[inds])
    for gi, g in enumerate(good_urls):
        good_urls[gi] = url_root + g
    if len(good_urls) == 0:
        raise ValueError('No MEGSIM dataset found with condition="%s",\n'
                         'data_format="%s", data_type="%s"'
                         % (condition, data_format, data_type))
    return good_urls

Example 19

Project: simpeg Source File: Optimization.py
    @Utils.count
    def bindingSet(self, x):
        """bindingSet(x)

            If we are on a bound and the negative gradient points away from the feasible set.

            Optimality condition. (Satisfies Kuhn-Tucker) MoreToraldo91

        """
        bind_up  = np.logical_and(x == self.lower, self.g >= 0)
        bind_low = np.logical_and(x == self.upper, self.g <= 0)
        return np.logical_or(bind_up, bind_low)

Example 20

Project: paegan Source File: c_grid.py
Function: get_xyind_from_bbox
    def get_xyind_from_bbox(self, var, bbox, **kwargs):
        grid = self.getgridobj(var)
        xbool = grid.get_xbool_from_bbox(bbox)
        ybool = grid.get_ybool_from_bbox(bbox)
        inds = np.where(np.logical_and(xbool, ybool))
        minrow = np.min(inds[0])
        mincol = np.min(inds[1])
        maxrow = np.max(inds[0])
        maxcol = np.max(inds[1])
        inds = (np.arange(minrow, maxrow+1), np.arange(mincol, maxcol+1))
        return inds, inds #xinds, yinds

Example 21

Project: qiime Source File: shared_phylotypes.py
def _calc_shared_phylotypes_multiple(otu_table, idxs):
    """Calculate shared otus between several samples indexed by values in idxes.

    otu_table: OTU table as a OTUtable subclass
    idxs: list of sample ids in the OTU table
    """

    if len(idxs) < 2:
        raise ValueError("calc_shared_phylotypes_multiple needs at least two "
                         "sampleIDs to comapre")

    shared_phylos = ones(len(otu_table.ids(axis='observation')))

    for id_ in idxs:
        shared_phylos = logical_and(shared_phylos, otu_table.data(id_, 'sample'))

    return shared_phylos.sum()

Example 22

Project: ldsc Source File: regressions.py
Function: init
    def __init__(self, z1, z2, x, w, N1, N2, M, hsq1, hsq2, intercept_hsq1, intercept_hsq2,
                 n_blocks=200, intercept_gencov=None, slow=False, twostep=None):
        self.intercept_hsq1 = intercept_hsq1
        self.intercept_hsq2 = intercept_hsq2
        self.hsq1 = hsq1
        self.hsq2 = hsq2
        self.N1 = N1
        self.N2 = N2
        y = z1 * z2
        step1_ii = None
        if twostep is not None:
            step1_ii = np.logical_and(z1**2 < twostep, z2**2 < twostep)

        LD_Score_Regression.__init__(self, y, x, w, np.sqrt(N1 * N2), M, n_blocks,
                                     intercept=intercept_gencov, slow=slow, step1_ii=step1_ii)
        self.p, self.z = p_z_norm(self.tot, self.tot_se)
        self.mean_z1z2 = np.mean(np.multiply(z1, z2))

Example 23

Project: paegan Source File: dataset.py
    def get_tind_from_bounds(self, var, bounds, convert=False, use_cache=True):
        assert var in self._current_variables
        time = self.gettimevar(var, use_cache)
        if convert:
            bounds = netCDF4.date2num(bounds, time._units + " since " + time.origin.isoformat())
        inds = np.where(np.logical_and(time.dates >= bounds[0].replace(tzinfo=time.dates[0].tzinfo), time.dates <= bounds[1].replace(tzinfo=time.dates[0].tzinfo)))
        return inds

Example 24

Project: rlpy Source File: PST.py
Function: is_terminal
    def isTerminal(self):
        sStruct = self.state2Struct(self.state)
        return (
            np.any(np.logical_and(sStruct.fuel <= 0,
                   sStruct.locations != UAVLocation.REFUEL))
        )

Example 25

Project: refinery Source File: TestFromSaved.py
Function: test_viable_init
  def test_viable_init(self):
    ''' Verify hmodel after init can be used to perform E-step
    '''
    initSavedParams = dict(initname='/tmp/', prefix='Test')
    self.hmodel.init_global_params(self.Data, **initSavedParams)
    assert self.hmodel.allocModel.K == self.K
    keysA = self.hmodel.allocModel.to_dict()
    keysB = self.modelB.allocModel.to_dict()
    assert len(keysA) == len(keysB)
    
    aLP = self.hmodel.calc_local_params(self.Data)
    assert np.all(np.logical_and(aLP['resp']>=0,aLP['resp']<=1.0))
    assert np.allclose(1.0, np.sum(aLP['resp'],axis=1))

Example 26

Project: refinery Source File: TestFromScratchGauss.py
Function: test_viable_init
  def test_viable_init(self, K=7):
    ''' Verify hmodel after init can be used to perform E-step
    '''
    for initname in ['randexamples', 'randsoftpartition', 'randexamplesbydist']:
      initParams = dict(initname=initname, seed=0, K=K)
      self.hmodel.init_global_params(self.Data, **initParams)
      LP = self.hmodel.calc_local_params(self.Data)
      resp = LP['resp']
      assert np.all(np.logical_and(resp >=0, resp <= 1.0))

Example 27

Project: hyperspy Source File: strategy.py
Function: refresh
    def refresh(self, overwrite, given_pixels=None):
        """Refreshes the database (i.e. constructs it again from scratch)
        """
        scale = self.samf._scale
        if overwrite:
            if given_pixels is None:
                good_pixels = self.samf.metadata.marker < 0
            else:
                good_pixels = given_pixels
            self.samf.metadata.marker[good_pixels] = -scale
        else:
            good_pixels = self.samf.metadata.marker < 0
            if given_pixels is not None:
                good_pixels = np.logical_and(good_pixels, given_pixels)
        self.samf.metadata.marker[~good_pixels] = scale
        self._update_database(None, 0)  # to force to update

Example 28

Project: GPy Source File: piecewise_linear.py
Function: gradients_x
    def gradients_X(self, dL_dF, X):
        x = X.flatten()

        #outside the range of the breakpoints, the function is just offset by a contant, so the partial derivative is 1.
        dL_dX = dL_dF.copy().flatten()

        #insude the breakpoints, the partial derivative is self.grads
        for low, up, g, v in zip(self. sorted_breaks[:-1], self.sorted_breaks[1:], self.grads, self.sorted_values[:-1]):
            i = np.logical_and(x>low, x<up)
            dL_dX[i] = dL_dF[i]*g

        return dL_dX.reshape(-1,1)

Example 29

Project: generativebot Source File: random.py
Function: random_points_in_rectangle
def random_points_in_rectangle(n, xx, yy, w, h):

  rnd = random((n,2))
  height_mask = logical_and(rnd[:,1]>yy-h*0.5, rnd[:,1]<yy+h*0.5)
  width_mask = logical_and(rnd[:,0]>xx-w*0.5, rnd[:,0]<xx+w*0.5)
  mask = logical_and(height_mask, width_mask)
  return rnd[mask,:]

Example 30

Project: scikit-tensor Source File: sptensor.py
Function: get_item
    def __getitem__(self, idx):
        # TODO check performance
        if len(idx) != self.ndim:
            raise ValueError('subscripts must be complete')
        sidx = ones(len(self.vals))
        for i in range(self.ndim):
            sidx = np.logical_and(self.subs[i] == idx[i], sidx)
        vals = self.vals[sidx]
        if len(vals) == 0:
            vals = 0
        elif len(vals) > 1:
            if self.accuemfun is None:
                raise ValueError('Duplicate entries without specified accuemulation function')
            vals = self.accuemfun(vals)
        return vals

Example 31

Project: research2epub Source File: findspaces.py
def rect_around_region(top, left, right, bottom):
	a = multi_and(
		indices[0] <= right, indices[0] >= left,
		indices[1] <= bottom, indices[1] >= top)
	b = multi_and(
		indices[0] < right, indices[0] > left,
		indices[1] < bottom, indices[1] > top)
	c =  numpy.logical_and(a, -b)
	#print a.sum(), b.sum()
	#print c.sum(), (top - bottom) * 2 + (right - left) * 2 - 2
	return c

Example 32

Project: SimpleCV Source File: demo_cv_thresh_sweep.py
def disp_thresh(lower, upper):
    depth, timestamp = freenect.sync_get_depth()
    depth = 255 * np.logical_and(depth > lower, depth < upper)
    depth = depth.astype(np.uint8)
    image = cv.CreateImageHeader((depth.shape[1], depth.shape[0]),
                                 cv.IPL_DEPTH_8U,
                                 1)
    cv.SetData(image, depth.tostring(),
               depth.dtype.itemsize * depth.shape[1])
    cv.ShowImage('Depth', image)
    cv.WaitKey(10)

Example 33

Project: PyEMMA Source File: fragmented_trajectory_reader.py
    def __get_ra_index_indices(self):
        """
        Returns a list containing indices of the ra_index array, which correspond to the separate trajectory fragments,
        i.e., ra_indices[fragment_indices[itraj]] are the ra indices for itraj (plus some offset by
        cuemulative length)
        """
        fragment_indices = []
        for idx, cuemlen in enumerate(self._cuemulative_lengths):
            cuemlen_prev = self._cuemulative_lengths[idx - 1] if idx > 0 else 0
            fragment_indices.append([np.argwhere(
                np.logical_and(self.ra_indices >= cuemlen_prev, self.ra_indices < cuemlen)
            )])
        return fragment_indices

Example 34

Project: sklearn-pmml Source File: utils.py
    @staticmethod
    def field_not_in_list(field, values):
        mv = pmml.MapValues(outputColumn='output', defaultValue=1)
        mv.append(pmml.FieldColumnPair(field=field, column='input'))
        it = pmml.InlineTable()
        for v in values:
            it.append(pmml_row(input=v, output=0))
        mv.append(it)
        return {
            DerivedFeatureTransformations.TRANSFORMATION: mv,
            DerivedFeatureTransformations.FUNCTION: lambda df: reduce(np.logical_and, [df[field] != _ for _ in values])
        }

Example 35

Project: TheCannon Source File: feh_alpha_cage.py
Function: load_data
def load_data():
    direc = "/Users/annaho/Data/LAMOST/Mass_And_Age"
    hdulist = pyfits.open("%s/catalog_paper.fits" %direc)
    tbdata = hdulist[1].data
    hdulist.close()
    snr = tbdata.field("snr")
    chisq = tbdata.field("chisq")
    teff = tbdata.field("cannon_teff")
    in_martig_range = tbdata.field("in_martig_range")
    choose = np.logical_and(in_martig_range, snr > 80)
    print(sum(choose))
    mh = tbdata.field("cannon_mh")[choose]
    afe = tbdata.field("cannon_am")[choose]
    age = tbdata.field("cannon_age")[choose]
    age_err = tbdata.field("cannon_age_err")[choose]
    return mh, afe, age, age_err

Example 36

Project: python-oceans Source File: datasets.py
Function: get_indices
def _get_indices(bbox, lons, lats):
    """Return the data indices for a lon, lat square."""
    lons = wrap_lon180(lons)

    idx_x = np.logical_and(lons >= bbox[0], lons <= bbox[1])
    idx_y = np.logical_and(lats >= bbox[2], lats <= bbox[3])
    if lons.ndim == 2 and lats.ndim == 2:
        inregion = np.logical_and(idx_x, idx_y)
        region_inds = np.where(inregion)
        imin, imax = _minmax(region_inds[0])
        jmin, jmax = _minmax(region_inds[1])
    elif lons.ndim == 1 and lats.ndim == 1:
        imin, imax = _minmax(np.where(idx_x))
        jmin, jmax = _minmax(np.where(idx_y))
    else:
        msg = 'Cannot understand input shapes lons {!r} and lats {!r}'.format
        raise ValueError(msg(lons.shape, lats.shape))
    return imin, imax+1, jmin, jmax+1

Example 37

Project: scikit-image Source File: test_random_noise.py
def test_salt_and_pepper():
    seed = 42
    cam = img_as_float(camera())
    cam_noisy = random_noise(cam, seed=seed, mode='s&p', amount=0.15,
                             salt_vs_pepper=0.25)
    saltmask = np.logical_and(cam != cam_noisy, cam_noisy == 1.)
    peppermask = np.logical_and(cam != cam_noisy, cam_noisy == 0.)

    # Ensure all changes are to 0. or 1.
    assert_allclose(cam_noisy[saltmask], np.ones(saltmask.sum()))
    assert_allclose(cam_noisy[peppermask], np.zeros(peppermask.sum()))

    # Ensure approximately correct amount of noise was added
    proportion = float(
        saltmask.sum() + peppermask.sum()) / (cam.shape[0] * cam.shape[1])
    assert 0.11 < proportion <= 0.18

    # Verify the relative amount of salt vs. pepper is close to expected
    assert 0.18 < saltmask.sum() / float(peppermask.sum()) < 0.33

Example 38

Project: GPy Source File: piecewise_linear.py
Function: f
    def f(self, X):
        x = X.flatten()
        y = x.copy()

        #first adjus the points below the first value
        y[x<self.sorted_breaks[0]]  = x[x<self.sorted_breaks[0]] + self.sorted_values[0] - self.sorted_breaks[0]

        #now all the points pas the last break
        y[x>self.sorted_breaks[-1]]  = x[x>self.sorted_breaks[-1]] + self.sorted_values[-1] - self.sorted_breaks[-1]

        #loop throught the pairs of points
        for low, up, g, v in zip(self. sorted_breaks[:-1], self.sorted_breaks[1:], self.grads, self.sorted_values[:-1]):
            i = np.logical_and(x>low, x<up)
            y[i] = v + (x[i]-low)*g

        return y.reshape(-1,1)

Example 39

Project: pycortex Source File: volume.py
Function: sample
    @staticmethod
    def _sample(pts, shape, norm):
        coords = pts.round().astype(int)[:,::-1]
        d1 = np.logical_and(0 <= coords[:,0], coords[:,0] < shape[0])
        d2 = np.logical_and(0 <= coords[:,1], coords[:,1] < shape[1])
        d3 = np.logical_and(0 <= coords[:,2], coords[:,2] < shape[2])
        valid = np.logical_and(d1, np.logical_and(d2, d3))
        if valid.any():
            idx = np.ravel_multi_index(coords[valid].T, shape)
            j, data = np.array(Counter(idx).items()).T
            return j, data / float(norm)

Example 40

Project: PRST Source File: wells_and_bc.py
def boundaryFaceIndices(G, side, I0, I1, I2):
    """
    Retrieve face indices belonging to a subset of global outer faces.

    Synopsis:
        ix = boundaryFaceIndices(G, side, I0, I1, I2)

    Arguments:
        G (Grid):
            pyrst.gridprocessing.Grid

        side (str):
            Global side from which to extract face indices. String. Must (case
            insensitively) match one of six alias groups:

                0) ["west",  "xmin", "left"  ]
                1) ["east",  "xmax", "right" ]
                2) ["south", "ymin", "back"  ]
                3) ["north", "ymax", "front" ]
                4) ["upper", "zmin", "top"   ]
                5) ["lower", "zmax", "bottom"]

            These groups correspond to the cardinal directions mentiond as the
            first alternative in each group.

        I0, I1 (list or ndarray):
            Index ranges for local (in-plane) axes one and two, respectively.
            No index range given (I1 or I2 is None) is interpreted as covering
            the entire corresponding local axis of `side` in the grid `G`. The
            local axes on a side in G are ordered `X` before `Y` before `Z`.

        I2 (list or ndarray):
            Index range for global axis perpendicular to `side`. The primary
            purpose of this parameter is to exclude faces *within* a reservoir
            from being added to the return value `ix`. Such faces typically
            occur in faulted reservoirs where a given face may be considered
            external by virtue of being connected to a single reservoir cell
            only.

    Returns:
        ix - Required face indices as a column array.

    Note:
        This function is mainly intended for internal use in this file. Its
        calling interface may change more frequently than functions in the
        BoundaryCondition class.

    See also:
        prst.params.wells_and_bc.BoundaryCondition
    """
    try:
        I0 = I0.ravel()
    except AttributeError:
        pass
    try:
        I1 = I1.ravel()
    except AttributeError:
        pass
    try:
        I2 = I2.ravel()
    except AttributeError:
        pass

    ## Extract all faces of cells within the given subset.
    cells, faceTag, isOutF = _boundaryCellsSubset(G, side, I0, I1, I2)

    fIX = G.cells.facePos;
    hfIX = mcolon(fIX[cells], fIX[cells+1])
    faces = G.cells.faces[hfIX, 0]
    tags = G.cells.faces[hfIX, 1]
    ix = faces[np.logical_and(isOutF[faces], tags == faceTag)]
    # return as column array
    return ix[:,np.newaxis]

Example 41

Project: python-qinfer Source File: multicos.py
Function: are_models_valid
    def are_models_valid(self, modelparams):
        return np.all(np.logical_and(modelparams > 0, modelparams <= 1), axis=1)

Example 42

Project: PyParticles Source File: pseudo_bubble.py
Function: update_force
    def update_force( self , pset ):
        
        D = self.__D
        
        D[:] = dist.squareform( dist.pdist( pset.X ) )
        ( n , m ) = np.where( np.logical_and( D <= self.__R , D != 0.0 ) )
        
        #print(np.where(b))
        
        self.__F[:] = 0.0
        
        self.__F[n,m] = ( -( self.__B / self.__R ) * D[n,m] + self.__B ) / D[n,m]
        
        for i in range( pset.dim ) :
            self.__V[:,:] = pset.X[:,i]
            self.__V[:,:] = ( self.__V[:,:].T - pset.X[:,i] ).T
            
            self.__A[:,i] = np.sum( self.__F * self.__V[:,:] / self.__M.T , 0 )
        
        return self.__A

Example 43

Project: scikit-beam Source File: arithmetic.py
def logical_nand(x1, x2, out=None):
    """Computes the truth value of NOT (x1 AND x2) element wise.

    This function enables the computation of the LOGICAL_NAND of two image or
    volume data sets. This function enables easy isolation of all data points
    NOT INCLUDED IN BOTH SOURCE DATA SETS. This function can be used for data
    comparison, material isolation, noise removal, or mask
    application/generation.

    Parameters
    ----------
    x1, x2 : array-like
        Input arrays. `x1` and `x2` must be of the same shape.

    output : array-like
        Boolean result with the same shape as `x1` and `x2` of the logical
        operation on corresponding elements of `x1` and `x2`.

    Returns
    -------
    output : {ndarray, bool}
        Boolean result with the same shape as `x1` and `x2` of the logical
        NAND operation on corresponding elements of `x1` and `x2`.


    Example
    -------
    >>> x1 = [[0,0,1,0,0], [2,1,1,1,2], [2,0,1,0,2]]
    >>> x2 = [[0,0,0,0,0], [2,1,1,1,2], [0,0,0,0,0]]
    >>> logical_nand(x1, x2)
    array([[ True,  True,  True,  True,  True],
           [False, False, False, False, False],
           [ True,  True,  True,  True,  True]], dtype=bool)

    """
    return logical_not(logical_and(x1, x2, out), out)

Example 44

Project: python-qinfer Source File: models.py
Function: are_models_valid
    def are_models_valid(self, modelparams):
        return np.logical_and(
            super(DiffusiveTomographyModel, self).are_models_valid(modelparams),
            modelparams[:, -1] > 0
        )

Example 45

Project: scikit-beam Source File: arithmetic.py
def logical_sub(x1, x2, out=None):
    """Compute truth value of x1 AND (NOT (x1 AND x2)) element wise.

    This function enables LOGICAL SUBTRACTION of one binary image or volume
    data set from another. This function can be used to remove phase
    information, interface boundaries, or noise, present in two data sets,
    without having to worry about mislabeling of pixels which would result
    from arithmetic subtraction. This function will evaluate as true for all
    "true" voxels present ONLY in Source Dataset 1. This function can be used
    for data cleanup, or boundary/interface analysis.

    Parameters
    ----------
    x1, x2 : array-like
        Input arrays. `x1` and `x2` must be of the same shape.

    output : array-like
        Boolean result with the same shape as `x1` and `x2` of the logical
        operation on corresponding elements of `x1` and `x2`.

    Returns
    -------
    output : {ndarray, bool}
        Boolean result with the same shape as `x1` and `x2` of the logical
        SUBTRACT operation on corresponding elements of `x1` and `x2`.

    Example
    -------
    >>> x1 = [[0,0,1,0,0], [2,1,1,1,2], [2,0,1,0,2]]
    >>> x2 = [[0,0,0,0,0], [2,1,1,1,2], [0,0,0,0,0]]
    >>> logical_sub(x1, x2)
    array([[False, False,  True, False, False],
           [False, False, False, False, False],
           [ True, False,  True, False,  True]], dtype=bool)
    """
    return logical_and(x1, logical_not(logical_and(x1, x2, out), out), out)

Example 46

Project: rlpy Source File: PST.py
Function: showdomain
    def showDomain(self, a=0):
        s = self.state
        if self.domain_fig is None:
            self.domain_fig = plt.figure(
                1, (UAVLocation.SIZE * self.dist_between_locations + 1, self.NUM_UAV + 1))
            plt.show()
        plt.clf()
         # Draw the environment
         # Allocate horizontal 'lanes' for UAVs to traverse

        # Formerly, we checked if this was the first time plotting; wedge shapes cannot be removed from
        # matplotlib environment, nor can their properties be changed, without clearing the figure
        # Thus, we must redraw the figure on each timestep
        #        if self.location_rect_vis is None:
        # Figure with x width corresponding to number of location states, UAVLocation.SIZE
        # and rows (lanes) set aside in y for each UAV (NUM_UAV total lanes).
        # Add buffer of 1
        self.subplot_axes = self.domain_fig.add_axes(
            [0, 0, 1, 1], frameon=False, aspect=1.)
        crashLocationX = 2 * \
            (self.dist_between_locations) * (UAVLocation.SIZE - 1)
        self.subplot_axes.set_xlim(0, 1 + crashLocationX + self.RECT_GAP)
        self.subplot_axes.set_ylim(0, 1 + self.NUM_UAV)
        self.subplot_axes.xaxis.set_visible(False)
        self.subplot_axes.yaxis.set_visible(False)

        # Assign coordinates of each possible uav location on figure
        self.location_coord = [0.5 + (self.LOCATION_WIDTH / 2) +
                               (self.dist_between_locations) * i for i in range(UAVLocation.SIZE - 1)]
        self.location_coord.append(crashLocationX + self.LOCATION_WIDTH / 2)

         # Create rectangular patches at each of those locations
        self.location_rect_vis = [mpatches.Rectangle(
            [0.5 + (self.dist_between_locations) * i,
             0],
            self.LOCATION_WIDTH,
            self.NUM_UAV * 2,
            fc='w') for i in range(UAVLocation.SIZE - 1)]
        self.location_rect_vis.append(
            mpatches.Rectangle([crashLocationX,
                                0],
                               self.LOCATION_WIDTH,
                               self.NUM_UAV * 2,
                               fc='w'))
        [self.subplot_axes.add_patch(self.location_rect_vis[i])
         for i in range(4)]
        self.comms_line = [lines.Line2D(
            [0.5 + self.LOCATION_WIDTH + (self.dist_between_locations) * i,
             0.5 + self.LOCATION_WIDTH + (
                 self.dist_between_locations) * i + self.RECT_GAP],
            [self.NUM_UAV * 0.5 + 0.5,
             self.NUM_UAV * 0.5 + 0.5],
            linewidth=3,
            color='black',
            visible=False) for i in range(UAVLocation.SIZE - 2)]
        self.comms_line.append(
            lines.Line2D(
                [0.5 + self.LOCATION_WIDTH + (self.dist_between_locations) * 2,
                 crashLocationX],
                [self.NUM_UAV * 0.5 + 0.5,
                 self.NUM_UAV * 0.5 + 0.5],
                linewidth=3,
                color='black',
                visible=False))

        # Create location text below rectangles
        locText = ["Base", "Refuel", "Communication", "Surveillance"]
        self.location_rect_txt = [plt.text(
            0.5 + self.dist_between_locations * i + 0.5 * self.LOCATION_WIDTH,
            -0.3,
            locText[i],
            ha='center') for i in range(UAVLocation.SIZE - 1)]
        self.location_rect_txt.append(
            plt.text(crashLocationX + 0.5 * self.LOCATION_WIDTH, -0.3,
                     locText[UAVLocation.SIZE - 1], ha='center'))

        # Initialize list of circle objects

        uav_x = self.location_coord[UAVLocation.BASE]

        # Update the member variables storing all the figure objects
        self.uav_circ_vis = [mpatches.Circle(
            (uav_x,
             1 + uav_id),
            self.UAV_RADIUS,
            fc="w") for uav_id in range(0,
                                        self.NUM_UAV)]
        self.uav_text_vis = [None for uav_id in range(0, self.NUM_UAV)]  # cuem
        self.uav_sensor_vis = [mpatches.Wedge(
            (uav_x + self.SENSOR_REL_X,
             1 + uav_id),
            self.SENSOR_LENGTH,
            -30,
            30) for uav_id in range(0,
                                    self.NUM_UAV)]
        self.uav_actuator_vis = [mpatches.Wedge(
            (uav_x,
             1 + uav_id + self.ACTUATOR_REL_Y),
            self.ACTUATOR_HEIGHT,
            60,
            120) for uav_id in range(0,
                                     self.NUM_UAV)]

        # The following was executed when we used to check if the environment needed re-drawing: see above.
                # Remove all UAV circle objects from visualization
        #        else:
        #            [self.uav_circ_vis[uav_id].remove() for uav_id in range(0,self.NUM_UAV)]
        #            [self.uav_text_vis[uav_id].remove() for uav_id in range(0,self.NUM_UAV)]
        #            [self.uav_sensor_vis[uav_id].remove() for uav_id in range(0,self.NUM_UAV)]

        # For each UAV:
        # Draw a circle, with text inside = amt fuel remaining
        # Triangle on top of UAV for comms, black = good, red = bad
        # Triangle in front of UAV for surveillance
        sStruct = self.state2Struct(s)

        for uav_id in range(0, self.NUM_UAV):
            # Assign all the variables corresponding to this UAV for this iteration;
            # this could alternately be done with a UAV class whose objects keep track
            # of these variables.  Elect to use lists here since ultimately the state
            # must be a vector anyway.
            # State index corresponding to the location of this uav
            uav_location = sStruct.locations[uav_id]
            uav_fuel = sStruct.fuel[uav_id]
            uav_sensor = sStruct.sensor[uav_id]
            uav_actuator = sStruct.actuator[uav_id]

            # Assign coordinates on figure where UAV should be drawn
            uav_x = self.location_coord[uav_location]
            uav_y = 1 + uav_id

            # Update plot wit this UAV
            self.uav_circ_vis[uav_id] = mpatches.Circle(
                (uav_x, uav_y), self.UAV_RADIUS, fc="w")
            self.uav_text_vis[uav_id] = plt.text(
                uav_x - 0.05,
                uav_y - 0.05,
                uav_fuel)
            if uav_sensor == SensorState.RUNNING:
                objColor = 'black'
            else:
                objColor = 'red'
            self.uav_sensor_vis[uav_id] = mpatches.Wedge(
                (uav_x + self.SENSOR_REL_X,
                 uav_y),
                self.SENSOR_LENGTH,
                -30,
                30,
                color=objColor)

            if uav_actuator == ActuatorState.RUNNING:
                objColor = 'black'
            else:
                objColor = 'red'
            self.uav_actuator_vis[uav_id] = mpatches.Wedge(
                (uav_x,
                 uav_y + self.ACTUATOR_REL_Y),
                self.ACTUATOR_HEIGHT,
                60,
                120,
                color=objColor)

            self.subplot_axes.add_patch(self.uav_circ_vis[uav_id])
            self.subplot_axes.add_patch(self.uav_sensor_vis[uav_id])
            self.subplot_axes.add_patch(self.uav_actuator_vis[uav_id])

        numHealthySurveil = np.sum(
            np.logical_and(
                sStruct.locations == UAVLocation.SURVEIL,
                sStruct.sensor))
        # We have comms coverage: draw a line between comms states to show this
        if (any(sStruct.locations == UAVLocation.COMMS)):
            for i in xrange(len(self.comms_line)):
                self.comms_line[i].set_visible(True)
                self.comms_line[i].set_color('black')
                self.subplot_axes.add_line(self.comms_line[i])
            # We also have UAVs in surveillance; color the comms line black
            if numHealthySurveil > 0:
                self.location_rect_vis[
                    len(self.location_rect_vis) - 1].set_color('green')
        plt.draw()
        sleep(0.5)

Example 47

Project: scipy Source File: bsplines.py
def _bspline_piecefunctions(order):
    """Returns the function defined over the left-side pieces for a bspline of
    a given order.

    The 0th piece is the first one less than 0.  The last piece is a function
    identical to 0 (returned as the constant 0).  (There are order//2 + 2 total
    pieces).

    Also returns the condition functions that when evaluated return boolean
    arrays for use with `numpy.piecewise`.
    """
    try:
        return _splinefunc_cache[order]
    except KeyError:
        pass

    def condfuncgen(num, val1, val2):
        if num == 0:
            return lambda x: logical_and(less_equal(x, val1),
                                         greater_equal(x, val2))
        elif num == 2:
            return lambda x: less_equal(x, val2)
        else:
            return lambda x: logical_and(less(x, val1),
                                         greater_equal(x, val2))

    last = order // 2 + 2
    if order % 2:
        startbound = -1.0
    else:
        startbound = -0.5
    condfuncs = [condfuncgen(0, 0, startbound)]
    bound = startbound
    for num in xrange(1, last - 1):
        condfuncs.append(condfuncgen(1, bound, bound - 1))
        bound = bound - 1
    condfuncs.append(condfuncgen(2, 0, -(order + 1) / 2.0))

    # final value of bound is used in piecefuncgen below

    # the functions to evaluate are taken from the left-hand-side
    #  in the general expression derived from the central difference
    #  operator (because they involve fewer terms).

    fval = factorial(order)

    def piecefuncgen(num):
        Mk = order // 2 - num
        if (Mk < 0):
            return 0  # final function is 0
        coeffs = [(1 - 2 * (k % 2)) * float(comb(order + 1, k, exact=1)) / fval
                  for k in xrange(Mk + 1)]
        shifts = [-bound - k for k in xrange(Mk + 1)]

        def thefunc(x):
            res = 0.0
            for k in range(Mk + 1):
                res += coeffs[k] * (x + shifts[k]) ** order
            return res
        return thefunc

    funclist = [piecefuncgen(k) for k in xrange(last)]

    _splinefunc_cache[order] = (funclist, condfuncs)

    return funclist, condfuncs

Example 48

Project: BioSPPy Source File: tools.py
def band_power(freqs=None, power=None, frequency=None, decibel=True):
    """Compute the avearge power in a frequency band.

    Parameters
    ----------
    freqs : array
        Array of frequencies (Hz) at which the power was computed.
    power : array
        Input power spectrum.
    frequency : list, array
        Pair of frequencies defining the band.
    decibel : bool, optional
        If True, input power is in decibels.

    Returns
    -------
    avg_power : float
        The average power in the band.

    """

    # check inputs
    if freqs is None:
        raise TypeError("Please specify the 'freqs' array.")

    if power is None:
        raise TypeError("Please specify the input power spectrum.")

    if len(freqs) != len(power):
        raise ValueError(
            "The input 'freqs' and 'power' arrays must have the same length.")

    if frequency is None:
        raise TypeError("Please specify the band frequencies.")

    try:
        f1, f2 = frequency
    except ValueError:
        raise ValueError("Input 'frequency' must be a pair of frequencies.")

    # make frequencies sane
    if f1 > f2:
        f1, f2 = f2, f1

    if f1 < freqs[0]:
        f1 = freqs[0]
    if f2 > freqs[-1]:
        f2 = freqs[-1]

    # average
    sel = np.nonzero(np.logical_and(f1 <= freqs, freqs <= f2))[0]

    if decibel:
        aux = 10 ** (power / 10.)
        avg = np.mean(aux[sel])
        avg = 10. * np.log10(avg)
    else:
        avg = np.mean(power[sel])

    return utils.ReturnTuple((avg,), ('avg_power',))

Example 49

Project: mpop Source File: tools.py
def estimate_cth(IR_108, cth_atm="standard"):

    '''
    Estimation of the cloud top height using the 10.8 micron channel
    limitations: this is the most simple approach
                 a simple fit of the ir108 to the temperature profile
                 * no correction for water vapour or any other trace gas
                 * no viewing angle dependency
                 * no correction for semi-transparent clouds

    optional input:
      cth_atm    * "standard", "tropics", "midlatitude summer", "midlatitude winter", "subarctic summer", "subarctic winter"
                  Matching the 10.8 micron temperature with atmosphere profile
                  (s)  AFGL atmospheric constituent profile. U.S. standard atmosphere 1976. (AFGL-TR-86-0110) 
                  (t)  AFGL atmospheric constituent profile. tropical.                      (AFGL-TR-86-0110)
                  (mw) AFGL atmospheric constituent profile. midlatitude summer.            (AFGL-TR-86-0110) 
                  (ms) AFGL atmospheric constituent profile. midlatitude winter.            (AFGL-TR-86-0110)
                  (ss) AFGL atmospheric constituent profile. subarctic summer.              (AFGL-TR-86-0110) 
                  (sw) AFGL atmospheric constituent profile. subarctic winter.              (AFGL-TR-86-0110)
                  Ulrich Hamann (MeteoSwiss)
                * "tropopause"
                  Assuming a fixed tropopause height and a fixed temperature gradient
                  Richard Mueller (DWD)
    output: 
      parallax corrected channel
                 the content of the channel will be parallax corrected.
                 The name of the new channel will be
                 *original_chan.name+'_PC'*, eg. "IR_108_PC". This name is
                 also stored to the info dictionary of the originating channel.

    Versions: 05.07.2016 initial version
              Ulrich Hamann (MeteoSwiss), Richard Mueller (DWD)
    '''

    print "*** estimating CTH using the 10.8 micro meter brightness temperature "

    if cth_atm.lower() != "tropopause":

        # define atmospheric temperature profile    
        import os
        from numpy import loadtxt, zeros, where, logical_and
        import mpop 

        mpop_dir = os.path.dirname(mpop.__file__)
        afgl_file = mpop_dir+"/afgl.dat"
        print "... assume ", cth_atm, " atmosphere for temperature profile"

        if cth_atm.lower()=="standard" or cth_atm.lower()=="s":
            z, T = loadtxt(afgl_file, usecols=(0, 1), unpack=True, comments="#")
        elif cth_atm.lower()=="tropics" or cth_atm.lower()=="t":
            z, T = loadtxt(afgl_file, usecols=(0, 2), unpack=True, comments="#")
        elif cth_atm.lower()=="midlatitude summer" or cth_atm.lower()=="ms":
            z, T = loadtxt(afgl_file, usecols=(0, 3), unpack=True, comments="#")
        elif cth_atm.lower()=="midlatitude winter" or cth_atm.lower()=="ws":
            z, T = loadtxt(afgl_file, usecols=(0, 4), unpack=True, comments="#")
        elif cth_atm.lower()=="subarctic summer" or cth_atm.lower()=="ss":
            z, T = loadtxt(afgl_file, usecols=(0, 5), unpack=True, comments="#")
        elif cth_atm.lower()=="subarctic winter" or cth_atm.lower()=="ss":
            z, T = loadtxt(afgl_file, usecols=(0, 6), unpack=True, comments="#")
        else:
            print "*** Error in estimate_cth (mpop/tools.py)"
            print "unknown temperature profiel for CTH estimation: cth_atm = ", cth_atm
            quit()

        height = zeros(IR_108.shape)
        # warmer than lowest level -> clear sky 
        height[where(IR_108 > T[-1])] = -1.
        print "     z0(km)   z1(km)   T0(K)   T1(K)  number of pixels"
        print "------------------------------------------------------"
        for i in range(z.size)[::-1]:

            # search for temperatures between layer i-1 and i
            ind =  np.where( logical_and( T[i-1]< IR_108, IR_108 < T[i]) )
            # interpolate CTH according to ir108 temperature
            height[ind] = z[i] + (IR_108[ind]-T[i])/(T[i-1]-T[i]) * (z[i-1]-z[i])
            # verbose output
            print " {0:8.1f} {1:8.1f} {2:8.1f} {3:8.1f} {4:8d}".format(z[i], z[i-1], T[i], T[i-1], len(ind[0]))

            # if temperature increases above 8km -> tropopause detected
            if z[i]>=8. and T[i] <= T[i-1]:
                # no cloud above tropopose
                break
            # no cloud heights above 20km
            if z[i]>=20.:
                break

        # if height is still 0 -> cloud colder than tropopause -> cth == tropopause height
        height[np.where( height == 0 )] = z[i]
        
    else:

        Htropo=11.0 # km
        # this is an assumption it should be optimized 
        # by making it dependent on region and season. 
        # It might be good to include the ITC in the  
        # region of interest, that would make a fixed Htropo 
        # value more reliable. 
        Tmin = np.amin(IR_108) 
        # for Tmin it might be better to use the 5th or 10th percentile 
        # else overshoting tops induces further uncertainties  
        # in the calculation of the cloud height. 
        # However numpy provides weird results for 5th percentile. 
        # Hence, for the working version the minima is used 

        print "... assume tropopause height ", Htropo, ", tropopause temperature ", Tmin, "K (", Tmin-273.16, "deg C)"
        print "    and constant temperature gradient 6.5 K/km"

        height = -(IR_108 - Tmin)/6.5 + Htropo 
        # calculation of the height, the temperature gradient 
        # 6.5 K/km is an assumption  
        # derived from USS and MPI standard profiles. It 
        # has to be improved as well 

    # convert to masked array
    # convert form km to meter
    height = np.ma.masked_where(height <= 0, height, copy=False) * 1000.

    if False:
        from trollimage.image import Image as trollimage
        from trollimage.colormap import rainbow
        from copy import deepcopy 
        # cloud top height
        prop = height
        min_data = prop.min()
        max_data = prop.max()
        print " estimated CTH(meter) (min/max): ", min_data, max_data
        min_data =     0
        max_data = 12000    
        colormap = deepcopy(rainbow)
        colormap.set_range(min_data, max_data)
        img = trollimage(prop, mode="L") #, fill_value=[0,0,0]
        img.colorize(colormap)
        img.show()

    # return cloud top height in meter
    return height

Example 50

Project: sfepy Source File: bspline.py
    def insert_knot(self, new):
        """
        Insert a new knot into the knot vector.

        Parameters
        ----------
        new : float
            The new knot value.
        """
        kn = self.knots
        dg = self.degree
        ncp = self.ncp
        cp = self.cp_coors
        idx = nm.where(nm.logical_and(kn[:-1] <= new, new < kn[1:]))[0]

        if len(idx) > 0:
            multi = len(nm.where(kn == new)[0])
            if multi < dg:
                # new knot
                newkn = nm.zeros((len(kn) + 1,), dtype=nm.float64)
                newkn[:(idx + 1)] = kn[:(idx + 1)]
                newkn[idx + 1] = new
                newkn[(idx + 2):] = kn[(idx + 1):]
                u1 = idx - dg + 1
                u2 = idx + 1

                # new control points
                newcp = nm.zeros((ncp + 1, cp.shape[1]), dtype=nm.float64)
                newcp[:u1,:] = cp[:u1,:]
                newcp[u2:,:] = cp[(u2 - 1):,:]

                for ii in range(u1, u2):
                    kn1 = kn[ii]
                    kn2 = kn[ii + dg]
                    dd = kn2 - kn1
                    newcp[ii,:] = (kn2 - new) / dd * cp[ii - 1] + \
                                  (new - kn1) / dd * cp[ii]

                self.knots = newkn
                self.cp_coors = newcp
                self.ncp = newcp.shape[0]

                # evaluate B-spline base functions for new configuration
                self.eval_basis()

            else:
                print('knot insertion failed: multiplicity = spline degree!')

        else:
            print('knot insertion failed: out of bounds!')
See More Examples - Go to Next Page
Page 1 Selected Page 2 Page 3