我们从Python开源项目中,提取了以下50个代码示例,用于说明如何使用numpy.complex64()。
def run_test_matmul_aa_correlator_kernel(self, ntime, nstand, nchan, misalign=0): x_shape = (ntime, nchan, nstand*2) perm = [1,0,2] x8 = ((np.random.random(size=x_shape+(2,))*2-1)*127).astype(np.int8) x = x8.astype(np.float32).view(np.complex64).reshape(x_shape) x = x.transpose(perm) x = x[..., misalign:] b_gold = np.matmul(H(x), x) triu = np.triu_indices(x.shape[-1], 1) b_gold[..., triu[0], triu[1]] = 0 x = x8.view(bf.DataType.ci8).reshape(x_shape) x = bf.asarray(x, space='cuda') x = x.transpose(perm) x = x[..., misalign:] b = bf.zeros_like(b_gold, space='cuda') self.linalg.matmul(1, None, x, 0, b) b = b.copy('system') np.testing.assert_allclose(b, b_gold, RTOL*10, ATOL)
def flip_code(code): if isinstance(code, (numpy.dtype,type)): # since several things map to complex64 we must carefully select # the opposite that is an exact match (ticket 1518) if code == numpy.int8: return gdalconst.GDT_Byte if code == numpy.complex64: return gdalconst.GDT_CFloat32 for key, value in codes.items(): if value == code: return key return None else: try: return codes[code] except KeyError: return None
def solver(self,gy, solver=None, *args, **kwargs): """ The solver of NUFFT :param gy: data, reikna array, (M,) size :param solver: could be 'cg', 'L1TVOLS', 'L1TVLAD' :param maxiter: the number of iterations :type gy: reikna array, dtype = numpy.complex64 :type solver: string :type maxiter: int :return: reikna array with size Nd """ import src._solver.solver_hsa try: return src._solver.solver_hsa.solver(self, gy, solver, *args, **kwargs) except: if numpy.ndarray==type(gy): print("input gy must be a reikna array with dtype = numpy.complex64") raise TypeError else: print("wrong") raise TypeError
def forward(self, gx): """ Forward NUFFT on the heterogeneous device :param gx: The input gpu array, with size=Nd :type: reikna gpu array with dtype =numpy.complex64 :return: gy: The output gpu array, with size=(M,) :rtype: reikna gpu array with dtype =numpy.complex64 """ self.x_Nd = self.thr.copy_array(gx) self._x2xx() self._xx2k() self._k2y() gy = self.thr.copy_array(self.y) return gy
def adjoint(self, gy): """ Adjoint NUFFT on the heterogeneous device :param gy: The output gpu array, with size=(M,) :type: reikna gpu array with dtype =numpy.complex64 :return: gx: The input gpu array, with size=Nd :rtype: reikna gpu array with dtype =numpy.complex64 """ self.y = self.thr.copy_array(gy) self._y2k() self._k2xx() self._xx2x() gx = self.thr.copy_array(self.x_Nd) return gx
def selfadjoint(self, gx): """ selfadjoint NUFFT (Teplitz) on the heterogeneous device :param gx: The input gpu array, with size=Nd :type: reikna gpu array with dtype =numpy.complex64 :return: gx: The input gpu array, with size=Nd :rtype: reikna gpu array with dtype =numpy.complex64 """ self.x_Nd = self.thr.copy_array(gx) self._x2xx() self._xx2k() self._k2y2k() self._k2xx() self._xx2x() gx2 = self.thr.copy_array(self.x_Nd) return gx2
def __init__(self): """ Constructor. :param None: :type None: Python NoneType :return: NUFFT: the pynufft_hsa.NUFFT instance :rtype: NUFFT: the pynufft_hsa.NUFFT class :Example: >>> import pynufft.pynufft >>> NufftObj = pynufft.pynufft.NUFFT_cpu() .. note:: requires plan() .. seealso:: :method:`plan()' .. todo:: test 3D case """ self.dtype=numpy.complex64 pass
def solve(self,gy, solver=None, *args, **kwargs): """ The solver of NUFFT_hsa :param gy: data, reikna array, (M,) size :param solver: could be 'cg', 'L1TVOLS', 'L1TVLAD' :param maxiter: the number of iterations :type gy: reikna array, dtype = numpy.complex64 :type solver: string :type maxiter: int :return: reikna array with size Nd """ from .._nonlin.solve_hsa import solve try: return solve(self, gy, solver, *args, **kwargs) except: if numpy.ndarray==type(gy): print("input gy must be a reikna array with dtype = numpy.complex64") raise TypeError else: print("wrong") raise TypeError
def selfadjoint(self, gx): """ selfadjoint NUFFT (Teplitz) on the heterogeneous device :param gx: The input gpu array, with size=Nd :type: reikna gpu array with dtype =numpy.complex64 :return: gx: The output gpu array, with size=Nd :rtype: reikna gpu array with dtype =numpy.complex64 """ self.x_Nd = self.thr.copy_array(gx) self._x2xx() self._xx2k() self._k2y2k() self._k2xx() self._xx2x() gx2 = self.thr.copy_array(self.x_Nd) return gx2
def guarantee_array(variable): ''' Guarantees that a varaible is a numpy ndarray and supports -, *, +, and other operators Args: variable (`number` or `numpy.ndarray`): variable to coalesce Returns: (type). Which supports * / and other operations with arrays ''' if type(variable) in [float, np.ndarray, np.int32, np.int64, np.float32, np.float64, np.complex64, np.complex128]: return variable elif type(variable) is int: return float(variable) elif type(variable) is list: return np.asarray(variable) else: raise ValueError(f'variable is of invalid type {type(variable)}')
def test_branch_cuts_complex64(self): # check branch cuts and continuity on them yield _check_branch_cut, np.log, -0.5, 1j, 1, -1, True, np.complex64 yield _check_branch_cut, np.log2, -0.5, 1j, 1, -1, True, np.complex64 yield _check_branch_cut, np.log10, -0.5, 1j, 1, -1, True, np.complex64 yield _check_branch_cut, np.log1p, -1.5, 1j, 1, -1, True, np.complex64 yield _check_branch_cut, np.sqrt, -0.5, 1j, 1, -1, True, np.complex64 yield _check_branch_cut, np.arcsin, [ -2, 2], [1j, 1j], 1, -1, True, np.complex64 yield _check_branch_cut, np.arccos, [ -2, 2], [1j, 1j], 1, -1, True, np.complex64 yield _check_branch_cut, np.arctan, [0-2j, 2j], [1, 1], -1, 1, True, np.complex64 yield _check_branch_cut, np.arcsinh, [0-2j, 2j], [1, 1], -1, 1, True, np.complex64 yield _check_branch_cut, np.arccosh, [ -1, 0.5], [1j, 1j], 1, -1, True, np.complex64 yield _check_branch_cut, np.arctanh, [ -2, 2], [1j, 1j], 1, -1, True, np.complex64 # check against bogus branch cuts: assert continuity between quadrants yield _check_branch_cut, np.arcsin, [0-2j, 2j], [ 1, 1], 1, 1, False, np.complex64 yield _check_branch_cut, np.arccos, [0-2j, 2j], [ 1, 1], 1, 1, False, np.complex64 yield _check_branch_cut, np.arctan, [ -2, 2], [1j, 1j], 1, 1, False, np.complex64 yield _check_branch_cut, np.arcsinh, [ -2, 2, 0], [1j, 1j, 1], 1, 1, False, np.complex64 yield _check_branch_cut, np.arccosh, [0-2j, 2j, 2], [1, 1, 1j], 1, 1, False, np.complex64 yield _check_branch_cut, np.arctanh, [0-2j, 2j, 0], [1, 1, 1j], 1, 1, False, np.complex64
def test_prod(self): ba = [1, 2, 10, 11, 6, 5, 4] ba2 = [[1, 2, 3, 4], [5, 6, 7, 9], [10, 3, 4, 5]] for ctype in [np.int16, np.uint16, np.int32, np.uint32, np.float32, np.float64, np.complex64, np.complex128]: a = np.array(ba, ctype) a2 = np.array(ba2, ctype) if ctype in ['1', 'b']: self.assertRaises(ArithmeticError, a.prod) self.assertRaises(ArithmeticError, a2.prod, axis=1) else: assert_equal(a.prod(axis=0), 26400) assert_array_equal(a2.prod(axis=0), np.array([50, 36, 84, 180], ctype)) assert_array_equal(a2.prod(axis=-1), np.array([24, 1890, 600], ctype))
def test_sum_complex(self): for dt in (np.complex64, np.complex128, np.clongdouble): for v in (0, 1, 2, 7, 8, 9, 15, 16, 19, 127, 128, 1024, 1235): tgt = dt(v * (v + 1) / 2) - dt((v * (v + 1) / 2) * 1j) d = np.empty(v, dtype=dt) d.real = np.arange(1, v + 1) d.imag = -np.arange(1, v + 1) assert_almost_equal(np.sum(d), tgt) assert_almost_equal(np.sum(d[::-1]), tgt) d = np.ones(500, dtype=dt) + 1j assert_almost_equal(np.sum(d[::2]), 250. + 250j) assert_almost_equal(np.sum(d[1::2]), 250. + 250j) assert_almost_equal(np.sum(d[::3]), 167. + 167j) assert_almost_equal(np.sum(d[1::3]), 167. + 167j) assert_almost_equal(np.sum(d[::-2]), 250. + 250j) assert_almost_equal(np.sum(d[-1::-2]), 250. + 250j) assert_almost_equal(np.sum(d[::-3]), 167. + 167j) assert_almost_equal(np.sum(d[-1::-3]), 167. + 167j) # sum with first reduction entry != 0 d = np.ones((1,), dtype=dt) + 1j d += d assert_almost_equal(d, 2. + 2j)
def test_zero_division(self): with np.errstate(all="ignore"): for t in [np.complex64, np.complex128]: a = t(0.0) b = t(1.0) assert_(np.isinf(b/a)) b = t(complex(np.inf, np.inf)) assert_(np.isinf(b/a)) b = t(complex(np.inf, np.nan)) assert_(np.isinf(b/a)) b = t(complex(np.nan, np.inf)) assert_(np.isinf(b/a)) b = t(complex(np.nan, np.nan)) assert_(np.isnan(b/a)) b = t(0.) assert_(np.isnan(b/a))
def test_basic(self): ba = [1, 2, 10, 11, 6, 5, 4] ba2 = [[1, 2, 3, 4], [5, 6, 7, 9], [10, 3, 4, 5]] for ctype in [np.int8, np.uint8, np.int16, np.uint16, np.int32, np.uint32, np.float32, np.float64, np.complex64, np.complex128]: a = np.array(ba, ctype) a2 = np.array(ba2, ctype) tgt = np.array([1, 3, 13, 24, 30, 35, 39], ctype) assert_array_equal(np.cumsum(a, axis=0), tgt) tgt = np.array( [[1, 2, 3, 4], [6, 8, 10, 13], [16, 11, 14, 18]], ctype) assert_array_equal(np.cumsum(a2, axis=0), tgt) tgt = np.array( [[1, 3, 6, 10], [5, 11, 18, 27], [10, 13, 17, 22]], ctype) assert_array_equal(np.cumsum(a2, axis=1), tgt)
def test_basic(self): ba = [1, 2, 10, 11, 6, 5, 4] ba2 = [[1, 2, 3, 4], [5, 6, 7, 9], [10, 3, 4, 5]] for ctype in [np.int16, np.uint16, np.int32, np.uint32, np.float32, np.float64, np.complex64, np.complex128]: a = np.array(ba, ctype) a2 = np.array(ba2, ctype) if ctype in ['1', 'b']: self.assertRaises(ArithmeticError, np.prod, a) self.assertRaises(ArithmeticError, np.prod, a2, 1) else: assert_equal(a.prod(axis=0), 26400) assert_array_equal(a2.prod(axis=0), np.array([50, 36, 84, 180], ctype)) assert_array_equal(a2.prod(axis=-1), np.array([24, 1890, 600], ctype))
def test_basic(self): ba = [1, 2, 10, 11, 6, 5, 4] ba2 = [[1, 2, 3, 4], [5, 6, 7, 9], [10, 3, 4, 5]] for ctype in [np.int16, np.uint16, np.int32, np.uint32, np.float32, np.float64, np.complex64, np.complex128]: a = np.array(ba, ctype) a2 = np.array(ba2, ctype) if ctype in ['1', 'b']: self.assertRaises(ArithmeticError, np.cumprod, a) self.assertRaises(ArithmeticError, np.cumprod, a2, 1) self.assertRaises(ArithmeticError, np.cumprod, a) else: assert_array_equal(np.cumprod(a, axis=-1), np.array([1, 2, 20, 220, 1320, 6600, 26400], ctype)) assert_array_equal(np.cumprod(a2, axis=0), np.array([[1, 2, 3, 4], [5, 12, 21, 36], [50, 36, 84, 180]], ctype)) assert_array_equal(np.cumprod(a2, axis=-1), np.array([[1, 2, 6, 24], [5, 30, 210, 1890], [10, 30, 120, 600]], ctype))
def test_shuffle(self): # Test lists, arrays (of various dtypes), and multidimensional versions # of both, c-contiguous or not: for conv in [lambda x: np.array([]), lambda x: x, lambda x: np.asarray(x).astype(np.int8), lambda x: np.asarray(x).astype(np.float32), lambda x: np.asarray(x).astype(np.complex64), lambda x: np.asarray(x).astype(object), lambda x: [(i, i) for i in x], lambda x: np.asarray([[i, i] for i in x]), lambda x: np.vstack([x, x]).T, # gh-4270 lambda x: np.asarray([(i, i) for i in x], [("a", object, 1), ("b", np.int32, 1)])]: np.random.seed(self.seed) alist = conv([1, 2, 3, 4, 5, 6, 7, 8, 9, 0]) np.random.shuffle(alist) actual = alist desired = conv([0, 1, 9, 6, 2, 4, 5, 8, 7, 3]) np.testing.assert_array_equal(actual, desired)
def waveform_to_file( station=0, clen=10000, oversample=10, filter_output=False, ): a = rep_seq( create_pseudo_random_code(clen=clen, seed=station), rep=oversample, ) if filter_output: w = numpy.zeros([oversample * clen], dtype=numpy.complex64) fl = (int(oversample + (0.1 * oversample))) w[0:fl] = scipy.signal.blackmanharris(fl) aa = numpy.fft.ifft(numpy.fft.fft(w) * numpy.fft.fft(a)) a = aa / numpy.max(numpy.abs(aa)) a = numpy.array(a, dtype=numpy.complex64) a.tofile('code-l%d-b%d-%06df.bin' % (clen, oversample, station)) else: a.tofile('code-l%d-b%d-%06d.bin' % (clen, oversample, station))
def for_all_dtypes_combination(names=('dtyes',), no_float16=False, no_bool=False, full=None, no_complex=False): """Decorator that checks the fixture with a product set of all dtypes. Args: names(list of str): Argument names to which dtypes are passed. no_float16(bool): If ``True``, ``numpy.float16`` is omitted from candidate dtypes. no_bool(bool): If ``True``, ``numpy.bool_`` is omitted from candidate dtypes. full(bool): If ``True``, then all combinations of dtypes will be tested. Otherwise, the subset of combinations will be tested (see description in :func:`cupy.testing.for_dtypes_combination`). no_complex(bool): If, True, ``numpy.complex64`` and ``numpy.complex128`` are omitted from candidate dtypes. .. seealso:: :func:`cupy.testing.for_dtypes_combination` """ types = _make_all_dtypes(no_float16, no_bool, no_complex) return for_dtypes_combination(types, names, full)
def addFilter(self,f,recenter=False): ''' f can be a function f(x,y) is defined over x = (-w/2, w/2] and y = (-h/2,h/2] and should be centered on the coord 0,0. TODO: At some point this function should be expanded to take filters represented by arrays. ''' if recenter == True: raise NotImplementedError if isinstance(f,GaborWavelet): filt = np.fft.fft2(f.mask(self.tile_size)) self.filters.append(filt) else: w,h = self.tile_size m = np.zeros((w,h),np.complex64) for x in range(-w/2,w/2): for y in range(-h/2,h/2): m[x,y] = f(x,y) filt = np.fft.fft2(m) self.filters.append(filt.conj())
def initSize(self,size): w,h = size w = int(round(w)) h = int(round(h)) if self.x.has_key((w,h)): return # Initializing translations self.x[(w,h)] = [] for x in range(w): mat = np.zeros((w,h),dtype=np.complex64) mat[x,0] = 1.0 filter = np.fft.fft2(mat) self.x[(w,h)].append(filter) self.y[(w,h)] = [] for y in range(h): mat = np.zeros((w,h),dtype=np.complex64) mat[0,y] = 1.0 filter = np.fft.fft2(mat) self.y[(w,h)].append(filter)
def log_power_spectrum_extractor(x, win_len, shift_len, win_type, is_log=False): samples = x.shape[0] frames = (samples - win_len) // shift_len stft = np.zeros((win_len, frames), dtype=np.complex64) spect = np.zeros((win_len // 2 + 1, frames), dtype=np.float64) if win_type == 'hanning': window = np.hanning(win_len) elif win_type == 'hamming': window = np.hamming(win_len) elif win_type == 'rectangle': window = np.ones(win_len) for i in range(frames): one_frame = x[i*shift_len: i*shift_len+win_len] windowed_frame = np.multiply(one_frame, window) stft[:, i] = np.fft.fft(windowed_frame, win_len) if is_log: spect[:, i] = np.log(np.power(np.abs(stft[0: win_len//2+1, i]), 2.)) else: spect[:, i] = np.power(np.abs(stft[0: win_len//2+1, i]), 2.) return spect
def stft_extractor(x, win_len, shift_len, win_type): samples = x.shape[0] frames = (samples - win_len) // shift_len stft = np.zeros((win_len, frames), dtype=np.complex64) spect = np.zeros((win_len // 2 + 1, frames), dtype=np.complex64) if win_type == 'hanning': window = np.hanning(win_len) elif win_type == 'hamming': window = np.hamming(win_len) elif win_type == 'rectangle': window = np.ones(win_len) for i in range(frames): one_frame = x[i*shift_len: i*shift_len+win_len] windowed_frame = np.multiply(one_frame, window) stft[:, i] = np.fft.fft(windowed_frame, win_len) spect[:, i] = stft[: win_len//2+1, i] return spect
def comp_ola_deconv(fs_gpu, ys_gpu, L_gpu, alpha, beta): """ Computes the division in Fourier space needed for direct deconvolution """ sfft = fs_gpu.shape block_size = (16,16,1) grid_size = (int(np.ceil(np.float32(sfft[0]*sfft[1])/block_size[0])), int(np.ceil(np.float32(sfft[2])/block_size[1]))) mod = cu.module_from_buffer(cubin) comp_ola_deconv_Kernel = mod.get_function("comp_ola_deconv_Kernel") z_gpu = cua.zeros(sfft, np.complex64) comp_ola_deconv_Kernel(z_gpu.gpudata, np.int32(sfft[0]), np.int32(sfft[1]), np.int32(sfft[2]), fs_gpu.gpudata, ys_gpu.gpudata, L_gpu.gpudata, np.float32(alpha), np.float32(beta), block=block_size, grid=grid_size) return z_gpu
def comp_ola_gdeconv(xx_gpu, xy_gpu, yx_gpu, yy_gpu, L_gpu, alpha, beta): """ Computes the division in Fourier space needed for gdirect deconvolution """ sfft = xx_gpu.shape block_size = (16,16,1) grid_size = (int(np.ceil(np.float32(sfft[0]*sfft[1])/block_size[0])), int(np.ceil(np.float32(sfft[2])/block_size[1]))) mod = cu.module_from_buffer(cubin) comp_ola_gdeconv_Kernel = mod.get_function("comp_ola_gdeconv_Kernel") z_gpu = cua.zeros(sfft, np.complex64) comp_ola_gdeconv_Kernel(z_gpu.gpudata, np.int32(sfft[0]), np.int32(sfft[1]), np.int32(sfft[2]), xx_gpu, xy_gpu, yx_gpu, yy_gpu, L_gpu.gpudata, np.float32(alpha), np.float32(beta), block=block_size, grid=grid_size) return z_gpu
def comp_ola_sdeconv(gx_gpu, gy_gpu, xx_gpu, xy_gpu, Ftpy_gpu, f_gpu, L_gpu, alpha, beta, gamma=0): """ Computes the division in Fourier space needed for sparse deconvolution """ sfft = xx_gpu.shape block_size = (16,16,1) grid_size = (int(np.ceil(np.float32(sfft[0]*sfft[1])/block_size[0])), int(np.ceil(np.float32(sfft[2])/block_size[1]))) mod = cu.module_from_buffer(cubin) comp_ola_sdeconv_Kernel = mod.get_function("comp_ola_sdeconv_Kernel") z_gpu = cua.zeros(sfft, np.complex64) comp_ola_sdeconv_Kernel(z_gpu.gpudata, np.int32(sfft[0]), np.int32(sfft[1]), np.int32(sfft[2]), gx_gpu.gpudata, gy_gpu.gpudata, xx_gpu.gpudata, xy_gpu.gpudata, Ftpy_gpu.gpudata, f_gpu.gpudata, L_gpu.gpudata, np.float32(alpha), np.float32(beta), np.float32(gamma), block=block_size, grid=grid_size) return z_gpu
def fst_delay_snd(fst, snd, samp_rate): # Verify argument shape. s1, s2 = fst.shape, snd.shape if len(s1) != 1 or len(s2) != 1 or s1[0] != s2[0]: raise Exception("Argument shape invalid, in 'fst_delay_snd' function") length = s1[0] half_len = int(length / 2) Xfst = numpy.fft.fft(fst) Xsnd_star = numpy.conj(numpy.fft.fft(snd)) Xall = numpy.zeros(length, dtype=numpy.complex64) for i in range(0, length): if Xsnd_star[i] == 0 or Xfst[i] == 0: Xall[i] = 0 else: Xall[i] = (Xsnd_star[i] * Xfst[i]) / abs(Xsnd_star[i]) / abs(Xfst[i]) R = numpy.fft.ifft(Xall) max_pos = numpy.argmax(R) if max_pos > half_len: return -(length - 1 - max_pos) / samp_rate else: return max_pos / samp_rate
def fst_delay_snd(fst, snd, samp_rate): # Verify argument shape. s1, s2 = fst.shape, snd.shape if len(s1) != 1 or len(s2) != 1 or s1[0] != s2[0]: raise Exception("Argument shape invalid, in 'fst_delay_snd' function") length = s1[0] half_len = int(length / 2) Xfst = numpy.fft.fft(fst) Xsnd = numpy.fft.fft(snd) Xsnd_star = numpy.conj(Xsnd) Xall = numpy.zeros(length, dtype=numpy.complex64) for i in range(0, length): Xall[i] = (Xsnd_star[i] * Xfst[i]) / abs(Xsnd_star[i]) / abs(Xfst[i]) R = numpy.fft.ifft(Xall) max_pos = numpy.argmax(R) if max_pos > half_len: delta = -(length - 1 - max_pos) / samp_rate else: delta = max_pos / samp_rate return delta, R, Xfst, Xsnd
def __init__(self, test_array, complex_numbers=False): """Figure out data settings from the test array. @param[in] test_array A list or numpy array containing test data""" super(TestingBlock, self).__init__() if isinstance(test_array, np.ndarray): if test_array.dtype == np.complex64: complex_numbers = True if complex_numbers: self.test_array = np.array(test_array).astype(np.complex64) header = { 'nbit': 64, 'dtype': 'complex64', 'shape': self.test_array.shape} self.dtype = np.complex64 else: self.test_array = np.array(test_array).astype(np.float32) header = { 'nbit': 32, 'dtype': 'float32', 'shape': self.test_array.shape} self.dtype = np.float32 self.output_header = json.dumps(header)
def main(self, input_rings, output_rings): """ @param[in] input_rings First ring in this list will be used for data @param[out] output_rings First ring in this list will be used for data output.""" data_accumulate = None for ispan in self.iterate_ring_read(input_rings[0]): if self.nbit < 8: unpacked_data = unpack(ispan.data_view(self.dtype), self.nbit) else: unpacked_data = ispan.data_view(self.dtype) if data_accumulate is not None: data_accumulate = np.concatenate((data_accumulate, unpacked_data[0])) else: data_accumulate = unpacked_data[0] if self.shape != [1, 1]: data_accumulate = np.reshape(data_accumulate, (self.shape[0], -1)) data_accumulate = data_accumulate.astype(np.complex64) self.out_gulp_size = data_accumulate.nbytes outspan_generator = self.iterate_ring_write(output_rings[0]) ospan = outspan_generator.next() result = np.fft.fft(data_accumulate).astype(np.complex64) ospan.data_view(np.complex64)[0] = result.ravel()
def main(self, input_rings, output_rings): """ @param[in] input_rings First ring in this list will be used for data input. @param[out] output_rings First ring in this list will be used for data output.""" data_accumulate = None for ispan in self.iterate_ring_read(input_rings[0]): if self.nbit < 8: unpacked_data = unpack(ispan.data_view(self.dtype), self.nbit) else: unpacked_data = ispan.data_view(self.dtype) if data_accumulate is not None: data_accumulate = np.concatenate((data_accumulate, unpacked_data[0])) else: data_accumulate = unpacked_data[0] data_accumulate = data_accumulate.astype(np.complex64) self.out_gulp_size = data_accumulate.nbytes outspan_generator = self.iterate_ring_write(output_rings[0]) ospan = outspan_generator.next() result = np.fft.ifft(data_accumulate) ospan.data_view(np.complex64)[0][:] = result[:]
def numpy2bifrost(dtype): if dtype == np.int8: return _bf.BF_DTYPE_I8 elif dtype == np.int16: return _bf.BF_DTYPE_I16 elif dtype == np.int32: return _bf.BF_DTYPE_I32 elif dtype == np.uint8: return _bf.BF_DTYPE_U8 elif dtype == np.uint16: return _bf.BF_DTYPE_U16 elif dtype == np.uint32: return _bf.BF_DTYPE_U32 elif dtype == np.float16: return _bf.BF_DTYPE_F16 elif dtype == np.float32: return _bf.BF_DTYPE_F32 elif dtype == np.float64: return _bf.BF_DTYPE_F64 elif dtype == np.float128: return _bf.BF_DTYPE_F128 elif dtype == ci8: return _bf.BF_DTYPE_CI8 elif dtype == ci16: return _bf.BF_DTYPE_CI16 elif dtype == ci32: return _bf.BF_DTYPE_CI32 elif dtype == cf16: return _bf.BF_DTYPE_CF16 elif dtype == np.complex64: return _bf.BF_DTYPE_CF32 elif dtype == np.complex128: return _bf.BF_DTYPE_CF64 elif dtype == np.complex256: return _bf.BF_DTYPE_CF128 else: raise ValueError("Unsupported dtype: " + str(dtype))
def numpy2string(dtype): if dtype == np.int8: return 'i8' elif dtype == np.int16: return 'i16' elif dtype == np.int32: return 'i32' elif dtype == np.int64: return 'i64' elif dtype == np.uint8: return 'u8' elif dtype == np.uint16: return 'u16' elif dtype == np.uint32: return 'u32' elif dtype == np.uint64: return 'u64' elif dtype == np.float16: return 'f16' elif dtype == np.float32: return 'f32' elif dtype == np.float64: return 'f64' elif dtype == np.float128: return 'f128' elif dtype == np.complex64: return 'cf32' elif dtype == np.complex128: return 'cf64' elif dtype == np.complex256: return 'cf128' else: raise TypeError("Unsupported dtype: " + str(dtype))
def run_test_matmul_aa_ci8_shape(self, shape, transpose=False): # **TODO: This currently never triggers the transpose path in the backend shape_complex = shape[:-1] + (shape[-1] * 2,) # Note: The xGPU-like correlation kernel does not support input values of -128 (only [-127:127]) a8 = ((np.random.random(size=shape_complex) * 2 - 1) * 127).astype(np.int8) a_gold = a8.astype(np.float32).view(np.complex64) if transpose: a_gold = H(a_gold) # Note: np.matmul seems to be slow and inaccurate when there are batch dims c_gold = np.matmul(a_gold, H(a_gold)) triu = np.triu_indices(shape[-2] if not transpose else shape[-1], 1) c_gold[..., triu[0], triu[1]] = 0 a = a8.view(bf.DataType.ci8) a = bf.asarray(a, space='cuda') if transpose: a = H(a) c = bf.zeros_like(c_gold, space='cuda') self.linalg.matmul(1, a, None, 0, c) c = c.copy('system') np.testing.assert_allclose(c, c_gold, RTOL, ATOL)
def run_test_matmul_ab_ci8_shape(self, shape, k, transpose=False): ashape_complex = shape[:-2] + (shape[-2], k * 2) bshape_complex = shape[:-2] + (k, shape[-1] * 2) a8 = (np.random.random(size=ashape_complex) * 255).astype(np.int8) b8 = (np.random.random(size=bshape_complex) * 255).astype(np.int8) a_gold = a8.astype(np.float32).view(np.complex64) b_gold = b8.astype(np.float32).view(np.complex64) if transpose: a_gold, b_gold = H(b_gold), H(a_gold) c_gold = np.matmul(a_gold, b_gold) a = a8.view(bf.DataType.ci8) b = b8.view(bf.DataType.ci8) a = bf.asarray(a, space='cuda') b = bf.asarray(b, space='cuda') if transpose: a, b = H(b), H(a) c = bf.zeros_like(c_gold, space='cuda') self.linalg.matmul(1, a, b, 0, c) c = c.copy('system') np.testing.assert_allclose(c, c_gold, RTOL, ATOL)
def run_benchmark_matmul_aa_correlator_kernel(self, ntime, nstand, nchan): x_shape = (ntime, nchan, nstand*2) perm = [1,0,2] x8 = ((np.random.random(size=x_shape+(2,))*2-1)*127).astype(np.int8) x = x8.astype(np.float32).view(np.complex64).reshape(x_shape) x = x.transpose(perm) b_gold = np.matmul(H(x[:,[0],:]), x[:,[0],:]) triu = np.triu_indices(x_shape[-1], 1) b_gold[..., triu[0], triu[1]] = 0 x = x8.view(bf.DataType.ci8).reshape(x_shape) x = bf.asarray(x, space='cuda') x = x.transpose(perm) b = bf.zeros_like(b_gold, space='cuda') bf.device.stream_synchronize(); t0 = time.time() nrep = 200 for _ in xrange(nrep): self.linalg.matmul(1, None, x, 0, b) bf.device.stream_synchronize(); dt = time.time() - t0 nflop = nrep * nchan * ntime * nstand*(nstand+1)/2 * 2*2 * 8 print nstand, '\t', nflop / dt / 1e9, 'GFLOP/s' print '\t\t', nrep*ntime*nchan / dt / 1e6, 'MHz'
def run_test_c2c_impl(self, shape, axes, inverse=False, fftshift=False): shape = list(shape) shape[-1] *= 2 # For complex known_data = np.random.normal(size=shape).astype(np.float32).view(np.complex64) idata = bf.ndarray(known_data, space='cuda') odata = bf.empty_like(idata) fft = Fft() fft.init(idata, odata, axes=axes, apply_fftshift=fftshift) fft.execute(idata, odata, inverse) if inverse: if fftshift: known_data = np.fft.ifftshift(known_data, axes=axes) # Note: Numpy applies normalization while CUFFT does not norm = reduce(lambda a, b: a * b, [known_data.shape[d] for d in axes]) known_result = gold_ifftn(known_data, axes=axes) * norm else: known_result = gold_fftn(known_data, axes=axes) if fftshift: known_result = np.fft.fftshift(known_result, axes=axes) x = (np.abs(odata.copy('system') - known_result) / known_result > RTOL).astype(np.int32) a = odata.copy('system') b = known_result compare(odata.copy('system'), known_result)
def run_test_c2r_impl(self, shape, axes, fftshift=False): ishape = list(shape) oshape = list(shape) ishape[axes[-1]] = shape[axes[-1]] // 2 + 1 oshape[axes[-1]] = (ishape[axes[-1]] - 1) * 2 ishape[-1] *= 2 # For complex known_data = np.random.normal(size=ishape).astype(np.float32).view(np.complex64) idata = bf.ndarray(known_data, space='cuda') odata = bf.ndarray(shape=oshape, dtype='f32', space='cuda') fft = Fft() fft.init(idata, odata, axes=axes, apply_fftshift=fftshift) fft.execute(idata, odata) # Note: Numpy applies normalization while CUFFT does not norm = reduce(lambda a, b: a * b, [shape[d] for d in axes]) if fftshift: known_data = np.fft.ifftshift(known_data, axes=axes) known_result = gold_irfftn(known_data, axes=axes) * norm compare(odata.copy('system'), known_result)
def test_data_sizes(self): """Test that different number of bits give correct throughput size""" for iterate in range(5): nbit = 2**iterate if nbit == 8: continue self.blocks[0] = ( SigprocReadBlock( './data/2chan' + str(nbit) + 'bitNoDM.fil'), [], [0]) open(self.logfile, 'w').close() Pipeline(self.blocks).main() number_fftd = np.loadtxt(self.logfile).astype(np.float32).view(np.complex64).size # Compare with simple copy self.blocks[1] = (CopyBlock(), [0], [1]) open(self.logfile, 'w').close() Pipeline(self.blocks).main() number_copied = np.loadtxt(self.logfile).size self.assertEqual(number_fftd, number_copied) # Go back to FFT self.blocks[1] = (FFTBlock(gulp_size=4096 * 8 * 8 * 8), [0], [1])
def test_equivalent_data_to_copy(self): """Test that the data coming out of this pipeline is equivalent the initial read data""" self.logfile = '.log.txt' self.blocks = [] self.blocks.append(( SigprocReadBlock( './data/1chan8bitNoDM.fil'), [], [0])) self.blocks.append((FFTBlock(gulp_size=4096 * 8 * 8 * 8 * 8), [0], [1])) self.blocks.append((IFFTBlock(gulp_size=4096 * 8 * 8 * 8 * 8), [1], [2])) self.blocks.append((WriteAsciiBlock(self.logfile), [2], [])) open(self.logfile, 'w').close() Pipeline(self.blocks).main() unfft_result = np.loadtxt(self.logfile).astype(np.float32).view(np.complex64) self.blocks[1] = (CopyBlock(), [0], [1]) self.blocks[2] = (WriteAsciiBlock(self.logfile), [1], []) del self.blocks[3] open(self.logfile, 'w').close() Pipeline(self.blocks).main() untouched_result = np.loadtxt(self.logfile).astype(np.float32) np.testing.assert_almost_equal(unfft_result, untouched_result, 2)
def default(self, obj): # convert dates and numpy objects in a json serializable format if isinstance(obj, datetime): return obj.strftime('%Y-%m-%dT%H:%M:%SZ') elif isinstance(obj, date): return obj.strftime('%Y-%m-%d') elif type(obj) in (np.int_, np.intc, np.intp, np.int8, np.int16, np.int32, np.int64, np.uint8, np.uint16, np.uint32, np.uint64): return int(obj) elif type(obj) in (np.bool_,): return bool(obj) elif type(obj) in (np.float_, np.float16, np.float32, np.float64, np.complex_, np.complex64, np.complex128): return float(obj) # Let the base class default method raise the TypeError return json.JSONEncoder.default(self, obj)
def testTwoSessions(self): optimizer = ctf.train.CplxAdamOptimizer() g = tf.Graph() with g.as_default(): with tf.Session(): var0 = tf.Variable(np.array([1.0+1.0j, 2.0+2.0j], dtype=np.complex64), name="v0") grads0 = tf.constant(np.array([0.1+0.1j, 0.1+0.1j], dtype=np.complex64)) optimizer.apply_gradients([(grads0, var0)]) gg = tf.Graph() with gg.as_default(): with tf.Session(): var0 = tf.Variable(np.array([1.0+1.0j, 2.0+2.0j], dtype=np.complex64), name="v0") grads0 = tf.constant(np.array([0.1+0.1j, 0.1+0.1j], dtype=np.complex64)) # If the optimizer saves any state not keyed by graph the following line # fails. optimizer.apply_gradients([(grads0, var0)])
def _testTypes(self, vals): for dtype in [np.complex64]: x = np.zeros(vals.shape).astype(dtype) y = vals.astype(dtype) var_value, op_value = self._initAssignFetch(x, y, use_gpu=False) self.assertAllEqual(y, var_value) self.assertAllEqual(y, op_value) var_value, op_value = self._initAssignAddFetch(x, y, use_gpu=False) self.assertAllEqual(x + y, var_value) self.assertAllEqual(x + y, op_value) var_value, op_value = self._initAssignSubFetch(x, y, use_gpu=False) self.assertAllEqual(x - y, var_value) self.assertAllEqual(x - y, op_value) if tf.test.is_built_with_cuda() and dtype in [np.float32, np.float64]: var_value, op_value = self._initAssignFetch(x, y, use_gpu=True) self.assertAllEqual(y, var_value) self.assertAllEqual(y, op_value) var_value, op_value = self._initAssignAddFetch(x, y, use_gpu=True) self.assertAllEqual(x + y, var_value) self.assertAllEqual(x + y, op_value) var_value, op_value = self._initAssignSubFetch(x, y, use_gpu=False) self.assertAllEqual(x - y, var_value) self.assertAllEqual(x - y, op_value)
def testHStack(self): with self.test_session(force_gpu=True): p1 = array_ops.placeholder(dtypes.complex64, shape=[4, 4]) p2 = array_ops.placeholder(dtypes.complex64, shape=[4, 4]) c = array_ops.concat([p1, p2], 0) params = { p1: (np.random.rand(4, 4) + 1j*np.random.rand(4, 4)).astype(np.complex64), p2: (np.random.rand(4, 4) + 1j*np.random.rand(4, 4)).astype(np.complex64), } result = c.eval(feed_dict=params) self.assertEqual(result.shape, c.get_shape()) self.assertAllEqual(result[:4, :], params[p1]) self.assertAllEqual(result[4:, :], params[p2])
def testVStack(self): with self.test_session(force_gpu=True): p1 = array_ops.placeholder(dtypes.complex64, shape=[4, 4]) p2 = array_ops.placeholder(dtypes.complex64, shape=[4, 4]) c = array_ops.concat([p1, p2], 1) params = { p1: (np.random.rand(4, 4) + 1j*np.random.rand(4, 4)).astype(np.complex64), p2: (np.random.rand(4, 4) + 1j*np.random.rand(4, 4)).astype(np.complex64), } result = c.eval(feed_dict=params) self.assertEqual(result.shape, c.get_shape()) self.assertAllEqual(result[:, :4], params[p1]) self.assertAllEqual(result[:, 4:], params[p2])
def testGradientWithUnknownInputDim(self): with self.test_session(use_gpu=True): x = array_ops.placeholder(dtypes.complex64) y = array_ops.placeholder(dtypes.complex64) c = array_ops.concat([x, y], 2) output_shape = [10, 2, 9] grad_inp = (np.random.rand(*output_shape) + 1j*np.random.rand(*output_shape)).astype(np.complex64) grad_tensor = constant_op.constant( [inp for inp in grad_inp.flatten()], shape=output_shape) grad = gradients_impl.gradients([c], [x, y], [grad_tensor]) concated_grad = array_ops.concat(grad, 2) params = { x: (np.random.rand(10, 2, 3) + 1j*np.random.rand(10, 2, 3)).astype(np.complex64), y: (np.random.rand(10, 2, 6) + 1j*np.random.rand(10, 2, 6)).astype(np.complex64), } result = concated_grad.eval(feed_dict=params) self.assertAllEqual(result, grad_inp)
def testShapeWithUnknownConcatDim(self): p1 = array_ops.placeholder(dtypes.complex64) c1 = constant_op.constant(np.complex64(10.0+0j), shape=[4, 4, 4, 4]) p2 = array_ops.placeholder(dtypes.complex64) c2 = constant_op.constant(np.complex64(20.0+0j), shape=[4, 4, 4, 4]) dim = array_ops.placeholder(dtypes.int32) concat = array_ops.concat([p1, c1, p2, c2], dim) self.assertEqual(4, concat.get_shape().ndims) # All dimensions unknown. concat2 = array_ops.concat([p1, p2], dim) self.assertEqual(None, concat2.get_shape()) # Rank doesn't match. c3 = constant_op.constant(np.complex64(30.0+0j), shape=[4, 4, 4]) with self.assertRaises(ValueError): array_ops.concat([p1, c1, p2, c3], dim)