Here are the examples of the python api numpy.roots taken from open source projects. By voting up you can indicate which examples are most useful and appropriate.
13 Examples
3
Example 1
def process_frame(X, window, num_formants, new_sr):
X = X * window
A, e, k = lpc(X, num_formants*2)
rts = np.roots(A)
rts = rts[np.where(np.imag(rts) >= 0)]
angz = np.arctan2(np.imag(rts), np.real(rts))
frqs = angz * (new_sr / (2 * np.pi))
frq_inds = np.argsort(frqs)
frqs = frqs[frq_inds]
bw = -1 / 2 * (new_sr / (2 * np.pi)) * np.log(np.abs(rts[frq_inds]))
return frqs, bw
2
Example 2
Project: python-control Source File: statesp.py
def zero(self):
"""Compute the zeros of a state space system."""
if self.inputs > 1 or self.outputs > 1:
raise NotImplementedError("StateSpace.zeros is currently \
implemented only for SISO systems.")
den = poly1d(poly(self.A))
# Compute the numerator based on zeros
#! TODO: This is currently limited to SISO systems
num = poly1d(poly(self.A - dot(self.B, self.C)) + ((self.D[0, 0] - 1) *
den))
return roots(num)
2
Example 3
Project: tyssue Source File: sheet_isotropic_model.py
def find_grad_roots(mod_specs):
poly = isotropic_grad_poly(mod_specs)
roots = np.roots(poly)
good_roots = np.real([r for r in roots if np.abs(r) == r])
np.sort(good_roots)
if len(good_roots) == 1:
return good_roots
elif len(good_roots) > 1:
return good_roots[0]
else:
return np.nan
2
Example 4
Project: rayopt Source File: gaussian_trace.py
@property
def eigenmodes(self):
n, m = self.system.paraxial_matrix(self.wavelength)
# FIXME only know how to do this for simple astigmatic matrices
# otherwise, solve qi*b*qi + qi*a - d*qi - c = 0
assert self.is_simple_astigmatic(m)
q = []
for axis in (0, 1):
a, b, c, d = m[axis::2, axis::2].flat
q.append(np.roots((c, d - a, -b)))
q = np.eye(2)[None, :]/np.array(q).T[:, :, None] # mode, axis
return q
0
Example 5
Project: pyomo Source File: colloc.py
def calc_cp(alpha,beta,k):
gamma = []
for i in range(k+1):
num = factorial(alpha+k)*factorial(alpha+beta+k+i)
denom = factorial(alpha+i)*factorial(k-i)*factorial(i)
gamma.insert(i,num/denom)
poly = []
for i in range(k+1):
if i == 0:
poly.insert(i,gamma[i])
else:
prod = [1]
j=1
while j<=i:
prod=conv(prod,[1,-1])
j=j+1
while len(poly)<len(prod):
poly.insert(0,0)
prod = [gamma[i]*t for t in prod]
poly = [sum(pair) for pair in zip(poly,prod)]
cp = numpy.roots(poly)
return cp
0
Example 6
Project: audiolazy Source File: lazy_lpc.py
def lsf(fir_filt):
"""
Find the Line Spectral Frequencies (LSF) from a given FIR filter.
Parameters
----------
filt :
A LTI FIR filter as a LinearFilter object.
Returns
-------
A tuple with all LSFs in rad/sample, alternating from the forward prediction
and backward prediction filters, starting with the lowest LSF value.
"""
den = fir_filt.denominator
if len(den) != 1:
raise ValueError("Filter has feedback")
elif den[0] != 1: # So we don't have to worry with the denominator anymore
fir_filt /= den[0]
from numpy import roots
rev_filt = ZFilter(fir_filt.numerator[::-1]) * z ** -1
P = fir_filt + rev_filt
Q = fir_filt - rev_filt
roots_p = roots(P.numerator[::-1])
roots_q = roots(Q.numerator[::-1])
lsf_p = sorted(phase(roots_p))
lsf_q = sorted(phase(roots_q))
return reduce(operator.concat, xzip(*sorted([lsf_p, lsf_q])), tuple())
0
Example 7
def irr(values):
"""
Return the Internal Rate of Return (IRR).
This is the rate of return that gives a net present value of 0.0.
Parameters
----------
values : array_like, shape(N,)
Input cash flows per time period. At least the first value would be
negative to represent the investment in the project.
Returns
-------
out : float
Internal Rate of Return for periodic input values.
Examples
--------
>>> np.irr([-100, 39, 59, 55, 20])
0.2809484211599611
"""
res = np.roots(values[::-1])
# Find the root(s) between 0 and 1
mask = (res.imag == 0) & (res.real > 0) & (res.real <= 1)
res = res[mask].real
if res.size == 0:
return np.nan
rate = 1.0/res - 1
if rate.size == 1:
rate = rate.item()
return rate
0
Example 8
def test_roots(self):
assert_array_equal(np.roots([1,0,0]), [0,0])
0
Example 9
Project: Mycodo Source File: method.py
def bezier_curve_y_out(shift_angle, P0, P1, P2, P3, second_of_day=None):
'''
For a cubic Bezier segment described by the 2-tuples P0, ..., P3, return
the y-value associated with the given x-value.
Ex: getYfromXforBezSegment((10,0), (5,-5), (5,5), (0,0), 3.2)
'''
seconds_per_day = 24*60*60
# Check if the second of the day is provided.
# If provided, calculate y of the Bezier curve with that x
# Otherwise, use the current second of the day
if second_of_day is None:
now = datetime.datetime.now()
dt = datetime.timedelta(hours=now.hour, minutes=now.minute, seconds=now.second)
seconds = dt.total_seconds()
else:
seconds = second_of_day
# Shift the entire graph using 0 - 360 to determine the degree
if shift_angle:
percent_angle = shift_angle/360
angle_seconds = percent_angle*seconds_per_day
if seconds+angle_seconds > seconds_per_day:
seconds_shifted = seconds+angle_seconds-seconds_per_day
else:
seconds_shifted = seconds+angle_seconds
percent_of_day = seconds_shifted/seconds_per_day
else:
percent_of_day = seconds/seconds_per_day
x = percent_of_day*(P0[0]-P3[0])
# First, get the t-value associated with x-value, where t is the
# parameterization of the Bezier curve and ranges from 0 to 1.
# We need the coefficients of the polynomial describing cubic Bezier
# (cubic polynomial in t)
coefficients = [-P0[0] + 3*P1[0] - 3*P2[0] + P3[0],
3*P0[0] - 6*P1[0] + 3*P2[0],
-3*P0[0] + 3*P1[0],
P0[0] - x]
# Find roots of the polynomial to determine the parameter t
roots = np.roots(coefficients)
# Find the root which is between 0 and 1, and is also real
correct_root = None
for root in roots:
if np.isreal(root) and 0 <= root <= 1:
correct_root = root
# Check a valid root was found
if correct_root is None:
print('Error, no valid root found. Are you sure your Bezier curve '
'represents a valid function when projected into the xy-plane?')
param_t = correct_root
# From the value for the t parameter, find the corresponding y-value
# using the formula for cubic Bezier curves
y = (1-param_t)**3*P0[1] + 3*(1-param_t)**2*param_t*P1[1] + 3*(1-param_t)*param_t**2*P2[1] + param_t**3*P3[1]
assert np.isreal(y)
# Typecast y from np.complex128 to float64
y = y.real
return y
0
Example 10
def irr(values):
"""
Return the Internal Rate of Return (IRR).
This is the "average" periodically compounded rate of return
that gives a net present value of 0.0; for a more complete explanation,
see Notes below.
Parameters
----------
values : array_like, shape(N,)
Input cash flows per time period. By convention, net "deposits"
are negative and net "withdrawals" are positive. Thus, for
example, at least the first element of `values`, which represents
the initial investment, will typically be negative.
Returns
-------
out : float
Internal Rate of Return for periodic input values.
Notes
-----
The IRR is perhaps best understood through an example (illustrated
using np.irr in the Examples section below). Suppose one invests 100
units and then makes the following withdrawals at regular (fixed)
intervals: 39, 59, 55, 20. Assuming the ending value is 0, one's 100
unit investment yields 173 units; however, due to the combination of
compounding and the periodic withdrawals, the "average" rate of return
is neither simply 0.73/4 nor (1.73)^0.25-1. Rather, it is the solution
(for :math:`r`) of the equation:
.. math:: -100 + \\frac{39}{1+r} + \\frac{59}{(1+r)^2}
+ \\frac{55}{(1+r)^3} + \\frac{20}{(1+r)^4} = 0
In general, for `values` :math:`= [v_0, v_1, ... v_M]`,
irr is the solution of the equation: [G]_
.. math:: \\sum_{t=0}^M{\\frac{v_t}{(1+irr)^{t}}} = 0
References
----------
.. [G] L. J. Gitman, "Principles of Managerial Finance, Brief," 3rd ed.,
Addison-Wesley, 2003, pg. 348.
Examples
--------
>>> round(irr([-100, 39, 59, 55, 20]), 5)
0.28095
>>> round(irr([-100, 0, 0, 74]), 5)
-0.0955
>>> round(irr([-100, 100, 0, -7]), 5)
-0.0833
>>> round(irr([-100, 100, 0, 7]), 5)
0.06206
>>> round(irr([-5, 10.5, 1, -8, 1]), 5)
0.0886
(Compare with the Example given for numpy.lib.financial.npv)
"""
res = np.roots(values[::-1])
mask = (res.imag == 0) & (res.real > 0)
if res.size == 0:
return np.nan
res = res[mask].real
# NPV(rate) = 0 can have more than one solution so we return
# only the solution closest to zero.
rate = 1.0/res - 1
rate = rate.item(np.argmin(np.abs(rate)))
return rate
0
Example 11
def test_roots(self):
assert_array_equal(np.roots([1, 0, 0]), [0, 0])
0
Example 12
Project: numba Source File: polynomial.py
@overload(np.roots)
def roots_impl(p):
# cast int vectors to float cf. numpy, this is a bit dicey as
# the roots could be complex which will fail anyway
ty = getattr(p, 'dtype', p)
if isinstance(ty, types.Integer):
cast_t = np.float64
else:
cast_t = np_support.as_dtype(ty)
def roots_impl(p):
# impl based on numpy:
# https://github.com/numpy/numpy/blob/master/numpy/lib/polynomial.py
if len(p.shape) != 1:
raise ValueError("Input must be a 1d array.")
non_zero = np.nonzero(p)[0]
if len(non_zero) == 0:
return np.zeros(0, dtype=cast_t)
tz = len(p) - non_zero[-1] - 1
# pull out the coeffs selecting between possible zero pads
p = p[int(non_zero[0]):int(non_zero[-1]) + 1]
n = len(p)
if n > 1:
# construct companion matrix, ensure fortran order
# to give to eigvals, write to upper diag and then
# transpose.
A = np.diag(np.ones((n - 2,), cast_t), 1).T
A[0, :] = -p[1:] / p[0] # normalize
roots = np.linalg.eigvals(A)
else:
roots = np.zeros(0, dtype=cast_t)
# add in additional zeros on the end if needed
if tz > 0:
return np.hstack((roots, np.zeros(tz, dtype=cast_t)))
else:
return roots
return roots_impl
0
Example 13
Project: numba Source File: test_polynomial.py
def roots_fn(p):
return np.roots(p)