Merge pull request #94 from dkazanc/methodscupy

Fixing backward incompatibility of regul toolkit

Demos/Python/DemoFISTA_2D.py

@@ -1,9 +1,7 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
-"""
-Script to generate 2D analytical phantoms and their sinograms with added noise 
+"""Script to generate 2D analytical phantoms and their sinograms with added noise
 and then reconstruct using the regularised FISTA algorithm.
-
 """
 import numpy as np
 import matplotlib.pyplot as plt
@@ -224,59 +222,3 @@
 print("RMSE for FISTA-OS is {}".format(RMSE_FISTA_os))
 print("RMSE for regularised FISTA-OS is {}".format(RMSE_FISTA_os_reg))
 # %%
-print("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
-print("Reconstructing with FISTA-KL-OS method")
-print("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
-from tomobar.methodsIR import RecToolsIR
-
-# set parameters and initiate a class object
-Rectools = RecToolsIR(
-    DetectorsDimH=P,  # Horizontal detector dimension
-    DetectorsDimV=None,  # Vertical detector dimension (3D case)
-    CenterRotOffset=None,  # Center of Rotation scalar
-    AnglesVec=angles_rad,  # A vector of projection angles in radians
-    ObjSize=N_size,  # Reconstructed object dimensions (scalar)
-    datafidelity="KL",  # Data fidelity, choose from LS, KL, PWLS
-    device_projector="gpu",
-)
-
-# prepare dictionaries with parameters:
-_data_ = {"projection_norm_data": noisy_sino, "OS_number": 10}  # data dictionary
-lc = Rectools.powermethod(
-    _data_
-)  # calculate Lipschitz constant (run once to initialise)
-
-# Run FISTA-OS reconstrucion algorithm without regularisation
-_algorithm_ = {"iterations": 50, "lipschitz_const": lc * 0.3}
-RecFISTA_os_kl = Rectools.FISTA(_data_, _algorithm_)
-
-# adding regularisation
-_regularisation_ = {
-    "method": "PD_TV",
-    "regul_param": 0.00003,
-    "iterations": 80,
-    "device_regulariser": "gpu",
-}
-
-# adding regularisation using the CCPi regularisation toolkit
-RecFISTA_os_kl_reg = Rectools.FISTA(_data_, _algorithm_, _regularisation_)
-
-plt.figure()
-plt.subplot(121)
-plt.imshow(RecFISTA_os_kl, vmin=0, vmax=1, cmap="gray")
-plt.colorbar(ticks=[0, 0.5, 1], orientation="vertical")
-plt.title("FISTA-KL-OS reconstruction")
-plt.subplot(122)
-plt.imshow(RecFISTA_os_kl_reg, vmin=0, vmax=1, cmap="gray")
-plt.colorbar(ticks=[0, 0.5, 1], orientation="vertical")
-plt.title("Regularised FISTA-KL-OS reconstruction")
-plt.show()
-
-# calculate errors
-Qtools = QualityTools(phantom_2D[indicesROI], RecFISTA_os_kl[indicesROI])
-RMSE_FISTA_os = Qtools.rmse()
-Qtools = QualityTools(phantom_2D[indicesROI], RecFISTA_os_kl_reg[indicesROI])
-RMSE_FISTA_os_reg = Qtools.rmse()
-print("RMSE for FISTA-KL-OS is {}".format(RMSE_FISTA_os))
-print("RMSE for regularised FISTA-KL-OS is {}".format(RMSE_FISTA_os_reg))
-# %%

Demos/Python/DemoFISTA_NLTV_2D.py

@@ -1,12 +1,9 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
-"""
-Script to generate 2D analytical phantoms and their sinograms with added noise 
+"""Script to generate 2D analytical phantoms and their sinograms with added noise
 and then reconstruct using Non-local Total variation (NLTV) regularised FISTA algorithm.
-
 NLTV method is quite different to the generic structure of other regularisers, hence
 a separate implementation
-
 """
 import numpy as np
 import timeit
@@ -114,13 +111,12 @@
     pars["patchwindow"],
     pars["neighbours"],
     pars["edge_parameter"],
-    "gpu",
 )
-"""
-plt.figure()
-plt.imshow(Weights[0,:,:], vmin=0, vmax=1, cmap="gray")
-plt.colorbar(ticks=[0, 0.5, 1], orientation='vertical')
-"""
+
+# plt.figure()
+# plt.imshow(Weights[0,:,:], vmin=0, vmax=1, cmap="gray")
+# plt.colorbar(ticks=[0, 0.5, 1], orientation='vertical')
+
 # %%
 print("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
 print("%%%%%%%%%%%Reconstructing with FISTA-OS method%%%%%%%%%%%%%%")
@@ -155,10 +151,11 @@
 _regularisation_ = {
     "method": "NLTV",
     "regul_param": 0.0025,
-    "iterations": 5,
     "NLTV_H_i": H_i,
     "NLTV_H_j": H_j,
     "NLTV_Weights": Weights,
+    "NumNeighb": pars["neighbours"],
+    "IterNumb": 5,
     "device_regulariser": "gpu",
 }
 

conda-recipe/environment/environment.yml

@@ -3,7 +3,6 @@ channels:
   - conda-forge
   - httomo
 dependencies:
-  - httomo::ccpi-regularisation-cupy
   - httomo::ccpi-regulariser
   - conda-forge::astra-toolbox
   - conda-forge::cupy

docs/source/howto/installation.rst

@@ -2,9 +2,8 @@
 
 Installation Guide
 ------------------
-ToMoBAR is a Python package with several :ref:`ref_dependencies` to ensure its full functionality.
-It mostly relies on the GPU-enabled computations and therefore we suggest having a decent NVIDIA
-graphics card to support it.
+ToMoBAR is a Python package with several :ref:`ref_dependencies`. To ensure its full functionality it is recommended to install them.
+It mostly relies on the GPU-enabled computations and therefore we suggest using a decent NVIDIA graphics card to support it.
 
 .. _ref_python:
 
@@ -18,24 +17,31 @@ Install ToMoBAR as a pre-built conda Python package:
 
    $ conda install -c httomo tomobar
 
-or install with the main dependencies into a new environment:
+or install with the dependencies into a new environment:
 
 .. code-block:: console
 
    $ conda install -c httomo -c conda-forge tomophantom tomobar astra-toolbox ccpi-regulariser pypwt
 
+one can also try this installation, especially for other than Linux OSs:
+
+.. code-block:: console
+
+   $ conda install -c httomo -c ccpi -c conda-forge tomophantom tomobar astra-toolbox ccpi-regulariser
+
 Install ToMoBAR from PyPi:
 +++++++++++++++++++++++++++
 One can install ToMoBAR from PyPi, however, not all dependencies might be at PyPi yet.
 
 .. code-block:: console
 
    $ pip install tomobar
+   $ pip install ccpi-regularisation-cupy
 
-Create conda environment:
+Using conda environment:
 +++++++++++++++++++++++++
-Sometimes the one-liner above doesn't work. In that case we suggest creating a new conda environment,
-and **pip** install ToMoBAR into it as shown bellow.
+One can also create a new conda environment by using environment yaml file,
+and then **pip** install ToMoBAR into it.
 
 .. code-block:: console
 

tests/test_RecToolsIR.py

@@ -5,7 +5,7 @@
 from tomobar.methodsIR import RecToolsIR
 from tomobar.supp.suppTools import normaliser
 
-eps = 1e-06
+eps = 1e-05
 
 
 def test_SIRT2D(data, angles):
@@ -607,8 +607,8 @@ def test_ADMM3D_reg(data, angles):
 
     Iter_rec = RecTools.ADMM(_data_, _algorithm_, _regularisation_)
 
-    assert_allclose(np.min(Iter_rec), -0.009140163, rtol=eps)
-    assert_allclose(np.max(Iter_rec), 0.020371437, rtol=eps)
+    assert_allclose(np.min(Iter_rec), -0.008855395, rtol=0, atol=eps)
+    assert_allclose(np.max(Iter_rec), 0.020371437, rtol=0, atol=eps)
     assert Iter_rec.dtype == np.float32
     assert Iter_rec.shape == (128, 160, 160)
 
@@ -684,8 +684,8 @@ def test_FISTA3D(data, angles):
 
     Iter_rec = RecTools.FISTA(_data_, _algorithm_)
 
-    assert_allclose(np.min(Iter_rec), -0.0021881335, rtol=eps)
-    assert_allclose(np.max(Iter_rec), 0.024684845, rtol=eps)
+    assert_allclose(np.min(Iter_rec), -0.0021881335, rtol=0, atol=eps)
+    assert_allclose(np.max(Iter_rec), 0.024684845, rtol=0, atol=eps)
     assert Iter_rec.dtype == np.float32
     assert Iter_rec.shape == (128, 160, 160)
 
@@ -751,8 +751,8 @@ def test_FISTA_OS_3D(data, angles):
 
     Iter_rec = RecTools.FISTA(_data_, _algorithm_)
 
-    assert_allclose(np.min(Iter_rec), -0.008425578, rtol=eps)
-    assert_allclose(np.max(Iter_rec), 0.032162726, rtol=eps)
+    assert_allclose(np.min(Iter_rec), -0.008425578, rtol=0, atol=eps)
+    assert_allclose(np.max(Iter_rec), 0.032162726, rtol=0, atol=eps)
     assert Iter_rec.dtype == np.float32
     assert Iter_rec.shape == (128, 160, 160)
 
@@ -790,7 +790,7 @@ def test_FISTA_OS_regul_3D(data, angles):
 
     Iter_rec = RecTools.FISTA(_data_, _algorithm_, _regularisation_)
 
-    assert_allclose(np.min(Iter_rec), -0.000175, rtol=1e-04)
-    assert_allclose(np.max(Iter_rec), 0.021823, rtol=1e-04)
+    assert_allclose(np.min(Iter_rec), -0.000174, rtol=0, atol=eps)
+    assert_allclose(np.max(Iter_rec), 0.021823, rtol=0, atol=eps)
     assert Iter_rec.dtype == np.float32
     assert Iter_rec.shape == (128, 160, 160)

tomobar/astra_wrappers/astra_base.py

@@ -224,7 +224,7 @@ def _set_gpu_projection2d_parallel_geometry(self):
         self.A_optomo = astra.OpTomo(self.proj_id)
 
     def _set_gpu_projection3d_parallel_geometry(self):
-        """the classical 3D projection geometry (cpu)"""
+        """the classical 3D projection geometry"""
         vectors = _vec_geom_init3D(self.angles_vec, 1.0, 1.0, self.centre_of_rotation)
         self.proj_geom = astra.create_proj_geom(
             "parallel3d_vec", self.detectors_y, self.detectors_x, vectors

tomobar/methodsDIR_CuPy.py

@@ -130,7 +130,7 @@ def FBP(self, data: xp.ndarray, **kwargs) -> xp.ndarray:
         # filter the data on the GPU and keep the result there
         data = _filtersinc3D_cupy(data, cutoff=cutoff_freq)
         data = xp.ascontiguousarray(xp.swapaxes(data, 0, 1))
-        xp._default_memory_pool.free_all_blocks() # free everything related to the filtering before starting Astra
+        xp._default_memory_pool.free_all_blocks()  # free everything related to the filtering before starting Astra
         reconstruction = self.Atools._backprojCuPy(data)  # 3d backprojecting
         cache = xp.fft.config.get_plan_cache()
         cache.clear()  # flush FFT cache here before performing ifft to save the memory