Miscellaneous

ducc0.misc

Various unsorted utilities

Notes

The functionality in this module is not considered to have a stable interface and also may be moved to other modules in the future. If you use it, be prepared to adjust your code at some point in the future!

ducc0.misc.GL_thetas(nlat: int) numpy.ndarray
ducc0.misc.GL_weights(nlat: int, nlon: int) numpy.ndarray
class ducc0.misc.OofaNoise

Class for computing noise with a power spectrum that has a given slope between a minimum frequency f_min and a knee frequency f_knee, and is white outside this region.

Original implementation by Stephane Plaszczynski; for details see https://arxiv.org/abs/astro-ph/0510081.

filterGaussian(self: ducc0.misc.OofaNoise, rnd: numpy.ndarray) numpy.ndarray

Apply the noise filter to input Gaussian noise

Parameters:

rnd (numpy.ndarray((nsamples,), dtype=numpy.float64)) – input Gaussian random numbers with mean=0 and sigma=1

Returns:

the filtered noise samples with the requested spectral shape.

Return type:

numpy.ndarray((nsamples,), dtype=numpy.float64)

Notes

Subsequent calls to this method will continue the same noise stream; i.e. it is possible to generate a very long noise time stream chunk by chunk. To generate multiple independent noise streams, use different OofaNoise objects (and supply them with independent Gaussian noise streams)!

ducc0.misc.empty_noncritical(shape: list[int], dtype: object) numpy.ndarray

Creates an uninitialized array of the requested shape and data type, with a memory layout that avoids critical strides.

As an example, the array generated by np.zeros((1024,1024)) has a critical stride of 8192 bytes for axis 0. Accessing the array elements in the order [0, 0], [1, 0], [2,0] etc. is very slow in this case, but exactly this kind of access is required for, e.g., computing the FFT over axis 0 of the array. If this kind of critical stride is detected, empty_noncritical() embeds the array into a slightly larger array, in this case with the shape [1024,1027], which improves access times significantly.

This routine considers every stride which is a multiple of 4096 bytes as critical, which should be a decent heuristic.

Parameters:
  • shape (sequence of int) – the array shape

  • dtype – the requested data type

Returns:

An uninitialized numpy array with the requested properties

Return type:

numpy.ndarray (shape, dtype=dtype)

ducc0.misc.get_correction(beta: float, e0: float, W: int, npoints: int, dx: float) numpy.ndarray
ducc0.misc.get_deflected_angles(theta: numpy.ndarray, phi0: numpy.ndarray, nphi: numpy.ndarray, ringstart: numpy.ndarray, deflect: numpy.ndarray, calc_rotation: bool = False, res: object = None, nthreads: int = 1, dphi: object = None) numpy.ndarray

Obtains new pointing angles on the sky according to a deflection field on a set of isolatitude rings

Parameters:
  • theta (numpy.ndarray((nrings,), dtype=numpy.float64)) – colatitudes of the rings (nrings any number)

  • phi0 (numpy.ndarray((nrings,), dtype=numpy.float64)) – longitude of the first pixel in each ring

  • ringstart (numpy.ndarray((nrings,), dtype=numpy.uint64)) – index of the first pixel of each ring in output map

  • deflect (numpy.ndarray((npix, 2), dtype=numpy.float32 or numpy.float64)) – Spin-1 deflection field, with real and imaginary comp in first and second entry (typically, the output of a spin-1 alm2map_spin transform) The array layout and npix must be consistent with the given geometry

  • calc_rotation(optional) (boolean) – If set, also returns the phase correction (gamma in astro-ph/0502469v3)

  • nthreads(optional) (int) – Number of threads to use. Defaults to 1

  • res(optional) (numpy.ndarray((npix, 3 if calc_rotation is set or 2), same dtype as deflect)) – output array, containing new co-latitudes, new longitudes, and rotation gammma if required

  • dphi(optional) (numpy.ndarray((nrings,), dtype=numpy.float64)) – azimuthal distance between pixels in each ring (in radians)

Returns:

numpy.ndarray

Return type:

identical to res

ducc0.misc.get_kernel(beta: float, e0: float, W: int, npoints: int) numpy.ndarray
ducc0.misc.l2error(a: object, b: object) float

Compute the L2 error between two arrays or scalars. More specifically, compute sqrt(sum_i(|a_i - b_i|^2) / max(sum_i(|a_i|^2), sum_i(|b_i|^2))), where i goes over all array elements.

Parameters:
  • a (scalar or numpy.ndarray of any shape; dtype must be a float or complex type) –

  • b (scalar or numpy.ndarray of the same shape as a; dtype must be a float or complex type) –

Returns:

the L2 error between the two objects.

Return type:

float

Notes

The accumulations are performed in long double precision for good accuracy.

ducc0.misc.lensing_rotate(values: numpy.ndarray, gamma: numpy.ndarray, spin: int, nthreads: int = 1) None

Rotates complex values depending on given angles and spin

Parameters:
  • values (numpy.ndarray(<any shape>, dtype=complex)) – values to rotate; operation is applied in place

  • gamma (numpy.ndarray(same shape as values, dtype=float with same prcision as values)) – rotation angles

  • spin (int) – will be multiplied to gamma before rotation

  • nthreads(optional) (int) – Number of threads to use. Defaults to 1

ducc0.misc.make_noncritical(in: numpy.ndarray) numpy.ndarray

Returns a copy of the input array with a memory layout that avoids critical strides.

As an example, the array generated by np.zeros((1024,1024)) has a critical stride of 8192 bytes for axis 0. Accessing the array elements in the order [0, 0], [1, 0], [2,0] etc. is very slow in this case, but exactly this kind of access is required for, e.g., computing the FFT over axis 0 of the array. If this kind of critical stride is detected, make_noncritical() embeds the array into a slightly larger array, in this case with the shape [1024,1027], which improves access times significantly.

This routine considers every stride which is a multiple of 4096 bytes as critical, which should be a decent heuristic.

Parameters:

in (numpy.ndarray (float or integer dtype)) – the input array

Returns:

A copy of the array with noncritical strides

Return type:

numpy.ndarray (same dtype and content as in)

ducc0.misc.preallocate_memory(gbytes: float) bool
ducc0.misc.roll_resize_roll(inp: numpy.ndarray, out: numpy.ndarray, roll_inp: list[int], roll_out: list[int], nthreads: int = 1) numpy.ndarray

Performs operations equivalent to

tmp = np.roll(inp, roll_inp, axis=tuple(range(inp.ndim))) tmp2 = np.zeros(out.shape, dtype=inp.dtype) slices = tuple(slice(0, min(s1, s2)) for s1, s2 in zip(inp.shape, out.shape)) tmp2[slices] = tmp[slices] out[()] = np.roll(tmp2, roll_out, axis=tuple(range(out.ndim))) return out

Parameters:
  • inp (numpy.ndarray(any shape, dtype=float or complex)) – input array

  • out (numpy.ndarray(any shape, same dimensionality and dtype as in)) – output array

  • roll_inp (tuple(int), length=inp.ndim) – amount of rolling for the input array

  • roll_out (tuple(int), length=out.ndim) – amount of rolling for the output array

  • nthreads (int) – Number of threads to use. If 0, use the system default (typically the number of hardware threads on the compute node).

Returns:

numpy.ndarray

Return type:

identical to out

Notes

inp and out must not overlap in memory.

ducc0.misc.transpose(in: numpy.ndarray, out: numpy.ndarray, nthreads: int = 1) numpy.ndarray
ducc0.misc.vdot(a: object, b: object) object

Compute the scalar product of two arrays or scalars., i.e. sum_i(conj(a_i)*b_i) over all array elements.

Parameters:
  • a (scalar or numpy.ndarray of any shape; dtype must be a float or complex type) –

  • b (scalar or numpy.ndarray of the same shape as a; dtype must be a float or complex type) –

Returns:

the scalar product. If the result can be represented by a float value, it will be returned as float, otherwise as complex.

Return type:

float or complex

Notes

The accumulation is performed in long double precision for good accuracy.

ducc0.misc.wigner3j_int(l2: int, l3: int, m2: int, m3: int) object

Computes Wigner 3j symbols according to the algorithm of Schulten & Gordon: J. Math. Phys. 16, p. 10 (1975)

This special case only takes integer quantum numbers.

Parameters:
  • l2 (integer) – fixed quantum numbers

  • l3 (integer) – fixed quantum numbers

  • m2 (integer) – fixed quantum numbers

  • m3 (integer) – fixed quantum numbers

Returns:

  • int (the l1 quantum number of the first value in the returned array)

  • numpy.ndarray(dtype=numpy.float64) (3j symbols in order of increasing l1)