Updated example

etomo/tests/fourier_reconstruction.py

@@ -1,5 +1,6 @@
 """
-Example of a single image reconstruction
+Example of a single image reconstruction with fourier operators. See
+density_compensation for similar reconstruction with better results.
 """
 import numpy as np
 import matplotlib.pyplot as plt
@@ -11,7 +12,7 @@
 from modopt.opt.proximity import SparseThreshold
 
 from etomo.operators.utils import generate_locations_etomo_2D
-from etomo.operators import NFFT, WaveletPywt, HOTV
+from etomo.operators import NUFFT2, WaveletPywt, HOTV
 from etomo.reconstructors.forwardtomo import TomoReconstructor
 from etomo.operators.fourier.utils import estimate_density_compensation
 
@@ -23,14 +24,12 @@
 # Create radon operator and simulate data
 theta = np.arange(0., 180., 3.)
 kspace_loc = generate_locations_etomo_2D(img_size, theta)
-fourier_op = NFFT(samples=kspace_loc, shape=image.shape, normalized=True)
+fourier_op = NUFFT2(samples=kspace_loc, shape=image.shape, normalized=True)
 data = fourier_op.op(image.data)
 
-TV = HOTV(img_shape=image.shape, order=1)
-wavelet = WaveletPywt(wavelet_name='sym8', nb_scale=3)
-linear_op = TV
+linear_op = WaveletPywt(wavelet_name='bior4.4', nb_scale=3)
 
-regularizer_op = SparseThreshold(linear=Identity(), weights=3e-6)
+regularizer_op = SparseThreshold(linear=Identity(), weights=1e-5)
 reconstructor = TomoReconstructor(
     data_op=fourier_op,
     linear_op=linear_op,
@@ -43,7 +42,7 @@
 x_final, cost, *_ = reconstructor.reconstruct(
     data=data,
     optimization_alg='condatvu',
-    num_iterations=200,
+    num_iterations=300,
     cost_op_kwargs={'cost_interval': 5}
 )
 

etomo/tests/fourier_test.py

@@ -2,7 +2,7 @@
 import numpy as np
 from pysap.data import get_sample_data
 
-from etomo.operators import NUFFT2, NFFT, gpuNUFFT
+from etomo.operators import NUFFT2, gpuNUFFT
 from etomo.operators.utils import generate_locations_etomo_2D
 
 
@@ -28,36 +28,18 @@ def test2D(self):
         kspace_loc = generate_locations_etomo_2D(self.img_size, theta)
 
         fourier_pynufft = NUFFT2(2 * np.pi * kspace_loc, (self.img_size,) * 2)
-        fourier_pynfft = NFFT(kspace_loc, (self.img_size,) * 2)
         fourier_gpu = gpuNUFFT(kspace_loc, (self.img_size,) * 2)
 
         fake_data = np.random.rand(self.img_size, self.img_size)
         fake_adjoint_data = np.random.rand(self.nb_proj * int(
             np.sqrt(2) * self.img_size))
 
-        for ope in [fourier_pynfft, fourier_pynufft, fourier_gpu]:
+        for ope in [fourier_pynufft, fourier_gpu]:
             Fxy = np.sum(ope.op(fake_data) * np.conj(fake_adjoint_data))
             xFty = np.sum(fake_data * np.conj(ope.adj_op(fake_adjoint_data)))
             print(Fxy, xFty)
             self.assertTrue(np.allclose(Fxy, xFty, rtol=1e-6))
 
-    def test_results(self):
-        """
-        Compares results from pynfft and pynufft
-        """
-        image = get_sample_data('2d-mri').data
-        theta = np.arange(0., 180., 180. / self.nb_proj)
-        kspace_loc = generate_locations_etomo_2D(image.shape[0], theta)
-
-        fourier_pynufft = NUFFT2(2 * np.pi * kspace_loc, image.shape)
-        fourier_pynfft = NFFT(kspace_loc, image.shape)
-
-        pynfft_data = fourier_pynfft.op(image)
-        pynufft_data = fourier_pynufft.op(image)
-        print(pynufft_data)
-        print(pynfft_data)
-        self.assertTrue(np.allclose(pynufft_data, pynfft_data, rtol=1e-6))
-
 
 if __name__ == '__main__':
     unittest.main()