quadratic interpolation
509517
private static int interpolate(
        double y_minus, // = f(-1)
        double y_0, // = f(0)
        double y_plus, // = f(+1)
        double[] result // xe, ye, root1, root2
    ) {

        double a = 0.5 * (y_plus + y_minus) - y_0;
        double b = 0.5 * (y_plus - y_minus);

        double xe = -b / (2.0 * a);
        double ye = (a * xe + b) * xe + y_0;
        double dis = b * b - 4 * a * y_0;

        double root1 = Double.NaN;
        double root2 = Double.NaN;
        int count = 0;

        if (dis >= 0) {
            double dx = 0.5 * Math.sqrt(dis) / Math.abs(a);
            if (Math.abs(xe - dx) <= 1.0) {
                root1 = xe - dx;
                count++;
            }
            if (Math.abs(xe + dx) <= 1.0) {
                root2 = xe + dx;
                count++;
            }
        }

        result[0] = xe;
        result[1] = ye;
        result[2] = root1;
        result[3] = root2;

        return count;

    }
92772
private Vertex interpolate(Vertex start, Vertex end, double t)
	{
		double x = interpolate(start.getX(), end.getX(), t);
		double y = interpolate(start.getY(), end.getY(), t);
		double z = interpolate(start.getZ(), end.getZ(), t);

		return new Vertex(x, y, z);
	}
1275000
def interpolate_1d(xvalues, yvalues, method='linear', limit=None,
                   limit_direction='forward', limit_area=None, fill_value=None,
                   bounds_error=False, order=None, **kwargs):
    
    # Treat the original, non-scipy methods first.

    invalid = isna(yvalues)
    valid = ~invalid

    if not valid.any():
        # have to call np.asarray(xvalues) since xvalues could be an Index
        # which can't be mutated
        result = np.empty_like(np.asarray(xvalues), dtype=np.float64)
        result.fill(np.nan)
        return result

    if valid.all():
        return yvalues

    if method == 'time':
        if not getattr(xvalues, 'is_all_dates', None):
            # if not issubclass(xvalues.dtype.type, np.datetime64):
            raise ValueError('time-weighted interpolation only works '
                             'on Series or DataFrames with a '
                             'DatetimeIndex')
        method = 'values'

    valid_limit_directions = ['forward', 'backward', 'both']
    limit_direction = limit_direction.lower()
    if limit_direction not in valid_limit_directions:
        msg = ('Invalid limit_direction: expecting one of {valid!r}, '
               'got {invalid!r}.')
        raise ValueError(msg.format(valid=valid_limit_directions,
                                    invalid=limit_direction))

    if limit_area is not None:
        valid_limit_areas = ['inside', 'outside']
        limit_area = limit_area.lower()
        if limit_area not in valid_limit_areas:
            raise ValueError('Invalid limit_area: expecting one of {}, got '
                             '{}.'.format(valid_limit_areas, limit_area))

    # default limit is unlimited GH #16282
    if limit is None:
        # limit = len(xvalues)
        pass
    elif not is_integer(limit):
        raise ValueError('Limit must be an integer')
    elif limit < 1:
        raise ValueError('Limit must be greater than 0')

    from pandas import Series
    ys = Series(yvalues)

    # These are sets of index pointers to invalid values... i.e. {0, 1, etc...
    all_nans = set(np.flatnonzero(invalid))
    start_nans = set(range(ys.first_valid_index()))
    end_nans = set(range(1 + ys.last_valid_index(), len(valid)))
    mid_nans = all_nans - start_nans - end_nans

    # Like the sets above, preserve_nans contains indices of invalid values,
    # but in this case, it is the final set of indices that need to be
    # preserved as NaN after the interpolation.

    # For example if limit_direction='forward' then preserve_nans will
    # contain indices of NaNs at the beginning of the series, and NaNs that
    # are more than'limit' away from the prior non-NaN.

    # set preserve_nans based on direction using _interp_limit
    if limit_direction == 'forward':
        preserve_nans = start_nans | set(_interp_limit(invalid, limit, 0))
    elif limit_direction == 'backward':
        preserve_nans = end_nans | set(_interp_limit(invalid, 0, limit))
    else:
        # both directions... just use _interp_limit
        preserve_nans = set(_interp_limit(invalid, limit, limit))

    # if limit_area is set, add either mid or outside indices
    # to preserve_nans GH #16284
    if limit_area == 'inside':
        # preserve NaNs on the outside
        preserve_nans |= start_nans | end_nans
    elif limit_area == 'outside':
        # preserve NaNs on the inside
        preserve_nans |= mid_nans

    # sort preserve_nans and covert to list
    preserve_nans = sorted(preserve_nans)

    xvalues = getattr(xvalues, 'values', xvalues)
    yvalues = getattr(yvalues, 'values', yvalues)
    result = yvalues.copy()

    if method in ['linear', 'time', 'index', 'values']:
        if method in ('values', 'index'):
            inds = np.asarray(xvalues)
            # hack for DatetimeIndex, #1646
            if needs_i8_conversion(inds.dtype.type):
                inds = inds.view(np.int64)
            if inds.dtype == np.object_:
                inds = lib.maybe_convert_objects(inds)
        else:
            inds = xvalues
        result[invalid] = np.interp(inds[invalid], inds[valid], yvalues[valid])
        result[preserve_nans] = np.nan
        return result

    sp_methods = ['nearest', 'zero', 'slinear', 'quadratic', 'cubic',
                  'barycentric', 'krogh', 'spline', 'polynomial',
                  'from_derivatives', 'piecewise_polynomial', 'pchip', 'akima']

    if method in sp_methods:
        inds = np.asarray(xvalues)
        # hack for DatetimeIndex, #1646
        if issubclass(inds.dtype.type, np.datetime64):
            inds = inds.view(np.int64)
        result[invalid] = _interpolate_scipy_wrapper(inds[valid],
                                                     yvalues[valid],
                                                     inds[invalid],
                                                     method=method,
                                                     fill_value=fill_value,
                                                     bounds_error=bounds_error,
                                                     order=order, **kwargs)
        result[preserve_nans] = np.nan
        return result
790988
def interpolate(self, gmvs):
        
        # gmvs are clipped to max(iml)
        gmvs_curve = numpy.piecewise(
            gmvs, [gmvs > self.imls[-1]], [self.imls[-1], lambda x: x])
        idxs = gmvs_curve >= self.imls[0]  # indices over the minimum
        gmvs_curve = gmvs_curve[idxs]
        return self._mlr_i1d(gmvs_curve), self._cov_for(gmvs_curve), idxs
712185
def interpolate_xml_array(data, low_res_coords, shape, chunks):
        """"""
        xpoints, ypoints = low_res_coords

        return interpolate_xarray_linear(xpoints, ypoints, data, shape, chunks=chunks)
636848
function QuaternionLinearInterpolant( parameterPositions, sampleValues, sampleSize, resultBuffer ) {

	Interpolant.call( this, parameterPositions, sampleValues, sampleSize, resultBuffer );

}
854983
def interpolation_points(self, N):
        
        if N == 1:
            return np.array([0.])
        return np.cos(np.arange(N)*np.pi/(N-1))
1202257
def setup_interpolant(self):
        """"""
        # for computing nearby (z < 1) redshifts
        zs = numpy.linspace(0., 1., num=self.numpoints)
        ds = self.cosmology.luminosity_distance(zs).value
        self.nearby_d2z = interpolate.interp1d(ds, zs, kind='linear',
                                                bounds_error=False)
        # for computing far away (z > 1) redshifts
        zs = numpy.logspace(0, numpy.log10(self.default_maxz),
                            num=self.numpoints)
        ds = self.cosmology.luminosity_distance(zs).value
        self.faraway_d2z = interpolate.interp1d(ds, zs, kind='linear',
                                                 bounds_error=False)
        # store the default maximum distance
        self.default_maxdist = ds.max()
549536
protected function interpolate($x, $y, $nw, $ne, $sw, $se)
    {
        list($r0, $g0, $b0) = $this->getRGB($nw);
        list($r1, $g1, $b1) = $this->getRGB($ne);
        list($r2, $g2, $b2) = $this->getRGB($sw);
        list($r3, $g3, $b3) = $this->getRGB($se);

        $cx = 1.0 - $x;
        $cy = 1.0 - $y;

        $m0 = $cx * $r0 + $x * $r1;
        $m1 = $cx * $r2 + $x * $r3;
        $r  = (int) ($cy * $m0 + $y * $m1);

        $m0 = $cx * $g0 + $x * $g1;
        $m1 = $cx * $g2 + $x * $g3;
        $g  = (int) ($cy * $m0 + $y * $m1);

        $m0 = $cx * $b0 + $x * $b1;
        $m1 = $cx * $b2 + $x * $b3;
        $b  = (int) ($cy * $m0 + $y * $m1);

        return ($r << 16) | ($g << 8) | $b;
    }
335512
def quadratic_2d(data):
    
    arg_data_max = np.argmax(data)
    i, j = np.unravel_index(arg_data_max, data.shape)
    z_ = data[i-1:i+2, j-1:j+2]
    # our quadratic function is defined as
    # f(x, y | a, b, c, d, e, f) := a + b * x + c * y + d * x^2 + e * xy + f * y^2
    # therefore, the best fit coeffiecients are given as
    # note that they are unique and the uncertainty in each of them (#TODO) can be
    # computed following the derivations done by Vakili & Hogg (2016) and
    # Teague & Foreman-Mackey (2018)
    try:
        a = (-z_[0,0] + 2*z_[0,1] - z_[0,2] + 2*z_[1,0] + 5*z_[1,1] + 2*z_[1,2] -
             z_[2,0] + 2*z_[2,1] - z_[2,2]) / 9
        b = (-z_[0,0] - z_[0,1] - z_[0,2] + z_[2,0] + z_[2,1] + z_[2,2]) / 6
        c = (-z_[0,0] + z_[0,2] - z_[1,0] + z_[1,2] - z_[2,0] + z_[2,2]) / 6
        d = (z_[0,0] + z_[0,1] + z_[0,2] - z_[1,0]*2 - z_[1,1]*2 - z_[1,2]*2 +
             z_[2,0] + z_[2,1] + z_[2,2])/6
        e = (z_[0,0] - z_[0,2] - z_[2,0] + z_[2,2]) * .25
        f = (z_[0,0] - 2 * z_[0,1] + z_[0,2] + z_[1,0] - 2 * z_[1,1] + z_[1,2] +
             z_[2,0] - 2 * z_[2,1] + z_[2,2]) / 6
    except IndexError:
        return (i, j)

    # see https://en.wikipedia.org/wiki/Quadratic_function
    det = 4 * d * f - e ** 2
    xm = - (2 * f * b - c * e) / det
    ym = - (2 * d * c - b * e) / det
    return (i+xm, j+ym)
174016
def interpolate_jagged(xyz,nseg):
    
    
    # Spherical coordinates specifying the angles of all line
    # segments that make up the section path
    (r,theta,phi) = sequential_spherical(xyz)
    
    # cumulative length of section path at each coordinate
    rcum = np.append(0,np.cumsum(r))

    # breakpoints for segment paths along section path
    breakpoints = np.linspace(0,rcum[-1],nseg+1)
    np.delete(breakpoints,0)
    
    # Find segment paths
    seg_paths = []
    for a in range(nseg):
        path = []
        
        # find (x,y,z) starting coordinate of path
        if a == 0:
            start_coord = xyz[0,:]
        else:
            start_coord = end_coord # start at end of last path
        path.append(start_coord)

        # find all coordinates between the start and end points
        start_length = breakpoints[a]
        end_length = breakpoints[a+1]
        mid_boolean = (rcum > start_length) & (rcum < end_length)
        mid_indices = np.nonzero(mid_boolean)[0]
        for mi in mid_indices:
            path.append(xyz[mi,:])

        # find (x,y,z) ending coordinate of path
        end_coord = find_coord(end_length,xyz,rcum,theta,phi)
        path.append(end_coord)

        # Append path to list of segment paths
        seg_paths.append(np.array(path))
    
    # Return all segment paths
    return seg_paths
1160952
def _call(self, x, out=None):
        """"""
        # TODO: pass reasonable options on to the interpolator
        def linear(arg, out=None):
            """Interpolating function with vectorization."""
            if is_valid_input_meshgrid(arg, self.grid.ndim):
                input_type = 'meshgrid'
            else:
                input_type = 'array'

            interpolator = _LinearInterpolator(
                self.grid.coord_vectors, x, input_type=input_type)

            return interpolator(arg, out=out)

        return self.range.element(linear, vectorized=True)
36395
def linear_interpolate(p1, p2, fraction):
    ''''''
    return np.array((p1[0] + fraction * (p2[0] - p1[0]),
                     p1[1] + fraction * (p2[1] - p1[1]),
                     p1[2] + fraction * (p2[2] - p1[2])))
892141
def _linearInterpolationTransformMatrix(matrix1, matrix2, value):
    
    return tuple(_interpolateValue(matrix1[i], matrix2[i], value) for i in range(len(matrix1)))
1303067
def extrap1d(interpolator):
    
    
    xs = interpolator.x
    ys = interpolator.y

    def pointwise(x):
        if x < xs[0]:
            return ys[0] # +(x-xs[0])*(ys[1]-ys[0])/(xs[1]-xs[0])
        elif x > xs[-1]:
            return ys[-1] # +(x-xs[-1])*(ys[-1]-ys[-2])/(xs[-1]-xs[-2])
        else:
            return interpolator(x)

    def ufunclike(xs):
        return N.array(map(pointwise, N.array(xs)))

    return ufunclike
919547
def interpolate2dStructuredCrossAvg(grid, mask, kernel=15, power=2):
    '''
    

    '''

    vals = np.empty(shape=4, dtype=grid.dtype)
    dist = np.empty(shape=4, dtype=np.uint16)
    weights = np.empty(shape=4, dtype=np.float32)
    valid = np.empty(shape=4, dtype=bool)

    return _calc(grid, mask, power, kernel, vals, dist, weights, valid)
1066294
def _get_interpolator_nd(method, **kwargs):
    '''
    '''
    valid_methods = ['linear', 'nearest']

    try:
        from scipy import interpolate
    except ImportError:
        raise ImportError(
            'Interpolation with method `%s` requires scipy' % method)

    if method in valid_methods:
        kwargs.update(method=method)
        interp_class = interpolate.interpn
    else:
        raise ValueError('%s is not a valid interpolator for interpolating '
                         'over multiple dimensions.' % method)

    return interp_class, kwargs