```
%------------------- MATLAB
```

% read image in, convert to greyscale and double-precision

image=double(rgb2gray(imread('myimage.jpg')));

% crop to smallest dimension

[m,n]=size(im);

im=im(1:min(m,n),1:min(m,n));

% find image mean

imMean=mean(im(:));

% and subtract mean from the image

demean=im-imMean;

% find the 2D fourier transform

ft=fft2(demean);

% power spectrum is the element-by-element square of the fourier transform

spectrum=abs(ft).^2;

Seven lines of commands. Couldn't be simpler. That's the beauty of MATLAB. The snag is that it costs several hundreds to thousands of dollars for the program, which makes it difficult to share your code with other people. So recently I've thought about coding up some things I want to share with other people in Python.

Python is a cross-platform, open-source programming language and environment. It is not too dissimilar to MATLAB in many respects, especially when you install the right libraries. In order to do the same in Python as above, as easy as possible, I needed to install Python Imaging Library (PIL), NumPy and SciPy. Then the sequence of operations is as follows:

```
# ----------------PYTHON
```

#-----------------------

# declare libraries to use

import Image

import scipy.fftpack

from pylab import*

import cmath

from numpy import*

# load image, convert to greyscale

im=Image.open('myimage.jpg').convert("L")

# crop to smallest dimension

box=(0,0,im.ny,im.nx)

region=im.crop(box)

# find mean

imMean=mean(asarray(region))

# 2D fourier transform of demeaned image

ft=fft2(region-imMean)

# power spectrum

spectrum=pow(abs(ft),2)

Same answer. Not bad - not harder, just different. And completely free. Oh, and this'll work on UNIX, Linux, Mac, and PC when Python is installed (with a couple of extra libraries). Great!

Hmm, I thought, I wonder if it is possible to do the same in a piece of non-numerical software. I tried ImageMagick because its free and ubiquitous on Linux and Mac, and also available for PC, and it has a massive user-community. It is incredibly powerful as a command-line graphics package - even more powerful than PhotoShop and Illustrator, and like all UNIX-based applications, infinitely flexible. It turns out it is a little more tricky to calculate the power spectrum in ImageMagick, but here's how.

First, you need the right tools - ImageMagick is usually already present on most Linux distributions - if not, if a debian/Ubuntu user use sudo apt-get install imagemagick, or download from the webiste hyperlinked above. Then you'll need an open-source C-program for calculation of FFTs: this is provided by FTW - download fftw-3.2.1.tar.gz, untar, cd to that directory, then install the binary by issuing the usual commands:

This C-program is no good on its own, so the last ingredient is a wrapper for ImageMagick. Fortunately - nay, crucially, this has been provided by Sean Burke and Fred Weinhaus - introducing IMFFT! The program can be downloaded and installed from here. And here is an in-depth tutorial of how to use it (I'm only just getting to grips with it myself). For some reason the compiled program is called 'demo'. I move it to my /bin directory so my computer can see it wherever I am working, and in the process change the name to something more suitable - fft:

Right, it's a little bit more involved, but here is the bash shell script for the power spectrum in ImageMagick:

Voila! The spectrum is the image 'ggc_dm_F_Q16_spec.tiff'. I haven't yet worked out how to keep all the data inside the program, to eliminate the need (and cumulative errors associated with) writing the outputs at every stage to file, and eading them back in again. This means I can't use the spectrum for any scientific application. But in theory it should be possible. Any suggestions welcome!

```
./configure
```

sudo make && make install

This C-program is no good on its own, so the last ingredient is a wrapper for ImageMagick. Fortunately - nay, crucially, this has been provided by Sean Burke and Fred Weinhaus - introducing IMFFT! The program can be downloaded and installed from here. And here is an in-depth tutorial of how to use it (I'm only just getting to grips with it myself). For some reason the compiled program is called 'demo'. I move it to my /bin directory so my computer can see it wherever I am working, and in the process change the name to something more suitable - fft:

```
cp demo ./fft && sudo mv ./fft /bin
```

Right, it's a little bit more involved, but here is the bash shell script for the power spectrum in ImageMagick:

```
#----------------------------- IMAGEMAGICK
```

#-------------------------------------------------

#!/bin/bash

# read image and convert to greyscale

convert myimage.tif -colorspace gray input_grey.tiff

# find min dimension for cropping

min=`convert input_grey.tiff -format "%[fx: min(w,h)]" info:`

# crop image

convert input_grey.tiff -background white -gravity center -extent ${min}x${min} input_grey_crop.tiff

# de-mean image

# get data

data=`convert input_grey_crop.tiff -verbose info:`

# find mean as proportion

mean=`echo "$data" | sed -n '/^.*[Mm]ean:.*[(]\([0-9.]*\).*$/{ s//\1/; p; q; }'`

# convert mean to percent

m=`echo "scale=1; $mean *100 / 1" | bc`

# de-mean image

convert input_grey_crop.tiff -evaluate subtract ${m}% ggc_dm.tiff

# pad to power of 2

p2=`convert ggc_dm.tiff -format "%[fx:2^(ceil(log(max(w,h))/log(2)))]" info:`

convert ggc_dm.tiff -background white -gravity center -extent ${p2}x${p2} ggc_dm_pad.tiff

# forward fast fourier transform

fft f ggc_dm_pad.tiff

# creates ggc_dm_F_Q16.tiff, on which get spectrum

convert ggc_dm_F_Q16.tiff -contrast-stretch 0% -gamma .5 -fill "rgb(1,1,1)" -opaque black ggc_dm_F_Q16_spec.tiff

Voila! The spectrum is the image 'ggc_dm_F_Q16_spec.tiff'. I haven't yet worked out how to keep all the data inside the program, to eliminate the need (and cumulative errors associated with) writing the outputs at every stage to file, and eading them back in again. This means I can't use the spectrum for any scientific application. But in theory it should be possible. Any suggestions welcome!

## 1 comment:

:-) cool!

Post a Comment