The following are code examples for showing how to use . They are extracted from open source Python projects. You can vote up the examples you like or vote down the exmaples you don’t like. You can also save this page to your account.
Example 1
def correlate(self, imgfft):
#Very much related to the convolution theorem, the cross-correlation
#theorem states that the Fourier transform of the cross-correlation of
#two functions is equal to the product of the individual Fourier
#transforms, where one of them has been complex conjugated:
if self.imgfft is not 0 or imgfft.imgfft is not 0:
imgcj = np.conjugate(self.imgfft)
imgft = imgfft.imgfft
prod = deepcopy(imgcj)
for x in range(imgcj.shape[0]):
for y in range(imgcj.shape[0]):
prod[x][y] = imgcj[x][y] * imgft[x][y]
cc = Corr( np.real(fft.ifft2(fft.fftshift(prod)))) # real image of the correlation
# adjust to center
cc.data = np.roll(cc.data, int(cc.data.shape[0] / 2), axis = 0)
cc.data = np.roll(cc.data, int(cc.data.shape[1] / 2), axis = 1)
else:
raise FFTnotInit()
return cc Example 2
def genSpectra(time,dipole,signal):
fw, frequency = pade(time,dipole)
fw_sig, frequency = pade(time,signal,alternate=True)
fw_re = np.real(fw)
fw_im = np.imag(fw)
fw_abs = fw_re**2 + fw_im**2
#spectra = (fw_re*17.32)/(np.pi*field*damp_const)
#spectra = (fw_re*17.32*514.220652)/(np.pi*field*damp_const)
#numerator = np.imag((fw*np.conjugate(fw_sig)))
numerator = np.imag(fw_abs*np.conjugate(fw_sig))
#numerator = np.abs((fw*np.conjugate(fw_sig)))
#numerator = np.abs(fw)
denominator = np.real(np.conjugate(fw_sig)*fw_sig)
#denominator = 1.0
spectra = ((4.0*27.21138602*2*frequency*np.pi*(numerator))/(3.0*137.036*denominator))
spectra *= 1.0/100.0
#plt.plot(frequency*27.2114,fourier)
#plt.show()
return frequency, spectra Example 3
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 4
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 5
def test_spherical_harmonic(scheme):
'''Assert the norm of the spherical harmonic
Y_1^1(phi, theta) = -1/2 sqrt(3/2/pi) * exp(i*phi) * sin(theta)
is indeed 1, i.e.,
int_0^2pi int_0^pi
Y_1^1(phi, theta) * conj(Y_1^1(phi, theta)) * sin(theta)
dphi dtheta = 1.
'''
def spherical_harmonic_11(azimuthal, polar):
# y00 = 1.0 / numpy.sqrt(4*numpy.pi)
y11 = -0.5 * numpy.sqrt(3.0/2.0/numpy.pi) \
* numpy.exp(1j*azimuthal) * numpy.sin(polar)
return y11 * numpy.conjugate(y11)
val = quadpy.sphere.integrate_spherical(
spherical_harmonic_11,
rule=scheme
)
assert abs(val - 1.0) < 1.0e-15
return Example 6
def __init__(self, n):
self.degree = n
self.points = -numpy.cos(numpy.pi * (numpy.arange(n) + 0.5) / n)
# n -= 1
N = numpy.arange(1, n, 2)
length = len(N)
m = n - length
K = numpy.arange(m)
v0 = numpy.concatenate([
2 * numpy.exp(1j*numpy.pi*K/n) / (1 - 4*K**2),
numpy.zeros(length+1)
])
v1 = v0[:-1] + numpy.conjugate(v0[:0:-1])
w = numpy.fft.ifft(v1)
assert max(w.imag) < 1.0e-15
self.weights = w.real
return Example 7
def CoreShellScatteringFunction(mCore,mShell,wavelength,dCore,dShell,minAngle=0, maxAngle=180, angularResolution=0.5, normed=False):
# http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellScatteringFunction
xCore = np.pi*dCore/wavelength
xShell = np.pi*dShell/wavelength
theta = np.linspace(minAngle,maxAngle,int((maxAngle-minAngle)/angularResolution))*np.pi/180
thetaSteps = len(theta)
SL = np.zeros(thetaSteps)
SR = np.zeros(thetaSteps)
SU = np.zeros(thetaSteps)
for j in range(thetaSteps):
u = np.cos(theta[j])
S1,S2 = CoreShellS1S2(mCore,mShell,xCore,xShell,u)
SL[j] = (np.sum((np.conjugate(S1)*S1))).real
SR[j] = (np.sum((np.conjugate(S2)*S2))).real
SU[j] = (SR[j]+SL[j])/2
if normed:
SL /= np.max(SL)
SR /= np.max(SR)
SU /= np.max(SU)
return theta,SL,SR,SU Example 8
def CoreShellScatteringFunction(mCore,mShell,wavelength,dCore,dShell,minAngle=0, maxAngle=180, angularResolution=0.5, normed=False):
# http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellScatteringFunction
xCore = np.pi*dCore/wavelength
xShell = np.pi*dShell/wavelength
theta = np.linspace(minAngle,maxAngle,int((maxAngle-minAngle)/angularResolution))*np.pi/180
thetaSteps = len(theta)
SL = np.zeros(thetaSteps)
SR = np.zeros(thetaSteps)
SU = np.zeros(thetaSteps)
for j in range(thetaSteps):
u = np.cos(theta[j])
S1,S2 = CoreShellS1S2(mCore,mShell,xCore,xShell,u)
SL[j] = (np.sum((np.conjugate(S1)*S1))).real
SR[j] = (np.sum((np.conjugate(S2)*S2))).real
SU[j] = (SR[j]+SL[j])/2
if normed:
SL /= np.max(SL)
SR /= np.max(SR)
SU /= np.max(SU)
return theta,SL,SR,SU Example 9
def gain_substitution_scalar(gain, x, xwt):
nants, nchan, nrec, _ = gain.shape
newgain = numpy.ones_like(gain, dtype='complex')
gwt = numpy.zeros_like(gain, dtype='float')
# We are going to work with Jones 2x2 matrix formalism so everything has to be
# converted to that format
x = x.reshape(nants, nants, nchan, nrec, nrec)
xwt = xwt.reshape(nants, nants, nchan, nrec, nrec)
for ant1 in range(nants):
for chan in range(nchan):
# Loop over e.g. 'RR', 'LL, or 'xx', 'YY' ignoring cross terms
top = numpy.sum(x[:, ant1, chan, 0, 0] *
gain[:, chan, 0, 0] * xwt[:, ant1, chan, 0, 0], axis=0)
bot = numpy.sum((gain[:, chan, 0, 0] * numpy.conjugate(gain[:, chan, 0, 0]) *
xwt[:, ant1, chan, 0, 0]).real, axis=0)
if bot > 0.0:
newgain[ant1, chan, 0, 0] = top / bot
gwt[ant1, chan, 0, 0] = bot
else:
newgain[ant1, chan, 0, 0] = 0.0
gwt[ant1, chan, 0, 0] = 0.0
return newgain, gwt Example 10
def convolve_scalestack(scalestack, img):
"""Convolve img by the specified scalestack, returning the resulting stack
:param scalestack: stack containing the scales
:param img: Image to be convolved
:return: stack
"""
convolved = numpy.zeros(scalestack.shape)
ximg = numpy.fft.fftshift(numpy.fft.fft2(numpy.fft.fftshift(img)))
nscales = scalestack.shape[0]
for iscale in range(nscales):
xscale = numpy.fft.fftshift(numpy.fft.fft2(numpy.fft.fftshift(scalestack[iscale, :, :])))
xmult = ximg * numpy.conjugate(xscale)
convolved[iscale, :, :] = numpy.real(numpy.fft.ifftshift(numpy.fft.ifft2(numpy.fft.ifftshift(xmult))))
return convolved Example 11
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 12
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 13
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 14
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 15
def cumsum_snr(freqs, detector, data, hpf, hxf, theta, phi, psi, zeroFreq=False, normalizeTemplate=True ):
"""
returns the cumulative sum of the snr as a function of frequency
does NOT maximize over the phase at coalescence
"""
template = detector.project(freqs, hpf, hxf, theta, phi, psi, zeroFreq=zeroFreq)
PSD = detector.PSD(freqs)
deltaF = freqs[1]-freqs[0]
ans = 2*np.cumsum(deltaF*np.conjugate(data)*template/PSD).real
if normalizeTemplate:
ans /= np.sum(deltaF*np.conjugate(template)*template/PSD).real**0.5
return ans
#------------------------ Example 16
def __newlambdagammaC(self, theta, l):
""" Apply Eqns. 25-27 in Vidal 2003 to update lambda^C and gamma^C
(lambda and gamma of this qbit).
"""
gamma_ket = self.coefs[l+1].lam
gamma_bra = np.conjugate(gamma_ket)
Gamma_star = np.conjugate(self.coefs[l+1].gamma)
inputs = [Gamma_star, theta, gamma_bra, gamma_ket]
Gamma_star_idx = [1, -3, -2]
theta_idx = [-1, 1, -4, -5]
gamma_bra_idx = [-6]
gamma_ket_idx = [-7]
idx = [Gamma_star_idx, theta_idx, gamma_bra_idx, gamma_ket_idx]
contract_me = scon(inputs, idx)
svd_me = np.einsum('agibggg', contract_me)
evals, evecs = la.eigh(svd_me)
return evals, evecs Example 17
def paulidouble(i, j, tensor=True):
pauli_i = pauli(i)
pauli_j = pauli(j)
outer = np.zeros((4, 4), dtype=np.complex)
outer[:2, :2] = pauli_i
outer[2:, 2:] = pauli_j
if tensor:
outer.shape = (2, 2, 2, 2)
return outer
# def paulitwo_left(i):
# return np.kron(pauli(i), pauli(0))
# def paulitwo_right(i):
# return np.kron(pauli(0), pauli(i))
# def newrho_DK(Lket, theta_ij):
# Lbra = np.conjugate(L_before)
# theta_star = np.conjugate(theta_ij)
# in_bracket = scon(Lbra, Lket, theta_ij, theta_star,
# [1], [1], [2, 3, 1, Example 18
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 19
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 20
def qisunitary(self,op): mat = op[1] (r,c) = mat.shape if r != c: return False invmat = np.conjugate(np.transpose(mat)) pmat = np.asarray(mat * invmat) for i in range(r): for j in range(c): if i != j: if np.absolute(pmat[i][j]) > self.maxerr: return False else: if np.absolute(pmat[i][j]-1.0) > self.maxerr: return False return True
Example 21
def __init__(self, gain):
self.gain = gain * np.pi # pi term puts angle output to pi range
# components
self.conjugate = Conjugate()
self.complex_mult = ComplexMultiply()
self.angle = Angle()
self.out = Sfix()
# specify component delay
self._delay = self.conjugate._delay + \
self.complex_mult._delay + \
self.angle._delay + 1
# constants
self.gain_sfix = Const(Sfix(self.gain, 3, -14)) Example 22
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 23
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 24
def test_poly(self):
assert_array_almost_equal(np.poly([3, -np.sqrt(2), np.sqrt(2)]),
[1, -3, -2, 6])
# From matlab docs
A = [[1, 2, 3], [4, 5, 6], [7, 8, 0]]
assert_array_almost_equal(np.poly(A), [1, -6, -72, -27])
# Should produce real output for perfect conjugates
assert_(np.isrealobj(np.poly([+1.082j, +2.613j, -2.613j, -1.082j])))
assert_(np.isrealobj(np.poly([0+1j, -0+-1j, 1+2j,
1-2j, 1.+3.5j, 1-3.5j])))
assert_(np.isrealobj(np.poly([1j, -1j, 1+2j, 1-2j, 1+3j, 1-3.j])))
assert_(np.isrealobj(np.poly([1j, -1j, 1+2j, 1-2j])))
assert_(np.isrealobj(np.poly([1j, -1j, 2j, -2j])))
assert_(np.isrealobj(np.poly([1j, -1j])))
assert_(np.isrealobj(np.poly([1, -1])))
assert_(np.iscomplexobj(np.poly([1j, -1.0000001j])))
np.random.seed(42)
a = np.random.randn(100) + 1j*np.random.randn(100)
assert_(np.isrealobj(np.poly(np.concatenate((a, np.conjugate(a)))))) Example 25
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 26
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 27
def get_pmf_xi(self):
"""Return the values of the variable ``xi``.
The components ``xi`` make up the probability mass function, i.e.
:math:`\\xi(k) = pmf(k) = Pr(X = k)`.
"""
chi = np.empty(self.number_trials + 1, dtype=complex)
chi[0] = 1
half_number_trials = int(
self.number_trials / 2 + self.number_trials % 2)
# set first half of chis:
chi[1:half_number_trials + 1] = self.get_chi(
np.arange(1, half_number_trials + 1))
# set second half of chis:
chi[half_number_trials + 1:self.number_trials + 1] = np.conjugate(
chi[1:self.number_trials - half_number_trials + 1] [::-1])
chi /= self.number_trials + 1
xi = np.fft.fft(chi)
if self.check_xi_are_real(xi):
xi = xi.real
else:
raise TypeError("pmf / xi values have to be real.")
xi += np.finfo(type(xi[0])).eps
return xi Example 28
def test_generic_methods(self):
# Tests some MaskedArray methods.
a = array([1, 3, 2])
assert_equal(a.any(), a._data.any())
assert_equal(a.all(), a._data.all())
assert_equal(a.argmax(), a._data.argmax())
assert_equal(a.argmin(), a._data.argmin())
assert_equal(a.choose(0, 1, 2, 3, 4), a._data.choose(0, 1, 2, 3, 4))
assert_equal(a.compress([1, 0, 1]), a._data.compress([1, 0, 1]))
assert_equal(a.conj(), a._data.conj())
assert_equal(a.conjugate(), a._data.conjugate())
m = array([[1, 2], [3, 4]])
assert_equal(m.diagonal(), m._data.diagonal())
assert_equal(a.sum(), a._data.sum())
assert_equal(a.take([1, 2]), a._data.take([1, 2]))
assert_equal(m.transpose(), m._data.transpose()) Example 29
def test_testArrayMethods(self):
a = array([1, 3, 2])
self.assertTrue(eq(a.any(), a._data.any()))
self.assertTrue(eq(a.all(), a._data.all()))
self.assertTrue(eq(a.argmax(), a._data.argmax()))
self.assertTrue(eq(a.argmin(), a._data.argmin()))
self.assertTrue(eq(a.choose(0, 1, 2, 3, 4),
a._data.choose(0, 1, 2, 3, 4)))
self.assertTrue(eq(a.compress([1, 0, 1]), a._data.compress([1, 0, 1])))
self.assertTrue(eq(a.conj(), a._data.conj()))
self.assertTrue(eq(a.conjugate(), a._data.conjugate()))
m = array([[1, 2], [3, 4]])
self.assertTrue(eq(m.diagonal(), m._data.diagonal()))
self.assertTrue(eq(a.sum(), a._data.sum()))
self.assertTrue(eq(a.take([1, 2]), a._data.take([1, 2])))
self.assertTrue(eq(m.transpose(), m._data.transpose())) Example 30
def _solve_lens_equation(self, complex_value):
"""
Solve the lens equation for the given point (in complex coordinates).
"""
complex_conjugate = np.conjugate(complex_value)
return complex_value - (1. / (1. + self.q)) * (
(1./complex_conjugate) + (self.q / (complex_conjugate - self.s))) Example 31
def _jacobian_determinant_ok_WM95(self, source_x, source_y):
"""determinants of lens equation Jacobian for verified roots"""
roots_ok_bar = np.conjugate(self._polynomial_roots_ok_WM95(
source_x=source_x, source_y=source_y))
# Variable X_bar is conjugate of variable X.
denominator_1 = self._position_z1_WM95 - roots_ok_bar
add_1 = self.mass_1 / denominator_1**2
denominator_2 = self._position_z2_WM95 - roots_ok_bar
add_2 = self.mass_2 / denominator_2**2
derivative = add_1 + add_2
# Can we use Utils.complex_fsum()-like function here?
return 1.-derivative*np.conjugate(derivative)
# Can we use Utils.complex_fsum()-like function here? Example 32
def adj(self,x):
"""Returns Hermitian adjoint of a matrix"""
assert x.shape[0] == x.shape[1]
return np.conjugate(x).T Example 33
def test_conjugate(self):
a = np.array([1-1j, 1+1j, 23+23.0j])
ac = a.conj()
assert_equal(a.real, ac.real)
assert_equal(a.imag, -ac.imag)
assert_equal(ac, a.conjugate())
assert_equal(ac, np.conjugate(a))
a = np.array([1-1j, 1+1j, 23+23.0j], 'F')
ac = a.conj()
assert_equal(a.real, ac.real)
assert_equal(a.imag, -ac.imag)
assert_equal(ac, a.conjugate())
assert_equal(ac, np.conjugate(a))
a = np.array([1, 2, 3])
ac = a.conj()
assert_equal(a, ac)
assert_equal(ac, a.conjugate())
assert_equal(ac, np.conjugate(a))
a = np.array([1.0, 2.0, 3.0])
ac = a.conj()
assert_equal(a, ac)
assert_equal(ac, a.conjugate())
assert_equal(ac, np.conjugate(a))
a = np.array([1-1j, 1+1j, 1, 2.0], object)
ac = a.conj()
assert_equal(ac, [k.conjugate() for k in a])
assert_equal(ac, a.conjugate())
assert_equal(ac, np.conjugate(a))
a = np.array([1-1j, 1, 2.0, 'f'], object)
assert_raises(AttributeError, lambda: a.conj())
assert_raises(AttributeError, lambda: a.conjugate()) Example 34
def test_var_values(self):
for mat in [self.rmat, self.cmat, self.omat]:
for axis in [0, 1, None]:
msqr = _mean(mat * mat.conj(), axis=axis)
mean = _mean(mat, axis=axis)
tgt = msqr - mean * mean.conjugate()
res = _var(mat, axis=axis)
assert_almost_equal(res, tgt) Example 35
def test_oddfeatures_1(self):
# Test of other odd features
x = arange(20)
x = x.reshape(4, 5)
x.flat[5] = 12
assert_(x[1, 0] == 12)
z = x + 10j * x
assert_equal(z.real, x)
assert_equal(z.imag, 10 * x)
assert_equal((z * conjugate(z)).real, 101 * x * x)
z.imag[...] = 0.0
x = arange(10)
x[3] = masked
assert_(str(x[3]) == str(masked))
c = x >= 8
assert_(count(where(c, masked, masked)) == 0)
assert_(shape(where(c, masked, masked)) == c.shape)
z = masked_where(c, x)
assert_(z.dtype is x.dtype)
assert_(z[3] is masked)
assert_(z[4] is not masked)
assert_(z[7] is not masked)
assert_(z[8] is masked)
assert_(z[9] is masked)
assert_equal(x, z) Example 36
def test_basic_ufuncs(self):
# Test various functions such as sin, cos.
(x, y, a10, m1, m2, xm, ym, z, zm, xf) = self.d
assert_equal(np.cos(x), cos(xm))
assert_equal(np.cosh(x), cosh(xm))
assert_equal(np.sin(x), sin(xm))
assert_equal(np.sinh(x), sinh(xm))
assert_equal(np.tan(x), tan(xm))
assert_equal(np.tanh(x), tanh(xm))
assert_equal(np.sqrt(abs(x)), sqrt(xm))
assert_equal(np.log(abs(x)), log(xm))
assert_equal(np.log10(abs(x)), log10(xm))
assert_equal(np.exp(x), exp(xm))
assert_equal(np.arcsin(z), arcsin(zm))
assert_equal(np.arccos(z), arccos(zm))
assert_equal(np.arctan(z), arctan(zm))
assert_equal(np.arctan2(x, y), arctan2(xm, ym))
assert_equal(np.absolute(x), absolute(xm))
assert_equal(np.angle(x + 1j*y), angle(xm + 1j*ym))
assert_equal(np.angle(x + 1j*y, deg=True), angle(xm + 1j*ym, deg=True))
assert_equal(np.equal(x, y), equal(xm, ym))
assert_equal(np.not_equal(x, y), not_equal(xm, ym))
assert_equal(np.less(x, y), less(xm, ym))
assert_equal(np.greater(x, y), greater(xm, ym))
assert_equal(np.less_equal(x, y), less_equal(xm, ym))
assert_equal(np.greater_equal(x, y), greater_equal(xm, ym))
assert_equal(np.conjugate(x), conjugate(xm)) Example 37
def test_testUfuncs1(self):
# Test various functions such as sin, cos.
(x, y, a10, m1, m2, xm, ym, z, zm, xf, s) = self.d
self.assertTrue(eq(np.cos(x), cos(xm)))
self.assertTrue(eq(np.cosh(x), cosh(xm)))
self.assertTrue(eq(np.sin(x), sin(xm)))
self.assertTrue(eq(np.sinh(x), sinh(xm)))
self.assertTrue(eq(np.tan(x), tan(xm)))
self.assertTrue(eq(np.tanh(x), tanh(xm)))
with np.errstate(divide='ignore', invalid='ignore'):
self.assertTrue(eq(np.sqrt(abs(x)), sqrt(xm)))
self.assertTrue(eq(np.log(abs(x)), log(xm)))
self.assertTrue(eq(np.log10(abs(x)), log10(xm)))
self.assertTrue(eq(np.exp(x), exp(xm)))
self.assertTrue(eq(np.arcsin(z), arcsin(zm)))
self.assertTrue(eq(np.arccos(z), arccos(zm)))
self.assertTrue(eq(np.arctan(z), arctan(zm)))
self.assertTrue(eq(np.arctan2(x, y), arctan2(xm, ym)))
self.assertTrue(eq(np.absolute(x), absolute(xm)))
self.assertTrue(eq(np.equal(x, y), equal(xm, ym)))
self.assertTrue(eq(np.not_equal(x, y), not_equal(xm, ym)))
self.assertTrue(eq(np.less(x, y), less(xm, ym)))
self.assertTrue(eq(np.greater(x, y), greater(xm, ym)))
self.assertTrue(eq(np.less_equal(x, y), less_equal(xm, ym)))
self.assertTrue(eq(np.greater_equal(x, y), greater_equal(xm, ym)))
self.assertTrue(eq(np.conjugate(x), conjugate(xm)))
self.assertTrue(eq(np.concatenate((x, y)), concatenate((xm, ym))))
self.assertTrue(eq(np.concatenate((x, y)), concatenate((x, y))))
self.assertTrue(eq(np.concatenate((x, y)), concatenate((xm, y))))
self.assertTrue(eq(np.concatenate((x, y, x)), concatenate((x, ym, x)))) Example 38
def test_testUfuncRegression(self):
f_invalid_ignore = [
'sqrt', 'arctanh', 'arcsin', 'arccos',
'arccosh', 'arctanh', 'log', 'log10', 'divide',
'true_divide', 'floor_divide', 'remainder', 'fmod']
for f in ['sqrt', 'log', 'log10', 'exp', 'conjugate',
'sin', 'cos', 'tan',
'arcsin', 'arccos', 'arctan',
'sinh', 'cosh', 'tanh',
'arcsinh',
'arccosh',
'arctanh',
'absolute', 'fabs', 'negative',
'floor', 'ceil',
'logical_not',
'add', 'subtract', 'multiply',
'divide', 'true_divide', 'floor_divide',
'remainder', 'fmod', 'hypot', 'arctan2',
'equal', 'not_equal', 'less_equal', 'greater_equal',
'less', 'greater',
'logical_and', 'logical_or', 'logical_xor']:
try:
uf = getattr(umath, f)
except AttributeError:
uf = getattr(fromnumeric, f)
mf = getattr(np.ma, f)
args = self.d[:uf.nin]
with np.errstate():
if f in f_invalid_ignore:
np.seterr(invalid='ignore')
if f in ['arctanh', 'log', 'log10']:
np.seterr(divide='ignore')
ur = uf(*args)
mr = mf(*args)
self.assertTrue(eq(ur.filled(0), mr.filled(0), f))
self.assertTrue(eqmask(ur.mask, mr.mask)) Example 39
def power_process(data, sfreq, toffset, modulus, integration, log_scale, zscale, title):
""" Break power by modulus and display each block. Integration here acts
as a pure average on the power level data.
"""
if modulus:
block = 0
block_size = integration * modulus
block_toffset = toffset
while block < len(data) / block_size:
vblock = data[block * block_size:block * block_size + modulus]
pblock = (vblock * numpy.conjugate(vblock)).real
# complete integration
for idx in range(1, integration):
vblock = data[block * block_size + idx *
modulus:block * block_size + idx * modulus + modulus]
pblock += (vblock * numpy.conjugate(vblock)).real
pblock /= integration
yield power_plot(pblock, sfreq, block_toffset, log_scale, zscale, title)
block += 1
block_toffset += block_size / sfreq
else:
pdata = (data * numpy.conjugate(data)).real
yield power_plot(pdata, sfreq, toffset, log_scale, zscale, title) Example 40
def MatrixElements(m,wavelength,diameter,mu): # http://pymiescatt.readthedocs.io/en/latest/forward.html#MatrixElements x = np.pi*diameter/wavelength # B&H eqs. 4.77, where mu=cos(theta) S1,S2 = MieS1S2(m,x,mu) S11 = 0.5*(np.abs(S2)**2+np.abs(S1)**2) S12 = 0.5*(np.abs(S2)**2-np.abs(S1)**2) S33 = 0.5*(np.conjugate(S2)*S1+S2*np.conjugate(S1)) S34 = 0.5j*(S1*np.conjugate(S2)-S2*np.conjugate(S1)) return S11, S12, S33, S34
Example 41
def CoreShellS1S2(mCore,mShell,xCore,xShell,mu): # http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellS1S2 nmax = np.round(2+xShell+4*(xShell**(1/3))) an,bn = CoreShell_ab(mCore,mShell,xCore,xShell) pin,taun = MiePiTau(mu,nmax) n = np.arange(1,int(nmax)+1) n2 = (2*n+1)/(n*(n+1)) pin *= n2 taun *= n2 S1=np.sum(an*np.conjugate(pin))+np.sum(bn*np.conjugate(taun)) S2=np.sum(an*np.conjugate(taun))+np.sum(bn*np.conjugate(pin)) return S1,S2
Example 42
def CoreShellMatrixElements(mCore,mShell,xCore,xShell,mu): # http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellMatrixElements S1,S2 = CoreShellS1S2(mCore,mShell,xCore,xShell,mu) S11 = 0.5*(np.abs(S2)**2+np.abs(S1)**2) S12 = 0.5*(np.abs(S2)**2-np.abs(S1)**2) S33 = 0.5*(np.conjugate(S2)*S1+S2*np.conjugate(S1)) S34 = 0.5j*(S1*np.conjugate(S2)-S2*np.conjugate(S1)) return S11, S12, S33, S34
Example 43
def MatrixElements(m,wavelength,diameter,mu): # http://pymiescatt.readthedocs.io/en/latest/forward.html#MatrixElements x = np.pi*diameter/wavelength # B&H eqs. 4.77, where mu=cos(theta) S1,S2 = MieS1S2(m,x,mu) S11 = 0.5*(np.abs(S2)**2+np.abs(S1)**2) S12 = 0.5*(np.abs(S2)**2-np.abs(S1)**2) S33 = 0.5*(np.conjugate(S2)*S1+S2*np.conjugate(S1)) S34 = 0.5j*(S1*np.conjugate(S2)-S2*np.conjugate(S1)) return S11, S12, S33, S34
Example 44
def CoreShellS1S2(mCore,mShell,xCore,xShell,mu): # http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellS1S2 nmax = np.round(2+xShell+4*(xShell**(1/3))) an,bn = CoreShell_ab(mCore,mShell,xCore,xShell) pin,taun = MiePiTau(mu,nmax) n = np.arange(1,int(nmax)+1) n2 = (2*n+1)/(n*(n+1)) pin *= n2 taun *= n2 S1=np.sum(an*np.conjugate(pin))+np.sum(bn*np.conjugate(taun)) S2=np.sum(an*np.conjugate(taun))+np.sum(bn*np.conjugate(pin)) return S1,S2
Example 45
def CoreShellMatrixElements(mCore,mShell,xCore,xShell,mu): # http://pymiescatt.readthedocs.io/en/latest/forwardCS.html#CoreShellMatrixElements S1,S2 = CoreShellS1S2(mCore,mShell,xCore,xShell,mu) S11 = 0.5*(np.abs(S2)**2+np.abs(S1)**2) S12 = 0.5*(np.abs(S2)**2-np.abs(S1)**2) S33 = 0.5*(np.conjugate(S2)*S1+S2*np.conjugate(S1)) S34 = 0.5j*(S1*np.conjugate(S2)-S2*np.conjugate(S1)) return S11, S12, S33, S34
Example 46
def sum_visibility(vis: Visibility, direction: SkyCoord) -> numpy.array:
""" Direct Fourier summation in a given direction
:param vis: Visibility to be summed
:param direction: Direction of summation
:return: flux[nch,npol], weight[nch,pol]
"""
# TODO: Convert to Visibility or remove?
assert isinstance(vis, Visibility) or isinstance(vis, BlockVisibility), vis
svis = copy_visibility(vis)
l, m, n = skycoord_to_lmn(direction, svis.phasecentre)
phasor = numpy.conjugate(simulate_point(svis.uvw, l, m))
# Need to put correct mapping here
_, frequency = get_frequency_map(svis, None)
frequency = list(frequency)
nchan = max(frequency) + 1
npol = svis.polarisation_frame.npol
flux = numpy.zeros([nchan, npol])
weight = numpy.zeros([nchan, npol])
coords = svis.vis, svis.weight, phasor, list(frequency)
for v, wt, p, ic in zip(*coords):
for pol in range(npol):
flux[ic, pol] += numpy.real(wt[pol] * v[pol] * p)
weight[ic, pol] += wt[pol]
flux[weight > 0.0] = flux[weight > 0.0] / weight[weight > 0.0]
flux[weight <= 0.0] = 0.0
return flux, weight Example 47
def gain_substitution_vector(gain, x, xwt):
nants, nchan, nrec, _ = gain.shape
newgain = numpy.ones_like(gain, dtype='complex')
if nrec > 0:
newgain[..., 0, 1] = 0.0
newgain[..., 1, 0] = 0.0
gwt = numpy.zeros_like(gain, dtype='float')
# We are going to work with Jones 2x2 matrix formalism so everything has to be
# converted to that format
x = x.reshape(nants, nants, nchan, nrec, nrec)
xwt = xwt.reshape(nants, nants, nchan, nrec, nrec)
if nrec > 0:
gain[..., 0, 1] = 0.0
gain[..., 1, 0] = 0.0
for ant1 in range(nants):
for chan in range(nchan):
# Loop over e.g. 'RR', 'LL, or 'xx', 'YY' ignoring cross terms
for rec in range(nrec):
top = numpy.sum(x[:, ant1, chan, rec, rec] *
gain[:, chan, rec, rec] * xwt[:, ant1, chan, rec, rec], axis=0)
bot = numpy.sum((gain[:, chan, rec, rec] * numpy.conjugate(gain[:, chan, rec, rec]) *
xwt[:, ant1, chan, rec, rec]).real, axis=0)
if bot > 0.0:
newgain[ant1, chan, rec, rec] = top / bot
gwt[ant1, chan, rec, rec] = bot
else:
newgain[ant1, chan, rec, rec] = 0.0
gwt[ant1, chan, rec, rec] = 0.0
return newgain, gwt Example 48
def gain_substitution_matrix(gain, x, xwt):
nants, nchan, nrec, _ = gain.shape
newgain = numpy.ones_like(gain, dtype='complex')
gwt = numpy.zeros_like(gain, dtype='float')
# We are going to work with Jones 2x2 matrix formalism so everything has to be
# converted to that format
x = x.reshape(nants, nants, nchan, nrec, nrec)
xwt = xwt.reshape(nants, nants, nchan, nrec, nrec)
# Write these loops out explicitly. Derivation of these vector equations is tedious but they are
# structurally identical to the scalar case with the following changes
# Vis -> 2x2 coherency vector, g-> 2x2 Jones matrix, *-> matmul, conjugate->Hermitean transpose (.H)
for ant1 in range(nants):
for chan in range(nchan):
top = 0.0
bot = 0.0
for ant2 in range(nants):
if ant1 != ant2:
xmat = x[ant2, ant1, chan]
xwtmat = xwt[ant2, ant1, chan]
g2 = gain[ant2, chan]
top += xmat * xwtmat * g2
bot += numpy.conjugate(g2) * xwtmat * g2
newgain[ant1, chan][bot > 0.0] = top[bot > 0.0] / bot[bot > 0.0]
newgain[ant1, chan][bot <= 0.0] = 0.0
gwt[ant1, chan] = bot.real
return newgain, gwt Example 49
def solution_residual_scalar(gain, x, xwt):
"""Calculate residual across all baselines of gain for point source equivalent visibilities
:param gain: gain [nant, ...]
:param x: Point source equivalent visibility [nant, ...]
:param xwt: Point source equivalent weight [nant, ...]
:return: residual[...]
"""
nants, nchan, nrec, _ = gain.shape
x = x.reshape(nants, nants, nchan, nrec, nrec)
xwt = xwt.reshape(nants, nants, nchan, nrec, nrec)
residual = numpy.zeros([nchan, nrec, nrec])
sumwt = numpy.zeros([nchan, nrec, nrec])
for ant1 in range(nants):
for ant2 in range(nants):
for chan in range(nchan):
error = x[ant2, ant1, chan, 0, 0] - \
gain[ant1, chan, 0, 0] * numpy.conjugate(gain[ant2, chan, 0, 0])
residual += (error * xwt[ant2, ant1, chan, 0, 0] * numpy.conjugate(error)).real
sumwt += xwt[ant2, ant1, chan, 0, 0]
residual[sumwt > 0.0] = numpy.sqrt(residual[sumwt > 0.0] / sumwt[sumwt > 0.0])
residual[sumwt <= 0.0] = 0.0
return residual Example 50
def solution_residual_vector(gain, x, xwt):
"""Calculate residual across all baselines of gain for point source equivalent visibilities
Vector case i.e. off-diagonals of gains are zero
:param gain: gain [nant, ...]
:param x: Point source equivalent visibility [nant, ...]
:param xwt: Point source equivalent weight [nant, ...]
:return: residual[...]
"""
nants, nchan, nrec, _ = gain.shape
x = x.reshape(nants, nants, nchan, nrec, nrec)
x[..., 1, 0] = 0.0
x[..., 0, 1] = 0.0
xwt = xwt.reshape(nants, nants, nchan, nrec, nrec)
xwt[..., 1, 0] = 0.0
xwt[..., 0, 1] = 0.0
residual = numpy.zeros([nchan, nrec, nrec])
sumwt = numpy.zeros([nchan, nrec, nrec])
for ant1 in range(nants):
for ant2 in range(nants):
for chan in range(nchan):
for rec in range(nrec):
error = x[ant2, ant1, chan, rec, rec] - \
gain[ant1, chan, rec, rec] * numpy.conjugate(gain[ant2, chan, rec, rec])
residual += (error * xwt[ant2, ant1, chan, rec, rec] * numpy.conjugate(error)).real
sumwt += xwt[ant2, ant1, chan, rec, rec]
residual[sumwt > 0.0] = numpy.sqrt(residual[sumwt > 0.0] / sumwt[sumwt > 0.0])
residual[sumwt <= 0.0] = 0.0
return residual 