Visualizing 32-bit Integer PGMs (eg. from FTL) or other ASCII Image Formats

In this post I will post some Python code for visualising ASCII image formats such as PGMs (esp. from the Finite Transform Library (FTL)). These formats are commonly used because of their simplicity in implementation.

EDIT: The sMILX Viewer (v1.0 Alpha and above), part of the SMILI project for scientific visualisation now supports these PGMs. Just drag and drop the file into the viewer window to visualise them!

The ASCII PGM format for example has a simple header made up of a format string (made of two characters) such as 'P2' on the first line, followed by a comment line, the dimensions of the image (width and height) and the bit depth. Finally the data is next. An example is given below:

P2

# Generated PGM.

101 101

255

184 180 188 199 202 203 200 195 195 202 198 186 181 156 ........

Since most PGM and image viewers always expect a bit depth of 8-bit, the result isn't always shown correctly. I have written a simple Python module to do this type of reading and also utilise 32-bit PNGs as well. The full module can be found in this gist. The section for PGMs is as follows:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32

def readPGM(name): ''' Read a PGM image, where the PGM format is that of the FTL library for NTTs and FRTs. This is a slightly modified format where values are not limited to 8-bit values. Returns array of image and bit depth (image, depth) ''' inFile = open(name,"r") #read header formatLine = inFile.readline() commentLine = inFile.readline() if not "P2" in formatLine: print "Error: PGM not in correct format (P2)" print "Comment:", commentLine width, height = [int(x) for x in inFile.readline().split()] # read dimensions print "PGM Size:", width, "x", height bitDepth = [int(x) for x in inFile.readline().split()] # read bit Depth imageList = [] for line in inFile: # read remaining lines valueList = [int(x) for x in line.split()] #read integers on each line for value in valueList: imageList.append(value) #append as 1D list # print imageList #store as array image = np.array(imageList).reshape(height, width) return image, bitDepth

Then to load results from FTL, simply use a script as:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28

# -*- coding: utf-8 -*- """ Test Read PGM member Created on Sun Jan 25 21:52:38 2015 @author: shakes """ import imageio image, depth = imageio.readPGM("lena.pgm") frtSpace, depth = imageio.readPGM("frtSpace.pgm") #Plot import matplotlib.pyplot as plt #plot images fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(16, 5)) plt.gray() ax[0].imshow(image) ax[0].axis('off') ax[0].set_title('Image') ax[1].imshow(frtSpace) ax[1].axis('off') ax[1].set_title('FRT of Lena') plt.show()

The result is a Matplotlib plot of the images, both of Lena and her FRT space:

You can adapt it for your image format by just changing what/how the header is read. Feel free to fork the gist as I will put the latest version there. This module is part of a pure Python version of FTL coming soon, so watch this space..... better yet... Follow!

PS: A FTL plugin for my scientific visualisation software SMILI is on the way too. More soon in the next version.

UPDATE: Scientific visualisation (open-source) software SMILI supports 32-bit PGMs natively now via the FTL plugin. Get the binaries from GitHub or Sourceforge.

Cheers Shakes - L3mming

Robust Digital Image Reconstruction Example

In this post, we discuss how to employ the digital image reconstruction technique of Chandra et al. (2014):

Robust digital image reconstruction via the discrete Fourier slice theorem

S Chandra, N Normand, A Kingston, JP Guedon, I Svalbe

IEEE Sig. Proc. Lett. (2014)

using the FTL (implemented in C, available via LGPL license).

This method takes a sufficient set of discrete (rational angle) projections assuming the Dirac pixel model, i.e. digital image sampling where lines have said to have sampled a pixel iff the line goes through the centre of the pixel, and reconstructs them in O(nlogn), where n=N^2. Sufficiency is classified as those projections meeting the Katz criterion, i.e. basically all bins are sampled at least once and no unambiguous  solution (i.e. a ghost) cannot fit within the image. See also:

Fast Mojette transform for discrete tomography

SS Chandra, N Normand, A Kingston, J Guédon, I Svalbe

arXiv preprint arXiv:1006.1965

Once you have FTL built, you should have four binaries for this method:

• fmt_angles - Select the angle type and generate the rational angles for given n and N, the image and reconstruction (FFT space) sizes respectively.
• mt - Compute the discrete (rational angle) projections of the image, also known as the Mojette transform (MT).
• mt2frt - Convert the projections to those of the FRT/DRT, which are the inverse FFTs of the slices of the 2D FFT.
• ifrt - To reconstruct the resulting FRT projections in O(nlogn), no relation to the n before.

To illustrate the process we given a tutorial here of the whole process.

1. First, we crop Lena image to a 128x128 image of Lena from the centre:

./crop lena512.pgm 128 128 0 lena128.pgm

2. We create the angle set, we choose the L1 minimal set since it has a nice symmetry:

./fmt_angles 128 256 1 mt_angles_128_in_256.txt

3. Next we compute the MT:

./mt lena128.pgm mt_angles_128_in_256.txt mt_lena128.pgm

    Note that if you already have projections, such as those of a sinogram, then see this Google Groups     discussion. You can find the publication by my colleague Andrew Kingston on how to do this here.
4. Convert the MT projections into FRT space

 ./mt2frt mt_lena128.pgm mt_angles_128_in_256.txt 128 256 1.0 frt_lena128.pgm

5. Invert the FRT projections in O(nlogn) using the discrete Fourier slice theorem

./ifrt frt_lena128.pgm recon_lena128.pgm

This gives our nxn result reconstructed and padded into the NxN space.

EDIT: See my other post about visualising these and FTL (or PGM files) results in general.

HTH
Cheers Shakes - L3mming