 

THESE

ZOlO

 

LIPRARY
MiChI‘U state
University

 

 

 

This is to certify that the
dissertation entitled

A SEQUENTIAL MONTE CARLO BASED RECURSIVE
TECHNIQUE FOR SOLVING NDE INVERSE PROBLEMS

presented by

TARIQ MAIRAJ RASOOL KHAN

has been accepted towards fulfillment
of the requirements for the

PhD degree in Electrical Engineerig

 

 

Kill/~44

Major Professor's Signature

DG—cemgc-A {04 looq

 

Date

MSU is an Afﬁrmative Action/Equal Opportunity Employer

 

.—.._.-._.-.---o-n-I-o-n-.-.-a-.--—.-.-.-.-.—

 

PLACE IN RETURN BOX to remove this checkout from your record.
TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

 

DATE DUE DATE DUE DATE DUE

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5I08 K:IProi/Aoc&Pros/ClRC/DateDue.indd

A SEQUENTIAL MONTE CARLO BASED RECURSIVE TECHNIQUE
FOR SOLVING NDE INVERSE PROBLEMS
By

Tariq Mairaj Rasool Khan

A DISSERTATION

Submitted to
Michigan State University
in partial fulﬁllment of the requirements
for the degree of

DOCTOR OF PHILOSOPHY
Electrical Engineering

2009

ABSTRACT
A SEQUENTIAL MONTE CARLO BASED RECURSIVE TECHNIQUE FOR
SOLVING NDE INVERSE PROBLEMS
By
Tariq Mairaj Rasool Khan
Flaw proﬁle estimation from measurements is a typical inverse problem in
electromagnetic nondestructive evaluation (NDE). The NDE inverse problem is ill-
posed like most other inverse problems. Major issues with conventional solutions of
inverse problems include inaccurate solutions in the presence of noise and a high
computational cost. This research proposes a computationally efﬁcient and robust
approach for solving inverse problems, with a focus on NDE applications. In this
research the inverse problem is formulated in terms of a statistical inverse problem in
which posterior densities of unknown parameter (such as ﬂaw depth) are computed
recursively. This formulation also resembles a target tracking problem with state
transition and measurement models. The formulation can be extended to ﬂaw proﬁling
in any number of dimensions. This reformulation facilitates the application of nonlinear
ﬁltering tools based on sequential Markov Chain Monte Carlo methods known as
particle ﬁlters (PF). The recursive nature of the solution addresses the computational
issues inherent in conventional solutions. The formulation also allows considerable
ﬂexibility in the choice of measurement models, and an assessment of the best
measurement model (out of a given set of potential models) in terms of accuracy and
computational efﬁciency is also conducted. The proposed inversion framework has also

been modiﬁed to fuse data from complementary measurement modes. Principal

Component Analysis (PCA) is used with the modiﬁed framework to further improve
solution accuracy. While the focus of this dissertation was on NDE inverse problems,
the proposed fusion technique can be applied for fusing information in any sampling
importance resampling algorithm based particle ﬁltering application.

The proposed inversion algorithm is applied to a diverse set of simulated and
experimental NDE measurement data. The performance of the algorithm on these data
sets is characterized, and studies conducted to determine the effect of the algorithm
parameters. The performance of the proposed inversion algorithm is also characterized

using conﬁdence metrics such as Cramer Rao lower bounds and conﬁdence intervals.

ACKNOWLEDGMENTS

I wish to express my sincere gratitude to my advisor Dr. Pradeep Ramuhalli and
my committee members Dr. Lalita Udpa, Dr. Sarat Dass, and Dr. Shanker
Balasubramaniam for their unfailing support in my PhD study. Their advice on
technical matters, encouragement, motivation and constructive criticism were greatly
acknowledged.

My studies were partly supported by National University of Science and
Technology (NUST), Pakistan. This support is gratefully acknowledged. I would also
like to thank Department of Electrical and Computer Engineering, Michigan State
University for their support in time of need. I am also grateful to College of
Engineering, Michigan State University which has provided the best research resources
and facilities. I am also grateful to MSU Graduate School for awarding me dissertation
completion fellowship for fall 2009. I wish to thank the entire staff of Electrical and
Computer engineering, Division of Engineering Computing Services and Graduate
School for their valuable timely assistance. My special thanks to nondestructive
evaluation lab members for their support.

Finally, I wish to express my heartfelt gratitude to my parents and my in laws for
the constant support and encouragement. Words cannot express how thankful I am to
my dearest wife for her love and care. I would also like to thank my daughters for love
and understanding throughout my PhD studies. Without their support this dissertation

would not have been possible.

iv

TABLE OF CONTENTS

 

 

 

 

LIST OF TABLES ....................................................................................................... viii
LIST OF FIGURES ........................................................................................................ ix
CHAPTER 1 1
INTRODUCTION ........................................................................................................... 1
1.1 Inverse Problem .......................................................................................................... l
1.2 Overview of Current Solution Approaches to Inverse Problem ................................. 2
1.2.1 Classical Regularization Methods .................................................................... 2

1.2.2 Statistical Approaches ...................................................................................... 3

1.3 Research Objectives .................................................................................................. 5
1.3.1 NDE Inverse Problems ..................................................................................... 5

1.3.2 Objectives ......................................................................................................... 6

1.4 Dissertation layout ..................................................................................................... 7
CHAPTER 2 9
NDE INVERSE PROBLEMS ......................................................................................... 9
2.1 Electromagnetic NDE ................................................................................................ 9
2.1.1 Eddy Current NDE ............................................................ _ ............................... 9

2.1.2 Magnetic Flux Leakage NDE ......................................................................... 10

2.2 Solutions to Inverse Problems in NDE .................................................................... 11
2.2.1 Direct Approaches .......................................................................................... 12

2.2.2 Model-based Approaches ............................................................................... 12

2.2.3 Data Fusion in NDE ....................................................................................... 15
CHAPTER 3 17
A NOVEL FORMULATION OF THE NDE INVERSE PROBLEM .......................... 17
3.1 Problem Formulation ............................................................................................... 17
3.1.1 Statistical Inverse Problem ............................................................................. 17

3.1.2 Formulation in different dimensions .............................................................. 19
CHAPTER 4 22
APPLICATION OF PARTICLE FILTERS TO SOLVE NDE PROBLEMS .............. 22
4.1 The Inverse Problem as a Tracking Problem ......................................................... 22
4.2 Solutions to the Tracking Problem ......................................................................... 23
4.3 Functional Description of Particle Filters ............................................................... 26
4.3.1 Resampling .......................................................... , ......................................... 29

4.3.2 Pictorial Description of particle ﬁltering ....................................................... 30

4.4 Particle ﬁlter Variants .............................................................................................. 31
4.4.1 Sampling importance resampling (SIR) ﬁlter ................................................ 31

4.4.2 Auxiliary SIR (ASIR) Filter ........................................................................... 31

 

 

 

4.4.3 Particle Filters with Improved sample diversity ............................................. 31
4.4.4 Local Linearization Particle Filter .................................................................. 32

4.5 Particle ﬁltering for NDE inverse problems ............................................................ 32
4.5.1 Choice of Importance density ......................................................................... 32
4.5.2 Implementation scheme for NDE inverse problems ...................................... 33

4.6 Models ..................................................................................................................... 37
4.6.1 Measurement Model ....................................................................................... 37
4.6.1.1 Polynomial model ............................................................................... 37

4.6.1.2 Neural Network based models ............................................................ 37

4.6.2 State transition model ......................................................................... ‘ ............ 3 9
CHAPTER 5 40
PARTICLE FILTER BASED DATA FUSION ............................................................ 40
5.1 NDE Data Fusion .................................................................................................... 40
5.2 Data Fusion Using the Particle Filter ...................................................................... 40
5.3 Principal Component Analysis ............................................................................. 43
5.3.1 Theory ............................................................................................................. 43
5.3.2 Implementation ............................................................................................... 46
CHAPTER 6 48
RESULTS ...................................................................................................................... 48
6.1 Evaluation Metrics ................................................................................................... 49
6.2 Eddy Current NDE Data Description ...................................................................... 51
6.2.1 Simulation Data Description .......................................................................... 51
6.2.2 Experimental Data Description ...................................................................... 54

6.3 Comparison of Measurement Models ..................................................................... 57
6.4 Results of Simulated Data ....................................................................................... 62
6.4.1 Simulated Data without Noise ........................................................................ 62
6.4.2 Simulated Data with Noise ............................................................................. 69

6.5 Results with Experimental Data .............................................................................. 74
6.6 2—D Flaw Proﬁling ................................................................................................... 83
6.6.1 Data Description ............................................................................................. 83
6.6.2 Results ............................................................................................................ 85

6.7 Magnetic ﬂux leakage problem ............................................................................... 88
6.7.1 Data Description ............................................................................................. 88
6.7.2 Results ............................................................................................................ 88

6.8 Discussions .............................................................................................................. 91
CHAPTER 7 95
CONFIDENCE METRICS ........................................................................................... 95
7.1 Introduction ............................................................................................................. 95
7.2 Cramer—Rao lower bound (CRLB) ......................................................................... 95
7.2.1 Computation of CRLB ................................................................................... 96

7.3 Credible intervals ..................................................................................................... 97
7.3.1 Computation of credible interval .................................................................... 98

vi

7.4 Results ..................................................................................................................... 98

 

 

7.5 Discussion .............................................................................................................. 104
CHAPTER 8 105
SUMMARY AND FUTURE WORK ......................................................................... 105
8.1 Research Contribution ........................................................................................... 105
7.2 Future Work ........................................................................................................... 106
APPENDIX 107
REFERENCES ................................................................................................................ l 1 1

vii

LIST OF TABLES

Table 6.1 Flaw proﬁles used in comparison of measurement models ................... 58
Table 6.2 Simulated Flaw proﬁles ............................................................ 67
Table 6-3 Inversion results for simulated non-noisy data .................................. 69
Table 6-4 Inversion results for measurement data at SNR=5 ............................. 71
Table 6-5 Inversion results for measurement data at SNR =2 ............................. 73

Table 6-6 MSE between actual proﬁle and PME/MAP estimates using different

measurement modes as shown in ﬁgure 6-20 ................................................ 74
Table 6-7 Experimental ﬂaw proﬁles .......................................................... 76
Table 6-8 Inversion results for experimental data ........................................... 82
Table 6-9 Comparison of NN and particle ﬁlter based inversion techniques ............ 83
Table 6-10 Inversion results (Section C) ...................................................... 87
Table 6-11. Inversion results (Section D) ..................................................... 87
Table 6-12. Flaws used for MFL signal inversion ........................................... 88

Table 6-13. MSE of the predicted proﬁles using MAP estimates as a function of SNR
and L for MFL test data ......................................................................... 91

viii

LIST OF FIGURES

Images in this thesis/dissertation are presented in color.

Figure 1-1. Schematic representation of inverse problem .................................... 1
Figure 1-2. NDE inverse problem ............................................................... 6
Figure 2-1 Phenomological methods ........................................................... 14

Figure 3-1 Representation of NDE inverse problem l-D ﬂaw proﬁle reconstruction
(assumed L=0) .................................................................................... 20

Figure 3-2 Representation of NDE inverse problem 2-D ﬂaw proﬁle reconstruction

(assumed L=O) ..................................................................................... 21
Figure 4—1. Graphical representation of particle ﬁlter (with only 10 samples) ........... 31
Figure 4-2 Particle ﬁlter implementation for ﬂaw proﬁling ................................ 34
Figure 4-3. Implementation of particle ﬁlter ................................................. 36
Figure 5-1. Overall particle ﬁlter based data fusion algorithm ............................ 42
Figure 5-2 Pictorial representation of PCA analysis on toy data 45

Figure 6-1 (a) Top view of posterior PDFs (b) 3-D view of posterior PDF (0) Posterior

PDF at one position index k (d) Posterior PDF estimates ................................. 50
Figure 6-2 Eddy current measurements for a simulated ﬂaw .............................. 52
Figure 6-3. Experimental ECT ................................................................. 54
Figure 6-4 Eddy current measurements at 300 kHz for an experimental ﬂaw .......... 56
Figure 6-5 Illustration of Data registration ................................................... 57

ix

Figure 6-6. Comparison of measurement models ........................................... 59

Figure 6-7. MSE between PME and actual proﬁle versus SNR of measurement data

with different measurement models ............................................................. 60
Figure 6-8. Mean and standard deviations of MSE (between PME and true proﬁle)

versus SNR using different measurement models ............................................. 61
Figure 6-9. Proposed inversion technique ...................................................................... 63

Figure 6-10. Proposed inversion technique (a) Measurement data at multiple
frequencies (b) 3 D view of posterior PDF (c) top view of posterior PDF ((1) computed
estimates ............................................................................................ 65

Figure 6-11. Computational cost versus state vector length parameter L (size of
contributing neighborhood) ....................................................................... 66

Figure 6-12. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' at different measurement modes ......................................... 68

Figure 6-13. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L at different measurement modes ...................................................... 69

Figure 6-14. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' at different measurement modes with additive noise at SNR=5. . .70

Figure 6-15. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L at different measurement modes with additive noise at SNR=5 .................. 71

Figure 6-16. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' at different measurement modes with additive noise at SNR=2....72

Figure 6-17. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L at different measurement modes with additive noise at SNR=2 .................. 73

Figure 6-18 Proposed inversion technique using magnitude at 300 kHz ................... 75

Figure 6-19. Proposed inversion technique (a) Phase (b) Magnitude (c) 3 D view of
posterior PDF ((1) top view of posterior PDF (e) computed estimates ..................... 76

Figure 6-20. Inversion results using measurements at single and multiple frequencies.77

Figure 6-21. Two experimental ﬂaws .......................................................... 80

Figure 6—22. MSE between MAP of posterior PDF and true proﬁles versus state vector
length parameter 'L' at different measurement modes for experimental data .............. 81

Figure 6-23. Mean and standard deviations of MSE (between MAP and true proﬁle)

versus L at different measurement modes ....................................................... 82
Figure 6-24. Aircraft skin specimen .............................................................. 84
Figure 6-25. Inversion results —Section C (Single measurement modes) .................. 86
Figure 6-26. Inversion results-Section C (Data Fusion) ..................................... 87
Figure 6-27. Inversion results-Section D (Data Fusion) ..................................... 87

Figure 6—28. Posterior PDF estimates (L=2) for ﬂaw proﬁle from MFL test data (a) 2-D
view of PDFs (b) 3-D view of PDFs (c) Estimates evaluated from PDFs. . ................89

Figure 6-29. MAP estimates evaluated with MFL test data at different SNR levels. . . ..90

Figure 7-1. Credible intervals on posterior PDF ............................................... 98

Figure 7-2. Credible intervals .................................................................... 99

Figure 7-3. Covariance of state and CRLB for single and multiple measurement

modes ............... _ ............................................................................... 100
Figure 7-4. Credible intervals ................................................................... 100
Figure 7—5. Covariance of state and CRLB for single and multiple modes .............. 101

xi

Figure 7-6. Credible intervals ................................................................... 102

Figure 7-7. Covariance of state and CRLB for single and multiple modes .............. 102

Figure 7-8. Credible intervals .................................................................. 103

Figure 7-9. Covariance of state and CRLB for single and multiple modes. . . . .......103

xii

CHAPTER 1
INTRODUCTION

1.1 Inverse Problem

An inverse problem typically involves determination of model parameters from
observed data (measurements) [1]. This generic problem occurs in many branches of
science and engineering, such as geophysics, medical imaging, remote sensing, ocean
acoustic tomography, acoustic holography, nondestructive testing and astronomy. The
typical inverse problem may schematically be represented in Figure 1-1.

In accordance with convention [2], the collection of values that are to be

reconstructed is termed as the image. The image is usually denoted by f . The set of all

images is called the image space. The forward problem or the direct problem is the
mapping from the image to the measured quantities. In practice, exact measurements
cannot be acquired and the measured data d is a corrupted version of the exact
quantities, i.e., the measured data comprises of error-free data and noise n. Therefore,

the forward process is a mapping from the image to the measured data.

Forward Problem

 

 

Model Parameters: f Measurement Data: d

 

 

 

 

 

 

Inverse Problem

Figure 1-1. Schematic representation of inverse problem

If the mapping from image to actual data is assumed as a linear process, then the
mapping is given by the following relation [2]:

d = Af + n (1.1)

Here A is called the system matrix. The inverse problem is then the problem of

ﬁnding f given the data d and knowledge of the forward problem. Evaluation of

f given dis a well posed inverse problem if a unique solution exists for each

measurement [2]. In general, this is not the case, and the inverse problem is usually an
ill-posed problem. Hadamard [3] deﬁned the ill—posed problem as that where one or

more of the following conditions exist:
(i) An inverse does not exist because the data is outside the range of the
system matrix.
(ii) The inverse is not unique because more than one solution is mapped
to the same data.
(iii) An arbitrarily small change in the data can cause an arbitrarily large

change in the solution [3].
1.2 Overview of Current Solution Approaches to Inverse Problem
1.2.1 Classical Regularization Methods

The most common techniques to solve inverse problems are regularization
methods. Regularization methods constitute mathematical techniques which have been

developed to address the ill posedness [2, 4]. The key idea behind regularization is to

seek a new problem which is very much similar to the actual problem, such that the new
problem is uniquely solvable and robust in the sense that small errors in the data do not
corrupt the solution. Regularization involves introducing additional information in order
to solve an ill-posed problem. This information is usually in the form of a penalty, such
as restrictions on smoothness or bounds on the norm of the solution. The least-squares
method [4] can be viewed as a very simple form of regularization. Another simple form
of regularization applied to integral equations is termed as Tikhonov regularization [2,
4]. Tikhonov regularization is a trade-off between ﬁtting the data and reducing the norm
of the solution at the same time. A non-linear regularization method, such as total
variation (TV), is used if data is discontinuous or unsmoothed [4]. Truncated iterative
regularization techniques [2, 4] using Landweber-Fridman, gradient descent and
conjugate gradient algorithms transform the inverse problem into an optimization
problem [2]. The iterative algorithms attempt to solve the inverse problem by reducing
the error functional between the actual and computed solution. Genetic algorithm (GA)
based techniques [5] have also been applied to solve inverse problems. However, most
of these classical methods give inaccurate results when the measurement data is noisy
or the system matrix is ill-conditioned. Moreover, iterative techniques and GA based

inversion techniques offer computationally expensive solutions to the inverse problem.

1.2.2 Statistical Approaches

In a typical inverse problem, there are observable parameters (data) as well as
unknown parameters. The variables are related to each other using models. In most
inverse problems, some prior information about the unknown parameter is also

available. The objective of statistical approaches is to evaluate the available information

while accounting for uncertainty about the variables of interest. The objective is
achieved using all available knowledge of the models (measurement process) and prior
information about the unknowns. The solution to inverse problems using statistical
methods produces a probability distribution function, in contrast to regularization
methods which produce a single estimate of the solution. The resultant probability
distribution is known as the posterior probability [7]. The statistical inversion process is
generally divided into three subtasks [6]. The ﬁrst subtask is to evaluate the prior
probability of the unknown. The next subtask is to determine the likelihood function
which describes the relationship between the measurement and the unknown quantity.
The ﬁnal subtask is to explore the posterior probability distribution of the unknown
quantity. To date various statistical approaches have been adopted to solve NDE
inverse problems [2]. However, the posterior probability distribution as a solution is
generally useless unless ways of exploring this distribution are developed. Therefore,
single estimator evaluation approaches, such as maximum a posteriori (MAP) estimates
and conditional mean/co-variance estimates are used. Evaluation of MAP estimates of
the posterior distribution leads to possible non-existence and non—uniqueness problems.
Conditional mean and conditional covariance techniques lead to integration problems in
multi-dimensional spaces which can be computationally very expensive [6]. Another
major issue in analyzing inverse problems probabilistically is that the relationship
between measurement data and model parameters is generally nonlinear [7]. Non-
linearity makes the description of posterior probability difﬁcult due to possible

multimodality or non-availability of some moments.

To address these drawbacks an alternative approach in the paradigm of
statistical inverse problems can also be used. In this alternative approach, the posterior
density can be used to generate a set of points which supports the distributionfMarkov
Chain Monte Carlo (MCMC) [6—8] methods offer the required sets of points (samples).
MCMC methods are based on constructing a Markov chain [9] that has the desired
distribution as its equilibrium distribution. The state of the chain after a large number of
steps is then used as a sample from the desired distribution, i.e., the solution to the
inverse problem. The accuracy of the solution to inverse problems using this method
improves as the number of steps increases. Therefore, Monte Carlo based methods in
this form also offer a computationally expensive solution if highly accurate inversion is
required [6, 8]. Moreover, in case of high dimensionality of unknown variables, the
complexities of inverse problem as well as the computational cost increase

considerably.
1.3 Research Objectives

1.3.1 NDE Inverse Problems

Nondestructive Evaluation (NDE) is an interdisciplinary ﬁeld of study which is
concerned with the development of measurement technologies and analysis techniques
for the quantitative characterization of materials and structures by noninvasive means
[10]. Ultrasonic, radiographic, thermo-graphic, electromagnetic, and optic methods
have been employed to probe the interior rnicrostructure of materials for subsequent
characterization of subsurface features [10]. Applications of NDE inverse problems are

in non-invasive medical diagnosis, intelligent robotics, security screening, and on-line

manufacturing process control, as well as the traditional NDE areas of ﬂaw detection,

structural health monitoring, and materials characterization.

-- Lilli”; 'il'lfi. ..--~
' .i ii."
'l'qi My“.

Input Space (Measurement Signal) Output Space (Defect Proﬁle)

   

Figure 1-2. NDE inverse problem

Determination of ﬂaw from NDE measurements, as shown in Figure 1-2, is the
most common inverse problem in NDE. This inverse problem is ill-posed due to non-
uniqueness and instability of the estimated solution from the noisy measurement data.
Signiﬁcant research work, using both classical and statistical approaches, has been
undertaken to date to combat the ill-posedness present in these inverse problems [10].
However, novel inversion techniques are necessary to address the drawbacks of existing

approaches.

1.3.2 Objectives
This PhD dissertation proposes a novel approach for solving NDE inverse
problems that addresses the drawbacks of current techniques. The objectives of this

PhD research are:

(i)

(ii)

(iii)

(W)

(V)

Develop a recursive framework for solving NDE inverse problems.
Such an approach would address the computational issues inherent in
conventional solutions.

Apply nonlinear ﬁltering tools based on Markov Chain Monte Carlo
methods to obtain a low cost, high accuracy solution to NDE inverse
problems.

Extend the recursive nonlinear ﬁltering solution technique to fuse
information from multiple NDE measurement modes, and quantify
any resulting improvements in solution accuracy.

Characterize performance of the proposed algorithm using conﬁdence
metrics such as Cramer Rao Lower Bounds (CRLB) and conﬁdence
intervals. Evaluate the effect of data fusion on the conﬁdence metrics.
Evaluate the feasibility of the proposed algorithm using a wide
variety of NDE data, including simulated data with and without added
noise, as well as experimental data. Using the results on the different
data sets, quantify the accuracy, robustness to noise and

computational efﬁciency of the proposed algorithm.

1.4 Dissertation layout

Chapter 2 of this PhD dissertation describes the electromagnetic NDE inverse

problem followed by a literature survey of various techniques that have been applied to

solve the NDE inverse problem. Chapter 3 proposes a new formulation of inverse

problems that admits a recursive solution algorithm. Chapter 4 introduces the particle

ﬁltering technique and discusses its application for ﬂaw proﬁle reconstruction, using the

proposed problem formulation. Chapter 5 discusses the extension of the proposed
solution technique to fusing multimodal measurement data. Results of applying the
proposed solution technique to simulation and experimental NDE data are presented in
Chapter 6. Chapter 7 presents a study quantifying the efﬁciency and accuracy of the
proposed technique with and without data fusion. Finally, conclusions and

recommendations for future work are presented in Chapter 8.

CHAPTER 2
NDE INVERSE PROBLEMS

2.1 Electromagnetic NDE

All electromagnetic methods of non-destructive evaluation rely on the
interaction between a source ﬁeld and the material under test [10, ll]. Electromagnetic
NDE methods range from the magnetic ﬂux leakage technique to higher frequency
techniques such as eddy current and microwave NDE [11]. However, the focus of this

dissertation is restricted to magnetostatic and eddy current methods of NDE.

2.1.1 Eddy Current NDE

Eddy current testing is a popular and extensively used technique for anomaly
detection in conducting materials. Eddy current testing techniques are frequently used
in nuclear aerospace, marine, high pressure, high temperature and high speed
engineering systems where frequent inspection of systems is mandatory to avoid
catastrophic damages to systems and human lives [10]. These methods are suitable for
inspection of engines, nuclear reactor systems and machine parts. The technique is
based on the analysis of variation of the magnetic ﬁeld linked to an eddy current
perturbation [11, 12]. Eddy currents are induced within the specimen by a time-
harrnonic exciting source (coil) displaced along the length of the specimen. According
to Lenz’s law, a magnetic ﬂux is created in reaction to the ﬂow of eddy currents, and is

observed as a change of the coil impedance. When a defect is present, the eddy current

pattern is perturbed, which results in further change of the coil impedance. Changes in
impedance of the coil, as it scans the surface of the specimen, are used as a basis for
detecting the presence of discontinuities in the specimen. For eddy current NDE, the

governing equation is derived in terms of the magnetic vector potential A as [12]

c ” _. (2.1)
l-VX(VXA) = —oa—A+ Js
,u a:

where J s is the source current density and 0’ is the conductivity of the material. Under

sinusoidal steady state conditions the magnetic vector potential is a phasor quantity and

its derivative with respect to time can be written as

_ g (2.2)
arm,
at

to give the governing equation for eddy current NDE as [12]:

1

_. _. _. (2.3)
—VX(VXA) = -ijA + J
,u

where [J is the permeability, J is the source current density and 0' is the conductivity

of the material and A is the magnetic vector potential.

2.1.2 Magnetic Flux Leakage NDE

The magnetic ﬂux leakage method, frequently used in the inspection of
ferromagnetic materials, employs permanent magnets or currents to magnetize the
sample and a set of ﬂux-sensitive sensors to record the leakage ﬂux for analysis.
Leakage ﬂux arises because the presence of a defect causes an increase in the magnetic

ﬂux density in the vicinity of the ﬂaw. This causes a shift in the operating point on the

10

hysteresis curve and a corresponding decrease of local permeability, resulting in a
leakage of ﬂux into the surrounding medium (air) [12]. The leakage ﬂux is typically
recorded using Hall probes.

The governing equations for MFL can be derived from the static form of
Maxwell’s equations. Considering only source-free regions of the material and

assuming the region is homogeneous and isotropic, we can show that [12]

- V.(,u.V¢) = —yV2¢ = 0 (2.4)

where (p is the magnetic scalar potential, ,u is the permeability and w is related to the

magnetic ﬁeld strength H as

H = —V¢ (25)

The solution to (2.5) for realistic inspection geometries with complex defect shapes
necessitates the use of numerical techniques such as a ﬁnite difference or ﬁnite element
model (FEM). Considerable work has been done in the development of these models
[13- 14] for predicting all three components of the magnetic leakage ﬂux density B as a

function of spatial location.

2.2 Solutions to Inverse Problems in NDE

A signiﬁcant amount of research work has been carried out till date to address
the NDE inverse problem. The solution approaches to NDE inverse problems can be

divided into two broad categories.

11

2.2.1 Direct Approaches

Direct, or non-phenomenological, approaches typically rely on the use of signal
processing techniques to establish a relationship (mapping) between speciﬁc
characteristics of the signal and the geometry of the defect, ignoring the underlying
physical process. These methods typically pose the inverse problem as determining a
mapping from the measurement space to the material property space [15-16] where the
set of unknown parameters that deﬁne the mapping are determined from measurements.
Direct approaches ranging from calibration methods to more recent procedures based on
neural networks [17-18] have all been proposed. Another study [19] relates the
mechanical stresses with electromagnetic properties of deﬁned materials to reconstruct
electromagnetic maps using support vector regression machines (SVRM). However,
these techniques are limited in that they require data from known ﬂaws (training data)
in order to determine the mapping parameters, and can only be used 'when
measurements are acquired from ﬂaws similar to those used in the training process. The
advantages of this approach are their simplicity and speed; however, the approaches are
also very sensitive to the analytical models that are used as well as noise in the

measurements.

2.2.2 Model-based Approaches

Model based, or phenomenological, methods usually rely on a physical model to
accurately simulate the underlying physical phenomenon and predict the probe response
[20]. The model is used to estimate the measurement given the ﬂaw proﬁle, which is
iteratively derived by minimizing the difference between the estimated and actual

measurements (ﬁgure 2-1). The optimization schemes generally used in these studies

12

utilize gradient search methods to take advantage of their rapid convergence rates.
However, such methods often fail to converge to a global optimum due to the presence
of multiple local optima in the objective function. Minimization may be through
conventional techniques such as conjugate gradient, though other techniques such as
simulated annealing [21] or genetic algorithms [22] have also been proposed. Model-
based approaches usually involve signiﬁcant computational effort, since the physical
model needs to be solved repeatedly unless a database containing pre-computed
defect/signal pairs is available.

Key to this approach is a forward model that accurately represents the physics of
the imaging process. Numerical models such as a ﬁnite element model or integral
equation models [23-36] have been proposed for electromagnetic NDE signal inversion,
and though accurate, these tend to be computationally expensive since the models must
be solved during each iteration. A typical approach uses a ﬁnite-element forward model
and a gradient-based optimization algorithm, to invert electromagnetic data is given in
[23]. Hoole et al. [24] applied a similar approach for estimating the geometry of cracks
from electromagnetic ﬁeld measurements. Bowler [35-36] presented a fast integral
equation based model to calculate eddy-current probe responses together with gradient-
based optimization algorithms to reconstruct conductivity distributions. One of the ﬁrst
3D reconstructions is presented in reference [30], where the rrrinirnization procedure
works by solving a sequence of nrinimizations in a continuous space. The optimization
schemes generally used in these studies utilize gradient search methods to take
advantage of their rapid convergence rates. However, such methods often fail to

converge to a global optimum due to the presence of multiple local optima in the

13

objective function. To overcome this limitation, optimization methods based on genetic
algorithms [37-38] and other evolutionary strategies [39] have been proposed. A
different approach that is based on an analogy between crack proﬁle estimation in NDE
and digital communications over noisy channels has been presented in [40]. The
problem is formulated as a shortest path search on a trellis using the Viterbi algorithm
(VA), and achieves the global minimum at the price of computational complexity.

The primary drawback to these techniques is a very high computational cost due
to the need to solve the direct problem several times, particularly when numerical
models are used to solve the direct problem. Solutions that use a database containing
pre-computed defect/signal pairs have also been proposed; however they tend to be not
as accurate or applicable in general. The computational cost is exacerbated when
genetic algorithms are used, since these techniques require several direct solutions in
each iteration. Alternative forward-model based inversion techniques that are not

iterative have been proposed primarily for inverse microwave and acoustic scattering.

Experimental
Input Signal

Yes

   

 

 

 

 

 

 

 

 

 

 

 

 

 

Initial defect.
Pmﬁle ' Forward
#7 Model — < a (tolerance
V
N Desired
Update defect O l Defect Proﬁle

 

 

 

proﬁle

 

Figure 2-1. Phenomenological methods for solving inverse problems

14

These include the tomographic reconstruction using the Fourier Slice theorem [26], [41-
42], point source technique [43] as well as linear sampling [44-45], which uses far ﬁeld
measurements to solve the inverse scattering problem. Constraints on material
properties in the form of Markov random ﬁelds [46], as well as other regularization
methods have also been attempted [44, 47- 48]. These methods are not recursive (i.e.,
attempt to solve the inverse problem over the entire problem domain using all of the
data) and are generally computationally expensive. A Bayesian technique was also
proposed for defect signal analysis in NDE images [49]. The authors in [49-50]
designed a hierarchical Bayesian framework for detecting and estimating NDE defect

signals from noisy measurements.

2.2.3 Data Fusion in NDE

All the approaches discussed above use data from a single measurement mode to
solve the inverse problem. The use of measurement data from multiple NDE
measurement modes to solve the inverse problem is known as data fusion. Many
approaches to fusing the information in multiple NDE measurements have been
presented [51-52]. Commonly proposed solutions include neural networks [53-56],
Bayesian analysis based on Dempster-Shafer evidence theory [57-58], wavelet and
other multi-resolution algorithms [67], and image fusion [58-60] in time and frequency
domain. These methods have been applied to fuse NDE data from a range of sources,
including multifrequency eddy current measurements (ECT) [53,54,59], ECT data and
ultrasound measurements [52-54,60-61], ultrasound, x-ray and acoustic emission
measurements [59], and other measurement types (such as pulsed eddy current

measurements) [63-64]. In addition to these conventional fusion techniques, other

15

methods such as the Q-transform based technique [65] have also been recently
investigated.

The major drawback associated with most of the NDE data fusion techniques is
that these techniques provide a fused signal (or image) instead of ﬂaw proﬁle estimates.
Moreover, available fusion algorithms are generally based on processing signals or
images without regard to physics of the measurement process. In view of the
drawbacks, a new technique is necessary.

The discussion above highlights the need for an accurate, robust, and
computationally efﬁcient solution method for electromagnetic NDE inverse problems.
A new inversion technique is proposed in this dissertation which offers the desired

characteristics.

16

CHAPTER 3

A NOVEL FORMULATION OF THE NDE INVERSE PROBLEM

As discussed earlier, there is a need for an accurate, robust, and computationally
efﬁcient solution method for electromagnetic NDE inverse problems. In this chapter,
we discuss a novel formulation of the NDE inverse problem that admits a recursive
solution framework. The advantage of this recursive formulation is a simpler and
computationally efﬁcient implementation of the solution. AS in the previous studies, we
will assume that the goal of this formulation is to determine, from NDE measurements,

the ﬂaw depth (or some similar characteristic of the ﬂaw) at all locations.
3.1 Problem Formulation

3.1.1 Statistical Inverse Problem
Assume that the specimen consists of K discrete locations. The depth proﬁle of
ﬂaws in the specimen is expressed in terms of a

set}? ={x1,x2, ..... xk_1,xk,xk+1,....xK_1,xK}, where the elements of the set are

states (flaw depths) at discrete locations. The corresponding measurement set is
denoted byZ={Z1,z2, ..... ,zk_1,zk,zk+1,...zK_1,zK}, where each element zk is
the measurement at the discrete locationk. It is assumed that N k locations are in the
neighborhood of any position indexk. Two functions are deﬁned: h(xk , Vk) known

as the measurement function relates the state xk to the corresponding measurement Zk ,

and f (x j I 11.1sz ,7”) known as state transition function relates state xk to
j¢k

l7

states in the neighborhood i.e. x j I 151sz . Vk and 77k are terms that represent the
j¢k
error (or noise) in the state transition and measurement functions respectively.

The evaluation of the sequence of states X given the measurements Z is the
inverse problem. This inverse problem can be formulated in terms of a statistical

estimation problem due to the noise in the state transition and measurement processes.

—9

In this case, the unknown parameter is the set of states X while the corresponding

observation is the set of measurementsZ . The posterior PDF p()? I Z) of the states

X is computed in statistical inverse problems [6]. The states X are estimated from the

posterior PDF. The desired posterior PDF is given by [6]:
190? like p(Z|)?).p(X') (3.1)

Here p(Z I I? ) is the likelihood function and p0? ) is the prior information of the
states (set of ﬂaw depths). The estimation of pH? I Z) is a computationally expensive

process due to the high dimensionality of the state vector)? . Therefore, the posterior
PDF needs to be evaluated sequentially (i.e., evaluated at each position by stepping
through the position index k) to keep the computation cost and complexity of the
inverse problem low. Three conditions are assumed to enable the sequential solution to
the inverse problem [8, 51]. Firstly a locally independent Markov ﬁeld is adopted to

solve the problem. Therefore, the relationship offered by (3.1) holds locally and the

equation can also be written for a neighborhood around the position indexk , which

runs across the whole region of interest (1 I K ).

p(xk lzk) °c p(Zk lxk).p(xk) (3.2)

18

The second assumption is that the observations (measurements) are mutually

independent within the neighborhood of the state given true values of the unknown

states. Thus,

Nk (3-3)
pm. lxk)= II P(Zj Ixj)

i=1

Finally, the prior density of the unknown state is assumed to be a product of exponential

densities centered on the value of the states in the neighborhood.
(3.4)

1” ”xi 'xk IIP
_ z , ,

_ ’=1 P
P(kuijj=1:Nk)—e J .
j¢k

 

where p is a scalar, with value selected to be around 1 to ensure a low amount of

variability in the state value within the neighborhood.

3.1.2 Formulation in different dimensions

The formulation can be easily extended to ﬂaw proﬁling in any number of
dimensions. The problem domain (NDE specimen) can be divided into K discrete
locations with a neighborhood deﬁned for each spatial position. In two spatial
dimensions, the neighborhood may be deﬁned using a lexicographic arrangement of
pixels commonly used in the image processing literature [66]. Examples of the
formulation for 1-D (crack-like) and 2-D (volumetric ﬂaws such as corrosion, wear or

wall thickness loss) proﬁling are Shown in ﬁgures 3-1 and 3-2 respectively.

19

Measurements
zk -1 Zk zk+1

     
 

     
    
 

o
o
.0

h(Xk+1, Vk+1) 1] ﬂ Uh(xk+lvvk+l)
I ’ ~ \ \l I ‘ \ \

¢—~‘ 4_~
I
3‘

/\ I

\

\‘~_—” \‘~_—’)"~_.—’”
Region k-l k k+l
State xk—l xk xk+1

Figure 3-1. Representation of NDE inverse problem for 1-D ﬂaw proﬁle reconstruction

(assuming L=0)

20

 

 

Figure 3-2. Representation of NDE inverse problem for 2-D ﬂaw proﬁle reconstruction
(assuming L=0)

As shown in Figure 3-1, the neighborhood deﬁned for a 1-D proﬁle construction
is a 1-D window expressed as xk—(2L+l):k—l,k+l:k+(2L+1)- L is a scalar
parameter which describes the size of the neighborhood. The neighborhood deﬁned for
2-D proﬁle construction is a 2-D window as shown in Figure 3-2, and deﬁned as

xp,q Ip=m—-(2L+1):m—l,m+l:m+(2L+1), . The 2-D problem formulation is similar
q=n—(2L+l):n—1,n+1:n+(2L+1)

to the 2D Ising model [67].

21

CHAPTER 4

APPLICATION OF PARTICLE FILTERS TO SOLVE NDE INVERSE
PROBLEMS

4.1 The Inverse Problem as a Tracking Problem

As discussed earlier, the inverse problem in NDE is to determine the best ﬂaw
characteristics that match the measurements. Given the representation in Chapter 3, the
problem of ﬂaw characterization from measurements may also be formulated as a
problem of determining the best sequence of states from the measurements. In
particular, this problem may be modeled by using two different equations [68] - a state
transition equation that is used to model the evolution of the states with respect to
spatial position and a measurement model which relates the measurements to the states
at a given position:

xk =fk(xj |J=1.Nk 9Vk) (4.1)

j¢k
Zk = hk(xk,uk)

where f k and hk are functions that model the state transition process and the

measurement process respectively, and Vk and ,uk represent the process noise and
measurement noise respectively. Note that the state transition and measurement
equations are applied to the local neighborhood of state xk. The problem described by

equation (4.1) is often referred to as a tracking problem, and is used frequently in target

22

tracking applications [68]. In these applications, the state transition function models the

motion of a target from a known position x j , while the measurement function

describes some function of the target position. Process and measurement noise densities
represent the uncertainty in the state and measurement models. The tracking problem as
deﬁned in equation (4.1) is a dynamic state estimation problem [68]. The Bayesian
approach to this problem is to construct a posterior PDF of the state, based on the
sequence of measurements.

The problem described by equation (4.1) can be shown to be equivalent to the
statistical inverse problem deﬁned by equation (3.2) by the following arguments. The
state transition function (when taken with the associated process noise distribution) is
equivalent to the local prior PDF used in the statistical inverse problem. Similarly the
measurement function is equivalent to the local likelihood PDF in the statistical inverse
problem when a locally independent Markov ﬁeld is assumed (sec. 3.1.1). The solution
to the reformulated NDE inverse problem as described in section 3.1 is also based on
computation of the posterior PDF of the ﬂaw characteristic at a spatial location given
the likelihood and prior PDF based on the local neighborhood of the spatial location.
The functions are applied on the ﬂaw proﬁle characteristics within the local
neighborhood of each spatial location. Therefore, the reformulated NDE inverse
problem, as deﬁned by equation (3.1), is equivalent to a sequential tracking problem,
with corresponding state transition and measurement functions. The advantage of this
equivalence is the potential applicability of solution techniques for tracking problems to

solve the NDE inverse problem.

23

Note also that the formulation of the inverse problem as statistical inverse
problem (3.2) or equivalently, as a tracking problem (4.1), implicitly assume that the

NDE measurement process is a localized process, i.e., the measurement at a particular

location k is only affected by the state of the sample in a neighborhood of k . This
Markovian assumption is valid for low-frequency electromagnetic NDE (magnetostatic

and eddy current NDE) [40].

4.2 Solutions to the Tracking Problem

Kalman ﬁltering [68-69] provides an optimal solution to the tracking problem if
the following two conditions are met:

(i) The function which relates the states in a neighborhood (i.e., the state
transition function) and the function which relates the state to
measurement (i.e., the measurement function) are linear.

(ii) The likelihood and prior PDFS are Gaussian.

In general, for the NDE inverse problem, these functions can be non-linear and
the PDFS can be non Gaussian (multi-modal). Therefore, the Kalman ﬁlter cannot
provide an optimal solution and suboptimal algorithms may be necessary to evaluate the
ﬂaw depth.

Suboptimal algorithms [69] include the extended Kalman ﬁlter (EKF) and
Particle ﬁlters. In EKF, the state transition and measurement equations (4.1) are locally
linearized. However the sate transition and measurement noise PDFS are assumed to be
Gaussian. An extension of EKF is the unscented Kalman ﬁlter (UKF) [69]. UKF

considers a set of points that are deterrrrinistically selected from the Gaussian

24

approximation to posterior PDF p(xk I z k ). These selected points are all propagated

through the true nonlinearity. The UKF has been shown to give better performance than
a standard EKF for some problems, as it better approximates the nonlinearity. However,
in UKF the PDFS are again assumed to be Gaussian. If the true density is non-Gaussian
(e.g., if it is bimodal or heavily skewed), then a Gaussian can never describe it well. So
the EKF or UKF cannot solve this tracking problem accurately. Therefore, a more
generalized ﬁltering technique is required to solve this inverse problem. Particle ﬁlters
offer such a generalized ﬁltering technique.

Particle ﬁlters are sequential Monte Carlo methods based on point mass (or
“particle”) representations of probability densities that can be applied to any state-space
model and which generalize traditional Kalman ﬁltering methods [68-69]. In this
approach to dynamic state estimation, the posterior probability density function (PDF)
of the state is constructed based on all available information, including the set of
received measurements. Since this PDF embodies all available statistical information, it

may be said to be the complete solution to the estimation problem. In principle, an

optimal estimate of the state may be obtained from the PDF. The PDF of the state xk

conditioned on all measurements up to (and includingzk), p(xk,zlzk)may be

obtained recursively in two stages: prediction and update. The prediction stage uses the

system model to predict the PDF forward from one measurement location to the next.
Suppose that the required PDF p(xk I 21; k— 1 ) at location k -1 is available. The

prediction stage involves using the system model to obtain the prediction density of the

state at k via the Chapman-Kolmogorov equation [68]:

25

P(xk I lek—l) = IP(xk I xk-1)P(xk—1 I Z1;k—1)dxk_1. (4-1)

If we consider a Markov process of order one then
p(xk ka_1,z1:k_1)= p(xk ka_1). Since the state is usually subject to
unknown disturbances (modeled as random noise), the prediction step generally
translates, deforms, and otherwise distorts the PDF. The update operation uses the latest
measurement to modify the prediction PDF. This is achieved using Bayes’ theorem as
follows:

P(Zk karzlzk—1)P(xk IZ1;k_1) (4-2)
P(Zk erzk—r)

 

P(xk IZk) = PW: | zkrzlzk—1)=

= P(Zk ka)P(xk I Zl:k-1)
P(Zk I21:k—1)

 

9

where the normalizing constant is given by:

( I )=l( IXMX’I - )dx <43)
[9 zk z1:k-1 P Zk P zl.k—1 k -
The normalization constant depends on the likelihood function p(z k I xk ) , deﬁned

by the measurement model.

The key idea behind the particle ﬁlter is the representation of the required
posterior density function p(xk Izk) by a weighted set of samples and the

computation of estimates based on these samples and weights. As the number of
samples becomes very large, this Monte Carlo characterization becomes an equivalent

representation to the usual functional description of the posterior PDF and the particle

26

ﬁlter approaches the optimal Bayesian estimator [68-69]. A functional description of the

technique is given below.

4.3 Functional Description of Particle Filters

A simple algorithm for implementing the particle ﬁlter is known as the
sequential importance sampling (SIS) algorithm [68-71]. In order to apply particle
ﬁltering, the desired posterior PDF is represented in terms of samples and associated
weights at each location. Let the posterior PDF at position k be given by (4.4)

N, l, l, (4.4)
p(xk IZk) z Zwk§(xk -xk).

I:
Here, x2 are the samples (or particles) used to represent the posterior density,
i = 1 I N s is the total number of samples used and WIC is the weight associated with

sample x2. The samples are drawn from prior distribution given by (3.4). As discussed

above, the prior distribution depends upon the value of the parameter in the

neighborhood of. locationk. Normalized weights are chosen using the principle of

importance sampling [68-74]. If the samples x2 were drawn from an importance

density, q(x1: k I 2:1: k ), the weights are given by

' (4.5)
. (x'. z . )

(“xi/C I zlzk)

 

Suppose at position k -l we have samples constituting an approximation

to p(x1:k_1 I lek-1)- With the reception of measurement z k at positionk , we wish

to approximate p(x1: k I 21; k )with a new set of samples. If the importance density is

chosen to factorize such that

27

61061:]: | lek) = q(xk IXIzk-1,Z1:k)q(x1:k—1 I lek—l)’ (46)

we can obtain samples xi, k =3 q(x1:k Izlzk)by augmenting each of the existing
samples xI‘k—l z q(x1:k_1 I z1:k—1) with the new
statex;c z q(xk |x1;k_1, lek)- To evaluate the weight update equation, the PDF
p(x1:k I lek) is ﬁrst expressed in terms of, p(x1:k_1 I lek—1)r p(zk ka) and

P(xk ka—1)=
p(xl:k I zlzk) °c p(Zk ka)P(xk ka_1)p(x1;k_1 I lek—l) (4-7)

Given the set of weights Wk—l at position k — 1, the weights at position k may be

computed recursively using the weight update equation derived from the principle of
importance sampling as

' ' ' 4.8
i i P(Zklx2)p(x;(|x2_l) ( )

W = W
k k— - -
q(x;C Ix;(_1,zk)

 

Now, if q(xk Ix1:k_1,zk)=q(xk ka_1,zk), then the importance density

becomes only dependent on xk_1 and zk. This is particularly useful when

p(xk I lek) is required at each position index. In such scenarios, only x2 need be

stored, and the path x1:k_1 and the history of observations Z1: k—l can be discarded.

The modiﬁed weight is then

- - ° 4.9
i i P(Zk|xI)P(xICIXIC_1) ( )

W 0‘ W
k k-l ' ‘
q(x;( Ix;(_1, Zk)

 

28

4.3.1 Resampling

The weights are assigned to each sample at each step in the recursive estimation
procedure. It can be shown that the variance of the importance weights can only
increase at each step [73] for the importance function of the form given by (4.6). After
several steps, all but one sample (particle) will have negligible normalized weights. A
lot of computational effort is then wasted working with samples with negligible
weights. This phenomenon is known as degeneracy. It is impossible to avoid

degeneracy in the SIS framework [73]. A measure of degeneracy is effective sample
size N cf?“ [75]:

_ 1 (4.10)
2.41M.) '

 

N eﬂ
where 15 N eﬁ S. N S . If the weights are uniform then N eff = N s . However, if for

i=1: NS,W]‘(i =1 and WI; Ii¢j= O for some value ofj, then Neﬁr =1. Thus, a
small value of N eff indicates severe degeneracy.

In order to avoid degeneracy (where all samples except one have negligible
weight), re-sampling [66, 70] of ka ’Wk III-=1st is performed. Resampling

introduces additional samples with higher weights and discards samples with low
weights. In a direct implementation of resampling, N s i.i.d variables are generated

from a uniform distribution. These variables are sorted in ascending order and then

compared with the cumulative sum of normalized weights (CSW) iteratively for

29

i = 1 2 N s . Only those samples x2 are assigned weights for which the i.i.d variable is

less than CSW. This comparison ensures a larger number of samples of normalized
weights where the weight assigned to the sample is large. Similarly fewer samples are
produced where the weight assigned to the sample is small. Moreover, samples with

negligible weights are ignored.

4.3.2 Pictorial Description of partlcle fllterlng

A pictorial representation of the operation of a particle ﬁlter is shown in ﬁgure

4-1 [68]. Uniformly weighted samples initially approximate the density

p(xk IZk—l) and the importance weights are calculated using the likelihood

function p(z k I xk ). The result is a random measure of samples and their weights.

The random measure is Ix; , WICI for i th sample, which is an approximation

of p(xk I z k ) . AS shown in the ﬁgure, some samples have large weights while others

have negligible weights. If the weights are assigned in this manner to the sample at each
recursive step, all but one particle will have negligible normalized weights. As
discussed above, in order to avoid degeneracy (where all samples except one have

negligible weight), re-sampling [73, 75] is performed. The resampling step results in a
'4: _

measure, Ix; , N s 1 I which represents a new set of samples with normalized weights.

The last step is the prediction step, which introduces variety (due to process noise) and

results in a measure Ix;c +1,Ns_lI that approximates p(xk+1 I 2k). This measure

approximates the prediction density at the next cycle of the particle ﬁlter.

30

I‘i’NTII OO 0 00000

 

 

 

 

 

 

 

 

P(Zkak)I— / \gJ//\
{XIWII CO 000 °°°°O
O
O
O

i* i 0 O O
k’wkI 0

OO
O

it..~.—1Ic'> 6 0/80 666% c'>

Figure 4-1. Graphical representation of particle ﬁlter (with only 10 samples)

4.4 Particle filter Variants

There are several variants of particle ﬁlters [68]. These variants are based on the
SIS algorithm by appropriate choice of the importance density and/or modiﬁcation of

the resampling step.

4.4.1 Sampling Importance resampling (SIR) filter

In the SIR ﬁlter [68], the transitional prior is used as the importance density and
resampling is performed at each step. The advantage of the SIR ﬁlter is that importance

weights can be easily evaluated and the importance density can be easily sampled.

4.4.2 Auxiliary SIR (ASIR) Filter

The resampling step is performed at k - 1(using the available measurement at
positionk ), before the samples (particles) are propagated to position k [68]. The

advantage of the ASIR ﬁlter is that it generates samples from the state atk -1,

31

conditioned on the measurement atk . Therefore it is more likely to be in the region of

high likelihood. However, when state transition noise is large, the approximation of

p(xk ka_1) is poor.

4.4.3 Particle Filters with Improved sample diversity

In general, resampling is performed to address degeneracy. However, resampling
introduces other problems such as loss of diversity among theiparticles. Two techniques
to improve sample diversity include [68]:

(i) Regularized particle ﬁlter
(ii) MCMC Move step

These can be used when process noise is small.

4.4.4 Local LInearlzatIon Particle Filter
The optimal importance density can be approximated by using the most current
measurement z k via a bank of extended or unscented Kalman ﬁlters [68]. The idea is to

use for each particle a separate extended Kalman ﬁlter (EKF) or unscented Kalman
ﬁlter (UKF) to generate and propagate a Gaussian importance distribution. The

implementation of this variant of the particle ﬁlter is computationally very expenSive.
4.5 Particle filtering for NDE inverse problems

4.5.1 Choice of Importance density

The most commonly used variant of particle ﬁlter is Sampling Important re-

sampling (SIR) algorithm [67-72] due to its simple implementation. The SIR ﬁlter is

32

therefore used in this research. The importance density is chosen as the transitional prior

in the implementation of SIR algorithm:

405;; IXI;_1,Zk) = p(xIc ka-l) (4-11)
Therefore from (4.10) and (4.11),
wk °‘ WI-1P(Zk IXI) (4.12)
We can also write (4.11) as

(4.13)

W}, °c p(Zk IxI)

where p(z k ka) is the likelihood function. The likelihood function can be

represented by a noise PDF centered at the computed measurement for state x}: using a

measurement model. The noise PDF represents the difference between the actual

measurement and the computed measurement.

4.5.2 Implementatlon scheme for NDE inverse problems

In the NDE inverse problem, as discussed in chapters 2 and 3, the state is the
ﬂaw depth and the measurements corresponds to the NDE measurement. The posterior
PDF of the unknown parameter (i.e., ﬂaw depth) is evaluated at each position index
k =1: K in the problem domain using (3.2) iteratively. The sequence of steps for

evaluation of the posterior PDF of state (ﬂaw depth) is shown in ﬁgure 4-2.

33

83:58
5‘33. :32

9.8m

 

meal
“cw—Smog

 

858% use ocean .8 55303 .3 05mm

 

 

 

 

 

 

>29me

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 
 

 

 

 

 

 

 

 

 

 

 

 

 

Sun— and.

 

 

 

 

 

 

 

 

 

38 cocomeoEoO AI
AN _ NCQ 3%
E 02 I
mam—named _
8‘ new ,
~R. a. oz
. ON 2 II. ~
3 as t... s _
.. .. s Ems; wages: I
T . I .
H+~n~ TITS; «2T: SN _
s. _ ESE 30895982
«in: . x e

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

34

As discussed above, the posterior PDF is represented by samples and associated

weights. N 3 samples are sampled from the prior PDF at each position index. The prior

PDF as given in equation (3.2) relates the state (ﬂaw depth) at position kto state (ﬂaw

depth) in the neighborhood of position k deﬁned as x;- I 131er Therefore, the
j ¢k

sampling step can be represented in terms of the state transition model. The next step is
the assignment of weights to each sample. Weights are assigned using the likelihood
PDF (3.3) which relates state (flaw depth) to NDE measurement. The likelihood PDF is

represented by error in the measurement model. The posterior PDFS (samples and

associated weights) are estimated at all locations k = 1: K . Single point estimates of
the unknown state from the evaluated posterior PDF are also computed at each position
index (these estimates are discussed further in Chapter 6). The posterior PDFS at every
location are computed sequentially. These PDFS are then used as initial values and are
updated iteratively till the error between the single point ﬂaw proﬁle estimates at two

consecutive iterations is less than some preset convergence criterion. If each iteration is

denoted by I , and the estimate evaluated at I is estimate(I ) , then the convergence

criterion is given by

(estimate(l +1) — estimate(l))2 < T (4.13)

K _

 

Here, I is the preset convergence threshold. The posterior PDF is computed iteratively
till convergence (as deﬁned by equation (4.13)). At convergence, estimate(l +1) is

assumed to be the predicted ﬂaw proﬁle. The complete particle ﬁltering algorithm is

summarized in ﬁgure 4-3.

35

5:9:qu 2&on Bea H8 mains—m 20:58 me Bomto>0 .mé 053n—

ewaﬁ ESE

 

wouﬂEﬁ—mo
“2:5. 5.02 Ease
e
”Em

 

 

_ 838mm— wEEEH I—

 

 

 

 

eotoumom J
Egan seem? 5E 225m

:N— «5%

 

 

L

682 3082382

 

 

 

 

 

«2”; as

   

  

Etc—2
com—55¢. 38m

 

 

 

  

 

 

 

 

r
a «a
3% 60.5 a N
ooeosvom . s. o uooﬁonnwaz
“5805302 «2.7k— LR H «a a .

 

 

 

 

36

4.6 Models

As shown in ﬁgures 4-2 and 4-3, state transition and measurement models are used in

the implementation of particle ﬁltering to solve the NDE inverse problem.

4.6.1 Measurement Model

The measurement model typically relates the state to measurements and the
results of inversion may be expected to depend on the choice of measurement model.
Any type of measurement model can be used. In this study, both polynomial and neural
network based models are used due to their mathematically simple form and ease of
use. A comparative study of using different measurement models, in terms of accuracy
of inversion and computational load is also carried out [76]. The measurement model
parameters are derived from a training database of ﬂaw proﬁles and the corresponding

measurements .

4.6.1.1 Polynomial model

The measurements (response) are approximated from the state through an R-

order polynomial. The model can be expressed as follows:

R r (4.14)
Zk = Z Crxk
r=0

The coefﬁcients of the polynomial C are determined from a training database of known

states and corresponding measurements.

4.6.1.2 Neural Network based models

An artiﬁcial neural network is a mathematical model that mimics the structure

37

and function of biological neural networks. It consists of a group of artiﬁcial neurons
with weighted interconnections [77-78]. Each neuron transforms its inputs using a
transfer function and presents the result to all of its interconnected neurons. In most
cases, these neurons are connected in a layered fashion, with the outputs of neurons in a
layer applied as input to neurons in the subsequent layer. The interconnection weights
are updated using a learning phase based on training data. Their ability to model
complex functional mappings is utilized in this study, to create a functional mapping
between inputs (states) and outputs (NDE measurements). Two separate neural network
structures are investigated for use as measurement models - multilayer perceptrons
(MLP) and generalized regression neural networks (GRNN). MLP networks often use
the log-sigmoid, tan-sigmoid and/or linear transfer functions in each neuron, and
frequently use the error back propagation technique [78] for updating the weights given
a training data set. The use of radially symmetric Gaussian transfer functions results in

the generalized regression neural network (GRNN) [77]. The GRNN estimates the joint
probability density p(xk , z k ) of the input and output from a training database.

Specht [77] shows that the density function is estimated optimally when

S (4.15)
Zk = Z Wk,s¢s (Ika _CSII)

s=1

where Wk,s is the target weight corresponding to input training state x k and

output Z k , (Us is the radial basis function centered at C s- S is the total number of

radial basis functions used to train the network. Assigning the weights (i.e. training of
network) is not computationally expensive, however after training, these networks are

generally slower to use, requiring more computation to perform a function

38

approximation. Both types of neural networks have been shown to be universal

approximators, i.e., they are capable of approximating arbitrary functions [78].

4.6.2 State transition model

The state transition model relates the State xk (ﬂaw depth) at position index k

to the states x j I j=11Nk in the neighborhood. As given in [79], crack growth in most
j atk

rrricrostructures is a stochastic phenomenon. The state transition model is therefore
expressed in terms of a PDF with contributions from states in the neighborhood. This

prior PDF, as discussed in chapter 3, is given by:

 

2k Ila—yr ““6
p(Xk Ilej=1:Nk)=e 1:1

j¢k

39

CHAPTER 5

PARTICLE FILTER BASED DATA FUSION

5.1NDE Data Fusion

Frequently, NDE measurements from multiple measurement modes, conducted
on the same object or ﬂaw, are available. Often, multiple measurement modes contain
complementary information, and fusing data from multiple measurement modes can
potentially improve the inversion accuracy and reliability. In this chapter, the particle
ﬁlter based inversion algorithm is extended to make inversion using multimodal

measurements possible.
5.2 Data Fusion Using the Particle Filter

Assume that the measurements 21:1 Iq=le from Q measurement modes are

available at each position indexk. The objective is to evaluate the posterior PDF

p(kuz,1c,zi,...,sz) of ﬂaw depths at position k given the

measurements Z21 Iq=l:Q° As discussed in Chapter 4, the particle ﬁlter represents the

posterior PDF by means of samples and their associated weights. Therefore, as given in
equation (4.4), the posterior PDF p(xk I zk) at position index k is a weighted sum
of N s samples. The weight assignment to each sample is a function of the likelihood

PDF as given in equation (4.8). When multiple measurement modes are available, the

likelihood PDFS corresponding to each measurement mode need to be considered in the

40

weight assignment process. If we assume that qu is the weight of the sample i at

position index k assigned by measurement mode q , then for every sample at each

position index, Q weights can be assigned using the Q likelihood densities. Then the

likelihood function corresponding to the qth measurement mode is given by

'. ' (5.1)
Wi" °<= PM? Iii)-

If measurement processes are assumed to be independent to each other, then the joint

likelihood due to measurement modes q = 1,2, ..... , Q is the product of individual
likelihoods:
o 1 . 2 o o (5.2)

p(zk IxIC) = p(zk Ix;C ).p(zk Ix; ).,......,.p(sz IxIC).
Therefore, from (5.1) and (5.2), we have

u 1 o 2 o u (5.3)
w]: oc p(zk Ix;c ).p(zk IxIc) ........ p(sz IxIC).
Using (5.1) and (5.3),

1' i1 i,2 i.Q (5.4)

Therefore, the ﬁnal weight assigned to each sample is the product of the weights
assigned to the sample by the individual measurement modes. The overall algorithm for
particle based data fusion is shown in ﬁgure 5-1. The algorithm is executed at each
position index and the ﬁnal weight assigned to each sample is computed using equation

(5.4). Finally the posterior PDF of the state is computed at each position index.

41

Barrow—w .553 See Ewan LEE 225m =80>O .Tm BEE

Ewe: a5.
m N L

woumEm—mo

ease as .r sac mean:
is. 50:

      
   

.N. sac maﬁa;

 

.0. see mesa;

 

 

 

 

 

 

”HE L . O. 388 3082332
Hotcuwom
.N. BUOE aaOEOuu—wmoz

 

 

 

 

 

 

.7 SEE 80803302 TI

 

 

 

 

 

 

 

 

 

 

 

 

 

 

«2“; as
A «R _ SN 3 a ,
a; _ as: 8355 3st
a a a s
3% sub 3% eat 3% sub
.2. ooeosvom . I .N. boson—wow .7 ooeoscom a:
EoEoSmmoE 3088332 “3805302 «ZWH: NR | «R we eoosuonnwmoz

 

 

 

 

42

In using equation (5.4), it is assumed that the measurement processes are independent of
each other. However, this assumption may not be valid in many cases. To address this
issue, the principal component analysis technique is applied on the data from different

measurement modes.
5.3 Principal Component Analysis

5.3.1 Theory

Principal component analysis (PCA) is a powerful tool for analyzing data [80].
It involves a mathematical procedure which transforms possibly correlated variables
into uncorrelated variables called principal components. The central idea behind PCA is
to reduce the dimensionality of data, retaining only those uncorrelated dimensions of
data which can explain all the variations present in the data. However, in this study the
dimensionality of data is not changed; correlated data is transformed into uncorrelated
data only. The dimensionality of the data is retained in order to compare the inversion
results of particle ﬁlter based data fusion with and without principal component analysis
of measurement data.

PCA is mathematically deﬁned as an orthogonal linear transformation that
transforms the data to a new coordinate system. Figure 5-2 represents this
transformation of coordinate system for a sample two dimensional data set. In the new
co-ordinate system, the greatest variance by any projection of the data lies on the ﬁrst
coordinate (called the ﬁrst principal component), the second greatest variance on the
second coordinate, and so on. Given a set of points in Euclidean space, the ﬁrst
principal component (the eigenvector with the largest eigenvalue) corresponds to a line

that passes through the mean and minimizes sum squared error with those points. The

43

second principal component corresponds to the same concept after all correlation with
the ﬁrst principal component has been subtracted out from the points. Each eigenvalue
indicates the portion of the variance that is correlated with each eigenvector. Thus, the
sum of all the eigenvalues is equal to the sum squared distance of the points with their
mean divided by the number of dimensions. PCA essentially rotates the set of points

around their mean in order to align with the principal components.

«tum mama—nae Loam 852. San SV 83 353.5 A3 8% Bum—=86 :o mum—ace <8 he sown—5358 3.585 Nd Raw?—

e 3

 

 

 

 

 

 

 

 

 

 

 

 

45

5.3.2 Implementation

PCA can also be applied on data from multiple measurement modes to evaluate

the principal components. At each position indexk , the measurements ZIEQ I 1:1sz

in the neighborhood of k from all measurement modes (1: Q) are considered. The

resulting Q-dimensional data vector is the input to the PCA technique. The following

steps are carried out to evaluate the principal components [80-81]:

Step 1: Data from dimensions zisQ I j=12 Nk is stored as a Q-dimensional vector

§ q |q=1:Q -
Step 2: Each measurement vector is adjusted by subtracting out the mean of the data.

‘_' "_° —' 5.5
gq=gq—mean(gq). ( )

Step 3: The adjusted data vectors are arranged as rows of a matrix. This newly formed

matrix is termed as g.

= —'—’ —" (5.6)
g=[g1,g2, ...... ,gQ].

Step 4: The covariance matrix of E is computed.

= = (5.7)
0‘1 = cov(gq) .
Step 5: The eigenvectors /Iq of E are then evaluated.

(5.8)

2‘1 = eig(0'q).

Step 6: The matrix /I of computed eigenvectors is formed.

46

(5.9)

_.

Tia/TIE, ...... ,AQ].

Step 7: The product of /I and the 5‘ matrix results in the principal components:

= 2;. (5.10)

ﬁll

The rows of VI are the principal (uncorrelated) components in the data.
(5.11)

—.

The resultant Q-dimensional components W q are independent. These
components are then treated as data from separate measurement modes, and applied to

the particle ﬁlter based data fusion technique as described in section 5 .2.

47

CHAPTER 6

RESULTS

The proposed particle ﬁlter based technique was applied to simulated and
experimental NDE data. Evaluation metrics were used to evaluate the performance of
proposed inversion algorithm. Details of these metrics are given in section 6.1. Three
different databases were used in the evaluation of the proposed technique. Simulated
and experimental eddy current NDE data from steam generator tubing were used to
initially evaluate the algorithm and choose an appropriate measurement model. A
description of this data is given in section 6.2. A comparative study of the performance
with the measurement models discussed in Chapter 4 is presented in section 6.3. This is
followed by the inversion results on simulated data (with and without additive noise),
and experimental data. Section 6.6 presents the results with the second database (eddy
current measurements from aircraft lap joints) while Section 6.7 presents results

obtained for magnetic ﬂux leakage (MFL) NDE.

6.1 Evaluation Metrics

As discussed in the previous chapters, the posterior PDF of the unknown
quantity is evaluated at discrete positions using the proposed technique. Two types of
single-point estimates can be computed from the posterior PDF at each location —
Posterior mean and maximum a posteriori (MAP) estimates — which are assumed to be

the predicted ﬂaw proﬁles in this study. The posterior mean estimate [82] is given by

48

AMMSE (6.1)
We = Elxk lzkl= ka-P(xk law/.-

while the maximum a posteriori (MAP) estimate [82] is given as

AMAP (6.2)

kuk = arg max p(xk I zk).
xk

Figure 6-1 shows an example of the posterior PDF of ﬂaw depths evaluated throughout
the length of the NDE specimen. Figures 6-1 (a) and (b) show the top view and the 3-D
view, respectively, of the posterior PDF. Figure 6-1 (c) shows the posterior PDF at one
position index k. The MAP and MMSE estimates are indicated. In ﬁgure 6-1 ((1) the
predicted ﬂaw proﬁles, using both the MAP and MSE, are plotted with the true ﬂaw
depth at each position index. The resulting predicted proﬁles using the posterior mean
estimates and the MAP estimates as shown in ﬁgure 6.1 (d) were compared with the

true proﬁles using a mean-square error (MSE) metric [83]:

 

K (6.3)
2 (xk( predicted) T xk(actual))
MSE = "=1
K

49

88:88 2&on 3am AB SE 59: acumen one “a ma...— uotoumom 8v mam 85qu mo 33> Dim A8 can 303 nob AS To oEwE

6C
593 3?.
E

       

 

I .

 

_
ease s _ Amaze seam
n 582 Breech
_
_
H _ ex SE
. _
m _ A ~_ 3%
ﬁe
woman oEEmw c3
59.3 .32 . awed a5
0
M
8m 5
9
00¢ m.
w com W
w
, 89 cor

ecc—

 

50

6.2 Eddy Current NDE Data Description

6.2.1 Simulation Data Description

A 2-D ﬁnite element model (FEM) [84] is used to simulate the eddy current
inspection of steam generating tubing (SGT). The model generates the measurement
signals (eddy current data) from known crack proﬁles at three different frequencies (200
kHz, 400 kHz and 800 kHz). The eddy current data consists of magnitude, phase, vertical
(imaginary) component and horizontal (real) component of the complex eddy current
signal at different frequencies. An example of the measurements for a simulated ﬂaw is
shown in Figure 6-2. The tube conductivity and relative permeability for this example
were 106 S/m and 1, respectively. The tube thickness and axial length was 2.5 mm and

3.75 mm, respectively.

51

Ea Baa—58mm a e8 BeoEoSmeoE 39:8 batm— .N.© oSmE

 

 

 

 

 

 

 

 

 

 

 

2%.; 3E became: a> as.
v m n _ o .8”—
OS
% -
m .
m.
M r
m
M” l
3
W1 1
e
T
eves; 095:3:
4 m N _ o a. m N
. . m2 .
- . o2 -
m
m. T
- - m
- . m2 .
. 4 r . . . . an; t

 

 

 

 

 

 

 

Amurﬁeurl

52

The parameters of the measurement models used in this study (polynomial or
neural network based models) were derived from a training database consisting of the
ﬂaw depth proﬁle and the corresponding measurement. This database was constructed
using the FEM model described above. To evaluate the performance of particle ﬁlter
inversion algorithm, ﬂaw proﬁles were predicted from measurement sequences on test
data that were not present in the training database. Artiﬁcial noise at several signal-to-
noise ratio (SNR) levels was also added to the measurements to assess the performance
of the algorithm both in the presence and absence of noise. In this study, the SNR was

deﬁned as the ratio of measured signal power and added white Gaussian noise power:

 

P 2 (6.4)
signal 0'
SNR = =
P noise N
rms

where Psignal is the power of measurement signal, Pnoise is the power of added

noise, and 0' and N rm s are the RMS amplitudes of measured signal and added noise

respectively.

53

6.2.2 Experimental Data Description

The experimental setup for eddy current testing is shown in ﬁgure 6-3 [85].

I Nircumferential

Axial Phase agnitude

 

 

Impedance Plane
Inconel tube

ﬂ HF Pancake coil
Pancake coil + Point coil

 

MRPC Probe

EC image

Figure 6-3. Experimental ECT

The experimental data was collected using a motorized rotating probe coil
(MRPC). The MRPC probe contains three coils: plus (‘+’) point coil, pancake coil and
high frequency (HF) pancake coil. The data was acquired through axial rotation and
simultaneous longitudinal motion of probe, resulting in an image as shown in ﬁgure 6-
3. A rectangular region of interest was selected on each image and the maximum
magnitude and corresponding phase at each position along the tube axis within the
region were computed. The data was collected from both laboratory-induced cracks and

cracks in ﬁeld tubes. Metallographic results for all ﬂaws were used as “ground truth”

54

for the ﬂaw proﬁles. For this study, plus (‘+’) point probe magnitude and phase at 100
kHz, 200 kHz and 300 kHz were used as measurements. The magnitude and phase of
eddy current signal at 300 kHz for an experimental ﬂaw of 7 mm width and maximum

thickness 100% of wall thickness are shown in ﬁgure 6-4.

Training data (for estimating the parameters of the measurement models)
consisted of a subset of ﬂaws from this database, and the algorithm performance was
veriﬁed usingia separate test database of ﬂaws. As with the simulated data, results were
evaluated by a direct comparison between the actual and predicted proﬁles (Posterior

mean and maximum a posterior estimates from posterior PDFS).

55

Volts

 

 

 

Radians

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5 0 5
Magnitude Phase
a
o
5
.o
S
'8
3
.8
B
is"
100 1 L . . . . e 1 .
0 5
Flaw Proﬁle

Figure 6-4. Eddy current measurements at 300 kHz for an experimental ﬂaw

56

 

 
  

Measurement

 

 

  

 

AAAAAAAAAA

 

-------—-_---
l- I

 

 

I
l'
p.
I:
p
P
l-
I-
p
II

 

Proﬁle

Figure 6-5. Illustration of the data registration problem

Data registration was an important issue in working with the experimental data
(ﬁgure 6-5). Metallographic information was available at a higher spatial resolution,
and the measurement data was interpolated using linear interpolation to address the

issue of registration.

6.3 Comparison of Measurement Models

Measurements from six different simulated ﬂaw proﬁles were used to test the

57

algorithm performance using the three different measurement models discussed in
Chapter 4 (polynomial, MLP neural network and generalized regression neural
network). The MLP neural network was a feed forward network with three layers and
sigmoid transfer functions. The back propagation algorithm was used to train the neural
network. The details regarding the ﬂaw proﬁles used in this study are tabulated in Table
6-1. The resulting predicted proﬁles were compared with the true proﬁles using mean-
square error (MSE) metric (under both noisy and non-noisy measurement conditions).
The associated computational cost for each measurement model was also evaluated. In
all cases, the results were evaluated for the state vector length parameter L=l. The

algorithm was tested on an 8-core dual Intel Xeon E5345 2.33GHz 64 bit processor,

with 8 GB RAM. The convergence criterion used for evaluation was 2' = 10_3.

Table 6-1 Flaw proﬁles used in comparison of measurement models

 

 

 

 

 

 

 

 

 

 

S.No Width @m) Depth (%)
1 1.27, 0.635 20,40
2 1.9,0.381 60,40
3 0.635, 1.27,0.3175 20,50,30
4 0.381,0.635,0.635 15,40,25
5 0.635,1.27 30,50
6 0.317,1.27 75,35

 

58

 

meme»; 39:84“ 3v $38.8 mo “page SEQ, @9808 as “woo Reouﬁaﬁou A3 ﬂeece “582388 mo commaaﬁoo 6-0 oSwE

coo H

moaﬁem

3v

com

 

 

 

 

“52 I
boa ImI

 

 

..

 

 

(asw) 301

338.8

 

 

 

3
338mm
.. com I
. ‘\.
ZZMO I
n32 I
boa I¢I

 

 

 

ow

ow

 

 

spuooos

59

Figure 6-6(a) displays the computational cost involved in ﬂaw proﬁling as a
function of the number of samples used at each spatial position in the inversion
technique for the different measurement models. Figure 6-6(b) displays the accuracy
in terms of the logarithm of the average MSE between the predicted and true proﬁle
for all six ﬂaws as a function of number of samples. Figure 6-7 presents the MSE
between the predicted and true proﬁle as a function of SNR for all six ﬂaws. The
different models are indicated by different colors in the ﬁgure. A summary of this
information (in the form of mean and standard deviation of the MSE) using different

measurement models are shown in ﬁgure 6-8.

 

a Polvnomial

 

 

 

 

Flaw ID "'1 SNR

Figure 6-7. MSE between PME and actual proﬁle, versus SNR for different
measurement models

60

0.07 r 1 I I

 

 

Polynomial
—- MLP
"- - GRNN

 

0E5 -

 

 

 

___ . . I_____._.‘
l

I
I
I

T

0.05

0.04 ~

 

0.03

I

MSE
I-——- _.
I——»\\-—'»..-.

0.02 - _

0.01 L _

 

 

 

Int 10 5 2
SNR

Figure 6-8. Mean and standard deviations of MSE (between PME and true proﬁle)
versus SNR using different measurement models

The results presented above indicate that the polynomial measurement
model offers better accuracy when compared to the neural network models on this
limited set of ﬂaws. The computational cost using the polynomial based model is
also seen to be the lowest. The accuracy of polynomial based measurement model is
also seen to be better in the presence of noise. It is possible that the use of other
neural network architectures (or different training databases) may improve the
results. However, the selection of optimal neural network architecture 5 is a difﬁcult
task, and will increase the complexity associated with solving the inverse problem.
For this reason, the polynomial model was selected as the measurement model for

subsequent studies of the proposed inversion technique.

61

The choice of number of samplest , length of state vector (governed by

parameter L), and type of measurement model used affect the computational cost
and accuracy of the inversion results. The parameters for ﬂaw proﬁling using the
simulated eddy current tubing inspection data were selected using the parametric
studies described in this section. Similar studies were conducted for the other

databases, and used in the selection of these parameters.
6.4 Results on Simulated Data

6.4.1 Simulated Data without Noise

The magnitude of the complex eddy current signal is used as measurement

data to further evaluate the performance of the proposed algorithm. The number of

samples used at each position index is 1000, while a 2nd -order polynomial model
is used as the measurement model. The effect of changing the state vector length L

was also investigated in this study, with L set to 0, 1, and 2. The convergence

. . . -3
crrterron as drscussed above was 2' = 10 .

62

83:53 339:8 AB "HE
Sigma mo 32> Q-m A8 ”HE Horsemen mo 32> mom. A8 See E250 have noun—:86 mEms 05$an .5652: @8395 db oEmE

8v
ewe: E3

oiod m<2
85.0 m2."—

MwE Egon S

 

BC
0223 oEEmm 3

59.3 .32 Ewe: 35..

Ausuap Jouarsod

cow

coeds ojdums

 

coo—

 

63

An example of ﬂaw depth proﬁling from eddy current data for a multi-level
rectangular ﬂaw (lengths: 0.9, 0.7, 0.7 mm, depths 20%, 40%, 60% of tube wall
thickness) in the absence of added measurement noise is shown in ﬁgure 6-9. Figure 6;
9(b) shows the (3-D view of) estimated PDF of ﬂaw depths as a function of position
while ﬁgure 6-9(a) shows the same information as a top view (2-D view). Note that the '
PDF at any location may be multimodal (as shown in ﬁgures 6-9 (a) and (b)). The state
vector length parameter L was set to 1 and measurement polynomial order was set to 2
while evaluating the PDF estimates. These PDF estimates are the output of the particle
ﬁlter, and are used to compute posterior mean and MAP estimates of the ﬂaw depth at
each location. Figure 6-9(c) compares the estimated ﬂaw proﬁles (using the posterior
mean and MAP estimates) with the true proﬁle. MSEs between PME/MAP estimates
and actual proﬁle are also given.

Figure 6-10 shows the inversion results using magnitude measurements at
multiple frequencies for a ﬂaw of width [0.75mm,0.635mm,0.7mm] and
depth=[40%,80%,60%] . Figure 6-10(a) shows the magnitude measurements, ﬁgure 6-
10(b) & 6-10(c) show the 3-D and top views of posterior PDFS computed using the

proposed data fusion technique without PCA.

838:8 AB “am 85309 mo 32> 98 A3 32> Q m 39 3655on 03258 “a 8% mDZ A8 038508 c2802: .26 oSmE

:8 6»

:53 :62 €33 E3

 

 

 

 

 

 

 

 

  

[usuap JOMOJSOd

 

 

 

 

65

A study was also carried out to determine the effect of state vector length
parameter L on the computational cost. The computational cost for evaluation of
posterior PDF at each location (spatial position index) is plotted versus different values
of L in ﬁgure 6-11. The number of samples at each position index is kept constant (1000
samples for all values of L). As expected, the cost associated with the inversion

procedure increases with L. However, the increase is seen to be linear.

 

Seconds

 

 

 

 

Figure 6-11. Computational cost versus state vector length parameter L (size of
contributing neighborhood). The data presented is the average computational cost over
20 ﬂaws, without data fusion.

Twenty different multi-level rectangular defect proﬁles were selected from the
test database for further investigating the performance of the proposed inversion
technique. The ﬂaw proﬁles are tabulated in Table 6-2. As indicated in the Table, the

width of the ﬂaws increases from 0.75 mm to 5.50 mm. Figure 6-12 presents mean

square errors between predicted and true ﬂaw proﬁles for different measurement modes

66

as a function of state vector length parameterL = 0,1,2. The MSE between predicted and

true ﬂaw proﬁle using measurement data at three individual frequencies (200,400 and
800 kHz) are shown in blue, red and cyan colors in the ﬁgure. The MSE for data fusion
(i.e., combining the information in all the frequencies) without PCA and with PCA are
shown in green and magenta colors, respectively, in the ﬁgure. A summary of the
inversion results for these twenty ﬂaws is presented in ﬁgure 6-13, which shows the
mean and standard deviation of the MSE (averaged over all 20 ﬂaws). The mean error

(over all 20 ﬂaws) is also tabulated in Table 6-3.

Table 6-2 Simulated ﬂaw proﬁles used in particle ﬁlter evaluation

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Flaw Width Depth

ID (mm) (% tube wall thickness)
1 0.75 20
2 1.00 20
3 1.25 20
4 1.50 20
5 1.75 20
6 2.00 40
7 2.25 40
8 2.50 40
9 2.75 40
10 3.00 40
11 3.25 60
12 3.50 60
13 3.75 60
14 4.00 60
15 4.25 60
16 4.50 80
17 4.75 80
18 5.00 80
19 5.25 80
20 5.50 80

 

67

 

 

 

 

.........
...........

 

 

 

   

012
01201201
2012 20 FlawlD

L

Figure 6-12. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' using different measurement modes

68

0.018 I I T l 17

 

 

  

 

 

 

 

 

 

 

. ——200kHz
0015 ~ 1 —-—400kHz ~
i ._ ---~-— 800kHz
0.014 - N ._._ DF
1 ‘ - 1 ——DF PCA
0.012 - 1 -
0.01 1 _
UJ
‘9
0.003 - .
0.005 - -
0.004 . Lu. -
\‘K‘
0.002 ~ ---1 ————1 -
D A l l I l
-05 0 0.5 1 1.5 2 2.5

Figure 6-13. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L using different measurement modes

Table 6-3. Inversion results for simulated non-noisy data

 

 

 

 

 

 

 

 

 

 

State vector length parameter ‘L’
Measurement 0 l 2
mode PME MAP PME MAP PME MAP

ZOOld-Iz 0.0140 0.0] 71 0.0 l 29 0.0340 0.0093 0.029

400kHz 0.0921 0.1017 0.0070 0.08l0 0.0058 0.0063

800kHZ 0.0 l 51 0.0216 0.0129 0.0093 0.0089 0.0114
Data fusion 0.0036 0.0028 0.0021 0.0026 0.0018 0.0034
Data fusion 0.0028 0.0039 0.0020 0.0045 0.0020 0.0049
(with PCA)

 

 

 

 

 

 

6.4.2 Simulated Data with Noise

 

To verify the robustness of the algorithm, the inversion was also carried out with

additive white noise added to the complex measurements at various SNR levels. The

69

magnitude of the complex eddy current signal is used as measurement data to evaluate
the inversion results. The inversion results are evaluated for two different SNR levels
(SNR=5 and 2). The ﬂaw proﬁles used for computing results for noisy data are the same
as those listed in Table 6-2. Figures 6-14 and 6-16 present the inversion results for
different measurement modes, with SNR equal to 5 and 2 respectively. A summary of
the inversion results is shown in ﬁgures 6-15 and 6-17 respectively. Mean errors (over

all twenty ﬂaws) are tabulated in Tables 6-4 and 6-5 respectively.

 

— 200 kHz
— 400 kHz
800 kHz

 

 

 

 

 

 

 

 

 

 

10

 

012012012012012 20
FIawID

L

Figure 6-14. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' using different measurement modes with additive noise (SNR=5)

70

0.04

0.035

0.03

0.025

0.02

MSE

0.015

0.01

0.005

 

I

I

I

I

I

I

 

 

 

 

 

———2El]kHz

-—4EDI<HZ .

 

 

 

 

 

0.5

I—dr-I-H

1.5

M-rH

2.5

Figure 6-15. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L using different measurement modes with additive noise (SNR=5)

Table 6-4. Inversion results for measurement data at SNR=5

 

State vector length parameter ‘L’

 

 

 

 

 

 

 

 

Measurement mode 1 2
PME MAP PME MAP PME MAP
ZOOkI-Iz 0.0302 0.0512 0.0160 0.0321 0.0108 1 0.0164
400ml 0.0251 0.0312 0.0133 0.0165 0.0090 0.0150
8001(1-12 0.0343 0.0452 0.0 l 78 0.0329 0.011 1 0.0179
Data fusion 0.0032 0.0051 0.0019 0.0017 0.0017 0.0039
Data fusion 0.0030 0.0078 0.0024 0.0060 0.0018 0.0032
(with PCA)

 

 

 

 

 

 

 

71

 

 

 

 

 

 

 

 

 

 

 

 

 

 

    

10

01201201201201; 20

Flaw ID
L

Figure 6-16. MSE between PME of posterior PDF and true proﬁles versus state vector
length parameter 'L' using different measurement modes with additive noise (SNR=2)

72

0.12

 

 

 

 

 

 

 

 

 

 

 

2CDkHz
-—4mkHz
0.1 1 1l-.._..._ BEDkHz
1 —— DF
1 —— DF PCA
0.08 -
0 005 ~
2
0.04 -
0.02 ~
H : ———-——-__—-:§
0 1 1 I I l
-05 0 0.5 1 1.5 2
L

2.5

Figure 6-17. Mean and standard deviations of MSE (between PME and true proﬁle)
versus L using different measurement modes with additive noise (SNR=2)

Table 6-5. Inversion results for measurement data at SNR =2

 

State vector length parameter ‘L’

 

 

 

 

 

 

 

 

 

Measurement mode 1 2

PME MAP PME MAP PME MAP

ZOOkHz 0.0705 0.0736 0.0385 0.0302 0.0252 0.0263

400kHz 0.0604 0.0580 0.0331 0.0346 0.0216 0.0179

800kI-Iz 0.0930 0.1202 0.0479 0.0394 0.0302 1 0.0301

Data ﬁision 0.0116 0.0186 0.0076 0.0140 0.0072 0.0189

Data fusion 0.0078 0.0150 0.0050 0.0150 0.0038 0.0123
(with PCA)

 

 

 

 

 

 

 

6.5 Results with Experlmental Data

Measurement data used for inversion of experimental data (discussed in section
6.2) is the magnitude and phase of the complex eddy current signal at 100 kHz, 200 kHz
and 300 kHz. The number of samples used at each position index is 2000, while a 3rd

order polynomial model is used as the measurement model. The state vector length L

was varied between 0, 1 and 2. The convergence criterionT =10_3. Figure 6-18 shows
the inversion results using magnitude only measurements at 300 kHz for a ﬂaw of width
7 mm and maximum depth of 100% of tube thickness. Figures 6-18 (a) and (b) show
top and 3-D views of the computed posterior PDF respectively. Figure 6-18 (c) shows
the estimated ﬂaw proﬁle. The results for test data, i.e., magnitude & phase
measurements (shown in Figure 6-19 (a) and (b)) are shown in ﬁgure 6-19 (c).
Inversion results using magnitude measurement at each of three frequencies (100, 200
and 300 kHz) are shown in Figure 6—20 (a), (c) and (e). Results using magnitude and
phase measurements at each frequency are shown in Figure 6-20 (b), (d) and (f). Finally,
the inversion results using magnitude and phase measurements at three frequencies are
shown in the ﬁgure 6-20 (g). Table 6-6 tabulates the MSE between true and predicted

proﬁles for this ﬂaw.

74

 

 

 

 

 

 

 

- 5 I t : It;

 

 

............

 

 

 

 

 

 

 

 

 

835qu 33988 A8 gags—2 3V omega A3 059508 gangs .36 Rama...“

A8

 

 

 

 

 

w 3 o
o - . . o
. f L m.
-om - ._m
2 H .
s . . a

 

 

 

 

 

76

3h omen so .358 memo—Hogan 025 98 xenon—wot 2mg 3 3580532: mam: 338.— :2835 .86 nun—mi

NE oom©omaam 0% ocaEwn—z A3

 

NE oom©ovBEuaE A8

 

 

 

»

. 11.1.13 IT???

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

  

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

77

88:03on 8:: can maggot 2w?“ 3 33805308 m5? 8.3.8 comp—2,5 dab onE
N5 838q§ @665 a 83282 3

 

EM 03 @823 a. wean—Ma: S 2.3 OS ©ovBEm~2 A3

@422 IT
MSE .. .1
out.

 

 

78

Table 6-6. MSE between actual proﬁle and PME/MAP estimates using different
measurement modes as shown in ﬁgure 6-20

 

 

 

 

 

 

 

 

 

Frequency Magnitude Magnitude & Phase
(KHz) MAP PME MAP PME

300 0.0141 0.0566 0.0112 0.0340

200 0.0412 0.0562 0.0391 0.0411

100 0.0457 0.0639 0.0242 0.0481

100,200 & 300 - - 0.0105 0.0240

 

 

 

The proposed technique was applied to twelve additional ﬂaws selected from an

industry test database. Summary information about the ﬂaw proﬁles (maximum width

and depth) are tabulated in Table 6-7. Figure 6-21 presents some examples of the

proﬁles in this test database.

Figure 6-22 presents mean square errors between MAP and true ﬂaw proﬁles for

these twelve ﬂaws, when using one or more inspection frequencies, as a function of L.

Table 6-7 Experimental ﬂaw proﬁles

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Flaw Width Max Depth
ID (mm) (% tube wall thickness
1 2.75 20
2 7.00 100
3 5.50 100
4 3.00 85
5 2.75 90
6 2.25 70
7 5.25 78
8 l 1.00 60
9 2.50 58
10 3.00 71
l 1 3.25 39
12 3.50 76

 

 

79

 

 

 

 

% Depth

 

100*

 

 

 

 

 

 

(a) Flaw ID 5 (b) Flaw ID 8

Figure 6-21. Two experimental ﬂaws
MSEs between MAP and true ﬂaw proﬁle using measurement data (Magnitude and
phase) at three individual frequencies (100,200 and 300 kHz) are shown in blue, red and
cyan colors. MSEs between PME and true proﬁles using data fusion of measurement
modes without PCA and with PCA are shown in green and magenta colors in Figure 6-
22. A summary of the inversion results for the twelve ﬂaws tabulated in Table 6-7 is
given in Figure 6-23. Average and standard deviation of MSEs for each measurement
mode are shown in Figure 6-23 in the respective colors. Averages of MSEs between

PME/MAP estimates and true ﬂaw proﬁles are tabulated in Table 6-8.

80

 

 

 

 

 

Figure 6-22. MSE between MAP of posterior PDF and true proﬁles versus state vector
length parameter 'L' using different measurement modes for experimental data

81

 

 

 

 

 

 

 

 

 

 

 

 

 

0.05 . . . 1
— 10011112
0055 - ——200kHz
~» -------3001<Hz 1
0.05 - ——-—— 01:
—-—DF PCA .
0.045 - \
0.04 - [ ~
”01 0.035 - 1. -
2 ..
0.03 - ~1 1 .
1r _ 1 ' ’ t :-
0025 - 1" 1 t -
0.02 ~ ~-L~X I —
0.015 1 3‘ 1 ﬂ -
001 I l J i J
-05 0 0 5 1 1 5 2 2 5

Figure 6-23. Mean and standard deviations of MSE (between MAP and true proﬁle)
versus L using different measurement modes

Table 6-8 Inversion results for experimental data

 

 

 

 

 

 

 

 

 

 

State vector length parameter ‘L’
Measurement mode 0 1 2

PME MAP PME MAP PME MAP

1 OOkHz 0.0464 0.0360 0.0312 0.0313 0.0332 0.0305

ZOOkHz 0.0676 0.0512 0.0646 0.0447 0.0579 0.0413

300kHZ 0.0418 0.0291 0.0694 0.0237 0.0301 0.0215

Data fusion 0.0311 0.0216 0.0790 0.0176 0.0389 0.0174

Data fusion 0.0441 0.0207 0.0590 0.0171 0.0463 0.0156
(with PCA)

 

 

 

 

 

 

82

 

Results on some of the experimental data were compared to those from a
previous study that used multilayer perceptron (MLP) neural network based ﬂaw
proﬁling (described in [86]). The comparative results are given in Table 6—9. The
particle ﬁlter (PF) based inversion results were obtained using fusion of data at 100, 200
kHz and 300 kHz, with L equal to 2. The higher accuracy result is underlined for each

ﬂaw in the table.

Table 6-9 Comparison of neural network and particle ﬁlter based inversion techniques

 

 

 

 

 

 

 

 

Flaw ID Flaw Proﬁle MSE
Width Max. Depth NN PF
(min) (%)
3 5.75 100 0.0314 0.0461
6 2.50 80 0.0131 0.008;
7 5.25 78 0.0291 0.0085
9 3.5 70 0.0218 0.0169
1 1 3.75 50 0.0387 0.015;
12 4.50 85 0.0191 0.0230

 

 

 

 

 

 

 

6.6 2-D Flaw Profiling

6.6.1 Data Description

Quantiﬁcation of wall thickness loss due to corrosion in aircraft lap joints was
selected as a test problem to investigate the applicability of the proposed inversion
technique to 2-D ﬂaw proﬁle reconstruction. The NDE measurement data was acquired

from a 30-year old service-retired Boeing 727 aircraft.

83

 

Figure 6-24. Aircraft lap joint specimen

Two specimens (C and D), each consisting of a two-layer aluminum 2024-T3
lap joint cut out from below the cargo ﬂoor (ﬁgure 6—24), were inspected using the
multifrequency eddy current method. The thickness of each layer of the lap joint was
0.045 inches. The specimens were disassembled and cleaned of all corrosion products
after the inspection. X-rays mappings were then acquired to assess the true plate
thickness. The inversion results were compared to the true thickness to determine the
accuracy of the proposed inversion technique. Measurements from multiple frequency
eddy current inspections (5.5 kHz, 17 kHz and 30 kHz) were available, and used in the

inversion. The number of samples, N s in the particle ﬁlter was 2000. The resulting

predicted proﬁles using posterior mean estimates and the MAP estimates are compared

with the true proﬁles using a root mean-square error (RMSE) metric [87].

(6.4)
lg; lg; (xpredicted _ x )2 }6
m,n m,n
RMSE= m=1n=1

 

MN

84

6.6.2 Results

True and predicted plate thickness (from single measurement modes) at all
positions for specimen C is plotted in ﬁgure 6—25. The inversion result evaluated by
fusing multiple measurement modes for section C and section D are shown in ﬁgure 6-
26 and ﬁgure 6-27 respectively. A summary of the results for both specimens are
tabulated in Tables 6—10 and 6-11 respectively, and show the potential for improved

accuracy when fusing multisensor data. In Table 6-10 the inversion results on section C

are given. RMSE when using data fusion mode was found to be 1.86X10—3. In

contrast, the RMSE for data fusion using a wavelet neural network technique (WNN),
applied on the same data, was 3.017X10_3, while that when using a generalized

additive method (GAM) was 0.89X10—3 [86]. In Table 6-11 the inversion results on

section D are given. Again, RMSE when using data fusion was 2.21 X103. In contrast,
the RMSE for WNN was 3.53X10_3 , while that for a generalized additive method

GAM was 1.40X10‘3 [37].

85

 
  

Amoco:— EoEoSmon 2mva 0 cowoowl 20:8me Ema. 9: Ho.“ 3:52 .5652: .26 PEER

~04 0mm A3

Na :Hm— 3v

 

Good

Bod

086

956

 

NEV— om (rm 3

 

03H 3

 

 

 

86

3065 Sue Q concom- 538%.: Ema. as com 3:52 86.435 nab oBmE
5me San— A8

03H. A3

omod

 

‘ . 'l-u' I“, 53”..”

it”?

9.06

 

Evan 85 a:

 

cap. A8

299

 

Omod

 

 

mvod

 

 

 

87

Table 6-10. Inversion results for lap joint inspection (Section C)

 

 

 

 

 

 

 

 

S.No Measurement mode RMSE (X103)
1. 17 kHz Eddy Current 3.19
2. 30Khz Eddy Current 2.96
3. P-ET LOI 3.89
4 Data Fusion (l7kHzET+30kHzET+P-ET LOI) 1.86

 

 

Table 6-11. Inversion results for lap joint inspection (Section D)

 

 

 

 

 

 

 

 

 

S.No Measurement mode RMSE (X108)
1 17 kHz Eddy Current 4.31
2 30Khz Eddy Current 3.82
3. P-ET LOI 5.13
4 Data Fusion (l7kHZET+30kHZET+P—ET L01) 2.21

 

6.7 Magnetic flux leakage problem

6.7.1 Data Descriptlon

A ﬁnite element model for MFL was used to simulate the tangential component
of the leakage ﬂux due to ﬂaws in steel bars, and a range of ﬂaws were simulated to
generate large training and test databases. White Gaussian noise at different SNR levels
is also added to verify the robustness. The state transition model used here was similar
to the model used for eddy current data inversion (Section 6.2). Polynomial based
model of order 2 is used as measurement model. 1000 samples were used at each

position index to evaluate the posterior PDF of ﬂaw depth at that spatial position.

6.7.2 Results

Flaw proﬁling results from MFL measurements are shown in Figure 6-28. The

ﬁgure shows the proﬁling results for an irregular ﬂaw proﬁle with non-noisy test data.

88

The value of L is set to 1 and measurement model order N is set to 2. A summary of the
MSE between the predicted (PME) and actual proﬁles for ﬁve different ﬂaws without
and with different levels of measurement noise added, are presented in Figure 6-29 and
Table 6-12. Five different ﬂaw proﬁles are tabulated in Table 6-12. The results are
presented for values of L from 0 to 2 at different values of SNR in Figure 6-26 and
Table 6-13. Convergence criterion parameter 1 value was set to 10'3 for evaluation of
these results.

Table 6-12. Flaws used for MFL signal inversion

 

 

 

 

 

 

 

Flaw ID Width Depth
(mm) (% of tube wall thickness)
1 15 65
2 23.5 40
3 35 25
4 47.5 70
5 54.5 40

 

 

 

 

 

— No noise
_ SNR=5

 

 

 

 

2 Flaw 10

Figure 6-28. MAP estimates evaluated with MFL test data at different SNR levels

89

. \\ \ﬁ
\\\ “EN?
\\\1 \\§ \\‘\’~:\1 \

 

 

 

 

 

 

 

 

 

 

 

i

' E
"M. V
e
01
l:
, .3
if) E

Indoo %

(a)

Axial Length (mm)

 

coeds eldureg

90

(C)

Figure 6-29. Posterior PDF estimates (L=2) for ﬂaw proﬁle from MFL test data (a) 2-D view of PDFS (b) 3-D view of PDFS (c)

Estimates evaluated from PDFs

Table 6-13. MSE of the predicted proﬁles using MAP estimates as a function of SNR

 

 

 

 

 

 

 

 

and L for MFL test data
Flaw SNR=oo SNR=5 SNR=2
ID L=0 L=l L=2 L=0 L=l L=2 L=0 L=1 L=2

1 0.0390 0.0276 0.0246 0.0708 0.0501 0.0580 0.014 0.0796 0.046
2 0.0308 0.0298 0.0195 0.0589 0.0483 0.041 1 0.0868 0.0735 0.0390
3 0.0279 0.0281 0.0168 0.0291 0.0285 0.0288 0.0339 0.0296 0.0177
4 0.0166 0.0166 0.0099 0.0181 0.0182 0.0174 0.0179 0.0163 0.0093
5 0.0123 0.0104 0.0075 0.0097 0.0103 0.0105 0.0098 0.0089 0.0070

 

 

 

 

 

 

 

 

 

 

 

6.8 Discussion

Flaw depth proﬁling results from the different data sets demonstrate the
feasibility of the proposed inversion algorithm, even in the presence of noise. The
parameters that impact the performance of the proposed technique include the state
vector length L and the number of samples. Increasing L results in a larger

neighborhood and provides improved ability to track the ﬂaw proﬁle (as indicated by
the improvement in accuracy with increasing L). However, larger values of L increase
the computational cost.

Results on simulated eddy current data seem to indicate that the error in the
solution is somewhat dependent on the ﬂaw length (error reduces with increasing ﬂaw
length) A plausible explanation for this is that the NDE measurement is less sensitive to
shorter ﬂaw lengths [87]. Moreover, the size of the neighborhood will affeCt the results
using larger L values more than those with smaller L values. As described above, the
prior density of the ﬂaw depth depends on ﬂaw depths in a neighborhood deﬁned by L.
Therefore, sampling from the prior density would also give better results for longer
ﬂaws. The results also indicate that the use of data fusion improves the accuracy, with

the use of PCA before fusion further improving the results.

91

The proposed technique is seen to be robust to additive noise. While the
accuracy of the result degrades in the presence of noise, the degradation is seen to be
low. Average degradation in inversion accuracy with simulated eddy current data
(Table 6.2) at SNR = 2 with respect to non-noisy data is 30%, with inversion accuracy
in this case measured in terms of the MSE between PME estimates and true proﬁles.
The impact of increasing L, and fusing information from other measurement modes, is
particularly felt when dealing with noisy measurements, with these steps contributing to
a signiﬁcant increase in solution accuracy.

Errors are seen to be higher when dealing with experimental data (when
compared to those from simulated data). One reason for this increase is the complexity
of the shape of natural cracks, with abrupt variations in depth. This fact makes the
tracking of the ﬂaw proﬁle by the proposed algorithm a challenging proposition. The
use of a larger neighborhood improves inversion accuracy, possibly by providing a
regularizing effect on the inverse problem. However, the improvement is not as
pronounced as it is for simulated data. Again, the effect of including information from
different modes is seen to be beneﬁcial.

The inversion results on steam generator tube inspection experimental data are
also compared with the results of neural network based ﬂaw proﬁling technique
(Section 6.5). As shown in Table 6-9, the inversion results using the proposed technique
are better in most cases. It is expected that the inversion results using the particle ﬁlter
may improve further if a larger training database is used to determine the measurement
model parameters. The application of the proposed technique on additional data sets

(including a two-dimensional proﬁling problem) shows the capability of the technique

92

to handle a diverse set of problems and measurements. In particular, the inversion
results using data fusion for corrosion quantiﬁcation are comparable to (or better than)
other contemporary techniques reported in the literature [87].

The particle ﬁlter algorithm may also be applied in a non-recursive manner to
determine the complete ﬂaw depth proﬁle from the complete set of measurements. In
the non-recursive approach, each particle corresponds to a potential ﬂaw proﬁle, and the
measurement model may be used to estimate the weight assigned to each particle. The
overall PDF obtained from this problem setup can then be used to compute the most
likely proﬁle given the measurements. This problem formulation may also be
encapsulated by the recursive formulation proposed with the state vector length set to K
(where K is the number of positions at which the depth must be evaluated). The
drawback to this approach is the resulting high dimensionality of the state vector, and
the corresponding high dimensional PDF. This in turn will require (in general) a large
number of samples (typically larger than in a recursive approach) to accurately
represent the PDF. Furthermore, the computational cost of this technique is expected to
be higher (in part due to the necessity for a highly accurate measurement model in high
dimensional space). For this reason, the use of particle ﬁltering in a recursive approach
is preferred. 1

The recursive approach proposed is computationally more demanding than a
simple direct inversion approach (such as those using neural networks). On the other
hand, the approach is signiﬁcantly faster than iterative methods that use numerical
models. For instance, the ﬁnite element models that were used to generate the data used

in this study took (on average) about 31.28 seconds (using a 2.4 GHz processor, with 1

93

GB RAM) to compute the measurements given the ﬂaw depth proﬁle. Using these (or
similar) models in an iterative fashion will increase the computational load (since these
models will need to be solved during each iteration). On the other hand, the proposed
technique (with polynomial measurement models, and simulated on the same processor
in a MATLAB environment) took 23.77 seconds for evaluating the same ﬂaw depth
proﬁle (value of L used was 0). These results indicate the lower computational cost
associated with the proposed algorithm. The results further indicate that particle ﬁlter-
based inversion using polynomial based measurement models are faster than particle
ﬁlter-based inversion using neural network based measurement models. A very high
computational cost of employing a numerical model as measurement models precluded

their use in this research.

94

CHAPTER 7

CONFIDENCE METRICS

7.1 Introduction

A sequential Monte Carlo based technique is presented to solve the NDE inverse
problem in the preceding chapters. In this chapter, we discuss conﬁdence metrics to
better evaluate the accuracy and efﬁciency of the proposed technique. The metrics can
be used as a benchmark for comparison to other inversion techniques. These conﬁdence

metrics include the Cramer-Rao lower bound and credible intervals.

7.2 Cramer-Rao lower bound (CRLB)

The CRLB provides a lower bound on the variance of estimators of a
deterministic parameter [68]. The NDE inverse problem is formulated in this
dissertation in terms of a non-linear ﬁltering problem. A closed form analytical solution
of the nonlinear ﬁltering cannot be obtained for this problem [68]. However, the best
achievable error performance can be determined in terms of the CRLB. Traditionally,
the bound has been used a bench mark for comparison of suboptimal ﬁltering

algorithms [68].

7.2.1 Computation of CRLB

The CRLB states that the variance of any estimator is at least as high as the

inverse of the Fisher information matrix [68, 89]. Any unbiased estimator which

95

achieves this lower bound is said to be efﬁcient. The non-linear problem formulated in

chapter 3 is repeated here:

j¢k
Zk =hk(xk)+uk

A
If xkl k is an unbiased estimator of the state xk based on measurement sequence and

A
prior density, the covariance matrix of xklk , denoted by Pklk , has CRLB given by

[68,89]

A A T (7.2)
Pklsz (Xka-Xk](Xk|k—Xk] 2.ng.

Matrix J k is the ﬁltering information matrix and its inverse is the CRLB. Tichavsky

[89] provided an elegant method to compute the information matrix J k recursively:

22 21 11-1 12 7.3
Jk+1=Dk —Dk (Jk+Dk) Dk. ( )

where

J0 = El [on 10g ”’81)"on log p(x0)]T }, (7.4)

96

D1161 : —E{ka [ka 10g p(xk+1 I xk )F} (7.5)

13,31 = —E{ka [ka+1 10g p(xk+1 IMF}

D1152 = _E{ka+1[vxk 10g p(xk+l ka)]T }= [D1311]: .
0132 = _E{ka+1[vxk+1 log P(xk+1 lxk )IT}

— E{ka+1[vxk+l log P(Zk+l lxk+1)]T}

As discussed in the preceding chapters, the prior (state transition) density as
given by equation (3.4), being an exponential density, can be used in the analytical
expression (7.5). The likelihood density is a non-Gaussian density and depends upon the
measurement model used and measurement noise. The non-availability of an analytical
expression for the likelihood density is a major issue in solving equation (7.5). An
alternative computational approach is described here. First an ensemble of states (ﬂaw
depth) is created using Monte Carlo approximation by sampling from the prior density.
Then, the likelihood density is evaluated using a histogram of differences between
computed measurements (when using the measurement model) and the true
measurements for each state. The likelihood PDF is centered at the true measurement
(test data). The likelihood density is then transformed in terms of a mixture of Gaussian
densities using K—means clustering and the expected maximization algorithm [90].
Details of the expected maximization algorithm for computation of parameters of

Gaussian densities in Gaussian mixture are provided in appendix A.

7.3 Credible intervals

The credible interval (CI) is a posterior probability interval which is used for interval

estimation [8]. Bayesian credible intervals are similar to conﬁdence intervals except

97

that the credible intervals incorporate information from the prior distribution into the

estimate while conﬁdence intervals are based solely on the data.

7.3.1 Computation of credible interval
The (l — a) credible interval on posterior PDF p(X k I Z k ) is given by [8]:
U3 (7.6)

p(LB < Xk <UB)= ]p(Xk |zk)ka =(1—a).
LB

In this study, the percentile method [8] is used to evaluate the credible interval of

posterior PDFs computed from both single and multiple measurement modes. A
bootstrap [8] 1- a conﬁdence interval based on percentile method could be

constructed by ﬁnding ((1- a / 2)100)th and ((0! / 2)100)th empirical

percentiles in the histogram denoted by LL and UL respectively in ﬁgure 7-1.

 

 

O 1

4— Possible ﬂaw depths ——>

Figure 7-1. Credible intervals on posterior PDF
7.4 Results

Conﬁdence metrics, as discussed above, were computed for a selected set of

simulated and experimental eddy current data, from steam generator tubing (sections

98

6.4 and 6.5). The metrics were computed for inversion results using single as well as
multiple measurement mode data. Figure 7-2 shows the credible interval for a ﬂaw of
width 1.9 mm and depth 40% of tube thickness. Figure 7-2 (a) shows the credible
intervals when data from single measurement mode is used. Figure 7.2 (b) shows the
credible intervals when data from multiple measurement modes is used. Figure 7.3
shows the expected value of the covariance of the state as well as the CRLB for single
and multiple measurement modes for the same ﬂaw. Similar results for a simulated ﬂaw
of width (0.7 mm, 0.635 mm, 0.7 mm) and depth (40%, 80%, 60 %) are shown in
ﬁgures 7—4 and 7-5. Results for ﬂaws from the experimental database are shown in
Figures 7-6 and 7-7 (width 7 mm and maximum depth 100% of tube wall thickness),

and in ﬁgures 7-8 and 7-9 (width 4 mm and maximum depth 45%).

 

 

 

True

 

 

True

 

 

 

 

 

 

 

 

 

 

 

 

1 r l l 1 ‘ ‘ ‘ 4

 

 

(a) Single measurement mode (b) Multiple measurement modes
(Magnitude @400kHz) (200,400 & 800 kHz)

Figure 7-2. Credible intervals

99

 

I l

10°

 

—q.— Single Err-Var
—+_ Single CRLB

__..._ Fusion Err-Var
__ _¢_ Fusion CRLB I

 

 

 

10.2 E—

 

104'

 

 

] l l

Axial Length

 

Figure 7-3. Covariance of state and CRLB using single and multiple measurement
modes

 

True . . .
_ PME Tme

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(a) Single measurement mode (b) Multiple measurement mode
(Magnitude @400kHz) (Magnitude @200,400 & 800 kHz)

Figure 7—4. Credible intervals

100

 

__..._ Single Err-Var
. f 4, Single CRLB

100 H... Fusion Err-Var
; __,_ Fusron CRLB

 

 

 

Ill]

10-2:

10-4»

 

 

l 1

Axial Length

 

—n
N

Figure 7-5. Covariance of state and CRLB using single and multiple measurement
modes

101

 

 

 

 

 

 

 

 

  
 

True
PM E

 

True
PME

 

     

 

 

 

 

 

 

 

 

 

 

.__*._ UL _*__ UL
__.*._ LL 1 I I L + LL 1 1 1 1
(a) Single measurement mode (b) Multiple measurement modes

Figure 7-6. Credible intervals

 

+ Single Err-Var
—+—— Single CRLB
——IiI— Fusion Err-Var
—+— Fusion CRLB

 

 

 

 

 

 

 

 

Figure 7-7. Covariance of state and CRLB using single and multiple measurement
modes

102

 

 

xk/ \:_\

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

WA. f.
1
True
True
, PME , _ STE
+UL —*_
—-1-—LL +LL
(a) Single measurement mode (1)) Multiple measurement mode
(Magnitude @300kHz) (Magnitude and Phase 100,200 & 300 kHz)
Figure 7-8. Credible intervals
100 : I I r l l I
:r + Single Err-Var
j __,__ Single CRLB
-2- _,.,_ Fusion Err-Var
10 r I Fusion CRLB
W1 1 1 1 1 1 1 —1——1 ';
1 0.4 h J J 1 l L 1
2 . 5
Axral Length (mm)

Figure 7-9. Covariance of state and CRLB using single and multiple measurement
modes

103

7.5 Discussion

As seen from the results, the range of the credible intervals decreases when
measurement data from multiple measurement modes is fused. Taken with the results in
the previous chapter, this “tightening” of the intervals may be taken as an indication of
improved accuracy. Similarly, data fusion is seen to also improve the covariance of
error and the CRLB. This behavior is similar for both simulated and experimental data.
However, from this limited study, it is clear that the variance of the estimation errors do
not reach CRLB in all the cases. This suggests room for improvement in the inversion

process.

104

CHAPTER 8

CONCLUSIONS AND FUTURE WORK

8.1 Research Contributions

In this PhD dissertation, a sequential Monte Carlo based technique was proposed
to solve NDE inverse problems. A novel formulation for NDE inverse problems was
described to facilitate the application of a recursive Bayesian inversion technique. Prior
and likelihood densities required for the statistical inversion were also formulated in
terms of state transition and measurement models. The posterior density of the unknown
parameter (which was taken to be the ﬂaw depth in this dissertation, but could be any
other characteristic of the material) was then evaluated at each spatial location
recursively. Finally the estimates were computed from the posterior PDFs. The
proposed technique was applied to the problem of ﬂaw proﬁle estimation from

simulated and experimental eddy current NDE data and MFL data.

A second contribution of this PhD dissertation is a novel sequential Montecarlo
based data fusion technique along with principal component analysis. While the focus
of this dissertation was on NDE inverse problems, the proposed fusion technique can be
applied for fusing information in any SIR algorithm based particle ﬁltering application.
The proposed data fusion technique improved solution accuracy for NDE inverse
problems. The use of principal component analysis prior to fusion was seen to further

improve the inversion results in terms of accuracy. Conﬁdence metrics based on the

105

CRLB and credible intervals were also proposed to quantify the accuracy and efﬁciency

of the proposed techniques.

The proposed technique was applied to invert simulated and experimental NDE
data from a wide range of applications. Results indicated that the proposed technique is
computationally efﬁcient and produces accurate inversion results. Results on simulated
data with noise, as well as experimental data, also showed that the technique is robust to
noise. The use of the proposed conﬁdence metrics was seen to provide a mechanism to

evaluate the accuracy of the result.
8.2 Future Work

Future work will focus on further developing the data fusion technique when the
measurement data are not independent. As discussed in chapter 5, particle ﬁlter based
data fusion is used with the assumption that different measurements modes are
independent of each other. In practice the measurement modes may be correlated, and
the correlation information can potentially be used to improve the results. In the future
work, covariance intersection/covariance union algorithms and copula models [91] will
be considered for incorporating the correlation information. The impact of algorithm

parameters (such as L, number of samples Ns, measurement model order R as well as
the model itself, and better techniques for computing p(z k + 1 I xk+1)) on the

conﬁdence metrics will be studied. The application of this technique on other data sets
is another focus area for future work. These data sets include X-ray tomographic data

and ultrasonic data.

106

APPENDIX

EXPECTATION-MAXIMIZATION (EM) ALGORITHM FOR
CALCULATION OF GAUSSIAN MIXTURE DENSITIES
PARAMETERS

A-1.1 EM Algorithm

The EM algorithm [82] is used for ﬁnding maximum likelihood estimates of
parameters in probabilistic models, where the model depends on unobserved latent
variables. EM is an iterative method which alternates between performing an
expectation (E) step and a maximization (M) step. The expectation step computes an
expectation of the log likelihood with respect to the current estimate of the distribution
for the latent variables. The maximization step computes the parameters which
maximize the expected log likelihood found in the E step. These parameters are then

used to determine the distribution of the latent variables for the next E step.

A-1.2 Finding Gaussian Mixture Densities Parameters

Let X be a sample of independent observations from a mixture of m = l I M

N

Gaussian densities with mean and variance 6m = (,um , 2m) , and let Y = {yi } i=1

be the latent variables (unobserved data items) that determine the component from

which the observation originates. It is assumed for eachi , yi = m if the ith sample

107

was generated by the mth mixture component. To begin, assume the following

probabilistic model [89]:

M (A-1.1)
p<x|®>= Zampmmﬁr).

m=1

where 6,. =(ﬂpz1), and p,(x|6,.) is a Gaussian density. The aim is to estimate the
unknown parameters 9 = (a1, ..... ,aMﬂl, ...... ,6M)representing the "mixing" value a

between the Gaussians and the means and standard deviations of each density:

9 =(aj=1:M’6j=1:M)- (A'l-z)

The log likelihood function is deﬁned as [88]:

N (A-l.3)
log(L((~) | X,Y) = log(p(X,Y | O) = Zlog(ayi pyi (x,- I6Yi )).

I:

We assume an initial guess for the unknown parameters of the mixture densities as
98 .Given 98 we can easily compute p j (x,- I 615 )for each i and j . The mixing
parameters a j can be thought of as prior probabilities. Therefore using Bayes rule [89],

a .p .(x- 6’ .) a .p .(x- 19 .) (A-1-4)
P(Yilxi’@g)= Yr Yr 1' Yr = Yr Yr 1' Yr .

. g M
”“119 ) 2a§pk<xilef>
k=1

 

 

108

and

N (A-1.5)
p(YIX,(")g)= Hpm Ixzﬂg).
i=1
E—step
The E-step results in the function:
Q(9 I 9‘”) = E110g<L<®I X.Y»p<Y I Mag 1 (“6)
Equation (A—1.6) can be simpliﬁed and written as:
(t) M N (A-l.7)
Q(@|9 )= Z Zlog(arpr(xr191))P(l|xi,@g)
l=li=l
M N
= Z Zlog(ar)p(l 16L))P(l 1311,93)
l=li=1
M N
+ 2 2103001 (xi |91))P(l | 161,98)
l=1i=1
M-step

We can maximize the term containing a, and the term containing 01 independently in
(A-1.7) since these terms are not related. Introducing Lagrange multipliers and a

constraint that Z a, = 1, the mixing parameters can be evaluated as:
l

109

alrrew ___

[V12

1 g (A-1.8)
— l x-,® .
Ni p(| r )

l

The update equations for Gaussian density parameters 61 = (111,01) can also be

evaluated using [89]:

 

N (A-1.9)
inp(l|xr,@g)
new _ i=1
#1 — N ’
ZIP(1 I xi19g)
i=1
and
N g new new T (A-1.10)
2100196119 )(xi ‘1“, X961 “/11 )
Znew = i=1
1

 

N
Zpalxrﬂg)
i=1

110

REFERENCES

[1]

[2]

[31

[4]

[5]

[6]

[7]

[8]

[9]

[101

[11]

[12]

[13]

[14]

Tarantola, A. (1987). Inverse Problem Theory: Methods for Model Parameter
Estimation Elsevier.

 

Bertero, M. and Boccacci, P. (1998). Introduction to Inverse Problems in
Imaging CRC Press.

Hadamard, J. (1902). Sur les problemes aux derivees partielles et leur
signiﬁcation Physique, Princeton Universig Bulletin, 13.

Vogel, C. (2002). Computational Methods For Inverse Problem;
Philadelphia,PA: SIAM.

Weile, D. and Michielssen, E. (1997). Genetic algorithm optimization applied to
electromagnetics: A review,” IEEE Trans. Antennas Propagat,. 45, 343-353.

Kiapio, J. and Somersalo, E. (2004). Statistical and Computational inverse
problems, Springer.

Mossegard, K. (2006). Monte Carlo Analysis of Inverse Problems. J & R
Frydenberg A/S — Copenhagen N. Denmark.

Givens, G. and Hoeting, J. (2005). Computational Statistics, Wiley Series.

Norris. J. (1998). Markov Chains Cambridge University Press.

 

McMaster, R. (1985). The origins of electromagnetic testing. Materials
evaluation. Handbook of Non-destructive Testing, 43, 8, Columbus, OH:
American Society of Non destructive testing.

Collins, R. (1995). Nondestructive Testing Handbook of Materials, Netherlands:
103 Press.

 

Ida, N. (1983). Three dimensional ﬁnite element modeling of electromagnetic
non-destructive testing phenomerg, Ph.D. dissertation, Colorado State Univ.,

Fort Collins.

 

Sylvester, P. P. and Ferrari, R. L. (1990) Finite elements for electrical
engineers, Cambridge, University Press, Cambridge, England.

Zhou, P. (1993) Numerical analysis of electromagnetic ﬁelda, Springer-Verlag,
Berlin.

111-

[15]

[16]

[17]

[18]

[19]

[20]

[21]

[22]

[23]

[24]

[25]

[261

Udpa, L. and Udpa, S. (1997). Neural application of signal processing and
pattern recognition techniques to inverse problems in NDE, Int. J. Appl.
Electromagn. Mech, 8, 99—117.

(2003). A review of structural health monitoring literature: 1996-2001, Lg
Mos National Laboratory. Tech. Rep. LA-l3976-MS.

Ramuhalli, P. Udpa, L. and Udpa, S. (2002).E1ectromagnetic NDE signal
inversion using function approximation neural networks, IEEE Trans. Magn, 38,
6, 3633—3642.

Ramuhalli, P. Udpa, L. and Udpa, S. (2001). Neural network algorithm for
electromagnetic NDE signal inversion, Electromagnetic Nondestruc.Eval. (\Q.
121—128,Amsterdam, The Netherlands: IOS Press.

Cacciola, M. Forresta, F. Morabito, F. and Versaci, M. (2007). Advanced use of
soft computing and eddy current test to evaluate mechanical integrity of metallic
plates, NDT & E international, 40.

Yan, M. Afzal, M. Udpa, S. Mandayam, S. Sun, Y. Udpa, L. and Sacks,
P.(1998) Iterative algorithms for electromagnetic NDE signal inversion in

ENDE (II), Series: Studies in Applied Electromagaetic and Mechanics, 14, 287-
296 Eds. Amsterdam, The Netherlands: 108 Press.

Kirkpatrick, S, Gelatt, C and Vecchi, M (1983) Optimization by simulated
annealing, Science, 220, 4598, 671—680.

 

Haupt, R. An introduction to genetic algorithms for electromagnetic, IEEE
Trans. Antennas Propag., 37, 2, 7—15.

Li, Y. Udpa, L and Udpa, S. (2004). Three dimensional defect reconstruction
from eddy-current NDE signals using a genetic local search algorithm, IEEE
Trans. Magn., 27, 2, 410—417.

Hoole, S. Subramanian, S. Saldanha, R. and Coulomb, J. (1991). Inverse
problem methodology and ﬁnite elements in the identiﬁcation of cracks,
sources, materials, and their geometry in inaccessible locations, IEEE Trans.
Magn., 27, 3, 3433—3443.

Arunachalam, K. Melapudi, V. Rothwell, E. Udpa, L. and Udpa, S. (2005),
Microwave NDE for reinforced concrete in Proc. Ouapt. Nondestructive. Eva].

Zeng, Z, Lu, C, Shanker, B and Udpa, L (2006). Element-free Galerkin’s

method in modeling microwave inspection of civil structures, in Proc.12th
Biennial IEEE Conf. Electromagn. Field Computation, 27, 268—1268.

112

[27]

[28]

[291

[30]

[31]

[321

[33]

[34]

[351

[36]

[371

[38]

Arunachalam, K (2006), Investigation of a deformable mirror microwaLe
imaging and therapy technique for breast cancer, Ph.D. dissertation, Michigan
State Univ., East Lansing.

Kim, Y, Jofre, L, Flaviis, F and Feng, M (2003) Microwave reﬂection
tomographic array for damage detection of civil structures, IEEE Trans.
Antennas Propagation , 51,11, 3022—3032.

Li, Y. (2001), An adjoint equation based method for 3D eddy current NDE
signal inversion, in Electromagnetic Nondestructive. Eval. (V).Amsterdam, The
Netherlands: IOS Press, 89—96.

Monebhurrun, V. Duchene, B. and Lesselier, D. (1998). 3D inversion of eddy
current data for nondestructive evaluation of steam generator tubes, Inverse
Problems, 14, 707—724.

Morozov, M. Rubinnaci, G. Tamburrino, A. and Ventre, S. (2006) “Numerical
models of volumetric insulating cracks in eddy-current testing with experimental
validation, IEEE Trans. Magn., 42, 5, 1568—1576.

Balasubramaniam, S. Shanker, B. and Udpa, L. (2001). A fast integral equation
based scheme for computing magnetostatic ﬁelds and its application to NDE
problems, Proc. Quant. Nondestructive Evaluation, 20, 331—337.

 

Schiﬁni, R. and Bruno, A. (2005). Experimental veriﬁcation of a ﬁnite element
model used in a magnetic ﬂux leakage inverse problem, J. Appl.ths.. 38,
1875—1880.

Caorsi, S. Gragnani, G. Medicina, S. Pastorino, M. and Zunino, G.
(1994). Microwave imaging based on a Markov random ﬁeld model, IEEE
Trans. Antennas Propag., 42, 2, 293—303.

Norton, S. and Bowler, J. (1994). Theory of ECT inversion, J. Appl. Phys, 73,
501—512. -

J. Bowler. (1994). Inversion of eddy current data to reconstruct ﬂaws by an
optimization process, Int. J. Appl. Elecuonargn. Mat_e_;, 4, 277-284.

Chen, Z. and Miya, K. (1998). ECT inversion using a knowledge-based forward
solver, J. Nondestruc. Eval_., 17, 167—175.

Zaoui, F. Marchand, C. and Pavo, J. (2001) Stochastic crack inversion by an

integral approach, in lﬂectromagnetic Nondestruc. Eval. (V). Amsterdam, The
Netherlands: IOS Press, 129-136.

113

[39]

[401

[41]

[42]

[43]

[44]

[45]

[46]

[47]

[48]

[49]

[50]

[51]

[521

loan, D. Mihai, M. and Duca, A. (2001). Use of evolution agents to solve ENDE
inverse problems, in Electromagnetic Nondestruct. Eval.(V). Amsterdam, The
Netherlands: 108 Press, 59—66.

Tamburrino, A. (2000). A communication theory approach for electromagnetic
inverse problems, IEEE Trans. Maga, 36, 4, 1136—1139.

Barrett, H. and Myes, K. (2004). Foundation_s of Image Science. New York:
Wiley-Interscience.

Kim, Y. Jofre, L. Flaviis, F. and Feng, Q. (2003). Microwave reﬂection using
tomographic array for damage detection of civil structures, IEEETrans.
Anteraras Propag._, 51, 11, 3022—3032, 2003.

Potthast, R. (2001). Point sources and multipoles in inverse scattering theory, in
Point Sources and Multipoles in Inverse Scattering Theory. Boca Raton, FL:
Chapman & Hall/CRC.

Colton, D. and Kirsch, A. (1996). A simple method for solving inverse
scattering problems in the resonance region, Inverse Probl. , 12, 383-393.

Colton, D. and Kress,R. (1998). Inverse Acoustic and Electromagnetic
Scattering Theog. New York: Springer-Verlag.

Nikalova, M. Idier, M. and Djafari, MA. (1998). Inversion of large support ill-
posed linear operators using a piecewise Gaussian MRF, IEEE Trans. Image
Process. 7, 4, 571—585.

 

Vogel, C. (2002). Computational Methoda for Inverse Problema, SIAM.

Colton, D. and Kress, R. (1998), Inverse agoustic and electromagnetic scattering
theog, Springer-Verlag.

Dogandzic, A. (2007). Bayesian NDE defect signal analysis, IEEE Trans. Sig.
Proc., 55,1,372-378.

Dogandzic, A. and Zhang, B. (2007). Markov chain Montecarlo defect
identiﬁcation in NDE images.

Liu, Z .et al. (2007). Survey: State of the art in NDE data fusion techniques,”
IEEE Trans. Instrum. Meas., 56, 6, 2435-2451.

Gross, X.E. (1997). NDT Data Fusion, Butterworth-Heinemann, Boston.

114

[53]

[54]

[55]

[56]

[57]

[581

[591

[60]

[61]

162]

[63]

[64]

[65]

[66]

Yim, J. Udpa, S. Mina, M. and Udpa, L. (1996). Optimum ﬁlter based
techniques for data fusion. Review of Progress in ONDE. 15, D. O. Thompson
and D. E. Chimenti, Eds. New York: Plenum, 773-780.

Yim, J. Udpa, S. Udpa, L and Lord, W. (1995). Neural network approaches to
data fusion. Review of Progress in ONDE. 14, D. O. Thompson and D. E.
Chimenti, Eds. New York: Plenum, 819-826.

Simone,G. and Morabito, RC. (2001). NDT image fusion using eddy current
and ultrasonic data. Int. J. Comput. Math. Elect. ﬂectron. Eng” 20, 3, 857—868.

Dion, J. Kumar, M and Ramuhalli, P. Multi-sensor data fusion for high-
resolution material characterization, Review of Progress in ONDE. AIP
Conference Proceedings , 894, 1189-1196.

Gros, X. E. Strachan, P. and Lowden, D.W. (1995). Theory and implementation
of NDT data fusion. Res. Nondestruct. Eval., 6, 4, 227—236.

Gros, X. B. Liu, Z. and Hanasaki, K. (2000).Experimenting with pixel level
NDT data fusion techniques. IEEE Trans. Instrum. Meas., 49, 5, 1083—1090.

Kaftandjian, V. and Francois, N. (2002). Use of data fusion methods to improve
reliability of inspection: Synthesis of the work done in the frame of a European
thematic network, Proc. 8thECNDT , Spain,8. NDT. net, 2.

Liu, Z. et a1. (1998). One-dimensional eddy current multi-frequency data fusion:
A multi-resolution analysis approach. Insight, 40, 4, 286—289.

Y. W. Song and S. S. Udpa, A new morphological algorithm for fusing
ultrasonic and eddy current images. Proc. IEjE_E Ultrason. Supp, 649—652.

Liu, Z. et a]. (1999). 3D visualization of ultrasonic inspection data by using
AVS. Proc. 5th Far-East Conf. Nondestructive Testing, Taiwan, 549—554.

Forsyth, D.S. et a1. (2000). The role of data fusion in NDE for aging aircraft.
Proc. SPIE, 3394, 47-58.

Sun, K. et a1. (1996). Registration issues in the fusion of eddy current and
ultrasonic NDE data using Q-transforms. Review of Progress in ONDE. 15, D.
0. Thompson and D. E. Chimenti, Eds. New York: Plenum, 813—820.

Hall, L. and Llina, J. (2001). _H_andbool_< of Multisensor Data Fusion. CRC Press.
Nguyen, N. Milanfar, P. and Golub, G. (2001). A computationally efﬁcient

super resolution image reconstruction algorithm, IEEE Trans. Image Processing,
10, 573-583.

115

[67]

[68]

[69]

[70]

[72]

[73]

[741

[75]

[76]

[77]

[78]

[79]

[80]

[81]

Alves, N. (1990). Ising model Monte Carlo simulations: density of states and
mass gap, Phys. Rev. B 41 383.

Ristic, B. Arulampalam, S. and Gordon, 1. (2004). Beyond the Kalman Filter.
London, UK: Artech House.

S. Arulampalam, S. Maskell, N. Gordon, and T. Clapp, “A tutorial on particle
ﬁlters for on-line non-linear/non-Gaussian Bayesian tracking,”IEEE Trans.
Signal Process, vol. 50, pp. 174—189, 2002.

(2002). Special issue on Monte Carlo methods for statistical signal processing,
IEEE Trans. Signal Proc., 50.

Gordon, N. J. Salmond, DJ. and Smith, A. F. M. (1993). Novel approach to
nonlinear/non-Gaussian Bayesian state estimation, IEEE Proc.-FzRadar and
Signal Proc,l40, 107-113.

Doucet, A. Godsill, S. and Andrieu, C. (2000). On sequential Monte Carlo
sampling methods for Bayesian ﬁltering, Saatistics and Computing, 10, 197-
208.

Carpenter, J. Clifford, P. and Fearenheard, P.(1999). Improved particle ﬁlter for
non-linear problems,” IEEE Proc. Part F: Radar and Sonar navigation, 146, 2-7.

Kong, A. Liu, J. and Wong, W. (1994). Sequential imputations and Bayesian
missing data problems. Journal of American Sgtistical Association, 89,425,278-
288.

 

Khan, T. and Ramuhalli, P. (2009). Sequential Monte Carlo Methods for
Electromagnetic NDE Inverse Problems—Evaluation and Comparison of
Measurement Models, IEEE Trans. Mag.,45, 3, 1566-1569.

Specht, D. F. (1991). A general regression neural network, IEEE Transact.
NeuraLlNetworlgs, 2, 568-576.

Haykin, S. (1994.) Neural Networks : A comprehensive foundation. Macmillan
College Publishing.

Sobczyk, K. Kirkner, D. (2001). Stochastic modeling of Microstructures,
Birkhauser, Boston.

Jolliffe, LT. (2002). Principal component analysis, Springer.

Morris, C. N. (1983). Parametric empirical Bayes inference: theory and
applications. J. Amer. Statist. Assoc. 78, 47-65.

116

[82]

[83]

[84]

[85]

[86]

[87]

[88]

[89]

[90]

[91]

Little, R. J. A., and Rubin, D. B. (1987). Statistical Analysis with Missing Data.
John Wiley, New York.

Allen, D. M. (1970). Mean square error of prediction as a criterion for selecting
variables. Technometrics 13,469-475.

Xuan, L. Zeng, Z. Balasubramaniam, S. and Udpa, L. (1999). Element-Free
Galerkin Method for Static and Quasi-Static Electromagnetic Field
Computation, IEEE Trans. Magn., 40,12-20.

(2003). Rotating Probe Eddy Current Data Analysis, EPRI, Palo Alto, CA.

Udpa, L. Ramuhalli, P. Benson, J and Udpa, S. (2004). Automated analysis of
eddy current signals in steam generator tube inspection, 16th World Conference
on Nondestructive Evaluation. Canada.

Liu, Z. et al. (2008). Combining multiple nondestructive inspection images with
a generalized additive model. Miasurement Science and Technology, 19,8.

Rummel, W. D. (1998). Probability of Detection as a Quantitative Measure of
Nondestructive Testing End-to—End Process Capabilities. Materials Evaluation,
56, 1, 29-35.

Tichavsky, P. Muravchik, C. Nehorai, A. (1998). Posterior Cramer—Rao bounds
for discrete-time nonlinear ﬁltering, IEEE Trans. Signal Processing ,46
(5),]386—1396.

Bilems, A. J. (1998). A gentle tutorial on the EM algorithm and its application to
parameter estimation for Gaussian mixture and hidden Markov models.

Technical Report ICSI-TR-97-02. University of Berkelgr.

Herubiniand, U. Luciano, E. and Vecchiato, W. (2004). Copula Methods in
Finance, John Wiley & Sons, New York.

117

       

MIC

11111111111111111111111“