summaryrefslogtreecommitdiffstats
path: root/Wrappers/Python/demos
diff options
context:
space:
mode:
Diffstat (limited to 'Wrappers/Python/demos')
-rw-r--r--Wrappers/Python/demos/demo_cpu_inpainters.py192
-rw-r--r--Wrappers/Python/demos/demo_cpu_regularisers.py43
-rw-r--r--Wrappers/Python/demos/demo_gpu_regularisers.py18
3 files changed, 237 insertions, 16 deletions
diff --git a/Wrappers/Python/demos/demo_cpu_inpainters.py b/Wrappers/Python/demos/demo_cpu_inpainters.py
new file mode 100644
index 0000000..3b4191b
--- /dev/null
+++ b/Wrappers/Python/demos/demo_cpu_inpainters.py
@@ -0,0 +1,192 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Demonstration of CPU inpainters
+@authors: Daniil Kazantsev, Edoardo Pasca
+"""
+
+import matplotlib.pyplot as plt
+import numpy as np
+import os
+import timeit
+from scipy import io
+from ccpi.filters.regularisers import NDF_INP, NVM_INP
+from qualitymetrics import rmse
+###############################################################################
+def printParametersToString(pars):
+ txt = r''
+ for key, value in pars.items():
+ if key== 'algorithm' :
+ txt += "{0} = {1}".format(key, value.__name__)
+ elif key == 'input':
+ txt += "{0} = {1}".format(key, np.shape(value))
+ elif key == 'maskData':
+ txt += "{0} = {1}".format(key, np.shape(value))
+ else:
+ txt += "{0} = {1}".format(key, value)
+ txt += '\n'
+ return txt
+###############################################################################
+
+# read sinogram and the mask
+filename = os.path.join(".." , ".." , ".." , "data" ,"SinoInpaint.mat")
+sino = io.loadmat(filename)
+sino_full = sino.get('Sinogram')
+Mask = sino.get('Mask')
+[angles_dim,detectors_dim] = sino_full.shape
+sino_full = sino_full/np.max(sino_full)
+#apply mask to sinogram
+sino_cut = sino_full*(1-Mask)
+#sino_cut_new = np.zeros((angles_dim,detectors_dim),'float32')
+#sino_cut_new = sino_cut.copy(order='c')
+#sino_cut_new[:] = sino_cut[:]
+sino_cut_new = np.ascontiguousarray(sino_cut, dtype=np.float32);
+#mask = np.zeros((angles_dim,detectors_dim),'uint8')
+#mask =Mask.copy(order='c')
+#mask[:] = Mask[:]
+mask = np.ascontiguousarray(Mask, dtype=np.uint8);
+
+plt.figure(1)
+plt.subplot(121)
+plt.imshow(sino_cut_new,vmin=0.0, vmax=1)
+plt.title('Missing Data sinogram')
+plt.subplot(122)
+plt.imshow(mask)
+plt.title('Mask')
+plt.show()
+#%%
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+print ("___Inpainting using linear diffusion (2D)__")
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+
+## plot
+fig = plt.figure(2)
+plt.suptitle('Performance of linear inpainting using the CPU')
+a=fig.add_subplot(1,2,1)
+a.set_title('Missing data sinogram')
+imgplot = plt.imshow(sino_cut_new,cmap="gray")
+
+# set parameters
+pars = {'algorithm' : NDF_INP, \
+ 'input' : sino_cut_new,\
+ 'maskData' : mask,\
+ 'regularisation_parameter':5000,\
+ 'edge_parameter':0,\
+ 'number_of_iterations' :5000 ,\
+ 'time_marching_parameter':0.000075,\
+ 'penalty_type':0
+ }
+
+start_time = timeit.default_timer()
+ndf_inp_linear = NDF_INP(pars['input'],
+ pars['maskData'],
+ pars['regularisation_parameter'],
+ pars['edge_parameter'],
+ pars['number_of_iterations'],
+ pars['time_marching_parameter'],
+ pars['penalty_type'])
+
+rms = rmse(sino_full, ndf_inp_linear)
+pars['rmse'] = rms
+
+txtstr = printParametersToString(pars)
+txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
+print (txtstr)
+a=fig.add_subplot(1,2,2)
+
+# these are matplotlib.patch.Patch properties
+props = dict(boxstyle='round', facecolor='wheat', alpha=0.75)
+# place a text box in upper left in axes coords
+a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
+ verticalalignment='top', bbox=props)
+imgplot = plt.imshow(ndf_inp_linear, cmap="gray")
+plt.title('{}'.format('Linear diffusion inpainting results'))
+#%%
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+print ("_Inpainting using nonlinear diffusion (2D)_")
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+
+## plot
+fig = plt.figure(3)
+plt.suptitle('Performance of nonlinear diffusion inpainting using the CPU')
+a=fig.add_subplot(1,2,1)
+a.set_title('Missing data sinogram')
+imgplot = plt.imshow(sino_cut_new,cmap="gray")
+
+# set parameters
+pars = {'algorithm' : NDF_INP, \
+ 'input' : sino_cut_new,\
+ 'maskData' : mask,\
+ 'regularisation_parameter':80,\
+ 'edge_parameter':0.00009,\
+ 'number_of_iterations' :1500 ,\
+ 'time_marching_parameter':0.000008,\
+ 'penalty_type':1
+ }
+
+start_time = timeit.default_timer()
+ndf_inp_nonlinear = NDF_INP(pars['input'],
+ pars['maskData'],
+ pars['regularisation_parameter'],
+ pars['edge_parameter'],
+ pars['number_of_iterations'],
+ pars['time_marching_parameter'],
+ pars['penalty_type'])
+
+rms = rmse(sino_full, ndf_inp_nonlinear)
+pars['rmse'] = rms
+
+txtstr = printParametersToString(pars)
+txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
+print (txtstr)
+a=fig.add_subplot(1,2,2)
+
+# these are matplotlib.patch.Patch properties
+props = dict(boxstyle='round', facecolor='wheat', alpha=0.75)
+# place a text box in upper left in axes coords
+a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
+ verticalalignment='top', bbox=props)
+imgplot = plt.imshow(ndf_inp_nonlinear, cmap="gray")
+plt.title('{}'.format('Nonlinear diffusion inpainting results'))
+#%%
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+print ("Inpainting using nonlocal vertical marching")
+print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
+
+## plot
+fig = plt.figure(4)
+plt.suptitle('Performance of NVM inpainting using the CPU')
+a=fig.add_subplot(1,2,1)
+a.set_title('Missing data sinogram')
+imgplot = plt.imshow(sino_cut,cmap="gray")
+
+# set parameters
+pars = {'algorithm' : NVM_INP, \
+ 'input' : sino_cut_new,\
+ 'maskData' : mask,\
+ 'SW_increment': 1,\
+ 'number_of_iterations' : 150
+ }
+
+start_time = timeit.default_timer()
+(nvm_inp, mask_upd) = NVM_INP(pars['input'],
+ pars['maskData'],
+ pars['SW_increment'],
+ pars['number_of_iterations'])
+
+rms = rmse(sino_full, nvm_inp)
+pars['rmse'] = rms
+
+txtstr = printParametersToString(pars)
+txtstr += "%s = %.3fs" % ('elapsed time',timeit.default_timer() - start_time)
+print (txtstr)
+a=fig.add_subplot(1,2,2)
+
+# these are matplotlib.patch.Patch properties
+props = dict(boxstyle='round', facecolor='wheat', alpha=0.75)
+# place a text box in upper left in axes coords
+a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
+ verticalalignment='top', bbox=props)
+imgplot = plt.imshow(nvm_inp, cmap="gray")
+plt.title('{}'.format('Nonlocal Vertical Marching inpainting results'))
+#%%
diff --git a/Wrappers/Python/demos/demo_cpu_regularisers.py b/Wrappers/Python/demos/demo_cpu_regularisers.py
index 3567f91..986e3e9 100644
--- a/Wrappers/Python/demos/demo_cpu_regularisers.py
+++ b/Wrappers/Python/demos/demo_cpu_regularisers.py
@@ -44,13 +44,31 @@ u0 = Im + np.random.normal(loc = 0 ,
u_ref = Im + np.random.normal(loc = 0 ,
scale = 0.01 * Im ,
size = np.shape(Im))
-
+(N,M) = np.shape(u0)
# map the u0 u0->u0>0
# f = np.frompyfunc(lambda x: 0 if x < 0 else x, 1,1)
u0 = u0.astype('float32')
u_ref = u_ref.astype('float32')
-
+# change dims to check that modules work with non-squared images
+"""
+M = M-100
+u_ref2 = np.zeros([N,M],dtype='float32')
+u_ref2[:,0:M] = u_ref[:,0:M]
+u_ref = u_ref2
+del u_ref2
+
+u02 = np.zeros([N,M],dtype='float32')
+u02[:,0:M] = u0[:,0:M]
+u0 = u02
+del u02
+
+Im2 = np.zeros([N,M],dtype='float32')
+Im2[:,0:M] = Im[:,0:M]
+Im = Im2
+del Im2
+"""
+#%%
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("_______________ROF-TV (2D)_________________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
@@ -288,7 +306,6 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
imgplot = plt.imshow(fgp_dtv_cpu, cmap="gray")
plt.title('{}'.format('CPU results'))
-
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("__________Total nuclear Variation__________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
@@ -301,9 +318,8 @@ a.set_title('Noisy Image')
imgplot = plt.imshow(u0,cmap="gray")
channelsNo = 5
-N = 512
-noisyVol = np.zeros((channelsNo,N,N),dtype='float32')
-idealVol = np.zeros((channelsNo,N,N),dtype='float32')
+noisyVol = np.zeros((channelsNo,N,M),dtype='float32')
+idealVol = np.zeros((channelsNo,N,M),dtype='float32')
for i in range (channelsNo):
noisyVol[i,:,:] = Im + np.random.normal(loc = 0 , scale = perc * Im , size = np.shape(Im))
@@ -344,25 +360,19 @@ plt.title('{}'.format('CPU results'))
# Uncomment to test 3D regularisation performance
#%%
"""
-N = 512
slices = 20
-
-filename = os.path.join(".." , ".." , ".." , "data" ,"lena_gray_512.tif")
-Im = plt.imread(filename)
-Im = np.asarray(Im, dtype='float32')
-
-Im = Im/255
perc = 0.05
-noisyVol = np.zeros((slices,N,N),dtype='float32')
-noisyRef = np.zeros((slices,N,N),dtype='float32')
-idealVol = np.zeros((slices,N,N),dtype='float32')
+noisyVol = np.zeros((slices,N,M),dtype='float32')
+noisyRef = np.zeros((slices,N,M),dtype='float32')
+idealVol = np.zeros((slices,N,M),dtype='float32')
for i in range (slices):
noisyVol[i,:,:] = Im + np.random.normal(loc = 0 , scale = perc * Im , size = np.shape(Im))
noisyRef[i,:,:] = Im + np.random.normal(loc = 0 , scale = 0.01 * Im , size = np.shape(Im))
idealVol[i,:,:] = Im
+
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("_______________ROF-TV (3D)_________________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
@@ -403,6 +413,7 @@ a.text(0.15, 0.25, txtstr, transform=a.transAxes, fontsize=14,
imgplot = plt.imshow(rof_cpu3D[10,:,:], cmap="gray")
plt.title('{}'.format('Recovered volume on the CPU using ROF-TV'))
+
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("_______________FGP-TV (3D)__________________")
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
diff --git a/Wrappers/Python/demos/demo_gpu_regularisers.py b/Wrappers/Python/demos/demo_gpu_regularisers.py
index b873700..f3ed50c 100644
--- a/Wrappers/Python/demos/demo_gpu_regularisers.py
+++ b/Wrappers/Python/demos/demo_gpu_regularisers.py
@@ -44,10 +44,28 @@ u0 = Im + np.random.normal(loc = 0 ,
u_ref = Im + np.random.normal(loc = 0 ,
scale = 0.01 * Im ,
size = np.shape(Im))
+(N,M) = np.shape(u0)
# map the u0 u0->u0>0
# f = np.frompyfunc(lambda x: 0 if x < 0 else x, 1,1)
u0 = u0.astype('float32')
u_ref = u_ref.astype('float32')
+"""
+M = M-100
+u_ref2 = np.zeros([N,M],dtype='float32')
+u_ref2[:,0:M] = u_ref[:,0:M]
+u_ref = u_ref2
+del u_ref2
+
+u02 = np.zeros([N,M],dtype='float32')
+u02[:,0:M] = u0[:,0:M]
+u0 = u02
+del u02
+
+Im2 = np.zeros([N,M],dtype='float32')
+Im2[:,0:M] = Im[:,0:M]
+Im = Im2
+del Im2
+"""
print ("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
print ("____________ROF-TV regulariser_____________")