Elastic stiffness characterization using three-dimensional full-field deformation obtained with optical coherence tomography and digital volume correlation

Jiawei Fu; Fabrice Pierron; Pablo D. Ruiz

doi:10.1117/1.JBO.18.12.121512

17 December 2013 Elastic stiffness characterization using three-dimensional full-field deformation obtained with optical coherence tomography and digital volume correlation

Jiawei Fu, Fabrice Pierron, Pablo D. Ruiz

Author Affiliations +

Journal of Biomedical Optics, Vol. 18, Issue 12, 121512 (December 2013). https://doi.org/10.1117/1.JBO.18.12.121512

Abstract

This paper presents a methodology for stiffness identification from depth-resolved three-dimensional (3-D) full-field deformation fields. These were obtained by performing digital volume correlation on optical coherence tomography volume reconstructions of silicone rubber phantoms. The effect of noise and reconstruction uncertainties on the performance of the correlation algorithm was first evaluated through stationary and rigid body translation tests to give an indication of the minimum strain that can be reliably measured. The phantoms were then tested under tension, and the 3-D deformation fields were used to identify the elastic constitutive parameters using a 3-D manually defined virtual fields method. The identification results for the cases of uniform and heterogeneous strain fields were compared with those calculated analytically through the constant uniaxial stress assumption, showing good agreement.

1. Introduction

The measurement of material deformation is necessary to assess the material mechanical properties. A wide variety of techniques have been developed to measure the deformation of materials under load, ranging from point-wise sensors, such as the resistive strain gauge¹ and optical fiber Bragg gratings,² to two-dimensional (2-D) full-field measurements, including digital image correlation (DIC),³ the grid method,⁴ speckle interferometry,⁵ and Moiré interferometry.⁶ For homogeneous and isotropic materials, these techniques usually provide enough information to investigate their mechanical behavior. However, for materials with more complex structures, such as biological tissues and composites, surface measurements are much less adequate to address their complete mechanical behavior since the deformation may vary significantly between the bulk and the surface. In this case, a depth-resolved three-dimensional (3-D) full-field measurement of the deformation would be highly desirable.

Thanks to the development of the various tomographic techniques, digital volume correlation (DVC) has become a popular measurement technique for depth-resolved 3-D deformation fields. It is effectively a 3-D extension of DIC, widely applied to determine the surface deformation fields. DIC was first introduced in the early 1980s³ and has been widely applied in many disciplines, such as mechanical engineering, material science, medical science, etc.⁷^–¹¹ It basically determines the deformation field by tracking a speckle pattern imprinted on the surface of the sample between a reference and a deformed state. Based on the same principle, DVC was developed to measure the internal 3-D deformation behavior of materials by tracking internal features that resemble 3-D speckle patterns contained in the reconstructed volumes.¹² It requires sufficient speckle contrast in the reconstructed volumes to ensure the correlation algorithm runs successfully. DVC has made its way into clinical and industrial applications,¹³^–²⁰ where it is mainly applied on x-ray computed tomography data of materials such as composites, metals, foams, and trabecular bones.¹³^–¹⁶ In all these cases, the pattern contrast of the reconstructed volume is provided by differences in x-ray absorption of the different constituents or phases of the material. However, for soft biological tissues dominated by collagen, such as cornea, arteries, or skin, optical coherence tomography (OCT) is a more suitable technique to reconstruct the material microstructure.

OCT is a noninvasive imaging technique that can acquire micrometer resolution of cross-sectional images from within semitransparent, light scattering materials. It is an interferometric technique that uses a broadband light source with short coherence length to provide depth-resolved information of the object microstructure. The contrast of the images encodes refractive index changes in the material (between fibers and a matrix into which they are embedded, for example). Indeed, this rapidly developing imaging technique has already been applied in ophthalmology, cardiology, gastroenterology, and dermatology,²¹^–²⁴ and many commercially available OCT systems have been developed for clinical and research purposes.²⁵ With the aid of OCT, the volumetric data that represent the whole microstructure of the scanned sample can be reconstructed. The reconstructed data can then be utilized for structure analysis or deformation analysis using DVC. To the best of our knowledge, this is the first time that measurements of depth-resolved deformations have been undertaken combining OCT and DVC.

Given a constitutive model, constitutive parameters are constants that describe the mechanical behavior of a material. Thanks to the development of full-field deformation measurement techniques, several methods have been developed to identify the material constitutive parameters. Finite element model updating (FEMU) is a widely applied approach among these methods. It compares the experimental measurements with their numerical counterparts obtained from an FE model and builds up a cost function using the difference between numerical and experimental values in terms of displacement or strain. By minimizing the cost function with respect to the sought constitutive parameters, the solution of the problem can be reached iteratively. FEMU has already been shown to be feasible to identify constitutive parameters in many cases, such as elasticity, hyperelasticity, etc.²⁶^–²⁸ This approach, however, exhibits some shortcomings. The initial values to start the iteration procedure, the numerical model, and the quality of the cost function, for instance, can all affect the convergence time and the quality of the results. The virtual fields method (VFM) is an alternative option for solving the identification problem,²⁹ which is based on the principle of virtual work and retrieves the constitutive parameters by utilizing the full-field deformation measurements. This method is more effective than FEMU in terms of computation time since for the latter, an FE model needs to be created and updated iteratively. So far, the VFM has successfully been applied to the identification of constitutive parameters for linear elastic materials such as composites,³⁰^–³³ elasto-plasticity for metals³⁴^–³⁶ as well as hyperelasticity for soft and biological materials such as artery walls.³⁷^,³⁸

The aim of this study is to develop an effective methodology to investigate the internal deformation behavior of silicone rubber phantoms under load and to identify their elastic stiffness parameters. This methodology will eventually be applied to study the mechanical properties of biological tissues, such as the vertebrate eye cornea.³⁹ First, a brief description of the material, the experimental setup, the experimental methods as well as the identification theory is presented. Then, the performance of DVC coupled with OCT is evaluated through stationary and rigid body translation tests. After that, DVC results of the tensile tests are examined and critically discussed. Finally, the depth-resolved full-field deformation measurement results are used to identify the material elastic stiffness parameters using VFM and the identification results are also discussed.

2. Materials and Methods

2.1.

Materials

In the present study, two rectangular flat phantom strips were fabricated using silicone gel (MM240-TV), one of them with a notch. The material comprises two parts, a rubber base and a hardener. They were mixed to a ratio of $10 ∶ 1$ . The nominal Young’s modulus of the silicone rubber is 1.88 MPa, which can change with the proportion of hardener. Since the silicone gel does not present suitable speckle contrast for the application of DVC, copper particles (with an average particle size of $\sim 1 μ m$ ) were seeded into the silicone gel mixture to provide the speckle contrast. The mixture was then put into molds and left to cure at room temperature for $\sim 24 h$ . The flat strips were cut from a larger piece using a scalpel to $60 \times 1.4 \times 10 {mm}^{3}$ ( $length \times thickness \times width$ ). Figure 1 shows the rectangular and the notched phantom strips and their loading and observation configurations.

Fig. 1

Schematic of the fabricated silicone rubber phantom strips and the loading configuration.

2.2.

Experimental Setup and Image Acquisition

The experimental setup consists of a tensile test fixture and a swept-source OCT system (SS-OCT, Thorlabs OCS1300SS, Ely, United Kingdom, lateral resolution 25 μm, depth resolution 12 μm in air). For the test, the phantom strip was mounted with one end fixed to the fixture and the other end loaded by a dead weight. At first, a dead weight of 50 g was applied as a preload to insure the phantom strip was taut and vertical. This is a necessary step since the silicone rubber is rather compliant. The preload state was taken as the reference state. A 10 g dead weight was then added to the preload and considered as the deformed state, referred to as load step 1 hereafter.

For both reference and deformed states, a 3-D volume image sequence of the specimen was acquired using the SS-OCT system. It uses a rapidly tuned narrowband light source with central wavelength of 1325 nm and spectral bandwidth of 100 nm and records the information with a single photodetector. The frequency of the interference signal is proportional to the optical path difference between the sample and reference arms of an interferometer. Depth profiles (A-scans) of the sample are obtained by evaluating the Fourier transform of the signal for each wavelength scan. The 5- to 6-mm coherence length of the laser enables $\sim 3 mm$ depth measurement. Adjacent A-scans are then synthesized to create an image in the $x y$ plane. Multiple adjacent 2-D images in the $z$ direction then form the reconstructed 3-D volume. In the present work, a region with dimensions of $11 \times 3 \times 11 {mm}^{3}$ was scanned, corresponding to $x$ , $y$ , and $z$ directions, respectively. A 3-D data volume of $1024 \times 512 \times 1024 voxels$ was obtained. The acquisition time for each 3-D volume is $\sim 5 \min$ . It should be pointed out that each time before acquiring the volume images, the specimen was left for 10 min under load to accommodate significant short-term viscoelastic creep. Along the lateral scanning directions $x$ and $z$ , the image voxel size was determined by dividing the 11 mm dimension by the number of corresponding 1024 voxels, which is equal to 10.7 μm. Regarding the through-thickness $y$ direction, which corresponds to the optical path, the voxel size depends on the refractive index of the medium. For the silicone rubber phantom, the voxel size inside the medium along the $y$ direction is 4.1 μm, which is determined by dividing its thickness, 1.4 mm, by the number of corresponding voxels, here 345. The reconstructed volume and a typical central transverse $z$ -slice ( $1024 \times 512 voxels$ ) of the specimen are given in Figs. 2(a) and 2(b), respectively. Due to light saturation, the voxels at the top surface exhibit very high intensity values, as can be seen from the white line at the top of the phantom in Fig. 2(b). It should be pointed out that all the regions outside the phantom strip and the saturated voxels at the top surface were masked out in Fig. 2(a) for a better visualization of the 3-D reconstructed volume. The reconstructed volumes for the reference and deformed states were recorded and DVC was then used to compute the 3-D displacement and strain fields.

Fig. 2

(a) Three-dimensional (3-D) view of the reconstructed volume. (b) Two-dimensional view of the central transverse $z$ slice of the rectangular phantom strip generated through the swept-source optical coherence tomography system.

2.3.

Digital Volume Correlation

DVC is the 3-D extension of the widely applied DIC used to measure surface deformations. During the DVC procedure, the correlated volumes are first divided into subvolumes. The displacement vector of each subvolume is determined by tracking and matching the voxels of the subvolumes in the reference and deformed states. This is performed by maximizing the correlation coefficient, which measures the degree of similarity of the gray-level distributions in the subvolumes in the reference and deformed states. The best prediction of the displacement leads to the highest degree of similarity of the gray-level distributions and thus the maximal correlation coefficient (a correlation coefficient value close to 1 indicates a perfect match). In the present study, the displacement fields were calculated using the DaVis® (LaVision) software package based on a fast Fourier transform (FFT) algorithm. It evaluates a normalized cross-correlation coefficient (NCC) defined as

Eq. (1)

C = \frac{\sum \sum \sum F (x, y, z) G (\hat{x}, \hat{y}, \hat{z})}{\sqrt{\sum \sum \sum {[F (x, y, z)]}^{2} \sum \sum \sum {[G (\hat{x}, \hat{y}, \hat{z})]}^{2}}},

where

F (x, y, z)

represents the gray level at a voxel

(x, y, z)

in the subvolume of the reference state, while

G (\hat{x}, \hat{y}, \hat{z})

represents the gray level at a point

(\hat{x}, \hat{y}, \hat{z})

in the subvolume of the deformed state. The coordinates

(x, y, z)

and

(\hat{x}, \hat{y}, \hat{z})

stand for the same material point in the reference and deformed states, respectively, and are related by the 3-D affine transformation in the form of rigid body motion combined with displacement gradients.

Eq. (2)

\hat{x} = x + U_{x} + \frac{\partial U_{x}}{\partial x} Δ x + \frac{\partial U_{x}}{\partial y} Δ y + \frac{\partial U_{x}}{\partial z} Δ z \hat{y} = y + U_{y} + \frac{\partial U_{y}}{\partial x} Δ x + \frac{\partial U_{y}}{\partial y} Δ y + \frac{\partial U_{y}}{\partial z} Δ z \hat{z} = z + U_{z} + \frac{\partial U_{z}}{\partial x} Δ x + \frac{\partial U_{z}}{\partial y} Δ y + \frac{\partial U_{z}}{\partial z} Δ z,

where

U_{x}

,

U_{y}

, and

U_{z}

are the rigid body displacement components of the subset center in the

x

,

y

, and

z

directions, respectively.

Δ x

,

Δ y

, and

Δ z

represent the distance between the subset center and the point

(x, y, z)

. A double-pass approach was used whereby large subvolumes were initially used to capture large displacements. Subsequent to this, these initially calculated displacements were used to displace smaller subvolumes and, thus, ensure the pattern was followed and signal-to-noise ratio increased. Gaussian curve-fitting of the correlation function peak was used to detect the position of the displacement with subvoxel resolution. The strains were then determined from the centered finite difference of the calculated displacement fields, without any additional smoothing.

2.4.

Evaluation of Measurement Performance

The DVC algorithm requires sufficient contrast in a subvolume in order to determine a displacement vector. The size of the subvolumes influences the value of the correlation coefficient and, thus, the displacement and strain uncertainties. In the present study, four different subvolume sizes were selected and compared to determine an optimal size, considering the displacement and strain spatial variation as well as the spatial resolution. A double-pass approach used initial subvolume sizes of $24^{3}$ , $48^{3}$ , $72^{3}$ , and $96^{3}$ in the first pass, followed by 12³, $24^{3}$ , $36^{3}$ , and $48^{3}$ in the second pass, and each has 50% overlap with the six adjacent neighbors. Thus, the distance between each subvolume center and its immediate neighbors is 6, 12, 18, and 24 voxels, respectively.

Prior to the tensile tests, it is necessary to evaluate the errors caused by all sources of noise and reconstruction uncertainties. This can be done by performing DVC on two subsequent reconstructed volumes of the stationary phantom strip. Since the stationary specimen was not subjected to any applied force, the correlation results should show zero displacement and strain fields over the whole field of view. This is not the case in practice due to the influence of electronic noise in the SS-OCT system, environmental vibration, the volume reconstruction algorithm, etc. Therefore, any nonzero results should be attributed to the contribution of noise and other reconstruction uncertainties. The standard deviations of the spurious strains were calculated to evaluate the resolution (uncertainty) of the strain measurement.

Then, two reconstructed volumes were recorded after introducing a rigid body translation of 40 μm ( $\sim 10 voxels$ ) in the $y$ direction between the two volumes. This not only tests the effect of all sources of noise, but also tests the performance of DVC subvoxel interpolation (tri-linear in the present study) and accuracy of the correlation algorithm in determining the displacement fields for a translated specimen. The same procedure of data processing as for the stationary test was applied to the rigid body translation test and the strain measurement resolution was computed. The above tests give an overall idea of the resolution of the whole setup so that the significance of the tensile test results can be better analyzed.

2.5.

Virtual Fields Method

VFM, which is based on the principle of virtual work, was adopted to identify the constitutive parameters of the silicone rubber phantom. This equation, provided below, is the integral form of the equilibrium equation for the standard deformable continuous solid model.

Eq. (3)

- \int_{V} σ : ε^{*} d V + \int_{\partial V} \bar{T} \cdot U^{*} d S + \int_{V} b \cdot U^{*} d V = \int_{V} ρ a \cdot U^{*} d V .

In this equilibrium equation, $σ$ is the actual stress tensor, $ε^{*}$ is the virtual strain tensor, $U^{*}$ is the virtual displacement vector, $\bar{T}$ is the traction vector acting on the boundary $\partial V$ of the solid volume $V$ , $b$ is the volume force vector, $a$ is the acceleration vector, and $ρ$ is the mass per unit volume; “·” indicates vector dot product and “:” is the contracted product for second-order tensors or the matrix dot product. In this study, the phantom strip was loaded statically and the body forces can be neglected. Therefore, their contribution to the virtual work can be canceled out. Thus, the equilibrium equation becomes

Eq. (4)

- \int_{V} σ : ε^{*} d V + \int_{\partial V} \bar{T} \cdot U^{*} d S = 0 .

This equilibrium equation is valid for any continuous and differentiable virtual displacement field. The actual strain field and the load information are provided by the experiment. The stress components can be expressed through the material constitutive parameters and the strain components through an appropriate constitutive equation. Here, the stress and strain components are written in columns as

Eq. (5)

σ : {\begin{matrix} σ_{x x} & σ_{y y} & σ_{z z} & σ_{y z} & σ_{x z} & σ_{x y} \end{matrix}}^{T} ε : {\begin{matrix} ε_{x x} & ε_{y y} & ε_{z z} & 2 ε_{y z} & 2 ε_{x z} & 2 ε_{x y} \end{matrix}}^{T},

where

T

is the transpose operator.

For the silicone rubber phantom, isotropic elasticity was assumed. Thus, the stress-strain relation can be expressed as

Eq. (6)

{\begin{matrix} σ_{x x} \\ σ_{y y} \\ σ_{z z} \\ σ_{y z} \\ σ_{x z} \\ σ_{x y} \end{matrix}} = [\begin{matrix} Q_{x x} & Q_{x y} & Q_{x y} & 0 & 0 & 0 \\ Q_{x y} & Q_{x x} & Q_{x y} & 0 & 0 & 0 \\ Q_{x y} & Q_{x y} & Q_{x x} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{Q_{x x} - Q_{x y}}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{Q_{x x} - Q_{x y}}{2} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{Q_{x x} - Q_{x y}}{2} \end{matrix}] {\begin{matrix} ε_{x x} \\ ε_{y y} \\ ε_{z z} \\ 2 ε_{y z} \\ 2 ε_{x z} \\ 2 ε_{x y} \end{matrix}},

where

Q_{x x}

and

Q_{x y}

are the two stiffness components to identify, relating to the elastic modulus

E

and the Poisson’s ratio

ν

through

Eq. (7)

{\begin{matrix} Q_{x x} = \frac{(1 - ν) E}{(1 + ν) (1 - 2 ν)} \\ Q_{x y} = \frac{ν E}{(1 + ν) (1 - 2 ν)} \end{matrix} .

After introducing Eqs. (5) and (6) into Eq. (4), the equilibrium equation becomes

Eq. (8)

\int_{V} Q_{x x} (ε_{x x} ε_{x x}^{*} + ε_{y y} ε_{y y}^{*} + ε_{z z} ε_{z z}^{*} + 2 ε_{y z} ε_{y z}^{*} + 2 ε_{x z} ε_{x z}^{*} + 2 ε_{x y} ε_{x y}^{*}) d V + \int_{V} Q_{x y} (ε_{x x} ε_{z z}^{*} + ε_{z z} ε_{x x}^{*} + ε_{y y} ε_{z z}^{*} + ε_{z z} ε_{y y}^{*} + ε_{x x} ε_{y y}^{*} + ε_{y y} ε_{x x}^{*} - 2 ε_{y z} ε_{y z}^{*} - 2 ε_{x z} ε_{x z}^{*} - 2 ε_{x y} ε_{x y}^{*}) d V = \int_{\partial V} \bar{T} \cdot U^{*} d S,

where

ε_{i j}^{*}

is the virtual strain component deduced from the virtual displacement. In the present study, the material is macroscopically homogeneous. Therefore,

Q_{x x}

and

Q_{x y}

can be taken out of the integrals and the equilibrium equation [Eq. (8)] becomes

Eq. (9)

Q_{x x} \int_{V} (ε_{x x} ε_{x x}^{*} + ε_{y y} ε_{y y}^{*} + ε_{z z} ε_{z z}^{*} + 2 ε_{y z} ε_{y z}^{*} + 2 ε_{x z} ε_{x z}^{*} + 2 ε_{x y} ε_{x y}^{*}) d V + Q_{x y} \int_{V} (ε_{x x} ε_{z z}^{*} + ε_{z z} ε_{x x}^{*} + ε_{y y} ε_{z z}^{*} + ε_{z z} ε_{y y}^{*} + ε_{x x} ε_{y y}^{*} + ε_{y y} ε_{x x}^{*} - 2 ε_{y z} ε_{y z}^{*} - 2 ε_{x z} ε_{x z}^{*} - 2 ε_{x y} ε_{x y}^{*}) d V = \int_{\partial V} \bar{T} \cdot U^{*} d S .

Any new virtual field in the equilibrium equation leads to an equation involving the stiffness components. Therefore, a proper choice of the virtual displacement fields enables the identification of the unknown stiffness components $Q_{x x}$ and $Q_{x y}$ . There are different methods for choosing virtual fields. In the present work, two manually defined polynomial virtual fields were employed for the two sought constitutive parameters. More details about other choices of virtual fields, such as virtual fields minimizing noise effects, piecewise virtual fields, etc., can be found in Ref. 29. The first virtual displacement field and the corresponding virtual strain field are

Eq. (10)

U_{x}^{* (1)} = x - L; U_{y}^{* (1)} = 0; U_{z}^{* (1)} = 0; ε_{x x}^{* (1)} = 1; ε_{y y}^{* (1)} = 0; ε_{z z}^{* (1)} = 0; ε_{y z}^{* (1)} = 0; ε_{x z}^{* (1)} = 0; ε_{x y}^{* (1)} = 0 .

For the second virtual field,

Eq. (11)

U_{x}^{* (2)} = 0; U_{y}^{* (2)} = x y (x - L); U_{z}^{* (2)} = 0; ε_{x x}^{* (2)} = 0; ε_{y y}^{* (2)} = x (x - L); ε_{z z}^{* (2)} = 0; ε_{y z}^{* (2)} = 0; ε_{x z}^{* (2)} = 0; ε_{x y}^{* (2)} = x y - \frac{1}{2} y L,

where

L

is the length (dimension along the

x

direction) of the phantom strips.

After introducing the above two virtual fields [Eqs. (10) and (11)] into the equilibrium equation [Eq. (9)], a linear equation system [Eq. (12)] can be formed to directly determine the sought stiffness components $Q_{x x}$ and $Q_{x y}$ .

Eq. (12)

A Q = B,

A : [\begin{array}{l} \int_{V} ε_{x x} d V & \int_{V} (ε_{y y} + ε_{z z}) d V \\ \int_{V} [x (x - L) ε_{y y} + (2 x y - y L) ε_{x y}] d V & \int_{V} [x (x - L) (ε_{x x} + ε_{z z}) - (2 x y - y L) ε_{x y}] d V \end{array}] Q : {\begin{matrix} Q_{x x} \\ Q_{x y} \end{matrix}} B : {\begin{matrix} - F L \\ 0 \end{matrix}},

where

F

is the tension load at the position

x = 0

,

- L

represents the virtual displacement at the same position calculated from Eq. (10) for the first virtual field, and

- F L

represents the external virtual work done by the tension load.

It should be pointed out that due to the discrete nature of the deformation measurement, the integrals above such as in Eq. (9) must be approximated by discrete sums. For instance, an integral $\int_{V} ε_{x x} ε_{x x}^{*} d V$ can be approximated by $\sum_{i = 1}^{n} ε_{x x}^{(i)} ε_{x x}^{* (i)} v^{(i)}$ , where $n$ is the number of data points, $ε_{x x}^{(i)}$ is the actual strain $ε_{x x}$ at data point $(i)$ , $ε_{x x}^{* (i)}$ is the virtual strain $ε_{x x}^{*}$ at data point $(i)$ , and $v^{(i)}$ is the volume of each measurement point. Due to the uniaxial tension configuration, the strain field should be dominated by low spatial frequencies. In addition, a relatively small subvolume size was chosen to enable a high strain spatial resolution. Therefore, the error involved due to the discretization is small and considered as acceptable for the current tests. Thus, matrix $A$ in the linear equation system [Eq. (12)] becomes

Eq. (13)

A : [\begin{array}{l} \sum_{i = 1}^{n} ε_{x x}^{(i)} v^{(i)} & \sum_{i = 1}^{n} (ε_{y y}^{(i)} + ε_{z z}^{(i)}) v^{(i)} \\ \sum_{i = 1}^{n} [x^{(i)} (x^{(i)} - L) ε_{y y}^{(i)} & \sum_{i = 1}^{n} [x^{(i)} (x^{(i)} - L) (ε_{x x}^{(i)} + ε_{z z}^{(i)}) \\ + (2 x^{(i)} y^{(i)} - y^{(i)} L) ε_{x y}^{(i)}] v^{(i)} & - (2 x^{(i)} y^{(i)} - y^{(i)} L) ε_{x y}^{(i)}] v^{(i)} \end{array}] .

For these tests, $v^{(i)}$ of each data point is constant because of the constant subvolume size, which is equal to $V / n$ . Therefore, $v^{(i)}$ can be taken out of the sum and a sum in matrix $A$ , for instance, becomes

Eq. (14)

\sum_{i = 1}^{n} ε_{x x}^{(i)} v^{(i)} = \frac{V}{n} \sum_{i = 1}^{n} ε_{x x}^{(i)} = V \bar{ε_{x x}},

where

\bar{ε_{x x}}

denotes the arithmetic spatial average of

ε_{x x}

and

V

is the total volume of the processed data. Finally, the linear equation system [Eq. (12)] to solve can be written as follows:

Eq. (15)

[\begin{array}{l} \bar{ε_{x x}} & \bar{ε_{y y} + ε_{z z}} \\ \bar{x (x - L) ε_{y y}} + \bar{(2 x y - y L) ε_{x y}} & \bar{x (x - L) (ε_{x x} + ε_{z z})} - \bar{(2 x y - y L) ε_{x y}} \end{array}] {\begin{matrix} Q_{x x} \\ Q_{x y} \end{matrix}} = {\begin{matrix} \frac{- F L}{V} \\ 0 \end{matrix}} .

3. Results and Discussions

3.1.

Strain Noise Analysis

A fundamental condition in image correlation techniques is that the changes in the intensity pattern are in one-to-one correspondence with the displacements of the surface. DVC results from OCT volume reconstructions of the reference and deformed states are likely to be affected by a variety of factors such as (1) electronic noise of the detectors; (2) light source stability; (3) reconstruction algorithms; (4) contrast reduction through the sample thickness due to material absorption, scattering, dispersion, defocusing, and spectral roll-off; and (5) strain-induced speckle decorrelation. We studied the combined effect of (1) to (4) by performing a stationary test and a rigid translation of the phantoms. The effect of strain in speckle decorrelation was explored through a numerical simulation described below.

3.1.1.

Stationary test

From the DVC results of the stationary test, the influence of subvolume size was analyzed quantitatively by comparing the standard deviations of the strain components for a central $z$ -slice $(x, y)$ for the four final subvolume sizes. All six strain components were derived from the centered finite differences of the spurious displacements as follows, without any smoothing:

Eq. (16)

ε_{i j} = \frac{1}{2} (U_{i, j} + U_{j, i}),

where the commas stand for the partial derivatives. To study the effect of subvolume size, the standard deviations over the whole field of view are compared in Fig. 3(a) for

ε_{x x}

,

ε_{x y}

, and

ε_{y y}

only, for the sake of legibility. The standard deviations of the strain components drop from

\sim 3 \times 10^{- 3}

to

\sim 5 \times 10^{- 4}

when increasing the subvolume size from

12^{3}

to

24^{3}

. When further increasing the subvolume size to

36^{3}

and

48^{3}

[see Fig. 3(b)], the strain fluctuations are further reduced to

\sim 2 \times 10^{- 4}

. This is so because the larger number of features contained in the

48^{3}

subvolume assists the convergence of the volume correlation algorithm and enables more accurate tracking of the subvolume deformation. A smaller number of features in the

12^{3}

subvolume leads to bigger tracking errors. In the case of thin specimens such as corneas,³⁹ or the phantoms studied in this work,

36^{3}

subvolumes were found to be a good compromise in terms of strain and spatial resolution. Depending on the OCT spatial sampling rate, the speckle field may be undersampled, leading to aliasing and interpolation errors.⁴⁰ This issue will be discussed later on in the paper.

Fig. 3

(a) Comparison of the strain standard deviations for different subvolume sizes, calculated for a stationary rectangular phantom strip. (b) 3-D views of the subvolumes with different sizes.

Based on the volume strain fields of the stationary test, the standard deviation of each strain component was calculated for each $z$ -slice and the results are plotted in Fig. 4(a). It can be observed that the standard deviations of all the strain components generally remain stable along the different $z$ -slices, with a slight increase toward the ends. This is expected as the DVC results are usually noisier near edges due to the lack of data. Although with fluctuation, all the standard deviations are generally between $1.5 \times 10^{- 4}$ and $2.5 \times 10^{- 4}$ without any smoothing, which is considered as satisfactory compared with the strain levels in the tensile tests later on.

Fig. 4

Strain standard deviations obtained with $36^{3}$ voxels subvolume and 50% overlap for (a) stationary test and (b) rigid body translation test.

3.1.2.

Translation test

Slightly larger values of strain standard deviation, between $2 \times 10^{- 4}$ and $3.5 \times 10^{- 4}$ , were obtained for the rigid body translation test, see Fig. 4(b). This is so as extra sources of error add to those present in the stationary test, mainly subvoxel interpolation error as the DVC algorithm tracks the subvolumes between the reference and displaced states. It was observed that lateral translations (in the $x z$ plane) lead to strain standard deviation values between those obtained for the stationary and the axial translation tests. These levels are still low compared to the strain levels of $\sim 1 %$ in the tensile tests and were thus considered satisfactory.

3.1.3.

Strain-induced speckle decorrelation

Due to the backscatter illumination/observation configuration, the complex 3-D point spread function (PSF) of the OCT system has $\sim 18$ fringes across it along the axial direction [ratio between the depth resolution, 8.3 μm, and the half wavelength of the light source in the medium of refractive index 1.45, i.e., $1.325 μm / (2 \times 1.45) = 0.457 μ m]$ . The magnitude of the measured OCT signal corresponds to the convolution between the 3-D PSF and the scattering particles within the phantom. This magnitude, which determines the brightness of the 3-D speckle grain at any particular position within the sample, does not change with rigid body motion of the sample as relative phase differences between scatterers within the 3-D PSF remain constant. In case of strain, however, there is a limit within which the magnitude of the speckle remains nearly unchanged and beyond which an incremental DVC approach would be required.

In order to estimate the level of OCT speckle decorrelation due to strain, we performed a 2-D (on the $x y$ plane) numerical simulation involving the following steps:

1. Generate a 2-D random distribution of scatterers such that there are many of them ( $\sim 100$ ) inside the PSF of the OCT system.
2. Evaluate the 2-D speckle field due to the spatial distribution of scatterers considering the numerical aperture of the system, the central wavelength and bandwidth of the source, and the refractive index of the medium. This was done using linear systems theory⁴¹^,⁴² by first calculating the transfer function of the OCT system and then evaluating the complex PSF and convolving it with the scatterers. The speckle field was oversampled to recover phase information within the PSF. In this way, the correlation coefficient evaluation is free from undersampling effects. We used images of $1024 \times 1024 pixels$ with a speckle size given by the dimensions of the PSF ( $\sim 1024 / 4 pixels$ in the axial direction $y$ and $\sim 1024 / 2 pixels$ in the lateral direction $x$ ). Using the Rayleigh resolution criteria, this leads to $\sim 8 \times 4 = 32$ speckles in the simulated subset, which compares well with the number of speckles observed on the $x y$ face of the $36^{3}$ subvolume shown in Fig. 3(b).
3. Deform the spatial distribution of scatterers with a horizontal strain $ε_{x x}$ from 0 to 0.5 in steps of 0.001 from 0 to 0.02 and steps of 0.05 thereafter. Poisson’s contraction in the vertical direction was also considered, using $ν = 0.42$ as estimated for our phantoms in Sec. 3.3. For each deformation state, the intensity (magnitude) of the corresponding speckle field was calculated.
4. Finally, the normalized cross-correlation as defined in Eq. (1) was evaluated between the first speckle field for $ε_{x x} = 0$ and all others in the sequence.

Figure 5 shows that NCC drops to $\sim 0.987$ for $ε_{x x} = 1 %$ and to $\sim 0.961$ for $ε_{x x} = 2 %$ . This latter strain corresponds to a total maximum through-thickness relative displacement of the scatterers in the PSF of $\sim 70 nm$ . This is equivalent to $\sim 1 / 6$ th of the fringe period inside the PSF and is inversely proportional to the Poisson’s ratio. NCC drops to $\sim 0.9$ for $ε_{x x}$ $\sim 4 %$ , which corresponds to $ε_{y y}$ $\sim 1.68 %$ using $ν = 0.42$ . This strain level is probably a good estimate to the maximum strain that we can measure with OCT and DVC without using an incremental approach. Above this level, the correlation coefficient is too low to guarantee a good estimate of displacements. When the changes in the magnitude of the interference of light scattered by particles within the PSF are large, the speckle is said to boil, i.e., its structure changes while only the average speckle size is preserved. Even though these results correspond to a 2-D case (B-scan), a 3-D simulation is expected to render results similar to those reported here. Zaitsev et al.⁴³ have found similar results performing a numerical simulation and evaluating the zero mean normalized cross-correlation coefficient as a function of strain. They report that speckle boiling fully decorrelates the speckle for axial strain levels of $\sim 1.5 %$ , which compares well with our figure of $\sim 1.68 %$ mentioned above.

Fig. 5

Normalized cross-correlation of a speckle field in a simulated B-scan as a function of longitudinal strain $ε_{x x}$ .

3.2.

DVC Results for Tensile Tests

3.2.1.

Correlation coefficient maps

The 3-D deformation field was measured under tension after performing DVC on the OCT reconstructed volume data for the rectangular and notched phantom strips. As DVC was performed using the subvolume size of $36^{3}$ voxel and 50% overlap, the reconstructed volume of interest thus contained $55 \times 18 \times 50$ measurement points, corresponding to the dimensions of $10.5 \times 1.4 \times 10 {mm}^{3}$ . The reliability of the deformation measurements was assessed in terms of the 3-D correlation coefficient maps shown in Figs. 6(a) and 6(b) for the rectangular and notched phantom strips, respectively. In order to see the correlation coefficients within the specimens, subvolumes of the whole fields are represented here and are obtained by cutting the volumes in the $x y$ plane at $z$ -slice 25. For $y$ -slice 18 at the top of the rectangular specimen, the mean value of the correlation coefficient is 0.95, while it is 0.92 for $y$ -slice 1 at the bottom. Similar results were obtained for the notched specimen: 0.95 and 0.84, respectively. This decrease in correlation coefficient through the thickness can be attributed to a depth-dependent speckle contrast reduction as a result of signal attenuation due to light scattering, defocusing, and spectral roll-off.

Fig. 6

3-D views of the correlation coefficient maps for (a) rectangular phantom strip and (b) notched phantom strip, cut at the position $z = 25$ , load step 1.

In our experiments, the maximum $ε_{x x}$ was $\sim 1.4 %$ (see results for load step 2 in Sec. 3.2.2) and therefore speckle boiling was not an issue. Without further analysis it would seem prudent to consider NCC values $< 0.9$ as inappropriate to measure strains with OCT and DVC.

3.2.2.

DVC results for rectangular phantom strip

$U_{x}$ displacement maps for the rectangular phantom strip, which denote the displacement along the tension direction, are shown in Fig. 7. In Fig. 7, one can see the evolution of the $U_{x}$ displacement in cross-sections cut at different $z$ -slices. The absolute value of $U_{x}$ increases continuously along the $x$ direction from the fixed side to the other, ranging from 0 to 0.07 mm, as expected from the loading configuration. Figure 8(a) shows that the mean value of $U_{x}$ evolves linearly along the $x$ direction. Nevertheless, a sinusoidal oscillation is apparent when the difference between the actual values and a linear fit is plotted as a function of $x$ in Fig. 8(b). An analysis of these displacement oscillations is provided below.

Fig. 7

$U_{x}$ displacement distribution obtained for a rectangular phantom strip under tension (load step 1), showing different subregions of the data volume.

Fig. 8

(a) $U_{x}$ displacement averaged within $x$ slices along the $x$ -direction and a linear data fit obtained for a rectangular phantom strip under tension (load step 1). (b) Difference between the averaged $U_{x}$ and the linear fit shown in (a).

The strain maps were derived using the same procedure as for the noise analysis (centered finite difference of the displacement data, without any smoothing). All six strain components for central $z$ -slice 25 and central $y$ -slice 10 are shown in Figs. 9 and 10, respectively. As expected from the loading configuration, the $ε_{x x}$ strain maps for both $z$ - and $y$ - slices show positive values around 0.007. Strain maps $ε_{y y}$ and $ε_{z z}$ show negative values, indicating a Poisson’s contraction along the corresponding directions. It is interesting to note that $ε_{z z}$ is very small (close to zero) at the right-hand side, where the grip prevents Poisson’s contraction in the $z$ direction. Regarding $ε_{y y}$ , it is not zero at the right-hand side because the constraint from the grip acts only at the surface. The reason why $ε_{y y}$ is actually larger in magnitude may be from some material nonlinearity due to the compression in the grip. Figure 11 shows that the measured strain components are significantly larger than the corresponding strain noise obtained for the stationary test. Since this is a pure tensile test for the rectangular strip, all the shear strain components should be close to noise level, as confirmed in Figs. 9 and 10. Some irregularities, however, can be seen from these strain maps. The fringes observed in the $ε_{x x}$ strain component in Figs. 9 and 10 could be due to either a spatial variation of the elastic modulus or interpolation bias.⁴⁰ In the case of the FFT-based DIC algorithm used here, the period of the oscillation due to interpolation bias corresponds to a displacement equivalent to 2 voxels.⁴⁴^–⁴⁶ This is consistent with 2.5 fringes observed in Fig. 8(b) for a total deformation corresponding to 5 voxels. A bias in displacement directly leads to bias in strain, proportional to the slope of the displacement bias. When aliasing arises due to spatial undersampling by the SS-OCT system, the displacement values obtained when comparing subvolumes between reference and deformed states are likely to suffer from larger interpolation errors, which are further amplified when strain is calculated. It has been shown that aliasing typically shows up as a Moiré-like fringe pattern in the displacement and strain fields and that it is more obvious in the latter. In order to confirm the nature of the observed fringes, a second load step was performed on the same rectangular phantom strip with an extra 10 g dead weight (we refer to this case as load step 2). The strain maps showed twice as many fringes in Fig. 12, confirming that these are due to interpolation bias in this elastic material. One way of reducing this bias is to perform presmoothing on both the reference and deformed volume images using a Gaussian low-pass filter before correlation to reduce high-spatial frequency components⁴⁷^,⁴⁸ or to increase the sampling density of the OCT reconstruction.

Fig. 9

Spatial distributions of the strain components obtained for a rectangular phantom strip under tension (load step 1). The central $z$ slice 25 is shown.

Fig. 10

Spatial distributions of the strain components obtained for a rectangular phantom strip under tension (load step 1). The central $y$ slice 10 is shown.

Fig. 11

Normal strain components obtained at $z$ slice 25 for a rectangular phantom strip under tension (load step 1) and the corresponding strain noise obtained in the stationary test.

Fig. 12

Spatial distributions of the normal strain $ε_{x x}$ obtained for a rectangular phantom strip under tension (load step 2): (a) central $z$ slice 25 and (b) central $y$ slice 10.

3.2.3.

Results for notched phantom strip

$U_{x}$ displacement maps for the notched phantom strip are shown in Fig. 13. One shows the internal $U_{x}$ displacement in a cross-section cut at central $z$ -slice 25, while the other shows the whole volume. It can be seen that in each $z$ -slice, the absolute value of $U_{x}$ increases continuously along the $x$ direction from the fixed side to the other, as expected from the loading configuration and consistent with the displacement maps for the rectangular phantom strip as shown in Fig. 7. In Fig. 13, a slight bending of the strip can be observed from the larger $U_{x}$ displacement values at the top half of the strip compared to those at the bottom half. This is probably due to the slight geometrical asymmetry between the two notches. The geometrical asymmetry was induced during the manufacturing process when cutting the strip to a notched shape from a larger piece.

Fig. 13

$U_{x}$ displacement distribution obtained for a notched phantom strip under tension (load step 1), showing different subregions of the data volume.

For the notched phantom strip, all six strain components for central $z$ -slice 25 and central $y$ -slice 10 are shown in Figs. 14 and 15, respectively. The $ε_{x x}$ strain map shows positive values, while the $ε_{y y}$ and $ε_{z z}$ strain maps show negative values, as expected. Larger deformation is expected in the notched region, which can be observed in the maps of the normal strain components in Figs. 14 and 15. In addition, strain concentration is observed near the top notch tip of the strip for the normal strain components in Fig. 15. This is consistent with the larger $U_{x}$ displacement found near the top notch in Fig. 13. The explanation for this local strain concentration is the geometrical asymmetry of the notched strip, which has already been stated earlier. In Fig. 15, the $ε_{x z}$ strain map shows an antisymmetric shear strain distribution. Regarding the other two shear strain components, they are close to zero. These results are consistent for this type of test. It should be pointed out, however, that these results also suffer from the interpolation errors due to aliasing, especially evident in $ε_{x x}$ .

Fig. 14

Spatial distributions of the strain components obtained for a notched phantom strip under tension (load step 1). The central $z$ slice 25 is shown.

Fig. 15

Spatial distributions of the strain components obtained for a notched phantom strip under tension (load step 1). The central $y$ slice 10 is shown.

3.2.4.

Bias reduction using Gaussian presmoothing

Presmoothing using a Gaussian filter has proved effective in reducing interpolation bias.⁴⁷^,⁴⁸ A Gaussian filter with a kernel size of $7 \times 7 \times 7 voxels$ and a standard deviation of 1.5 voxel was applied to the OCT volume reconstructions of the tensile tests prior to correlation. In Fig. 16, for the rectangular phantom strip, the plots of the difference between the measured and the fitted mean $U_{x}$ displacement of each $x$ -slice as a function of $x$ position with and without presmoothing are compared. The difference is substantially reduced after presmoothing, which is expected according to Refs. 47 and 48. Gaussian smoothing was also applied to the volume images of the stationary and rigid body translation tests in order to check its influence in the strain resolution. The results show an increase in the strain standard deviations with presmoothing, generally ranging from $4 \times 10^{- 4}$ to $8 \times 10^{- 4}$ , as shown in Figs. 17(a) and 17(b) for the stationary and rigid body translation tests, respectively. The reason for this increase has been explained in Ref. 48, which states that the sum of squares of subset intensity gradient (SSSIG) is reduced after smoothing, and the standard deviation error is inversely proportional to the SSSIG value. In any case, these noise levels may still be considered as acceptable compared with the strain levels of the tensile tests (about one order of magnitude). Figure 18 shows the $ε_{x x}$ strain maps with and without presmoothing for the rectangular and notched phantom strips. The fringes due to interpolation bias are eliminated after presmoothing.

Fig. 16

Difference between the $U_{x}$ displacement averaged within $x$ slices along the $x$ -direction and a linear data fit obtained for a rectangular phantom strip under tension (load step 1), with and without presmoothing.

Fig. 17

Strain standard deviations obtained with $36^{3}$ voxels subvolume, 50% overlap, and Gaussian premoothing for (a) stationary test and (b) rigid body translation test.

Fig. 18

Spatial distributions of the normal strain $ε_{x x}$ obtained at central $y$ slice 10 for phantom strips under tension (load step 1). (a) Rectangular, no presmoothing. (b) Rectangular, with presmoothing. (c) Notched, no presmoothing. (d) Notched, with presmoothing.

3.3.

Identification Results

Using the experimental strain results in the linear equation system [Eq. (15)], the material stiffness components $Q_{x x}$ and $Q_{x y}$ can then be directly determined for the two tensile tests. The identification results with and without presmoothing are listed in Tables 1 and 2, respectively. From these results, Young’s modulus $E$ and Poisson’s ratio $ν$ can be calculated through the relations stated in Eq. (7). In order to provide a reference for the stiffness parameters obtained through the VFM, Young’s modulus and Poisson’s ratio were also calculated for the rectangular strip based on the assumption of constant uniaxial stress through the relation.

Eq. (17)

{\begin{matrix} E = \frac{σ_{x}}{ε_{x}} \\ ν = - \frac{ε_{y}}{ε_{x}} \end{matrix},

where uniform stress

σ_{x}

in the

y z

cross-section of the rectangular phantom strip was determined through the equation

σ_{x} = F / A

.

F

is the tension load and

A

is the

y z

cross-sectional area.

ε_{x}

and

ε_{y}

are the average values of the corresponding strain components over the whole field of view. Thus, both Young’s modulus and Poisson’s ratio can be derived and the results are listed in Tables 1 and 2. In both tensile tests, the results obtained from the VFM are consistent with each other and also with those calculated from the constant uniaxial stress assumption. This indicates that the 3-D VFM is an effective tool to identify the constitutive parameters of materials with nonuniform stress/strain states when the constant uniaxial stress assumption is no longer applicable. Moreover, the identification results with presmoothing are in agreement with those without presmoothing when comparing Table 1 to Table 2.

Table 1

Identified material elastic stiffness parameters without presmoothing.

	$Q_{x x}$ (MPa)	$Q_{x y}$ (MPa)	$E$ (MPa)	$ν$
Virtual fields method (VFM) (rectangular, load step 1)	3.45	2.44	1.44	0.41
VFM (rectangular, load step 2)	3.48	2.44	1.47	0.41
Constant stress (rectangular, load step 1)	3.65	2.64	1.43	0.42
VFM (notched, load step 1)	3.24	2.18	1.49	0.40

Table 2

Identified material elastic stiffness parameters with presmoothing.

	$Q_{x x}$ (MPa)	$Q_{x y}$ (MPa)	$E$ (MPa)	$ν$
VFM (rectangular, load step 1)	3.39	2.40	1.40	0.41
VFM (rectangular, load step 2)	3.26	2.23	1.45	0.41
Constant stress (rectangular, load step 1)	3.55	2.57	1.39	0.42
VFM (notched, load step 1)	3.18	2.09	1.51	0.40

It is interesting to note that when a larger preload of 115 g was used in the tensile test, both parameters increased ( $E$ from 1.44 to 1.82 MPa and $ν$ from 0.42 to 0.48). In both cases, a load of 10 g was used.

4. Conclusions

We have shown that DVC can provide, by means of a single channel OCT system, multicomponent displacement fields from which all the strain components required by inversion methods such as the 3-D VFM can be evaluated. OCT+DVC has low displacement sensitivity compared to phase-sensitive OCT elastography and seems appropriate for strain as large as $\sim 1.68 %$ (in the axial direction) at which point an incremental approach should be used to avoid speckle decorrelation. A strain uncertainty in the order of $\sim 4 \times 10^{- 4}$ to $8 \times 10^{- 4}$ was observed when Gaussian presmoothing is used to reduce bias noise. Strain below this uncertainty level would require an alternative approach, such as phase-sensitive OCT, capable of detecting subwavelength displacements with low noise. In the cases studied in this work, one fringe across a 10 mm field of view would correspond to a relative axial displacement equal to 0.45 μm and an average strain of $4.5 \times 10^{- 5}$ , an order of magnitude better than the uncertainty we report for OCT+DVC. Even though most phase-sensitive OCT elastography systems have so far focused on phase measurements with only axial sensitivity, a new system with sensitivity to all displacement components has been recently developed based on a wavelength scanning OCT system using multiple illumination directions and a single observation direction.⁴⁹ In this work, strain was evaluated with a centered finite difference operator applied to the displacement field. No displacement smoothing was used before strain calculation in order to achieve maximal spatial resolution with a view to further studies on thin biological tissues such as the vertebrate eye cornea. In cases where strain accuracy is paramount, a weighted-least squares strain estimator would be more appropriate.⁵⁰

Uniform and nonuniform 3-D strain fields measured with OCT+DVC were used to identify the elastic stiffness components of rectangular and notched silicone rubber phantoms using the VFM with 3-D manually defined virtual fields. The material moduli extracted from this approach are consistent with those calculated from the constant uniaxial stress approach. In order to test the proposed identification methodology, simple uniaxial tensile tests and isotropic materials were used.

This is the first time that volume strain data are derived by performing DVC on OCT reconstructed volumes. Future work will be aimed at applying this methodology to measuring the internal 3-D full-field deformation of biological tissues under more complex loading conditions and also identifying spatial distributions of constitutive parameters.

Acknowledgments

The authors are thankful to the China Scholarship Council and the Wolfson School of Mechanical Engineering, Loughborough University, for their financial support. They are also thankful to the reviewers for their constructive comments and insight.

References

1.

H. E. Couleset al., “Residual strain measurement for arc welding and localised high-pressure rolling using resistance strain gauges and neutron diffraction,” J. Strain Anal. Eng. Des., 47 (8), 576 –586 (2012). http://dx.doi.org/10.1177/0309324712453408 JSADDZ 0309-3247 Google Scholar

2.

C. Fernández-ValdivielsoI. R. MatıasF. J. Arregui, “Simultaneous measurement of strain and temperature using a fiber Bragg grating and a thermochromic material,” Sens. Actuators A Phys., 101 (1–2), 107 –116 (2002). http://dx.doi.org/10.1016/S0924-4247(02)00188-7 SAAPEB 0924-4247 Google Scholar

3.

T. C. Chuet al., “Applications of digital-image-correlation techniques to experimental mechanics,” Exp. Mech., 25 (3), 232 –244 (1985). http://dx.doi.org/10.1007/BF02325092 EXMCAZ 0014-4851 Google Scholar

4.

R. Moulartet al., “On the realization of microscopic grids for local strain measurement by direct interferometric photolithography,” Opt. Laser Eng., 45 (12), 1131 –1147 (2007). http://dx.doi.org/10.1016/j.optlaseng.2007.06.009 OLENDN 0143-8166 Google Scholar

5.

W. AnT. E. Carlsson, “Speckle interferometry for measurement of continuous deformations,” Opt. Laser Eng., 40 (5–6), 529 –541 (2003). http://dx.doi.org/10.1016/S0143-8166(02)00085-4 OLENDN 0143-8166 Google Scholar

6.

D. PostW. A. Baracat, “High-sensitivity Moiré interferometry—a simplified approach,” Exp. Mech., 21 (3), 100 –104 (1981). http://dx.doi.org/10.1007/BF02326365 EXMCAZ 0014-4851 Google Scholar

7.

B. L. Boyceet al., “Full-field deformation of bovine cornea under constrained inflation conditions,” Biomaterials, 29 (28), 3896 –3904 (2008). http://dx.doi.org/10.1016/j.biomaterials.2008.06.011 BIMADU 0142-9612 Google Scholar

8.

V. LibertiauxF. PasconS. Cescotto, “Experimental verification of brain tissue incompressibility using digital image correlation,” J. Mech. Behav. Biomed., 4 (7), 1177 –1185 (2011). http://dx.doi.org/10.1016/j.jmbbm.2011.03.028 JMBBCP 1878-0180 Google Scholar

9.

Q. Z. Fanget al., “Rate-dependent large deformation behavior of PC/ABS,” Polymer, 50 (1), 296 –304 (2009). http://dx.doi.org/10.1016/j.polymer.2008.10.042 POLMAG 0032-3861 Google Scholar

10.

Y. WangA. M. Cuitino, “Full-field measurements of heterogeneous deformation patterns on polymeric foams using digital image correlation,” Int. J. Solids Struct., 39 (13–14), 3777 –3796 (2002). http://dx.doi.org/10.1016/S0020-7683(02)00176-2 IJSOAD 0020-7683 Google Scholar

11.

E. M. Parsonset al., “Three-dimensional large-strain tensile deformation of neat and calcium carbonate-filled high-density polyethylene,” Polymer, 46 (7), 2257 –2265 (2005). http://dx.doi.org/10.1016/j.polymer.2005.01.045 POLMAG 0032-3861 Google Scholar

12.

B. K. Bayet al., “Digital volume correlation: three-dimensional strain mapping using x-ray tomography,” Exp. Mech., 39 (3), 217 –226 (1999). http://dx.doi.org/10.1007/BF02323555 EXMCAZ 0014-4851 Google Scholar

13.

F. Forsberget al., “3-D micro-scale deformations of wood in bending: synchrotron radiation μCT data analyzed with digital volume correlation,” J. Struct. Biol., 164 (3), 255 –262 (2008). http://dx.doi.org/10.1016/j.jsb.2008.08.004 JSBIEM 1047-8477 Google Scholar

14.

I. JandejsekO. JirousekD. Vavrik, “Precise strain measurement in complex materials using digital volumetric correlation and time lapse micro-CT data,” Procedia Eng., 10 1730 –1735 (2011). http://dx.doi.org/10.1016/j.proeng.2011.04.288 PERNBE 1877-7058 Google Scholar

15.

F. ForsbergC. R. Siviour, “3-D deformation and strain analysis in compacted sugar using x-ray microtomography and digital volume correlation,” Meas. Sci. Technol., 20 (9), 1 –8 (2009). http://dx.doi.org/10.1088/0957-0233/20/9/095703 MSTCEP 0957-0233 Google Scholar

16.

N. Limodinet al., “Influence of closure on the 3-D propagation of fatigue cracks in a nodular cast iron investigated by x-ray tomography and 3-D volume correlation,” Acta Mater., 58 (8), 2957 –2967 (2010). http://dx.doi.org/10.1016/j.actamat.2010.01.024 ACMAFD 1359-6454 Google Scholar

17.

A. Benoitet al., “3-D analysis from micro-MRI during in situ compression on cancellous bone,” J. Biomech., 42 (14), 2381 –2386 (2009). http://dx.doi.org/10.1016/j.jbiomech.2009.06.034 JBMCB5 0021-9290 Google Scholar

18.

C. Francket al., “Three-dimensional full-field measurements of large deformations in soft materials using confocal microscopy and digital volume correlation,” Exp. Mech., 47 (3), 427 –438 (2007). http://dx.doi.org/10.1007/s11340-007-9037-9 EXMCAZ 0014-4851 Google Scholar

19.

A. Germaneauet al., “3-D mechanical analysis of aeronautical plain bearings: validation of a finite element model from measurement of displacement fields by digital volume correlation and optical scanning tomography,” Opt. Laser Eng., 48 (6), 676 –683 (2010). http://dx.doi.org/10.1016/j.optlaseng.2010.01.010 OLENDN 0143-8166 Google Scholar

20.

S. A. Maskarinecet al., “Quantifying cellular traction forces in three dimensions,” Proc. Natl. Acad. Sci. U. S. A., 106 (52), 22108 –22113 (2009). http://dx.doi.org/10.1073/pnas.0904565106 PNASA6 0027-8424 Google Scholar

21.

T. GambichlerV. JaedickeS. Terras, “Optical coherence tomography in dermatology: technical and clinical aspects,” Arch. Dermatol. Res., 303 (7), 457 –473 (2011). http://dx.doi.org/10.1007/s00403-011-1152-x ADREDL 0340-3696 Google Scholar

22.

B. J. Kaluzyet al., “Spectral optical coherence tomography: a novel technique for cornea imaging,” Cornea, 25 (8), 960 –965 (2006). http://dx.doi.org/10.1097/01.ico.0000224644.81719.59 CORNDB 0277-3740 Google Scholar

23.

N. Hutchingset al., “Swelling of the human cornea revealed by high-speed, ultrahigh-resolution optical coherence tomography,” Invest. Ophthalmol. Vis. Sci., 51 (9), 4579 –4584 (2010). http://dx.doi.org/10.1167/iovs.09-4676 IOVSDA 0146-0404 Google Scholar

24.

A. F. Fercher, “Optical coherence tomography—development, principles, applications,” Med. Phys., 20 (4), 251 –276 (2010). http://dx.doi.org/10.1016/j.zemedi.2009.11.002 MPHYA6 0094-2405 Google Scholar

25.

D. F. KiernanW. F. MielerS. M. Hariprasad, “Spectral-domain optical coherence tomography: a comparison of modern high-resolution retinal imaging systems,” Am. J. Ophthalmol., 149 (1), 18 –31 (2010). http://dx.doi.org/10.1016/j.ajo.2009.08.037 AJOPAA 0002-9394 Google Scholar

26.

M. R. BryantP. J. McDonnell, “Constitutive laws for biomechanical modeling of refractive surgery,” J. Biomech. Eng., 118 (4), 473 –481 (1996). http://dx.doi.org/10.1115/1.2796033 JBENDY 0148-0731 Google Scholar

27.

K. AndersonA. E. SheikhT. Newson, “Application of structural analysis to the mechanical behaviour of the cornea,” J. R. Soc. Interface, 1 (1), 3 –15 (2004). http://dx.doi.org/10.1098/rsif.2004.0002 1742-5689 Google Scholar

28.

T. D. NguyenB. L. Boyce, “An inverse finite element method for determining the anisotropic properties of the cornea,” Biomech. Model Mechanobiol., 10 (3), 323 –337 (2011). http://dx.doi.org/10.1007/s10237-010-0237-3 1617-7959 Google Scholar

29.

F. PierronM. Grédiac, The Virtual Fields Method, Springer, New York (2012). Google Scholar

30.

F. PierronS. ZhavoronokM. Grédiac, “Identification of the through-thickness properties of thick laminated tubes using the virtual fields method,” Int. J. Solids Struct., 37 (32), 4437 –4453 (2000). http://dx.doi.org/10.1016/S0020-7683(99)00149-3 IJSOAD 0020-7683 Google Scholar

31.

M. GrédiacF. PierronY. Surrel, “Novel procedure for complete in-plane composite characterization using a single T-shaped specimen,” Exp. Mech., 39 (2), 142 –149 (1999). http://dx.doi.org/10.1007/BF02331118 EXMCAZ 0014-4851 Google Scholar

32.

M. GrédiacF. PierronA. Vautrin, “The Iosipescu in-plane shear test applied to composites: a new approach based on displacement field processing,” Compos. Sci. Technol., 51 (3), 409 –417 (1994). http://dx.doi.org/10.1016/0266-3538(94)90109-0 CSTCEH 0266-3538 Google Scholar

33.

M. Grédiac, “On the direct determination of invariant parameters governing the bending of anisotropic plates,” Int. J. Solids Struct., 33 (27), 3969 –3982 (1996). http://dx.doi.org/10.1016/0020-7683(95)00220-0 IJSOAD 0020-7683 Google Scholar

34.

M. GrédiacF. Pierron, “Applying the virtual fields method to the identification of elasto-plastic constitutive parameters,” Int. J. Plasticity, 22 (4), 602 –627 (2006). http://dx.doi.org/10.1016/j.ijplas.2005.04.007 IJPLER 0749-6419 Google Scholar

35.

Y. Pannieret al., “Identification of elasto-plastic constitutive parameters from statically undetermined tests using the virtual fields method,” Exp. Mech., 46 (6), 735 –755 (2006). http://dx.doi.org/10.1007/s11340-006-9822-x EXMCAZ 0014-4851 Google Scholar

36.

S. Avrilet al., “Stress reconstruction and constitutive parameter identification in plane-stress elasto-plasticity problems using surface measurements of deformation fields,” Exp. Mech., 48 (4), 403 –419 (2008). http://dx.doi.org/10.1007/s11340-007-9084-2 EXMCAZ 0014-4851 Google Scholar

37.

S. AvrilP. BadelA. Duprey, “Anisotropic and hyperelastic identification of in vitro human arteries from full-field optical measurements,” J. Biomech., 43 (15), 2978 –2985 (2010). http://dx.doi.org/10.1016/j.jbiomech.2010.07.004 JBMCB5 0021-9290 Google Scholar

38.

J. H. Kimet al., “Experimental characterization of rupture in human aortic aneurysms using a full-field measurement technique,” Biomech. Model Mechanobiol., 11 (6), 841 –853 (2012). http://dx.doi.org/10.1007/s10237-011-0356-5 BMMICD 1617-7959 Google Scholar

39.

J. Fuet al., “Assessment of corneal deformation using optical coherence tomography and digital volume correlation,” in Conf. Proc. SEM, 155 –160 (2013). Google Scholar

40.

M. A. SuttonJ. J. OrteuH. W. Schreier, Image Correlation for Shape, Motion and Deformation Measurements, Springer, New York (2009). Google Scholar

41.

J. CouplandJ. Lobera, “Optical tomography and digital holography,” Meas. Sci. Technol., 19 (7), 070101 (2008). http://dx.doi.org/10.1088/0957-0233/19/7/070101 MSTCEP 0957-0233 Google Scholar

42.

P. D. RuizJ. M. HuntleyJ. M. Coupland, “Depth-resolved imaging and displacement measurement techniques viewed as linear filtering operations,” Exp. Mech., 51 (4), 453 –465 (2011). http://dx.doi.org/10.1007/s11340-010-9393-8 EXMCAZ 0014-4851 Google Scholar

43.

V. Y. Zaitsevet al., “Elastographic mapping in optical coherence tomography using an unconventional approach based on correlation stability,” J. Biomed. Opt., 19 (2), 021107 (2014). http://dx.doi.org/10.1117/1.JBO.19.2.021107 JBOPFO 1083-3668 Google Scholar

44.

K. Madiet al., “Computation of full-field displacements in a scaffold implant using digital volume correlation and finite element analysis,” Med. Eng. Phys., 35 (9), 1298 –1312 (2013). http://dx.doi.org/10.1016/j.medengphy.2013.02.001 MEPHEO 1350-4533 Google Scholar

45.

T. AstaritaG. Cardone, “Analysis of interpolation schemes for image deformation methods in PIV,” Exp. Fluids, 38 (2), 233 –243 (2005). http://dx.doi.org/10.1007/s00348-004-0902-3 EXFLDU 0723-4864 Google Scholar

46.

G. BesnardF. HildS. Roux, “Finite-element displacement fields analysis from digital images: application to Portevin-Le Châtelier bands,” Exp. Mech., 46 (6), 789 –803 (2006). http://dx.doi.org/10.1007/s11340-006-9824-8 EXMCAZ 0014-4851 Google Scholar

47.

P. Lavaet al., “Assessment of measuring errors in DIC using deformation fields generated by plastic FEA,” Opt. Laser Eng., 47 (7–8), 747 –753 (2009). http://dx.doi.org/10.1016/j.optlaseng.2009.03.007 OLENDN 0143-8166 Google Scholar

48.

B. Pan, “Bias error reduction of digital image correlation using Gaussian pre-filtering,” Opt. Laser Eng., 51 (10), 1161 –1167 (2013). http://dx.doi.org/10.1016/j.optlaseng.2013.04.009 OLENDN 0143-8166 Google Scholar

49.

S. ChakrabortyP. D. Ruiz, “Measurement of all orthogonal components of displacement in the volume of scattering materials using wavelength scanning interferometry,” J. Opt. Soc. Am. A, 29 (9), 1776 –1785 (2012). http://dx.doi.org/10.1364/JOSAA.29.001776 JOAOD6 0740-3232 Google Scholar

50.

B. F. Kennedyet al., “Strain estimation in phase-sensitive optical coherence elastography,” Biomed. Opt. Express, 3 (8), 1865 –1879 (2012). http://dx.doi.org/10.1364/BOE.3.001865 BOEICL 2156-7085 Google Scholar

Citation Download Citation

Jiawei Fu, Fabrice Pierron, and Pablo D. Ruiz "Elastic stiffness characterization using three-dimensional full-field deformation obtained with optical coherence tomography and digital volume correlation," Journal of Biomedical Optics 18(12), 121512 (17 December 2013). https://doi.org/10.1117/1.JBO.18.12.121512

Published: 17 December 2013

Access the abstract

JOURNAL ARTICLE
16 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

CITATIONS

Cited by 54 scholarly publications.

Explore citations on Lens.org

RIGHTS & PERMISSIONS

Get copyright permission Get copyright permission on Copyright Marketplace

KEYWORDS

Optical coherence tomography

Speckle

Silicon

Point spread functions

Reconstruction algorithms

Digital image correlation

Light scattering

1.

Introduction

2.

Materials and Methods

2.1.

Materials

Fig. 1

2.2.

Experimental Setup and Image Acquisition

Fig. 2

2.3.

Digital Volume Correlation

Eq. (1)

Eq. (2)

2.4.

Evaluation of Measurement Performance

2.5.

Virtual Fields Method

Eq. (3)

Eq. (4)

Eq. (5)

Eq. (6)

Eq. (7)

Eq. (8)

Eq. (9)

Eq. (10)

Eq. (11)

Eq. (12)

Eq. (13)

Eq. (14)

Eq. (15)

3.

Results and Discussions

3.1.

Strain Noise Analysis

3.1.1.

Stationary test

Eq. (16)

Fig. 3

Fig. 4

3.1.2.

Translation test

3.1.3.

Strain-induced speckle decorrelation

Fig. 5

3.2.

DVC Results for Tensile Tests

3.2.1.

Correlation coefficient maps

Fig. 6

3.2.2.

DVC results for rectangular phantom strip

Fig. 7

Fig. 8

Fig. 9

Fig. 10

Fig. 11

Fig. 12

3.2.3.

Results for notched phantom strip

Fig. 13

Fig. 14

Fig. 15

3.2.4.

Bias reduction using Gaussian presmoothing

Fig. 16

Fig. 17

Fig. 18

3.3.

Identification Results

Eq. (17)

Table 1

Table 2

4.

Conclusions

Acknowledgments

References

Show All Keywords

Keywords/Phrases

Search In: