Imaging of defects in girth welds using inverse wave field extrapolation of ultrasonic data

(1)

Inverse Wave Field Extrapolation of

Ultrasonic Data

(2)

Imaging of Defects in Girth Welds using

Inverse Wave Field Extrapolation of

Ultrasonic Data

PROEFSCHRIFT

ter verkrijging van de graad van doctor

aan de Technische Universiteit Delft,

op gezag van de Rector Magnificus prof. dr. ir. J.T. Fokkema,

voorzitter van het College voor Promoties,

in het openbaar te verdedigen

op dinsdag 6 november 2007 om 10:00 uur

door

Niels P ¨

ORTZGEN

(3)

Toegevoegd promotor: Dr. ir. D.J. Verschuur

Samenstelling promotiecommissie:

Rector Magnificus, voorzitter

Prof. dr. ir. A. Gisolf, Technische Universiteit Delft, promotor

Dr. ir. D.J. Verschuur, Technische Universiteit Delft,

toegevoegd promotor

Prof. dr. ir. C.P.A. Wapenaar, Technische Universiteit Delft

Prof. dr. H.P. Urbach, Technische Universiteit Delft

Prof. dr. A. Erhard, BAM Berlin

Dr. ir. G. Blacqui`ere, TNO

Dr. ir. M. Lorenz Shell Global Solutions

ISBN 978-90-9022387-2

Copyright c 2007, by N. P¨ortzgen, Laboratory of Acoustical Imaging and Sound Control, Faculty of Applied Sciences, Delft University of Technology, Delft, The Netherlands. All rights reserved. No part of this publication may be reproduced, stored in a retrieval system or transmitted in any form or by any means, electronic, mechanical, photocopying, recording or otherwise, without the prior written permission of the author N. P¨ortzgen, Faculty of Applied Sciences, Delft University of Technology, P.O. Box 5046, 2600 GA, Delft, The Netherlands.

SUPPORT

The research for this thesis was financially supported by RTD bv. and by a subsidy within the Dutch ’Technologische Samenwerking’ program.

Typesetting system: LA_TEX.

(4)

To my sons Thijs and Daan

(5)

(6)

1

Introduction

1.1 Non-destructive inspection techniques

Non-destructive inspection (NDI) methods are used in a wide range of application fields. In general, the objective of NDI methods is to inspect the inside of an object without damaging it. NDI makes a substantial contribution to safety, the economy and the en-vironment. In many cases, imperfections like small cracks, corrosion or imperfections in welds caused near accidents or accidents which lead to sometimes environmental disasters, economic failure and even fatalities. Therefore, governments developed regulations to set a standard for the quality of newly constructed objects in the petrochemical industry such as pressure vessels, pipe lines and storage tanks, or in the aviation industry such as aero-plane wings. In addition, regulations were made for the maintenance of those objects.

A simple NDI method is visual inspection. Just by carefully looking at an object something can be said about the quality of the manufacturing of it while the shape and condition of the object remains unchanged. A shortcoming of visual inspection is, of course, that in many cases the interior of an object cannot be inspected. Other non-destructive methods have been developed to measure the condition of the interior of the object from the outside. Those NDI methods are usually based on:

• Radiography

(11)

All these methods can be used to measure parameters that are related to the physical prop-erties of materials. The applications for which a certain NDI method is suitable depend on the physical principles involved. In general, all NDI methods have got their weak and strong points concerning inspection speed, resolution, safety, accuracy and probability of detection of inhomogeneities. A useful table that summarizes the advantages and disad-vantages of the major NDI methods can be found in Lorenz [1993].

In this thesis, the application of interest is girth weld inspection. The non-destructive technique for the inspection of girth welds is based on ultrasonic wave propagation. An ultrasonic transducer transmits a wave field that is scattered by inhomogeneities in the ob-ject of interest, after which the scattered wave field is measured by one or more receivers, which often are the same transducer.

To judge the condition of an object from these measurements, the raw data has to be pro-cessed interpreted. After data analysis it can be concluded if an inhomogeneity is present. An inhomogeneity may not necessarily be unacceptable. Throughout this thesis, the term ’defect’ will be used for all unacceptable inhomogeneities. A thorough understanding of the physical principles of the used NDI method will improve the results of the interpretation.

1.2 Conventional defect sizing

In ultrasonic NDI a wave field usually is generated with the use of piezo electric crystals or composites. When such a crystal is exposed to a mechanical vibration, an electric potential is generated. Vise versa, a mechanical vibration is generated when the crystal is subjected to an electric potential. Crystals with these characteristics are called transducers. In practice, the crystal is exposed to a short potential pulse causing the crystal to vibrate with a frequency bandwidth and directivity pattern characteristic for the crystal design. Directly after the pulse, the potential over the crystal is measured during a limited time. When waves reflect or diffract at medium boundaries, a response signal is recorded. This signal is the basis of any ultrasonic NDI measurement. In NDI the recorded signal is called an A-scan, where in seismic exploration it is often called a trace. The physical parameter that is measured directly is a voltage caused by the strain in the piezo-electric material. The amplitude of the signal usually is related to a reference value and presented as a percentage of the full screen height (FSH).

When ultrasonic responses are measured the results have to be evaluated. A response can be caused by an inhomogeneity or by the boundaries defined by the geometry of the object (see figure 1.1). In case of ultrasonic inspection techniques, responses caused by the ir-regular geometry of boundaries can usually be identified. Before the inspection starts, the geometry of the weld and the materials the weld and the pipe consist of, are known. Loca-tions that cause responses are, for example, the weld reinforcements (cap and root), clad layers or buffer zones. Since the positions of those locations are known, the concomitant signals can be identified by their travel times. Sometimes the travel time from a defect to the receiver is almost the same as the travel time from a boundary to the receiver. In this case the defect’s responses are masked.

(12)

1.2 Conventional defect sizing 7

cap geometry

lack of fusion defect

Figure 1.1: Cross section of a narrow gap weld. The cap of the weld can cause ultrasonic responses. This weld contains a lack of fusion defect.

an inhomogeneity. The acceptance or rejection of the weld on the basis of the detection of an inhomogeneity depends on established codes that prescribe acceptance criteria. Most codes for weld inspection are based on radiographic inspection because this was the first NDI method that was used. The acceptance criteria in earlier codes are based on good workmanship, later codes are based on fracture mechanics or engineering critical assess-ment (ECA). Good workmanship codes are mainly meant as a guidance to contracting parties in the industry and are usually conservative. Codes based on ECA contain fitness for purpose criteria. As a consequence, defect characterization is dictated by the result of the assessment and must be as accurate as possible.

(13)

Then the dimensions of the reference defect with the corresponding amplitude obtained from the curve, are use as defect size. Amplitude based sizing has got some disadvantages: • By using the dimensions of a corresponding reference reflector, the assumption is made that the shape of the defect is identical to the shape of the reference defect. • The amplitude of a reflected signal is highly dependent on the orientation of the

defect. The consequence can be that a large defect under a certain orientation is accepted because the amount of received energy is much lower than the total reflected energy.

Because the assumptions made in amplitude based sizing do no always hold in practice and can result in large defect size variations, as described by Lozev [2005], acceptance criteria must be conservative.

Defect sizing with temporal information is based on the travel times from waves that are diffracted at defect tips, the most common method is the time of flight diffraction method (ToFD), see e.g. Ravenscroft [1991] and Ogilvy [1983]. The transmitting and receiving transducers are placed in a pitch-catch configuration. The travel time from source to de-fect tip to receiver contains the location information of the dede-fect. The ToFD technique is less dependent on defect orientation. When diffractions caused by the upper tip and lower tip are measured, reasonably accurate sizing is possible, depending on the frequency bandwidth of the signal. A disadvantage of the technique is the ’dead zone’ caused by the direct wave traveling just below the surface, also called the lateral wave. Cracks connected to the surface are obscured by the lateral wave (Erhard and Ewert [1999]). The ToFD technique is widely accepted and special standards are available.

In most practical situations a combination of pulse-echo techniques and ToFD techniques is used to increase the probability of detection of defects and to improve sizing by com-bining the results. Good results have been obtained with both the pulse-echo technique and ToFD techniques in controlled laboratory and field circumstances. Still, those results involve interpretation by experienced operators, because they do not directly show the de-fects location, orientation, shape and size. Pattern recognition and neural networks have been applied for uniform interpretation, such as described by Moura [2004]. However, the success of these techniques depends on training that does not involve the physics of wave propagation.

Advanced ultrasonic techniques have been studied and applied to make an image of a defect illustrating the defects characteristics. Phased array sectorial scans are used suc-cessfully in medical imaging, such as described by Wells [2000], Fenster [2001], Jesse and Smith [2002]. With the development and miniaturization of ultrasonic array equipment, sectorial scans have become popular for industrial applications (e.g. R/D-Tech [2004]). However, the same drawbacks regarding defect shape and orientation remain.

(14)

1.2 Conventional defect sizing 9

need for conservatism can be minimized by accurate failure analysis and by accurate sizing of detected defects. Minimizing the need for conservatism will lead to lower repair rates and thus substantial economic benefit. Defect sizing still dictates the extend of conservatism. In addition, with the development of more efficient and faster computers, numerical fracture mechanics analysis has improved. Ultrasonic imaging, as presented in this thesis, aims for accurate and unambiguous defect characterization and sizing. We will introduce in chapter 2 an imaging approach that has been developed in seismic exploration. With the develop-ment of ultrasonic arrays we will see that this approach becomes feasible for the imaging of defects from ultrasonic data.

One of the techniques that aims for better characterization of the inhomogeneities is the synthetic aperture focussing technique (SAFT), as described by Lorenz [1993]. The SAFT technique was originally developed for radar applications. With the SAFT technique multi-ple pulse-echo measurements are taken whereby the distances between the transmitter and receiver are varied. Two single element transducers are used for each measurement, the locations of the transducers are usually stored by using an encoder. If an inhomogeneity is present, signals caused by the inhomogeneity will arrive at different times depending on:

• the transmitter-receiver separation,

• the distance between the transducers and the inhomogeneity, • the wall thickness of the inspected object,

• the wave mode (transversal or longitudinal or indirect insonification via boundaries). The arrival times can be determined for each combination. For a point of interest all scans can be compensated for certain travel times by applying phase shifts. After the compen-sation of the travel times, the scans are added and the amplitude at zero time is used as contribution for that point in a defined image space. If the inhomogeneity was present at the chosen point, the signals caused by the inhomogeneity will all be shifted to zero time and the amplitudes will add constructively resulting in a large contribution. Usually, only one path is regarded in SAFT. Multi-SAFT, as presented by Lorenz [1993], also takes different modes into account to increase the spatial bandwidth. The orientation and size of the defect appear after the contributions of all points in the image space are calculated. The principle of filling an image space is unique in the ultrasonic application field. Despite its good potential, the SAFT technique is not regarded as a standard ultrasonic inspection technique, nor is it described in codes. In addition, the data collection and processing can be time consuming.

(15)

1 2 3 4 5 6 7 zone cap 3rd fill 2nd fill 1st fill hot pass lcp root

Figure 1.2: The automated ultrasonic girth weld inspection approach is based on zone dis-crimination. Each zone is inspected with a dedicated ultrasonic set-up.

1.3 Automated girth weld inspection based on zone

discrimi-nation

A well known application for NDI is girth weld inspection. Traditionally, girth welds were inspected with radiography. However, with advances in ultrasonic technology, ultrasonic inspection is now a well accepted inspection method. The first ultrasonic inspection meth-ods were based on moving an ultrasonic probe with a particular angle (usually 45◦_{, 60}◦

and 70◦_{) towards and from the weld in a so called ’Meander’ movement. This}

inspec-tion method has been automated resulting in a mechanized scanner with several ultrasonic probes to cover the entire weld volume. Mechanized ultrasonic inspection of girth welds can be considered as a reliable alternative for radiography (see de Raad and Dijkstra [1995] and de Raad and Dijkstra [1998]).

The inspection method for automated ultrasonic testing (AUT) is based on zone discrimi-nation for full coverage of the weld and heat affected zones (see Dub´e et al. [1998] or Findlay and van der Ent [2001]). A complete overview of AUT and zone discrimination is given by Ginzel [2006]. The weld is divided into zones of typically 1- 3 mm height. Each zone is then inspected with an dedicated ultrasonic beam, generated by a probe that is fixed in a frame called the probe pan. The probe pan is fixed to a carrier that can be mounted on a band attached to the pipe circumference. The carrier moves once around the entire cir-cumference during the weld inspection (see figure 1.3). The probes are optimized to reflect defects at the fusion face. Depending on the weld configuration, single crystal probes can be used for a perpendicular insonification, whereas steep fusion faces are inspected with a dual crystal probe in a so called tandem configuration. The near surface zones are usually not insonified by a configuration that is aimed in the direction of the fusion face. Here, the ultrasound will reflect as a result of the interaction with the defect and the surface (corner reflection).

(16)

1.3 Automated girth weld inspection based on zone discrimination 11

Figure 1.3: Automated ultrasonic testing system. The weld is inspected with zone discrim-ination, whereby each zone is examined by a dedicated ultrasonic beam. The ultrasonic probes are mounted into a frame. The frame is fixed to a carrier that moves around the circumference during inspection.

used. Once the position of the probes are optimized, the gain will be adjusted such that the amplitude reaches 80% full screen height (FSH). With the computer, a time gate will be placed over the calibration signal so that only the information in the time gate will be stored.

The data display of AUT is done with strip charts (see figure 1.5). Each column of a strip chart contains the recordings of the highest amplitude and the position within the gate of a certain probe. Some columns record the full digitized signal within the gate. These columns are referred to as mapping channels and aim to detect volumetric flaws such as porosity or inclusions. In addition, signals caused by the cap and root re-enforcements can be monitored. The top of the column represents the start position and the end of the column represents the end position around the circumference of the pipe. During an inspection the real-time view resembles a waterfall of data building up from the top of the screen to the bottom. In the presence of a defect, an amplitude will appear in the column that corresponds to the inspected zone.

(17)

flat bottom holes

notches

Figure 1.4: The sensitivity of the automated ultrasonic inspection system is set with the use of calibration reflectors in a special designed calibration block.

(18)

1.4 Ultrasonic arrays in NDI 13

Ginzel [2000], Gross et al. [2001], Morgan [2002] and Heck¨auser [2006]. The improvement of this methodology forms a substantial objective of this thesis.

1.4 Ultrasonic arrays in NDI

In section 1.2 several ultrasonic NDI techniques for detection and sizing of defects were discussed and in section 1.3 an application was presented. All techniques make use of ul-trasonic beams with a certain angle, beam spread, focal point, frequency and wave mode. When the pulse echo technique is applied, the beam characteristics are dictated by the defects of interest (a more detailed discussion is given in section 3.1). As discussed in the previous section, in an automated inspection setup for the inspection of girth welds, multiple probe configurations are use for the inspection of the entire volume of the weld. The weld configuration is usually different from project to project. This requires the man-ufacturing of new probes for almost each project. All the probes are fixed in a frame. The positions in the frame must be optimized for each probe, making the setup labor intensive. The number of probes and the setup time for automated ultrasonic testing systems can be reduced with the use of ultrasonic array technology, as described by P¨ortzgen et al. [2002], Moles et al. [2005] and Huang et al. [2004]. The characteristics of beams generated by ultrasonic arrays can be controlled with a computer, within certain limits. One array probe can even be used to replace all the single crystal probes in an automated inspection setup.

An array probe consists of typically 32 to 64 crystal elements which have a size in the order of half the wavelength of the center frequency of the crystal. The center frequency is typically between 1 and 10 MHz and the element size 0.3mm - 2mm. Each element is connected to the ultrasonic hardware and can be fired individually. The directivity pattern of a single element resembles the directivity pattern of a dipole source. When adjacent elements are fired simultaneously, the resulting wave front gets the shape of a wave front generated by a single crystal transducer with the same outer dimensions. This wave front is parallel to the crystal. When the elements are fired in a sequence of pre-determined delays, the direction and focal point of the resulting wave front can be controlled, as explained by Wooh and Shi [1999] and Lee and Choi [2000]. In figure 1.6 two examples of beam forming using delay times are illustrated.

For girth weld inspection, array probes have more elements then necessary for the gener-ation of a single beam. When a group of elements is used to generate a beam, this group can be shifted along the array( P¨ortzgen et al. [2002], Moles et al. [2003] and Moles and Labb´e [2005]). The position of the beam (or index point) can be varied while the position of the probe remains unchanged. The active group of elements may be different for the transmitting beam and the receiving beam. Tandem configurations can be constructed in this way.

(19)

Figure 1.6: Steering and focusing of an ultrasonic beam with a linear ultrasonic array. When the the active group of elements is displaced, the index point of the beam can be varied.

Figure 1.7: Sectorial scan of an unborn baby in the womb. and Harvey et al. [2002], see figure 1.7).

(20)

1.5 The use of arrays in seismic exploration 15

therefore not appear in the sectorial scan. Imaging by displaying multiple time-amplitude responses (or A-scans as they are called in the ultrasonic application field) with a color key that represents amplitude height might give the operator more insight. However, the information in the scans does not increase (P¨ortzgen [2006]).

Although ultrasonic arrays offer operational benefits and allow alternative inspection tech-niques like sectorial scans, the basic principles of ultrasonic inspection remains unchanged. In fact, the only difference with the use of phased arrays is the mechanism of transmitting and receiving the ultrasonic beams. Once the wave fronts have been formed, practically no difference can be observed between the wave front that was generated by a single crystal transducer and a phased array transducer. Consequently, the results will be similar and the same drawbacks and benefits concerning probability of detection, defect characterization, sizing and accuracy can be expected, see Armitt [2006].

Rather than using the elements of an ultrasonic array to form a wave front with certain characteristics, the individual elements can be used also as single sources and receivers. A-scans can be measured for each source receiver combination with the array in a fixed position. The data set created this way represents the most complete data set that can be gathered along the aperture of the array. All A-scans obtained from a beam created by a group of active elements can be reconstructed from this complete data set, as demonstrated by von Bernus et al. [2006]. This is possible because the A-scans from the data set can be phase shifted and added with phase shifts corresponding to any desired beam angle and focal spot. In the application field of seismic exploration, it is more common to make use of such data sets. In the next section, it is explained how images are constructed based on such data matrices in seismic exploration and how this application field can be linked with ultrasonic NDI.

1.5 The use of arrays in seismic exploration

In seismic exploration, arrays have been used for many years to image the earth’s interior. Seismic data can be processed into a representative image, with minimal assumptions and without the use of reference reflectors. The algorithms are based on fundamental wave theory and are usually referred to as migration techniques, see e.g. Schnneider [1978], Berkhout [1982] and Claerbout [1985]. Most of the imaging techniques make use of inverse wave field extrapolation with the Rayleigh integral (see e.g. Berkhout [1982]). This ap-proach will be discussed in detail in chapter 2.

(21)

Figure 1.8: 3D seismic image of a medium containing a salt dome (from the SEG/EAGE salt model presented by Burch and Burton [1984]).

the matrix contains the information of all source-receiver combination so that it is com-pletely filled. In addition, the aperture where the sources and receivers are located must be long compared to the area of interest. Images can be obtained by applying mathematical operations on the data matrix. Different implementations for the imaging techniques can be used. The most important differences with standard ultrasonic imaging techniques as used in NDI are:

• The A-scans are processed to an image after the measurements. The algorithms to produce an image are based on fundamental wave theory.

• A single point in the image contains information acquired from multiple A-scans. The A-scans were taken from different positions so that the medium is insonified from many directions.

• The images represent the reflectivity properties of the medium at the corresponding positions. As a consequence, interpretation of the final image is transparent and relatively unambiguous.

• The data contains all the medium properties (reflectivity and position of layers, reflec-tors and diffracreflec-tors, sound velocities and densities). In principle, all this information can be recovered without extensive a-priori knowledge.

• Ultrasonic techniques make use of reference block for system calibration. The re-sponse of a defect will be compared with the rere-sponse of a known calibration reflector. In seismic exploration, calibration is not standard practice, since all the variables of interest can be derived from the images.

(22)

1.6 Outline and objectives of this thesis 17

for the probes must be optimized. In seismic exploration the diversity of parameters is much less and the measurement set-up is much more standardized.

The most important similarities between seismic exploration and ultrasonic NDI are the use of sound and the use of array recorders. Although the scale is different, fundamental wave theory can be used in a similar way for both application fields. Hence, the appealing benefits of imaging that are common in seismic exploration (like unambiguous interpretation, no calibration procedure and standardized measurement set-up) can now also be explored in ultrasonic NDI.

1.6 Outline and objectives of this thesis

The main objective of this thesis is to apply the imaging concepts as developed in seismic exploration to the application field of ultrasonic NDI in 2D and to give a proof of concept in 3D. The imaging concepts will be verified with ultrasonic data obtained from simple geometrical reflectors and with reflectors that are representative for defects in girth welds, with the exception of transverse cracks.

As ultrasonic NDI is a broad application field, here we will concentrate on weld inspection of carbon steel girth welds that are common for newly constructed pipe lines. We will assume that carbon steel is a homogeneous and isotropic material. The geometry and the material constants (density and sound velocity) of the welds are assumed to be known and also typical for girth welds. However, the theory can be extended to more complex geometries and other materials like dissimilar welds.

As a consequence of the geometry and the nature of the material, defects will be insonified multiple times via boundaries and by mode converted wave fronts. It will be demonstrated that these different wave fronts have different arrival properties and can be imaged indi-vidually. Images from different arrived wave fronts may contain more information on the defect because the defect may be insonified better. For example, defects that are located near the upper surface are better insonified via the back wall.

Although the different wave fronts may lead to better images, they can also cause cer-tain artifacts in the images. It will be explained what causes these artifacts, what shape they have and where they are located. From that, a practical approach will be developed and demonstrated, to predict and suppress these artifacts.

(23)

can be used.

In chapter 3, the application of girth weld inspection and the defect types of interest will be discussed in more detail. The requirements of the imaging approach will be described in relation with girth weld inspection, based on existing codes and standards. Attention will be given to the design of the array transducer. This also involves an analysis of the reso-lution that can be achieved with the imaging approach. From this analysis, the expected performance will be compared with the requirements.

In chapter 4, results of 2D images will be presented. The images were constructed from measurements on test blocks with machined defects. The machined defects are representa-tive for defects that are common to girth welds such as embedded and surface breaking lack of fusion defects, porosity, cold lap defects and lack of cross penetration defects. Images from different insonification paths and scatter path will be presented and discussed. In chapter 5, results of 3D images constructed with the two pass method will be presented. The images were constructed from measurements on blocks with machined reflectors and with an actual weld with an embedded tungsten fragment. The machined reflectors are not very representative for real defects. However, the purpose of these reflectors is to demon-strate the characteristics of the two pass method. The real weld with the tungsten fragment can be considered to be a real defect.

In chapter 6, artifacts in images from the direct longitudinal insonification and the direct transversal scattering paths (the L-T path) caused by waves from the direct longitudinal insonification and longitudinal scatter path (the L-L path) will be discussed. An analytical example will be given to demonstrate the presence and shape of such an artifact correspond-ing to a point scatterer. The artifacts will be illustrated in images from real measurements. In addition, a procedure will be presented to suppress these artifacts.

In chapter 7 the conclusions of the research will be presented and discussed.

(24)

2

Wave theory and imaging

In this chapter, the principles of the wave equation based imaging approach are presented. The wave theory will be presented that is used for imaging based on inverse wave field extrapolation of acoustic waves. It also will be demonstrated that this theory can be applied to elastic waves. The imaging approach will be demonstrated with a numerical example. Finally, some practical implementation aspects for 2D and 3D imaging will be discussed with the aid of the matrix notation.

2.1 Introduction of the imaging philosophy

In ultrasonic NDI, the words image and imaging are often used for different treatments of ultrasonic data. Confusing is the difference in names that are used to address the imaging methods. For example, names that are used for imaging in NDI are synthetic aperture fo-cusing technique (SAFT), such as described by Mayer et al. [1990], Thomson [1984], Lorenz [1993], Burch and Burton [1984], Johnson and Barna [1983], inverse wave field extrapolation (IWEX) (see P¨ortzgen [2006] or P¨ortzgen [2004]) sampling phased array (see von Bernus et al. [2006] orChiao and Thomas [1994]) or total focussing method TFM (see Wilcox et al. [2006] or Holmes et al. [2005]).

Ultrasonic imaging techniques roughly can be divided into three categories:

(25)

whereby the amplitudes are represented with a color or gray scale. This was briefly discussed in sections 1.3 and 1.4. The mapping channels, as discussed in section 1.3, are referred to as B-scans.

[2] Techniques that remove propagation effects from source to defect and from defect to receiver. The already mentioned examples above are of this category, whereby the process to remove the propagation effects are carried out in different ways (like in the time-space domain or the frequency-space domain or by only using data with coinciding source-receiver position). The techniques of this category are all based on linear wave theory (the Born approximation). By linear wave theory we mean that the wave field is linear in the deviations of the media properties from the steel properties, of a particular weld configuration. As a consequence, waves that were caused by scattering resulting from insonification via other scatterers are not taken into account. The imaging approach as described in this thesis belongs to this cate-gory. Because our approach is based on inverse wave field extrapolation, we refer to this approach as the IWEX approach (Inverse Wave field EXtrapolation).

[3] Techniques based on non linear wave theory. These techniques aim to find the model that explains the data, see Marklein et al. [2006] and Gisolf and Verschuur [2005]. The imaging techniques of this category are usually performed with an iterative approach and are therefore calculation intensive and time consuming.

(26)

2.2 The Rayleigh II integral 21

That time equals the travel time from the source to this point. Thus, a general imaging process for a single analysis point in the image space consists of three steps:

[1] Removal of the propagation effects from an image point to the receivers. This step is performed in the space-frequency domain with the Rayleigh II integral for back propagation. After this step, for every source location the data consist of a single recording of the scattered field, by a virtual receiver located in the chosen image point.

[2] Removal of the propagation effects from the sources to the image point. When a well sampled source geometry is available, it can be demonstrated, using the argument of reciprocity, that this step also involves back propagation with the Rayleigh II integral. After this step, the full data set has been reduced to a single signal generated by a virtual source in the image point and detected by a single virtual receiver located in the same position.

[3] At the image point, the amplitude at t = 0 of the coinciding virtual source - virtual receiver recording is assigned to that point as image amplitude. In seismic explo-ration, this step is called applying the imaging condition.

The three steps must be repeated for each point in the image space. The L-L and L-T imaging modes are selected by using L-wave propagation velocity in steps 1 and 2, or by using transversal wave (T-wave) velocity in step 1 and longitudinal wave (L-wave) velocity in step 2, respectively.

In de following sections the theoretical principles involved in the imaging process as de-scribed above will be further examined.

2.2 The Rayleigh II integral

The basis of the imaging process is the possibility to extrapolate a wave field from known values at a certain surface to any location in space. Starting point of the imaging theory is the general Rayleigh II integral derived from the Kirchhoff integral (see e.g. Berkhout [1987]), formulated as, P (~rA, ω) = Z S P (~r, ω)∂G ∂ndS, (2.2.1)

where P (~r, ω) is the temporal Fourier transform of the measured pressure field p(~r, t), ~rA

is the position vector of a point A not on the observation surface S, ~r is the position vector of an observation point on S and ~n is the direction normal to the surface (see figure 2.1). Furthermore, ω is the angular frequency and G is the Green’s function. If we assume that we are dealing with a homogeneous medium, the Fourier transformed pressure P (~r, ω) in equation (2.2.1) obeys the Helmholtz equation:

∇2P (~r, ω) +ω

2

(27)

Figure 2.1: The Rayleigh II integral describes how the wave field at point A can be computed from the wave field at observation surface S.

where c is the sound propagation velocity in the medium. The Green’s function G describes the wave field due to a point source in A, which is defined in a homogeneous medium as,

G±= 1 2π∆re

∓jω

c∆r_, _(2.2.3)

where ∆r =| ~r − ~rA | is the distance between an observation point on S and the point

source in A. The Green’s function can be causal (G+) or anti-causal (G−_{). The causal}

Green’s function is the wave field of a causal point source, whereas the anti-causal Green’s function is the wave field of an anti causal point sink.

When the causal Green’s function is used in equation (2.2.1), the Rayleigh II integral can be used to extrapolate a measured wave field away from the source. On the other hand, the anti-causal Green’s function can be used to extrapolate the wave field towards the source. In our application the sources are the defects that act as reflectors and/or diffractors. If a defect is hit by an incoming wave front, it starts to act as a secondary source or distribution of sources, depending on the shape and size of the defect. By spatially back propagation the wave field due to such a source (distribution), we can find its location and produce an image of the defect. Hence, for the purpose of imaging, the measured wave field must be extrapolated in a direction closer to the source (distribution) and as a consequence we must use the anti-causal Green’s function. Furthermore, we choose the recording plane S with a normal vector in the negative z-direction. Substituting the anti-causal Green’s function, equation (2.2.3), in the general formulation of the Rayleigh II integral for inverse wave field extrapolation, equation (2.2.1), yields (Berkhout [1987]):

P (~rA, ω) = zA− z0 2π ∞ Z −∞ ∞ Z −∞ P (~r, ω)1 − j ω c∆r ∆r3 e jω c∆r_dxdy, _(2.2.4) with ∆r = p

(x − xA)2+ (y − yA)2+ (z0− zA)2and where z0indicates the location of the

(28)

2.3 Wave field extrapolation in elastic media 23 P (~rA, ω) = jω c(zA− z0) 2 ∞ Z −∞ P (~r, ω) 1 ∆rH (2) 1 (− ω c∆r)dx, (2.2.5) where H1(2)(ωc∆r) denotes the first-order Hankel function of the second kind and with

∆r =

p

(x − xA)2+ (z0− zA)2. Equation (2.2.5) represents a convolution integral over the

spatial coordinate x. A far field expression can be formulated for ω

c∆r >> 1 (Abramowitz and Stegun [1964]): P (~rA, ω) ≈ − r −j ω c 2π(zA− z0) ∞ Z −∞ P (~r, ω)√ 1 ∆r∆re jω c∆r_dx. _(2.2.6) In our application, we examine carbon steel girth welds. Carbon steel is an elastic medium, hence two wave modes, longitudinal and transversal, can exist. The different wave modes can convert to each other during transmission and reflection. In the next section, it is discussed that equations (2.2.4), (2.2.5) and (2.2.6)also can be used for imaging of wave fields in elastic media.

2.3 Wave field extrapolation in elastic media

The derivation of the Rayleigh II integrals for forward and inverse extrapolation was ob-tained with the use of the Green’s functions given by equation (2.2.3), describing the wave fields of an acoustic point source and a point sink. For girth weld inspection, we cannot assume acoustic wave fields, since the involved metals are elastic media. Shear waves as well as compressional waves will occur, which means that we need to consider the elastic wave equation. 2D elastic waves can be described in the Fourier domain by two scalar potential functions Φ and Ψ (Aki and Richards [2002]) obeying the scalar Helmholtz equations:

∇2Φ + ω 2 cp2 Φ = 0, (2.3.7) and: ∇2Ψ + ω 2 cs2 Ψ = 0. (2.3.8)

Equation (2.3.7) describes the wave field of compressional (longitudinal) waves with sound velocity cland equation (2.3.8) describes the wave field of shear (transversal) waves with

sound velocity ct. From Φ the displacement vector ~Ul of the compressional wave can be

found through ~Ul=_ρω12∇Φ, whereas the displacement vector ~Utof the shear vertical waves follows from ~Ut = _ρω12(−∂Ψ_∂z,∂Ψ_∂x). The significance of these results is that compressional and shear waves are uncoupled during propagation in a homogeneous medium and can be treated in the same way as acoustic waves. On reflection and transmission, compression waves can be converted to shear waves and vice versa. With the factor 1

ρω2 in the definition for ~Uland ~Ut, the wave potentials Φ and Ψ have the dimension of pressure. At the front

(29)

From equations (2.3.7) and (2.3.8) corresponding Green’s functions can be derived. Both compressional and shear waves can be extrapolated in space using the Rayleigh II integrals with appropriate velocities.

2.4 The imaging condition

We now define an image domain below a recording plane and we use a source located in the recording plane. Once the source is fired, for example one element of an ultrasonic array, the wave propagates downwards into the medium. The receivers record all the reflected and diffracted waves traveling upwards. The recorded wave field is now back propagated to a point A in the image domain. The back propagation removes all propagation effects from point A to the receivers. If a discontinuity is present in A, the back propagated wave field will have a non-zero amplitude at the moment that A became a secondary source, i.e. at the moment the incident wave field reached A. This principle is the so called imaging condition. Using the position of the original source location and forward propagation techniques such as ray tracing ˇCerven´y [2001], it is easy to find the moment when A became a secondary source if a scatterer was present. For homogeneous media, this travel time can be found analytically. For point A in the image domain, the amplitude of the inversely extrapolated recordings at the time that A became a secondary source is selected and used as image amplitude. This process is repeated for all grid points in the image domain. This process is called Kirchhoff migration and can be repeated for all source location, after which all sub-images can be summed.

The image amplitude is related to the scattering coefficient in point A. If the down-ward propagating wave field from the real source to point A is denoted by P↓(~r; ω) and if the scattered field propagating from the secondary source to the receivers is denoted by P↑_{(~r; ω), the scattering coefficient R(~r}

A; ω) is defined by

P↑(~rA; ω) = R(~rA; ω)P↓(~rA; ω), (2.4.9)

and the image amplitude I(~rA) is defined by

I(~rA) ≡ ∞

Z

−∞

P↑(~rA, ω)ejωtSAdω, (2.4.10)

where tSAis the travel time from the source S located at the surface to point A.

(30)

2.5 An example of imaging with simulated acoustic data 25 lenght in (mm) depth in (mm) 2 4 6 8 10 12 1 2 3 4 5 6

Figure 2.2: Model to demonstrate the imaging process. On the top of the model 128 receivers are placed. The black dot represents a circular diffractor in a homogeneous media.

2.5 An example of imaging with simulated acoustic data

To demonstrate the imaging process, we use an acoustic material with a longitudinal sound velocity cl= 2475 m/s. The thickness of the material is assumed infinite, or all boundaries

are absorbing. We place an array of 128 transducer elements at level z = 0. For this example, a 2D acoustic finite difference simulation code was used, such as described by Reynolds [1978] and Youzwischen and Margrave [1999]. The element width is 0.05 mm and the spacing between two adjacent elements is also 0.05 mm, hence the heart-to-heart distance between two adjacent elements is 0.1 mm and the total aperture is 12.8 mm. The wavelet used in the simulations has a center frequency of 4 MHz and a 50% bandwidth. The widths of the elements are small compared to the dominant wavelength so that the wave front due to one element will not have a strong directivity and can be assumed cylindrical. At level z = 3 mm a discontinuity is present below the center of the array aperture. The model is illustrated in figure 2.2.

We use element number 20 as a source, and we follow the wave field while it propagates downwards, is scattered by the discontinuity and propagates upwards back to the receivers. The sequence in figure 2.3 shows the result of the finite difference simulation. At t = 2.4 µs we can see that the wave front hits the discontinuity. The discontinuity becomes a secondary point source at this time. At t = 3.8 µs, the scattered wave field has reached the receivers of the transducer array that are above the diffractor. In figure 2.4 the recordings are displayed as a function of time.

The data recorded on all 128 receivers due to the scattering of the secondary source has got the characteristic shape of a hyperbola, the linear and stronger recordings above the hyperbola are due to the direct arrival.

(31)

(a) t = 1.0 µs (b) t = 2.4 µs, time of imaging condition

(c) t = 3.8 µs

Figure 2.3: Results of finite difference simulation of a source located at position 20. The wave front propagates downwards and reflects at the diffractor at t = 2.401 µs.

receiver number time in µ s 20 40 60 80 100 120 0 1 2 3 4 5 6 7 8

(32)

2.5 An example of imaging with simulated acoustic data 27 receiver number time in µ s 20 40 60 80 100 120 1 2 3 4 5 6 7 8

(a) Virtual receivers at depth = 1 mm receiver number time in µ s 20 40 60 80 100 120 1 2 3 4 5 6 7 8 (b) Virtual receivers at depth = 2 mm receiver number time in µ s 20 40 60 80 100 120 1 2 3 4 5 6 7 8 (c) Virtual receivers at depth = 3 mm

Figure 2.5: Recordings of figure 2.4, inverse extrapolated to different depth levels.

The sequence in figure 2.5 illustrates the back propagating process. The recorded data is back propagated to a number of virtual receiver levels with an increasing depth. At a depth of 3 mm, we see that the data is focussed at the position of the discontinuity at a virtual recording time of 2.4µs. To find the contribution to the image in that point, we must apply the imaging condition, i.e. we calculate the time at which the secondary source was activated by the source wave field from source at location number 20. For the defect at 3 mm depth this time is t = 2.4 µs. In figure 2.5c, we can see that the data focusses exactly at that time, thus giving a large image amplitude at that location.

If we repeat the imaging process for each point in the image space, we obtain an im-age based on a single source (element number = 20) and all 128 receivers. The result is presented in figure 2.6. The maximum image amplitude is found at the position of the discontinuity. In the image of figure 2.6 we can see two weak curves crossing the maximum point. These curves are caused by the finite aperture of the receiver array. These effects can be reduced strongly when we combine the images obtained from several source positions, and by spatial tapering of the data records before inverse wave field extrapolation. As an example, we repeat the imaging process with elements 10, 20, 30,..., 120 as source elements and add the resulting images. The result is presented in figure 2.7.

(33)

distance in (mm) depth in (mm) 2 4 6 8 10 12 1 2 3 4 5

Figure 2.6: Image obtained with a single source and all 128 receivers. The image shows highest energy at the location of the circular reflector from figure 2.2. Also two ‘tails’ are visible because of a limited receiver aperture. The asymmetric shape of the artifacts is due to the oblique illumination from the single source at element number 20 (at distance 2 mm).

distance in (mm) distance in (mm) 0 2 4 6 8 10 12 1 2 3 4 5 6

Figure 2.7: Stacked image obtained with images from sources 10, 20, 30,..., 120 and all 128 receivers. The ‘tails’ visible in figure 2.6 have largely disappeared.

receiver record with sources still at the acquisition surface, is inversely propagated to the same point in the image domain.

If we compare the image in figure 2.7 with the original model in figure 2.2, we can conclude that imaging by inverse wave field extrapolation is possible for this single point defect. A discontinuity with an irregular shape like a crack or lack of fusion defect, will behave as a distribution of point defects, all of which will be imaged in their correct positions. This allows a more accurate sizing and determination of the orientation of the defect.

(34)

2.6 2D imaging implementation using the matrix notation 29 z x n=1 n=N m=M m=1 ∆rm,n (xm, zA) zA (xn,0) ∆x A

Figure 2.8: The Rayleigh II integral for back propagation can be formulated in discrete matrix notation.

2.6 2D imaging implementation using the matrix notation

According to the Rayleigh II integral for inverse propagation equation (2.2.4) zero spaced pressure measurements are required over an infinitely wide area. In addition, the source positions should be varied over the same area with the same spacing to acquire total insonification of the defects. In practice, recordings will be taken along a finite aperture on a limited number of locations and with a discrete set of source locations. For our ultrasonic measurements, linear arrays with typically 64 elements are used. The one-dimensional Rayleigh integral of equation 2.2.6 can be written as a sum over all the available receiver elements (figure 2.8). With z0= 0 this yields:

P (xm, zA; ω) ≈ − r −jωc 2π N X n=1 P (xn, 0; ω)p zA ∆rm,n∆rm,n ejωc∆rm,n_∆x, _(2.6.11) with P (xm, zA; ω) being the Fourier transformed pressure at an image point A located at

(xm, zA), P (xn, 0; ω) the pressure recorded at element number n, N the total amount of

recorder elements, ∆x the heart-to-heart distance of the recorder elements and ∆rm,n =

p

(xn− xm)2+ zA2, with xn= (n−1)∆x (figure 2.8). Note that in principle the summation

in equation (2.6.11) should be infinite. Through the limited length of the receiver array, artifacts may arise as discussed earlier. We now write equation (2.6.11) as,

P (xm, zA; ω) = N

X

n=1

Qm,nP (xn, 0; ω), (2.6.12)

where we have defined Qm,nas,

(35)

Equation (2.6.12) can be formulated as a vector product: P (xm, zA; ω) = ~QTm(zA; 0; ω) ~P (0; ω), (2.6.14) or similarly, P (xm, zA; ω) = Qm,1 Qm,n . . . Qm,N 0 B B B B P (x1, 0; ω) P (xn, 0; ω) .. . P (xN, 0; ω) 1 C C C C A . (2.6.15)

Equation (2.6.14) formulates the discrete Rayleigh II integral for back propagation in vector notation for a single point that was insonified by a single source. The notation can be extended to a matrix vector operation so that it accounts for more points at depth level zA. The result yields,

~ P (zA; ω) = Q(zA; 0; ω) ~P (0; ω), (2.6.16) or similarly, 0 B B B B P (x1, zA; ω) P (xm, zA; ω) .. . P (xM, zA; ω) 1 C C C C A = 0 B B B B Q1,1 Q1,n . . . Q1,N Qm,1 Qm,n . . . .. . ... . .. QM,1 QM,N 1 C C C C A 0 B B B B P (x1, 0; ω) P (xn, 0; ω) .. . P (xN, 0; ω) 1 C C C C A , (2.6.17)

where m = 1, 2, ...M with M the total number of image points and with Q(zA; 0; ω) the so

called inverse propagation operator. It can be demonstrated for homogeneous media that Q(zA; 0; ω) takes the form of a Toeplitz structure, meaning that elements along diagonal

directions are identical.

The matrix Q(zA; 0; ω) removes the propagation effects from each of the M points in

the image space to all the N receivers. In fact, for a single point located in (xm, zA), all

the receiver positions are redatumed to a virtual point located in (xm, zA). Later, this

virtual point will be used together with the argument of reciprocity, to remove also the propagation effects from all the sources to this virtual receiver point.

Equation (2.6.16) can be used for the back propagation of a wave field that was caused by scatterers that were insonified by a single source element. Like the extended formulation from one image point to several image points, we can extend equation (2.6.16) for more source elements, with a maximum of N elements. Then, the vectors ~P (0; ω) and ~P (zA; ω)

become matrices, which yields:

(36)

2.6 2D imaging implementation using the matrix notation 31 and P(0; 0; ω) = 0 B B B B P (x1, 0; x1, 0; ω) P (x1, 0; xn, 0; ω) . . . P (x1, 0; xN, 0; ω) P (xn, 0; x1, 0; ω) P (xn, 0; xn, 0; ω) . . . P (xn, 0; xN, 0; ω) .. . ... . .. P (xN, 0; x1, 0; ω) P (xN, 0; xN, 0; ω) 1 C C C C A . (2.6.20)

Matrix P(0; 0; ω) contains the Fourier transformed data of all source-receiver combinations for one temporal frequency. Each column refers to a source position and each row refers to a receiver position. Matrix P(zA; 0; ω) contains the pressure recordings of upward traveling

waves at virtual receivers located in (xm, zA) for one frequency ω, i.e. this is one complex

number.

When all the recorded wave fields are back propagated to all the M points, the imag-ing condition must be applied for all the N sources to find the imagimag-ing amplitudes. In fact, we must remove the propagation effects from the sources to the imaging points located on depth level zA. We can do this efficiently with the argument of reciprocity. This means

that we may switch receiver locations with source locations. Hence, all the virtual receivers located in the points (xA, zm) become sources and all the sources located at the surface

be-come receivers. To remove the propagation effects from the new virtual source locations to the new receiver locations, we can use the back propagation matrix Q(zA; 0; ω) again. We

must transpose this matrix since we have switched the sources and receivers. Multiplying both sides of equation (2.6.18) with QT(0; zA; ω) yields,

P(zA; 0; ω)QT(0; zA; ω) = Q(zA; 0; ω)P(0; 0; ω)QT(0; zA; ω). (2.6.21)

The result of equation (2.6.21) can be interpreted as pressure recordings of upward traveling waves by virtual receivers located at depth level zA, caused by virtual sources also located

at depth level zA, hence

P↑(zA; zA; ω) = Q(zA; 0; ω)P(0; 0; ω)QT(0; zA; ω). (2.6.22)

Similar to equation (2.4.9), we can define the upward traveling wave field as the downward traveling wave field multiplied by a matrix that contains scattering coefficient,

P↑(zA; zA; ω) ≡ R(zA; ω)P↓(zA; zA; ω). (2.6.23)

The physical interpretation R(zA; ω) is difficult. The diagonal elements are the scattering

strengths of defect scatterers along depth zA. Interpretation of the diagonal elements is

meaningful only for plane reflection interfaces. In seismic exploration, this is often a real-istic condition, but in weld inspection this condition is only met for planar defects that are large compared to the wavelength.

Because the virtual sources and receivers are located at the same depth level zAin equation

2.6.23, we can find the imaging amplitudes for image points located at this level by apply-ing the imagapply-ing condition at t = 0. This yields a summation over all relevant frequency components. For the imaging amplitudes ~I(zA) we can define,

(37)

where ωk are the discrete frequency components and K is the total number of frequency

components. Note that we have selected the elements on the diagonal of the matrix. This diagonal represents the pressure recordings from virtual receivers caused by virtual sources located at coinciding locations. In case scatterers are present at depth level zA, the

ele-ments on the diagonal of R(zA; ωk) will have a non-zero value. From equation (2.6.24) we

can see that the image amplitude is also determined by the downward traveling wave field P↓_(z

A; zA; ωk). Usually, this wave field is caused by a band limited source. The signature

of this source is recognized in the image.

For the implementation as described above, we assume that the data of all source-receiver combinations are available and the matrix P(0; 0; ω) is completely filled. Due to practical limitations, it is possible that only the data can be obtained for the source-receiver combi-nations that have the same position. This is called zero offset data and in this case only the diagonal of matrix P(0; 0; ω) is available. For 3D imaging, we will use linear arrays to measure the full data matrix in the direction parallel to the array, this is called the in-line direction. This can be done for many positions along the surface, so that zero offset data is obtained in the direction perpendicular the the array, the cross-line direction. In the next section, it will be discussed how zero-offset data can be processed efficiently.

2.7 2D Zero offset imaging

Here, we will discuss the case when recordings are taken with coinciding transmitter and receiver positions. Such a data set is usually referred to as zero offset, meaning that only the diagonal of the data matrix P(0; 0; ω) is filled with data. Note that for primary reflec-tions in zero offset data the path from source to scatterer is the same as the path from scatterer to receiver. The data matrix contains the recorded data of N experiments. If only the diagonal of this matrix is filled, we can interpret the diagonal as a vector that contains the recordings on all the N receiver positions as though they were recorded from a single experiment, whereby the scatterers in the subsurface are considered point sources that ignite at once, see figure 2.9. To account for the travel path from the surface to the scatters and back in reality, the waves caused by buried sources must travel with half the medium sound velocity. This interpretation is the well known exploding-reflector model, as used by Claerbout [1985], Loewenthal et al. [1976], Schnneider [1978] and MÆland [1988]. Spherical spreading should be corrected before the exploding reflector model can be ap-plied.

We can write the vector with zero offset recordings as,

~

Pzo(0; ω) = diag[P(0; 0; ω)], (2.7.25)

with ~Pzo(0; ω) the vector that contains the frequency components of the Fourier

(38)

2.7 2D Zero offset imaging 33 z x n=1 n=N m=M m=1 ∆rm,n (xm, zA) zA (xn,0) ∆x

Figure 2.9: Zero offset imaging can be done with data of coinciding transmitter and receiver location. Zero offset imaging can be interpreted by back propagating the wave field of an exploding source located in A with half the medium velocity.

for back propagation (equation (2.2.6)), using half the medium sound velocity, hence

P (xm, zA; ω) = − s −j2ωc 2π ∞ Z −∞ Pzo(0; ω)√zA ∆r∆re j2ω c∆rdx. (2.7.26)

Since it is assumed that the exploding reflectors are located in the subsurface, a second back propagation step to compensate for the propagation effects from the real sources (located at the surface) to the scatterers, can be omitted.

Similar to equation (2.6.12), equation (2.7.26) can be written in a discrete notation for a limited number of sources and receivers,

Pzo(xm, zA; ω) = N X n=1 ˆ Qm,nPzo(xn, 0; ω), (2.7.27)

with ˆQm,n defined by equation (2.6.13) with half the medium sound velocity. Following

the same steps of the previous section, we may write for more points on depth level zA in

matrix-vector notation, ~

Pzo(zA; ω) = ˆQ(zA; 0; ω) ~Pzo(0; ω), (2.7.28)

where ~Pzo(zA; ω) contains the Fourier transformed pressures at depth level zA. The imaging

condition must be applied to ~Pzo(zA; ω), since the reflectors ’exploded’ at t = 0. We may

(39)

Equation (2.7.28) implies only a matrix vector multiplication. This operation is more efficient to evaluate with a computer then equation (2.6.24) for all source-receiver combi-nations, derived in the previous section. Therefore, it will require less computer operations and hence, less calculation time. Note, however, that zero offset measurements do not carry the same information as full offset measurements.

In the next section, we will discuss 3D imaging by using both 2D imaging approaches in combination. Given the nature of the experiments, we will see that the first 2D imag-ing step will involve all source receiver combinations of equation (2.6.24), and the second 2D imaging step will involve zero offset imaging of equation (2.7.28) in the perpendicular direction.

2.8 3D imaging with linear arrays, the two pass method

In theory, it is possible to construct a 3D image from ultrasonic data with the use of the 2D Rayleigh integral for back propagation given by equation (2.2.4). To evaluate the integral, recordings are necessary over a surface area S. For a good reconstruction of a defect, also sources are required over the same surface area for insonification in all direction. With a linear array, all combinations of transmitters and receivers can be measured over the line where the array is located. To obtain all the combinations of transmitters and receivers over a surface area, a 2D matrix array is required. Ultrasonic matrix arrays have already been studied, manufactured and used in the medical application field (Whittingham [1999], Fenster [2001]) and for NDI (Jesse and Smith [2002], Mahaut et al. [2004]).

Although the technology is available, acquiring the data can be quite time consuming and requires a huge data storage capacity. In addition, processing the data into an image is also time consuming and, therefore, not yet practical in real time. The same issues can be found in seismic exploration, where a full areal source and receiver grid is not feasible in practice. In this application field, these issues have been studied and published. An efficient approach to 3D imaging was developed, that consists of two cascading 2D imaging processes in orthogonal lateral directions whereby a 3D scheme is obtained (Jakubowicz and Levin [1983] and Gibson et al. [1983]). This 3D scheme, called the two pass approach, is exact for homogeneous media as proved by Jakubowicz and Levin [1983]. In this section, we will describe the two pass approach that is used to obtain 3D images with ultrasonic linear array data.

(40)

2.8 3D imaging with linear arrays, the two pass method 35

in-line direction (say parallel to the x-coordinate) must be constructed for each cross-line position (say parallel to the y-coordinate). With equation (2.2.6) we can back propagate the measured wave field of a single source S located at (xS, y, 0) to a depth level zA for

each cross-line location y, this yields P′(xS, y, 0; xA, y, zA′; ω) = − r −jω 2πc ∞ Z −∞ P (xS, y, 0; xR, y, 0; ω) z ′ A √ ∆rR∆rR ejωc∆rR_dx R, (2.8.30) where P (xS, y, 0; xR, y, 0; ω) are the Fourier transformed pressure recordings at z = 0,

P′_(x

S, y, 0; xA, y, z′A; ω) are pseudo-pressure recordings of virtual receivers located at zA′,

∆rR =

p

(xR− xA)2+ z′2A where xR is the x-coordinate of the receivers and where z′A is

the (generally incorrect) position to which the scatterers have been imaged. All source locations can be taken into account by applying equation (2.2.6) again. This yields,

P′(xA, y, zA′; xA, y, zA′; ω) =− jω 2πc ∞ Z −∞ ∞ Z −∞ P (xS, y, 0; xR, y, 0; ω) z ′2 Aej ω c(∆rR+∆rS) √ ∆rR∆rR √ ∆rS∆rS dxRdxS, (2.8.31) where ∆rS= p

(xS− xA)2+ zA′2and where P′(xS, y, zA′; xA, y, zA′; ω) represents the Fourier

transformed pseudo-pressure recordings of upward traveling waves measured by virtual re-ceivers located at z′

A, caused by scattered downward traveling waves generated by virtual

sources also located at z′

A. In order to find the imaging amplitudes I(xA, y, z′A) of the 2D

images for all cross-line positions, we must apply the imaging condition at t = 0. This yields integration of equation (2.8.31) over all frequencies:

I(xA, y, z′A) ≡ p′(xA, y, z′A; xA, y, z′A; 0) = ∞

Z

−∞

P′(xA, y, zA′; xR, y, zA′; ω)dω. (2.8.32)

In the second pass of the process, we convert the 2D images obtained in the first pass with equation (2.8.32), to pseudo-pressure recordings. Because only zero offset recordings are available, we will use the zero offset imaging approach as was described in section 2.7 for the second pass. We assume that the medium is homogeneous so that we can convert the space dimension of zA′ to time with the sound velocity, we define

p′′(xA, y, 0; t) ≡ p′(xA, y, zA′; xA, y, zA′; 0), (2.8.33) with t ≡2z ′ A c . (2.8.34) We can consider p′′_(x

A, y, 0; t) as a zero offset data set. The factor 2 in equation (2.8.34)

(41)

- ? y t x cross-line direction in-line direction s s s (x0, y0, T0) (x, y, t) (x0, y, t0) K i ?_t

-Figure 2.10: The amplitude at an arbitrary point (x, y, t) can be mapped to the apex (x0, y0, T0) either directly or in a two-step process.

where P′′_(x

A, yA, zA; ω) is the Fourier transformed pressure in point A = (xA, yA, zA),

P′′_(x

A, y, 0; ω) is the Fourier transformed zero offset data, and r =

p

(y − yA)2+ zA. The

image amplitudes for the final 3D image can be obtained by applying the imaging condition at t = 0 to equation (2.8.35), hence I(xA, yA, zA) ≡ p′′(xA, yA, zA; 0) = ∞ Z −∞ P′′(xA, yA, zA; ω)dω. (2.8.36)

An intuitive proof of the two pass process for zero offset data was given by Gibson et al. [1983], where 3D imaging is regarded as a summation over trace (A-scan) amplitudes along the surface of hyperboloids in the space time domain, corresponding to diffraction paths and placing each sum at the apex or minimum-time position of the associated hy-perboloid. For a medium of constant velocity c, the surface of the hyperboloids are given by, t2= T02+4(x − x0 )2 c2 + 4(y − y0)2 c2 , (2.8.37)

where x and y are orthogonal coordinates on the recording surface, and t is the two-way time at position (x, y) on the hyperboloid. The hyperboloid’s apex is at position (x0, y0)

and time T0, as illustrated in figure 2.10. Summation to obtain the output amplitude at a

fixed apex position (x0, y0, T0) can be viewed as a mapping of reflection amplitudes from the

(42)

2.8 3D imaging with linear arrays, the two pass method 37

same velocity c. The total summation can be done in two steps as illustrated in figure 2.10. First a temporary sum over a hyperbolic path in the x-direction is done for each value of y. The temporary sums will be placed at a local apex position at (x0, y, t0). Then in the second

step, these temporary sums are themselves summed along the hyperbola in the x-direction. The resultant amplitude at the apex point (x0, y0, T0) is the image amplitude at that point.

(43)

(44)

3

Application of imaging techniques

for weld inspection

In this chapter, an overview of the most common weld properties is given, together with the practical consequences for imaging1. It will be illustrated that most defects and their characteristics are related to the weld properties and to the welding process. Not all defects will cause a failure and, therefore, acceptance criteria have been developed. In this chapter, the requirements for the quality of images of defects will be derived based on these criteria. Finally, given the practical restrictions and requirements, some aspects of the resolution will be discussed.

3.1 Properties of girth welds in pipelines and defects

Pipelines can be found all over the globe, both on-shore and off-shore, for the transportation of various products like oil and gas. The diameter of pipelines vary from 10” to 50” and the wall thickness from 5 mm to 35 mm. The pipeline and the weld properties depend on the product that is transported by the pipeline and of the location of the pipeline. The type of defects that can occur in newly constructed girth welds are, therefore, a result of the properties related to the pipe and the weld. Properties that influence the defect characteristics are, for example:

1_{The application of imaging based on inverse wave field extrapolation to girth weld inspection}

Imaging of defects in girth welds using inverse wave field extrapolation of ultrasonic data

Inverse Wave Field Extrapolation of

Ultrasonic Data

Imaging of Defects in Girth Welds using

Inverse Wave Field Extrapolation of

Ultrasonic Data

PROEFSCHRIFT

ter verkrijging van de graad van doctor

aan de Technische Universiteit Delft,

op gezag van de Rector Magnificus prof. dr. ir. J.T. Fokkema,

voorzitter van het College voor Promoties,

in het openbaar te verdedigen

op dinsdag 6 november 2007 om 10:00 uur

door

Niels P ¨

ORTZGEN

Contents

1

Introduction

1.1

Non-destructive inspection techniques

1.2

Conventional defect sizing

1.3

Automated girth weld inspection based on zone

discrimi-nation

flat bottom holes

notches

1.4

Ultrasonic arrays in NDI

1.5

The use of arrays in seismic exploration

1.6

Outline and objectives of this thesis

2

Wave theory and imaging

2.1

Introduction of the imaging philosophy

2.2

The Rayleigh II integral

2.3

Wave field extrapolation in elastic media

2.4

The imaging condition

2.5

An example of imaging with simulated acoustic data

2.6

2D imaging implementation using the matrix notation

2.7

2D Zero offset imaging

2.8

3D imaging with linear arrays, the two pass method

3

Application of imaging techniques

for weld inspection

3.1

Properties of girth welds in pipelines and defects