diff options
Diffstat (limited to 'Wrappers/Python')
3 files changed, 25 insertions, 18 deletions
diff --git a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Gaussian.py b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Gaussian.py index cba5bcb..9d00ee1 100644 --- a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Gaussian.py +++ b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Gaussian.py @@ -54,7 +54,7 @@ import os, sys loader = TestData(data_dir=os.path.join(sys.prefix, 'share','ccpi')) # Load Data -N = 256 +N = 200 M = 300 @@ -65,7 +65,7 @@ M = 300 # TestData.CAMERA 2D # TestData.RESOLUTION_CHART 2D # TestData.SIMPLE_PHANTOM_2D 2D -data = loader.load(TestData.PEPPERS, size=(N,M), scale=(0,1)) +data = loader.load(TestData.BOAT, size=(N,M), scale=(0,1)) ig = data.geometry ag = ig @@ -212,9 +212,9 @@ if cvx_not_installable: plt.colorbar() plt.show() - plt.plot(np.linspace(0,N,N), pdhg.get_output().as_array()[int(N/2),:], label = 'PDHG') - plt.plot(np.linspace(0,N,N), u.value[int(N/2),:], label = 'CVX') - plt.plot(np.linspace(0,N,N), data.as_array()[int(N/2),:], label = 'Truth') + plt.plot(np.linspace(0,N,M), pdhg.get_output().as_array()[int(N/2),:], label = 'PDHG') + plt.plot(np.linspace(0,N,M), u.value[int(N/2),:], label = 'CVX') + plt.plot(np.linspace(0,N,M), data.as_array()[int(N/2),:], label = 'Truth') plt.legend() plt.title('Middle Line Profiles') diff --git a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Poisson.py b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Poisson.py index 4d53635..1d887c1 100644 --- a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Poisson.py +++ b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Denoising_Poisson.py @@ -49,7 +49,7 @@ import matplotlib.pyplot as plt from ccpi.optimisation.algorithms import PDHG from ccpi.optimisation.operators import BlockOperator, Identity, Gradient -from ccpi.optimisation.functions import ZeroFunction, KullbackLeibler, \ +from ccpi.optimisation.functions import ZeroFunction, KullbackLeibler, IndicatorBox, \ MixedL21Norm, BlockFunction from skimage.util import random_noise @@ -70,7 +70,7 @@ noisy_data = ImageData(n1) # Show Ground Truth and Noisy Data -plt.figure(figsize=(15,15)) +plt.figure(figsize=(10,10)) plt.subplot(2,1,1) plt.imshow(data.as_array()) plt.title('Ground Truth') @@ -81,7 +81,6 @@ plt.title('Noisy Data') plt.colorbar() plt.show() -#%% # Regularisation Parameter alpha = 2 @@ -101,9 +100,8 @@ if method == '0': f1 = alpha * MixedL21Norm() f2 = KullbackLeibler(noisy_data) - f = BlockFunction(f1, f2) - - g = ZeroFunction() + f = BlockFunction(f1, f2) + g = IndicatorBox(lower=0) else: @@ -124,7 +122,7 @@ tau = 1/(sigma*normK**2) pdhg = PDHG(f=f,g=g,operator=operator, tau=tau, sigma=sigma) pdhg.max_iteration = 3000 pdhg.update_objective_interval = 200 -pdhg.run(3000, verbose=False) +pdhg.run(3000, verbose=True) plt.figure(figsize=(15,15)) @@ -165,16 +163,17 @@ if cvx_not_installable: ##Construct problem u1 = Variable(ig.shape) q = Variable() + w = Variable() DY = SparseFiniteDiff(ig, direction=0, bnd_cond='Neumann') DX = SparseFiniteDiff(ig, direction=1, bnd_cond='Neumann') # Define Total Variation as a regulariser - regulariser = alpha * sum(norm(vstack([DX.matrix() * vec(u1), DY.matrix() * vec(u1)]), 2, axis = 0)) + regulariser = alpha * sum(norm(vstack([DX.matrix() * vec(u1), DY.matrix() * vec(u1)]), 2, axis = 0)) + fidelity = sum(kl_div(noisy_data.as_array(), u1)) + + constraints = [q>=fidelity, u1>=0] - fidelity = sum( u1 - multiply(noisy_data.as_array(), log(u1)) ) - constraints = [q>= fidelity, u1>=0] - solver = ECOS obj = Minimize( regulariser + q) prob = Problem(obj, constraints) diff --git a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Tomo2D_poisson.py b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Tomo2D_poisson.py index 830ea00..8f39f86 100644 --- a/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Tomo2D_poisson.py +++ b/Wrappers/Python/demos/PDHG_examples/PDHG_TV_Tomo2D_poisson.py @@ -27,8 +27,8 @@ import matplotlib.pyplot as plt from ccpi.optimisation.algorithms import PDHG -from ccpi.optimisation.operators import BlockOperator, Gradient, Identity -from ccpi.optimisation.functions import ZeroFunction, KullbackLeibler, \ +from ccpi.optimisation.operators import BlockOperator, Gradient +from ccpi.optimisation.functions import KullbackLeibler, \ MixedL21Norm, BlockFunction, IndicatorBox from ccpi.astra.ops import AstraProjectorSimple @@ -68,6 +68,14 @@ detectors = N angles = np.linspace(0, np.pi, N) ag = AcquisitionGeometry('parallel','2D',angles, detectors) + +device = input('Available device: GPU==1 / CPU==0 ') + +if device=='1': + dev = 'gpu' +else: + dev = 'cpu' + Aop = AstraProjectorSimple(ig, ag, 'cpu') sin = Aop.direct(data) |