Journal of Magnetic Resonance 239 (2014) 16–22
Contents lists available at ScienceDirect
Journal of Magnetic Resonance journal homepage: www.elsevier.com/locate/jmr
Spatially resolved measurements of mean spin–spin relaxation time constants Ruben Emanuel Nechifor, Konstantin Romanenko, Florea Marica, Bruce J. Balcom ⇑ MRI Research Centre, Department of Physics, University of New Brunswick, 8 Bailey Drive, Fredericton, New Brunswick E3B 5A3, Canada
a r t i c l e
i n f o
Article history: Received 4 October 2013 Revised 20 November 2013 Available online 3 December 2013 Keywords: MRI CPMG prepared Centric Scan SPRITE Relaxation time mean Relaxation time distribution T2 mapping, SE-SPI
a b s t r a c t Magnetic Resonance measurements of the T2 distribution have become very common and they are a powerful way to probe microporous fluid bearing solids. While the structure of the T2 distribution, and changes in the structure, are often very informative, it is common to reduce the T2 distribution to a mean numeric quantity in order to provide a quantitative interpretation of the distribution. Magnetic Resonance Imaging measurements of the T2 distribution have recently been introduced, but they are time consuming, especially for 2 and 3 spatial dimensions. In this paper we explore a direct MRI measurement of the arithmetic mean of 1/T2, characterizing the distribution by using the initial slope of the spatially resolved T2 decay in a CPMG prepared Centric Scan SPRITE experiment. The methodology is explored with a test phantom sample and realistic petroleum reservoir core plug samples. The arithmetic mean of 1/T2 is related to the harmonic mean of T2. The mean obtained from the early decay is explored through measurements of uniform saturated core plug samples and by comparison to other means determined from the complete T2 distribution. Complementary data were obtained using SE-SPI T2 distribution MRI measurements. The utility of the arithmetic mean 1/T2 is explored through measurements of centrifuged core plug samples where the T2 distribution varies spatially. The harmonic mean T2 obtained from the early decay was employed to estimate the irreducible water saturation for core plug samples. Ó 2013 Elsevier Inc. All rights reserved.
1. Introduction T2 distribution measurements have become very important and common in investigations of fluids in microporous samples. The structure of the T2 distribution is often informative, but it is common to reduce the distribution to a mean numeric quantity, in order to provide a quantitative interpretation of the T2 distribution. T2 distribution measurements are most frequently bulk measurements, but the distribution may be spatially resolved [1]. Spatially resolved MRI measurements of the T2 distribution are relatively new [2,3]. Spatially resolved T2 distribution measurements are time consuming measurements in 2 and 3 spatial dimensions. The T2 distribution is most commonly determined by inverse Laplace transform (ILT) of the T2 decay. The ILT process is an illconditioned problem [4–9]. Noisy data hinders determination of the T2 distribution, therefore a better representation of the relaxation time distribution may be possible using means independently obtained [10]. This paper presents a direct way to determine the spatially resolved T2 relaxation time based on the arithmetic mean
of the relaxation rate. The arithmetic mean obtained in this way is independent of the inverse Laplace transform process. The arithmetic mean of the relaxation rate is obtained from the initial slope of the CPMG decay and is used to directly determine the harmonic mean T2. The accuracy of this harmonic mean is tested by comparing it with different relaxation time means determined from the T2 distribution. Spatially resolved measurements are undertaken with two pure phase encode MRI methods. Uniform phantoms and brine saturated porous rock core plugs are employed as test samples. Direct measurement of the mean T2 relaxation time may be advantageous for core analysis, or studies of different types of porous media. In this work the harmonic mean T2 relaxation time is used to estimate the irreducible water saturation for centrifuged core plugs. The irreducible water saturation obtained based on the T2 harmonic mean is compared with the irreducible water saturation obtained by the classic cutoff bulk volume of irreducible water (CBVI) method [11,12] and gravimetric measurements.
2. Mean T2 parameters and MRI methods ⇑ Corresponding author. Fax: +1 (506) 453 4581. E-mail address: [email protected] (B.J. Balcom). 1090-7807/$ - see front matter Ó 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.jmr.2013.11.012
The T2 CPMG echo decay experiment [13,14] is often fit to a single or multiple exponential decay model, with Eq. (1):
17
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
AðtÞ ¼ f ðtÞ ¼
n X t Ai e T 2i
ð1Þ
f ð0Þ ¼
i¼1
1 T2
The T2 distribution may be characterized by different means, for example the arithmetic mean, harmonic mean and logarithmic mean. These means are outlined below. The arithmetic mean (AM) [15,16] of the probability distribution of the relaxation times is given by Eq. (2):
ð2Þ
where T2i are the discrete values of the relaxation times and Pi(T2i) are the corresponding probabilities. The arithmetic mean is known to be sensitive to large extremes of the discrete values of the distribution [10,15]. The harmonic mean (HM) is defined [15,16] in terms of the T2 distribution by Eq. (3):
T 2HM
" #1 n n X X 1 ¼ Pi ðT 2i Þ P i ðT 2i Þ T 2i i¼1 i¼1
ð3Þ
where T2i are the discrete values of the relaxation times and Pi ðT 2i Þ T1 are the corresponding probabilities of the relaxation 2i rate. The harmonic mean is known to be sensitive to the small extremes of the discrete values of the distribution [10,15]. The logarithmic mean (LM) is defined in two ways [2,17], but in terms of the T2 distribution the logarithmic mean is calculated by Eq. (4):
T 2LM ¼ exp
Pn i¼1 ½P i lnðT 2i Þ Pn i¼1 P i ðT 2i Þ
ð4Þ
The logarithmic mean computed using Eq. (4) is identical to the geometric mean. Both are relatively insensitive to the extremes of the discrete values of the distribution [10,15]. 2.2. Arithmetic mean of the relaxation rate In this work we explore direct measurement of a related mean, the mean of the relaxation rate. Differentiating Eq. (1) which describes the CPMG decay curve, with respect to time we obtain:
f 0 ðtÞ ¼
n X 1 t exp Ai T 2i T 2i i¼1
ð5Þ
and for time t = 0 this (6) becomes:
f 0 ð0Þ ¼
n X 1 Ai T 2i i¼1
ð6Þ
f0 (0)represents the slope of the tangent line to the CPMG decay at maximum amplitude when t = 0. Eq. (6) can be rewritten:
f 0 ð0Þ ¼
Df Dt
Eq. (1) at t = 0 equals:
Experimentally f(0) represents the maximum signal amplitude of the CPMG decay. If this maximum amplitude cannot be directly acquired experimentally, it can be determined by extrapolation of the early decay. Dividing Eq. (7) by Eq. (8), in accord with Eq. (2), the arithmetic mean of the relaxation rate will be obtained:
2.1. Mean T2s
Pn i¼1 T 2i P i ðT 2i Þ P n i¼1 P i ðT 2i Þ
ð8Þ
i¼1
where A(t) is the time dependent signal observed, Ai is the amplitude of the ith component of the signal, T2i is the ith component of the relaxation time and t is the time. If there is a distribution of the relaxation time it is common to use the ILT in order to determine the T2 distribution. From the T2 distribution it is possible to determine a variety of different mean T2s.
T 2AM ¼
n X Ai ¼ Amax
ð7Þ
Df f 0 ð0Þ ¼ ¼ Dt ¼ f ð0Þ f ð0Þ
Pn
1 i¼1 Ai T 2i Pn i¼1 Ai
¼
Initial Slope Amax
ð9Þ
Based on Eq. (9) the initial slope and the maximum amplitude can be used to obtain the arithmetic mean of the relaxation rate. The brackets h i in Eq. (9) indicate the arithmetic mean of the relaxation rate. As suggested by Eq. (3) the arithmetic mean of the relaxation rate is related to the harmonic mean of the relaxation time, as Eq. (10) describes:
Pn Ai Amax T 2HM ¼ P i¼1 ¼ n 1 Initial Slope i¼1 Ai T 2i
ð10Þ
The arithmetic mean of the relaxation rate can be used to determine the harmonic mean of the relaxation time, Eq. (11):
1 T2
¼
1 T 2HM
ð11Þ
The relationship of Eq. (11) has been verified through simulation of model data sets with a variety of T2 distributions. The early decay approach has been used successfully in MR measurements of motion [18]. It has recently been employed in MRI, for velocity imaging [19,20]. 2.3. NMR and MRI pulse sequences Experimental data were acquired using bulk CMPG and three MRI techniques: CPMG-prepared Centric Scan SPRITE (Single Point Ramped Imaging with T1 Enhancement) [3], SPRITE [21] and SE-SPI (Spin-Echo Single Point Imaging) [2]. Both pure phase-encode MRI techniques were employed for one dimensional T2 mapping experiments. Fig. 1(a) presents the CPMG-prepared Centric Scan SPRITE pulse sequence and (b) presents the SE-SPI pulse sequence. The CPMGprepared Centric Scan SPRITE pulse sequence commences with a CPMG preparation followed by a 90° RF pulse for Z-storage, and a spoiler gradient. The preparation process can have a variable numbers of echoes. The last part of the pulse sequence consists of a 1D Centric Scan SPRITE image for spatial encoding. The pure phase encoding SE-SPI pulse sequence, Fig. 1(b), contains a phase encode gradient between the 90° and 180° RF pulses. The phase encoded magnetization is then read out through multiple refocusing. From each scan a single k-space point is acquired for all echoes. These k-space data sets are then Fourier transformed to generate a series of T2-weighted profiles. Based on these profiles the local CPMG decay can be extracted. The SE-SPI experiment has limitations on the duration of the first echo, therefore this experiment will not permit spatial resolution of shorter T2 components. A full CPMG decay with the CPMGprepared Centric Scan SPRITE technique is a time consuming experiment. Nevertheless the CPMG-prepared SPRITE technique [3] can be applied successfully, in a time efficient manner, if we reduce considerably the number of profiles which comprise the CPMG decay. In this case the harmonic mean T2 is measured only from the initial slope, as suggested by Eqs. (9) and (10).
18
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
Fig. 1. (a) CPMG-prepared Centric Scan SPRITE method employed to measure T2weighted profiles. These profiles are used to generate local CPMG decays. The harmonic mean T2 is determined based on the initial slope. (b) SE-SPI pulse sequence employed to measure spatially resolved T2 decay. The T2 distribution may be determined by inverse Laplace Transform and early echo decays can be used to determine the harmonic mean T2.
3. Results and discussion
Fig. 2. (a) T2-weighted profiles (1st, 20th, 40th, 70th, 100th, 130th, 210th and 320th profile) obtained with CPMG-prepared Centric Scan SPRITE measurements for H2O doped with GdCl3 using an echo time TE = 2 ms. (b) T2 values spatially resolved for the water phantom. Based on a single exponential fit of the CPMG-prepared SPRITE experiment (d) T2 varied spatially between 189 and 195 ms. The bulk CPMG value of T2 (–) 207 ms was assigned to each pixel as a reference. The harmonic mean obtained from the early decay (s) had T2 values between 224 and 244 ms. The T2 values with highest probability extracted from the T2 distribution (.) varied between 169 and 179 ms.
3.1. Phantom measurements Measurements commenced with a uniform doped water phantom, where we anticipate single exponential T2 decay. Bulk CPMG measurements were undertaken to test the T2 behavior and the data was analyzed in three ways: (i) single exponential fitting, (ii) inverse Laplace transformation and (iii) analysis of the early decay mean. Single exponential decay fitting gave a bulk T2 of 207 ms. The ILT procedure, as anticipated, gave a narrow distribution of T2s centered on 210 ms. The T2 distribution was used to calculate the means described in Section 2.1. As anticipated for a single exponential decay, the various means yielded very similar results. The arithmetic mean was 207 ms, the harmonic mean was 205 ms, and the logarithmic mean was 206 ms. Analysis of the harmonic mean T2 based on the early decay, Section 2.2, resulted in a harmonic mean T2 of 210 ms, which also closely agrees with the above results. The same phantom sample was imaged with the CPMG-prepared Centric Scan SPRITE pulse sequence. 320 T2-weighted profiles were acquired for the phantom sample. Fig. 2(a) shows select T2-weighted profiles. Fig. 2(b) shows the spatially resolved T2 values. The solid line is the bulk CPMG single exponential result of T2 = 207 ms, which was assigned as a reference value for each pixel. One data set in Fig. 2(b) is based on single exponential decay fitting, for which the mean over all pixels
was T2 = 194 ms. Another data set shows the most probable T2 determined by ILT of the CPMG-prepared Centric Scan SPRITE data, for which the mean over all pixels was T2 = 178 ms. The figure also shows the harmonic mean T2 determined based on the early slope of the CPMG decay. The mean over the displayed pixels is 235 ms. The spatially resolved mean T2 values do not agree as well as for the bulk measurement, presumably due to reduced SNR and complications of the imaging process. Analogous measurements with the SE-SPI pulse sequence show better agreement. Fig. 3(a) shows select T2-weighted profiles acquired using the SE-SPI pulse sequence. Fig. 3(b) shows the spatially resolved T2 values. Three sets of data are also reported from SE-SPI data for each pixel: (i) the single exponential decay fitting, (ii) the most probable T2 in each pixel determined by ILT and (iii) the harmonic mean T2 determined based on the early slope of the CPMG decay in each pixel. The T2 pixel mean for each data set reported are: (i) T2 = 208 ms, (ii) T2 = 209 ms and (iii) T2 = 217 ms respectively. The solid line again represents the bulk CPMG single exponential result of T2 = 207 ms, which we assign to each pixel as a reference. The simple water phantom results show that the early decay yields a high quality harmonic mean T2 consistent with the bulk measurement. CPMG-prepared SPRITE measurement and SE-SPI measurement give similar results.
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
19
Fig. 4. The T2 distribution generated by ILT of a bulk CPMG measurement of a brine saturated Berea core plug. The short lifetime extreme of the distribution is unreliable, T2 values reported are much shorter than TE. This will bias the harmonic mean T2 calculated from the T2 distribution.
Fig. 3. (a) T2-weighted profiles (1st, 20th, 40th, 70th, 100th, 130th, 210th and 512th profile) obtained with SE-SPI measurements for H2O doped with GdCl3 using an echo time TE = 2 ms. (b) T2 values spatially resolved in the water phantom. Based on a single exponential fit of the SE-SPI experiment (d) T2 varied spatially between 206 and 209 ms. The bulk CPMG value of T2 (–) 207 ms was assigned to each pixel as a reference. The harmonic mean obtained from the early decay (s) had T2 values between 211 and 222 ms. The T2 values with highest probability extracted from the T2 distribution (.) varied between 206 and 210 ms.
3.2. Uniform core plug measurements We can anticipate that the various means described in Sections 2.1 and 2.2 will not yield the same result in the case of samples with T2 distributions that are not unimodal and which extend over orders of magnitude in T2. The various means have different sensitivities to long and short T2 components in the distribution [10,15]. The T2 distribution of water saturated reservoir core plugs are of intense experimental interest and they are frequently bimodal with a T2 range of three to four orders of magnitude. Bulk and spatially resolved T2 measurements were undertaken on a brine saturated Berea rock core plug. We assume for a fully saturated core plug sample that although the sample is microscopically heterogeneous (there is a pore size distribution) it is macroscopically homogeneous (the pore size distribution is the same throughout the sample). The fully saturated Berea core plug presents a bi-modal distribution when the bulk decay curve is fit through an ILT, Fig. 4. The T2 means described in Sections 2.1 and 2.2 were determined from the bulk measurements. The logarithmic mean T2 computed from the T2 distribution was 28 ms. The arithmetic mean was 53 ms, while the harmonic mean was 10 ms. The T2 harmonic mean based on the early CPMG decay was 30 ms. The harmonic mean calculated from the bulk T2 distribution is surprisingly different from the early decay harmonic mean T2. The divergent experimental results were explored through simulations with model data sets. Simulations revealed that the T2 distribution will generate the same harmonic mean as the
corresponding T2 early decay under specific conditions imposed by Eqs. (9) and (10). Simulations showed high quality Amax and high quality early decay rate data yield the same results for the harmonic mean from the early decay and from the T2 distributions. A high quality result requires that the echo time is 5 times shorter than the shortest T2 value and noise should not be more than 10% of maximum signal. A short echo time is required in order to obtain at least 3–5 collinear points in the early decay. These conditions are satisfied for our Berea core plug data. The difference between the harmonic means computed based on the early decay and based on the T2 distribution results from the fact that one HM is computed based on the early decay and the other HM is based on the T2 distribution. The T2 distribution of Fig. 4 shows T2 values significantly shorter than the echo time (2 ms). These values do not make physical sense. The harmonic mean T2 computed based on the T2 distribution is strongly affected by these short components [10,15], generating a harmonic mean of 10 ms instead of 30 ms, the value obtained from the early decay data. These results, and simulations run for different data sets, prove that the early CPMG decay can be used successfully to compute the harmonic mean. The harmonic mean is particularly sensitive to short lifetime components from the T2 distribution. The harmonic mean directly calculated from the early decay avoids difficulties associated with the ILT procedure. Fig. 5(a) presents select T2-weighted profiles obtained for the fully saturated Berea sample with the SE-SPI pulse sequence. Based on these profiles, local T2 decays were extracted. Fig. 5(b) shows various T2 means. The harmonic mean obtained based on the early decay was T2 = 39 ms when averaged over all pixels displayed. Other means based on the T2 distribution, once more reported as means over all pixels, were arithmetic mean T2 68 ms, harmonic mean T2 18 ms and logarithmic mean T2 41 ms. Fig. 5(b) shows that with a real distribution there will be differences between different means even if these means are computed based on the same T2 distribution. The difference between these mean T2s results from the known sensitivity of these means to extremes of the distribution. Extremes of the T2 distribution are associated with large or small pores yielding long or short lifetime T2s. If long lifetime components are present in the T2 distribution then the arithmetic mean T2 will be significantly influenced by those components, while the short lifetime components will affect the harmonic mean T2 [10,15]. The logarithmic and geometric means are not affected by extremes in the T2 distribution. The new initial slope procedure used to compute the harmonic mean T2 avoids issues associated with the ILT procedure but the mean will still be strongly influenced by short lifetime T2 components.
20
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
Fig. 5. (a) The T2-weighted profiles (1st, 4th, 8th, 16th, 32nd, 56th, 80th, 128th and 256th profiles) obtained with SE-SPI measurements of a Berea core plug fully saturated with brine using an echo time TE = 2 ms. (b) Harmonic mean T2 obtained from the initial CPMG decay (d) with spatial variation between 37 and 41 ms. The arithmetic mean computed from the T2 distribution (r) with T2 varying between 66 and 70 ms. The harmonic mean computed from the T2 distribution (j) with T2 varying between 16 and 19 ms. The logarithmic mean computed from the T2 distribution (N) with T2 varying between 39 and 43 ms.
Different T2 means have been related to petro-physical parameters in the literature. The geometric mean T2 is a good parameter for predicting sandstone permeability [20] and irreducible watersaturation in sandstones [10,20]. Sandstone permeability and hydrocarbon viscosity have been correlated with the logarithmic mean T2. The pore size distribution can be related to the width of the T2 distribution, while in hydrocarbon studies the mean chain length is related to the relaxation times [5,6,22]. 4. Centrifuged core plugs measurements 4.1. Heterogeneous core plugs measurements Having explored various T2 means and the early decay harmonic mean, we turn to macroscopically inhomogeneous samples where we anticipate spatially resolved differences in the T2 distribution. Based on previous work [2] we consider centrifuged core plugs where spatial variation of saturation is associated with differences in the T2 distribution. For this class of sample, the CPMG-prepared SPRITE technique was used. Fig. 6(a) shows the first 10 T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a brine saturated Berea core plug. The early decay harmonic mean T2 was determined for different positions through the sample, in Fig. 6(b). The early decay harmonic mean T2 varies from 24 ms to 38 ms within the field of view. The harmonic mean T2 was computed based on the first 10 collinear decay points, in order to ensure a high quality result.
Fig. 6. (a) The first 10 T2-weighted profiles of the Berea sandstone, acquired using the CPMG-prepared SPRITE technique after centrifugation of the core plug. (b) The harmonic mean T2 obtained based on the early slope of the CPMG extracted from the T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a centrifuged Berea core plug. The harmonic mean T2 varies spatially from 24 ms to 38 ms for positions between 1.6 cm and 4.7 cm.
Fig. 6(b) shows less variation of the mean T2 compared with previous results obtained by Li [2]. In Fig. 6(b) the early decay harmonic mean T2 is computed, while Li computed the logarithmic mean T2. The shallow slope of Fig. 6(b) is due to the relative insensitivity of the early decay harmonic mean T2 to changes produced by emptying large pores which are associated with long T2 components. A brine saturated Buff Berea core plug was also examined. For this sample the first 10 T2-weighted profiles were measured with a short echo time of TE = 500 ls, as can be observed in Fig. 7(a). A short echo time measurement results in less visible profile attenuation but the early decay mean is still readily calculated. Fig. 7(b) shows the harmonic mean plotted for different positions through the sample obtained based on the T2-weighted profiles. The harmonic mean of the relaxation time shown in Fig. 7(b) varies spatially as was anticipated, for positions between 2.8 and 5.8 cm. The harmonic mean T2 varied from 10 ms to 30 ms. Data at the edge of the core plug was not analyzed, Fig. 7(a) shows an obvious artefact at the object left edge.
4.2. Harmonic mean and irreducible water saturation determination The irreducible water saturation is a critical parameter for core plug analysis. It has been shown that the irreducible water saturation, Swirr, can be obtained using mean T2 values [2,12]. A linear relationship between the logarithmic mean T2 and the residual water saturation has been established [2,12]. We now consider if the harmonic mean is similarly related to local saturation.
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
21
Fig. 8. The harmonic mean T2, computed based on the early CPMG decay using the CPMG-prepared SPRITE technique, versus saturation of a centrifuged Buff Berea core plug. The irreducible water saturation was estimated to be 21% based on the analysis of the data with Eq. (12).
Fig. 7. (a) The first 10 T2-weighted profiles of the Buff Berea sandstone, acquired using the CPMG-prepared SPRITE technique. (b) The harmonic mean T2 obtained based on the slope of the CPMG extracted from the T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a centrifuged brine saturated Buff Berea core plug. The harmonic mean T2 varies spatially from 10 ms to 30 ms for positions between 2.8 cm and 5.8 cm. The edges of the sample were not measured.
Li et al. [2] proposed based on Coates et al. [11,12] a linear relationship between the T2LM and residual water saturation, as described by the equation:
SWðT 2LM Þ ¼ m T 2LM Swirr þ Swirr
ð12Þ
where m is a parameter which describes the pore geometry, SWðT 2LM Þ is the residual water saturation, T2LM is the logarithmic mean T2 associated with the ith position of the sample, and Swirr is the irreducible water saturation. In this paper we seek to determine if the relation of Eq. (12) remains true when we employ a different mean T2. In particular we wish to employ the T2HM determined from the early echo decay. Knowing the T2HM for different saturations of the sample, we obtained the irreducible water saturation by plotting harmonic mean T2 versus the residual water saturation. The local saturation is determined based on simple 1D DHK Centric Scan SPRITE profiles of fully saturated and centrifuged core plugs. The spin density measured with SPRITE pulse sequence for a fully brine saturated core plug represents SW = 1 (100% saturation). Measuring the spin density for a centrifuged core plug and keeping constant the experimental parameters for SPRITE pulse sequence, the variation of saturation across the sample was established. Fig. 8 shows the harmonic mean T2 obtained using CPMGprepared SPRITE measurement for a centrifuged Buff Berea core plug. Fitting the decay of water saturation versus the T2 mean with Eq. (12) it was determined that the irreducible water saturation for Buff Berea was 21%. Using the classical cutoff BVI techinique [12], the irreducible water saturation was estimated to be 25%. From gravimetric measurements an irreducible water saturation of 22%
was obtained. These results show that for Buff Berea, or for core plugs which present similar pore structure, the harmonic mean T2 of centrifuged samples could be used to determine the irreducible water saturation. The irreducible water saturation was also calculated for the Berea core plug. Using Eq. (12) it was determined that the irreducible water saturation for the centrifuged Berea was 20%. The irreducible water saturation was computed considering only the positions between 2.1 cm and 4.6 cm, neglecting the edges of the sample. The result in this case is not very reliable because of the scatter in the T2 mean for the short lifetime component. From gravimetric measurements the irreducible water saturation was known to be 18%. The conventional T2 cutoff value of 33 ms used in classical cutoff BVI technique for sandstone clearly gives an erroneous irreducible water saturation of 50%.
5. Experimental The experiments were performed on an 8.51 MHz Oxford Maran DRX HF instrument (Oxford Instruments, Witney, UK). All measurements were carried out at ambient magnet temperature of 25 °C. The T2 distribution was determined employing WinDXP software (Oxford Instruments, Witney, UK). The profiles obtained from CPMG-prepared SPRITE experiments were computed employing various routines written with Interactive Data Language (Exelis, Boulder, USA). Matlab Math Works Natick, MA, USA was employed for simulation. SE-SPI measurements were applied to a uniform brine saturated Berea. Trapezoidal phase encode gradients were employed with ramp times of 150 ls and a plateau time of 100 ls. The time between the first two RF pulses was 2 ms and 1024 echoes were acquired. In order to increase the signal-to-noise ratio, 20 points from each echo were acquired. The field of view (FOV) was 9.6 cm with two scans acquired. To maintain the phase shift upon refocusing, an XY-16 phase cycle was applied [23] and to compensate RF field heterogeneity, composite pulses were employed. The CPMG-prepared Centric Scan SPRITE measurements were undertaken with a phantom sample, doped with GdCl3 and brine saturated core plugs. GdCl3 (0.1 mM) was added to reduce T1 to 220 ms. For the phantom sample and brine saturated Berea the CPMG-prepared Centric Scan SPRITE technique had an encoding time of 400 ls, FOV = 60 mm, NS = 32, 4-step phase cycle, gradient spoil duration of 500 ls at 2 G/cm, P90 = 32 ls, P180 = 64 ls and 64 k-space points were acquired. The flip angle was 4°, the sweep
22
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
width was 50 kHz, the probe ring-down time was 5 ls, and the filter stabilization time was 284 ls. The echo time was TE = 2 ms and 1024 echoes were acquired. The CPMG-prepared Centric Scan SPRITE experiment employs longer echo time than necessary in order to permit a direct comparison with SE-SPI experiments. For the phantom sample, 320 profiles were acquired in 6 h for the same TE and different numbers of echoes in order to record a full CPMG decay. For the brine saturated Berea core plug the flip angle was 9°, and only 10 profiles were acquired. For each profile 128 averages were acquired in 8 min. For the brine saturated Buff Berea core plug the CPMG-prepared Centric Scan SPRITE technique had an encoding time of 150 ls, FOV = 90 mm, NS = 512, 4-step phase cycle, gradient spoil duration of 500 ls at 2 G/cm, P90 = 43 ls, P180 = 86 ls and 64 k-space points were acquired. The flip angle was 9°, the sweep width was 125 kHz, the probe ring-down time was 30 ls, and the filter stabilization time was 26 ls with an echo time of TE = 500 ls. Only 10 profiles were acquired and for each profile 512 averages were acquired in 22 min. The gradient stabilization time was 2 ms in all cases. Two different core plugs were used: Berea and Buff Berea (Kocurek Industries, Caldwell, TX, USA). Both were cylindrical core plugs, 3.7 cm in diameter and 5 cm in length. The porosities were 21% for Berea and 22% for Buff Berea as determined gravimetrically. The core plugs were saturated with brine (1.5% wt. NaCl in H2O) under vacuum for more than 48 h. To ensure partial saturation the samples were centrifuged for 60 min at 2000 rpm at room temperature. A Beckman Model J2-12M centrifuge was employed. The variation of saturation across the samples was computed employing Centric Scan SPRITE measurements of fully and partially saturated samples. The SPRITE pulse sequence parameters were the same as the parameters used for CPMG-prepared SPRITE. 6. Conclusions The harmonic mean T2 obtained from the initial slope of the CPMG decay is a powerful way to characterize the mean T2. The accuracy of the harmonic mean T2 computed based on the early CPMG decay is tested in this paper in various ways: through simulations, on phantom sample and on uniform and partial saturated core plugs. The harmonic mean was also compared with different means computed based on the T2 distribution. The experimental data were obtained using bulk CPMG technique, CPMG-prepared Centric Scan SPRITE technique and SE-SPI technique [2,24–27]. The advantage of the CPMG-prepared Centric Scan SPRITE technique is that it can have a short echo time consistent with bulk CPMG measurements. Although spatial encoding has been explored with phase encoding, frequency encoding with short echo time is also possible. The harmonic mean computed from the early decay does not require intermediate data manipulation, as with the inverse Laplace Transform approach. Therefore the harmonic mean is not affected by different parameters used in processing the data. Regardless of the method chosen to determine the harmonic mean, it was shown that it will still be biased by the short lifetime T2 component. Using this new technique the harmonic mean of the relaxation time or the arithmetic mean of the relaxation rate can be obtained. The early decay harmonic mean T2 can be used to study the variation of saturation for porous media and the harmonic mean T2 based on early decay for different positions in the sample can be obtained. The harmonic mean T2 can be used also in other studies, such as studies of fluid flow process through the sample.
Acknowledgments B.J.B. thanks NSERC for a Discovery grant and an Idea to Innovation grant. B.J.B. thanks the Canada Chairs program for a Research Chair in MRI of Materials. REN thanks R. MacGregor for his assistance with the experiments and B. MacMillan for helpful discussions. References [1] R. Kimmich, NMR Tomography, Diffusometry, Relaxometry, Springer-Verlag, Berlin, 1997. [2] L. Li, H. Han, B.J. Balcom, Spin echo SPI methods for quantitative analysis of fluids in porous media, J. Magn. Reson. 198 (2009) 252–260. [3] O.V. Petrov, G. Ersland, B.J. Balcom, T2 distribution mapping profiles with phase-encode MRI, J. Magn. Reson. 209 (2011) 39–46. [4] L. Venkataramanan, Y.-Q. Song, M.D. Hürlimann, Solving Fredholm integrals of the first kind with tensor product structure in 2 and 2.5 dimensions, IEEE Trans. Signal Process. 50 (2002) 1017–1026. [5] L. Venkataramanan, F.K. Gruber, T.M. Habashy, D.E. Freed, Mellin transform of CPMG data, J. Magn. Reson. 206 (2010) 20–31. [6] L. Venkataramanan, F.K. Gruber, T.M. Habashy, D.E. Freed, Continuous moment estimation of CPMG data using Mellin transform, J. Magn. Reson. 216 (2012) 43–52. [7] M. Simina, R. Nechifor, I. Ardelean, Saturation-dependent nuclear magnetic resonance relaxation of fluids confined inside porous media with micrometersized pores, Magn. Reson. Chem. 49 (2011) 314–319. [8] J.G. McWhirter, E.R. Pike, On the numerical inversion of the Laplace transform and similar Fredholm integrals of the first kind, J. Phys. A: Math. Gen. 11 (1978) 1729–1745. [9] C.L. Epstein, J. Schotland, The bad truth about Laplace’s transform, SIAM Rev. 50 (2008) 504–520. [10] G.C. Borgia, R.J.S. Brown, P. Fantazzini, Different ‘‘average’’ nuclear magnetic resonance relaxation times for correlation with fluid-flow permeability and irreducible water saturation in water-saturated sandstone, J. Appl. Phys. 82 (9) (1997) 4197–4204. [11] G.R. Coates, D. Marschall, D. Mardon, J. Galford, A new characterization of bulkvolume irreducible using magnetic resonance, Log Anal. 39 (1998) 51–63. [12] G.R. Coates, L. Xiao, M.G. Prammer, NMR Logging – Principles and Applications, Halliburton Energy Service, Houston, 1999. [13] S. Meiboom, D. Gill, Modified spin echo method for measuring nuclear relaxation times, Rev. Sci. Instrum. 29 (1958) 688–691. [14] P.T. Callaghan, Principles of Nuclear Magnetic Resonance Microscopy, Clarendon Press, Oxford, 1991. pp. 339–340. [15] G.H. Hardy, J.E. Littlewood, G. Pólya, Inequalities, Cambridge University Press, London, 1934. [16] M. Abramowitz, I.A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs and Mathematical Tables, National Bureau of Standards, Washington, DC, 1964. [17] K.B. Stolarsky, Generalizations of the logarithmic mean, Math. Mag. 48 (2) (1975) 87–92. [18] P. Dyverfeldt, A. Sigfridsson, H. Knutsson, T. Ebbers, A novel MRI framework for the quantification of any moment of arbitrary velocity distributions, Magn. Reson. Med. 65 (2011) 725–731. [19] K.V. Romanenko, B.J. Balcom, Permeability mapping in porous media by magnetization prepared centric-scan SPRITE, Exp. Fluids 50 (2011) 301–312. [20] G.C. Borgia, R.J.S. Brown, P. Fantazzini, Estimates of permeability and irreducible water saturation by means of a new robust computation of fractional power average relaxation times, Magn. Reson. Imaging 16 (1998) 613–615. [21] B.J. Balcom, R.P. MacGregor, S.D. Beyea, D.P. Green, R.L. Armstrong, T.W. Bremner, Single point ramped imaging with T1 enhancement (SPRITE), J. Magn. Reson. A 123 (1996) 131–134. [22] D.E. Freed, Dependence on chain length of NMR relaxation times in mixtures of alkanes, J. Chem. Phys. 126 (2007) 174502-1–174502-10. [23] T. Gullion, The effect of amplitude imbalance on compensated Carr-Purcell sequences, J. Magn. Reson. Ser. A 101 (1993) 320–323. [24] B.J. Balcom, SPRITE imaging of short relaxation time nuclei, in: B. Blumich, W. Kuhn (Eds.), Spatially Resolved Magnetic Resonance: Methods, Materials, Medicine, Biology, Rheology, Geology, Ecology, Hardware, Wiley-VCH, Weinheim, 1998, pp. 75–86. [25] A.V. Ouriadov, R.P. MacGregor, B.J. Balcom, Thin film MRI – high resolution depth imaging with a local surface coil and spin echo SPI, J. Magn. Reson. 169 (2004) 174–186. [26] S. Gravina, D.G. Cory, Sensitivity and resolution of constant-time imaging, J. Magn. Reson. Ser. B 104 (1994) 53–61. [27] O.V. Petrov, B.J. Balcom, Two-dimensional T2 distribution mapping in porous solids with phase encode MRI, J. Magn. Reson. 212 (2011) 102–108.
Contents lists available at ScienceDirect
Journal of Magnetic Resonance journal homepage: www.elsevier.com/locate/jmr
Spatially resolved measurements of mean spin–spin relaxation time constants Ruben Emanuel Nechifor, Konstantin Romanenko, Florea Marica, Bruce J. Balcom ⇑ MRI Research Centre, Department of Physics, University of New Brunswick, 8 Bailey Drive, Fredericton, New Brunswick E3B 5A3, Canada
a r t i c l e
i n f o
Article history: Received 4 October 2013 Revised 20 November 2013 Available online 3 December 2013 Keywords: MRI CPMG prepared Centric Scan SPRITE Relaxation time mean Relaxation time distribution T2 mapping, SE-SPI
a b s t r a c t Magnetic Resonance measurements of the T2 distribution have become very common and they are a powerful way to probe microporous fluid bearing solids. While the structure of the T2 distribution, and changes in the structure, are often very informative, it is common to reduce the T2 distribution to a mean numeric quantity in order to provide a quantitative interpretation of the distribution. Magnetic Resonance Imaging measurements of the T2 distribution have recently been introduced, but they are time consuming, especially for 2 and 3 spatial dimensions. In this paper we explore a direct MRI measurement of the arithmetic mean of 1/T2, characterizing the distribution by using the initial slope of the spatially resolved T2 decay in a CPMG prepared Centric Scan SPRITE experiment. The methodology is explored with a test phantom sample and realistic petroleum reservoir core plug samples. The arithmetic mean of 1/T2 is related to the harmonic mean of T2. The mean obtained from the early decay is explored through measurements of uniform saturated core plug samples and by comparison to other means determined from the complete T2 distribution. Complementary data were obtained using SE-SPI T2 distribution MRI measurements. The utility of the arithmetic mean 1/T2 is explored through measurements of centrifuged core plug samples where the T2 distribution varies spatially. The harmonic mean T2 obtained from the early decay was employed to estimate the irreducible water saturation for core plug samples. Ó 2013 Elsevier Inc. All rights reserved.
1. Introduction T2 distribution measurements have become very important and common in investigations of fluids in microporous samples. The structure of the T2 distribution is often informative, but it is common to reduce the distribution to a mean numeric quantity, in order to provide a quantitative interpretation of the T2 distribution. T2 distribution measurements are most frequently bulk measurements, but the distribution may be spatially resolved [1]. Spatially resolved MRI measurements of the T2 distribution are relatively new [2,3]. Spatially resolved T2 distribution measurements are time consuming measurements in 2 and 3 spatial dimensions. The T2 distribution is most commonly determined by inverse Laplace transform (ILT) of the T2 decay. The ILT process is an illconditioned problem [4–9]. Noisy data hinders determination of the T2 distribution, therefore a better representation of the relaxation time distribution may be possible using means independently obtained [10]. This paper presents a direct way to determine the spatially resolved T2 relaxation time based on the arithmetic mean
of the relaxation rate. The arithmetic mean obtained in this way is independent of the inverse Laplace transform process. The arithmetic mean of the relaxation rate is obtained from the initial slope of the CPMG decay and is used to directly determine the harmonic mean T2. The accuracy of this harmonic mean is tested by comparing it with different relaxation time means determined from the T2 distribution. Spatially resolved measurements are undertaken with two pure phase encode MRI methods. Uniform phantoms and brine saturated porous rock core plugs are employed as test samples. Direct measurement of the mean T2 relaxation time may be advantageous for core analysis, or studies of different types of porous media. In this work the harmonic mean T2 relaxation time is used to estimate the irreducible water saturation for centrifuged core plugs. The irreducible water saturation obtained based on the T2 harmonic mean is compared with the irreducible water saturation obtained by the classic cutoff bulk volume of irreducible water (CBVI) method [11,12] and gravimetric measurements.
2. Mean T2 parameters and MRI methods ⇑ Corresponding author. Fax: +1 (506) 453 4581. E-mail address: [email protected] (B.J. Balcom). 1090-7807/$ - see front matter Ó 2013 Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.jmr.2013.11.012
The T2 CPMG echo decay experiment [13,14] is often fit to a single or multiple exponential decay model, with Eq. (1):
17
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
AðtÞ ¼ f ðtÞ ¼
n X t Ai e T 2i
ð1Þ
f ð0Þ ¼
i¼1
1 T2
The T2 distribution may be characterized by different means, for example the arithmetic mean, harmonic mean and logarithmic mean. These means are outlined below. The arithmetic mean (AM) [15,16] of the probability distribution of the relaxation times is given by Eq. (2):
ð2Þ
where T2i are the discrete values of the relaxation times and Pi(T2i) are the corresponding probabilities. The arithmetic mean is known to be sensitive to large extremes of the discrete values of the distribution [10,15]. The harmonic mean (HM) is defined [15,16] in terms of the T2 distribution by Eq. (3):
T 2HM
" #1 n n X X 1 ¼ Pi ðT 2i Þ P i ðT 2i Þ T 2i i¼1 i¼1
ð3Þ
where T2i are the discrete values of the relaxation times and Pi ðT 2i Þ T1 are the corresponding probabilities of the relaxation 2i rate. The harmonic mean is known to be sensitive to the small extremes of the discrete values of the distribution [10,15]. The logarithmic mean (LM) is defined in two ways [2,17], but in terms of the T2 distribution the logarithmic mean is calculated by Eq. (4):
T 2LM ¼ exp
Pn i¼1 ½P i lnðT 2i Þ Pn i¼1 P i ðT 2i Þ
ð4Þ
The logarithmic mean computed using Eq. (4) is identical to the geometric mean. Both are relatively insensitive to the extremes of the discrete values of the distribution [10,15]. 2.2. Arithmetic mean of the relaxation rate In this work we explore direct measurement of a related mean, the mean of the relaxation rate. Differentiating Eq. (1) which describes the CPMG decay curve, with respect to time we obtain:
f 0 ðtÞ ¼
n X 1 t exp Ai T 2i T 2i i¼1
ð5Þ
and for time t = 0 this (6) becomes:
f 0 ð0Þ ¼
n X 1 Ai T 2i i¼1
ð6Þ
f0 (0)represents the slope of the tangent line to the CPMG decay at maximum amplitude when t = 0. Eq. (6) can be rewritten:
f 0 ð0Þ ¼
Df Dt
Eq. (1) at t = 0 equals:
Experimentally f(0) represents the maximum signal amplitude of the CPMG decay. If this maximum amplitude cannot be directly acquired experimentally, it can be determined by extrapolation of the early decay. Dividing Eq. (7) by Eq. (8), in accord with Eq. (2), the arithmetic mean of the relaxation rate will be obtained:
2.1. Mean T2s
Pn i¼1 T 2i P i ðT 2i Þ P n i¼1 P i ðT 2i Þ
ð8Þ
i¼1
where A(t) is the time dependent signal observed, Ai is the amplitude of the ith component of the signal, T2i is the ith component of the relaxation time and t is the time. If there is a distribution of the relaxation time it is common to use the ILT in order to determine the T2 distribution. From the T2 distribution it is possible to determine a variety of different mean T2s.
T 2AM ¼
n X Ai ¼ Amax
ð7Þ
Df f 0 ð0Þ ¼ ¼ Dt ¼ f ð0Þ f ð0Þ
Pn
1 i¼1 Ai T 2i Pn i¼1 Ai
¼
Initial Slope Amax
ð9Þ
Based on Eq. (9) the initial slope and the maximum amplitude can be used to obtain the arithmetic mean of the relaxation rate. The brackets h i in Eq. (9) indicate the arithmetic mean of the relaxation rate. As suggested by Eq. (3) the arithmetic mean of the relaxation rate is related to the harmonic mean of the relaxation time, as Eq. (10) describes:
Pn Ai Amax T 2HM ¼ P i¼1 ¼ n 1 Initial Slope i¼1 Ai T 2i
ð10Þ
The arithmetic mean of the relaxation rate can be used to determine the harmonic mean of the relaxation time, Eq. (11):
1 T2
¼
1 T 2HM
ð11Þ
The relationship of Eq. (11) has been verified through simulation of model data sets with a variety of T2 distributions. The early decay approach has been used successfully in MR measurements of motion [18]. It has recently been employed in MRI, for velocity imaging [19,20]. 2.3. NMR and MRI pulse sequences Experimental data were acquired using bulk CMPG and three MRI techniques: CPMG-prepared Centric Scan SPRITE (Single Point Ramped Imaging with T1 Enhancement) [3], SPRITE [21] and SE-SPI (Spin-Echo Single Point Imaging) [2]. Both pure phase-encode MRI techniques were employed for one dimensional T2 mapping experiments. Fig. 1(a) presents the CPMG-prepared Centric Scan SPRITE pulse sequence and (b) presents the SE-SPI pulse sequence. The CPMGprepared Centric Scan SPRITE pulse sequence commences with a CPMG preparation followed by a 90° RF pulse for Z-storage, and a spoiler gradient. The preparation process can have a variable numbers of echoes. The last part of the pulse sequence consists of a 1D Centric Scan SPRITE image for spatial encoding. The pure phase encoding SE-SPI pulse sequence, Fig. 1(b), contains a phase encode gradient between the 90° and 180° RF pulses. The phase encoded magnetization is then read out through multiple refocusing. From each scan a single k-space point is acquired for all echoes. These k-space data sets are then Fourier transformed to generate a series of T2-weighted profiles. Based on these profiles the local CPMG decay can be extracted. The SE-SPI experiment has limitations on the duration of the first echo, therefore this experiment will not permit spatial resolution of shorter T2 components. A full CPMG decay with the CPMGprepared Centric Scan SPRITE technique is a time consuming experiment. Nevertheless the CPMG-prepared SPRITE technique [3] can be applied successfully, in a time efficient manner, if we reduce considerably the number of profiles which comprise the CPMG decay. In this case the harmonic mean T2 is measured only from the initial slope, as suggested by Eqs. (9) and (10).
18
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
Fig. 1. (a) CPMG-prepared Centric Scan SPRITE method employed to measure T2weighted profiles. These profiles are used to generate local CPMG decays. The harmonic mean T2 is determined based on the initial slope. (b) SE-SPI pulse sequence employed to measure spatially resolved T2 decay. The T2 distribution may be determined by inverse Laplace Transform and early echo decays can be used to determine the harmonic mean T2.
3. Results and discussion
Fig. 2. (a) T2-weighted profiles (1st, 20th, 40th, 70th, 100th, 130th, 210th and 320th profile) obtained with CPMG-prepared Centric Scan SPRITE measurements for H2O doped with GdCl3 using an echo time TE = 2 ms. (b) T2 values spatially resolved for the water phantom. Based on a single exponential fit of the CPMG-prepared SPRITE experiment (d) T2 varied spatially between 189 and 195 ms. The bulk CPMG value of T2 (–) 207 ms was assigned to each pixel as a reference. The harmonic mean obtained from the early decay (s) had T2 values between 224 and 244 ms. The T2 values with highest probability extracted from the T2 distribution (.) varied between 169 and 179 ms.
3.1. Phantom measurements Measurements commenced with a uniform doped water phantom, where we anticipate single exponential T2 decay. Bulk CPMG measurements were undertaken to test the T2 behavior and the data was analyzed in three ways: (i) single exponential fitting, (ii) inverse Laplace transformation and (iii) analysis of the early decay mean. Single exponential decay fitting gave a bulk T2 of 207 ms. The ILT procedure, as anticipated, gave a narrow distribution of T2s centered on 210 ms. The T2 distribution was used to calculate the means described in Section 2.1. As anticipated for a single exponential decay, the various means yielded very similar results. The arithmetic mean was 207 ms, the harmonic mean was 205 ms, and the logarithmic mean was 206 ms. Analysis of the harmonic mean T2 based on the early decay, Section 2.2, resulted in a harmonic mean T2 of 210 ms, which also closely agrees with the above results. The same phantom sample was imaged with the CPMG-prepared Centric Scan SPRITE pulse sequence. 320 T2-weighted profiles were acquired for the phantom sample. Fig. 2(a) shows select T2-weighted profiles. Fig. 2(b) shows the spatially resolved T2 values. The solid line is the bulk CPMG single exponential result of T2 = 207 ms, which was assigned as a reference value for each pixel. One data set in Fig. 2(b) is based on single exponential decay fitting, for which the mean over all pixels
was T2 = 194 ms. Another data set shows the most probable T2 determined by ILT of the CPMG-prepared Centric Scan SPRITE data, for which the mean over all pixels was T2 = 178 ms. The figure also shows the harmonic mean T2 determined based on the early slope of the CPMG decay. The mean over the displayed pixels is 235 ms. The spatially resolved mean T2 values do not agree as well as for the bulk measurement, presumably due to reduced SNR and complications of the imaging process. Analogous measurements with the SE-SPI pulse sequence show better agreement. Fig. 3(a) shows select T2-weighted profiles acquired using the SE-SPI pulse sequence. Fig. 3(b) shows the spatially resolved T2 values. Three sets of data are also reported from SE-SPI data for each pixel: (i) the single exponential decay fitting, (ii) the most probable T2 in each pixel determined by ILT and (iii) the harmonic mean T2 determined based on the early slope of the CPMG decay in each pixel. The T2 pixel mean for each data set reported are: (i) T2 = 208 ms, (ii) T2 = 209 ms and (iii) T2 = 217 ms respectively. The solid line again represents the bulk CPMG single exponential result of T2 = 207 ms, which we assign to each pixel as a reference. The simple water phantom results show that the early decay yields a high quality harmonic mean T2 consistent with the bulk measurement. CPMG-prepared SPRITE measurement and SE-SPI measurement give similar results.
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
19
Fig. 4. The T2 distribution generated by ILT of a bulk CPMG measurement of a brine saturated Berea core plug. The short lifetime extreme of the distribution is unreliable, T2 values reported are much shorter than TE. This will bias the harmonic mean T2 calculated from the T2 distribution.
Fig. 3. (a) T2-weighted profiles (1st, 20th, 40th, 70th, 100th, 130th, 210th and 512th profile) obtained with SE-SPI measurements for H2O doped with GdCl3 using an echo time TE = 2 ms. (b) T2 values spatially resolved in the water phantom. Based on a single exponential fit of the SE-SPI experiment (d) T2 varied spatially between 206 and 209 ms. The bulk CPMG value of T2 (–) 207 ms was assigned to each pixel as a reference. The harmonic mean obtained from the early decay (s) had T2 values between 211 and 222 ms. The T2 values with highest probability extracted from the T2 distribution (.) varied between 206 and 210 ms.
3.2. Uniform core plug measurements We can anticipate that the various means described in Sections 2.1 and 2.2 will not yield the same result in the case of samples with T2 distributions that are not unimodal and which extend over orders of magnitude in T2. The various means have different sensitivities to long and short T2 components in the distribution [10,15]. The T2 distribution of water saturated reservoir core plugs are of intense experimental interest and they are frequently bimodal with a T2 range of three to four orders of magnitude. Bulk and spatially resolved T2 measurements were undertaken on a brine saturated Berea rock core plug. We assume for a fully saturated core plug sample that although the sample is microscopically heterogeneous (there is a pore size distribution) it is macroscopically homogeneous (the pore size distribution is the same throughout the sample). The fully saturated Berea core plug presents a bi-modal distribution when the bulk decay curve is fit through an ILT, Fig. 4. The T2 means described in Sections 2.1 and 2.2 were determined from the bulk measurements. The logarithmic mean T2 computed from the T2 distribution was 28 ms. The arithmetic mean was 53 ms, while the harmonic mean was 10 ms. The T2 harmonic mean based on the early CPMG decay was 30 ms. The harmonic mean calculated from the bulk T2 distribution is surprisingly different from the early decay harmonic mean T2. The divergent experimental results were explored through simulations with model data sets. Simulations revealed that the T2 distribution will generate the same harmonic mean as the
corresponding T2 early decay under specific conditions imposed by Eqs. (9) and (10). Simulations showed high quality Amax and high quality early decay rate data yield the same results for the harmonic mean from the early decay and from the T2 distributions. A high quality result requires that the echo time is 5 times shorter than the shortest T2 value and noise should not be more than 10% of maximum signal. A short echo time is required in order to obtain at least 3–5 collinear points in the early decay. These conditions are satisfied for our Berea core plug data. The difference between the harmonic means computed based on the early decay and based on the T2 distribution results from the fact that one HM is computed based on the early decay and the other HM is based on the T2 distribution. The T2 distribution of Fig. 4 shows T2 values significantly shorter than the echo time (2 ms). These values do not make physical sense. The harmonic mean T2 computed based on the T2 distribution is strongly affected by these short components [10,15], generating a harmonic mean of 10 ms instead of 30 ms, the value obtained from the early decay data. These results, and simulations run for different data sets, prove that the early CPMG decay can be used successfully to compute the harmonic mean. The harmonic mean is particularly sensitive to short lifetime components from the T2 distribution. The harmonic mean directly calculated from the early decay avoids difficulties associated with the ILT procedure. Fig. 5(a) presents select T2-weighted profiles obtained for the fully saturated Berea sample with the SE-SPI pulse sequence. Based on these profiles, local T2 decays were extracted. Fig. 5(b) shows various T2 means. The harmonic mean obtained based on the early decay was T2 = 39 ms when averaged over all pixels displayed. Other means based on the T2 distribution, once more reported as means over all pixels, were arithmetic mean T2 68 ms, harmonic mean T2 18 ms and logarithmic mean T2 41 ms. Fig. 5(b) shows that with a real distribution there will be differences between different means even if these means are computed based on the same T2 distribution. The difference between these mean T2s results from the known sensitivity of these means to extremes of the distribution. Extremes of the T2 distribution are associated with large or small pores yielding long or short lifetime T2s. If long lifetime components are present in the T2 distribution then the arithmetic mean T2 will be significantly influenced by those components, while the short lifetime components will affect the harmonic mean T2 [10,15]. The logarithmic and geometric means are not affected by extremes in the T2 distribution. The new initial slope procedure used to compute the harmonic mean T2 avoids issues associated with the ILT procedure but the mean will still be strongly influenced by short lifetime T2 components.
20
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
Fig. 5. (a) The T2-weighted profiles (1st, 4th, 8th, 16th, 32nd, 56th, 80th, 128th and 256th profiles) obtained with SE-SPI measurements of a Berea core plug fully saturated with brine using an echo time TE = 2 ms. (b) Harmonic mean T2 obtained from the initial CPMG decay (d) with spatial variation between 37 and 41 ms. The arithmetic mean computed from the T2 distribution (r) with T2 varying between 66 and 70 ms. The harmonic mean computed from the T2 distribution (j) with T2 varying between 16 and 19 ms. The logarithmic mean computed from the T2 distribution (N) with T2 varying between 39 and 43 ms.
Different T2 means have been related to petro-physical parameters in the literature. The geometric mean T2 is a good parameter for predicting sandstone permeability [20] and irreducible watersaturation in sandstones [10,20]. Sandstone permeability and hydrocarbon viscosity have been correlated with the logarithmic mean T2. The pore size distribution can be related to the width of the T2 distribution, while in hydrocarbon studies the mean chain length is related to the relaxation times [5,6,22]. 4. Centrifuged core plugs measurements 4.1. Heterogeneous core plugs measurements Having explored various T2 means and the early decay harmonic mean, we turn to macroscopically inhomogeneous samples where we anticipate spatially resolved differences in the T2 distribution. Based on previous work [2] we consider centrifuged core plugs where spatial variation of saturation is associated with differences in the T2 distribution. For this class of sample, the CPMG-prepared SPRITE technique was used. Fig. 6(a) shows the first 10 T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a brine saturated Berea core plug. The early decay harmonic mean T2 was determined for different positions through the sample, in Fig. 6(b). The early decay harmonic mean T2 varies from 24 ms to 38 ms within the field of view. The harmonic mean T2 was computed based on the first 10 collinear decay points, in order to ensure a high quality result.
Fig. 6. (a) The first 10 T2-weighted profiles of the Berea sandstone, acquired using the CPMG-prepared SPRITE technique after centrifugation of the core plug. (b) The harmonic mean T2 obtained based on the early slope of the CPMG extracted from the T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a centrifuged Berea core plug. The harmonic mean T2 varies spatially from 24 ms to 38 ms for positions between 1.6 cm and 4.7 cm.
Fig. 6(b) shows less variation of the mean T2 compared with previous results obtained by Li [2]. In Fig. 6(b) the early decay harmonic mean T2 is computed, while Li computed the logarithmic mean T2. The shallow slope of Fig. 6(b) is due to the relative insensitivity of the early decay harmonic mean T2 to changes produced by emptying large pores which are associated with long T2 components. A brine saturated Buff Berea core plug was also examined. For this sample the first 10 T2-weighted profiles were measured with a short echo time of TE = 500 ls, as can be observed in Fig. 7(a). A short echo time measurement results in less visible profile attenuation but the early decay mean is still readily calculated. Fig. 7(b) shows the harmonic mean plotted for different positions through the sample obtained based on the T2-weighted profiles. The harmonic mean of the relaxation time shown in Fig. 7(b) varies spatially as was anticipated, for positions between 2.8 and 5.8 cm. The harmonic mean T2 varied from 10 ms to 30 ms. Data at the edge of the core plug was not analyzed, Fig. 7(a) shows an obvious artefact at the object left edge.
4.2. Harmonic mean and irreducible water saturation determination The irreducible water saturation is a critical parameter for core plug analysis. It has been shown that the irreducible water saturation, Swirr, can be obtained using mean T2 values [2,12]. A linear relationship between the logarithmic mean T2 and the residual water saturation has been established [2,12]. We now consider if the harmonic mean is similarly related to local saturation.
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
21
Fig. 8. The harmonic mean T2, computed based on the early CPMG decay using the CPMG-prepared SPRITE technique, versus saturation of a centrifuged Buff Berea core plug. The irreducible water saturation was estimated to be 21% based on the analysis of the data with Eq. (12).
Fig. 7. (a) The first 10 T2-weighted profiles of the Buff Berea sandstone, acquired using the CPMG-prepared SPRITE technique. (b) The harmonic mean T2 obtained based on the slope of the CPMG extracted from the T2-weighted profiles acquired using the CPMG-prepared SPRITE technique for a centrifuged brine saturated Buff Berea core plug. The harmonic mean T2 varies spatially from 10 ms to 30 ms for positions between 2.8 cm and 5.8 cm. The edges of the sample were not measured.
Li et al. [2] proposed based on Coates et al. [11,12] a linear relationship between the T2LM and residual water saturation, as described by the equation:
SWðT 2LM Þ ¼ m T 2LM Swirr þ Swirr
ð12Þ
where m is a parameter which describes the pore geometry, SWðT 2LM Þ is the residual water saturation, T2LM is the logarithmic mean T2 associated with the ith position of the sample, and Swirr is the irreducible water saturation. In this paper we seek to determine if the relation of Eq. (12) remains true when we employ a different mean T2. In particular we wish to employ the T2HM determined from the early echo decay. Knowing the T2HM for different saturations of the sample, we obtained the irreducible water saturation by plotting harmonic mean T2 versus the residual water saturation. The local saturation is determined based on simple 1D DHK Centric Scan SPRITE profiles of fully saturated and centrifuged core plugs. The spin density measured with SPRITE pulse sequence for a fully brine saturated core plug represents SW = 1 (100% saturation). Measuring the spin density for a centrifuged core plug and keeping constant the experimental parameters for SPRITE pulse sequence, the variation of saturation across the sample was established. Fig. 8 shows the harmonic mean T2 obtained using CPMGprepared SPRITE measurement for a centrifuged Buff Berea core plug. Fitting the decay of water saturation versus the T2 mean with Eq. (12) it was determined that the irreducible water saturation for Buff Berea was 21%. Using the classical cutoff BVI techinique [12], the irreducible water saturation was estimated to be 25%. From gravimetric measurements an irreducible water saturation of 22%
was obtained. These results show that for Buff Berea, or for core plugs which present similar pore structure, the harmonic mean T2 of centrifuged samples could be used to determine the irreducible water saturation. The irreducible water saturation was also calculated for the Berea core plug. Using Eq. (12) it was determined that the irreducible water saturation for the centrifuged Berea was 20%. The irreducible water saturation was computed considering only the positions between 2.1 cm and 4.6 cm, neglecting the edges of the sample. The result in this case is not very reliable because of the scatter in the T2 mean for the short lifetime component. From gravimetric measurements the irreducible water saturation was known to be 18%. The conventional T2 cutoff value of 33 ms used in classical cutoff BVI technique for sandstone clearly gives an erroneous irreducible water saturation of 50%.
5. Experimental The experiments were performed on an 8.51 MHz Oxford Maran DRX HF instrument (Oxford Instruments, Witney, UK). All measurements were carried out at ambient magnet temperature of 25 °C. The T2 distribution was determined employing WinDXP software (Oxford Instruments, Witney, UK). The profiles obtained from CPMG-prepared SPRITE experiments were computed employing various routines written with Interactive Data Language (Exelis, Boulder, USA). Matlab Math Works Natick, MA, USA was employed for simulation. SE-SPI measurements were applied to a uniform brine saturated Berea. Trapezoidal phase encode gradients were employed with ramp times of 150 ls and a plateau time of 100 ls. The time between the first two RF pulses was 2 ms and 1024 echoes were acquired. In order to increase the signal-to-noise ratio, 20 points from each echo were acquired. The field of view (FOV) was 9.6 cm with two scans acquired. To maintain the phase shift upon refocusing, an XY-16 phase cycle was applied [23] and to compensate RF field heterogeneity, composite pulses were employed. The CPMG-prepared Centric Scan SPRITE measurements were undertaken with a phantom sample, doped with GdCl3 and brine saturated core plugs. GdCl3 (0.1 mM) was added to reduce T1 to 220 ms. For the phantom sample and brine saturated Berea the CPMG-prepared Centric Scan SPRITE technique had an encoding time of 400 ls, FOV = 60 mm, NS = 32, 4-step phase cycle, gradient spoil duration of 500 ls at 2 G/cm, P90 = 32 ls, P180 = 64 ls and 64 k-space points were acquired. The flip angle was 4°, the sweep
22
R.E. Nechifor et al. / Journal of Magnetic Resonance 239 (2014) 16–22
width was 50 kHz, the probe ring-down time was 5 ls, and the filter stabilization time was 284 ls. The echo time was TE = 2 ms and 1024 echoes were acquired. The CPMG-prepared Centric Scan SPRITE experiment employs longer echo time than necessary in order to permit a direct comparison with SE-SPI experiments. For the phantom sample, 320 profiles were acquired in 6 h for the same TE and different numbers of echoes in order to record a full CPMG decay. For the brine saturated Berea core plug the flip angle was 9°, and only 10 profiles were acquired. For each profile 128 averages were acquired in 8 min. For the brine saturated Buff Berea core plug the CPMG-prepared Centric Scan SPRITE technique had an encoding time of 150 ls, FOV = 90 mm, NS = 512, 4-step phase cycle, gradient spoil duration of 500 ls at 2 G/cm, P90 = 43 ls, P180 = 86 ls and 64 k-space points were acquired. The flip angle was 9°, the sweep width was 125 kHz, the probe ring-down time was 30 ls, and the filter stabilization time was 26 ls with an echo time of TE = 500 ls. Only 10 profiles were acquired and for each profile 512 averages were acquired in 22 min. The gradient stabilization time was 2 ms in all cases. Two different core plugs were used: Berea and Buff Berea (Kocurek Industries, Caldwell, TX, USA). Both were cylindrical core plugs, 3.7 cm in diameter and 5 cm in length. The porosities were 21% for Berea and 22% for Buff Berea as determined gravimetrically. The core plugs were saturated with brine (1.5% wt. NaCl in H2O) under vacuum for more than 48 h. To ensure partial saturation the samples were centrifuged for 60 min at 2000 rpm at room temperature. A Beckman Model J2-12M centrifuge was employed. The variation of saturation across the samples was computed employing Centric Scan SPRITE measurements of fully and partially saturated samples. The SPRITE pulse sequence parameters were the same as the parameters used for CPMG-prepared SPRITE. 6. Conclusions The harmonic mean T2 obtained from the initial slope of the CPMG decay is a powerful way to characterize the mean T2. The accuracy of the harmonic mean T2 computed based on the early CPMG decay is tested in this paper in various ways: through simulations, on phantom sample and on uniform and partial saturated core plugs. The harmonic mean was also compared with different means computed based on the T2 distribution. The experimental data were obtained using bulk CPMG technique, CPMG-prepared Centric Scan SPRITE technique and SE-SPI technique [2,24–27]. The advantage of the CPMG-prepared Centric Scan SPRITE technique is that it can have a short echo time consistent with bulk CPMG measurements. Although spatial encoding has been explored with phase encoding, frequency encoding with short echo time is also possible. The harmonic mean computed from the early decay does not require intermediate data manipulation, as with the inverse Laplace Transform approach. Therefore the harmonic mean is not affected by different parameters used in processing the data. Regardless of the method chosen to determine the harmonic mean, it was shown that it will still be biased by the short lifetime T2 component. Using this new technique the harmonic mean of the relaxation time or the arithmetic mean of the relaxation rate can be obtained. The early decay harmonic mean T2 can be used to study the variation of saturation for porous media and the harmonic mean T2 based on early decay for different positions in the sample can be obtained. The harmonic mean T2 can be used also in other studies, such as studies of fluid flow process through the sample.
Acknowledgments B.J.B. thanks NSERC for a Discovery grant and an Idea to Innovation grant. B.J.B. thanks the Canada Chairs program for a Research Chair in MRI of Materials. REN thanks R. MacGregor for his assistance with the experiments and B. MacMillan for helpful discussions. References [1] R. Kimmich, NMR Tomography, Diffusometry, Relaxometry, Springer-Verlag, Berlin, 1997. [2] L. Li, H. Han, B.J. Balcom, Spin echo SPI methods for quantitative analysis of fluids in porous media, J. Magn. Reson. 198 (2009) 252–260. [3] O.V. Petrov, G. Ersland, B.J. Balcom, T2 distribution mapping profiles with phase-encode MRI, J. Magn. Reson. 209 (2011) 39–46. [4] L. Venkataramanan, Y.-Q. Song, M.D. Hürlimann, Solving Fredholm integrals of the first kind with tensor product structure in 2 and 2.5 dimensions, IEEE Trans. Signal Process. 50 (2002) 1017–1026. [5] L. Venkataramanan, F.K. Gruber, T.M. Habashy, D.E. Freed, Mellin transform of CPMG data, J. Magn. Reson. 206 (2010) 20–31. [6] L. Venkataramanan, F.K. Gruber, T.M. Habashy, D.E. Freed, Continuous moment estimation of CPMG data using Mellin transform, J. Magn. Reson. 216 (2012) 43–52. [7] M. Simina, R. Nechifor, I. Ardelean, Saturation-dependent nuclear magnetic resonance relaxation of fluids confined inside porous media with micrometersized pores, Magn. Reson. Chem. 49 (2011) 314–319. [8] J.G. McWhirter, E.R. Pike, On the numerical inversion of the Laplace transform and similar Fredholm integrals of the first kind, J. Phys. A: Math. Gen. 11 (1978) 1729–1745. [9] C.L. Epstein, J. Schotland, The bad truth about Laplace’s transform, SIAM Rev. 50 (2008) 504–520. [10] G.C. Borgia, R.J.S. Brown, P. Fantazzini, Different ‘‘average’’ nuclear magnetic resonance relaxation times for correlation with fluid-flow permeability and irreducible water saturation in water-saturated sandstone, J. Appl. Phys. 82 (9) (1997) 4197–4204. [11] G.R. Coates, D. Marschall, D. Mardon, J. Galford, A new characterization of bulkvolume irreducible using magnetic resonance, Log Anal. 39 (1998) 51–63. [12] G.R. Coates, L. Xiao, M.G. Prammer, NMR Logging – Principles and Applications, Halliburton Energy Service, Houston, 1999. [13] S. Meiboom, D. Gill, Modified spin echo method for measuring nuclear relaxation times, Rev. Sci. Instrum. 29 (1958) 688–691. [14] P.T. Callaghan, Principles of Nuclear Magnetic Resonance Microscopy, Clarendon Press, Oxford, 1991. pp. 339–340. [15] G.H. Hardy, J.E. Littlewood, G. Pólya, Inequalities, Cambridge University Press, London, 1934. [16] M. Abramowitz, I.A. Stegun, Handbook of Mathematical Functions with Formulas, Graphs and Mathematical Tables, National Bureau of Standards, Washington, DC, 1964. [17] K.B. Stolarsky, Generalizations of the logarithmic mean, Math. Mag. 48 (2) (1975) 87–92. [18] P. Dyverfeldt, A. Sigfridsson, H. Knutsson, T. Ebbers, A novel MRI framework for the quantification of any moment of arbitrary velocity distributions, Magn. Reson. Med. 65 (2011) 725–731. [19] K.V. Romanenko, B.J. Balcom, Permeability mapping in porous media by magnetization prepared centric-scan SPRITE, Exp. Fluids 50 (2011) 301–312. [20] G.C. Borgia, R.J.S. Brown, P. Fantazzini, Estimates of permeability and irreducible water saturation by means of a new robust computation of fractional power average relaxation times, Magn. Reson. Imaging 16 (1998) 613–615. [21] B.J. Balcom, R.P. MacGregor, S.D. Beyea, D.P. Green, R.L. Armstrong, T.W. Bremner, Single point ramped imaging with T1 enhancement (SPRITE), J. Magn. Reson. A 123 (1996) 131–134. [22] D.E. Freed, Dependence on chain length of NMR relaxation times in mixtures of alkanes, J. Chem. Phys. 126 (2007) 174502-1–174502-10. [23] T. Gullion, The effect of amplitude imbalance on compensated Carr-Purcell sequences, J. Magn. Reson. Ser. A 101 (1993) 320–323. [24] B.J. Balcom, SPRITE imaging of short relaxation time nuclei, in: B. Blumich, W. Kuhn (Eds.), Spatially Resolved Magnetic Resonance: Methods, Materials, Medicine, Biology, Rheology, Geology, Ecology, Hardware, Wiley-VCH, Weinheim, 1998, pp. 75–86. [25] A.V. Ouriadov, R.P. MacGregor, B.J. Balcom, Thin film MRI – high resolution depth imaging with a local surface coil and spin echo SPI, J. Magn. Reson. 169 (2004) 174–186. [26] S. Gravina, D.G. Cory, Sensitivity and resolution of constant-time imaging, J. Magn. Reson. Ser. B 104 (1994) 53–61. [27] O.V. Petrov, B.J. Balcom, Two-dimensional T2 distribution mapping in porous solids with phase encode MRI, J. Magn. Reson. 212 (2011) 102–108.
