GPU calculation : step 4

After multiple problems regarding the use of Skynet and  the installation of Jacket on it, I finally run my benchmark code. This bar chart sums up the results.

Fig .6 : Calculation time on CPUs and GPU

The horizontal axis legend indicates the size of the final image then the number of projections.

The benchmarks were run on Skynet's GPU and on two different CPUs. The speed-up factors are plotted in this graph.

Fig .5 : Speed-up factor GPU/CPU on two differents CPUs

GPU calculation : step 3

For running the benchmarks, I will use the Jacket plugin for Matlab. Jacket is a Matlab add-on which provide easy-to-use GPU enabled functions.

After running the tests, the results are not those expected, the GPU calculation is slower than the CPU calculation. After discussion with the Jacket customer support, it appears that our GPU (an nVidia Quadro 2000D) is not powerful enough.

Hopefully the University of Oxford just allows me to use one of their supercomputers. This computer called 'Skynet' is equipped with 4 nVidia Tesla C2020 GPUs.

Fig. 5 : Emerald, Oxford's new supercomputer

GPU calculation : step 2

GPU calculation becomes interesting when you have to run the same calculation a large number of times. Our aim is to split the iradon transform of one sinogram in many identical tasks.

If we manage to do that, we will be able to run these task in parallel on the GPU and, we hope, gain calculation time.

In one sinogram, each row corresponds to a given angle in the imaging plane. We will create a sinogram of a phantom for 18 projections angles and compute a reconstructed image using the inverse radon transform.

Fig .2 : Original phantom (left), reconstructed image (right)

Now we will try to do the same, but we will perform one inverse radon transform per projection angle.

Here are the results we get :

Fig .3 : Reconstructed images for different theta values

And after calculating a sum of the images :

Fig .4 : Sum of the images from Fig3

The next step will be to run in parallel the reconstruction for different theta values and see if this as an impact on our calculation time.

GPU calculation : step 1

Computing an image from a sinogram takes us an average of 25s per image. The CPU we use is an IntelCore i7-3960X and our images are 2500x2500.

Fig .1 : Intel Core i7-3960X, state of the art CPU

The bottleneck of our code is the inverse Radon transform. We would like to create a version of this function running on the GPU in order to parallelize the calculation.

Here are two interesting publications dealing with the subject :

• Real time Radon transform
• Computation of Inverse Radon Transform on Graphics Cards

Metal artifact : final step

The last time, we had created a clean sinogram. We will use an inverse radon transform to reconstruct the image from the sinogram.

Fig .7 : Corrected image, without metal replacing (left), with metal replacing (right)

Here we can see that the metal artifact is gone, the image is ready to be analyzed. We have the choice to replace or not the metallic part after removing the artifact.

One thing to do would be to find a way to segment the metallic part without using a threshold. The threshold has to be determined by the user for each image and some datasets contains more than 500 images.

The method we're currently using with quite good results is the mean-shift segmentation.

Finally, a photo of a real case, here the artifacts were messing with the Digital Image Correlation:

Fig .8 : CT scan of an Mg alloy, before correction (left), after correction (right)

It 's the end of a long series of posts, I hope you enjoyed it.

Metal artifact : step 4

Now we have separated the metallic part sinogram, we will suppress it from the original sinogram by just replacing the values corresponding to the metal by NaNs (Not a Number value).

The the missing datas (NaN) will be interpolated using a simple linear interpolation algorithm.

Fig .6 : (a) : original sinogram, (b) : metal sinogram, (c) : metal sinogram removed, (d) : interpolated sinogam

Metal artifacts : step 3

To begin, we will try to segment the metallic object using a simple threshold value. ImageJ is good for that.

Fig .4 : Threshold choice with ImageJ

Here we can see that all the pixels above the value of 119 belongs to the metallic part. We will know create an image containing only the metallic part and calculate its sinogram.

Fig .5 : Metallic part (left), Sinogram of the metallic part (right)

Now we have separated the metallic part sinogram, it will be easy to remove it from the global sinogram.