@@ -0,0 +1,1657 @@
+import odl%0Aimport numpy as np%0Afrom odl.contrib.pyshearlab import PyShearlabOperator%0Aimport scipy.misc%0A%0Aimg = scipy.misc.imread(%0A   '/home/hkohr/Downloads/Yv2Ki.png')%0Aimg = np.sum(img, axis=-1)%0Aimg = img%5B:, :356%5D%5B::2, ::2%5D%0A%0Aspace = odl.uniform_discr(%5B-1, -1%5D, %5B1, 1%5D, img.shape)%0Aimg = space.element(img) / np.max(img)%0Aphantom = odl.phantom.derenzo_sources(space) + img%0A%0Ageom = odl.tomo.parallel_beam_geometry(space)%0Aray_trafo = odl.tomo.RayTransform(space, geom,%0A                                  impl='astra_cpu')%0Adata = ray_trafo(phantom)%0Anoise = odl.phantom.white_noise(ray_trafo.range) * 0.01%0Anoisy_data = data + noise%0A%0Ashear_op = PyShearlabOperator(space, scales=2)%0Awave_op = odl.trafos.WaveletTransform(space, 'haar',%0A                                      nlevels=2)%0Awave_op = wave_op / 89%0A%0Asol_space = odl.ProductSpace(space, 2)%0A%0Al1norm_wave = odl.solvers.L1Norm(wave_op.range)%0Al1norm_shear = odl.solvers.L1Norm(shear_op.range)%0Adata_matching = odl.solvers.L2NormSquared(ray_trafo.range)%0Adata_matching = data_matching.translated(noisy_data)%0A%0Af = odl.solvers.ZeroFunctional(sol_space)%0Apenalizer = odl.solvers.SeparableSum(100 * l1norm_wave,%0A                                     l1norm_shear)%0Asum_op = odl.ReductionOperator(%0A    odl.IdentityOperator(space), 2)%0Acoeff_op = odl.DiagonalOperator(wave_op, shear_op)%0A%0Ag = %5Bdata_matching, penalizer%5D%0AL = %5Bray_trafo * sum_op, coeff_op%5D%0A%0Atau = 1%0Aopnorms = %5Bodl.power_method_opnorm(op) for op in L%5D%0Asigma = %5B1 / opnorm ** 2 for opnorm in opnorms%5D%0A%0Acallback = odl.solvers.CallbackShow()%0A%0Aimages = sol_space.zero()%0Aodl.solvers.douglas_rachford_pd(%0A    images, f, g, L, tau, sigma, niter=100,%0A    callback=callback)%0A

add separation example to shearlets (#2)