Here are the examples of the python api numpy.rollaxis taken from open source projects. By voting up you can indicate which examples are most useful and appropriate.
78 Examples
3
Example 1
def load_batch(fname):
with open('data/cifar-10-batches-py/%s'%fname, 'rb') as f:
data = pickle.load(f)
x = data['data'].reshape((-1,3,32,32)).astype('float32')/255
x = np.rollaxis(x, 1, 4)
y = np.array(data['labels'])
return x, y
3
Example 2
def log_normalize(v, axis=0):
"""Normalized probabilities from unnormalized log-probabilites"""
v = np.rollaxis(v, axis)
v = v.copy()
v -= v.max(axis=0)
out = logsumexp(v)
v = np.exp(v - out)
v += np.finfo(np.float32).eps
v /= np.sum(v, axis=0)
return np.swapaxes(v, 0, axis)
3
Example 3
@deprecated("The function log_normalize is deprecated in 0.18 and "
"will be removed in 0.20.")
def log_normalize(v, axis=0):
"""Normalized probabilities from unnormalized log-probabilites"""
v = np.rollaxis(v, axis)
v = v.copy()
v -= v.max(axis=0)
out = logsumexp(v)
v = np.exp(v - out)
v += np.finfo(np.float32).eps
v /= np.sum(v, axis=0)
return np.swapaxes(v, 0, axis)
3
Example 4
def center_and_norm(x, axis=-1):
""" Centers and norms x **in place**
Parameters
-----------
x: ndarray
Array with an axis of observations (statistical units) measured on
random variables.
axis: int, optional
Axis along which the mean and variance are calculated.
"""
x = np.rollaxis(x, axis)
x -= x.mean(axis=0)
x /= x.std(axis=0)
3
Example 5
Project: scipy Source File: _cubic.py
def __init__(self, x, y, axis=0, extrapolate=None):
x = _asarray_validated(x, check_finite=False, as_inexact=True)
y = _asarray_validated(y, check_finite=False, as_inexact=True)
axis = axis % y.ndim
xp = x.reshape((x.shape[0],) + (1,)*(y.ndim-1))
yp = np.rollaxis(y, axis)
dk = self._find_derivatives(xp, yp)
data = np.hstack((yp[:, None, ...], dk[:, None, ...]))
_b = BPoly.from_derivatives(x, data, orders=None)
super(PchipInterpolator, self).__init__(_b.c, _b.x,
extrapolate=extrapolate)
self.axis = axis
3
Example 6
Project: eofs Source File: xarray.py
def _map_and_coords(pcs, field, mapfunc, *args, **kwargs):
info = {}
for array in (field, pcs):
info[array.name] = _coord_info(array)
cmap_args = [np.rollaxis(array.values, info[array.name][0])
for array in (pcs, field)]
cmap_args += args
cmap = mapfunc(*cmap_args, **kwargs)
dim_args = [info[array.name][1] for array in (pcs, field)]
dims = covcor_dimensions(*dim_args)
return cmap, dims
3
Example 7
Project: menpofit Source File: gn.py
def ds_dp(self):
r"""
Calculates the shape jacobian. That is
.. math::
\frac{d\mathcal{S}}{d\mathbf{p}} = \mathbf{J}_S = \mathbf{U}_S
with size :math:`2 \times n \times n_S`.
:type: `ndarray`
"""
return np.rollaxis(self.transform.d_dp(None), -1)
3
Example 8
Project: hebel Source File: plotting.py
def show_filters(W, img_dims, columns=10, normalize=True, **kwargs):
import matplotlib.pyplot as plt
if isinstance(W, gpuarray.GPUArray): W = W.get()
D, N = W.shape
if normalize:
W = W - W.min() #[np.newaxis,:]
W = W / W.max() #[np.newaxis,:]
rows = int(ceil(N / columns))
fig = plt.figure(1, **kwargs)
plt.subplots_adjust(left=0., right=.51, wspace=.1, hspace=.01)
filters = np.rollaxis(W.reshape(img_dims + (N,)), 2)
filters = np.vstack([np.hstack(filters[i:i+columns]) for i in range(0, N, columns)])
plt.axis('off')
plt.imshow(filters, cmap=plt.cm.gray, interpolation='nearest', figure=fig)
3
Example 9
def _precompute(self):
# compute warp jacobian
dW_dp = np.rollaxis(
self.transform.d_dp(self.template.indices()), -1)
self.dW_dp = dW_dp.reshape(dW_dp.shape[:1] + self.template.shape +
dW_dp.shape[-1:])
3
Example 10
def fct0(self, u, fu, S):
"""Fast Cheb transform of x-direction. No FFT, just align data in x-direction and do fct."""
U_mpi2 = self.work_arrays[((self.num_processes, self.Np[0], self.Np[1], self.N[2]), self.float, 0)]
UT = self.work_arrays[((self.N[0], self.Np[1], self.N[2]), self.float, 0)]
U_mpi2[:] = rollaxis(u.reshape(self.Np[0], self.num_processes, self.Np[1], self.N[2]), 1)
self.comm.Alltoall([U_mpi2, self.mpitype], [UT, self.mpitype])
fu = S.fct(UT, fu)
return fu
3
Example 11
def _maxout(x, W, b):
W_r = numpy.rollaxis(W, 2)
y = numpy.tensordot(_as_mat(x), W_r, axes=1)
if b is not None:
y += b
return numpy.max(y, axis=2)
3
Example 12
def ols(self, axis):
y = np.rollaxis(self.y, 0, axis+1) ## time index is axis
X = self.X
m = glm(y, X, axis=axis)
m1 = glm(y, X, axis=axis, method='kalman')
b = m.beta
b1 = m1.beta
tcon = m.contrast([1,0])
tcon1 = m1.contrast([1,0])
z = tcon.zscore()
z1 = tcon1.zscore()
assert_almost_equal(b, b1)
3
Example 13
Project: dipy Source File: registration_example.py
def save_volumes_as_mosaic(fname,volume_list):
import Image
vols=[]
for vol in volume_list:
vol=np.rollaxis(vol,2,1)
sh=vol.shape
arr=vol.reshape(sh[0],sh[1]*sh[2])
arr=np.interp(arr,[arr.min(),arr.max()],[0,255])
arr=arr.astype('ubyte')
print 'arr.shape',arr.shape
vols.append(arr)
mosaic=np.concatenate(vols)
Image.fromarray(mosaic).save(fname)
3
Example 14
Project: windspharm Source File: tools.py
def __order_dims(d, inorder):
if 'x' not in inorder or 'y' not in inorder:
raise ValueError('a latitude-longitude grid is required')
lonpos = inorder.lower().find('x')
latpos = inorder.lower().find('y')
d = np.rollaxis(d, lonpos)
if latpos < lonpos:
latpos += 1
d = np.rollaxis(d, latpos)
outorder = inorder.replace('x', '')
outorder = outorder.replace('y', '')
outorder = 'yx' + outorder
return d, outorder
3
Example 15
Project: thunder Source File: blocks.py
def toseries(self):
"""
Converts blocks to series.
"""
from thunder.series.series import Series
if self.mode == 'spark':
values = self.values.values_to_keys(tuple(range(1, len(self.shape)))).unchunk()
if self.mode == 'local':
values = self.values.unchunk()
values = rollaxis(values, 0, values.ndim)
return Series(values)
3
Example 16
def resize(img, size):
"""resize image into size x size
Attributes:
img: The image to resize
size: The size to resize the image
"""
img = cv2.resize(img, (size, size))
img = np.rollaxis(img, 2)
return img
3
Example 17
Project: instaseis Source File: reciprocal_merged_instaseis_db.py
def _get_and_reorder_utemp(self, id_elem):
# We can now read it in a single go!
utemp = self.meshes.merged.f["MergedSnapshots"][id_elem]
# utemp is currently (nvars, jpol, ipol, npts)
# 1. Roll to (npts, nvar, jpol, ipol)
utemp = np.rollaxis(utemp, 3, 0)
# 2. Roll to (npts, jpol, nvar, ipol)
utemp = np.rollaxis(utemp, 2, 1)
# 3. Roll to (npts, jpol, ipol, nvar)
utemp = np.rollaxis(utemp, 3, 2)
return utemp
3
Example 18
def _precompute(self):
# compute warp jacobian
dW_dp = np.rollaxis(self.transform.d_dp(self.template.indices()), -1)
dW_dp = dW_dp.reshape(dW_dp.shape[:1] + self.template.shape +
dW_dp.shape[-1:])
# compute steepest descent images
self.filtered_J, J = self.residual.steepest_descent_images(
self.template, dW_dp)
# compute hessian
self.H = self.residual.hessian(self.filtered_J, sdi2=J)
3
Example 19
Project: opendr Source File: geometry.py
def compute_d1(self):
# To stay consistent with numpy, we must upgrade 1D arrays to 2D
mtx1r = row(self.mtx1.r) if len(self.mtx1.r.shape)<2 else self.mtx1.r
mtx2r = col(self.mtx2.r) if len(self.mtx2.r.shape)<2 else self.mtx2.r
if mtx1r.ndim <= 2:
return sp.kron(sp.eye(mtx1r.shape[0], mtx1r.shape[0]),mtx2r.T)
else:
mtx2f = mtx2r.reshape((-1, mtx2r.shape[-2], mtx2r.shape[-1]))
mtx2f = np.rollaxis(mtx2f, -1, -2) #transpose basically
result = sp.block_diag([np.kron(np.eye(mtx1r.shape[-2], mtx1r.shape[-2]),m2) for m2 in mtx2f])
assert(result.shape[0] == self.r.size)
return result
3
Example 20
def reconstruct(time_series, images, axis=0):
# Reconstruct data from remaining components
n_tps = time_series.shape[0]
images = np.rollaxis(images, axis)
ncomps = images.shape[0]
img_size = np.prod(images.shape[1:])
rarr = images.reshape((ncomps, img_size))
recond = np.dot(time_series, rarr)
recond = recond.reshape((n_tps,) + images.shape[1:])
if axis < 0:
axis = axis + images.ndim
recond = np.rollaxis(recond, 0, axis+1)
return recond
3
Example 21
Project: spectralDNS Source File: shentransform.py
def fft(self, u, fu):
"""Fast Fourier transform of y and z"""
# Intermediate work arrays
Uc_mpi = self.work_arrays[((self.num_processes, self.Np[0], self.Np[1], self.Nf), self.complex, 0)]
Uc_hatT = self.work_arrays[(self.complex_shape_T(), self.complex, 0)]
Uc_hatT = rfft2(u, Uc_hatT, axes=(1,2), threads=self.threads, planner_effort=self.planner_effort['rfft2'])
Uc_mpi[:] = rollaxis(Uc_hatT.reshape(self.Np[0], self.num_processes, self.Np[1], self.Nf), 1)
self.comm.Alltoall([Uc_mpi, self.mpitype], [fu, self.mpitype])
return fu
3
Example 22
def center_and_norm(x, axis=-1):
""" Centers and norms x **in place**
Parameters
-----------
x: ndarray
Array with an axis of observations (statistical units) measured on
random variables.
axis: int, optional
Axis along which the mean and variance are calculated.
"""
x = np.rollaxis(x, axis)
x -= x.mean(axis=0)
x /= x.std(axis=0)
3
Example 23
Project: chainer Source File: test_rollaxis.py
def check_forward(self, x_data):
x = chainer.Variable(x_data)
y = functions.rollaxis(x, self.axis, self.start)
expect = numpy.rollaxis(self.x, self.axis, self.start)
testing.assert_allclose(y.data, expect)
3
Example 24
Project: chainer Source File: test_permutate.py
def check_forward(self, x_data, ind_data):
x = chainer.Variable(x_data)
indices = chainer.Variable(ind_data)
y = functions.permutate(x, indices, axis=self.axis, inv=self.inv)
y_cpu = cuda.to_cpu(y.data)
y_cpu = numpy.rollaxis(y_cpu, axis=self.axis)
x_data = numpy.rollaxis(self.x, axis=self.axis)
for i, ind in enumerate(self.indices):
if self.inv:
numpy.testing.assert_array_equal(y_cpu[ind], x_data[i])
else:
numpy.testing.assert_array_equal(y_cpu[i], x_data[ind])
3
Example 25
Project: eofs Source File: iris.py
def _map_and_dims(pcs, field, mapfunc, *args, **kwargs):
"""
Compute a set of covariance/correlation maps and the resulting
dimensions.
"""
info = {}
for cube in (field, pcs):
info[cube.name()] = _time_coord_info(cube)
cmap_args = [np.rollaxis(cube.data, info[cube.name()][0])
for cube in (pcs, field)]
cmap_args += args
dim_args = [info[cube.name()][1] for cube in (pcs, field)]
cmap = mapfunc(*cmap_args, **kwargs)
dims = covcor_dimensions(*dim_args)
return cmap, dims
1
Example 26
Project: python-skyfield Source File: timelib.py
def __getattr__(self, name):
# Cache of several expensive functions of time.
if name == 'P':
self.P = P = compute_precession(self.tdb)
return P
if name == 'PT':
self.PT = PT = rollaxis(self.P, 1)
return PT
if name == 'N':
self.N = N = compute_nutation(self)
return N
if name == 'NT':
self.NT = NT = rollaxis(self.N, 1)
return NT
if name == 'M':
self.M = M = einsum('ij...,jk...,kl...->il...', self.N, self.P, B)
return M
if name == 'MT':
self.MT = MT = rollaxis(self.M, 1)
return MT
# Conversion between timescales.
if name == 'tai':
self.tai = tai = self.tt - tt_minus_tai
return tai
if name == 'utc':
utc = self._utc_tuple()
utc = array(utc) if self.shape else utc
self.utc = utc = utc
return utc
if name == 'tdb':
tt = self.tt
self.tdb = tdb = tt + tdb_minus_tt(tt) / DAY_S
return tdb
if name == 'ut1':
self.ut1 = ut1 = self.tt - self.delta_t / DAY_S
return ut1
if name == 'delta_t':
table = self.ts.delta_t_table
self.delta_t = delta_t = interpolate_delta_t(table, self.tt)
return delta_t
if name == 'gmst':
self.gmst = gmst = sidereal_time(self)
return gmst
if name == 'gast':
self.gast = gast = self.gmst + earth_tilt(self)[2] / 3600.0
return gast
raise AttributeError('no such attribute %r' % name)
0
Example 27
Project: py-sdm Source File: extract_image_features.py
def _find_working_imread(modes=IMREAD_MODES):
"Finds an image-reading mode that works; returns the name and a function."
if isinstance(modes, str_types):
modes = [modes]
for mode in modes:
try:
if mode.startswith('skimage-'):
from skimage.io import use_plugin, imread
use_plugin(mode[len('skimage-'):])
elif mode == 'cv2':
import cv2
def imread(f):
img = cv2.imread(f)
if img.ndim == 3:
b, g, r = np.rollaxis(img, axis=-1)
return np.dstack([r, g, b])
return img
elif mode == 'matplotlib':
import matplotlib.pyplot as mpl
imread = lambda f: mpl.imread(f)[::-1]
return mode, imread
except ImportError:
pass
else:
raise ImportError("couldn't import any of {}".format(', '.join(modes)))
0
Example 28
Project: Numpy_arraysetops_EP Source File: utility.py
def axis_as_object(arr, axis=-1):
"""cast the given axis of an array to a void object
if the axis to be cast is contiguous, a view is returned, otherwise a copy is made
this is useful for efficiently sorting by the content of an axis, for instance
Parameters
----------
arr : ndarray
array to view as void object type
axis : int
axis to view as a void object type
Returns
-------
ndarray
array with the given axis viewed as a void object
"""
shape = arr.shape
# make axis to be viewed as a void object as contiguous items
arr = np.ascontiguousarray(np.rollaxis(arr, axis, arr.ndim))
# number of bytes in each void object
nbytes = arr.dtype.itemsize * shape[axis]
# void type with the correct number of bytes
voidtype = np.dtype((np.void, nbytes))
# return the view as such, with the reduced shape
return arr.view(voidtype).reshape(np.delete(shape, axis))
0
Example 29
Project: scikit-image Source File: test_adapt_rgb.py
def test_each_channel():
filtered = edges_each(COLOR_IMAGE)
for i, channel in enumerate(np.rollaxis(filtered, axis=-1)):
expected = img_as_float(filters.sobel(COLOR_IMAGE[:, :, i]))
assert_allclose(channel, expected)
0
Example 30
def preprocess(net, img):
#print np.float32(img).shape
return np.float32(np.rollaxis(img, 2)[::-1]) - net.transformer.mean['data']
0
Example 31
Project: Numpy_arraysetops_EP Source File: utility.py
def object_as_axis(arr, dtype, axis=-1):
"""
cast an array of void objects to a typed axis
Parameters
----------
arr : ndarray, [ndim], void
array of type np.void
dtype : numpy dtype object
the output dtype to cast the input array to
axis : int
position to insert the newly formed axis into
Returns
-------
ndarray, [ndim+1], dtype
output array cast to given dtype
"""
# view the void objects as typed elements
arr = arr.view(dtype).reshape(arr.shape + (-1,))
# put the axis in the specified location
return np.rollaxis(arr, -1, axis)
0
Example 32
Project: pims Source File: base_frames.py
def _as_grey(self, as_grey, process_func):
# See skimage.color.colorconv in the scikit-image project.
# As noted there, the weights used in this conversion are calibrated
# for contemporary CRT phosphors. Any alpha channel is ignored."""
if as_grey:
if process_func is not None:
raise ValueError("The as_grey option cannot be used when "
"process_func is specified. Incorpate "
"greyscale conversion in the function "
"passed to process_func.")
shape = self.frame_shape
ndim = len(shape)
# Look for dimensions that look like color channels.
rgb_like = shape.count(3) == 1
rgba_like = shape.count(4) == 1
if ndim == 2:
# The image is already greyscale.
process_func = None
elif ndim == 3 and (rgb_like or rgba_like):
reduced_shape = list(shape)
if rgb_like:
color_axis_size = 3
calibration = [0.2125, 0.7154, 0.0721]
else:
color_axis_size = 4
calibration = [0.2125, 0.7154, 0.0721, 0]
reduced_shape.remove(color_axis_size)
self._im_sz = tuple(reduced_shape)
def convert_to_grey(img):
color_axis = img.shape.index(color_axis_size)
img = np.rollaxis(img, color_axis, 3)
grey = (img * calibration).sum(2)
return grey.astype(img.dtype) # coerce to original dtype
self.process_func = convert_to_grey
else:
raise NotImplementedError("I don't know how to convert an "
"image of shaped {0} to greyscale. "
"Write you own function and pass "
"it using the process_func "
"keyword argument.".format(shape))
0
Example 33
Project: scikit-image Source File: test_adapt_rgb.py
def test_each_channel_with_filter_argument():
filtered = smooth_each(COLOR_IMAGE, SIGMA)
for i, channel in enumerate(np.rollaxis(filtered, axis=-1)):
assert_allclose(channel, smooth(COLOR_IMAGE[:, :, i]))
0
Example 34
Project: pims Source File: display.py
def export_pyav(sequence, filename, rate=30, bitrate=None,
width=None, height=None, codec='mpeg4', format='yuv420p',
autoscale=True):
"""Export a sequence of images as a standard video file using PyAv.
N.B. If the quality and detail are insufficient, increase the
bitrate.
Parameters
----------
sequence : any iterator or array of array-like images
The images should have two dimensions plus an
optional third dimensions representing color.
filename : string
name of output file
rate : integer
frame rate of output file, 30 by default
bitrate : integer
Video bitrate is crudely guessed if None is given.
width : integer
By default, set the width of the images.
height : integer
By default, set the height of the images. If width is specified
and height is not, the height is autoscaled to maintain the aspect
ratio.
codec : string
a valid video encoding, 'mpeg4' by default
format: string
Video stream format, 'yuv420p' by default.
autoscale : boolean
Linearly rescale the brightness to use the full gamut of black to
white values. If the datatype of the images is not 'uint8', this must
be set to True, as it is by default.
"""
if av is None:
raise("This feature requires PyAV with FFmpeg or libav installed.")
output = av.open(filename, 'w')
stream = output.add_stream(bytes(codec), rate)
stream.pix_fmt = bytes(format)
ndim = None
for frame_no, img in enumerate(sequence):
if not frame_no:
# Inspect first frame to set up stream.
if bitrate is None:
bitrate = _estimate_bitrate(img.shape, rate)
stream.bit_rate = int(bitrate)
if width is None:
stream.height = img.shape[0]
stream.width = img.shape[1]
else:
stream.width = width
stream.height = (height or
width * img.shape[0] // img.shape[1])
ndim = img.ndim
if ndim == 3:
if img.shape.count(3) != 1:
raise ValueError("Images have the wrong shape.")
# This is a color image. Ensure that the color axis is axis 2.
color_axis = img.shape.index(3)
img = np.rollaxis(img, color_axis, 3)
elif ndim == 2:
# Expand into color to satisfy PyAV's expectation that images
# be in color. (Without this, an assert is tripped.)
img = np.repeat(np.expand_dims(img, 2), 3, axis=2)
else:
raise ValueError("Images have the wrong shape.")
# PyAV requires uint8.
if img.dtype is not np.uint8 and (not autoscale):
raise ValueError("Autoscaling must be turned on if the image "
"data type is not uint8. Convert the datatype "
"manually if you want to turn off autoscale.")
if autoscale:
normed = (img - img.min()) / (img.max() - img.min())
img = (255 * normed).astype('uint8')
frame = av.VideoFrame.from_ndarray(np.asarray(img), format=b'bgr24')
packet = stream.encode(frame)
output.mux(packet)
output.close()
0
Example 35
@staticmethod
def profile_line(img, src, dst, axes, linewidth=1,
order=1, mode='constant', cval=0.0):
"""Return the intensity profile of an image measured along a scan line.
Parameters
----------
img : numeric array, shape (M, N[, C])
The image, either grayscale (2D array) or multichannel
(3D array, where the final axis contains the channel
information).
src : 2-tuple of numeric scalar (float or int)
The start point of the scan line.
dst : 2-tuple of numeric scalar (float or int)
The end point of the scan line.
linewidth : int, optional
Width of the scan, perpendicular to the line
order : int in {0, 1, 2, 3, 4, 5}, optional
The order of the spline interpolation to compute image values at
non-integer coordinates. 0 means nearest-neighbor interpolation.
mode : string, one of {'constant', 'nearest', 'reflect', 'wrap'},
optional
How to compute any values falling outside of the image.
cval : float, optional
If `mode` is 'constant', what constant value to use outside the
image.
Returns
-------
return_value : array
The intensity profile along the scan line. The length of the
profile is the ceil of the computed length of the scan line.
Examples
--------
>>> x = np.array([[1, 1, 1, 2, 2, 2]])
>>> img = np.vstack([np.zeros_like(x), x, x, x, np.zeros_like(x)])
>>> img
array([[0, 0, 0, 0, 0, 0],
[1, 1, 1, 2, 2, 2],
[1, 1, 1, 2, 2, 2],
[1, 1, 1, 2, 2, 2],
[0, 0, 0, 0, 0, 0]])
>>> profile_line(img, (2, 1), (2, 4))
array([ 1., 1., 2., 2.])
Notes
-----
The destination point is included in the profile, in contrast to
standard numpy indexing.
"""
import scipy.ndimage as nd
p0 = ((src[0] - axes[0].offset) / axes[0].scale,
(src[1] - axes[1].offset) / axes[1].scale)
p1 = ((dst[0] - axes[0].offset) / axes[0].scale,
(dst[1] - axes[1].offset) / axes[1].scale)
perp_lines = Line2DROI._line_profile_coordinates(p0, p1,
linewidth=linewidth)
if img.ndim > 2:
idx = [ax.index_in_array for ax in axes]
if idx[0] < idx[1]:
img = np.rollaxis(img, idx[0], 0)
img = np.rollaxis(img, idx[1], 1)
else:
img = np.rollaxis(img, idx[1], 0)
img = np.rollaxis(img, idx[0], 0)
orig_shape = img.shape
img = np.reshape(img, orig_shape[0:2] +
(np.product(orig_shape[2:]),))
pixels = [nd.map_coordinates(img[..., i], perp_lines,
order=order, mode=mode, cval=cval)
for i in range(img.shape[2])]
i0 = min(axes[0].index_in_array, axes[1].index_in_array)
pixels = np.transpose(np.asarray(pixels), (1, 2, 0))
intensities = pixels.mean(axis=1)
intensities = np.rollaxis(
np.reshape(intensities,
intensities.shape[0:1] + orig_shape[2:]),
0, i0 + 1)
else:
pixels = nd.map_coordinates(img, perp_lines,
order=order, mode=mode, cval=cval)
intensities = pixels.mean(axis=1)
return intensities
0
Example 36
Project: sharedmem Source File: array.py
def reduce(self, a, axis=0, dtype=None, chunksize=1024 * 1024):
rt = [None]
if axis != 0:
a = numpy.rollaxis(a, axis)
with sharedmem.MapReduce() as pool:
def work(i):
sl = slice(i, i+chunksize)
if len(a[sl]) == 0:
return self.ufunc.identity
else:
return self.ufunc.reduce(a[sl], axis, dtype)
def reduce(r):
if rt[0] is None:
rt[0] = r
elif r is None:
pass
else:
rt[0] = self.ufunc(rt[0], r, dtype=dtype)
pool.map(work,
range(0, len(a), chunksize),
reduce=reduce)
return rt[0]
0
Example 37
Project: spectralDNS Source File: shentransform.py
def forward(self, u, fu, fun, dealias=None):
# Intermediate work arrays
Uc_hat = self.work_arrays[(self.complex_shape(), self.complex, 0, False)]
if self.num_processes == 1:
if not dealias == '3/2-rule':
assert u.shape == self.real_shape()
Uc_hat = rfft2(u, Uc_hat, axes=(1,2), threads=self.threads, planner_effort=self.planner_effort['rfft2'])
fu = fun(Uc_hat, fu)
else:
if not self.dealias_cheb:
Upad_hat = self.work_arrays[(self.complex_shape_padded(), self.complex, 0, False)]
Upad_hat_z = self.work_arrays[((self.N[0], int(self.padsize*self.N[1]), self.Nf), self.complex, 0, False)]
Upad_hat = rfft(u, Upad_hat, axis=2, threads=self.threads, planner_effort=self.planner_effort['rfft'])
Upad_hat_z = SlabShen_R2C.copy_from_padded(Upad_hat, Upad_hat_z, self.N, 2)
Upad_hat_z[:] = fft(Upad_hat_z, axis=1, overwrite_input=True, threads=self.threads, planner_effort=self.planner_effort['fft'])
Uc_hat = SlabShen_R2C.copy_from_padded(Upad_hat_z, Uc_hat, self.N, 1)
fu = fun(Uc_hat/self.padsize**2, fu)
else:
# Intermediate work arrays required for transform
Upad_hat = self.work_arrays[(self.complex_shape_padded_0(), self.complex, 0, False)]
Upad_hat0 = self.work_arrays[(self.complex_shape_padded_0(), self.complex, 1, False)]
Upad_hat2 = self.work_arrays[(self.complex_shape_padded_2(), self.complex, 0, False)]
Upad_hat3 = self.work_arrays[(self.complex_shape_padded_3(), self.complex, 0, False)]
# Do ffts and truncation in the padded y and z directions
Upad_hat3 = rfft(u, Upad_hat3, axis=2, threads=self.threads, planner_effort=self.planner_effort['rfft'])
Upad_hat2 = SlabShen_R2C.copy_from_padded(Upad_hat3, Upad_hat2, self.N, 2)
Upad_hat2[:] = fft(Upad_hat2, axis=1, threads=self.threads, planner_effort=self.planner_effort['fft'])
Upad_hat = SlabShen_R2C.copy_from_padded(Upad_hat2, Upad_hat, self.N, 1)
# Perform fst of data in x-direction
Upad_hat0 = fun(Upad_hat, Upad_hat0)
# Truncate to original complex shape
fu[:] = Upad_hat0[:self.N[0]]/self.padsize**2
return fu
if not dealias == '3/2-rule':
Uc_hatT = self.work_arrays[(self.complex_shape_T(), self.complex, 0, False)]
if self.communication == 'alltoall':
Uc_mpi = Uc_hat.reshape((self.num_processes, self.Np[0], self.Np[1], self.Nf))
Uc_hatT = rfft2(u, Uc_hatT, axes=(1,2), threads=self.threads, planner_effort=self.planner_effort['rfft2'])
Uc_mpi[:] = rollaxis(Uc_hatT.reshape(self.Np[0], self.num_processes, self.Np[1], self.Nf), 1)
self.comm.Alltoall(MPI.IN_PLACE, [Uc_hat, self.mpitype])
elif self.communication == 'Alltoallw':
if len(self._subarraysA) == 0:
self._subarraysA, self._subarraysB, self._counts_displs = self.get_subarrays()
# Do 2 ffts in y-z directions on owned data
Uc_hatT = rfft2(u, Uc_hatT, axes=(1,2), threads=self.threads, planner_effort=self.planner_effort['rfft2'])
self.comm.Alltoallw(
[Uc_hatT, self._counts_displs, self._subarraysB],
[Uc_hat, self._counts_displs, self._subarraysA])
fu = fun(Uc_hat, fu)
else:
Uc_hatT = self.work_arrays[(self.complex_shape_T(), self.complex, 0, False)]
Uc_mpi = Uc_hat.reshape((self.num_processes, self.Np[0], self.Np[1], self.Nf))
if not self.dealias_cheb:
Upad_hatT = self.work_arrays[(self.complex_shape_padded_T(), self.complex, 0, False)]
Upad_hat_z = self.work_arrays[((self.Np[0], int(self.padsize*self.N[1]), self.Nf), self.complex, 0, False)]
Upad_hatT = rfft(u, Upad_hatT, axis=2, threads=self.threads, planner_effort=self.planner_effort['rfft'])
Upad_hat_z = SlabShen_R2C.copy_from_padded(Upad_hatT, Upad_hat_z, self.N, 2)
Upad_hat_z[:] = fft(Upad_hat_z, axis=1, threads=self.threads, planner_effort=self.planner_effort['fft'])
Uc_hatT = SlabShen_R2C.copy_from_padded(Upad_hat_z, Uc_hatT, self.N, 1)
if self.communication == 'alltoall':
Uc_mpi[:] = rollaxis(Uc_hatT.reshape(self.Np[0], self.num_processes, self.Np[1], self.Nf), 1)
self.comm.Alltoall(MPI.IN_PLACE, [Uc_hat, self.mpitype])
elif self.communication == 'Alltoallw':
if len(self._subarraysA) == 0:
self._subarraysA, self._subarraysB, self._counts_displs = self.get_subarrays()
self.comm.Alltoallw(
[Uc_hatT, self._counts_displs, self._subarraysB],
[Uc_hat, self._counts_displs, self._subarraysA])
fu = fun(Uc_hat/self.padsize**2, fu)
else:
assert self.num_processes <= self.N[0]/2, "Number of processors cannot be larger than N[0]/2 for 3/2-rule"
assert u.shape == self.real_shape_padded()
# Intermediate work arrays required for transform
Upad_hat = self.work_arrays[(self.complex_shape_padded_0(), self.complex, 0, False)]
Upad_hat0 = self.work_arrays[(self.complex_shape_padded_0(), self.complex, 1, False)]
Upad_hat1 = self.work_arrays[(self.complex_shape_padded_1(), self.complex, 0, False)]
Upad_hat2 = self.work_arrays[(self.complex_shape_padded_2(), self.complex, 0, False)]
Upad_hat3 = self.work_arrays[(self.complex_shape_padded_3(), self.complex, 0, False)]
# Do ffts and truncation in the padded y and z directions
Upad_hat3 = rfft(u, Upad_hat3, axis=2, threads=self.threads, planner_effort=self.planner_effort['rfft'])
Upad_hat2 = SlabShen_R2C.copy_from_padded(Upad_hat3, Upad_hat2, self.N, 2)
Upad_hat2[:] = fft(Upad_hat2, axis=1, threads=self.threads, planner_effort=self.planner_effort['fft'])
Upad_hat1 = SlabShen_R2C.copy_from_padded(Upad_hat2, Upad_hat1, self.N, 1)
if self.communication == 'alltoall':
# Transpose and commuincate data
U_mpi = Upad_hat.reshape(self.complex_shape_padded_0_I())
U_mpi[:] = rollaxis(Upad_hat1.reshape(self.complex_shape_padded_I()), 1)
self.comm.Alltoall(MPI.IN_PLACE, [Upad_hat, self.mpitype])
elif self.communication == 'Alltoallw':
if len(self._subarraysA_pad) == 0:
self._subarraysA_pad, self._subarraysB_pad, self._counts_displs = self.get_subarrays(padsize=self.padsize)
self.comm.Alltoallw(
[Upad_hat1, self._counts_displs, self._subarraysB_pad],
[Upad_hat, self._counts_displs, self._subarraysA_pad])
# Perform fst of data in x-direction
Upad_hat0 = fun(Upad_hat, Upad_hat0)
# Truncate to original complex shape
fu[:] = Upad_hat0[:self.N[0]]/self.padsize**2
return fu
0
Example 38
def preprocess(net, img):
return np.float32(np.rollaxis(img, 2)[::-1]) - net.transformer.mean['data']
0
Example 39
def rgb2gray(rgb):
r, g, b = np.rollaxis(rgb[..., :3], axis= -1)
return 0.21 * r + 0.71 * g + 0.07 * b
0
Example 40
def logsumexp(arr, axis=0):
"""Computes the sum of arr assuming arr is in the log domain.
Returns log(sum(exp(arr))) while minimizing the possibility of
over/underflow.
Examples
--------
>>> import numpy as np
>>> from sklearn.utils.extmath import logsumexp
>>> a = np.arange(10)
>>> np.log(np.sum(np.exp(a)))
9.4586297444267107
>>> logsumexp(a)
9.4586297444267107
"""
arr = np.rollaxis(arr, axis)
# Use the max to normalize, as with the log this is what accuemulates
# the less errors
vmax = arr.max(axis=0)
out = np.log(np.sum(np.exp(arr - vmax), axis=0))
out += vmax
return out
0
Example 41
@lazyprop
def m1(self):
m = np.rollaxis(self.ubasis.jac_nodal_basis_at(self.upts), 2)
return m.reshape(self.nupts, -1)
0
Example 42
def deltaE_ciede94(lab1, lab2, kH=1, kC=1, kL=1, k1=0.045, k2=0.015):
"""Color difference according to CIEDE 94 standard
Accommodates perceptual non-uniformities through the use of application
specific scale factors (`kH`, `kC`, `kL`, `k1`, and `k2`).
Parameters
----------
lab1 : array_like
reference color (Lab colorspace)
lab2 : array_like
comparison color (Lab colorspace)
kH : float, optional
Hue scale
kC : float, optional
Chroma scale
kL : float, optional
Lightness scale
k1 : float, optional
first scale parameter
k2 : float, optional
second scale parameter
Returns
-------
dE : array_like
color difference between `lab1` and `lab2`
Notes
-----
deltaE_ciede94 is not symmetric with respect to lab1 and lab2. CIEDE94
defines the scales for the lightness, hue, and chroma in terms of the first
color. Consequently, the first color should be regarded as the "reference"
color.
`kL`, `k1`, `k2` depend on the application and default to the values
suggested for graphic arts
========== ============== ==========
Parameter Graphic Arts Textiles
========== ============== ==========
`kL` 1.000 2.000
`k1` 0.045 0.048
`k2` 0.015 0.014
========== ============== ==========
References
----------
.. [1] http://en.wikipedia.org/wiki/Color_difference
.. [2] http://www.brucelindbloom.com/index.html?Eqn_DeltaE_CIE94.html
"""
L1, C1 = np.rollaxis(lab2lch(lab1), -1)[:2]
L2, C2 = np.rollaxis(lab2lch(lab2), -1)[:2]
dL = L1 - L2
dC = C1 - C2
dH2 = get_dH2(lab1, lab2)
SL = 1
SC = 1 + k1 * C1
SH = 1 + k2 * C1
dE2 = (dL / (kL * SL)) ** 2
dE2 += (dC / (kC * SC)) ** 2
dE2 += dH2 / (kH * SH) ** 2
return np.sqrt(dE2)
0
Example 43
Project: PyFR Source File: elements.py
@lazyprop
def _smats_djacs_mpts(self):
# Metric basis with grid point (q<=p) or pseudo grid points (q>p)
mpts = self.basis.mpts
mbasis = self.basis.mbasis
# Dimensions, number of elements and number of mpts
ndims, neles, nmpts = self.ndims, self.neles, self.nmpts
# Physical locations of the pseudo grid points
x = self.ploc_at_np('mpts')
# Jacobian operator at these points
jacop = np.rollaxis(mbasis.jac_nodal_basis_at(mpts), 2)
jacop = jacop.reshape(-1, nmpts)
# Cast as a matrix multiply and apply to eles
jac = np.dot(jacop, x.reshape(nmpts, -1))
# Reshape (nmpts*ndims, neles*ndims) => (nmpts, ndims, neles, ndims)
jac = jac.reshape(nmpts, ndims, ndims, neles)
# Transpose to get (ndims, ndims, nmpts, neles)
jac = jac.transpose(1, 2, 0, 3)
smats = np.empty((ndims, nmpts, ndims, neles))
if ndims == 2:
a, b, c, d = jac[0, 0], jac[1, 0], jac[0, 1], jac[1, 1]
smats[0, :, 0], smats[0, :, 1] = d, -b
smats[1, :, 0], smats[1, :, 1] = -c, a
djacs = a*d - b*c
else:
dtt = []
for dx in jac:
# Compute x cross x_(chi)
tt = np.cross(x, dx, axisa=1, axisb=0, axisc=1)
# Jacobian of x cross x_(chi) at the pseudo grid points
dt = np.dot(jacop, tt.reshape(nmpts, -1))
dt = dt.reshape(nmpts, ndims, ndims, -1).swapaxes(0, 1)
dtt.append(dt)
# Kopriva's invariant form of smats; JSC 26(3), 301-327, Eq. (37)
smats[0] = 0.5*(dtt[2][1] - dtt[1][2])
smats[1] = 0.5*(dtt[0][2] - dtt[2][0])
smats[2] = 0.5*(dtt[1][0] - dtt[0][1])
# Exploit the fact that det(J) = x0 . (x1 ^ x2)
djacs = np.einsum('ij...,ji...->j...', jac[0], smats[0])
return smats.reshape(ndims, nmpts, -1), djacs
0
Example 44
def apply_transform(x, transform_matrix, channel_index=2, fill_mode='nearest', cval=0.):
"""Return transformed images by given transform_matrix from ``transform_matrix_offset_center``.
Parameters
----------
x : numpy array
Batch of images with dimension of 3, [batch_size, row, col, channel].
transform_matrix : numpy array
Transform matrix (offset center), can be generated by ``transform_matrix_offset_center``
channel_index : int
Index of channel, default 2.
fill_mode : string
Method to fill missing pixel, default ‘nearest’, more options ‘constant’, ‘reflect’ or ‘wrap’
- `scipy ndimage affine_transform <https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.ndimage.interpolation.affine_transform.html>`_
cval : scalar, optional
Value used for points outside the boundaries of the input if mode='constant'. Default is 0.0
- `scipy ndimage affine_transform <https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.ndimage.interpolation.affine_transform.html>`_
Examples
--------
- See ``rotation``, ``shift``, ``shear``, ``zoom``.
"""
x = np.rollaxis(x, channel_index, 0)
final_affine_matrix = transform_matrix[:2, :2]
final_offset = transform_matrix[:2, 2]
channel_images = [ndi.interpolation.affine_transform(x_channel, final_affine_matrix,
final_offset, order=0, mode=fill_mode, cval=cval) for x_channel in x]
x = np.stack(channel_images, axis=0)
x = np.rollaxis(x, 0, channel_index+1)
return x
0
Example 45
def __tree_reduce(array, axis=None, op=np.add):
# reduce all axes
if axis is None:
result = array
for axis in range(len(array.shape))[::-1]:
result = __tree_reduce(result, axis=axis, op=op)
return result
assert 0 <= axis and axis < len(array.shape)
result = np.rollaxis(array, axis)
while True:
l = result.shape[0]
# if axis is small enough to neglect rounding errors do a faster numpy-reduce
if l < 10000:
return op.reduce(result, axis=0)
if l % 2 == 0:
result = op(result[0:l/2], result[l/2:])
else:
tmp = result[-1]
result = op(result[0:l/2], result[l/2:-1])
result[0] = op(result[0], tmp)
0
Example 46
def __init__(self, interpolator):
"""
Given an already created :class:`scipy.interpolate.interp1d` instance, return a callable object
which supports linear extrapolation.
.. deprecated :: 1.10
"""
self._interpolator = interpolator
self.x = interpolator.x
# Store the y values given to the interpolator.
self.y = interpolator.y
"""
The y values given to the interpolator object.
.. note:: These are stored with the interpolator.axis last.
"""
# Roll interpolator.axis to the end if scipy no longer does it for us.
if not self.roll_y:
self.y = np.rollaxis(self.y, self._interpolator.axis, self.y.ndim)
0
Example 47
Project: polar2grid Source File: swath.py
def _swath_from_var(var_name, h5_var, tool=None):
"""
given a variable by name, and its hdf5 variable object,
return a normalized numpy masked_array with corrected indexing order
:param var_name: variable name, used to consult internal guidebook
:param h5_var: hdf5 object
:return: numpy masked_array with missing data properly masked and dimensions corrected to
[in-track, cross-track, layer] for 3D variables
"""
if tool is not None:
data = tool(h5_var)
else:
data = h5_var[:]
shape = data.shape
if len(shape) == 3:
# roll the layer axis to the back, eg (101, 84, 60) -> (84, 60, 101)
LOG.debug('rolling %s layer axis to last position' % var_name)
data = np.rollaxis(data, 0, 3)
if 'missing_value' in h5_var.attrs:
mv = float(h5_var.attrs['missing_value'][0])
LOG.debug('missing value for %s is %s' % (var_name, mv))
mask = np.abs(data - mv) < 0.5
data[mask] = np.nan # FUTURE: we'd rather just deal with masked_array properly in output layer
data = np.ma.masked_array(data, mask) # FUTURE: convince scientists to use NaN. also world peace
LOG.debug('min, max = %s, %s' % (np.min(data.flatten()), np.max(data.flatten())))
else:
LOG.warning('no missing_value attribute in %s' % var_name)
data = np.ma.masked_array(data)
return data
0
Example 48
def deltaE_cmc(lab1, lab2, kL=1, kC=1):
"""Color difference from the CMC l:c standard.
This color difference was developed by the Colour Measurement Committee
(CMC) of the Society of Dyers and Colourists (United Kingdom). It is
intended for use in the textile industry.
The scale factors `kL`, `kC` set the weight given to differences in
lightness and chroma relative to differences in hue. The usual values are
``kL=2``, ``kC=1`` for "acceptability" and ``kL=1``, ``kC=1`` for
"imperceptibility". Colors with ``dE > 1`` are "different" for the given
scale factors.
Parameters
----------
lab1 : array_like
reference color (Lab colorspace)
lab2 : array_like
comparison color (Lab colorspace)
Returns
-------
dE : array_like
distance between colors `lab1` and `lab2`
Notes
-----
deltaE_cmc the defines the scales for the lightness, hue, and chroma
in terms of the first color. Consequently
``deltaE_cmc(lab1, lab2) != deltaE_cmc(lab2, lab1)``
References
----------
.. [1] http://en.wikipedia.org/wiki/Color_difference
.. [2] http://www.brucelindbloom.com/index.html?Eqn_DeltaE_CIE94.html
.. [3] F. J. J. Clarke, R. McDonald, and B. Rigg, "Modification to the
JPC79 colour-difference formula," J. Soc. Dyers Colour. 100, 128-132
(1984).
"""
L1, C1, h1 = np.rollaxis(lab2lch(lab1), -1)[:3]
L2, C2, h2 = np.rollaxis(lab2lch(lab2), -1)[:3]
dC = C1 - C2
dL = L1 - L2
dH2 = get_dH2(lab1, lab2)
T = np.where(np.logical_and(np.rad2deg(h1) >= 164, np.rad2deg(h1) <= 345),
0.56 + 0.2 * np.abs(np.cos(h1 + np.deg2rad(168))),
0.36 + 0.4 * np.abs(np.cos(h1 + np.deg2rad(35)))
)
c1_4 = C1 ** 4
F = np.sqrt(c1_4 / (c1_4 + 1900))
SL = np.where(L1 < 16, 0.511, 0.040975 * L1 / (1. + 0.01765 * L1))
SC = 0.638 + 0.0638 * C1 / (1. + 0.0131 * C1)
SH = SC * (F * T + 1 - F)
dE2 = (dL / (kL * SL)) ** 2
dE2 += (dC / (kC * SC)) ** 2
dE2 += dH2 / (SH ** 2)
return np.sqrt(dE2)
0
Example 49
def apply_transform(x, transform_matrix, channel_index=2, fill_mode='nearest', cval=0.):
"""Return transformed images by given transform_matrix from ``transform_matrix_offset_center``.
Parameters
----------
x : numpy array
Batch of images with dimension of 3, [batch_size, row, col, channel].
transform_matrix : numpy array
Transform matrix (offset center), can be generated by ``transform_matrix_offset_center``
channel_index : int
Index of channel, default 2.
fill_mode : string
Method to fill missing pixel, default ‘nearest’, more options ‘constant’, ‘reflect’ or ‘wrap’
- `Scipy ndimage affine_transform <https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.ndimage.interpolation.affine_transform.html>`_
cval : scalar, optional
Value used for points outside the boundaries of the input if mode='constant'. Default is 0.0
- `Scipy ndimage affine_transform <https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.ndimage.interpolation.affine_transform.html>`_
Examples
--------
- See ``rotation``, ``shift``, ``shear``, ``zoom``.
"""
x = np.rollaxis(x, channel_index, 0)
final_affine_matrix = transform_matrix[:2, :2]
final_offset = transform_matrix[:2, 2]
channel_images = [ndi.interpolation.affine_transform(x_channel, final_affine_matrix,
final_offset, order=0, mode=fill_mode, cval=cval) for x_channel in x]
x = np.stack(channel_images, axis=0)
x = np.rollaxis(x, 0, channel_index+1)
return x
0
Example 50
Project: polar2grid Source File: tools.py
def retrieval_read(prod, name, single_level=False):
"""pick out retrieval products from a SND format file
reorganize geographically (120 FOVs -> 2x60FOVs, with 0,3 detectors and 1,2 detectors)
return array of in the form value[level][0..1][0..60]
optional scaling factor is applied as 10**(-sf)
see pdf_ten_980760-eps-iasi-l2.pdf from EUMETSAT
"""
prod = open_product(prod)
lines = prod._records_['mdr']
# read in as [scanline][ifov][level]
sil = array( [prod.get('%s.%s' % (line,name)) for line in lines], dtype=float64 ) # * 10**(-scaling_factor)
if single_level:
nlevels = 1
else:
nlevels = sil.shape[-1]
nscanlines = sil.shape[0]
# convert to [scanline][fieldofregard][detector][level] for ifov_remap use
# ASSUMPTION: this assumes that detectors ordered 0 1 2 3 similar to the L1B products!
# have not confirmed this in the docuementation
sfdl = reshape(sil, (nscanlines,30,4,nlevels))
# convert to [pseudo-scanline][fieldsofview][level] using ifov_remap
# two pseudoscanlines per scanline
# two fields of view per field of regard (four total)
pvl = ifov_remap(sfdl)
# roll level to the front
lpv = rollaxis(pvl, 2)
# assert(lpv is not None)
return lpv