Graef M. Introduction to conventional transmission electron microscopy

Подождите немного. Документ загружается.

108 Basic quantum mechanics

(Chapter 10), we are in the unfortunate situation that the Fourier transformed point

spread function for the TEM usually contains multiple zero crossings.

2.5.4 Numerical computation of Fourier transforms and convolutions

The Fourier transform, as deﬁned in equations (2.45) and (2.46), is a continuous

integral transform, suitable for analytical computations. Numerical computations

operate on discrete functions rather than continuous functions. The discrete di-

rect Fourier transform (DDFT) of a one-dimensional function f (x

) deﬁned on a

discrete grid of N equidistant points x

= jδ, with δ the grid stepsize, is deﬁned as

f (q

) =

N −1



j=0

f (x

−2πiij/N

. (2.57)

The “frequencies” q

are deﬁned by q

= j/N δ. A straightforward (but naive)

numerical implementation of this equation would require the computation of 2N

trigonometric functions (cos and sin). In terms of processor cycles, this would be a

very expensive computation.

It is customary to introduce a shorthand notation for the exponential factor:

≡ e

2π

, (2.58)

and hence

f (q

) =

N −1



j=0

f (x

)ω

−ij

. (2.59)

The discrete inverse Fourier transform (DIFT) is obtained by changing the sign of

the exponent of ω

f (x

) =

N −1



i=0

f (q

)ω

. (2.60)

Higher-dimensional discrete Fourier transforms are deﬁned in a similar way. The

two-dimensional DDFT is deﬁned by

f (q

, r

) =

N −1



k=0

M−1



l=0

f (x

, y

)ω

−ik

−jl

. (2.61)

The DDFT is essentially a matrix multiplication. If we denote the function values

f (x

)byX

, and f (q

)byQ

, then it is easy to see that for N = 4, equation (2.59)

2.5 Fourier transforms and convolutions 109

is explicitly given by



















1111

1 ω

−1

−2

−3

1 ω

−2

−4

−6

1 ω

−3

−6

−9



















, (2.62)

or, in matrix notation:

Q = MX. (2.63)

The matrix M is a symmetric matrix, and only four different numbers need be

computed in the case N = 4; since ω

−k

= ω

−k±N

,wehaveω

−6

= ω

−2

, ω

−9

= ω

−1

and ω

−4

= 1, and the matrix reduces to

M =







1111

1 ω

−1

−2

−3

1 ω

−2

1 ω

−2

1 ω

−3

−2

−1







. (2.64)

In general, the N × N matrix M contains only N different numbers, which are all

powers of a single complex number ω

. Note that ω

requires two evaluations of

trigonometric functions, using Euler’s formula e

= cos x + i sin x. The number

and its powers are known as twiddle factors, and their computation can be

performed in a very efﬁcient way using the so-called Singleton algorithm [Sin67],

which minimizes the total number of trigonometric function calls. For a description

of the Singleton algorithm we refer the reader to Chu and George [CG00, pp. 14–15].

After the computation of the twiddle factors, we could in principle construct the

matrix M and compute the product MX. This would again be inefﬁcient from

the computational point of view, and it would also require signiﬁcant computer

memory for large values of N , and for multi-dimensional arrays. In 1965 Cooley

and Tukey [CT65] published an algorithm for the fast computation of summations

of the type (2.59). Their algorithm became known as the fast Fourier transform

(FFT) algorithm. The FFT algorithm makes use of the fact that if N is a power of

2, then the original summation (2.59) can be split into two summations:

N/2−1



j=0

2 j

−i(2 j )

−i

N/2−1



j=0

2 j+1

−i(2 j )

. (2.65)

Since ω

= ω

N/2

and also ω

N/2

= ω

(i+N /2) j

N/2

,weﬁnd

N/2−1



j=0

2 j

−ij

N/2

−i

N/2−1



j=0

2 j+1

−ij

N/2

, (2.66)

110 Basic quantum mechanics

in which each summation is now a DDFT of size N/2 [CG00, pp. 22–23]. Each

of the new summations can again be split into two summations of size N/4, and

we repeat this k times until N/k = 2, which is a trivial sum. The computational

gains obtained from using the FFT algorithm are signiﬁcant. A conventional “brute

force” DDFT of a 1-D array with N = 2

entries takes of the order of N

ﬂoat-

ing point operations. The FFT algorithm takes only N log

N = nN operations,

which represents a signiﬁcant time saving factor: for N = 1024 we have n = 10

and hence the FFT takes of the order of 10 240 operations compared with 1048 576

for the naively implemented DDFT. For multi-dimensional discrete Fourier trans-

forms, the time saving factors are even larger. We will make extensive use of the

FFT in Chapters 7 and 10, when we talk about the multi-slice image simulation

method.

There are many different numerical implementations of this recursive algorithm.

The algorithm can be extended to arrays of arbitrary size N ; in such cases the trans-

form is decomposed into summations corresponding to N = 2

...,

in other words, N is decomposed into powers of prime numbers. While the FFT

algorithm is intuitively easy to understand, the reader is advised against attempting

to implement the algorithm him/herself; there are many commercial

and public do-

main implementations available, and most of them have been optimized for speed.

In fact, in a 1995 review article, Sorensen, Burrus, and Heideman [SBH95] list

more than 30 different implementations, and a total of about 3500 papers devoted

to FFTs, since their inception in 1965.

For a detailed description of the basic algorithm and its implementation we

refer to Numerical Recipes: The Art of Scientiﬁc Computing (FORTRAN Ver-

sion) [PFTV89], Section 12.2. There are several public domain implementations

available, among others the

fftpack library (www.netlib.org). The source code

in this book makes use of the

FFTW public domain library of FFT routines;

†

the

source code can be downloaded from

http://www.fftw.org. The library is easy

to install and test. Although all of the routines are written in C, the package provides

a set of wrapper-routines that convert the internal data to a column-based ordering,

so that Fortran programs can call all of the C routines. Instructions for installing the

FFTW library can be found on the website. The main reasons for using the FFTW

library are its speed, and the fact that the algorithm automatically selects the best

decomposition of N into prime numbers, so that the restriction N = 2

, common

to many implementations, is not present. This will become useful in later chapters,

where we will need to compute FFTs of arrays for which the dimensions are not

simple powers of 2.

†

FFTW was written by Matteo Frigo and Steven G. Johnson. The source code is included on the website under

the GNU General Public Licence.

2.6 The electrostatic lattice potential 111

Convolution products can be computed readily using the direct and inverse FFT

routines and the convolution theorem (equation 2.54):

f (r) ⊗ g(r) = F

−1

[

F[ f (r)]F[g(r)]

]

, (2.67)

and therefore do not require a specialized algorithm. The use of this equation

assumes that both functions have been deﬁned on the same discrete grid.

2.6 The electrostatic lattice potential

The Schr¨odinger equation (2.37) contains two important contributions: the kinetic

energy term, which is experimentally determined by the acceleration voltage of the

electron microscope, and the electrostatic lattice potential V (r), which describes

the interaction of the beam electrons with the sample. Finding the correct potential

for a given crystal structure is one of the more fundamental problems of solid-state

physics, and most so-called ﬁrst-principles methods produce this potential, often

as the result of a self-consistent iteration. Fortunately for electron microscopists, it

turns out that, to a very good approximation, we can consider a solid to consist of

spherical point scatterers, represented by atomic scattering factors f

(s). The true

lattice potential is then the potential due to spherical scatterers plus contributions

from interatomic bonds. For most TEM applications, we can safely ignore bonding

contributions and assume that the solid consists of isolated spherical point scatterers.

In the following sections, we will ﬁrst deﬁne the atomic scattering factors, and then

derive expressions for the Fourier coefﬁcients of the electrostatic lattice potential

for an inﬁnite and a ﬁnite crystal. We will only consider elastic scattering events;

i.e. scattering events for which only the direction of the wave vector changes. The

electron wavelength, or, equivalently, its energy, remains unchanged.

2.6.1 Elastic scattering of electrons by an individual atom

The scattering of electrons by a single atom is determined by the total charge

distribution of that atom:

ρ(r) =|e|[ρ

(r) − ρ

(r)],

where ρ

is the nuclear charge distribution and is ρ

the electron charge distribution.

The electrostatic atom potential V

(r) is related to the total charge distribution

through Poisson’s equation:

V

(r) =−

|e|



[ρ

(r) − ρ

(r)], (2.68)

112 Basic quantum mechanics

where  is the Laplacian differential operator. The probability amplitude for an

electron to be scattered in a certain direction can now be computed using the bra-

ket concept introduced earlier in this chapter. The probability amplitude P that an

incident plane wave with wave vector k will be scattered by the atomic potential

(r) into the direction k



is given by

P =

2πik



·r

(r)

2πik·r

. (2.69)

We deﬁne the atomic scattering amplitude, f

(k), to be equal to this probabil-

ity amplitude, or, equivalently, to the Fourier transform of the atomic potential

distribution V

(r):

(k) ≡

...

(r)e

−2πik·r

dr,

k∆

(2.70)

where k = k



− k is the momentum transfer vector. The atomic scattering am-

plitude thus describes how the momentum of an incident electron changes upon

elastic scattering by the atom. The atomic scattering factor expresses the probability

amplitude that a certain elastic momentum transfer will occur.

Using the fact that the x-ray scattering amplitude f

(k) is deﬁned as the Fourier

transform of the electron charge density ρ

(r),

†

[Shm96], and approximating the

nuclear charge density by a delta-function of weight Z (the atomic number) located

in the origin, we can introduce the inverse Fourier transforms of all quantities in

Poisson’s equation:





...

(k)e

2πik·r

dk



|e|



...

(k) − Z

2πik·r

dk.

(2.71)

Bringing the Laplacian operator (operating on r) inside the integral and dropping

all integrals we ﬁnd:

(k) =

|e|

4π



|k|

Z − f

(k)

. (2.72)

This relation is valid for an isolated atom. Since we will be studying crystals, with

regular arrangements of atoms, we will often rewrite the scattering vector k as

k = 2s; for the Bragg condition, k



= k + g or, equivalently, s = g/2. Using the

†

The contribution of the nuclear charge to x-ray photon scattering is several orders of magnitude smaller than

that of the orbital electrons because of the mass difference between electrons and nucleons.

2.6 The electrostatic lattice potential 113

Bragg equation we ﬁnd for the magnitude of s:

|s|=

hkl

sin θ

hkl

. (2.73)

This means that for a given crystal structure the atomic scattering factor is “sampled”

only at scattering vectors s corresponding to half the reciprocal lattice vectors.

Since s is independent of the wavelength of the electrons (1/2d is constant for a

given crystal structure!), the following expression for the atomic scattering factor

is independent of the experimental conditions:

(s) =

|e|

16π



|s|

Z − f

(s)

(Mott–Bethe formula). (2.74)

Fig. 2.9 shows the electron scattering factor for a copper atom; the vertical lines

indicate the reciprocal lattice spacings for the face-centered cubic Cu structure

(a = 0.361 47 nm). Scattering functions for other atoms can be displayed using the

ION-routine

scatfac.pro on the website.

2468100

f (s) for Cu

s (nm )

-1

111

200

220

311

222

400

331

420

422

333

300 K

1000 K

Fig. 2.9. Atomic scattering factor for Cu, with superimposed lines indicating the ﬁrst 10 re-

ciprocal lattice spacings |g|/2 for face-centered cubic (fcc) Cu. The scattering factor is

“sampled” at those values of s = sin θ/λ. The values in this ﬁgure must be multiplied by

0.047 878 01 to yield scattering factors in units of V nm

. The dashed lines correspond to the

scattering factor for fcc copper at T = 300K(B(300) = 0.005 75 nm

) and T = 1000 K

(B(1000) = 0.018 69 nm

) [GP99].

114 Basic quantum mechanics

The electron scattering factors have been calculated by many researchers and

are often parameterized in a power series or an exponential expansion. For the

low atomic number atoms, this computation is actually not very difﬁcult and one

can use programs such as the Herman–Skillman code (based on non-relativistic

Hartree–Fock–Slater equations) [HS80] to compute the atomic scattering factors.

For heavier atoms, relativistic corrections become important and the computation

of f

(s) becomes far more involved. There are several parameterizations of the

electron and/or x-ray scattering factors available in the literature:

Smith and Burge [SB62];

Doyle and Turner [DT68];

Weickenmeier and Kohl [WK91];

Rez, Rez and Grant [RRG94].

A more extensive list and discussion can be found in Appendix D of E.J. Kirkland’s

book on Advanced Computing in Electron Microscopy [Kir98, pp. 199–219]. The

source code available from the

website uses two sets of scattering factors: a combi-

nation of Doyle–Turner and Smith–Burge parameterizations, identical to that used

in the

EMS package [Sta87], and the Weickenmeier–Kohl

FSCATT routine.

†

Both Doyle–Turner and Smith–Burge parameterizations are exponential expan-

sions of the form:

(s) =

|e|

16π





Z − f

(s)



= 0.047 878 01



j=1

−b

, (2.75)

where N = 4 for Doyle–Turner parameters and N = 3 for Smith–Burge parame-

ters; s must be expressed in nm

−1

. The units of f

(s) are then V nm

. This expansion

is accurate for values of s up to about 20 nm

−1

. The parameterization can be found

in the module

diffraction.f90 on the website.

Example 2.2 Compute the electron scattering factor f

for the (111) planes of

copper, with lattice parameter a = 0.361 47 nm.

Answer: First, we must evaluate the scattering parameter s = sin θ/λ = 1/2d for

this particular situation. The value for d

111

in a cubic crystal is easily found using

hkl

√

+ k

from which we ﬁnd d

111

= 0.361 47/

√

3 = 0.208 695 nm. The scattering parameter

s is then equal to 2.3958 nm

−1

. Substitution of this value, and the parameters for a

†

The original source code of FSCATT was developed by Andreas Wieckenmeier and Helmut Kohl, then at the

Technische Hochschule in Darmstadt, Germany, and is included on the

website with the permission of the

authors.

2.6 The electrostatic lattice potential 115

and b

for copper (see diffraction.f90 ﬁle), into the atomic scattering factor equation

results in

(2.3958) = 0.047 878 01 ×

1.579e

−1.820×(2.3958)

+ 1.820e

−0.124 53×(2.3958)

+1.658e

−0.025 04×(2.3958)

+ 0.532e

−0.003 33×(2.3958)

;

= 0.047 878 01 × 2.891 086 = 0.138 42 Vnm

Note that this number is independent of the wavelength of the electrons being used

since the number 1/2d

hkl

is independent of the wavelength.

The International Tables for Crystallography provide tables of the electron scat-

tering factors for all elements and their most common ions [WP99, pp. 226–244].

†

The values in these tables must be multiplied by the factor of 0.047 878 01 to give

the numerical values above. The scattering factor for copper, shown in Fig. 2.9,

must also be multiplied by this constant to have units of V nm

. When comparing

scattering factors from various sources, the reader should study carefully how the

factors have been deﬁned in each source.

The Weickenmeier–Kohl expansion uses seven ﬁtting parameters for the elastic

scattering amplitudes of the neutral atoms. The ﬁtting function is of the form

(s) =



i=1

)

1 − e

−B

, (2.76)

where

= A

= 0.023 95

3(1 + V )

−1

; A

= A

= VA

The speciﬁc deﬁnition of A

guarantees that the asymptotic behavior of the scat-

tering amplitude for large s has the correct shape. The values of the six param-

eters B

and V are listed in Table 1 in [WK91], and are tabulated in the FSCATT

function in the others.f90 module. Since this function also returns a second param-

eter, the absorptive form factor, we will defer a more complete description until

Section 6.5.1.

†

The reader should be warned that the older edition of the International Tables for Crystallography [HL52]

uses a different deﬁnition of the electron scattering factors; they are commonly denoted by f

(B for Born

approximation) and are related to the current deﬁnition by

(s) =

2πγm

|e|

(s).

The disadvantage of this formulation is that these scattering factors are dependent upon the acceleration voltage

(through γ ). The factors f

are expressed in nm when s is in nm

−1

. For a more detailed discussion of the

relation between the two sets of factors see [SZ92, pp. 30–35].

116 Basic quantum mechanics

2.6.2 Elastic scattering by an inﬁnite crystal

An inﬁnite point lattice can be described by a set of unit-weight delta-functions,

located at the lattice positions:

T (r) =

+∞



u,v,w=−∞

δ(r − t

uvw

). (2.77)

Each delta-function corresponds to a single unit cell of the lattice. Every unit cell

contains N atoms located at positions r

and the potential of a single unit cell is

then described by

cell

(r) =



j=1

(r − r

The complete lattice potential is determined by “copying” the single unit cell po-

tential at every lattice site; this operation is described by a convolution product:

V (r) = V

cell

(r) ⊗ T (r).

The potential V (r) has, by construction, the periodicity of the underlying Bravais

lattice. This allows us to expand the potential as a discrete Fourier series; this is

essentially a plane wave expansion, where the “planes” are the actual lattice planes

described by the reciprocal lattice vectors g:

V (r) =



2πig·r

; (2.78)

the coefﬁcients V

are the Fourier coefﬁcients of the electrostatic lattice potential.

There is one coefﬁcient (in general, a complex number) for each set of planes in the

crystal. Since the potential V (r) is real, we must have V

−g

= V

∗

. If, in addition,

the origin of the reference frame is located at an inversion center, then the Fourier

coefﬁcients are all real numbers and we have V

−g

= V

We anticipate that, for an inﬁnite crystal, the Fourier transform of the lattice

potential will consist of discrete peaks at the reciprocal lattice positions. As we will

show below, for a ﬁnite crystal this is no longer true, and the Fourier coefﬁcients

become continuous functions of a reciprocal space vector q,or

V (q) =



...

V (r)e

−2πiq·r

dr. (2.79)

This can be rewritten as

V (q) =



F[V (r)],

2.6 The electrostatic lattice potential 117

and using the convolution theorem we have

V (q) =



F[V

cell

(r)]

12 3

F[T (r)].

12 3

(2.80)

The ﬁrst factor on the right-hand side of equation (2.80) describes the contribu-

tion of a single unit cell. Using the deﬁnition of the atomic scattering factors, we

arrive at

F[V

cell

(r)] =



j=1

−2πiq·r

. (2.81)

It is easy to show that the second factor in equation (2.80) can be rewritten as a sum

of delta-functions at the reciprocal lattice positions (e.g. [Cow81, p. 43]):

F[T (r)] =



δ(q − g). (2.82)

We deﬁne the Fourier coefﬁcients of the electrostatic lattice potential as

≡





j=1

−2πig·r

, (2.83)

and we ﬁnd for the Fourier transform of the potential

V (q) =



δ(q − g). (2.84)

As we anticipated earlier, the Fourier transform (for an inﬁnite crystal) is a discrete

function which is non-zero only at the reciprocal lattice points. The numerical

computation of the Fourier coefﬁcients V

will be discussed in Section 2.6.7.

Example 2.3 Compute the probability P that an electron will be scattered from a

state with wave vector k into a state with wave vector k



by a periodic electrostatic

potential V (r).

Answer: A periodic potential V (r) can be written as a Fourier series:

V (r) =



2πig·r