Fourier Image Analysis

Settings and Helper Functions

Avoid error messages for similar symbols (e.g. viewI, viewF, viewC, viewM)

Off[General::spell];Off[General::spell1];

Re-align images

recenter[array_]:=RotateLeft[array,Floor[Dimensions[array]/2]]

Create {image,fourier} from source image or fourier-transform image

SetOptions[Fourier,FourierParameters->{1,1}];
SetOptions[InverseFourier,FourierParameters->{1,1}];
ItoIF[imageI_]:={imageI,Chop[Fourier[imageI]]};
FtoIF[imageF_]:={Chop[InverseFourier[imageF]],imageF}

Import color image and convert to {grayscale,fourier}

loadIF[imagePath_]:=ItoIF[recenter[ImageData[ColorConvert[Import[imagePath],"Grayscale"]]]]

Given function of x and y, and range for x and y, create {grayscale,fourier}

createIF[fn_,range_]:=ItoIF[
recenter[Table[fn,{x,-range[[1]]/2,range[[1]]/2-1},{y,-range[[2]]/2,range[[2]]/2-1}]]];
createFI[fn_,range_]:=FtoIF[
recenter[Table[fn,{x,-range[[1]]/2,range[[1]]/2-1},{y,-range[[2]]/2,range[[2]]/2-1}]]];

Create {image,fourier} from function[x,y]

convolve[imageA_,imageB_] := FtoIF[imageA[[2]] imageB[[2]]];
multiply[imageA_,imageB_] := ItoIF[imageA[[1]] imageB[[1]]]

Functions for viewing image and Fourier transform data.
View image, with color values scaled to fill range. View fourier transform, scaled so mean+one standard deviation is white (otherwise, DC terms swamps other terms and you just get an image with one visible dot in the middle)

viewI[im_]:=ImageAdjust[Image[recenter[Abs[im[[1]]]]]];
viewF[im_]:=Block[
{a=Abs[im[[2]]],m, s},m=Mean[Flatten[a]];s=StandardDeviation[Flatten[a]];
ImageAdjust[Image[recenter[a]],1,{0,m+s}]];

Display image, transform and inverse(transform(image)) side-by-side

viewIF[im_]:=GraphicsRow[{viewI[im],viewF[im]}]

Layout with two source images pairs (image/fourier) and result image below.

view3[imageA_,imageB_,imageC_]:=GraphicsGrid[{
{viewIF[imageA],SpanFromLeft,viewIF[imageB],SpanFromLeft},
{"",viewIF[imageC],SpanFromLeft,""}}]

Basic transform/inverse transform

viewIF[createIF[Cos[((x+32)/2+y)/8*π],{64,64}]]

image=loadIF["/Users/olano/Class/demos/fourier/mathematica/maggy.jpg"];
viewIF[image]

square=createIF[If[Abs[x]<2.5&&Abs[y]<2.5,1,0],Dimensions[image[[1]]]];
viewIF[square]

view3[image,square,convolve[image,square]]

circle=createIF[If[(x^2)+(y^2)<(2.5^2),1,0],Dimensions[image[[1]]]];
viewIF[circle]

view3[image,circle,convolve[image,circle]]

fsquare=createFI[If[(x^2)+(y^2)<(20^2),1,0],Dimensions[image[[1]]]];
viewIF[fsquare]

view3[image,fsquare,convolve[image,fsquare]]

gaussian=createIF[N[(1/(5(Sqrt[2π])))Exp[-(((x^2)+(y^2))/(2*(5^2)))]],Dimensions[image[[1]]]];
viewIF[gaussian]

view3[image,gaussian,convolve[image,gaussian]]

dg=createIF[Derivative[0,1][(1/(5(Sqrt[2π])))Exp[-((((#1)^2)+((#2)^2))/(2*(5^2)))]&][x,y],Dimensions[image[[1]]]];
viewIF[dg]

view3[image,dg,convolve[image,dg]]

sample=createIF[If[Mod[x,5]==2 && Mod[y,5]==2,1,0],Dimensions[image[[1]]]];
viewIF[sample]

sampled=multiply[image,sample];
view3[image,sample,sampled]

view3[sampled,square,convolve[sampled,square]]

view3[sampled,gaussian,convolve[sampled,gaussian]]

lowpass=FtoIF[Block[{s,s2},s=Dimensions[image[[1]]];s2=s/2; 
Table[N[If[
Abs[Mod[x+s2[[1]],s[[1]]]-s2[[1]]]<s2[[1]]/5&&Abs[Mod[y+s2[[2]],s[[2]]]-s2[[2]]]<s2[[2]]/5,1,0]],
{x,0,s[[1]]-1},{y,0,s[[2]]-1}]]];
viewIF[lowpass]

lowpassImage=convolve[image,lowpass];
view3[image,lowpass,lowpassImage]

sampledLowpass=multiply[lowpassImage,sample];
view3[lowpassImage,sample,sampledLowpass]

view3[sampledLowpass,lowpass,convolve[sampledLowpass,lowpass]]

Pattern Matching

mask=createIF[N[If[(x^2)+(y^2)<(7^2),1,0]],Dimensions[image[[1]]]];
viewI[{(.5+mask[[1]])*RotateLeft[image[[1]],{-11,7}]}]

pattern=ItoIF[mask[[1]]*RotateLeft[image[[1]],{-11,7}]];
crossCorrelation=convolve[image,Conjugate[pattern]];view3[image,Conjugate[pattern],crossCorrelation]

patternSq=Total[((pattern[[1]])^2),2];
imageSq=convolve[ItoIF[((image[[1]])^2)], mask];
match=imageSq[[1]]-2 crossCorrelation[[1]]+patternSq;
GraphicsRow[{GraphicsGrid[{
{viewI[image],viewI[imageSq]},
{viewI[pattern],viewI[crossCorrelation]}}],
ImageApply[If[#1<0.1,{1,0,0},{#1,#1,#1}]&,viewI[{match}]]}]