python类conjugate()的实例源码

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?

def adj(self,x):
        """Returns Hermitian adjoint of a matrix"""
        assert x.shape[0] == x.shape[1]
        return np.conjugate(x).T

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())

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)

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)

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))

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))))

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))

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)

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

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

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

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

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

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

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

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

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

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

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