Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def test_to_pytorch_function_complex(self):
A = linop.FFT([3])
x = np.array([1 + 1j, 2 + 2j, 3 + 3j], np.complex)
y = np.ones([3], np.complex)
with self.subTest('forward'):
f = pytorch.to_pytorch_function(
A,
input_iscomplex=True,
output_iscomplex=True).apply
x_torch = pytorch.to_pytorch(x)
npt.assert_allclose(f(x_torch).detach().numpy().ravel(),
A(x).view(np.float))
with self.subTest('adjoint'):
y_torch = pytorch.to_pytorch(y)
loss = (f(x_torch) - y_torch).pow(2).sum() / 2
loss.backward()
self.check_linop_pickleable(A)
data_shape = [2, 3, 4]
filt = util.randn([4, 2, 2, 3], dtype=dtype)
A = linop.ConvolveData(
data_shape, filt,
mode=mode,
multi_channel=True)
self.check_linop_linear(A, dtype=dtype, device=device)
self.check_linop_adjoint(A, dtype=dtype, device=device)
self.check_linop_pickleable(A)
data_shape = [2, 3, 4]
filt = util.randn([4, 2, 2, 3], dtype=dtype)
strides = [2, 2]
A = linop.ConvolveData(
data_shape, filt,
mode=mode, strides=strides,
multi_channel=True)
self.check_linop_linear(A, dtype=dtype, device=device)
self.check_linop_adjoint(A, dtype=dtype, device=device)
self.check_linop_pickleable(A)
xp = device.xp
with device:
pulses = xp.zeros((Nc, Nt), xp.complex)
# set up the system matrix
if explicit:
A = rf.linop.PtxSpatialExplicit(sens, coord, dt,
target.shape, b0)
else:
A = sp.mri.linop.Sense(sens, coord, weights=None, tseg=tseg,
ishape=target.shape).H
# handle the Ns * Ns error weighting ROI matrix
W = sp.linop.Multiply(A.oshape, xp.ones(target.shape))
if roi is not None:
W = sp.linop.Multiply(A.oshape, roi)
# apply ROI
A = W * A
# Unconstrained, use conjugate gradient
if st is None:
I = sp.linop.Identity((Nc, coord.shape[0]))
b = A.H * W * target
alg_method = sp.alg.ConjugateGradient(A.H * A + alpha * I,
b, pulses, P=None,
max_iter=max_iter, tol=tol)
# Constrained case, use SDMM
else:
# vectorize target for SDMM
coord[:, 1])))
# add sensitivities to system matrix
AFullExplicit = xp.empty(AExplicit.shape)
for ii in range(nc):
if three_d:
tmp = xp.squeeze(sens[ii, :, :, :]).flatten()
else:
tmp = sens[ii, :, :].flatten()
D = xp.transpose(xp.tile(tmp, [coord.shape[0], 1]))
AFullExplicit = xp.concatenate((AFullExplicit, D * AExplicit),
axis=1)
# remove 1st empty AExplicit entries
AFullExplicit = AFullExplicit[:, coord.shape[0]:]
A = sp.linop.MatMul((coord.shape[0] * nc, 1), AFullExplicit)
# Finally, adjustment of input/output dimensions to be consistent with
# the existing Sense linop operator. [nc x nt] in, [dim x dim] out
Ro = sp.linop.Reshape(ishape=A.oshape, oshape=sens.shape[1:])
Ri = sp.linop.Reshape(ishape=(nc, coord.shape[0]),
oshape=(coord.shape[0] * nc, 1))
A = Ro * A * Ri
A.repr_str = 'pTx spatial explicit'
# output a sigpy linop or a numpy array
if ret_array:
return A.linops[1].mat
else:
return A
return gradh_x
if self.R is None:
gamma_primal = self.lamda + self.mu
else:
gamma_primal = self.mu
else:
gradh = None
gamma_primal = 0
if self.G is None:
proxfc = prox.L2Reg(y.shape, 1, y=y)
gamma_dual = 1
else:
A = linop.Vstack([A, self.G])
proxf1c = prox.L2Reg(self.y.shape, 1, y=y)
proxf2c = prox.Conj(self.proxg)
proxfc = prox.Stack([proxf1c, proxf2c])
proxg = prox.NoOp(self.x.shape)
gamma_dual = 0
if self.tau is None:
if self.sigma is None:
self.sigma = 1
S = linop.Multiply(A.oshape, self.sigma)
AHA = A.H * S * A
max_eig = MaxEig(
AHA,
dtype=self.x.dtype,
device=self.x_device,
def Amult(self, x, A, mu, lam):
M = sp.linop.Multiply(x.shape, np.ones(x.shape) * np.sqrt(1 / mu))
L = sp.linop.Multiply(x.shape, np.ones(x.shape) * np.sqrt(lam))
Y = sp.linop.Vstack((A, M, L))
Ry = sp.linop.Reshape((Y.oshape[0], 1), Y.oshape)
Y = Ry * Y
return Y
def __init__(self, y, L, lamda=0.005,
mode='full', multi_channel=False, device=sp.cpu_device, **kwargs):
self.y = sp.to_device(y, device)
self.L = sp.to_device(L, device)
self.lamda = lamda
self.mode = mode
self.multi_channel = multi_channel
self.device = device
self._get_params()
self.A_R = sp.linop.ConvolveInput(self.R_shape, self.L,
mode=self.mode,
input_multi_channel=True,
output_multi_channel=self.multi_channel)
proxg_R = sp.prox.L1Reg(self.R_shape, lamda)
super().__init__(self.A_R, self.y, proxg=proxg_R, **kwargs)
"""
device = sp.Device(device)
xp = device.xp
with device:
w = xp.ones(coord.shape[:-1], dtype=coord.dtype)
img_shape = sp.estimate_shape(coord)
# Get kernel
x = xp.arange(n, dtype=coord.dtype) / n
kernel = xp.i0(beta * (1 - x**2)**0.5).astype(coord.dtype)
kernel /= kernel.max()
G = sp.linop.Gridding(img_shape, coord, width, kernel)
with tqdm(total=max_iter, desc="PipeMenonDCF",
disable=not show_pbar) as pbar:
for it in range(max_iter):
GHGw = G.H * G * w
w /= xp.abs(GHGw)
resid = xp.abs(GHGw - 1).max().item()
pbar.set_postfix(resid='{0:.2E}'.format(resid))
pbar.update()
return w