Function File: filtered = rho_filter (proj, type, scaling)

Filters the parallel ray projections in the columns of proj, according to the filter type chosen by type. type can be chosen from

  • ’none’
  • ’Ram-Lak’ (default)
  • ’Shepp-Logan’
  • ’Cosine’
  • ’Hann’
  • ’Hamming’

If given, scaling determines the proportion of frequencies below the nyquist frequency that should be passed by the filter. The window function is compressed accordingly, to avoid an abrupt truncation of the frequency response.

Function File: [filtered, filter] = rho_filter (…)

This form also returns the frequency response of the filter in the vector filter.

Performs rho filtering on the parallel ray projections provided.

Rho filtering is performed as part of the filtered back-projection method of CT image reconstruction. It is the filtered part of the name. The simplest rho filter is the Ramachadran-Lakshminarayanan (Ram-Lak), which is simply |rho|, where rho is the radial component of spatial frequency. However, this can cause unwanted amplification of noise, which is what the other types attempt to minimise, by introducing roll-off into the response. The Hann and Hamming filters multiply the standard response by a Hann or Hamming window, respectively. The cosine filter is the standard response multiplied by a cosine shape, and the Shepp-Logan filter multiplies the response with a sinc shape. The ’none’ filter performs no filtering, and is included for completeness and to enable incorporating this function easily into scripts or functions that may offer the ability to choose to apply no filtering.

This function is designed to be used by the function iradon, but has been exposed to facilitate custom inverse radon transforms and to more clearly break down the process for educational purposes. The operations

    filtered = rho_filter (proj);
    reconstruction = iradon (filtered, 1, 'linear', 'none');

are exactly equivalent to

    reconstruction = iradon (proj, 1, 'linear', 'Ram-Lak');

Usage example:

  P = phantom ();
  projections = radon (P);
  filtered_projections = rho_filter (projections, 'Hamming');
  reconstruction = iradon (filtered_projections, 1, 'linear', 'none');
  figure, imshow (reconstruction, [])

Demonstration 1

The following code

 P = phantom ();
 projections = radon (P);
 filtered_projections = rho_filter (projections, 'Hamming');
 reconstruction = iradon (filtered_projections, 1, 'linear', 'none');
 figure, imshow (reconstruction, [])

Produces the following figure

Figure 1

Package: image