Fourier Theory, Fourier Transforms, Diffraction and Computational Simulation of Transmission Electron Microscopy

Acknowledgement: Prof. S. Allen graciously provided useful physical insight on diffraction interpretation in these notes.

Creating Sample Atomic Data and Lattices for Fourier Analysis and Diffaction Simulation.

Here we create an image such as might be seen in an transmission electron microscope. For our first set of  data, we use intererence of two sine waves to produce a simulation of the density of scattering centers in an atomic lattice.  Most of the physical aspects of atomic imaging and diffraction can be simulated with the two-dimensional techniques that are produced in these notes.

ISize=64;
AtomDensity= N[Table[(1+ Sin[4(x+y)2 Pi/ISize])(1 + Sin[2(x-2y)2 Pi/ISize])/4,{x,1,ISize},{y,1,ISize}]];

We show the simplest visulations of our data, we will improve on the coloring methods below.  Here we show our two dimensional crystal as a gray-intensity map and with a three dimensional plot of scattering intensity over an x-y plane.

GraphicsRow[{ArrayPlot[AtomDensity],ListPlot3D[AtomDensity]},ImageSize-> Full]

Improving Visualization of Data and Diffraction with Contrast Functions

To view the data better, we create two versions of a contrast function, the first (normalcontrast) is useful when we wish to view the entire range of intensities, the second version (highcontrast) when we wish to resolve differences at the low-end of the intensities.

In order to view the data better, we wish to develop contrast functions. Experience suggests that it is useful to construct two type of contrast, this one, normalcontrast will be used as a ColorFunction when we want to resolve the entire range of intensities.

normalcontrast= Function[intens,If[intens<0.00001,GrayLevel[0],Blend[{RGBColor[0,0,.2] ,Red, White},intens]]];

This second version, highcontrast, will be useful when we wish to resolve differences in the low-intensity region. Note how the colors change rapidly in the lower region.

highcontrast=Function[intens,Blend[{Black,GrayLevel[0.1],Lighter[Blue],Lighter[Green],Green,Darker[Green],Lighter[Red],White},Abs[intens]^(0.25)]];

ContrastGraphics=GraphicsColumn[{ArrayPlot[Table[i/256,{j,1,24},{i,1,256}],ColorFunction-> normalcontrast,ImageSize->Large,PlotLabel-> "normalcontrast"],
ArrayPlot[Table[i/256,{j,1,24},{i,1,256}],ColorFunction-> highcontrast,ImageSize->Large,PlotLabel-> "highcontrast"]},ImageSize->Full];

ContrastGraphics

Visualizing Atomic Densities

ImagePlot will be our "work-horse" function to dispay images, as well as Fourier data.  It takes an rectangular array of data as input, and has thee optional arguments for a plot-label, the colorfunction to use (i.e., highcontrast or normalcontrast), and an indication of the image size.

ImagePlot[data_?MatrixQ, label_:None,colfunc_:highcontrast, imagesize_:Medium] := Module[{absdata=Abs[data],min,max},ArrayPlot[absdata,ColorFunction-> colfunc,BaseStyle->{Small,FontFamily->"Helvetica"}, PlotLabel->label,ImageSize->imagesize,Ticks-> False]]

Here we show the normal contrast and the highcontrast images of the simulated atomic data side-by-side.

SimulatedGraphics=GraphicsRow[{ImagePlot[AtomDensity,"Simulated Atomic Density:normal contrast",normalcontrast],ImagePlot[AtomDensity,"Simulated Atomic Density: high contrast (Default)"]}];

SimulatedGraphics

Discrete Fourier Transforms on Simple 2D Lattices

We create a function that shows the original data, its Fourier transform, and then its inverse transform (hopefully) back to the original image.

FourierRow[data_]:= Module[{fourdat=Fourier[data]},
 GraphicsRow[{ImagePlot[data,"Original Data",normalcontrast],ImagePlot[fourdat,"Fourier Transform (ℱ)"],ImagePlot[InverseFourier[fourdat],"(ℱ^(-1))[ℱ[Original Data]],",normalcontrast]},ImageSize->Full ]]

We will see the peaks of the image are located near (kx,ky) = 2 π (a,b)/(Image Size), where 
(a,b) = {(0,0), (Image Size,0), (0,Image SIze), (Image Size, Image Size)}.  These correspond to the longest wavelength periodicities.

FourierRow[AtomDensity]

The magnitude of fourier the fourier transform corresponds to the constructive interference of wave - like objects scattered from a real - space set of data.  When the real-space data has translational symmetry, such as in a lattice, the fourier intensity results in bright spots corresponding to periodicities present in the lattice.The distance between spots is inversely proportional to the distance between periodic scattering features, and the orientation of the (kx, ky) vector corresponds to the normal to the groups of constructively interering scatterers.

Simulating Diffraction

Materials Scientists, Microscopists, and Crystallographers are used to seeing the (OverVector[k])=0 spot at the center of the diffraction image; so we write a function that takes the Fourier data and manipulates to so as to move spots to the center.

KZeroMiddle[fourmat_?MatrixQ]:= Module[{nrows,ncols},{nrows,ncols}=Dimensions[fourmat];RotateRight[fourmat,{Quotient[nrows,2],Quotient[ncols,2]}]]

And, we modify our FourierRow function to use the k-at-zero transformation.

FourierRowK0[data_]:= Module[{fourdat=Fourier[data]},
 GraphicsRow[{ImagePlot[data,"Original Data",normalcontrast],ImagePlot[KZeroMiddle[fourdat],"ℱ[data] (OverVector[k])=(OverVector[0]) centered (high contrast)",highcontrast],ImagePlot[InverseFourier[fourdat],"(ℱ^(-1))[ℱ[Original Data]]",normalcontrast]},ImageSize->Full]]

FourierRowK0[AtomDensity]

The middle image is a good simulation of what is observed for electron diffraction in the TEM.  The white spot (at the center) is called the "direct beam" or the "transmitted beam" and the other spots correspond to diffraction peaks from different (hkl) planes (or (hk) in 2D).  The pattern corresponds the the reciprocal lattice of the crystal object.  A vector from the origin (white spot) to a particular peak is the reciprocal lattice vector for that particular set of diffracting planes.

Alternative Diffraction Visualization

Here is an alternative visualiazation of the diffraction peaks using a third dimension to indicate the intensity of the peak.  We write it so we can specify the range optionally.

Spots3D[data_, range_:All] := ListPlot3D[KZeroMiddle[Abs[Fourier[data]]],ColorFunction->(highcontrast[#3]&),PlotRange-> range,ImageSize->Medium]

Spots3DRow[data_, range_] := Module[{plt},plt=Spots3D[data];
GraphicsRow[{plt,Show[plt,PlotRange-> range]},ImageSize->Full]]

Spots3DRow[AtomDensity,{0,7}]

Visualizing Diffraction Patterns of Crystals with a Defect

Algorithm to Create a Defect in a Simulated Atomic Density Image : 
Here, the algorithm will create new data which will be called 
AtomDensityWithDefect

Look for a maximum in the atomic density map, and replace a circle around it with a "hole"

MinMax[alist_List] := Module[{flatlist=Flatten[alist]},Return[{Min[flatlist],Max[flatlist]}]];
{minval,maxval}=MinMax[AtomDensity];
{randx,randy}=RandomChoice[Position[AtomDensity,maxval]];
MakeDefect[data_?MatrixQ,size_]:=Module[{defect = {},nx,ny,ssq=size^2, minval, maxval, rx,ry},{minval,maxval}=MinMax[data];
{rx,ry}=RandomChoice[Position[data,maxval]];
{nx,ny}= Dimensions[data];Table[If[(i^2 + j^2)≤ ssq,AppendTo[defect,{Mod[i+rx,nx,1],Mod[j+ry,ny,1]}]],{i,-8,8},{j,-8,8}];
Return[ReplacePart[data,defect->0]]]
AtomDensityWithDefect=MakeDefect[AtomDensity,4];

We will show some examples of the diffraction and image reconstruction of different kinds of simulated data, these data will be explored below in the Computational Microscopy sections.

FourierRowK0[AtomDensityWithDefect]

In the above, the defect causes diffuse intensity to appear in the pattern. Using methods produced below, we can investigate  how much of reciprocal space is needed to reproduce the defect in the inverse transform (i.e., the image of the scatters that is produced in a TEM).

Spots3DRow[AtomDensityWithDefect,{0,1/2}]

Visualizing Diffraction Patterns of Crystals Thermal Noise

Function to Create A Lattice of Squares : 
NoisyLattice[TotalSize, LatticeVector1, LatticeVector2, SquareSize, RandomDisplacements]

NoisyLattice[size_,{xa_,ya_},{xb_,yb_},edge_,{noisex_,noisey_}]:= Module[{data,px,py,maxa,maxb,rx,ry,vals},data=vals=ConstantArray[0,{size,size}];
maxa=If[xa≥  ya,Quotient[size,xa],Quotient[size,ya]];
maxb=If[xb≥  yb,Quotient[size,xb],Quotient[size,yb]];
 For[ib=-maxb,ib≤ maxb,ib++,
For[ia=-maxa,ia≤ maxa,ia++,ry=noisey*RandomInteger[];rx=noisex*RandomInteger[];
{px,py} = {Mod[ia*xa+ ib*xb,size,1],Mod[ia*ya+ ib*yb,size,1]};If[vals[[px,py]]==0,vals[[px,py]]=1;
Table[data[[Mod[px+rx +ax,size,1],Mod[py+ry+ay,size,1]]]=1,{ax,0,edge},{ay,0,edge}]]
]
];Return[data]
]

First look at an example with no noise, we create a lattice of squares with lattice vectors (8,4) and (16,16).  This function will create smaller unit cells if the last vectors are divisors of the total size. Here we set the edge of the square to 1 and indicate that we want 0,0 noise in the x,y directions.

NoNoise=NoisyLattice[64,{8,4},{16,16},1,{0,0}];

FourierRowK0[NoNoise]
GraphicsRow[{Spots3D[NoNoise],Spots3D[NoNoise,{0,2}]}]

We reconstruct the example from above, but introduce random unit displacements of the squares in the x and y directions.

SomeNoise=NoisyLattice[64,{8,4},{16,16},1,{1,1}];
FourierRowK0[SomeNoise]
GraphicsRow[{Spots3D[SomeNoise],Spots3D[SomeNoise,{0,2}]}]

The appearance of diffuse intensity (e.g. the intensity that does not correspond to the bright spots for the ideal lattice) corresponds to random shifts in the original figure. This is a simulation of isotropic thermal vibrations. In crystals, the amplitudes would be larger in elastically soft directions.  It is instructive  to use the last two arguments in the NoisyLattice function to make the "noise" anisotropic and observe its effect on the diffration pattern--the streaks that you will observe correspond to observations in the TEM.

Computational Microscopy

Transmission electron microscopy works by accelerating electrons with a voltage difference and sending them towards a target.  The electrons interact with the target and scatter. After the electrons have passed through the sample, they are focused with a magnetiic objective lens (or typically lenses). This lens produces a plane at which electrons scattered in the same direction arrive
at the same point---this is the diffraction pattern.  Because diffraction transforms periodic elements into points, it is closely related to the fourier transform.  An image of the target is created beyond this diffraction plane.  An operator of an electron microscope can toggle between looking at the diffraction-plane or the image-plane.

Specific periodic features can be imaged by blocking all the electrons except those that appear at a point in the diffrations plance. To put it another way, one puts an aperature around a particular spot in order to selectively image the features that from that spot's periodicity. We can do this with a Fourier transform by selectling only a circular subset of points to be subsequently inversely transformed into an image.

The "Bright-Field Image" consists of using a central aperature around the direct beam to block off all others from contributing to the image.

The "Dark-Field Image" consists of selecting a specific diffracted peak with the aperature and using that to form an image

A "Structure image," or a "lattice image," uses the direct beam and one or more diffracted beams to form the image.  In this case, the aperatures are typically much larger than for bright- or dark-field imaging. 

Aperature size is effectively limited because of spherical aberrations that become significant for beams that are "off-axis" by a significant amount.  So, in practice one can only use part of the Fourier spectrum (reciprocal space) to produce an image in TEM.  You always lose some structure information in the image formation.

Function to Create Two Circular Aperatures (i.e, to remove all data from a Fourier Transform except the regions inside two specified circles):

CircAps[Center1, Center2, Radius1, Radius2, FourierData]

(*CircAps[{xpos1_,ypos1_},{xpos2_, ypos2_},rad1_,rad2_,fourmat_?MatrixQ]:=
Module[{rad1sq,rad2sq,ncols,nrows,y1,y2},
{nrows,ncols}=Dimensions[fourmat];
mindim=Min[{nrows,ncols}=Dimensions[fourmat]];
y1=nrows + 1 -ypos1;
y2=nrows + 1 -ypos2;
rad1sq=rad1^2;
rad2sq=rad2^2;
Table[
If[
(Abs[j-xpos1] ≤ rad1 && Abs[i-y1] ≤ rad1 &&(xpos1-j)^2 + (y1-i)^2≤ rad1sq) ||
(Abs[j-xpos2] ≤ rad2 && Abs[i-y2] ≤ rad2 &&(xpos2-j)^2 + (y2-i)^2≤ rad2sq),
fourmat[[i,j]],0],{i,1,nrows},{j,1,ncols}]
]
*)

CircAps[{xpos1_,ypos1_},{xpos2_, ypos2_},rad1_,rad2_,fourmat_?MatrixQ]:=
Module[{ncols,nrows,y1,y2,x1,x2,aperatureImage},
{nrows,ncols}=Dimensions[fourmat];
aperatureImage=ImageData[
Image[
Graphics[{Disk[{xpos1,ypos1},rad1],Disk[{xpos2,ypos2},rad2]},PlotRange-> {{1,ncols},{1,nrows}}],"Bit",ImageSize-> {nrows,ncols},ColorSpace-> "Grayscale"]
];
fourmat*(1-aperatureImage)
]

Function to Perform Diffraction Spot Microscopy on Real-Space Images

DiffractionMicroscopy[RealData]

Creates a Interactive Structure with Four Windows Arranged in a Square.
	
	NorthWest: All Fourier Data from Real Image with Movable Aperatures, use a clicked mouse to move aperatures to various diffraction peaks.
	
	NorthEast: The original real space image.
	
	SouthWest: The Fourier image of the aperature filtered data. (Don't try to move these aperatures)
	
	SouthEast: The reconstructed image from the aperature filtered data.

DiffractionMicroscopy[data_] := 
DynamicModule[{fd,nrows,ncols,image,fourimage,rawimage,applot,gcirc1,gcirc2,apdata },
rawimage=ImagePlot[data,"Input\nData",normalcontrast,Scaled[0.25]];
fd=KZeroMiddle[Fourier[data]];
{nrows,ncols}=Dimensions[data];
fourimage= ImagePlot[fd,"Fourier Representation of Data\nClick on Aperatures to Move",highcontrast,Scaled[0.5]];
Manipulate[
gcirc1={Lighter[Purple],Thick,Circle[pt1,radius1 ]};
gcirc2={Yellow,Thick,Circle[pt2,radius2 ]};
apdataGlob = apdata=CircAps[pt1,pt2,radius1,radius2,fd];
Grid[
{{Show[fourimage,Epilog->{gcirc1,gcirc2},ImageSize->Scaled[0.25]],
rawimage
},
{
applot=Show[ImagePlot[apdata,"Selected Fourier Data",highcontrast],Epilog-> {gcirc1,gcirc2},ImageSize-> Scaled[0.25]];
Show[applot],
ImagePlot[InverseFourier[apdata],"Filtered Image",normalcontrast,Scaled[0.25]]
}}
],
{{radius1,ncols/20,"Purple\nRadius"},0,Round[0.707*ncols],1},
{{radius2,ncols/45,"Yellow\nRadius"},0,Quotient[ncols,20],1},
{{pt1,{nrows/4,ncols/6}},Locator,Appearance->None},{{pt2,{nrows/2,ncols/2}},Locator,Appearance->None},ControlPlacement->{Top,Bottom},ContinuousAction->False,SaveDefinitions->True]
]

Function to Create Simulated Grain Structure for Subsequent Computational Microscopy Demo

GrainStructure[TotalSize, GrainSize]

Creates an array of "square grains" with stripes oriented in somewhat random orientations

GrainStructure[size_,grainsize_]:= Module[{xgrains,ygrains, ngrains,rotmat,thetalist,xrot,yrot},
xgrains=size/grainsize;
ygrains = size/grainsize;
ngrains=xgrains*ygrains;
thetalist=RandomSample[Table[xc*Pi/ngrains,{xc,1,ngrains}]];
rotmat[i_,j_] :=  rotmat[i,j]=N[RotationMatrix[thetalist[[(i-1)*ygrains + Mod[j,xgrains,1]]]]];
Return[Table[{xrot,yrot}=rotmat[Quotient[xv-1,grainsize]+1,Quotient[yv-1,grainsize]+1].{xv,yv};(1+ Sin[xrot 2 Pi 5/grainsize])/2,{xv,1,size},{yv,1,size}]
]
]

Simulated Grain Structure (Polycrystal)

Here is a simple example of four grains with the same modulated density, but each with a different orientation.

Grains=GrainStructure[128,64];
ImagePlot[Grains,"Simulated Grain Structure",normalcontrast,Large]

Bright Field, Dark Field, Lattice Imaging of a Polycrystal

Diffraction Spot Microscopy on grain structure.  Notice that if you leave the small yellow aperature at the center, and move the green aperature over other diffraction spots, then the grain in that particular orientation will "light up" and its modulation will appear.  Now toggle to the "Dark-Field Imaging" mode by  removubg the yellow aperature (by shrinking its radius to zero); notice that only one grain "lights up" but without any modulation.

DiffractionMicroscopy[Grains]

Bright Field, Dark Field, Lattice Imaging of a Polycrystal

Diffraction Spot Microscopy on Idealized Lattice.  Notice that when you select pairs of spots, that particular atom planes are visualized. Use the purple aperature to select a spot on the same ray (i.e., with the same plane normals), the modulation doubles. You are visualizing (abc)-planes in one case and (abc)/2 planes in the second case.  When two apperatures are used, you are using "lattice image"  or "structural image" mode.

DiffractionMicroscopy[NoNoise]

Diffraction Spot Microscopy on Noisy Lattice.  Notice that when you select pairs of the brightest spots, you still get a fairly robust image of those planes, even though they are not perfectly alligned in the original image.  Most of the information is contained in the sharp peaks and those must represent perfect modulations at that point. The noise information is diffused over the entire Fourier spectrum.

DiffractionMicroscopy[SomeNoise]

Frivolous Example: Beautiful Data are in the Eye of the Beholder

Import some random intensities that we have collected off the web:

AnImage = Import["http://pruffle.mit.edu/3.016/MP.png","GrayLevels"]/.{a_,b_}-> a;

AnImage

Let' s image that data with ArrayPlot

ArrayPlot[AnImage]

Oops, it' s  inside out.

AnImage = Import["http://pruffle.mit.edu/3.016/MP.png","GrayLevels"]/.{a_,b_}-> 1-a;

ArrayPlot[AnImage]

Let's view the data in real space and diffract it.

DiffractionMicroscopy[AnImage]