Python numpy.polyfit方法代码示例

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def fit_loglog(x, y):
    """
    Fit a line to isotropic spectra in log-log space

    Parameters
    ----------
    x : `numpy.array`
        Coordinate of the data
    y : `numpy.array`
        data

    Returns
    -------
    y_fit : `numpy.array`
        The linear fit
    a : float64
        Slope of the fit
    b : float64
        Intercept of the fit
    """
    # fig log vs log
    p = np.polyfit(np.log2(x), np.log2(y), 1)
    y_fit = 2**(np.log2(x)*p[0] + p[1])

    return y_fit, p[0], p[1] 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def remove_linear_BG_XAS_preedge(
    xmcd_data, scanparams, process_parameters=None, process_number=-1
):
    """Should remove a linear bg based on the preedge average"""
    preedge_spectrum = get_preedge_spectrum(process_parameters, xmcd_data)

    preedge_poly = np.poly1d(
        np.polyfit(preedge_spectrum["Energy"], preedge_spectrum["XAS"], 1)
    )

    xas_bg = preedge_poly(xmcd_data["Energy"])

    for xas in ["XAS+", "XAS-", "XAS"]:
        xmcd_data[xas] -= xas_bg

    return (xmcd_data, {"xas_bg_poly_coeffs": " ".join(map(str, preedge_poly.coeffs))}) 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def stable_fit(xfit, yfit):
    p = np.polyfit(xfit, yfit, 2)
    steprange = np.max(xfit)
    minstep = np.min(xfit)
    est_min = -p[1] / (2 * p[0])
    if est_min > steprange and p[0] > 0:  # minimum past the search radius
        est_min = steprange
    if est_min < minstep and p[0] > 0:  # mimimum behind the search radius
        est_min = minstep
    if p[0] < 0:
        plin = np.polyfit(xfit, yfit, 1)
        if plin[0] < 0:
            est_min = steprange
        if plin[0] > 0:
            est_min = minstep
    # print("estimated minimum adjusted", est_min, flush=True)
    return est_min 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def quad_interpol_argmax(y, x=None):
        """
        Find argmax for quadratic interpolation around argmax of y.

        :param x: x corresponding to (a) peak in y, if not set, ``np.argmax(y)`` is used
        :param y: array
        :return: float (index) of interpolated max, strength
        """
        if x is None:
            x = np.argmax(y)
        if x == 0 or x == y.shape[0] - 1:
            return x, y[x]
        z = np.polyfit([x - 1, x, x + 1], [y[x - 1], y[x], y[x + 1]], 2)
        # find (float) x value for max
        argmax = -z[1] / (2. * z[0])
        height = z[2] - (z[1] ** 2.) / (4. * z[0])
        return argmax, height 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def horner(circulation_time, times, temp):
    '''
    horner_bht_temp (circulation_time, times, temp)
    *Input parameters:
    - circulation_time - hours from last circulation;
    - times - total time since circulation stopped at 1st Run, 2nd Run and so on ...
    - temp - a list o temperatures coresponding to 1st Run, 2nd Run and so on ...
    *Returns:
    - horner_temp - formation temperature estimated by Horner method (thermometer readings
    from different runs)
    *Exemple of usage:
    horner(6, (7.0,11.5,19.5), (100,105,108))
        where:
        circulation_time = 6           # time since circulation stopped (hours)
        times = (7.0,11.5,19.5)        # total time since circulation stopped at 1st, 2nd, 3rd RUN (hours)
        temp=(100,105,108)             # temperatures recorded at 1st, 2nd, 3rd RUN (Celsius degrees)
    '''
    horner_time = np.array(times) / (circulation_time + np.array(times))
    slope,intercept = np.polyfit (np.log(horner_time), temp, 1)
    horner_temp=round(slope*np.log(1) +intercept,2)
    return horner_temp 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def test_polyfit_build(self):
        # Ticket #628
        ref = [-1.06123820e-06, 5.70886914e-04, -1.13822012e-01,
               9.95368241e+00, -3.14526520e+02]
        x = [90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103,
             104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115,
             116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 129,
             130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141,
             146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157,
             158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169,
             170, 171, 172, 173, 174, 175, 176]
        y = [9.0, 3.0, 7.0, 4.0, 4.0, 8.0, 6.0, 11.0, 9.0, 8.0, 11.0, 5.0,
             6.0, 5.0, 9.0, 8.0, 6.0, 10.0, 6.0, 10.0, 7.0, 6.0, 6.0, 6.0,
             13.0, 4.0, 9.0, 11.0, 4.0, 5.0, 8.0, 5.0, 7.0, 7.0, 6.0, 12.0,
             7.0, 7.0, 9.0, 4.0, 12.0, 6.0, 6.0, 4.0, 3.0, 9.0, 8.0, 8.0,
             6.0, 7.0, 9.0, 10.0, 6.0, 8.0, 4.0, 7.0, 7.0, 10.0, 8.0, 8.0,
             6.0, 3.0, 8.0, 4.0, 5.0, 7.0, 8.0, 6.0, 6.0, 4.0, 12.0, 9.0,
             8.0, 8.0, 8.0, 6.0, 7.0, 4.0, 4.0, 5.0, 7.0]
        tested = np.polyfit(x, y, 4)
        assert_array_almost_equal(ref, tested) 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def testEncodeDecodeShift(self):
        x = np.linspace(-1, 1, 1000).astype(np.float32)
        with self.test_session() as sess:
            encoded = mu_law_encode(x, QUANT_LEVELS)
            decoded = mu_law_decode(encoded, QUANT_LEVELS)
            roundtripped = sess.run(decoded)

        # Detect non-unity scaling and non-zero shift in the roundtripped
        # signal by asserting that slope = 1 and y-intercept = 0 of line fit to
        # roundtripped vs x values.
        coeffs = np.polyfit(x, roundtripped, 1)
        slope = coeffs[0]
        y_intercept = coeffs[1]
        EPSILON = 1e-4
        self.assertNear(slope, 1.0, EPSILON)
        self.assertNear(y_intercept, 0.0, EPSILON) 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def data_analysis(e_ph, flux, method="least"):

    if method == "least":
        coeffs = np.polyfit(x=e_ph, y=flux, deg=11)
        polynom = np.poly1d(coeffs)


        x = np.linspace(e_ph[0], e_ph[-1], num=100)
        pd = np.polyder(polynom, m=1)
        indx = np.argmax(np.abs(pd(x)))
        eph_c = x[indx]

        pd2 = np.polyder(polynom, m=2)
        p2_roots = np.roots(pd2)
        p2_roots = p2_roots[p2_roots[:].imag == 0]
        p2_roots = np.real(p2_roots)
        Eph_fin = find_nearest(p2_roots,eph_c)
        return Eph_fin, polynom

    elif method == "new method":
        pass

        #plt.plot(Etotal, total, "ro")
        #plt.plot(x, polynom(x)) 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def quad_fit(x, y, yerr, name, unit):
    """ Fit a qudratic to the SNR to make a lookup error table """
    print("performing quad fit")
    qfit = np.polyfit(x, y, deg=2, w = 1 / yerr)
    print(qfit)
    plt.figure()
    #plt.scatter(x, y)
    plt.errorbar(x, y, yerr=yerr, fmt='.', c='k')
    xvals = np.linspace(min(x), max(x), 100)
    print(xvals)
    yvals = qfit[2] + qfit[1]*xvals + qfit[0]*xvals**2
    print(yvals)
    plt.plot(xvals, yvals, color='r', lw=2)
    plt.xlabel("%s" %snr_label, fontsize=16)
    plt.ylabel(r"$\sigma %s \mathrm{(%s)}$" %(name,unit), fontsize=16)
    plt.show() 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def test_polyfit_build(self):
        """Ticket #628"""
        ref = [-1.06123820e-06, 5.70886914e-04, -1.13822012e-01,
                9.95368241e+00, -3.14526520e+02]
        x = [90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103,
             104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115,
             116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 129,
             130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141,
             146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157,
             158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169,
             170, 171, 172, 173, 174, 175, 176]
        y = [9.0, 3.0, 7.0, 4.0, 4.0, 8.0, 6.0, 11.0, 9.0, 8.0, 11.0, 5.0,
             6.0, 5.0, 9.0, 8.0, 6.0, 10.0, 6.0, 10.0, 7.0, 6.0, 6.0, 6.0,
             13.0, 4.0, 9.0, 11.0, 4.0, 5.0, 8.0, 5.0, 7.0, 7.0, 6.0, 12.0,
             7.0, 7.0, 9.0, 4.0, 12.0, 6.0, 6.0, 4.0, 3.0, 9.0, 8.0, 8.0,
             6.0, 7.0, 9.0, 10.0, 6.0, 8.0, 4.0, 7.0, 7.0, 10.0, 8.0, 8.0,
             6.0, 3.0, 8.0, 4.0, 5.0, 7.0, 8.0, 6.0, 6.0, 4.0, 12.0, 9.0,
             8.0, 8.0, 8.0, 6.0, 7.0, 4.0, 4.0, 5.0, 7.0]
        tested = np.polyfit(x, y, 4)
        assert_array_almost_equal(ref, tested) 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def fit_1d(self, acc):
        nx_list = 10 ** np.linspace(1.75, 2., 10)
        Lx = np.pi

        log_err_list = []
        log_dx_list = []

        for nx in nx_list:
            x = np.linspace(0., Lx, nx)
            dx = x[1] - x[0]
            f = np.sin(x)
            d_dx = FinDiff(0, dx)
            fx = d_dx(f, acc=acc)
            fxe = np.cos(x)
            err = np.max(np.abs(fxe - fx))
            log_dx_list.append(log(dx))
            log_err_list.append(log(err))

        fit = np.polyfit(log_dx_list, log_err_list, deg=1)
        return fit[0] 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def fit_2d(self, acc):
        nx_list = [10, 30, 100, 300]
        ny_list = [10, 30, 100, 300]
        Lx, Ly = 3, 3

        log_err_list = []
        log_dx_list = []

        for nx, ny in zip(nx_list, ny_list):
            x = np.linspace(0, Lx, nx)
            y = np.linspace(0, Ly, ny)
            dx, dy = x[1] - x[0], y[1] - y[0]
            X, Y = np.meshgrid(x, y, indexing='ij')
            f = np.sin(X) * np.sin(Y)
            d_dx = FinDiff(0, dx)
            fx = d_dx(f, acc=acc)
            fxe = np.cos(X) * np.sin(Y)
            err = np.max(np.abs(fxe - fx))
            log_dx_list.append(log(dx))
            log_err_list.append(log(err))

        fit = np.polyfit(log_dx_list, log_err_list, deg=1)
        return fit[0] 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def calculate_rmsd(data, time=None, slope=False):
    assert isinstance(data, np.ndarray) and data.ndim == 1
    assert time is None or isinstance(time, np.ndarray) and time.ndim == 1

    avg = data.mean()

    if time is None:
        time = np.arange(data.size)

    fit = np.polyfit(time, data, 1)

    def f(x):
        return fit[0] * x + fit[1]

    if slope:
        rmsd = 0
        for t, d in zip(time, data):
            rmsd += (d - f(t)) ** 2
        rmsd = np.sqrt(rmsd / data.size)
    else:
        rmsd = data.std()

    return avg, rmsd, fit[0] 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def polyfit(x, y, degree):
    """
    Wrapper around np.polyfit

    :param x: x values
    :param y: y values
    :param degree: polynomial degree
    :return: results, polynomial has coefficients, determination has r-squared
    :rtype: dict
    """
    results = {}
    coeffs = np.polyfit(x, y, degree)
    results['polynomial'] = coeffs.tolist()
    # r squared
    p = np.poly1d(coeffs)
    yhat = p(x)
    ybar = np.sum(y) / len(y)
    ssreg = np.sum((yhat - ybar) ** 2)
    sstot = np.sum((y - ybar) ** 2)
    results['determination'] = ssreg / sstot
    return results 

# 需要导入模块: import numpy [as 别名]
# 或者: from numpy import polyfit [as 别名]
def vander(points, deg):
    """N-dim Vandermonde matrix for data `points` and a polynomial of degree
    `deg`.

    Parameters
    ----------
    points : see polyfit()
    deg : int
        Degree of the poly (e.g. 3 for cubic).

    Returns
    -------
    vander : 2d array (npoints, (deg+1)**ndim)
    """
    powers = poly_powers(points.shape[1], deg)
    # low memory version, slower
    ##npoints = points.shape[0]
    ##vand = np.empty((npoints, (deg+1)**ndim), dtype=float)
    ##for ipoint in range(npoints):
    ##    vand[ipoint,:] = (points[ipoint]**powers).prod(axis=1)
    tmp = (points[...,None] ** np.swapaxes(powers, 0, 1)[None,...])
    return tmp.prod(axis=1)