Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

The characterization of porous rocks with nuclear magnetic resonance and the confocal laser scanning.. 1992

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
ubc_1992_spring_chapman_alice.pdf
ubc_1992_spring_chapman_alice.pdf [ 6.11MB ]
Metadata
JSON: 1.0052496.json
JSON-LD: 1.0052496+ld.json
RDF/XML (Pretty): 1.0052496.xml
RDF/JSON: 1.0052496+rdf.json
Turtle: 1.0052496+rdf-turtle.txt
N-Triples: 1.0052496+rdf-ntriples.txt
Citation
1.0052496.ris

Full Text

THE CHARACTERIZATION OF POROUS ROCKS WITH NUCLEAR MAGNETIC RESONANCE AND THE CONFOCAL LASER SCANNING MICROSCOPE by ALICE ELIZABETH CHAPMAN B.Sc., The University of British Columbia, 1989 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE in THE FACULTY OF GRADUATE STUDIES (Department of Geological Sciences) We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA September, 1992 © Alice Elizabeth Chapman, 1991 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of ;77-L() (-7/C JfrL The University of British Columbia Vancouver, Canada Date Cc /3,, /9gQ DE-6 (2/88) ABSTRACT The characterization of porous media with the confocal laser scanning microscope (CLSM) and nuclear magnetic resonance (NMR) was investigated. Pore casts (epoxy replicas) of one sandstone, Berea 100, and two carbonates, S55 and S77, were imaged with the CLSM. Isolated two-dimensional images and serial sections, a series of two-dimensional images separated by a specific vertical distance, were collected. The isolated images were analysed and mean pore radii of 25.3, 17.5, and 90.2p.m were determined for Berea 100, S55 and S77 respectively. The serial sections were reconstructed to form three-dimensional images of the pore networks. However, the pores in the samples were large compared to the depth of the pore casts, resulting in images with insufficient depth to allow three-dimensional pore size analysis. Spin-lattice (T1)measurements of twenty-two saturated carbonates (samples S31 to S49, S54, S55, S77) and six partially saturated carbonates were analysed to determine pore sizes and fluid distributions respectively. Pore sizes were calculated using the “two-fractions in fast exchange” model and the slow diffusion model of Brownstein and Tarr (1979). Samples S31 to S49, with pore radii of less than 41.Lm, were within the fast diffusion regime and thus the “two- fractions in fast exchange model” produced reliable pore sizes. Samples S54, S55 and S77, with mean pore radii of 5.1, 19.3, and 85.3.tm, were not within the fast diffusion regime and so the accurate analysis of data for these samples required the slow diffusion model. The spin-lattice time constants calculated for the partially saturated samples were assumed to be proportional to the sizes of the saturated pores. Thus it was concluded that the larger pores within the samples were drained preferentially; however, large pores which were accessible through smaller pore throats remained saturated to low saturations. Further investigation into these methods could focus on the production of pore casts to enable reconstruction of an image with sufficient depth for three-dimensional analysis or on the development of a pore model which more closely approximate the shape of the pores. II TABLE OF CONTENTS ABSTRACT .ii TABLE OF CONTENTS m LIST OF TABLES vi LIST OF FIGURES vii ACKNOWLEDGEMENTS xii CHAPTER ONE Introduction 1 CHAPTER TWO Characterization of Porous Rocks with the Confocal Laser Scanning Microscope 4 INTRODUCrION 4 Characterizing Porous Media 4 Mercury Injection 5 Gas AdsorptionfDesorption 5 Electrical and Optical Microscopy 6 THE CONFOCAL LASER SCANNING MICROSCOPE (CLSM) 9 EXPERIMENTAL PROCEDURES 11 Sample Preparation 11 Materials 12 Pore Cast Production 14 Impregnation Techniques 15 Wet Impregnation 15 Vacuum Impregnation Immersion Technique 16 Vacuum Impregnation Lowering Technique 16 Etching 17 Double Pore Casts 18 DATA COLLECfION AND ANALYSIS 18 RESULTS 24 DISCUSSION AND CONCLUSIONS 48 I I I CHAPTER THREE Characterizing Porous Rocks with Nuclear Magnetic Resonance 54 INTRODUCTION 54 HISTORY OF NMR IN PETROLEUM GEOLOGY 55 Nuclear Magnetic Logging (NML) 55 Relationships Between NMR and Physical Properties 56 Pore Sizes 56 Fluid Distributions 57 THEORY 57 Longitudinal Relaxation Time (Ti) 60 Transverse Relaxation Time (T2) 61 The NMR Response of Fluids Confmed in Rocks 62 Slow Diffusion 63 Fast Diffusion 65 EXPERIMENTAL PROCEDURES 67 Sample Preparation 67 Sample Characterization 68 Porosity and Permeability 68 Mineralogy 69 Samples S31 to S49 72 Samples S54,S55, andS77 72 DATA COLLECHON 76 Sample Preparation 76 Fully Saturated Samples 76 Partially Saturated Samples 77 NMR Measurements 78 Spin-Lattice Decay Curves 78 Spin-Spin Decay Curves 79 DATA ANALYSIS 79 Data Inversion 80 Spin-Lattice Relaxation Time (Ti) Distributions 80 Pore Sizes - Slow Diffusion Regime 82 Porosity Determination 83 Pore Sizes - Fast Diffusion Regime 84 RESULTS 84 Fully Saturated Samples 84 Ti Distributions 84 Porosity 85 Pore Sizes 90 Fast Diffusion 90 Slow Diffusion 92 iv Partially Saturated Samples .102 Fluid Distributions 102 Fluid Saturations 102 DISCUSSION AND CONCLUSIONS 108 CHAPTER FOUR Summary 113 INTRODUCTION 113 CONFOCAL MICROSCOPY 113 NMR 116 Determination of Pore Sizes 116 Fluid Distributions 119 CONCLUSIONS 120 REFERENCES 121 APPENDIX 1 BioRad MRC-500 Image File Format 127 APPENDIX 2 Image Analysis Programs 129 APPENDIX 3 Berea 100 Confocal Pore Size Distributions 165 APPENDIX 4 S55 Confocal Pore Size Distributions 174 APPENDIX 5 S77 Confocal Pore Size Distributions 185 APPENDIX 6 Carbonate Samples NMR Ti and Pore Size Distributions 193 APPENDIX 7 NMR Ti distributions and M(0) vs. Sw graphs for Partially Saturated carbonate samples 238 V LIST OF TABLES 2.1 Spurr Low-Viscosity Embedding Media 13 2.2 Berea 100 Mean Pore Sizes 47 2.3 Carbonate (S55) Mean Pore Sizes 47 2.4 Carbonate (S77) Mean Pore Sizes 48 2.5 Mean Pore Sizes 50 3.1 Porosity and Permeability 70 3.2 Stains for Differentiating Between Calcite and Dolomite 71 vi LIST OF FIGURES 2.1.20 Pixel intensities in rows 40 and 250 of image ABR1. The intensity value 50 is indicated as possible separation value between pore pixels and solid pixels 2.2a 21 Binary versions of Berea 100 image ABR1 with pore/solid separation value of 40. Black = solid, and white = pore space. 2.2b 22 Binary versions of Berea 100 image ABR1 with pore/solid separation value of 50. Black = solid, and white pore space. 2.2c 23 Binary versions of Berea 100 image ABR1 with pore/solid separation value of 60. Black = solid, and white = pore space. 2.3 25 Original version of a confocal laser scanning microscope image of Berea 100. The image is called ABR1 and the pore space corresponds to the brightest areas of the image. 2.4a 26 CLSM image of vugs in the carbonate sample, S55. The image is called F78 and the brightest areas of the image represent porosity. 2.4b 27 CLSM image of vugs in the carbonate sample, S55. The image is called T78 and the brightest areas of the image represent porosity. 2.5 29 Original version of a confocal laser scanning microscope image of carbonate sample S77. The image is called Fl and the pore space corresponds to the brightest areas of the image. 2.6a 31 Binary version of carbonate S55 CLSM image F78 (Figure 2.4a) produced using a solid/pore separation value of 50. White = porosity, and Black = solid. 2.6b 32 Binary version of carbonate S55 CLSM image F78 (Figure 2.4a) produced using a solid/pore separation value of 100. White = porosity, and Black = solid. 2.7 33 Binary version of carbonate S77 CLSM image Fl (Figure 2.5) produced using a solid/pore separation value of 15. White = porosity, and Black = solid. vii 2.8.35 Binary versions of Berea 100 image ABR1 with pore throats separated from the pores and the centers from which the sizes of the pores and pore throats will be measured indicated by a contrasting colour within the void space. Black = solid, the White = porosity, and the dotted areas = pore throats. 2.9a 37 Examples of the Berea 100 pore size distributions. The upper plot represents the pore size distribution of image ABR 1, the lower of image HBR 1. 2.9b 38 Berea 100 cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 16 images. 2.lOa 39 Examples of the carbonate (S55) pore size distributions. The upper plot represents the pore size distribution of image F78, the lower of image N78. 2.lOb 40 Carbonate (S55) cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 20 images. 2.lla 42 Examples of the carbonate (S77) pore size distributions. The upper plot represents the pore size distribution of image Gi, the lower of image Ji. 2.llb 43 Carbonate (S77) cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 13 images. 2.12 44 Porosity (the ratio of pore pixels to the total number of pixels) of the Berea 100 images. 2.13 45 Porosity (the ratio of pore pixels to the total number of pixels) of the carbonate S55 images. The upper distribution represents the binary images produced using a low intensity separation value, and the lower distribution represents those images produced with a high intensity separation value. 2.14 46 Porosity (the ratio of pore pixels to the total number of pixels) of the carbonate S77 images. 2.15 51 Berea 300 capillary pressure pore size distributions. (from Crocker et al., 1983 viii 2.16.53 Consecutive optical sections of Berea 100. The spacing was 9.6 microns between the sections. From right to left, the images are the second, third and final (fourth) images collected. ABR1 was the first of this set (Figure 2.4). The bright areas in the images are the pores, and the darker areas represent the solid. 3.1 59 Nuclear Magnetic Resonance Reference Frame 3.2 73 Partial dolomitization in Elkton Member samples, S31 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. 3.3 74 Examples of crinoid stems and ossicles in Elkton Member samples, S31 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. Lower photo contains a cross section of a crinoid stem 3.4 75 Examples of bryozoan and brachiopods which occur in Elkton Member samples, S31 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. 3.5 86 Sample S39 Ti distributions.(a) discrete and (b) continuous. 3.6 87 Sample S31 Ti spectra.(a) discrete and (b) continuous. 3.7 88 Sample S77 Ti spectra. (a) discrete and (b) continuous, interpolate a coefficient of 2.25x10-9 m2/s at 40°C, the temperature of the sample chamber in the NMR spectrometer used to collect the data for this study. 3.8 89 Correlation of the integral over the the discrete and continuous Ti distributions (M(0)) with measured porosity. 3.9 91 100,000 ppm NaCl brine relaxation time. 3.10 93 (a) Sample S39 Ti and fast diffusion pore size distributions. (b) S3 i fast diffusion pore size distribution. 3.11 94 (a) Sample S77 Ti and fast diffusion pore size distributions. (b) S55 fast diffusion pore size distribution. ix 3.12.95 Slow diffusion method results of inverting theoretical data produced using a theoretical pore size spectrum similar to that of S77. Variation in the pore sizes is a result of using two time ranges to calculate the theoretical decay curve. The upper graph was produced with data over the same time range as the data for S77 was collected, the lower with an extended range. 3.13 96 Slow diffusion method results of inverting theoretical data based on a theoretical pore size spectrum similar to that of S77. The pore sizes resulting when 1 and 5% noise are added to the theoretical data are shown in the upper and lower graphs respectively. 3.14 99 Slow diffusion pore size distributions. (a) S43 (b) S31 3.15 100 Slow diffusion pore size distributions. (a) S55 (b) S77 3.16 101 Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the bulk fluid relaxation time, Tib, is varied. 3.17 103 Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the self diffusion coefficent of the bulk fluid, D, is varied. 3.18 104 Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the surface sink parameter, p, is varied. 3.19 105 Example of the decrease in amplitude and relaxation time with saturation for partially saturated samples. 3.20 106 Example of the increase of the amplitude of higher relaxation times at low saturations. 3.21 107 Examples of the correlation between fluid saturation and the total integral over the Ti distribution (M(0)). 3.22 iii Comparisons of the pore size spectra calculated from the NMR data with the fast and slow diffusion methods with the spectra determined using the confocal laser scanning microscope. x 4.1.114 Berea 100 pore size distributions. The upper distribution is a capillary pressure pore size distribution (after Crocker et aL, 1983) and the lower is the distribution determined from 16 confocal laser scanning microscope images. 4.2 115 Carbonate (S55 and S77) cumulative pore size distributions. The upper plot represents the summation of the distributions of 20 S55 images, the lower of 13 S77 images. 4.3 118 Ti and fast diffusion pore size disthbutions. (a) Sample S39 (b) Sample S31 xi ACKNOWLEDGEMENTS I would like to express my appreciation of the enthusiasm and support provided by my supervisor, Dr. Rosemary Knight. I would also like to thank Michael Knoll, Paulettte Tercier, Ana Abad, Dave Goertz and Dave Butler for their comments and assistance throughout In addition, I would like to thank Ken Whittall for his assistance in analysing the data. This project ran considerably smoother due to advice and assistance provided by the technicians, in particular Yvonne Dumas and Ray Rodway. Funding for this project was provided by Shell Canada Limited, Esso Resources Canada Limited, Petro-Canada Resources, and the Natural Sciences and Engineering Research Council of Canada. xii CHAPTER ONE Introduction The physical properties of porous rocks (i.e. elastic wave velocities, electrical resistivity, dielectric constant) are all governed, to varying extents, by the microgeometry of the pores and pore fluids. Empirical laws are commonly used to link physical properties to characteristics of these components; however, such laws are generally quite limited in application. In order to determine more fundamental relationships between physical properties, the pore space, and fluid distributions, analytical methods are required which can provide accurate characterization of the pore space and fluid distribution of rocks. The methods currently available for characterizing porous media include mercury injection porosimetry, gas adsorption/desorption methods, electron and optical microscopy, and nuclear magnetic resonance (NMR). Mercury injection porosimetry and gas adsorption/desorption methods operate on the same basic principle: fluids which are not attracted to the surface of rocks must be forced into a rock by increasing the pressure on the fluid. The size of the pores or pore throats through which the fluid enters or escapes from the rock is related to the fluid pressure in the system. The results of these two methods are generally presented as plots of the volume fraction of the total void space which is accessible through specific of pore throat sizes. The most accurate characterization of porous media is accomplished using microscopy to image the pore networks. Both electrical and optical microscopy are used to image the surfaces of thin or polished sections of rocks. Pore sizes may be determined from microscope images either through statistical extrapolation of three-dimensional pore sizes from two-dimensional images or direct measurement of pore sizes from a three-dimensional image of the pore network. High quality, two-dimensional images have been collected for many years; however, pore networks are inherently three-dimensional. Thus one of the goals of this thesis was to attempt to collect two- 1 dimensional images which could be used to reconstruct a three-dimensional image of a pore network. Three-dimensional data have been obtained through the collection of stereo pairs with the scanning electron microscope (SEM) and optical microscopes (Pitmann and Duschatko, 1970). In addition, three-dimensional images have been reconstructed from vertical series of two- dimensional images, or serial sections, produced through physical sectioning of the media (Straley and Minnis, 1983, and Lin and Hamasaki, 1983). The method used in this thesis avoids many of the problems associated with stereo pair analysis and physical sectioning of rocks by collecting a series of two-dimensional images of a pore cast through optical sectioning of the sample with the Confocal Laser Scanning Microscope (CLSM). The CLSM was used to collect independent two-dimensional images as well as vertical series of two-dimensional images of three samples: Berea 100, a quarry sandstone, and two carbonates, samples S55 and S77. The independent two-dimensional images were collected to determine the pore sizes of the samples. The serial sections were collected to determine if the CLSM could be used to produce three-dimensional images of a pore network. Nuclear Magnetic Resonance was the second method used to characterize porous media in this thesis. Advantages of NMR include speed, non-destructiveness, and the fact that samples may be wet at the time of the measurement. NMR may be used to determine pore sizes because the spin-lattice or longitudinal (T1) relaxation times are dependent on the geometry of the pore space. The relation between pore sizes and T1 relaxation times is based on the increase in the rate of relaxation in fluids near a surface compared to the relaxation rate in bulk fluids. This results in the spin-lattic relaxation of fluids in porous media being a function of the surface to volume ratio of the pores (SfVre). Prior investigations of NMR have focussed primarily on determining the porosity, permeability, and pore sizes of sandstones. The determination of fluid distributions and the investigation of pores in carbonates, which are irregularly shaped and often very large, have rarely been considered. Thus, the other objective of this thesis was to test the current theory of NMR on 2 carbonates, and to determine what information on fluid distributions may be obtained using this technique. NMR data on twenty-two carbonate samples, nineteen of which were from the same formation, was collected and analysed using two approaches: a fast diffusion method and a slow diffusion method. The fast diffusion method is theoretically applicable only to those samples which contain very small pores (<20 microns) while the slow diffusion method is applicable to samples containing pores of any size. The nineteen samples from the Mississippian Elkton Member of the Turner Valley Formation had very little visible porosity, and therefore it was expected that these samples would be within the fast diffusion regime. However, the porosity of the three samples from the Mississippian Turner Valley Formation, S54, S55 and S77, was dominated by visible porosity and thus the accurate analysis of the NMR data for these samples was expected to be possible only with the slow diffusion method. Determining fluid distributions and the pore sizes of carbonates from NMR data has rarely been attempted and the CLSM has never been used to image pore sizes in rocks. However, the advantages associated with these methods may result in more accurate characterization of the pore space and fluid distributions in rocks than has previously been possible. 3 Characterization of Porous Rocks with the Confocal Laser Scanning Microscope INTRODUCTION The physical properties of porous rocks (eg. elastic wave velocities, electrical resistivity, dielectric constant) are all governed, to a varying extent, by the microgeomeiry of the pore space. Empirical laws are commonly used to link physical properties to characteristics of the pore space, but such laws are generally quite limited in application. In order to determine more fundamental relationships between physical properties and the pore space, analytical methods are needed that can provide accurate characterization of the pore space in rocks. The objective of this study was to develop a means of obtaining and analyzing three-dimensional digitized images of the pore space using a confocal laser scanning microscope. Characterizing Porous Media The current methods used to determine pore network characteristics in geologic samples include mercury injection, gas adsorptionldesorption, optical microscopy, and electron microscopy. Each of these techniques, and their associated problems, is briefly described below. 4 Mercury Injection The mercury injection method (Wardlaw, 1976) produces a size distribution which is more representative of the pore throats than of the pores. The general procedure for mercury injection is to remove all fluids from the sample, and then, with the sample under vacuum, force the mercury into the evacuated sample. The mercury is forced into the sample by progressively increasing the pressure on the system. The volume of mercury entering the pores and the corresponding pressure increases are recorded and plotted as a pressure-volume or mercury injection capillary-pressure curve. This curve is used to calculate pore throat size distributions based on the assumptions that the pressure required to force mercury into a pore is related to the size of the pore throat. Thus, the volume of mercury which enters the sample between specific pressure limits is related to that fraction of the total pore volume which is connected by pore throats within specific size limits. The pore throats are frequently assumed to be tubular and the resulting equation relating pressure to the pore size is: D _4acos9 IC— where Pc is the pressure, a is the surface tension of the liquid (mercury), e is the contact angle of the liquid on the solid, and d is the diameter of the tube, or pore throat (Wardilaw, 1976). Gas Adsorption/Desorption The gas adsorption/desorption method described by Pierce (1953) is similar to the mercury injection method in that the pressure on the system is altered to either push the gas into the sample or allow the gas to escape. An adsorption or desorption isotherm of pressure versus volume of gas adsorbed or desorbed is plotted and one of the isotherms is used to calculate the pore sizes of the 5 sample. The isotherm selected to calculate the pore sizes is always the one which represents capillary equilibrium in the sample. This is checked by comparing the surface area which can be calculated from the isotherm with the surface area calculated using the BET method (Brunauer et al., 1938). For either isotherm, the volume at the end of the pressure step is assumed to be related to the size of the pore which is filling or emptying, and, using the Kelvin radius computation, the radius of the core within the pores which is not filled (rk) may be determined: r= 4.14 log(j) where P0 is the initial pressure and p is the current pressure. The film of gas on the walls of the pores is assumed to contain “n” statistical layers and thus the pore radius (rn), is equal to: r = rk+3.6n where 3.6 A is the thickness of an adsorbed layer of nitrogen. Electrical and Optical Microsconv The most accurate characterization of pore networks in three-dimensions is based, as all pore characterization was in the early twentieth century, on the use of microscopy to image the networks. Both electrical and optical microscopy are used to image the surfaces of thin or polished sections of rocks. Currently, most researchers are focussing on two objectives: producing high quality images which can easily be separated into the pore and grain components, and collecting and analyzing three-dimensional data. High quality, clear images are often produced by injecting a material into the pore space which will enhance the contrast between pore space and grains in the rock. Conventionally, rock 6 samples in which the voids are to be examined were injected with a blue-dyed epoxy prior to thin sectioning. Other methods of improving the contrast between the pore and solid include injecting the pore space with otherwise modified epoxy, woods metal, or wax (Straley and Minnis, 1983, Soeder, 1990, Winnsauer et al. 1952, and Nutting, 1929). These methods all make it possible to see pores which, due to their small size or proximity to a large grain, would otherwise be impossible to see. However, even using this technique, microporosity is often difficult to differentiate from the surrounding grains because many of the injected materials are translucent in the small quantities required to fill such porosity. During the past two decades, the imaging of microporosity has been achieved using both the fluorescent optical microscope and the scanning electron microscope (SEM). Fluorescent dyes have been employed to render microporosity visible to fluorescent optical microscopes (Soeder, 1990, Ruzyla and Jezek, 1987, Rassineux et al., 1987, and Gies, 1987) and epoxy pore casts have been used to create high contrast images using the SEM (Pitmann and Duschatko, 1970, Slraley and Minnis, 1983, Lin and Hamasaki, 1983). The pore casts examined with the SEM consist of untreated epoxy casts (Pitmann and Duschatko, 1970) and “double” pore casts, in which an untreated and a brominated epoxy were used. These “double” pore casts were produced by injecting epoxy into the pore space, dissolving the rock with hydrochloric and/or hydrofluoric acid and then injecting a second epoxy in place of the rock (Straley and Minnis, 1983, Lin and Hamasaki, 1983). Due to the large difference in atomic number between untreated and brominated epoxy, a high contrast back- scattered SEM image is obtained from such double pore casts, allowing very small pores to be imaged. The methods just described solve the problem associated with the low contrast between the pore space and grains in most geological samples; however, the lack of quantifiable three dimensional data is a problem which, until recently, was not addressed. Three dimensional data is important because pore networks are inherently three-dimensional while most images are only two dimensional. Three-dimensional data has been obtained both with the SEM and optical microscopes by collecting stereo pairs of pore networks (Pitmann and Duschatko, 1970); however, 7 these images are very difficult to interpret quantitatively. The use of automated computer analysis is making this approach more viable, as Lin and Hamasaki (1983) have shown in their work on three-dimensional reconstruction and analysis. The only other approach which has been attempted, in order to obtain three-dimensional images, was to physically section the sample to obtain a series of two-dimensional images, or serial sections. Two methods routinely used to produce serial sections are: successively grinding layers off the surface of a sample (Lin and Hamasaki, 1983), and producing a two-epoxy pore replica which is then microtomed to produce successive thin sections of the epoxy (Straley and Minnis, 1983, Lin and Hamasaki, 1983). Both of these methods allow a “reconstruction” of the series of two-dimensional sections to fomi a three-dimensional picture of the network. The drawbacks associated with these methods are primarily due to the distortion of the samples during the physical sectioning and the fact that large amounts of time and skill are required to produce the data. When a sample is ground, the exact amount removed between imaging must be the same each time and, more importantly, the successive surfaces must be parallel, If distortion of the slices is to be minimized, care is also required when slices are produced with a microtome. Even when a microtome is skillfully operated, a permanent distortion of five percent is usually expected and the force of the blade used to cut the sample may cause the sample to twist, thus introducing more error (Straley and Minnis, 1983). The method used in this thesis to image pore networks avoids many of the problems encountered with the methods described above. A pore cast of the sample can be used to produce a series of two-dimensional images without physically sectioning the sample, and the collection and analysis of the images is relatively quick and easy. Also, the pore cast is not altered in any way, therefore re-examination of the same areas in their original orientation is always possible. The production of serial sections without physical sectioning is made possible by the use of the confocal laser scanning microscope (CLSM), the history and capabilities of which are described below. 8 THE CONFOCAL LASER SCANNING MICROSCOPE (CLSM) The CLSM is one of many microscopes made possible by the introduction of the conventional scanning microscope (Young and Roberts, 1951). Conventional microscopes image over a large field of view, while scanning microscopes are only required to image one point at a time and the final image is constructed through scanning. Removing the requirement to image a large area allowed the introduction of modified optical systems for specialized use. One of these modifications produced the confocal scanning microscope. The commercial confocal scanning microscopes available now are generally reflection and/or fluorescence microscopes. The light source is a laser focussed through a pinhole by a condenser lens. The light from this point source is reflected through the objective lens by a dichroic mirror which focusses it on the sample. Reflected or fluorescent light from the sample is collected by the same objective lens, passed back through the dichroic mirror and focussed into a detector pinhole. The detector, often a photomultiplier tube, measures the intensity of the light reflected or emitted by the sample (Bacallao and Steizer, 1989). The confocal optical system described above was first described by Minsky (1957), who recognized its possibilities for enhanced resolution and depth discrimination. The enhanced resolution is a result of Lukosz’s principle, which states that resolution may be increased at the cost of the depth of field (Lukosz, 1966). The depth discrimination or optical sectioning property of the confocal microscope is produced by the combined use of a point source and point detector. Thus, light from either above or below the focal plane is focussed behind or in front of the detector pinhole and is not imaged. The accuracy of the depth discrimination has been discussed and calculated theoretically by many researchers: for a dry lens with 0.95 numerical aperture, the depth of field is approximately 0.24 .tm and for oil immersion lenses, even smaller depths of field are possible (Russ, 1990). The loss of a large depth of view makes the confocal microscope unsuitable for those applications where a three-dimensional view is required or the sample is very uneven, causing a 9 great deal of scattering and thus degrading the image. However, there are many applications in which the optical sectioning property is actually a requirement for producing high quality images. Applications of microscopy which have been made possible by the optical sectioning of the confocal microscope include many previously iniractable problems in the biological sciences (Wilson, 1990 and Shuman et al., 1989). The primaiy difficulty facing most biologists using conventional microscopes was the haze of out-of-focus fluorescence produced by overlapping tissues. This haze hopelessly degraded the images and was first addressed by microtoming cells which had been fixed in plastic or wax. Slices 2 lim to 10 p.m thick were removed with the microtome, transferred onto a slide, and examined with a microscope (Marshall et aL, 1990). Microtomoscopy, however, has similar difficulties in biology to those encountered in geology. An immense number of slides have to be imaged, it is difficult to re-align the images, and distortion in the original architecture of the cell are caused both by the fixation process and by microtoming. The majority of the difficulties encountered in mechanically sectioning cells are avoided by the use of the CLSM. Cells of any shape and size may be viewed, and most importantly, they may be viewed without fixation. In some cases, they may be imaged while functioning (Bacallao and Stelzer, 1989, Egger and Petran, 1967). Thus a cell process, from neuron activity to bacterial growth, may be monitored as it occurs (Masters et al., 1991). These applications are only a small portion of those possible with the confocal microscope. Recently, a CLSM has been used in metallurgy to count pores in sintered alumina (Russ et al., 1989), in engineering to do surface profile measurements (Hamilton and Wilson, 1982), and by physicists to image the surfaces of microchips (Shuman et al., 1989). The success of the scanning confocal microscope in these fields suggest that it could also be an effective technique for characterizing porous media in the earth sciences. A more complete treatment of possible applications and adaptations of the confocal microscope for specialty work is contained in Confocal Microscopy (Wilson,1990). 10 EXPERIMENTAL PROCEDURES Sample Preparation Two sedimentary rock samples, one high porosity, high permeability quartz sandstone, and a low porosity, low permeability carbonate, were used to develop and test the epoxy impregnation techniques employed to create pore casts for this study. (Crocker et al., 1983). The quartz sandstone, Berea 100, has frequently been used in studies of porous media, therefore it was selected to determine the accuracy of the imaging and image analysis techniques. The carbonate, referred to as samples S55, is Devonian in age and composed primarily of calcite with less then one percent dolomite and silicate minerals. The low permeability and low porosity of sample S55 increased the difficulty of thoroughly impregnating this sample with epoxy. Thus the efficiency of the impregnation techniques tested were judged by how thoroughly S55 was infiltrated. The samples were cut with a rock saw to expose two parallel surfaces and small particles were cleaned off using compressed air. The dimensions of sample S55 was approximately 3mm x 4mm x 4mm. The Berea 100 samples were disks 8.9 mm in diameter and 3 to 4 mm thick. The carbonate sample was smaller than the Berea 100 sample due to limitations in the available sample. A pore cast of the third sample considered in this study was prepared by Core Laboratories in Calgary; the others were not prepared by Core Laboratories because their techniques were not suited to impregnating samples with permeabilities as low (< lmd) as that of sample S55. This third sample is a high porosity, high permeability limestone referred to as sample S77. S77 is Devonian in age and composed primarily of calcite with approximately ten percent dolomite and less than one percent silicate minerals present. 11 Materials Impregnation of permeable media with epoxy and/or other fluids has been used in such varied fields as soil science, biology, and geology. Biologists use epoxy to fix tissues in place for examination with microscopes. Soil scientists employ epoxy and wax to stabilize soils so the inter relationship of voids and grains in the samples may be studied. Geologists use epoxy impregnation to increase the visibility of the void space in their samples; however, a geologist’s rocks are sturdier than soils and larger than tissue samples, thus different epoxies and techniques are used to impregnate these three types of samples. Both the characteristics of the sample and the purpose for which it will be impregnated must be considered when selecting an epoxy. Many rock samples considered have a very low permeability, thus a low viscosity epoxy is required. The formation of a pore cast requires the exposure of the epoxy to dilute hydrochloric and/or hydrofluoric acid, therefore it must not react with these fluids once it is cured. Finally, the low permeability of many rock samples requires that they be immersed in the epoxy for several hours to ensure complete impregnation. This means that the pot life, or the amount of time for which the epoxy retains its initial viscosity, must be relatively long, with epoxies which require heating in order to cure being preferable to those which cure rapidly at room temperature. There are at least two epoxies which match these requirements: Maraglas 655 resin with 555 hardener (Soeder, 1990) and Spurr low viscosity embedding media (Spurr, 1969). Maraglass 655 is a low viscosity epoxy resin with a pot life of several days. The Spurr low viscosity embedding media was introduced by Spurr (1969) and has a viscosity of 7.8 cP at 20°C, similar to that of water, which is 1.002 cP at 20°C. Several different mixtures of the hardening and epoxy components are possible to produce fluids with different pot lives and curing schedules for the Spurr mixture. These possibilites are presented in Table 2.1. The mixture selected produced an epoxy with a pot life of approximately seven days and is the longer pot life/lower viscosity version 12 T A B LE 2. 1 Sp ur r Lo w -V is co si ty E m be dd in g M ed ia * St an da rd m ed iu m Su gg es ted m o di fi ca ti on s o f th e m ed iu m A B C D E IN GR ED IE NT FI RM H AR D SO FT RA PI D CU RE LO NG ER PO T LI FE LO W ER V IS CO SI TY Vi ny lcy clo he xe ne ox id e (V CD ) 10 .0 g 10 .0 g 10 .0 g 10 .0 g 10 .0 g Di gy cid yl et he r o f po ly pr op yl en eg ly co l( DE R 73 6) 6.0 g 4. 0 g 7. 0 g 6. 0 g 6. 0 g N on en yl s u c c in ic an hy dr id e (N SA ) 26 .0 g 26 .0 g 26 .0 g 26 .0 g 26 .0 g D im eth yl am in oe th an ol (D MA E) 0. 4 g 0. 4 g 0. 4 g 1.0 g 0. 2 g Cu re Sc he du le (hr )a t 7 0° C 8 8 8 3 16 Po tL ife (da ys) 3- 4 3- 4 3- 4 2 7 * Sp uf f, 19 69 C ) listed in Table 2.1. Both the Spurr and Maraglas epoxies could have been used to produce pore casts; however, Spurr’s epoxy resin was easily available and had a longer pot life, therefore it was selected. The pore casts were to be imaged using the fluorescence mode of a CLSM, so the epoxy had to be fluorescent before it was used to impregnate the rock samples. This requirement was not difficult to meet as most epoxies, including Spurr low-viscosity embedding media, are naturally fluorescent in the 460 to 490 nm range (Cercone and Pedone, 1987). This natural fluorescence is relatively weak, therefore the fluorescence of the epoxy was enhanced by the addition of 0.1 g of Rhodamine B per 25 g of epoxy. The visibility of the pore cast to optical microscopes was also improved by the addition of approximately 0.1 g of oracet blue dye per 25 g of the epoxy. Pore Cast Production Three impregnation techniques were tested to determine which method was best suited to our samples. In all techniques, the samples were dehydrated by heating in an oven at 105°C for several hours prior to impregnation. The dehydration was necessary because Spurr’s epoxy resin is incompatible with water, therefore the sample must be completely dry to enabe the epoxy to fully penetrate the void space. The samples were impregnated and cured in a teflon sample holder. Teflon sample holders were used because teflon is resistant to the 70°C temperature required for curing and it does not react with either the epoxy or the two acids used to etch the samples. The holders were cylindrical and the sides tapered at approximately 5° so that a few taps on a counter would cause the sample to simply slide out of the holder. After curing the sample at 70°C for at least 16 hours, the impregnated sample was removed from the holder and surface ground to expose the rock. The surface grinding was accomplished with grinding wheels, first with 250 mesh carborundum, then with 600 mesh (Al203)grit. A final grinding was done by hand with 900 mesh (A1203)grit on glass to smooth the surface. Problems 14 associated with grinding the sample include smearing the epoxy, removing too much of the sample, and not removing enough. Smearing of the epoxy over the grains and failing to remove the irregularities of the original surface may result in the pores appearing much larger than their actual size. Alternatively, removing too much epoxy may result in the loss of the pore cast in cases where the epoxy did not fully penetrate the sample. The exposed surface was etched with acid to remove the rock, leaving a cast of the pore system. Two acids were used: hydrochloric acid (HC1) to dissolve carbonate rocks, and hydrofluoric acid (F1F) to dissolve silicates. Impregnation Techniques Wet Impregnation A variation of the ethanol replacement method used in soil science and biology (Marshall et aL, 1990 and Smart and Tovey, 1982) was attempted. In this technique, the dehydrated sample was immersed in ethanol immediately after its removal from the oven. Several hours were allowed for the ethanol to completely saturate the sample and the degree of saturation of the samples with ethanol was checked by comparing the weight of the sample before and after soaking in the ethanol. The volume of ethanol adsorbed was calculated and compared to the porosity of the sample; every sample was found to be at least 80% saturated with ethanol. The sample was returned to the ethanol for an hour before the process of replacing the ethanol with epoxy was initiated. Ethanol was replaced with epoxy by removing fluid (epoxy and/or ethanol) from the sample holder until it barely covered the sample, then epoxy was poured in to increase the fluid level to at least twice its original level. This was repeated every two to three hours throughout the day for two days, for a total of seven or eight changes in the fluid. The epoxy is soluble in ethanol, thus if the two fluids mix thoroughly between each addition of epoxy, the amount of ethanol remaining in the fluid should be less than one percent. To ensure that the mixing was as 15 complete as possible, a special sample holder was constructed in which the sample was suspended on a porous platform and a stirring rod below the platform agitated the fluid continuously to aid in the mixing of the ethanol and epoxy. This method allowed the epoxy to easily penetrate even low permeability samples; however, the small amount of ethanol remaining in the epoxy after it was set reacted with the hydrofluoric acid, distorting or completely destroying the pore cast. Therefore this technique is not suitable for samples which will be etched with acid. Vacuum Impregnation: Immersion Technique The simplest technique used required the dehydrated sample to be immersed in the epoxy immediately after removal from the oven. The sample, immersed in epoxy, was placed within a vacuum dessicator and a vacuum applied to remove the air from the voids and allow the epoxy to penetrate them. The sample remained within the vacuum dessicator for four to five hours, the first two to three hours under vacuum and the final two hours equilibrating to the surrounding air pressure. Unfortunately, a high pressure vacuum was not available to pull the air out of the pores, thus the low permeability samples were not impregnated beyond the surface layer of the sample. Vacuum Impregnation: Lowering Technique The most successful technique used a form of vacuum impregnation that required the sample and the epoxy to be degassed separately, then the epoxy either poured over the sample or the sample lowered into the epoxy. The epoxy may be poured over the sample either from a beaker or a separatory funnel. The beaker must be suspended over the sample and it must be possible to manipulate it without opening the vacuum chamber. Vacuum apparatus similar to this is commercially available from several companies; however, a simpler solution is to use a lowering apparatus. This requires a vacuum vessel sealed with a rubber stopper through which a rod may be extended into the vessel. The sample holder containing epoxy is placed on the base of the 16 vessel, and the sample is suspended over the sample holder by attaching it to the rod with string or netting. The vacuum is applied to both the sample and the epoxy for an hour to remove all the air from both. The sample is then lowered into the epoxy and the netting or thread is released. The vacuum is continued for approximately thirty minutes to remove any gas which is trapped in the voids of the sample. Etching Following the curing of the epoxy, the sample was removed from the sample holder and surface ground to expose the rock surface as described previously. The exposed rock and epoxy sample was immersed in five percent hydrochloric acid (HC1) for between five and sixty minutes depending on the depth to which etching was desired. This dissolved any carbonate rock present in the upper portion of the sample. Residual HC1 was removed by rinsing the sample repeatedly with deionized water and then any silicate minerals present in the sample were dissolved by immersing the sample in 5% hydrofluoric acid (HF) for between five and sixty minutes. The HC1 was used first because HF reacts with carbonate rock to form insoluble calcium fluoride (Grayson, 1956). Both HC1 and HF were tested to ensure they did not react with Spurr’s epoxy resin by Pittman and Duschatko (1970). Residual acids were removed by repeated rinsing with deionized water before any imaging was done. All of the impregnation techniques described above worked quite well for the Berea 100 sample; however, the microstructure of the pore network in the low permeability carbonate, S55, was too fragile and collapsed when the surrounding mineral was removed. This problem was overcome by producing “double” pore casts of this sample. 17 Double Pore Casts Double pore casts were produced by re-impregnating the sample with a different epoxy after the carbonate portion was dissolved with HC1. The purpose of this second epoxy was to support the pore cast. The second epoxy used was London Resin White (LRW), a transparent, low viscosity, non-fluorescing resin with a very long pot life. This epoxy was cured at 72°C for at least 24 hours. Impregnation of the sample with a second epoxy stabilized the microstructure of the void space enough to enable the subsequent removal of silicate minerals in the sample using 5% HF for between five and sixty minutes. The LRW was tested to ensure it did not react with either HC1 or HF by immersing a block of it in both fluids for several hours and subsequently examining the block microscopically. The examination revealed no alteration in the character of the LRW so it was assumed safe to use in the production of these “double” pore casts. DATA COLLECTION AND ANALYSIS Images were collected using an MRC BioRad 500, a highly versatile microscope which is commercially available and operates either in reflection or fluorescence mode. The images collected with the BioRad 500 were 768 by 512 pixels in size and each pixel had an intensity value which corresponded to the intensity of the fluorescence at that point in the sample, with 0 being no response and 255 being the most intense response. The portion of the sample being imaged was scanned 10 to 15 times and the images were stacked to improve the signal to noise ratio. This was accomplished interactively; when the clarity of the image ceased to improve, the scanning was halted. The images were stored in a binary format, thus the algorithm written in the C programming language and listed in Appendix 1 was required to transform them into a format (ASCII) readable by the FORTRAN programs written to analyse the images. The primary requirement for the collection of data is that enough images be collected to provide a statistical sampling of all pore and pore throat sizes. Forty two-dimensional images of 18 Berea 100, twenty-one of sample S55, and thirty-three of sample S77 were collected. Some of these images were part of a vertical series of contiguous images, and therefore only one of the series was analysed to determine pore sizes to avoid biasing the pore size distribution towards the pores in that portion of the samples. The final analyses were accomplished with sixteen Berea 100 images, twenty-one S55 images, and thirteen S77 images. These images were manipulated to produce pore and pore throat distributions using the process described below, which was implemented using the programs listed in Appendix 2. Initially, the image was altered from one with an intensity range of 256 to one which was binary. The process of altering an image with a large intensity range to a binary image is often accomplished by visual means using an image analysis program or by using the intensity contrast between adjacent pixels to choose the boundaries of the objects. The separation of these images was accomplished similarly by using an imaging program to select the intensity ranges which corresponded to epoxy and void (pores and grains). Also, several rows of each image were plotted and the separation value chosen to separate pore from solid was super-imposed on this plot. A value which lay below the peaks which represented the pores and above the lows which represented the solid regions was chosen. Examples of the value chosen super-imposed on a plot of the image intensity are shown in Figure 2.1 for two rows in the image referred to as ABR1. The resulting binary images when the separation value chosen is varied are shown in Figure 2.2. The small holes in the pores and solid which may have been due to either the presence of very small pores or to uneven fluorescence intensity in the epoxy were removed. This was accomplished using a dilation/erosion sequence followed by an erosion/dilation sequence. These are standard procedures in image analysis. To illustrate the procedure, the process of removing small pores within the solid is described. First the solid is dilated by altering any pore pixels in contact with a solid to a pseudo-solid, then any pixels added to the solid via dilation which were in contact with a pore are changed back into pore pixels. Thus, a pore which was less than 1.5 pixels in diameter or, for the images collected at the highest magnification, 3.0 tm in diameter, was changed to solid. Similarly, small solid areas completely within a pore were removed by first 19 IMAGE ABR1 - ROW 40 > z z PIXEL NUMBER IMAGE ABR1 - ROW 250 I I I 300 200 100 0 300 0 200 400 600 800 200 l00 0• LL LrJ - fiJFJ J IrU\%J% 0 200 400 600 800 PIXEL NUMBER Figure 2.1: Pixel intensities in rows 40 and 250 of image ABR 1. The intensity value 50 is indicated as possible separation value between pore pixels and solid pixels. 20 Figure 2.2a: Binary versions of Berea 100 image ABR1 with pore/solid separation value of 40. Black = solid, and white = pore space. 2OOn 21 2OO.zm Figure 2.2b: Binary versions of Berea 100 image ABR1 with pore/solid separation value of 50. Black = solid, and white = pore space. 22 2OOini Figure 2.2c: Binary versions of Berea 100 image ABR1 with pore/solid separation value of 60. Black = solid, and white = pore space. 23 dilating the pores, then removing (eroding) those pixels which remained in contact with the solid. The pores were separated and those points within the pores which were farthest from the edges of the pores were found using the watershed algorithm described by Russ (1990). The most central coordinates of each pore and pore throat were recorded so that the size of the pores and pore throats could be determined. The size of the pores and pore throats was determined using methods similar to those employed by Ehrlich et al. (1984) and Doyen (1988). The pixels within the pores and pore throats which are most central are used to compare the pore or pore throat to a previously constructed circle in order to determine the largest diameter circle which is completely enclosed by the pore or pore throat. The radii measured in the manner described above are presented as volume fraction distributions which are the ratio of the area of the circles described by each individual radius to the total area of all of the measured pores. This ratio is assumed to be proportional to the corresponding three-dimensional ratio of the volume of a pore to the total volume of the pore space because, if sufficient sections are randomly imaged from a sample, the sizes, numbers, and area fractions of the pores images should be equal to the sizes, numbers, and volume fractions of the three dimensional objects (Wicksell, 1926, and Russ, 1990). RESULTS The images produced with the confocal laser scanning microscope were clear, high-contrast images which clearly show a major variation in pore shapes and sizes between the three samples. The Berea 100 images (Figure 2.3) contained the expected network of larger pores connected by much smaller pore throats. The carbonate images were dominated by thin, plate-like pores which were often associated with pores greater than 100 urn in diameter (Figure 2.4a). S55 contained only a few of the larger pores (Figure 2.4b); however, S77 contained a large number of such pores (Figure 2.5). 24 2OOim Figure 2.3: Original version of a confocal laser scanning microscope image of Berea 100. The image is called ABR1. The pore space is the brightest areas in the image. 25 2OOim Figure 2.4a: CLSM image of vugs in the carbonate sample, S55. The image is called F78 and the brightest areas of the image represent porosity. 26 Figure 2.4b:CLSM image of vugs in the carbonate sample, S55. The image is called T78 and the brightest areas of the image represent porosity. 200qxn 27 Figure 2.5: Original version of a CLSM image of carbonate sample S77. The image is called Fl. The pore space is the bright areas in the image. 28 29 The choice of a “cut-off’ or pore/solid separation value to alter the original images to binary ones was different for each sample. The magnitude and overall distribution of porosity changed very little for the Berea 100 image, referred to as ABR1, for solid/pore separation values between 40 and 60 (Figure 2.2), therefore the final value selected was 50. This value was tested on the other Berea images and found to be appropriate for them as well. The S55 images were more difficult to separate into two phases. The fluorescence of the thin, plate-like pores was faint compared to the larger pores, thus when a high separation value was chosen, the thin pores were neglected (Figure 2.6a), while if a lower value was used, the resulting image contained porosity which, in the original, seemed to be an out-of focus haze (Figure 2.6b). Separation of the S77 images into pores and solid was relatively easy as the intensity of the pore and solid regions were well separated in these images (Figure 2.7). The pores and pore throats were measured from the coordinates which were farthest from the edges of the pores or pore throats (Figure 2.8). Pore size and pore throat distributions were obtained from the binary images of Berea and S77, and from those binary S55 images which were produced using a high separation value. The pore size distributions presented here are volume fraction distributions, therefore the percentage associated with each pore radius indicates the fraction of the total pore volume occupied by pores of that size. The Berea 100 pore radius volume fraction distributions all have a similar shape (Appendix 3), as shown by the examples in Figure 2.9a. The distributions extend from 1 to 60 $.tm with the maxima occurring primarily between 10 and 30 im with several isolated peaks above 30 I.Lm.The pore size distributions of several of the images combined also extends from 1 to 60 .tm and the maximum occurs between 10 and 20 Lm (Figure 2.9b). The pore size distributions for sample S55 were quite different (Figure 2.lOa) and a copy of the distribution for each image is contained in Appendix 4. These distributions were less uniform, with most maxima occurring between 7 and 20 Lm and a few between 20 and 60 I.tm. The combined pore size distributions of all of the images extends between land 50 Im and reaches 30 Figure 2.6a: Binary version of carbonate S55 CLSM image F78 (Figure 2.4a) produced using a solid/pore separation value of 50. White = porosity, and Black = solid. 2OOun 31 Figure 2.6b: Binary version of carbonate S55 CLSM image F78 (Figure 2.4a) produced using a solid/pore separation value of 100. White = porosity, and Black = solid. 200pm w ‘ 32 Figure 2.7: Binary version of carbonate S77 CLSM image Fl (Figure 2.5) produced using a solid/pore separation value of 15. White = porosity, and Black = solid. 33 C ) I Figure 2.8: Binary versions of Berea 100 image ABR1 with pore throats separated from the pores. The centers from which the sizes of the pores and pore throats will be measured are indicated by a contrasting colour within the void space. Black = solid, the White = porosity, and the dotted areas = pore throats. 35  Figure 2.9a: Examples of the Berea 100 pore size distributions. The upper plot represents the pore size distribution of image ABR1, the lower of image HBR1. Image ABR1 106 0 .— I. J 0 C) 0 10 Pore Radius (m) 0.12 0.08 0.04 0.00 0.08 0.06 0.04 - 0.02 - 0.00 10 Image HBR1 & I ii 106 10 -I Pore Radius (m) 37 0 .— 0 0.08 5 Images A 0 C.) J 0.06 0.04 0.02 10 Pore Radius (m) L i04 16 Images 0.08 0.06 0.04 0.02 0.00 Figure 2.9b: Berea 100 cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 16 images. Pore Radius (m) 38 0 .-4 0 0 0 0 Image F78 0.3 0.2 0.1• 0.0 106 0.2 0.1 10 Pore Radius (m) i04 Image N78 0.0 10 L 10 Pore Radius (m) Figure 2. 10a Examples of the carbonate (S55) pore size distributions. The upper plot represents the pore size distribution of image F78, the lower of image N78. 39 20 Images Figure 2. lOb:Carbonate (S55) cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 20 images. 5 Images 0 C) J 0 0 Pore Radius (m) i04 0.12 0.08 0.04 0.00 -i 106 ii 10 Pore Radius (m) i04 40 a maximum at approximately 7 pm (Figure 2. lOb). The S77 pore size distributions (Appendix 5) resembled discrete pore size spectra , with the largest pores present in the image dominating the distribution (Figure 2.11 a). However, the largest pores present in the image varied from 50 to 220 pm in radius, resulting in a cumulative distribution which extended primarily between 10 and 220 pm with several maxima between 70 and 150 p.m (Figure 2.1 lb). The mean values of the volume fraction pore size distributions were calculated and are collated in Tables 2.2, 2.3 and 2.4 for the Berea 100, S55 and S77 samples respectively. The mean pore size of the Berea images varied from 20.2 to 30.4 p.m, with an overall mean pore size of 25.3 microns and a standard deviation of 15 pm for all of the images considered together. The pore sizes in S55 are generally smaller, ranging from 8.9 to 33.0 pm, with an overall mean pore size of 17.5 p.m with a standard deviation of 11 pm. The largest pore sizes and largest variation in mean pore sizes occur in sample S77. The mean pore size of this sample ranges from 30.6 to 150.9 pm for the individual images with an overall mean pore size of 90.2 p.m and standard deviation of 54 p.m for the cumulative pore size distribution. The porosity of the Berea 100 images varied a great deal, from 0.26 to 0.41 (Figure 2.12), with the majority of the images having porosities between 0.29 and 0.34. There is a similar variation in the porosity of the images for the carbonates, S55 (Figure 2.13) and S77 (Figure 2. 14).The minimum and maximum porosities of the S55 images produced using a high separation value are 0.04 to 0.16, with a single outlier at 0.25. The porosity of these same images when a low separation value was selected was much higher, ranging from 0.11 to 0.39, with most of the images having porosites between 0.15 and 0.25. S77 image porosities varied from 0.02 to 0.18 with an average of 0.088. 41 0 .— a) E 0 0 .-4 C) 0 Image Ji Image Gi 4AA1l 106 io4 Pore Radius (m) 0.5 0.4 0.3- 0.1 0.0 0.4 _ 0.3 0.1 - 0.0- I 106 io’ io3 Pore Radius (m) Figure 2.1 la: Examples of the carbonate (S77) pore size distributions. The upper plot represents the pore size distribution of image G 1, the lower of image Ji. I 42 5 Images io-5 Pore Radius (m) 13 Images 0.2 0.1 0.0 1 0.08 0.06 0.04 0.02 0 Q C > 0 C-) 0 ) ‘.3 liii I 0.00 1O..6 10” Pore Radius (m) Figure 2.1 ib: Carbonate (S77) cumulative pore size distributions. The upper plot represents the summation of the distributions of 5 images, the lower of 13 images. 43 C’) LU LL 0 LU D z BEREA IMAGES - POROSITY POROSITY Figure 2.12: Porosity (the ratio of pore pixels to the total number of pixels) of the Berea 100 images. CON. CD O CD.- CJ CD L() CON. CD O — 0000 000000000 0 44 S55 IMAGES - POROSITY 2 _____ LOW SEPARA11ON VALUE MEAN POROSITY: 0.22 0 - . —v—v—. . . - . . - . .i— .. • — -C’ e Ifl %O N 00 O- ‘f ‘D N 00 O\ -C4 fl ‘.0 N 000’. POROSITY S55 IMAGES - POROSITY 3 HK3H SEPARATION VALUE MEAN POROSITY: 0.11 C,, 1 0 Il’. ‘.0 N 00 0’ — - (‘. ‘.0 N 00 0’. - C.l e - In 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 POROSITY Figure 2.13: Porosity (the ratio of pore pixels to the total number of pixels) of the carbonate S55 images. The upper distribution represents the binary images produced using a low intensity separation value, and the lower distribution represents those images produced with a high intensity separation value. 45 N U M B ER O F IM A G ES 0 • 1% ) 0. 02 0. 03 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ Ci ) 0. 04 0. 05 _ _ _ _ _ _ _ _ — 0. 06 Z 0.0 7 0.0 8 C- ) 0. 09 I. < 0 o .i CI ) C D ‘ — 0. 11 _ 0. 12 C - 0. 13 o _ _ _ _ _ _ _ _ _ 0. 14 CD 0. 15 Ci ) — 0. 16 0. 17 0. 18 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ CD 0 0) TABLE 2.2 Berea 100 Mean Pore Sizes Image Mean Pore Standard Number of Mean Pore Standard Radius Deviation Images Radius Deviation jim jim jim jim ABR1 24.8 13 1 24.8 13 BBR1 23.7 15 2 24.2 14 CBR1 22.8 13 3 23.7 14 DBR1 26.4 16 4 24.4 15 FBR1 20.2 9 5 23.7 14 GBR1 21.5 12 6 23.3 14 HBR1 26.4 17 7 23.8 14 JBR1 28.6 19 8 24.3 15 KBR1 26.1 17 9 24.5 15 LBR1 24.9 13 10 24.6 15 MBR1 30.4 17 11 25.2 15 NBR1 27.4 17 12 25.4 15 OBR1 22.0 12 13 25.2 15 PBR1 24.3 14 14 25.1 15 QBR1 22.6 10 15 24.9 15 RBR1 29.7 17 16 25.3 15 TABLE 2.3 Carbonate (S55) Mean Pore Sizes Image Mean Pore Standard Radius Deviation jim j.tm A78 12.2 6 B78 18.4 11 C78 17.2 10 D78 15.7 10 E78 15.7 10 F78 10.6 6 G78 14.4 8 H78 8.9 4 J78 19.0 10 K78 14.0 6 L78 12.2 6 M78 14.1 8 N78 17.1 13 078 12.0 7 P78 15.7 9 Q78 29.2 16 R78 12.1 7 S78 18.0 9 T78 21.0 9 U78 33.0 13 Number of Mean Pore Standard Images Radius Deviation j.tm jim 1 12.2 6 2 17.0 11 3 17.1 10 4 16.6 10 5 16.4 10 6 15.8 10 7 15.6 10 8 15.0 10 9 15.5 10 10 15.4 10 11 15.2 9 12 15.1 9 13 15.4 10 14 15.1 10 15 15.2 10 16 16.6 11 17 16.3 11 18 16.4 11 19 17.0 11 20 17.5 11 47 TABLE 2.4 Carbonate (S77) Mean Pore Sizes Image Mean Pore Standard Number of Mean Pore Standard Radius Deviation Images Radius Deviation j.im jim jim jim Fl 113.6 50 1 113.6 50 Gi 150.9 71 2 126.3 60 Hi 38.2 13 3 121.8 62 Ii 30.6 13 4 119.3 63 Ji 97.6 33 5 116.3 60 Ki 50.5 14 6 111.3 60 Li 99.2 48 7 108.1 58 Ml 57.5 25 8 103.3 57 Ni 35.7 19 9 100.7 58 01 75.1 53 10 97.3 58 P1 63.4 26 ii 93.6 56 Qi 90.2 43 12 93.2 55 Ri 56.1 26 13 90.2 54 DISCUSSION AND CONCLUSIONS The confocal laser scanning microscope produced very clear, easily analysed images when the porosity being imaged was relatively uniform in size. Unfortunately, the S55 carbonate sample considered in this project contained pores of widely varying sizes, and the enhanced fluorescence of the larger pores overshadowed that of the small, thin pores. Thus the S55 images, although qualitatively useful, were very difficult to convert to a binary format. The two-phase images produced by using widely separated cut-off values (Figures 2.6a,b), when compared to the original image (Figure 2.4a) illusirate the primary difficulty in analyzing these images: the thin pores are either neglected completely (Figure 2.6b) or they are represented as being much thicker than they actually are (Figure 2.6a). Sample S77, although it also contained pores of widely varying sizes, did not produce the same problems. This may be a result of two effects: the epoxy/dye mixture used by Core Laboratories to impregnate this sample produced a stronger flourescence, or, considering Core Laboratories lack of success in impregnating S55, the smallest pores may not have been penetrated using their technique. 48 A possible solution to imaging samples with two or more widely separated pore sizes is to image the samples at two magnifications. The higher magnification images may be used to determine the size of the smaller pores and the low magnification images to correlate the number or volume of small pores associated with the larger ones. In this way the final pore size distribution could be determined by estimating the number and size of the associated smaller pores when the size distribution of the larger pores is determined. The choice of which pore/solid separation value to use in further analysis of the S55 carbonate images was based on the porosity of the resulting images. The measured porosity of the sample was 0.048 (Chapter 3), therefore the largest value of the cut-off reproduced the correct porosity most accurately with an average porosity of 0.11, (Figure 2.10), and these images were used for all further analysis. The large porosities of the images, an average of 0.322 compared to the measured value of 0.21 for Berea 100 (Tercier, 1992) and an average of 0.088 compared to the measured value of 0.09 for S77 (Chapter 3), are primarily due to experimental bias. The images are collected by scanning the surface of the sample and “randomly” selecting a portion to image. The selection was accomplished interactively, thus areas with little or no porosity were not imaged. If this experiment were repeated, one possible improvement is to mechanize the imaging to ensure that the entire surface of the sample is imaged regardless of the porosity included in the image. Alternatively, the number of image areas containing no pores could be recorded and included in the calculation of the total porosity of the sample. An example of the locations of the centers used to measure the sizes of the pores and pore throats in the images is present in Figure 2.8. This figure shows how these centers were spaced to ensure a minimum amount of overlap between the areas included within the measured circles. This results in elongated pores containing three or more centers, thus illustrating the primary difficulty in determining pore sizes: their shapes are not easily described. The assumption of circular pores which was used here is frequently used, despite its inaccuracy, because it simplifies most calculations. A more representative shape for pores would be elliptical; however, measuring the 49 sizes would be complicated by the irregularity in the outline of the pores. A possible solution to this difficulty is to approximate the size of ellipse which may be contained within a single pore based on the circles already used to measure the pores. This could be accomplished by determining the length of the major and minor axes of the ellipse by summing the diameters of circles paralleling the major axis and averaging those paralleling the minor axis. The overall mean pore sizes of the Berea 100 and carbonate images are listed below (Table 2.5). The results for the Berea 100 images were very encouraging, showing good agreement with mercury injection pore size distribution of Berea 300 presented by Crocker et al. (1983). The Berea 300 pore size distribution may be compared to that of Berea 100 as they are both Berea sandstones, the only variation between the the two being that Berea 100 has a permeability of approximately 100 mD, while Berea 300 has a permeability of approximately 300 mD. The Berea 300 distribution (Figure 2.15) has a mean pore radius of 9.0 jim, compared to a cumulative mean void radius of 25.3 p.m determined for the Berea 100 images. The 9.0 p.m determined by Crocker et al. (1983) is substantially lower, as expected since capillary presssure mercury injection data reflects the size of pore throats, making the mean size lower than the actual average pore size. TABLE 2.5 Mean Pore Sizes Sample Mean Void Radius Standard Deviation jim jim Berea 100 25.3 15 S55 17.5 11 S77 90.2 54 50 I40 30 20 10 0 Po re D ia m et er (pm ) Fi gu re 2. 15 : B er ea 30 0 ca pi lla ry pr es su re po re siz e di str ib ut io ns . (fr om Cr oc ke r e ta L, 19 83 ) 10 0 B er ea M ea n— 18 .0 gm M od e= 10 .0 u rn 0. 2 0. 4 1 2 4 10 20 40 (31 The primary objective of this study was to develop a method of collection and analysis of three-dimensional images of pore networks. Despite the lack of complete three dimensional images in this study, the confocal laser scanning microscope is a powerful tool which could be used to accomplish this goal. Unfortunately, the fragility of the pore casts ensured that the series of two- dimensional images which were obtained (Figure 2.16) did not encompass a thickness greater than the average pore diameter and thus could not be used to directly determine three-dimensional pore sizes. However, the correspondence between the Berea 100 pore size distribution and that determined by Crocker et al. (1983) for Berea 300 indicate that CLSM imaging is definitely a viable technique for obtaining accurate images of pore casts. 52 C) ’ C ) Fi gu re 2. 16 : Co ns ec ut iv e o pt ic al se ct io ns o f B er ea 10 0. Th e sp ac in g w as 9. 6 m ic ro ns be tw ee n th e se ct io ns . Fr om rig ht to le ft, th e im ag es ar e th e se co n d, th ird an d fin al (fo urt h) im ag es co lle ct ed . A BR 1 w as th e fir st o ft hi ss et (F igu re 2. 3). Th e br ig ht ar ea s in th e im ag es ar e th e po re s, an d th e da rk er ar ea s re pr es en t t he so lid . CHAPTER THREE Characterizing Porous Rocks with Nuclear Magnetic Resonance INTRODUCTION The characterization of porous media using Nuclear Magnetic Resonance (NMR) has become a focus of research as scientists seek to use it to address various questions in many different fields of endeavor. The initial application of NMR to describe the characteristics of porous media occurred in the oil industry. The Nuclear Magnetic Logging tool (NML) was used to correlate the NMR response of petroleum reservoirs with the physical properties of the rock (e.g. porosity, permeability, free fluid index, irreducible water saturation). Recently, the more indirect approach of determining pore sizes from NMR data and then calculating the physical properties from the pore sizes has been used. NMR has also been applied to the determination of size distributions in such diverse media as oil/water emulsions (Parker and Rees, 1972), biological cells (Brownstein and Tarr, 1979), and cheese (Callaghan et al., 1983). The objective of this study is to extend the uses of NMR in geology to the characterization of carbonate reservoir rocks. The suitability of the NMR method for determining pore size distributions of carbonate rocks will be investigated using the “two-fractions in fast exchange” model and the more general theory of Brownstein and Tarr (1979). In addition, the potential use of NMR to determine the fluid distributions in carbonates will be examined. The significance of this research derives from the fact that the microgeometry of the pore space and the pore fluids controls the petrophysical properties (i.e. permeability, electrical conductivity, dielectric constant) of porous media. These properties are used in the petroleum industry to interpret well-logging data and to estimate the capacity and production potential of oil reservoirs. The environmental industry uses the physical properties of rocks to predict the position and movement of groundwater and contaminants through the groundwater system. 54 HISTORY OF NMR IN PETROLEUM GEOLOGY Nuclear Magnetic Logging (NML) The basic principles of NMR were first published by Bloch (1946). The possible applications of this technique to the determination of reservoir properties were quickly recognized by the oil industry, and the first patent for a well logging instrument based on the principles of NMR was issued in the 1950s to the California Company (Jackson, 1984). The first nuclear magnetic logging tools (NML) were used in the 1960s. The industry expectations of determining fluid types, permeabilities, and the producible fraction of oil in place, were not met with these tools. However, more recent NMR logging tools do provide these data (Neuman and Brown, 1982). The unreliability of most of the NMR logging tools was due to difficulties with the magnetic field. Producing a strong, homogeneous magnetic field in the area suffounding the borehole where the electro-magnet must fit within the fifteen to twenty centimeter diameter of a well-logging tool was one of the problems. This problem has been addressed by the new technique invented in 1978 at the Los Alamos National Laboratory (Jackson, 1984). This technique is based on the creation of a toroidal region of homogeneous radial magnetic field surrounding the welibore at a specific distance from the center. The most recent development in the the field of nuclear magnetic well logging is the presentation of a prototype magnetic resonance imaging log, MRIL (Coates et al., 1991). The data which can be determined with this tool include a “sand-size” porosity, a bulk volume of irreducible porosity and a link to the permeability of the rock.” (Coates et al., 1991) 55 Relationships Between NMR and Physical Properties Empirical relations between the physical properties of rocks and NMR data have been suggested since 1956 when Brown and Fatt (1956) related the fractional wettability of oiffleld rocks to the free induction decay (FID) curve measured with an NMR apparatus. Correlations with the porosity, permeability, irreducible water saturation, the movable fluid in the reservoir and the oil/water saturations were also found (Timur, 1969, Robinson, 1974, and Vinegar et al., 1989). Pore Sizes The indirect approach of using NMR to determine pore sizes and then estimating the physical properties of the rocks was first proposed by Loren and Robinson (1970). The majority of research on this technique for geological materials has taken place in the last decade using spin- lattice (T1) laboratory measurements of sandstones and calculating pore sizes using the “two- fractions in fast exchange” model (Brownstein and Tarr, 1977, 1979). The calculation of pore sizes from Ti data using the “two-fractions in fast exchange” model is based on the assumption that the fluid in the pores is separated into two phases: a surface phase with a short relaxation time and a bulk phase with a substantially longer relaxation time. The majority of the nuclei are assumed to relax at the surface of the solid, an assumption which is accurate only if the pores are small enough so that all nuclei in the fluid can diffuse to the surface of the pore in a time less than the average relaxation time of the pore. If these assumptions are correct, then there is a one-to-one correlation between pore sizes and relaxation times. Various techniques have been used to calculate average relaxation times or distributions of relaxation times which accurately reproduce the collected NN’ER data. These methods include matching the data using a set number of relaxation times, using a “stretched exponential” representation (Borgia et a!., 1990 and Kenyon et al., 1988), and using inversion methods to produce discrete (Schmidt, 1986, Munn and Smith, 1987, and Gallegos et al., 1987, Whittall and 56 MacKay, 1989) or Continuous spectra of relaxation times (Gallegos and Smith, 1988, Gallegos et al., 1988, and Kenyon et al., 1989, Whittall and MacKay, 1989). The most recent research of Davies and co-workers (1990) has applied the theory of Brownstein and Tarr (1979) in which the fast diffusion assumption is not required. This study focussed primarily on sandstones, as did most of the others mentioned above. This predeliction to the study of sandstones is present because the pore network in sandstones is relatively uniform, the majority of oil reservoirs in the Unites States are sandstones, and the surface effect is more pronounced in sandstones than in carbonates (Loren and Robinson, 1970 and Schmidt et al., 1986). Fluid Distributions The microscopic distribution of fluid in porous media has a pronounced effect on their physical properties. This effect has been observed in laboratory measurements of compressional wave velocities (Domenico, 1976, Domenico, 1977, and Knight and Nolen-Hoeksema, 1990, seismic attenuation (Bourbie and Zinzner, 1984), electrical resistivity (Longeron et al., 1989 and Knight, 1991), and dielectric constant (Knight and Nur, 1987). Although NMR is a logical method for determining fluid distributions in sedimentary rocks, only two studies of partially saturated rocks have been completed with this goal in mind (Schmidt et al., 1986, and Davies et al, 1990). THEORY An excellent review of the basic theory of Nuclear Magnetic Resonance may be found in Experimental Pulse NMR. A Nuts and Bolts Approach (Fukushima and Roeder, 1981). The portions of this theory which are relevant to this study are summarized below. 57 Various nuclei (i.e: ‘H, ‘3C) possess an intrinsic magnetic field, or magnetic moment (M). This is created by the interaction between the angular momentum, or spin, of the nucleus and the electric charge associated with the atom. Let us consider nuclei which do not interact with each other, and thus represent isolated magnetic moments. These magnetic moments will be randomly oriented; however, if an external magnetic field (a3) is imposed, they will equilibrate to the field by assuming two possible lowest- energy positions: precession around the directions parallel and anti-parallel to the external field. The parallel orientation is favoured slightly by thennodynaniics, thus a cumulative magnetic moment (M0)parallel to the external field is produced. The precession frequency, also referred to as the Larmor frequency, is a function of the type of nucleus and the strength of the magnetic field. M0 is much weaker than the constant field, H0, therefore the change in the overall field cannot be accurately measured. However, if a radiofrequency (if) pulse operating at the Larmor frequency is applied perpendicular to H0, the nuclei will absorb energy from this field and rotate away from H0 and the rf field. The angle through which M0 rotates depends on the amount of time for which the rf pulse is applied. Thus, a 90° or 180° pulse is a radiofrequency pulse which causes M0 to rotate through 90° or 180° respectively. Immediately following a 90° pulse, M0 is perpendicular to both the rf field and the external field and will produce an emf field which subsequently induces a measurable voltage in the if coil. Through this process of first producing a coherent magnetic moment, M0, by imposing an external field, and then shifting M0 90° from the external field by applying a radiofrequency pulse, the magnitude of M0 can be measured. Nuclei in condensed matter, solids or liquids, are coupled to the degrees of freedom, referred to as the lattice, of that medium. This spin-lattice coupling tends to cause the nuclear moments and the lattice to evolve towards thermal equilibrium concurrently. Thus when the external field, H0, is applied, the evolution towards M0 requires that both the nuclei and the lattice reach thermal equilibrium. The progression towards M0 from the original components, M, M and M, is assumed to be exponential with two time constants: T, for the longitudinal component M, and T2 for the transverse components, M and M (Figure 3.1). 58 MM Figure 3.1: Nuclear Magnetic Resonance Reference Frame H0 M 59 Longitudinal Relaxation Time (T1) The longitudinal relaxation time is the exponential time constant representing the evolution of the longitudinal magnetic moment (Mi) of the nuclei into the moment M3 when an external magnetic field, H0, is applied to the media. The equation representing this progression is: -Mz=*(Mz-Mo) [1] The evolution towards M0 from M is accompanied by an exchange of energy between the lattice and the spinning nuclei and so T1 is often referred to as the spin-lattice relaxation time. The magnitude of the spin-lattice relaxation time T1 reflects the efficiency of the energy transfer between the lattice and the spinning nuclei, thus T1 is a function of the type of atoms present and their interaction with the lattice. The spin-lattice relaxation time is usually measured using the inversion-recovery pulse sequence which consists of a 180° pulse followed by a 90° pulse. The 1800 pulse flips the magnetization around so that it is equal and opposite to the cumulative magnetic moment, M0. The 90° pulse is applied a time (t) after the initial pulse and shifts that component of the magnetic field which is in the z direction to the axis perpendicular to the applied and rf fields so that its magnitude can be determined. This is followed by several seconds in which no pulse is applied and then the two pulses are applied again with the time separation varying with each application until the full decay curve has been outlined. 60 Transverse Relaxation Time (T2) The transverse relaxation time, T2, represents the gradual reduction of the transverse components of the magnetic moment to zero as the nuclei seek equilibrium with the external field, H0. This process occurs through the dephasing of the spinning nuclei relative to each other and is not accompanied by an energy exchange with the lattice. Thus it is commonly referred to as the spin-spin relaxation. The equation representing the decay of the transverse component of magnetization to zero is: f M,y = - ;j!:— [2] The spinning nuclei dephase after the removal of the if field due to variations in the magnetic field experienced by the nuclei. The field experienced by each nucleus is a function of the external field and the surroundings of the nuclei. The surroundings of the nuclei vary widely, therefore the field and the precession frequencies of the nuclei also vary. The primary difficulty associated with measuring the exact spin-spin relaxation time is that the imposed field, Ho, is usually inhomogeneous and this has a larger effect on the dephasing of the spins than does their immediate surroundings. This effect is routinely removed using the CPMG, or the Carr-Purcell-Meiboom-Gill, pulse sequence. This pulse sequence Consists of one 90° pulse followed by a sequence of 1800 pulses separated by time t. The initial 900 pulse sets up the magnetization perpendicular to the applied field and then, due to inhomogeneities in the field, some of the nuclei precess faster, some slower and the dephasing is quite rapid. However, one 180° pulse applied at time ton the x axis will reverse the magnetization along the Y axis and have no effect on the magnetization along the X axis. This pulse moves the magnetic moment of those nuclei which precess faster to a position behind the average and those which precess slower to a position ahead of the average. Therefore, after another time interval t, the magnetic moments will 61 re-align and the dephasing, which was not a result of the inhomogeneous magnetic field, may be measured. Spin-spin relaxation times are shorter in liquids confmed in solids; however, the data collected is a function of the pulse separation time, t, and data sets may only be compared if the effects of the t variation may be accounted for (Hedberg, 1990, and Kleinberg and Horsfield, 1990). In addition, an investigation into the relation between the spin-spin relaxation data and pore sizes by Kleinberg and Horsfield (1990) indicates that this relation is more complex than the corresponding relation between spin-lattice relaxation data and pore sizes. Therefore, T2 data were not analysed in this study. The NMR Response of Fluids Confined in Rocks The relaxation time of fluids confmed in porous media are much shorter than the times measured for bulk samples of the same fluid. For rocks, this fact was first demonstrated by Brown and Fatt (1956) and has since been confirmed by other researchers. This phenomenon has also been documented in such diverse media as animal tissues, biological cells, cheese, and wood lumens. (Bratton et al., 1965, Menon et al., 1987, and Callaghan et al, 1983) Many explanations of the increased relaxation rate/decreased relaxation time have been presented. They include: 1) the nuclei at the surface may be “bound” to the surface and thus act more like nuclei in solids, which have a shorter relaxation time than those in liquids 2) paramagnetic impurities on the pore surface enhance relaxation 3) the contrast between the magnetic susceptibility of the fluid and the solid sets up a magnetic gradient within the pores The relaxation of the nuclei of interest in this study, protons, has been modelled using diffusion and exchange mechanisms (Knispel, 1981, Brownstein and Tarr, 1977 and 1979). The diffusion model is the most widely accepted model and assumes that an increased relaxation rate at 62 the surface and the natural diffusion of the nuclei account for the increase in the relaxation rate of nuclei confmed in porous media. Brownstein and Tarr (1979) divided pore sizes into three regimes: the fast, intermediate, and slow diffusion regimes. In the fast diffusion regime, the pores are small and all of the nuclei diffuse to the surface quickly enough to produce an average relaxation time for the pore. In the slow diffusion regime, the nuclei farther from the surface relax at rates closer to the bulk relaxation rate of the fluid. Slow Diffusion The theory of nuclear magnetic relaxation of fluids within pores in the slow diffusion regime was developed by Brownstein and Tarr (1979). Pores within the slow diffusion regime meet the condition that paiD> 10.0, where p is the strength of the surface to relax the protons and is referred to as the the surface sink parameter, a is the radius of the spherical pore, and D is the self diffusion coefficient of the fluid. Theoretically, pores large enough to be within the slow diffusion regime are associated with an infinite sum of exponential decays with relaxation times T and amplitudes I which , for a spherical geometry, are given by: 2 12(sinC - Ccos) ln= [3] - sin(2C)] T = a2/D [4] where Cn are the positive roots of: 1-cot=paiD [5] 63 The expression for the observed magnetization M(t,R) at time t for a pore of size R is M(t,R) = M(O) Ie4fTn [6] where M(O) is the initial uniform magnetization and is a function of the quantity of protons present in the pore, not a function of R. The relaxation times T and the amplitudes I are functions of the size of the pore, the diffusivity of the fluid, D, the surface sink parameter, p , and the relaxation time of the bulk fluid, Tb Therefore a rock with a distribution of pore sizes, in which the pores are assumed to be approximately spherical, the T1 or T2 relaxation curve is the integration of a set of infinite sums of the form [6] and is given by: y(t) = j c(R) M(t,R) dR [7] where c(R) is the number of spheres of radius R, where e is the signal density and converts from number of spheres to signal intensity. The equations presented above are generally applicable to any type of medium as they require very few assumptions: they may be applied to pores of any size, and modified to account for most geometrical shapes. The primary assumptions associated with these equations are (1) that the pores are isolated, i.e. very few or no nuclei from one pore diffuse into another pore during the relaxation process, (2) the pores may be approximated by a simple geometric shape such as a sphere or cylinder, and (3) the surface sink parameter, p. is constant throughout the medium. 64 Fast Diffusion The fast diffusion regime is a special case of the slow diffusion equations presented above. if the pores are small enough that pa/D << 1.0, diffusion will ensure that all of the nuclei reach the pore surface and relax more quickly than they would independently. This results in an average relaxation rate throughout such pores so that each pore is associated with one relaxation time and the expression for the decay curve of a rock simplifies to: y(t) = J s(T)eti1ldT [8] where y(t1) is the data measurement at time t1, and s(T) is the amplitude of the component at relaxation time T. s(T) is equivalent to M(0) in equation [1] because it is the initial magnetization associated with relaxation time T. The assumptions associated with the fast diffusion equation presented above are slightly more resthctive than those of the slow diffusion equations. The assumptions are: (1) the pores are small enough so that pR/D << 1, (2) the pores are isolated, and (3) the surface sink parameter is constant throughout the medium. However, when pore sizes are considered, the benefit of the fast diffusion regime is that a general surface area to volume ratio of the pore , S/Vpore, may be determined. If required, this ratio may then be used to calculate a pore size given an assumption of shape. The one-to-one correlation between pore sizes and relaxation times in the fast diffusion regime allows the direct calculation of pore sizes from spin-spin or spin-lattice decay curves. The approach usually taken to this problem is to use the “two fractions in fast exchange” model first presented by Zimmerman and Brittin (1957) and further developed by Brownstein and Tarr (1977). This model assumes the pores are within the fast diffusion regime and that the fluid within the pore consists of two distinct phases: a bulk phase with a relatively long relaxation time (Tib) 65 and a surface phase with a substantially shorter relaxation time (Tis). The expression which relates the two responses to the measured relaxation time for each pore is: [9] T1 Tib Tis (Munn and Smith, 1987, Gallegos et aL, 1986) where T1, Tib, and Ti are the relaxation times of the pore, the bulk fluid, and the surface layer of fluid respectively. fa and fb are the fraction of bulk and surface water. The modifications of this relation which have been utilized by researchers are all very similar, ranging from Munn and Smith’s (1987) equation: ....L_.L 2 10T1 Tib Tlsurfrp whereT1is actually the relaxation time of the surface fluid divided by the thickness of the surface layer, 8, and rp is a mean pore radius defined as two times the ratio of pore volume to surface area,2(V/S)pore, to the expression employed by Schmidt et al. (1986): 1(1-SSIV) + os/V [11] T1 Tib Tis which does not neglect the fraction of surface water as the previous equation, [10], does. The expressions described above are all forms of the general equation for a spherical pore: [12 T1 V Tib For example, in equation [11] the equivalent form would be: 1_ 6S1 1 1 . _,1 1 v—— )+ with P— 11 V ls ‘lb ‘lb ‘is ‘lb 66 This equivalence is direct as p. 6 and T1s are unknowns which may only be estimated from other data. EXPERIMENTAL PROCEDURES Twenty-two carbonate samples were selected for use in this study. These samples were chosen to provide an adequate test of the NMR and imaging techniques in this thesis as well as to enable the determination of what properties of the samples affect the NMR response. The carbonates were selected to provide a range of pore sizes, with six of the samples having pores greater than 0.5 mm in length and three of these with very visible (>1mm) pores. Nineteen of the carbonates, referred to as S3 1- S49, were from a single formation, the Elkton Member of the Turner Valley Formation. This group of carbonates was chosen to ensure that any variation in the NMR results for these samples would be due to variations in the pore network rather than differences in the interaction between the protons and their surroundings. The other three samples, referred to as S54, S55 and S77, were also from the Turner Valley Formation. These samples were selected because they possessed large, easily visible pores which were not present in the samples from the Elkton Member. Sample Preparation The carbonate samples were prepared and the porosity and permeability measurements completed by Shell Canada Ltd. These rocks were one inch diameter plugs which had previously been cleaned using methanol extraction followed by drying in a vacuum oven to remove all hydrocarbons. The samples were trimmed to one inch in length and the ends ground smooth and parallel to each other. The one inch length was required in order to ensure that the entire sample was within the area in the NMR sample chamber in which the magnetic field was constant. 67 Sample Characterization Porosity and Permeability The porosity of the carbonate samples was measured using a dry weight/saturated weight comparison. The dry weight of the sample was measured after methanol extraction and vacuum drying to ensure the samples were completely dry. The volume of the plugs was determined by measuring the length and diameter of the samples with calipers, and the saturated weight of the sample was measured. This method assumes that the naptha penetrates all of the “effective” or interconnected porosity in the sample. The equation used to calculate the porosity was: (Ms-Md) “ D(irr2h) where M5 is the mass of the saturated sample, Md is the mass of the dry sample, r is the radius of the sample, h is the length of the sample, and D is the density of the naptha. The air permeability of the carbonate samples was measured by determining the pressure drop across a sample for a constant flow rate of air. Darcy’s Law was used to calculate the permeability using the equation: — -kpA(h-h1) where Q is the total discharge of fluid (air) per unit time A is the cross-sectional area of the flow path L is the length of the flow path P is the density of the fluid 11 is the dynamic fluid viscosity 68 h2-h1 is the pressure drop across the flow path, k is the permeability of the sample. The porosity and permeability of the samples are summarized in Table 3.1. Minera1ov The mineralogy of the carbonate samples was determined by a combination of hand sample and thin section examination. The thin sections were examined twice: once before and once after they were stained to aid in the separation of the carbonate constituents. The stain used was made up of a mixture of hydrochloric acid, Alizarin Red S, and potassium ferricyanide as described by Lewis (1984) and Evamy (1969). Alizarin Red S and potassium ferricyanide solutions may be applied separately; however, calcite and dolomite cannot be distinguished without using both. Table 3.2 indicates the stains which result when solutions of these two chemicals are used alone and when they are combined. The process of staining required that the thin sections were clean and dry, the chemicals mixed with deionized water, and the thin sections immersed in the solutions for one to two minutes. The stained thin sections were dried using compressed air to ensure that the stain did not become streaky, that a residue of the stain was not left on the section, and that the section dried rapidly. The use of potassium ferricyanide requires several precautions as it is photo sensitive and will react with tap water. Therefore only deionized water was used and, if the solution was not going to be used within 48 hours, it was stored in a dark bottle away from the light. 69 TABLE 3.1 Porosity and Permeability S3 1 S32 S33 S34 S35 S36 S37 S38 S39 S40 S41 S42 S43 S44 S45 S46 S47 S48 S49 S54 S55 S77 0.038 0.054 0.107 0.076 0.094 0.068 0.041 0.028 0.055 0.053 0.019 0.063 0.06 1 0.060 0.059 0.049 0.066 0.049 0.042 0.057 0.048 0.082 0.050 0.070 0.550 0.410 0.130 0.060 0.040 0.050 0.060 0.050 0.040 0.050 0.050 0.050 0.050 0.050 0.050 0.050 0.050 0.070 0.060 5.43 SAMPLE POROSiTY PERMEABILITY mD 70 TA B LE 3. 2 St ai ns fo r D iff er en tia tin g B et w ee n C al ci te a n d D ol or ni te * S T A IN IN G C A L C IT E D O L O M IT E S IL IC A R E A G E N T S Co m po sit io ns ar e gi ve n in Fe ++ fre e F e po or F e ric h Fe ++ fre e Fe ++ < 1 F e > 1 w ei gh ts pe rc en t. Cn tic al M g M g Q U I{ ’1 so lu tio n st re ng th s ar e CH A LC ED O N Y u n de rli ne d. CA LC IT E D OL OM IT E FE R R O A N A N K ER IT E O PA L FE R RO A N CA LC IT E D O LO M IT E se n su st ic tO se n su st ric to 0.2 %h yd roc hlo ric ac id R ed R ed R ed N ot N ot N ot N ot St ai ne d 0. 2% A liz ar in Re d S St ai ne d St ai ne d St ai ne d Q, 2% hv dro ch lor ic ac id N ot L ig ht D ar k N ot L ig ht D ar k N ot St ai ne d 0. 2% po ia ss iu m fe rri cy an id e St ai ne d B lu e B lu e St ai ne d B lu e B lu e 0. 2% hy dr oc hl or ic ac id R ed M au ve V io le t N ot L ig ht D ar k N ot St ai ne d 0. 2% A liz ar in Re d S 0. 2% po ta ss iu m fe rri cy an id e S ta in ed B lu e B lu e *a fte r Ev am y, 19 69 Samples S31 to S49 These nineteen samples are from the Mississippian Elkton Member of the Turner Valley Formation. This formation has previously been described as consisting of three sub-members: a lower member of porous, coarsely crystalline dolomite which grades into less porous limestone in some areas, a middle member of dolomitic shale or argillaceous and silty dolomite, and an upper member of finely to coarsely crystalline dolomite containing varying amounts of chert (Penner, D.G., 1957). This indicates that there is a wide variation in the mineralogy and physical properties of rocks from this formation; however, from the literature it would be expected that the rocks would be at least partially dolomitized, the fossils present will be dominated by crinoids and that other fossils such as bryozoans, brachiopods, and relict texture suggesting pelletoids may be present (Penner, D.G,. 1957, Young, F.G., 1967, and Thomas and Glaister, 1960). The hand sample and thin section examinations revealed the samples were primarily matrix supported limestones which were partially dolomitized, the dolomite being present as small rhombs (<0.2 mm) within the matrix and in partially replaced fossils (Figure 3.2). A small component of shale (<10%) was also present and the samples were gray, with little to no visible porosity. What little porosity was visible was present in samples S31, S32 and S33. The fossils observed consisted primarily of crinoid ossicles and stems (Figure 3.3) with a small number of bryozoans and brachiopods (Figure 3.4). Throughout the nineteen samples, only one area with a relict texture of pelletoids was viewed. A general description of these samples is: grey, low porosity, partially dolomitized, matrix supported crinoid limestones. Samples S54, S55, and S77 These three limestones were from the Mississippian Turner Valley Formation. They were cream coloured and contained a great deal of visible porosity. Pores up to 1 mm in diameter were measured in samples S54 and S55 and 2 to 3 mm in diameter in sample S77. The samples were 72 Figure 3.2: Partial dolomitization in Ellcton Member samples, S31 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. 73 • ... ... - -.. .1’ 4 . •••: •.-. 7 A. •0•- I. - - S 0.0 0.25 mm • ••••-.;,• -.•..• • • •. 1.. ;* Figure 3.3: Examples of crinoid stems and ossicles in Elkton Member samples, S3 1 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. Lower photo contains a cross section of a crinoid stem; examples of crinoid particles in the upper photo are indicated by C ‘s. 74 Figure 3.4: Examples of bryozoan and brachiopods which occur in Elkton Member samples, S3 1 to S49. Dolomite is blue and usually rhombic, calcite is red or pink, and extremely small crystals consisting of carbonate and silicate minerals appear as dark areas. 75 primarily coarsely crystalline limestone and contained less than one percent silica. S54 and S55 contained little or no dolomite, while S77 contained approximately 10% dolomite. DATA COLLECTION Sample Preparation Fully Saturated Samples The samples were dried by heating them in a vacuum oven for several hours to remove the naptha after the porosity was measured. The samples were then immersed in 100,000 ppm NaCl for several hours at room temperature and pressure, a vacuum was applied for several hours, and finally, the samples and fluid were transferred to a pressure vessel and 600 psi of Nitrogen was applied for approximately fifteen minutes. Following saturation with the brine, the samples were immersed in 100,000 ppm NaCl within a sealed dessicator for storage. The dessicator was kept in a heat bath which was maintained at 40°C for at least one-half hour prior to any NMR measurements. The samples were heated so they would be in equilibrium with the NMR sample chamber which is maintained at 40°C. Spin-lattice relaxation times vary with temperature, soit was important to ensure the temperature of the sample was constant throughout the experiment. The sample to be measured was removed from the dessicator using tongs, surface dried to remove excess water, and then wrapped in 2 layers of cling film. Teflon tape was wrapped around the outside of the cling film in two alternating layers. The layering of cling film and teflon tape was used to ensure that a minimum of fluid loss from the sample during the measurements. A teflon filler surrounded the samples within the NIVIR sample holder. The teflon filler, which does not produce an NMR signal, was required to keep the samples centered within the holder which was approximately 3.8 centimeters in diameter, 1.27 centimeters larger than the samples. The 76 sample holder was a glass, square-ended test tube which was placed within the sample chamber of the NMR machine to collect both T1 and T2 decay curves. Partially Saturated Samples Six of the samples employed in the experiment on the fully saturated samples were cleaned to remove the brine and re-saturated with deionized water. The brine was removed by soaking the samples in deionized water for a period of several days, with the deionized water being replaced with fresh water daily. The saturation with de-ionized water was accomplished in the same way they were previously saturated with brine. The samples were saturated with deionized water because the variation in saturation was accomplished using evaporative drying, which would cause an increase in the salinity of the pore fluid if a brine were used. Evaporative drying was accomplished by placing samples in a dessicator with anhydrous CaSO4after the initial NMR measurements of the samples in a fully saturated state. Prior to the fully saturated measurement, the samples were stored in a sealed vacuum dessicator within a heat bath which was turned on for at least thirty minutes prior to the measurement of the NMR decay curves of a sample. The partially saturated samples were weighed before and after each NMR measurement to determine their saturation level. The preparation for the NMR measurement was the same as for the fully saturated samples: if necessary, the sample was surface dried, and then two layers of cling wrap and two of teflon tape were wrapped around the sample. The cling wrap and teflon tape reduced fluid loss from the samples both during the NMR measurement and as they were heated to 40°C. The temperature of the sample was increased to 40°C by placing it in a sealed NMR sample holder and placing the holder in a 40°C heat bath for at least twenty minutes prior to the NMR measurements. 77 NMR Measurements The NMR measurements were performed with a Bruker Minispec PC-20 NMR Spectrometer with a magnetic field strength of 0.47 T. The pulse sequences and parameters used to collect each data set were determined by two different EDMs (Experiment Definition Modules) and several calibration tests were performed regularly to ensure the NMR machine was operating properly. The homogeneity of the magnetic field, the diode offset and the phase of the phase sensitive detector were checked at least once a week. To ensure the highest possible signal to noise ratio, the built-in bandwidth filter was set to the lowest possible value of 150 kHz. Spin-Lattice Decay Curves The spin-lattice relaxation curves were collected using a Bruker EDM 510 which was programmed with the “inversion recovery” sequence described previously. A maximum of 160 measurement times were possible, with the time of each measurement calculated with the equation: t=CT1*CF CF varies from 0.04 to 5.12 and is determined by multiplying the previous value by: (n-i) — Vo.oi. where n is the number of points on the curve and the initial value of CF is 0.04. Although the exact time calculated in this manner is presented as the measurement time, the actual measurement times of the Minispec PC-20 occur at whole number millisecond times, thus if the Minispec claims to have measured the magnetization at 2.2 or 2.8 ms, the measurement was actually performed at 2 ms. CT1 is a constant which is an estimate of the average T1 of the sample, thus the time between 78 two “inversion recovery” pulse sequences is set to 5 times this value to ensure that the nuclei have time to relax completely between measurements. Between 3 and 45 measurements were taken at each time and stacked to produce the final result. The measurements were completed in approximately one hour for all samples since for high porosity samples, the data were clearer and the same data quality was obtained from a smaller number of measurements. Spin-Spin Decay Curves Experiment Defmition Modules 610 and 612 were used to collect the T2 decay curves of these samples. Both EDM’s employ the CPMG pulse sequence described in the theory, however, 610 measures the magnetic field strength after every tenth echo while 612 measures after every second echo. Thus, 610 is suited to samples with long T2s and 612 to samples with short T2s. Unfortunately, only 50 measurements of the curve may be taken with the 612 EDM while EDM 610 will collect up to 160 measurements. The nuclei must be allowed to reach equilibrium with the external field between each spin- echo experiment, thus the constant mentioned above, CT1, is used to ensure that 5 times the approximate spin-lattice relaxation time separates each spin-echo experiment. The time over which the curve is measured is set with the CT2 constant, which should be approximately equal to the actual T2 of the sample. The maximum time at which data is collected is greater than 2*CT2, ensuring that 90% of the decay curve with a time constant of CT2 would be collected. .CT2 was adjusted to collect the curve at its closest approach to zero; however, it was maintained at a high enough value so that noise effects did not cause the signal to decrease below zero volts. DATA ANALYSIS The spin-lattice relaxation data collected for this study were analysed to determine the porosity, pore sizes, and T1 distributions of the samples. Pore size distributions were calculated 79 using two methods: linear inversion of the relaxation data using the slow diffusion model, and direct calculation from the T1 distributions using a relation, equation [12], developed using the “two fractions in fast exchange” model. Data Inversion Spin-Lattice Relaxation Time (T1)Disthbutions The linear inversion methods presented by Whittall and MacKay (1989) were used to determine both continuous and discrete T1 distributions. These inversions assumed, as is the case in either fast or slow diffusion, that the decay curve produced by a fluid saturated rock will be multiexponential in form. Thus, the equation used to determine the appropriate distribution of Ti’s was equation [8]: y(t1) = J s(T)etiul’dT [13]a where y(t) is the datum at measurement time t1, and s(T) is the amplitude of the component at relaxation time T. The linear inversion method assumes a number of known relaxation times Tj and solves for the associated amplitudes s(T) or sj. In order to execute the inversion with a computer, the equation above was discretized assuming the model was a sum of M delta functions with areas sj at known relaxation times T and the resulting expressions are: s(T) = s(T-T) [14] 80 M y1=sjgjj i=1,...,N j= 1 where the kernel, gjj, is given by: fTj+l gjj = etiIJdT JT The non-negative least squares (NNLS) algorithm selected to invert the data would provide the closest possible fit to the data by calculating the solution which minimized the term: E I sjgjj - i=1 j=1 However, experimental data always contains a certain amount of noise, therefore the preferred solution should not fit the data exactly. Solutions to the inversion should misfit the data by an amount corresponding to the error in the data, which is represented by the standard deviation, a. 2This is accomplished usmg the X statistic N = (y - y)Ia [15] to quantify the misfit between the measured data, y, and the data predicted by the constructed model, y. The expected value of is equal to the number of data points, N, corresponding to a misfit of approximately one standard deviation for each data point. The NNLS algorithm employed the criteria as a measure of the appropriateness of the model, s(T), produced. Two approaches were implemented: one which assumed that the amplitudes were delta functions and produced a discrete spectrum as described above, and one 81 which incorporated additional constraints in order to produce a piece-wise continuous distribution. The piece-wise continuous distribution was obtained by calculating the solution which minimized the term: sjgjjyj2 +p.s? i=1 j=1 j=1 in which u is a trade-off parameter which, if 0, allows the NNLS solution to be calculated, and as it becomes larger, places more emphasis on a solution which minimizes the “energy” of the system. Minimizing the second term, referred to as the energy of the system, corresponds to minimizing the structure in the solution. Pore Sizes - Slow Diffusion Regime The NNLS inversion techniques described above were also employed to determine pore sizes based on the theory of Brownstein and Tarr (1979). The pores were assumed to be spherical in shape and the inversion equation, [7], was: y(t) = j c(R) M(t,R) dR The application of this equation was accomplished by discretizing the problem, providing a set range and number of pore radii as a basis for the inversion, and removing the constants from the equation and solving a simplified equation: y(t) = f c(R) R3 M(t,R) dR Jo 82 These radii and equations [3], [4], and [5] were used to calculate the distribution of amplitudes and relaxation times associated with the selected pore sizes. These amplitude/relaxation time distributions were then utilized to produce a best fit model to the NMR decay curve by determining the number of pores of each pore size present in the sample. The development and testing of the computer algorithm which implemented this inversion is described in further detail by Whittall (1991). Porosity Determination The porosity of the samples were determined based on a basic principle of NMR: the magnitude of the initial magnetization, M(O), is proportional to the number of protons present in the sample. Therefore, if a sample is fully saturated, the “effective” porosity can be determined by comparing the initial response of the sample to that of a fluid sample of the same volume. Similarly, the saturation of a partially saturated sample can be determined by comparing its response with that of the fully saturated state. The magnitude of M(O) can be determined from the T1 distributions of the samples by substituting t =0 into equation [13], producing the expression: y(O)=f s(T)dT [16] Jo This indicates that the magnitude of M(0) = y(O), and thus the relative quantity of fluid, may be calculated by integrating over the T1 distribution determined for a fully or partially saturated rock sample. 83 Pore Sizes - Fast Diffusion Regime The pores sizes of the samples were calculated using the “two fractions in fast exchange” model by using equation [12] _L T1 V T, Throughout the remainder of this thesis, the calculation of pore sizes from the T1 distribution utilizing this equation will be referred to as the fast diffusion method. RESULTS Fully Saturated Samples I1Distributions Continuous and discrete T1 distributions were produced for each fully saturated sample using the non-negative least squares inversion method described earlier. Graphs of these distributions are included in Appendix 6. The discrete relaxation spectra for samples S34 to S49 are very similar. As illustrated by the example shown in Figure 3.5a, most of the fluid in the pores had relaxation times just under 10 ms and the maximum relaxation time ranged from 60 to 200 ms. The inversion of data which produced a non-zero amplitude at relaxation times greater than 90 ms was difficult. Regardless of the range of relaxation times chosen to invert the data, a peak was always present at the maximum time chosen. 84 Samples S3l to S33 were slightly different than the others from the Elkton Member, as shown by the S3 1 spectra in Figure 3.6a. The maximum amount of fluid in the pores had a relaxation time between 40 and 60 ms and, for S31 and S32, a substantial amount of the fluid (>20%) possessed a relaxation time of 300 and 200 ms respectively. The continuous relaxation spectra for the samples from the Elkton Member are exemplefied by the S39 spectra in Figure 3.5b. These spectra all reached a maximum at 20 ms or less; however S31,S32 and S33 were represented by slightly larger T1 values and attained their maximum amplitude between 25 and 80 ms. The increase in relaxation times for S31, S32 and S33 relative to the other Elkton Member samples is illustrated by the continuous relaxation spectra of S31 shown in Figure 3.6b. The discrete and continuous spectra of S54,S55 and S77, of which the S77 spectra in Figure 3.7a and 3.7b are typical, are represented by relaxation times close to a full second. The discrete spectra were characterized by maximum amplitudes at the maximum relaxation times represented: 1750, 1000, and 1250 ms for S54, S55 and S77 respectively. The continuous spectra for these samples were similarly biased towards the longer relaxation times, with the maximum in the T1 curve occurring at 800, 750 and 1000 ms for S54, S55 and S77 respectively. Porosity The total integral of the T1 spectra, M(0), of each sample should be related to the porosity of that sample, as shown by equation [15]. This correlation is shown in Figure 3.8 for both the Continuous and discrete relaxation time spectra. The best fit straight line for the graph of porosity vs M(0) has a regression parameter of 0.728 for the continuous T1 spectra and 0.883 for the discrete spectra. 85 Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) 0.6 0.5 0.4 0.3 0.2 0.1 0.0= 1 10 100 1000 I t) Spin-Lattice Relaxation Times (ms) (a) Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.00 ...I 1 10 100 1000 (b) Figure 3.5: Sample S39 Ti distributions.(a) discrete and (b) continuous. 86 Discrete Relaxation Spectrum 0.4 j 0.1• 0.0 1 10 100 1000 Spin-Lattice Relaxation Times (ms) (a) Continuous Relaxation Spectrum 0.08 1 10 100 1000 Spin-Lattice Relaxation Times (ms) (b) Figure 3.6: Sample S31 Ti spectra.(a) discrete and (b) continuous. 87 Discrete Relaxation Spectrum (a) Continuous Relaxation Spectrum 0.08 10 100 1000 10000 (b) Figure 3.7: Sample S77 Ti spectra. (a) discrete and (b) continuous, interpolate a coefficient of 2.25x i0- m2/s at 40°C, the temperature of the sample chamber in the NMR spectrometer used to collect the data for this study. a) a) a) 0.8 0.6 0.4 0.2 0.0 1 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) 0.06’ 0.04’ 0.02 0.00 ‘I ...I Spin-Lattice Relaxation Times (ms) 88 Correlation of Porosity with M(O) Continuous Ti Spectra 12000 y= .-163.36+9.4015e+4x RA2=0728 8000 M(O) 60002000 EEEEEZEZEZ4000 0.00 0.02 0.04 0.06 0.08 0.10 0.12 Porosity Correlation of Porosity with M(O) Discrete Ti Spectra 12000- y = 834.74 + 1.0175e+5x R”2 0 883 10000 8000 M(O) 6000 4000 2000 0- • • I • I • I • 0.00 0.02 0.04 0.06 0.08 0.10 0.12 Porosity Figure 3.8: Coffelation of the integral over the the discrete and continuous Ti distributions (M(0)) with measured porosity. 89 Pore Sizes The pore sizes calculated using both the slow and fast diffusion methods require the prior knowledge or estimation of three parameters: D, the self-diffusion coefficient of the bulk liquid; p, the surface sink parameter; and Tib, the relaxation time of the bulk fluid. The surface sink parameter, p, may only be determined through comparison of the calculated pore size spectra with a previously determined pore size distribution; however, both D and Tib may by measured. The relaxation spectrum for the 100000 ppm NaC1 brine used to saturate all of the carbonate samples is shown in Figure 3.9. Tib was determined to be 2.0 s for water which had not been de gassed to remove the paramagnetic oxygen, a value which compares favourably to other measurements (Schmidt et al., 1986, Endom et al., 1967). The self-diffusion coefficient may only be measured with a magnetic gradient, which was not available, therefore the values determined by Endom et al. (1967) for water molecules in 100,000 ppm NaC1 at 0, 25, 50 and 80°C were used to interpolate a coefficient of 2.25xiO m2/s at 40°C. Fast Diffusion The fast diffusion method was applied to calculate those pore radii which corresponded to the relaxation times determined with the NNLS inversion (Appendix 6). The surface sink parameter, p. was calculated by determining the value required to equate the average T1 of the discrete and continuous S55 and S77 spectra with the average pore size determined with the confocal laser scanning microscope. The values calculated were lx10 m/s for three of the calculations, and 2x105mis in one case, therefore the most common value, lxl&5mis, was selected. When the T1 distributions were altered to pore size distributions, the discrete spectra were neglected because the porosity in most rocks is better represented by a continuous distribution. 90 Bulk Relaxation Time of 100000 ppm NaC1 1.0• 0.8 4 1.) 0.2 - 0.0- •__. 1 10 Spin-Lattice Relaxation Time - Ti (s) Figure 3.9: 100,000 ppm NaC1 brine relaxation time. 91 The two pore size ranges determined for samples S3 1 to S49 using the fast diffusion pore size calculations and the relation between the original T1 spectra and the final pore size spectra are illustrated by the T1 and pore size spectra of S39 and the pore size spectra of S31 shown in Figure 3.10. The most striking result for samples S31 to S49 was the direct alteration of relaxation times to pore sizes with no change in the overall shape of the distribution. The maximum amplitude of the disthbutions covered two size ranges: between 0.3 and 0.7 p.m for S33-S49 , and 1.0 to 3.0 p.m for S31 and S32. The T1 distributions of S54, S55 and S77 contained much larger relaxation times, therefore the resulting pore size distribution was “stretched” as the relaxation times in the original distribution began to approach the relaxation time of the bulk fluid (Figure 3.11). The maximum amplitude of these spectra occurred between 40 and 60 p.m. sizes ten to twenty times larger than the radii at which the maximum amplitudes of S31 and S32 occurred. Slow Diffusion The ability of the slow-diffusion inversion method to reproduce the correct pore size distributions was tested by producing a theoretical spin-lattice decay curve based on a theoretical pore size distribution similar to that of sample S77. The times at which the theoretical “measurements” were calculated were varied by matching the actual measurement times of the sample S77 data ( 38 to 4863 ms) and by expanding this to include a much larger range of measurement times, from 1 to 22000 ms. The inversion of the theoretical data matched the theoretical distribution more accurately when the time range of the data was expanded, as illustrated in Figure 3.12. The theoretical data was also altered by adding 1% and 5 % noise and the subsequent inversion produced pore size spectra which matched the actual distribution reasonably well. Figure 3.13 illustrates the loss of detail and accuracy with increasing noise through comparison plots of 92 Ti and Pore Size Distributions 0.08 Pore Sizes Ti 0.06 . 0.04 0.02 0.00• •..—.. ..‘ —I I ••• ... i08 io 10.6 i05 i0 i0 10.2 10i io° Ti (ms) and Pore Radius (m) (a) Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.00 108 106 Pore Radius (m) (b) Figure 3.10: (a) Sample S39 Ti and fast diffusion pore size distributions. (b) S31 fast diffusion pore size distribution. 93 Ti and Fast Diffusion Pore Size Spectra 0.08 0 Pore Sizes I. Ti 0.04’ I 0.02 o.oo —. i0 i0 io io’ io i02 10’ 100 iOi Ti (s) and Pore Radius (m) (a) Pore Size Spectra Fast Diffusion Method 0.06 0.05 C .— 0.04 tl) 0.03 C 0.01W 0.00• 106 i- Pore Radius (m) (b) Figure 3.11: (a) Sample S77 Ti and fast diffusion pore size distributions. (b) S55 fast diffusion pore size distribution. 94 Slow Diffusion and Confocal Pore Size Spectra 0.4 Theoretical Pore Size Spectrum Slow Diffusion Pore Size Spectrum Theoretical data range: 38-4863 ms 0.3 0 .— I. I 0.2 $ i k _ I’ II o Ij Ii I 0.1 _________ iLl . io-3 Pore Radius (m) Slow Diffusion and Confocal Pore Size Spectra 0.25 Theoretical Pore Size Spectra 0.20 Slow Diffusion Pore Size Spectra Theoretical data range: 1-22000 ms .. II 0.15- : II :: : 0.05 0.00 - — —f •—i—’.! 1 io5 io io Pore Radius (m) Figure 3.12: Slow diffusion method results of inverting theoretical data produced using a theoretical pore size spectrum similar to that of S77. Variation in the pore sizes is a result of using two time ranges to calculate the theoretical decay curve. The upper graph was produced with data over the same time range as the data for S77 was collected, the lower with an extended range. 95 Slow Diffusion and Confocal Pore Size Spectra 0.4 0.3 _____ Slow Diffusion Pore Size Spectrum Theoretical data range: 38-4863 ms Slow Diffusion Pore Size Spectrum 0.2 Theoretical data range: 38-4863ms Noise: 1% II I ________ __ io io3 Pore Radius (m) Slow Diffusion and Confocal Pore Size Spectra 0.6 Slow Diffusion Pore Size Spectrum Theoretical data range: 3 8-4863 ms 0.4 II Slow Diffusion Pore Size Spectrum Theoretical data range: 38-4863ms II Noise: 5% • 0.3 S 0.2 0.1 Ri i ______ ____ oo —- — io4 Pore Radius (m) Figure 3.13: Slow diffusion method results of inverting theoretical data based on a theoretical pore size spectrum similar to that of S77. The pore sizes resulting when 1 and 5% noise are added to the theoretical data are shown in the upper and lower graphs respectively. 96 the pore size spectra calculated from the perfect and noisy theoretical data. The pore size spectra determined from noisy data are characterized by the appearance or enhancement of the amplitude associated with the largest pore size and the representation of the pore size spectra by two or three pore sizes rather than a continuum of pore sizes. If the component at the largest pore size was included, the average pore size determined from the data with 1% noise was 153 pm, and from the data with 5% noise was 216 pm. This illustrates the sensitivity of the problem to additional noise, as 153 .tm is within 20% of the actual average of 122 p.m, while 216 p.m is almost double this expected value. The results when the final pore size is discounted are slightly more promising, with average pore sizes of 106.1 and 82.3 p.m for the data with 1% and 5% noise respectively. The variation in the X2 parameter with p was tested with the theoretical data. Whittall 2(1991) found that the minimum mX tends to underestimate the actual parameter, particularly for pore sizes in the slow diffusion regime. For the theoretical data produced for this study, the inversion of theoretical data without added noise produced a minimum X2 when p matched the value the data was produced with. When even 1% noise was added to this data, this minima was either altered from an obvious to a broad minimum, or, alternatively, the value decreased steadily as p was decreased. When applied to both the real and theoretical data, the slow diffusion method consistently produced spectra which more closely resembled discrete than continuous pore size distributions. Typical results for the theoretical data are shown in Figures 3.13, while examples of the spectra calculated from real data are shown in Figures 3.14 and 3.15. The surface sink parameter used for each real data set was varied and the situations described for the theoretical data, in which a broad minimum in or a constant decrease in with p occured, were also observed for the real data. In all cases, either the p value just above the minimum, a p value slightly smaller than that at which increased substantially or, if corroborating data was available, the p value which best fit the known pore size distribution was selected. The surface sink parameter for samples S3 1 to S49 ranged from 1.5 to 3.75x 10-5 mIs, with an average of 3.0 x105 ni/s. The minimum occured at 97 p values of 2x105 rn/s for S54 and S55, and 3.5x10 rn/s for S77. However, surface sink parameters of lx iO rn/s and 1x10 rn/s fit the confocal pore size distribution more accurately for S55 and S77 respectively and therefore these values were used for these two samples. The pore size spectra produced for all of the data sets possess a peak at the largest pore size included in the inversion range. As shown in S31 and S43 pore size spectra in Figures 3.14 and the S55 and S77 spectra in Figure 3.15, this peak represents from 15 to 25 percent of the integral over the complete spectrum for samples S3 1 to S49 and between 50 and 60 percent of the complete spectra for samples S54, S55 and S77. The pore size associated with the maximum amplitude of the S33 to S49 distributions ranged from 0.8 to 1.5 microns with several samples having a second peak at pore sizes between 1.5 and 6 microns, as illustrated by sample S43’s pore size spectra in Figure 3.14a. The spectra produced for S3 1 and S32 were represented by larger pores than the other samples from the Elkton Member. This increase in pore sizes relative to the S33 to S49 results is typified by the pore size spectra of S31 in Figure 3.14; the maximum amplitude and secondary peaks occur at 2.5 and 10 p.m for S31 and at 4 and 15 p.m for S32. Two peaks other than the one present at the largest pore size were produced in the distribution of each of the other three samples, occurring at 3 and 9 p.m for S54, at 1,5.2, 6.5, 37, and 43 J.tm for S55, and at 10 ,50, and 59 p.m for S77. These results are illustrated by the pore size spectra of S55 and S77 included in Figure 3.15. The pore size spectra described above are non-unique solutions to a complex problem, and, since the pore size and amplitude of the distribution are affected by the parameters required for the inversion, the variation in the spectra with D, p and Tib was determined. Several sets of data from all of the samples were tested and samples S3 1 to S49 were found to have similar responses, as did samples S54, S55 and S77. An example of the variation in the pore size spectra for the two types of responses to Tib is shown in Figure 3.16. The spin-lattice relaxation time has only a very small effect on the amplitude of the spectra for S3 1 to S49, but it causes an increase in pore sizes 98 Pore Size Spectrum Slow Diffusion Method 0.3 0 0.2 I 0.0 •, . A., i06 io5 Pore Radius (m) (a) Pore Size Spectrum Slow Diffusion Method 0 .— 0.2 I 0.1• 0.0• II i0 106 io5 1o Pore Radius (m) (b) Figure 3.14: Slow diffusion pore size distributions. (a) S43 (b) S31 99 Pore Size Spectrum Slow Diffusion Method 0.6 0.5 C 0.4 j ::: ::z-iA.. 106 io Pore Radius (m) (a) Pore Size Spectrum Slow Diffusion Method 0.8 C .— 0 a) 0.4 0 •1 Pore Radius (m) (b) Figure 3.15: Slow diffusion pore size distributions. (a) S55 (b) S77 100 0 C) 0 Figure 3.16: Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the bulk fluid relaxation time, Tib, is varied.as it is decreased for samples S54, S55 and S77. Pore Size Variation with the Bulk Relaxation Time Sample S39 2 0 I. 1• 4.Os •2.Os - 1.0s 0• 0.48 1.69 I 1.98 Pore Size 0.8 0.6 2.31 45.63 Pore Radius (m) Variation with the Bulk Relaxation Time Sample S55 4.Os 2.Os • 0.4 02 0.0- Pore Radius (m) 101 The effects of varying D are illustrated in Figure 3.17 by the variation in pore sizes for samples S39 and S55. For samples S31 to S49, the variation with D was not direct since the pore sizes did not change according to a pattern, and the overall average pore size remained relatively stable. However, as D was increased for samples S54, S55 and S77, the pore sizes also increased. The surface sink parameter has a defmite effect on the pore size spectra determined for all of the samples; as shown in Figure 3.18, a larger surface sink parameter produces a larger pore size. The amplitude of the peaks vary as well, with the amplitude at larger pore sizes increasing with the surface sink parameter. Partially Saturated Samples Fluid Distributions The Ti spectra corresponding to varying saturations of deionized water for six samples (S31,S41,S42,S44, S54 and S55) were determined and are included in Appendix 7. As shown in the S54 spectra in Figure 3.19, the T1 spectra of the six samples generally decreased in amplitude and shifted to the left towards smaller T1s as the saturation of the samples was decreased. The shift towards higher T1swhich was observed in varying degrees at the lowest saturations of samples S41, S42, S44, and S55, is illustrated by the T1 spectra of sample S44 in Figure 3.20. Fluid Saturations The integral over the T1 distributions (M(0)) was correlated with the saturation of the samples and, as shown in Figure 3.21 for samples S31 and S41, the relationship could be approximated by a straight line. Graphs of M(0) vs saturation for all six samples are in Appendix 7 and the regression parameter indicating how well the data fit a straight line ranged from 0.92 1 to 0.993. 1 02 C 0 a.) E C > Pore Radius (m) Pore Size Variation with the Self-Diffusion Coefficient Pore Radius (m) 00 C’ oOC r C — — C4 \O Figure 3.17: Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the self diffusion coefficent of the bulk fluid, D, is varied. 1 .2 1.0 -I 0.8-I 0.6-I Pore Size Variation with the Self-Diffusion Coefficient Sample S39 J 4.5e-9m/s 2.25e-9 m2/s R 1.125e-9m/s I ___ 0.4 -, 0.2-’ 0.48 1.69 1.98 2.31 45.63 Sample S55 C 0 a) — C 1.5 1.0 :: I—F 4.5e-9m21s 2.25e-9 m I 1.125e-9m/s 1 03 Pore Size Variation with the Surface Sink Parameter 0 j 0 0 E 0 1.0 Sample S39 0.0 -, 00 .: c cn c’i Radius (m) 7.5e-5 rn/s 3.75e-5 rn/s • 00 \C t o o o Pore Pore Size Variation with the Surface Sink Parameter Sample S55 2 4e-4 rn/s 1 - 2e-4 rn/s le-4 rn/s 0- ..-----.—--. • — OcCC00O00 Pore Radius (m) Figure 3.18: Examples of the variation in the pore sizes determined using the slow diffusion method when the value of the surface sink parameter, p, is varied. 1 04 Sa m pl e S5 4 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns Sa tu ra tio ns 0. 99 5 a 0. 92 0 — 0- — — 0. 88 6 — 1. — — — 0. 79 4 0. 73 0 — 0 -- — 0. 63 2 0. 55 7 0. 46 4 0. 36 2 - 1- 0. 23 7 0 (7 ’ Sp in -L at tic e Re la xa tio n Ti m e (m s) Fi gu re 3. 19 : Ex am pl eo f t he de cr ea se in am pl itu de an d re la xa tio n tim e w ith sa tu ra tio n fo rp ar tia lly sa tu ra te d sa m pl es 0. 08 0. 06 0. 04 0. 02 0. 00 10 10 0 10 00 10 00 0 Sa m pl e S4 4 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns 0. 06 0. 05 Sa tu ra tio n 0. 04 0. 91 4 0. 81 2 0. 03 - - - 0. 71 1 — 0 — — 0. 61 4 0. 02 0. 50 9 — 0 — — - 0. 32 4 0. 01 0. 00 10 00 Sp in -L at tic e Re la xa tio n Ti m e (m s) Fi gu re 3. 20 : Ex am pl e o ft he in cr ea se in lo ng T1 v al ue s at lo w sa tu ra tio ns re la tiv e to T1 v al ue s a th ig he rs at ur at io ns . 1 10 10 0 0) Sample S31 Correlation of Saturation with M(O) 2200 2000 y—41012+21238x R”2_0921 1800 M(O) 1600 1400 1200 • 0.6 0.7 0.8 0.9 1.0 Saturation Sample S41 Correlation of Saturation with M(O) 1300- y = 14.402 + 1203.lx R’2 = 0.993 1200 1100 M(O) 600’’ 900 800 700 0.5 0.6 0.7 0.8 0.9 1.0 Saturation Figure 3.21: Examples of the correlation between fluid saturation and the total integral over the Ti distribution (M(0)). 1 07 DISCUSSION AND CONCLUSIONS The T1 distributions determined for porous samples are expected to reflect the pore sizes of the samples, with larger T1 values indicating the presence of larger pores. This expectation was met for both the discrete and continuous T1 spectra. Considering only those samples from the Ellcton Member, the porosity of these samples was primarily microscopic porosity, with S31, S32, and S33 possessing minor amounts of visible porosity. Correspondingly, the continuous and discrete T1 spectra of samples S43 to S49 reached their maximum amplitudes at relaxation times of 12 ms or less, while samples S31 to S33 reached their maximum amplitudes at relaxation times between 25 and 80 ms. The correlation between pore sizes and relaxation times was further supported by the continuous and discrete T1 spectra of samples S54, S55 and S77 possessing maximum amplitudes at 750 to 1750 ms. The discrepancy between the magnitude of the S31 to S49 relaxation spectra and those of samples S54, S55 and S77 was expected because the majority of the porosity of the final three samples was visible. The average pore sizes determined for samples S55 and S77 using the confocal imaging technique discussed earlier were 25 and 122 jim, which compare qualitatively with the average relaxation times the discrete and continuous T1 spectra occurring at 580 (discrete)! 470 ms (continuous) for sample S55 and 909/1392 ms for sample S77 The correlation between porosity and the integral over the T1 distributions (M(0)) was not very accurate, with regression parameters for the best fit straight lines of 0.728 and 0.883 for the continuous and discrete spectra respectively. There are three possible causes for this disparity in the theoretical and experimental correlations between porosity and M(0). Firstly, the samples are assumed to be fully saturated when in reality they are not, and furthermore, each sample will not be saturated to the same extent, therefore a certain amount of scatter is expected. The variation in saturations and the expectation that the samples are not fully saturated was formed by measuring the dry and saturated weights of several of the samples, calculating the volume of fluid in the 108 sample, and comparing this to the expected volume of the porosity. In most cases, the porosity in the samples actual saturation ranged between 80 and 90%. The second possible source of error in the correlation between porosity and M(0) lies in the accuracy of the porosity measurements: with the method used, error in the porosity may be as large as 0.5% porosity. The final contribution to error may be attributed to the difficulty in determining the actual value of M(O). Many of the decay curves are very steep at short measurement times, thus small changes in the fit of the inversion data to this portion of the curve can result in substantially different values of M(O). The pore sizes calculated using both the slow and fast diffusion methods agreed qualitatively with the presence and/or absence of visible porosity in the samples, with pore sizes from smallest to largest being calculated for the samples in this order: S34 to S49, S33, S31, S32, S54 and S55, and S77. Quantitatively, the comparison between both fast and slow diffusion pore sizes and the confocal pore size spectra shows reasonable agreement between the three. This comparison is shown in Figure 3.22 for samples S55 and S77. The inaccuracy of the fast diffusion pore size spectra lies primarily in the stretching of the distribution at the highest pore sizes due to the fact that the equation used to calculate pore sizes from T1 values becomes undefmed as T1 approaches Tib. This occurs due to the term: (1 1 ‘i-i ‘T1 Tib’ in the equation used to calculate pore sizes from T1 values with the fast diffusion method. The primary difficulty involved in interpreting the inversion results for the samples lies in the determination of the surface sink parameter, p. The values of p used in the fast diffusion pore size calculations were an order of magnitude less than those which produced the best fit using the slow diffusion method, and, as indicated by the inversion of theoretical data with noise added to it, determining the correct p value is very difficult unless an independently determined pore size distribution is available. The problem just mentioned may, if the data are accurate, be turned into an alvantage since a minimum in will occur slightly below the actual value for the sample and 1 09 thus the surface sink parameter can be determined without the availability of an independent pore size disthbution. Thus the slow diffusion method can be more useful than the fast diffusion method. One of the objectives of this study was to determine whether the fast diffusion method was adequate to calculate pore size distributions for carbonate samples containing vugs. The results of this investigation indicate that, although relatively acurate results may be determined with the fast diffusion method, the pore sizes in such samples are not all within the fast diffusion regime; therefore the slow diffusion method is more reliable for these samples. Indications of pore sizes outside the fast diffusion regime include the stretching of the fast diffusion pore size spectra at the largest pore sizes for samples S54, S55 and S77, and the fact that both the bulk relaxation rate and self diffusion coefficients have an obvious effect on the pore sizes calculated with the slow diffusion method. The primary conclusion drawn from this investigation is that the inversion of NIvIR data to calculate pore sizes is a very unstable problem. Three possible solutions are to select the type of pore size distribution to be determined, to incorporate more constraints into the inversion, or to obtain more accurate data over the largest possible range of measurement times. Despite the instability of the inversion, the slow diffusion results are more accurate than those of the fast diffusion methods for vuggy carbonate samples. The results indicate the range of pore sizes present in the sample and at least approximate the average pore size, as shown by comparing the confocal and slow diffusion pore size spectra. The second objective of this investigation was to determine the amount of information about fluid distributions that can be extracted from NMR spin-lattice relaxation data on partially saturated samples. The shift towards lower relaxation rates and lower amplitudes as the saturation decreases indicates two things: lower amplitudes reflect the presence of less fluid in the sample, while lower relaxation times indicate those portions of the fluid with the longest relaxation times are being removed preferentially. The fluid with the longest relaxation times fill the largest pores when the pores are in the fast diffusion regime and are found in the center of the largest pores for 110 Carbonate (S55) Pore Size Spectra 0.4 Slow Diffusion Spectra Confocal Spectra 0.3 Fast Diffusion Spectra 0 .— C.) I, 3- I 0.2 J 0 I I I 0.1 I ••]•• .f: I0.0 iO-7 i06 io5 io4 io3 Pore Radius (m) Carbonate (S77) Pore Size Spectra 0.4 Slow Diffusion Spectra Confocal Spectra ____ I IO Fast Diffusion Spectra II I ii I II I0.2 I E ii - II I ii I o ii I > I I ii I I I I 0.0 ..-l,. 106 io lo io Pore Radius (m) Figure 3.22: Comparisons of the pore size spectra calculated from the NMR data with the fast and slow diffusion methods with the spectra determined using the confocal laser scanning microscope. 111 samples with pores in the slow diffusion regime. This fact, coupled with the overall decrease in relaxation times with saturation, indicates that although a range of pore sizes are being drained concurrently, the fluid is removed from the largest pores more easily. The anomalous increase in Ti spectra at low saturations for samples S41, S42, S44, and S55, is probably due to the method of de-saturation used. The evaporation (drainage) of fluid from a sample tends to empty the largest pores first. However, even the largest pores in the center of the sample will not drain until the fluid within the pore is in contact with air, necessitating the emptying of many pores in the outer area of the sample before pores in the center will begin to empty. Also, if the accessibility to a large pore is limited by narrow pore throats, the large pore will empty only after the pore throats. Therefore the increase in the relaxation times at low saturations could indicate either the trapping of fluid in large pores due to their distance from the edges of the sample or limited accessiblity due to very narrow pore throats connecting some of the larger pores to the remainder of the pore network. The correlation of saturations with M(O) was quite accurate, and this result, along with the expected emptying of large pores first indicates that NMR could be used to determine more detailed information on fluid distributions if theoretical representations of partially saturated pores were developed. 112 CHAPTER FOUR Summary INTRODUCTION The primary objective of this thesis was to explore two new methods of characterizing porous media: confocal laser scanning microscopy (CLSM) and nuclear magnetic resonance (NMR). The benefits of these two methods include their non-destructive nature and the relative speed of the measurements. The primary advantage of the confocal microscope over other microscopes is its ability to produce serial sections by optically sectioning a sample, rather than physically sectioning it. Several advantages of NMR include the possible use of wet samples, the direct relation between pore sizes and relaxation times, and the fact that it may be possible to determine fluid distributions from decay curves of partially saturated media. CONFOCAL MICROSCOPY Three samples, a sandstone (Berea 100) and two limestones (S55 and S77), were characterized using the CLSM. The results for Berea 100 compared favourably with a Berea 300 pore size distributiont from the literature, with an average pore size of 25.3 p.m determined using the CLSM compared to 18.0 p.m determined for Berea 300 employing capillary pressure methods (Crocker et a!., 1983). The pore size spectra determined using these two methods are also very similar in shape, as shown in Figure 4.1. Therefore the pore size spectra determined for samples S55 and S77 are assumed to be accurate representations of the actual pore sizes. The pore size spectra for samples S55 and S77 are shown in Figure 4.2. The average pore sizes determined for these two samples was 17.5 and 90.2 p.m for S55 and S77 respectively. 113 40— 10 0 Mean 18.0Mm Mode 10.0Mm .2 0.4 0 C) I a) E — 0 > 1 2 4 10 20 40 100 Pore diameter (jzml 16 Images I I. 20 0.08 _ 0.06 0.04 0.02 0.00. I 106 io5 io Pore Radius (m) Figure 4.1: Berea 100 pore size distributions. The upper distribution is a capillaiy pressure pore size distribution (after Crocker et aL, 1983) and the lower is the distribution determined from 16 confocal laser scanning microscope images. 114 20 Images _____ L ______ Figure 4.2: Carbonate (S55 and S77) cumulative pore size distributions. The upper plot represents the summation of the distributions of 20 S55 images, the lower of 13 S77 images. 0.12 0.08 0.04 0.00 106 0.08 0.06 0.04 0.02 C .— C) j C —4 C io5 Pore Radius (m) 10 13 Images 0.00• 106 io-5 Pore Radius (m) io3 115 The ability of the confocal microscope to optically section samples was used to produce serial sections of the samples in the hope that reconstruction would allow three-dimensional image analysis. However, the samples were fragile and therefore enough serial sections to adequately represent the pores could not be collected. Partial sets of serial sections were collected, an example of which is shown in Figure 2.16, and the quality of these images indicate that with further development of the methods used to produce the pore casts, confocal laser scanning microscopy could be a powerful tool for three-dimensional characterization of porous media. The primary difficulty is one of computing power. The images are 384 K in size, and for detailed work, 100 or more serial sections could be collected at spacings as small as 0.2 ELm. The storage and analysis of such data would require a powerful computer with a large memory to allow reconstruction and analysis of the images. NMR This thesis focusses on several approaches to NMR which have not previously been used: direct comparison of fast diffusion and slow diffusion results for samples with pore sizes in both regimes, focussing on carbonates to determine the accuracy with which pores in carbonates may be determined, and measuring partially saturated samples to investigate the possibility of determining fluid disthbutions with NMR. Determination of Pore Sizes The fast diffusion method has, until recently, been the only method used to calculate pore sizes from NMR data. Despite evidence that pores in some samples used for other studies were not in the fast diffusion regime, only Davies et al. (1990) have applied slow diffusion theory to this problem. Rather than utilizing the more complex slow diffusion relations, most investigations 116 either assumed or theoretically proved that the effect of the larger pores on the calculated pore size distributions was minimal. The pore sizes determined using the fast and slow diffusion methods are very similar to those determined using the CLSM , as shown in the pore size spectra of S55 and S77 in Figure 3.22. The inaccuracy of the fast diffusion pore size spectra lies primarily in the stretching of the distribution at the highest pore sizes due to the fact that the equation used to calculate pore sizes from T1 values becomes undefmed as T1 approaches Tib. This occurs due to the term: (1 1 -i ‘T1 T1b’ in the equation used to calculate pore sizes from T1 values with the fast diffusion method. The similarity between the T1 and fast diffusion pore size spectra for samples S31 to S49 is illustrated by the comparison plots in Figure 4.3 for samples S3 1 and S39. This correspondence between the two spectra indicates that these samples are within the fast diffusion regime. Further support for the assumption that samples S31 to S49 are in the fast diffusion regime and samples S54, S55 and S77 are in the slow diffusion regime is provided by the nature of the variation in the slow diffusion pore sizes when the surface sink parameter, the self diffusion coefficient, and the bulk relaxation rate are varied. Theoretically, small variations in the bulk relaxation rate and self diffusion coefficient should not affect the calculation of pore sizes for samples in the fast diffusion regime while similar variations in the surface sink parameter should have no effect for samples in the slow diffusion regime. Thus, the indifference of the results of samples S3 1 through S49 to variations in the bulk relaxation rate and the randomness of their response to variations in the self diffusion coefficient indicate these samples are within the fast diffusion regime. In comparison, the pore sizes calculated for samples S54, S55 and S77 are affected by all three parameters, indicating the presence of pore sizes in both the fast and slow diffusion regime. 117 Ti and Pore Size Distributions 0.08 Pore Sizes Ti 0.06 0 C) 1:1.4 0.04 a) 0.020 0.00 ••••• 108 io 106 i05 i03 10 10W’ 100 Ti (ms) and Pore Radius (m) (a) Ti and Pore Size Distributions 0.08 Pore Sizes Ti 0.06 0 .— C) 1:1.4 0.04 ..% j 0 0.02 0.00• 108 io7 106 io- i010 10 10_i 10° Ti (ms) and Pore Radius (m) (b) Figure 4.3 T1 and fast diffusion pore size distributions. (a) Sample S39 (b) Sample S31 118 The fact that samples S54, S55 and S77 contain pore sizes in both the fast and slow diffusion regime should ensure that the slow diffusion method produces the most accurate pore size spectra for these samples. The graphs of the confocal, slow diffusion and fast diffusion pore size spectra of these two samples in Figure 3.22 indicate that this is the case only if the peak at the maximum pore size is neglected. This peak represents 60% of the total distribution and, with it included, average pore sizes of 35 and 558 p.m are calculated for S77 and S55 respectively. However, the inversion of theoretical data indicated that the addition of even 1% noise to a theoretical decay curve will produce a peak at the largest pore size which represents approximately 30% of the total spectra. Therefore, the final peak is probably due to noise and may be neglected. The average pore sizes calculated neglecting this peak are 19.3 and 85.3 urn, both comparing very favourably to the confocal averages of 17.5 and 90.2 p.m for samples S55 and S77 respectively. Fluid Distributions The application of NIvIR to determine fluid distributions indicated that the fluid which evaporates from the sample most easily is that portion of the fluid with the longest relaxation time. This is shown by the decrease in the T1 values in the spectra calculated at lower saturations compared to those calculated at higher saturations. The portion of the fluid with the longest relaxation time is found in the largest pores when the pores are in the fast diffusion regime and in the center of the largest pores for samples with pores in the slow diffusion regime. This fact, coupled with the overall decrease in relaxation times with saturation, indicates that although a range of pore sizes are being drained concurrently, the fluid is removed from the largest pores more easily. The anomalous increase in Ti spectra at low saturations for samples S41, S42, S44, and S55, is probably due to the method of de-saturation used. The evaporation (drainage) of fluid from a sample tends to empty the largest pores first. However, even the largest pores in the center of the sample will not drain until the fluid within the pore is in contact with air, necessitating the 119 emptying of many pores in the outer area of the sample before pores in the center will begin to empty. Also, if the accessibility to a large pore is limited by narrow pore throats, the large pore wifi empty only after the pore throats do. Therefore the increase in the relaxation times at low saturations could indicate either the trapping of fluid in large pores due to their distance from the edges of the sample or limited accessiblity due to very narrow pore throats connecting some of the larger pores to the remainder of the pore network. CONCLUSIONS The ability to accurately characterize the pore space of geophysical materials is central to our understanding of the physical properties of these materials. In this study, it has been shown that scanning confocal laser microscopy and nuclear magnetic resonance can provide much information about the sizes and shapes of pores, and the nature of the fluids present. In order to further expand the application of these technologies to determine pore network characteristics, the actual geometrical shape of the pores must be approximated more closely, a more stable method of inverting NMR data in the slow diffusion regime must be developed. The collection of more accurate data over the largest possible measurement time range should also improve the results dramatically. However, if this method is ever to be applied to data collected in situ, it must be developed to deal with relatively noisy data which does not extend to measurement times below 20 ms (Neuman and Brown, 1982). The determination of fluid distributions could also be improved by determining the response of partially saturated pores and then using this information in the analysis of the data. 1 20 REFERENCES Bacallao, R., and Steizer, E.H.K. 1989: Preservation of Biological Specimens for Observation in a Confocal Fluorescence Microscope and Operational Principles of Confocal Fluorescence Microscopy; Methods in Cell Biology, v. 31, p. 437-453. Bloch, F 1946: Nuclear Induction; Physical Review, v. 70, p. 460-474. Borgia, G.C., Fantazzini, P., and E. Mesini 1990: Water ‘H Spin-Lattice Relaxation as a Fingerprint of PorousMedia; Magnetic Resonance Imaging, v. 8, p. 435-447. Bourbie, T., and B. Zinszner 1984: Saturation Methods and Attenuation Versus Saturation Relationships in Fontainebleau Sandstone; Fifty-Fourth International Meeting, Society of Exploration Geophysicists, Expanded Abstracts, p. 344-350. Brown, R.J.S., and I. Fatt 1956: Measurements of Fractional Wettability of Oilfield Rocks by the Nuclear Magnetic Relaxation Method; Petroleum Transactions, American Institue of Mining aned Metallurgical Engineers, v. 207, p. 262-264. Brownstein, K.R., and C.E. Tarr 1977: Spin-Lattice Relaxation in a System Governed by Diffusion; Journal of Magnetic Resonance, v. 26, p. 17-24. Brownstein, K.R., and C.E. Tarr 1979: Importance of Classical Diffusion in NMR Studies of Water in Biological Cells; Physical Review A, v. 19, p. 2446-2453. Brunauer, S., Deming, L.S., Deming, W.E., and Teller, E. 1938: Adsorption of Gases in Multimolecular Layers; American Chemical Society Journal, v. 60, p. 309-319. Callaghan, P.T., Jolley, K.W., and R.S. Humphrey 1983: Diffusion of Fat and Water in Cheese as Studied by Pulsed Field Gradient NMR; Journal of Colloid and Interface Science, v. 93, p. 52 1-529. Cercone, K.R., and Pedone, V.A. 1987: Fluorescence (Photoluminescence) of Carbonate Rocks: Instrumental and Analytical Sources of Observational Error; Journal of Sedimentary Petrology, v. 57(4), p. 780-782. Coates, G.R., Muller, M., and G. Henderson 1991: An Investigation of a New Magnetic Resonance Imaging Log; SPWLA Thirty- Second Annual Logging Symposium, Paper D. Crocker, M.E., Donaldson, E.C., and Marchin, L.M. 1983: Comparison and Analysis of Reservoir Rocks and Related Clays; Society of Petroleum Engineers Paper 11973. 121 Davies, S., and K.J. Packer 1990: Pore-size distributions from nuclear magnetic resonance spin-lattice relaxation measurements of fluid-saturated porous solids. I. Theory and simulation; Journal of Applied Physics, v. 67(6), p. 3163-3170. Davies, S., Kalam, M.Z., Packer, K.J., and F.O. Zelaya 1990: Pore-size distributions from nuclear magnetic resonance spin-lattice relaxation measurements of fluid-saturated porous solids. II. Applications to reservoir core samples; Journal of Applied Physics, v. 67(6), p. 3171-3176. Domenico, S.N 1976: Effect of Brine-Gas Mixture on Velocity in an Unconsolidated Sand Reservoir; Geophysics, v. 41, p. 882-894. Domenico, S.N. 1977: Elastic Properties of Unconsolidated Porous Sand Reservoirs; Geophysics, v. 42, p. 1339-1368. Doyen, P.M. 1988: Permeability, Conductivity and Pore Geometry of Sandstone; Journal of Geophysical Research: v. 93, p. 7729-7740. Egger, M.D., and Petran, M. 1967: New Reflected-Light Microscope for Viewing Unstained Brain and Ganglion Cells; Science, v. 157, p. 305-307. Ehrlich, R., Kennedy, S.K., Crabtree, S.J., and Cannon, R.L. 1984: Petrographic Image Analysis, I. Analysis of Reservoir Pore Complexes; Journal of Sedimentary Petrology, v. 54(4), p. 1365-1378. Evamy, B.D. 1969: The Precipitational Environment and Correlation of Some Calcite Cements Deduced From Artificial Staining; Journal of Sedimentary Petrology, v. 39(2), p. 787-821. Fukushima, E. and S.B.W. Roeder 1981: Experimental Pulse NMR: A Nuts and Bolts Approach; Addison-Wesley Publishing Company, Reading, Massachusetts, 539 p. Gallegos, D.P., Munn, K., Smith, D.M., and D.L. Stermer 1987: A NMR Technique for the Analysis of Pore Structure: Application to Materials with Well-Defined Pore Structure; Journal of Colloid and Interface Science, v. 119(1), p. 127-140. Gallegos, D.P., and D.M. Smith 1988: A NMR Technique for the Analysis of Pore Structure: Determination of Continuous Pore Size Distributions; Journal of Colloid and Interface Science, v. 122(1), p. 143- 153. Gallegos, D.P., Smith, D.M., and C.J. Brinker 1988: A NMR Technique for the Analysis of Pore Structure: Application to Mesopores and Micropores; Journal of Colloid and Interface Science, v. 124(1), p. 186-198. 1 22 Gies, R.M. 1987: An Improved Method for Viewing Micropore Systems in Rocks with the Polarizing Microscope; Society of Petroleum Engineers Formation Evaluation, v. 2(2), p. 209-214. Grayson, J.F. 1956: The Conversion of Calcite to Fluorite; Micropaleontology, v. 2, p. 7 1-78. Hamilton, D.K., and Wilson, T. 1982: Surface Profile Measurements using the Confocal Microscope; Journal of Applied Physics, v. 7, p. 5320-5322. Hedberg, S.A. 1990: The NMR Response of Sandstones and its application to the detection of Hydrocarbon Contamination; B.A.Sc. Thesis, University of British Columbia, Vancouver, British Columbia, Canada. Jackson, J.A. 1984: Nuclear Magnetic Resonance Well Logging; The Log Analyst, September-October, p. 16-30. Kenyon, W.E., Day, P.1., Straley, C., and J.F. Willemsen 1988: A Three-Part Study of NMR Longitudinal Relaxation Properties of Water-Saturated Sandstones; SPE Formation Evaluation, p. 622-633. Kenyon, W.E., Howard, J.J., Sezginer, A., Straley, C., Matteson, A., Horkowitz, K., and R.Ehrlich 1989: Pore-size Distribution and NMR in Microporous Cherty Sandstones; Society of Professional Well Log Analysts Thirtieth Annual Logging Symposium, Paper LL. Kleinberg, R.L., and Horsfield, M.A. 1990: Transverse Relaxation Processes in Porous Sedimentary Rock; Journal of Magnetic Resonance, v. 88, p. 9-19. Knight, R. 1991: Hysteresis in the Electrical Resistivity of Partially Saturated Sandstones; Geophysics, v. 56(12), p. 2139-2147. Knight, R., and R. Nolen-Hoeksema 1990: A Laboratory Study of the Dependence of Elastic Wave Velocities on Pore Scale Fluid Distribution; Geophysical Research Letters, v. 17(10), p. 1529-1532. Knight, R., and A. Nur 1987: Geometrical Effects in the Dielectric Response of Partially Saturated Sandstones; The Log Analyst, v. 28, p. 513-5 19. Knispel, R.R. 1981: Review of Two Models for Water Proton Relaxation in Tissue; Bulletin of Magnetic Resonance, v. 3, p. 104-108. Lin, C., and Hamasaki, J. 1983: Pore Geometry: A New System for Quantitative Analysis and 3-D display; Journal of Sedimentary Petrology, v.53, p. 670-672. 1 23 Longeron, D.G., Arguad, M.J., and J.P. Feraud 1989: Effect of Overburden Pressure and the Nature and Microscopic Distribution of Fluids on Electrical Properties of Rock Samples; Society of Petroleum Engineers Formation Evaluation, v. 4, p. 194-202. Loren, J.D., and J.D. Robinson 1970: Relations Between Pore Size Fluid and Matrix Properties, and NIvIL Measurements; Society of Petroleum Engineers Journal, v. 249 , p. 268-278. Lukosz, W. 1966: Optical Systems with Resolving Powers Exceeding the Classical Limit; Journal of the Optical Society of America, v. 56, p. 1463-1472. Marshall, G.E., Konstas, A.G., and Lee, W.R. 1990: Ultrastructural distribution of collagen yypes 1-VI in aging human retinal vessels; British Journal of Opthamology, v. 74(4), p. 228-232. Masters, B.R., Kriete, A., and Kukulies, J. 1991: Fluorescence Confocal Microscopy of the Living Eye; Scanning, v. 13, p. 1-16. Menon, R.E., MacKay, A.L., Hailey, J.R.T., Bloom, M., Burgess, A.E., and J.S. Swanson 1987: An NMR determination of the physiological water distribution in wood during drying; Journal of Applied Polymer Science, v. 33, p. 1141-1155. Minsky, M. 1961: Microscopy Apparatus; US Patent 3013467, Dec. 19 (Filed Nov. 7, 1957). Munn, K., and D.M. Smith 1987: A NMR Technique for the Analysis of Pore Structure: Numerical Inversion of Relaxation Measurements; Journal of Colloid and Interface Science, v. 119(1), p. 117- 126. Neuman, C.H., and R.J.S. Brown 1982: Applications of Nuclear Magnetism Logging to Formation Evaluation; Journal of Petroleum Technology, v. 34(12), p. 2853-2862. Nutting, P.G. 1929: Some physical problems in oil recovery; Oil and Gas Journal, v. 28, p. 44-45. Parker, K.J. and C. Rees 1972: Pulsed NMR Studies of Restricted Diffusion I. Droplet Size Distributions in Emulsions; Journal of Colloid and Interface Science, v. 40(2), p. 206-2 18. Penner, D.G. 1957: Facies and porosity relationships in the Mississippian Elkton carbonate cycle of southwestern Alberta; Alberta Society of Petroleum Geologists Journal, v. 5(5), p. 101-104. Pierce, C. 1953: Computation of Pore Sizes from Physical Adsorption Data; Journal of Physical Chemistry, v. 57, p. 149-152. 1 24 Pittman, E.D., and Duschatko, R.W. 1970: Use of pore casts and scanning electron microscopy to study pore geometry; Journal of Sedimentary Petrology, v. 40, p. 1153-1157. Rassineux, F., Beaufort, D., Meunier, A., and Bouchet, A. 1987: A Method of Coloration by Fluorescein Aqueous Solution for Thin-Section Microscopic Observation; Journal of Sedimentary Petrology, v. 57, p. 782-783. Robinson, J.D., Loren, J.D., Vajnar, E.A., and D.E. Hartman 1974: Determining Residual Oil with the Nuclear Magnetism Log; Journal of Petroleum Technology, v. 257, p. 226-236. Russ, J.C. 1990: Computer Assisted Microscopy: The Measurement and Analysis of Images, Plenum Press, 453 p. Russ, J.L., Palmour III, H., and Hare, T.M. 1989: Direct 3-D pore location measurement in alumina Journal of Microscopy, v. 155, p. RP1-RP2. Ruzyla, K., and Jezek, D.I. 1987: Staining Method for Recognition of Pore Space in Thin and Polished Sections; Journal of Sedimentary Petrology, v. 57, p. 777-778. Schmidt, E.J., Velasco, K.K., and A.M. Nur 1986: Quantifying solid-fluid interfacial phenomena in porous rocks with proton nuclear magnetic resonance; Journal of Applied Physics, v. 59(8), p. 2788-2797. Shuman, H., Murray, J.M., and DiLullo, C. 1989: Confocal Microscopy: An Overview; Biotechniques, v. 7, p. 154-163. Smart, P., and Tovey, N.K. 1982: Electron microscopy of soils and sediments: techniques; Clarendon Press, Oxford, p. 61-80. Soeder, D.J. 1990: Applications of Fluorescence Microscopy to Study of Pores in Tight Rocks; American Association of Petroleum Geologists Bulletin, v. 74, p. 30-40. Spurr, A.R. 1969: A Low-Viscosity Epoxy Resin Embedding Medium for Electron Microscopy; Journal of Ultrastructure Research, v. 26, p. 3 1-43. Straley, C., and Minnis, M.M. 1983: Epoxy Rock Replicas for Microtoming; Journal of Sedimentary Petrology, v. 53, p. 667-669. Straley, C., Morriss, C.E., Kenyon, W.E., and J.J. Howard 1991: NMR in Partially Saturated Rocks: Laboratory Insights on Free Fluid Index and Comparison with Borehole Logs; Society of Petroleum Well Log Analysts Thirty Second Annual Logging Symposium, Paper CC. Tercier, P. 1992: Electrical Resistivity of Partially Saturated Sandstones; M.Sc. Thesis, University of British Columbia, Vancouver, British Columbia, Canada. 1 25 Thomas, G.E. and R.P. Glaister 1960: Facies and Porosity Relationships in some Mississippian Carbonate Cycles of the Western Canada Basin; Bulletin of the American Association of Petroleum Geologists, v. 44(5), p. 569-588. Timur, A. 1969: Pulsed Nuclear Magnetic Resonance Studies of Porosity, Movable Fluid, and Permeability of Sandstones; Journal of Petroleum Technology, v. 246, P. 775-786. Vinegar, H.J., Tutunjian, P.N., Crabtree, P.T., Taffaldi, R.J., DiFoggio, R., and W.A. Edeistein 1989: NMR Spectroscopy of Tight, Gypsum-Bearing Carbonates; Society of Professional Well Log Analysts Thirtieth Annual Logging Symposium, Paper HH. Wardlaw, N.C. 1976: Pore Geometry of Carbonate Rocks as Revealed by Pore Casts and Capillary Pressure; American Association of Petroleum Geology Bulletin, v. 60(2), p. 245- 257. Whittall, K.P. and A.L. MacKay 1989: Quantitative Interpretation of NMR Relaxation Data; Journal of Magnetic Resonance, v. 84, p. 134-152. Whittall, K.P. 1991: Recovering Compartment Sizes from NMR Relaxation Data; Journal of Magnetic Resonance, v. 94, p. 486-492. Wicksell, S.D. 1926: The corpuscle problem. Second memoir. Case of ellipsoidal corpuscles; Biometrika, v. 18, p. 152-172. Wilson, T. (ed.) 1990: Confocal Microscopy, Academic Press, Great Britain, 426 p. Winsauer, W.O., Shearin, H.M. Jr., Masson, P.H., and Williams, M. 1952: Resistivity of brine-saturated sands in relation to pore geometry; American Association of Petroleum Geology Bulletin, v. 36, p. 253-277. Young, F.G. 1967: Elkton Reservoir of Edson Gas Field, Alberta; Bulletin of Canadian Petroleum Geology, v. 15(1), p. 50-64. Young J.Z., and Roberts, F. 1951: A Flying-Spot Microscope; Nature, v. 167, p. 231. Zimmerman, J.R. and W.E. Brittin 1957: Nuclear Magnetic Resonance Studies in Multiple Phase Systems: Lifetime of a Water Molecule in an Adsorbing Phase on Silica Gel;Journal of Physical Chemistry, v. 61(6), p. 1328-1333. 1 26 APPENDIX 1 BioRad MRC-500 Image File Format These picture files are standard MS-DOS files which consist of a 76 byte header, a pixel array of dimensions NX x NY x Npic, and a series of note records of 96 bytes each. Header Structure Contents 0-1 NX - X pixel dimension of the image 2-3 NY - Y pixel dimension of the image 4-5 Npic - number of image frames in the file 6-7 Contrast ramp 0 (pixel value for full black) 8-9 Contrast ramp 255 (pixel value for peak white) 10-13 Ignore 14-15 1 if byte pixels, 0 if word pixels 16-49 Ignore 50-51 1 if merged image, 0 if unmerged image 52-53 Low byte used for selecting RGB colours (7=white) 54-55 Set to 12345, used to check for valid picture file 56-57 Contrast ramp 0 for second picture of merged pair 58-59 Contrast ramp 255 for second picture of merged pair 60-6 1 Low byte used for RGB colour of second picture of merged pair 64-65 lens power 66-69 Magnification factor, real variable - see note below 70-75 Ignore NOTE: the length of a line of pixels in microns is computed as: ((# pixels)x(Scale factor))/((lens mag.)x(mag. factor)) where the scale factor is the pixel to micron ratio which is set up in the setup menu, and the magnification factor is the value stored in the file header. Pixel Array: The pixel array contains Npic images with X and Y pixel dimensions of NX and NY. This array is either NX x NY x Npic bytes in size or, if word (16 bit) pixels rather than the usual byte (8 bit) pixels are used, the array is 2 x (NX x NY x Npic) bytes in size. 1 27 - L I\) I . — J c I 9 ? - ‘ .0 ‘ Ii I - C — t o ; 10 I 0 II CD - . I— - CD CD 0 — ‘ . + f CD - 4 • E ‘ . c J - — . ic. ‘-‘D• 00 -f- CD E Ii r.crQ O 0 ri CD APPENDIX 2 Image Analysis Programs BINREP.F C C FEBUARY, 1992 ALICE CHAPMAN C ROCK PHYSICS LAB, UBC C C THIS PROGRAM WILL READ IN CONFOCAL IMAGE FILES WRITTEN C IN THE FORMAT 10013 AND AN ARRAY SIZE OF 768 BY 512 C PIXELS. THE INTENSITY DISTRIBUTIONS OF THE IMAGES ARE DE C TERMINED AND WR11TEN TO AN OUTPUT FILE, A PORE/SOLID C DIVISION IS DONE TO PRODUCE BINARY IMAGES, AND THE SMALL C HOLES IN THE PORES ARE CLOSED USING A DILATION/EROSION C SEQUENCE. FINALLY, SMALL PORES WHICH ARE ASSUMED TO BE C ARTIFACTS OF THE CHOSEN PORE/SOLD CUTOFF VALUE ARE C REMOVED FROM THE IMAGES USING AN EROSION/DILATION SEQUENCE. C C MODIFIED FROM: BIN.F C C*************************** C C THE VARIABLES ARE DECLARED. C INTEGER*2 IMG(530,780) INTEGER POR,SOL,INTNS(256),DIL,ERO,NREP CHARACTER*20 INP(20),OUTP(20) C C*************************** C C VARIABLE DESCRIPTION: C C IMG(530,780): THE CONFOCAL IMAGE IS PLACED IN THIS ARRAY C WITH THREE ROWS OF SOLIDS SURROUNDING IT C POR: THE NUMBER WHICH PIXELS WITHIN THE PORES OF THE IMAGE C WILL BE SET TO (255) C SOL: THE NUMBER WHICH PIXELS WiTHIN THE SOLID PORTION OF THE C 11VIAGE WILL BE SET TO (0) C INTNS(256): A RECORD OF THE NUMBER OF PIXELS IN THE ORIGINAL C IMAGE WHICH HAD A FLUORESCENT INTENSiTY OF 0 TO 255 C DIL: THE NUMBER WHICH SOLID PIXELS IN CONTACT WITH THE PORES C WILL BE SET TO DURING THE DILATION SEQUENCE C ERO: THE NUMBER WHICH PORE PIXELS IN CONTACT WITH THE SOLD) C WILL BE SET TO DURING THE EROSION SEQUENCE C NREP: THE NUMBER OF IMAGES WHICH WILL BE ALTERED WITH THIS C PROGRAM C INP: THE NAMES OF THE INPUT FILES (USUALLY THE NAMES OF THE C MODELS) C OUTP: THE NAMES OF THE OUTPUT FILES C C*************************** C 1 29 C ThE WIDTH (NX) AND LENGTH (NY) OF THE CONFOCAL IMAGE ARRAY C AND THE WORKING ARRAY (NTX) AND (NTY) ARE SET. C NX = 512 NY = 768 NTX = NX + 6 NTY = NY +6 C C THE INTENSITY OF FLUORESCENT LIGHT RESPONSE WHICH CORRESPONDS C TO THE CUT OFF BETWEEN PORE AND SOLID IN THE ORIGINAL IMAGE IS C SET BY THE USER. C WpJTE(9,*) “WHAT CUT OFF WOULD YOU LIKE FOR THE SOLUYr’ READ(9,*) NSOL C C THE NUMBER OF IMAGES WHICH WILL BE ALTERED IS INiTIALIZED TO ZERO, C THEN THE USER IS REQUESTED TO ENTER THE DESIRED NUMBER. C 11 NREP=0 WPJTE(9,*)”PLEASE ENTER THE NUMBER OF IMAGES & YOU WISH TO HAVE ALTERED (0 <INT < 100)” READ(9,*)NREP C C THE NUMBER OF IMAGES TO BE ALTERED IS CHECKED TO ENSURE THAT C iT IS LESS THAN 100 AND IS NOT NEGATIVE. IF IT DOES NOT MEET C THESE REQUIREMENTS, A MESSAGE IS WRITTEN FOR THE USER, AND C THEY ARE ASKED TO ENTER AN ACCEPTABLE VALUE. C IF ((NREP.LT.1).OR.(NREP.GT.99)) THEN IF (NREP.LT.1) THEN wPJTE(9,*)”THIs NUMBER IS NEGATIVE!” ELSE WRITE(9,*)”THIS NUMBER IS TOO LARGE!” ENDIF GOTO 11 ENDIF C C THE NAME OF THE INPUT AND OUTPUT FILES ARE REQUESTED FROM THE C USER. C DO 5 I=1,NREP wR1TE(9,*)”EITER THE IMAGE FILE NAME” READ(9,*)INP(I) WRITE(9,*) “ENTER THE OUTPUT FILE NAME” READ(9,*)OUTP(I) 5 CONTINUE C C THE LOOP WHICH WILL ALTER AND WRiTE OUT EACH IMAGE IS C BEGUN. C DO 200 NINC = 1,NREP C C THE NINC (FIRST,SECOND...) INPUT AND OUTPUT FILES ARE OPENED. 1 30 C OPEN(5 ,FILE=INP(NINC)) OPEN(6,FILE=OUTP(NINC)) C C THE TWO ARRAYS IN THE PROGRAM ARE IMTIALIZED TO ZERO. C DO 10I=1,NTX DO 10 J=1,NTY IMG(I,J) = 0 10 CONTINUE DO 20 1=1,256 1NTNS(I) =0 20 CONTINUE C C THE IMAGE IS READ IN FROM THE IMAGE FILE. C READ(5,91 )NTEST,NTEST2 91 FORM.AT(/////////18/18) DO 30 I=4,NX÷3 WRITE(9,*)”I = 11,1 READ(5,92)(IMG(I,J),J=4,NY÷3) 30 CONTINUE 92 FORMAT(10013) C C THE INTENSITY OF EACH PIXEL IN THE ORIGINAL IMAGE IS C RECORDED IN THE ARRAY INTNS AND WRHTEN TO THE FILE C “INTNS”INP WHICH WAS OPENED FOR THAT PURPOSE. C OUTP(NINC) = “INTNS”/IINP(NINC) OPEN(7,FILE = OUTP(NTNC)) C DO 40 I=4,NX+3 DO 40 J=4,NY+3 ITMP = IMG(I,J) + 1 INTNS(ITMP) = INTNS(1TMP) + 1 40 CONTINUE C WRITE(7,93) 93 FORMAT(”INTENSITY DISTRIBUTION”) C DO 501=1,256 WRITE(7,94)I,INTNS(I) 50 CONTINUE 94 FORMAT(2110) C CLOSE(7) C C THE PORE/SOLID DIVISION IS PERFORMED TO PRODUCE A BINARY C IMAGE OF ZERO’S (SOLID) AND 255’S (PORE). C SOL=0 POR = 255 C 131 DO 60 1=4,NX+3 DO 60 J=4,NY+3 IF (IMG(I,J).GT.NSOL) THEN IMG(I,J) = POR ELSE IMG(I,J) = SOL ENDIF 60 CONTINUE C C NOW AN EROSION/DILATION SEQUENCE IS PERFORMED TO REMOVE ANY C SMALL PORES WHICH ARE SIMPLY ARTIFACTS OF THE PROCESS OF C CHANGING THE IMAGE TO A BINARY IMAGE. C ERO = 245 C C ANY PORE PIXELS IN CONTACT WITH A SOLID ARE RE-NUMBERED TO C THE ERO VALUE OF 245. C DO 70 I=4,NX+3 DO 70 J=4,NY÷3 IF (IMG(I,J).EQ.SOL) THEN IF (IMG(I÷1,J).EQ.POR) IMG(I-i-1,J) = ERO IF (IMG(I-1,J).EQ.POR) IMG(I-1,J) = ERO IF (IMG(I,J+1).EQ.POR) IMG(I,J-i-1) = ERO IF (IMG(I,J-1).EQ.POR) IMG(I,J-1) = ERO 14 ENDIF 70 CONTINUE C C PIXELS WHICH WERE CHANGED TO ERO AN]) ARE STILL IN CONTACT C WITH THE PORE ARE RE-SET TO POR(255). EACH PIXEL IN THE C IMAGE IS CHECKED AND THE IMAGE IS GONE THROUGH ONLY ONCE TO C ENSURE THAT SMALL ROUGHNESSES ON THE PORES EXTENDING INTO C THE SOLID ARE SMOOTHED. C DO 80 I=4,NX+3 DO 80 J=4,NY+3 IF (IMG(I,J).EQ.POR) THEN IF (IMG(I+1,J).EQ.ERO) IMG(I+1,J) = POR IF (IMG(I-1,J).EQ.ERO) IMG(I-1,J) = POR IF (IMG(I,J+1).EQ.ERO) IMG(I,J+1) = POR IF (IMG(I,J-1).EQ.ERO) IMG(I,J-1) = POR 15 ENDIF 80 CONTINUE C C ANY ALTERED PIXELS WHICH ARE EITHER COMPLETELY SURROUNDED BY C SOLID PIXELS OR ARE PART OF EXTENWE ROUGHNESS ON THE PORE C SURFACE ARE STILL SET TO ERO AN]) ARE NOW IRREVOCABLY C RE-SET TO SOLID PIXELS (0). C 16 DO 90 I=4,NX÷3 DO 90 J=4,NY+3 IF (IMG(I,J).EQ.ERO) THEN IMG(I,J) = SOL ENDIF 1 32 90 CONTINUE C C NOW A DILATION/EROSION SEQUENCE IS PERFORMED TO CLOSE ANY C SMALL HOLES IN THE PORES. C DIL = 250 C C ANY SOLID PIXELS IN CONTACT WITH A PORE ARE RE-NUMBERED TO C THE DIL VALUE OF 250. C DO 100 I=4,NX+3 DO 100 J=4,NY÷3 IF (IMG(I,J).EQ.POR) THEN IF (IMG(I+1,J).EQ.SOL) IMG(I+1,J) = DIL IF (IMG(I-1,J).EQ.SOL) IMG(I-1,J) = DIL IF (IMG(I,J+1).EQ.SOL) IMG(I,J+1) = DIL IF (IMG(I,J-1).EQ.SOL) IMG(I,J-1) = DIL ENDIF 100 CONTINUE C C PIXELS WHICH WERE CHANGED TO DIL AND ARE STILL IN CONTACT C WITH THE SOLID ARE RE-SET TO SOL(0). EACH PIXEL IN THE C IMAGE IS CHECKED AND THE IMAGE IS GONE THROUGH ONLY ONCE, C THUS FILLING SMALL CHANNELS INTO THE PORE SPACE AN]) GENERALLY C ‘SMOOTHING’ THE BINARY IMAGE. C DO 110 I=2,NX+5 DO 110 J=2,NY+5 IF (IMG(I,J).EQ.SOL) THEN IF (IMG(I+1,J).EQ.DIL) IMG(I+1,J) = SOL IF (IMG(I-1,J).EQ.DIL) IMG(I-1,J) = SOL IF (IMG(I,J+1).EQ.DIL) IMG(I,J+1) = SOL IF (IMG(I,J-1).EQ.DIL) IMG(I,J-1) = SOL 12 ENDIF 110 CONTINUE C C ANY ALTERED PIXELS WHICH ARE COMPLETELY SURROUNDED BY PORE C PIXELS ARE STILL SET TO DIL AND ARE NOW IRREVOCABLY RE-SET C TO PORE PIXELS (255). C 13 DO 120 I=4,NX+3 DO 120 J=4,NY+3 IF (IMG(I,J).EQ.D1L) THEN IMG(I,J) = POR ENDIF 120 CONTINUE C C THE MODIFIED IMAGE IS NOW WREN INTO THE OUTPUT FILE. C DO 130 I=4,NX--3 WRITE(6,95)(IMG(I,J),J=4,NY+3) 130 CONTINUE 95 FORMAT(10014) C 1 33 C THE INPUT AND OUTPUT FILES ARE CLOSED. C CLOSE(5) CLOSE(6) C C THE NEXT STAThMENT IS THE END OF THE LOOP WHICH C IS REPEAThD FOR EACH IMAGE. C 200 CONTINUE C STOP END 1 34 IMAGANAL.F C C MARCH, 1992 ALICE CHAPMAN C ROCK PHYSICS LAB, UBC C C THIS PROGRAM READS iN BINARY CONFOCAL IMAGE FILES C PRODUCED BY THE PROGRAM BINIMGR.F. THESE IMAGES ARE C MANIPULATED TO ALTER ALL PIXELS WiTHIN PORE ThROATS C TO VALUES OF SEPR (20). A WATERSHED ALGORiTHM IS C THEN USED TO DETERMINE THE POINTS WITHIN THE PORES C AND PORE THROATS WHICH ARE FARTHEST FROM THE EDGE OF C THE PORE OR PORE THROATS. THESE PORE AND PORE THROAT C CENTERS ARE THEN USED TO DETERMINE THE SIZE OF THE C PORES AND PORE THROATS. FINALLY, THE OVERALL SURFACE C AREA AND VOLUME OF THE PORES iN THE IMAGE IS DETERMINED C AN]) THE POROSITY (RATIO OF THE VOLUME OF PORES TO C TOTAL VOLUME OF THE IMAGE) AND SURFACE AREA TO VOLUME C RATIO OF THE PORES IS CALCULATED. C C C*************************** C C THE VARIABLES ARE DECLARED. C INTEGER*2 IMG(530,780) INTEGER POR,SOL,SEPR,BDRY,PRDIST( 101) INTEGER NPR,NTHR,IMAX,THDIST( 101),MXRAD REAL POROS,PSV,TSV,THSV CHARACTER*20 INP(20),OUTP(20),MDL(20),TMP C C*************************** C C VARIABLE DESCRIPTION: C C IMG(530,780): THE CONFOCAL IMAGE IS PLACED IN THIS ARRAY C WITH THREE ROWS OF SOLIDS SURROUNDING IT C POR: THE NUMBER WHICH PIXELS WITHIN THE PORES OF THE C IMAGE WILL BE SET TO (255) C SOL: THE NUMBER WHICH PIXELS WITHIN THE SOLD) PORTION OF C THEIMAGEWILLBESETTO(0) C SEPR: THE NUMBER WHICH PIXELS IN THE SMALLEST PART OF THE C PORE THROATS WILL BE SET TO TO SEPARATE THE PORES C BDRY: THE INITIAL NUMBER WHICH PORE PIXELS BORDERING THE SOLD) C ARE SET TO IN THE WATERSHED ALGORiTHM. PORE PIXELS WHICH C BORDER THOSE WHICH ARE SET TO BDRY ARE THEN SET TO BDRY+1 C ETC, ETC. C CrR: THE NUMBER TO WHICH THE PIXEL AT THE CENTRE OF A PORE C OR PORE THROAT IS SET C NREP: THE NUMBER OF IMAGES WHICH WILL BE ALTERED WITH THIS C PROGRAM C PRDIST: A DISTRIBUTION OF THE PORE SIZES IN THE IMAGE. PRDIST(I) C WILL BE A NUMBER (J) INDICATING J PORES OF RADI1JS 1-0.5 C WERE PRESENT IN THE IMAGE C THDIST: A DISTRIBUTION OF THE PORE THROAT SIZES IN THE IMAGE. 1 35 C THDIST(I) WILL BE A NUMBER (J) INDICATING J PORE THROATS C OF RADIUS 1-0.5 WERE PRESENT IN THE IMAGE C NPR: THE TOTAL NUMBER OF PORES WHOSE SIZE WAS MEASURED C NTHR:THE TOTAL NUMBER OF PORE THROATS WHOSE SIZE WAS MEASURED C MXRAD:THE INTEGER REPRESENTING THE SIZE OF THE LARGEST PORE SIZE C MEASURED. MXRAD = RADIIJS+0.5 C POROS: THE POROSITY OF THE IMAGE: IE THE RATIO OF THE NUMBER C OF PORE AND PORE THROAT PIXELS IN THE IMAGE TO THE TOTAL C NUMBER OF PIXELS IN THE IMAGE. C PSV: THE SURFACE TO VOLUME RATIO (PSA/PVOL) OF THE PORES C THSV: THE SURFACE TO VOLUME RATIO (THSA1HVOL) OF THE PORE C THROATS C TSV:THE OVERALL SURFACE TO VOLUME RATIO OF THE VOII) SPACE IN THE C IMAGE: (PSA + THSA)/(PVOL + THVOL) C INP: THE NAMES OF THE INPUT FILES C OUTP: THE NAMES OF THE OUTPUT FILES C MDL: THE NAME OF THE MODELS, WHICH ARE USED TO NAME SOME C OUTPUT FILES C C*************************** C C THE WIDTH (NX) AND LENGTH (NY) OF THE CONFOCAL IMAGE ARRAY C AND THE WORKING ARRAY (NTX) AND (NTY) ARE SET. C NX = 512 NY = 768 NTX = NX + 6 NTY = NY +6 C C THE NUMBER OF IMAGES WHICH WILL BE ALTERED IS INITIALIZED TO C ZERO, THEN THE USER IS REQUESTED TO ENTER THE DESIRED NUMBER. C 15 NREP=0 WRITE(9,*)”PLEASE ENTER THE NUMBER OF IMAGES & YOU WISH TO HAVE ALTERED (0 <INT < 100)” READ(9,*)NREP C C THE VALUE OF NREP IS CHECKED TO ENSURE THAT IT IS LESS THAN 100 C AND IS NOT NEGATIVE. IF IT DOES NOT MEET THESE REQUIREMENTS, C A MESSAGE IS WRITTEN FOR THE USER, AND THEY ARE ASKED TO C ENTER A REASONABLE VALUE. C IF ((NREP.LT.1).OR.(NREP.GT.99)) THEN IF (NREP.LT.1) THEN WRITE(9,*)”TBIS NUMBER IS NEGATIVE!” ELSE WRITE(9,*)”THIS NUMBER IS TOO LARGE!” ENDIF GOTO 15 ENDIF C C THE NAME OF THE MODEL AND THE INPUT AND OUTPUT FILES ARE C REQUESTED FROM THE USER. C 1 36 DOS I=1,NREP wPJTE(9,*)”ENTER THE IMAGE FILE NAME” READ(9,*)INP(I) WPJTE(9*)”ENTER THE OUTPUT FILE NAME” pEAD(9,*)OUTp(I) WRITE(9,*)”ENTER THE NAME OF THE MODEL” PEAD(9,*)MDL(I) 5 CONTINUE C C THE LOOP WHICH WILL ALTER AND WRITE OUT EACH IMAGE IS C BEGUN. C DO 380 NINC = 1,NREP C C THE NINC (FIRST,SECOND...) INPUT FILE IS OPENED. C OPEN(5 ,FILE=INP(NINC)) C C THE ARRAYS IN THE PROGRAM ARE IN1TIALIZED TO ZERO. C DO 10 I=1,NTX DO 1OJ=1,NTY IMG(I,J) = 0 10 CONTINUE C C THE IMAGE IS NOW READ IN FROM THE INPUT FILE. C WRITE(9,*)”IMAGE “,NINC,” IS BEING ANALYSED” DO 30 I=4,NX+3 READ(5,9 1)(IMG(I,J),J=4,NY+3) 30 CONTINUE 91 FORMAT(100I4) C C THE PORE/SOLID, SEPARATION, CTh, AND BOUNDARY VALUES FOR THE C IMAGE FILE ARE INITIALIZED. C POR = 255 SOL=0 BDRY =50 SEPR=30 IMAX =63 CTR=200 C C THE “WATERSHED ALGORITHM IS NOW APPLIED. THIS SETS THE VALUES C OF THE PIXELS WITHIN THE PORES TO PROGRESSIVELY LOWER OR HIGHER C VALUES DEPENDING ON THEIR PROXIMITY TO THE SOLID. C DO 40 I=50,TMAX ITMP = I DO 50 J= 3, NX+4 DO 50 K= 3, NY÷4 1 37 IF (I.EQ.50) THEN IF (IMG(J,K).LT.SEPR) THEN IF (IMG(J+i,K).EQ.POR) IMG(J÷1,K) = ITMP IF (IMG(J-i,K).EQ.POR) IMG(J-1,K) = ITMP IF (IMG(J,K÷1).EQ.POR) IMG(J,K+i) = ITMP IF (IMG(J,K-1).EQ.POR) IMG(J,K-1) = ITMP ENDIF ELSE IF (IMG(J,K).EQ.(ITMP-1)) THEN IF (IMG(J+l,K).EQ.POR) IMG(J÷l,K) = ITMP IF (IMG(J-1,K).EQ.POR) IMG(J-1,K) = ITMP IF (IMG(J,K+1).EQ.POR) IMG(J,K-i-1) = ITMP IF (IMG(J,K-1).EQ.POR) IMG(J,K-1) = ITMP ENDIF ENDIF 50 CONTINUE 40 CONTINUE DO 360 MAX = 51,IMAX,4 SEPR =30 C C ALL PIXELS SET TO VALUES BETWEEN 50 AND MAX-i ARE RE-SET C TO SEPR. C DO 130 I=3,NX+4 DO 130 J=3,NY+4 IF (IMG(I,J).GT.SEPR) THEN IF (IMG(I,J).LT.MAX) THEN IMG(I,J) = SEPR ENDIF ENDIF 130 CONTINUE C C NOW THE WATERSHED ALGORITHM IS RE-APPLIED TO BUILD THE PORES C UP BASED ON THE REMAINING LAYER OF MAX’S IN THE PORES. C NTMP =48 DO 140 I=MAX-1,NTMP,-1 NCHK=I+i 1TMP I DO 150J=3,NX+4 DO 150 K= 3, NY-1-4 IF (IMG(J,K).EQ.NCHK) THEN IF (IMG(J+1,K).EQ.SEPR) IMG(J+1,K) = ITMP iF (IMG(J-1,K).EQ.SEPR) IMG(J-1,K) = ITMP IF (IMG(J,K+1).EQ.SEPR) IMG(J,K+1) = ITMP IF (IMG(J,K-1).EQ.SEPR) IMG(J,K-1) = ITMP ENDIF 150 CONTINUE 140 CONTINUE C C THE PIXELS NEWLY ALTERED TO THE PREVIOUS VALUE OF SEPR (30) C ARE CHECKED TO SEE WHAT THEIR SURROUNDINGS ARE. IF THERE IS 1 38 C APOREONONESIDEANDEITHERAPORETHROATORSOLIDONTHE C OTHER, THE PIXEL IS SET TO 15. IF ThERE WAS NOT A PORE THROAT C OR SOLID PRESENT, iT IS SET TO THE CURRENT SEPR VALUE, 20. C ALSO, IF THERE IS ONLY EiTHER SOLID OR PORE THROAT CLOSE BY, C THE PIXEL IS SET TO 20. C SEPR =20 DO 160 I=3,NX-i-4 DO 160 J=3,NY+4 IF (IMG(I,J).EQ.30) THEN NSOL=0 NPOR=0 Nil =0 DO 170 K=I-3,I+3 IF (IMG(K,J).GT.47) THEN NPOR= 1 ELSE IF (IMG(K,J).LT.1 1) THEN NSOL=i ELSE IF(IMG(K,J).EQ.11) Nil = 1 ENDIF ENDIF 170 CONTINUE IF (NPOR.EQ.1) THEN IF (NSOL.NE.1).OR.(Nll.EQ.1) THEN IMG(I,J) = SEPR GOTO 160 ELSE IMG(I,J) = 15 ENDIF ELSE IF (NSOL.EQ.1).OR.(Nll.EQ.1) THEN IMG(I,J) = SEPR ENDIF ENDIF NPOR=0 DO 180 K=J-3,J+3 IF (IMG(I,K).GT.47) THEN NPOR=l ELSE IF (IMG(I,K).LT.ll) THEN NSOL=NSOL + 1 ELSE IF(IMG(K,J).EQ.1l) Nil = Ni 1+1 ENDIF ENDIF 180 CONTINUE IF (NPOR.EQ.1) THEN IF (NSOL.NE.i).OR.(Ni1.NE.l) THEN IMG(I,J) = SEPR ELSE IMG(I,J) = 15 ENDIF 1 39 ELSE IF (NSOL.GT.1) THEN IMG(I,J) = 15 ENDIF IF (N11.GT.1) THEN IMG(I,J) = SEPR ENDIF ENDIF 12 ENDIF 160 CONTINUE C C PIXELS SET TO 30 ARE NOW CHECKED SO THAT ANY IN CONTACT C WITH PiXELS SET TO SEPR ARE ALSO SET TO SEPR. AS THE C SURROUNDING OF EACH SEPR PIXEL IS CHECKED, THE PIXEL IS C SET TO 25 SO THAT iT WILL ONLY BE CHECKED ONCE. C DO 190 L=1,200 NSEP0 DO 200 I=3,NX+4 DO 200 J=3,NY-i-4 IF (IMG(I,J).EQ.SEPR) THEN IMG(I,J) = 25 DO 210 K=I-1,I-i-1,2 IF (IMG(K,J).EQ.30).OR.(IMG(K,J).EQ.15) THEN IMG(K,J) = SEPR NSEP = 1 ENDIF 210 CONTINUE DO 220 K=J-1,J+1,2 IF (IMG(I,K).EQ.30).OR.(IMG(I,K).EQ. 15) THEN IMG(I,K) = SEPR NSEP = 1 ENDIF 220 CONTINUE ENDIF 200 CONTINUE IF (NSEP.EQ.0) GOTO 185 190 CONTINUE 185 CONTINUE C C PIXELS CHANGED TO 25 TO SPEED THE PREVIOUS LOOP ARE C RE-SET TO 20. C DO 195 I=4,NX÷3 DO 195 J=4,NY+3 IF (IMG(I,J).EQ.25) IMG(I,J)=20 195 CONTINUE C C AFTER EACH LOOP IN WHICH MAX CHANGES, THE PIXELS WHICH WERE C PREVIOUSLY SET TO BE PORE THROATS ARE RE-SET TO 11. IF MAX C IS GREATER THAN 60, ThEN ANY NEWLY CREATED PORE-THROATS IN C CONTACT WITH AN OLD PORE THROAT IS RE-SET TO 10, IE: THEY C ARE NOT ALLOWED TO REMAIN PORE THROATS. ALSO, ANY PIXELS C IN CONTACT WITH PIXELS SET TO 10 ALREADY ARE SET TO 10. 140 C IF (MAX.GT.60) THEN M1=2 M2= 1 INC=-1 ELSE M1=1 M2=2 INC= 1 ENDIF DO 485 M=M1,M2,INC rRfl’E(9,*)”M = N=9+M N3 = N N2= 10+2*M IF (MAX.GT.60) THEN N3=10 ENDIF DO 490 L=1,200 NSEP =0 DO 500 I=3,NX+4 DO 500 J=3,NY÷4 IF (IMG(I,J).EQ.N) THEN IMG(I,J) = N2 DO 510 K=I-1,I+1,2 IF (IMG(K,J).GT. 15) .AND. (IMG(K,J).LT. 31) THEN IMG(K,J) N3 NSEP = 1 ENDIF 510 CONTINUE DO 520 K=J-1,J÷1,2 IF (IMG(I,K).GT. 15).AND.(IMG(I,K).LT.3 1) THEN IMG(I,K) = N3 NSEP = 1 ENDIF 520 CONTINUE ENDIF 500 CONTINUE IF (NSEP.EQ.0) GOTO 485 490 CONTINUE 485 CONTINUE C C IN THE PREVIOUS LOOP, AFTER EACH PIXEL OF VALUE 10 OR 11 C WAS CHECKED, IT WAS RE-SET TO 12 OR 14(10 - 12 AND 11-14) C SO THAT IT WOULD NOT BE CHECKED TWICE. NOW, THE PIXELS C WHICH WERE SET TO 12 AND 14 ARE RE-SET TO 10 AND 11. C DO 550 I=4,NX+3 DO 550 J=4,NY+3 IF (IMG(I,J).EQ.12) THEN IMG(I,J) = 10 ELSE IF (IMG(I,J).EQ.14) THEN 141 IMG(I,J) = 11 ENDIF ENDIF 550 CONTINUE C C PIXELS STILL SET TO 30 OR 15 ARE SET TO 10, THOSE SET TO C SEPRARESETTO11. C DO 240 I=3,NXi-4 DO 240 J=3,NY+4 IF (IMG(I,J).EQ.20) THEN IMG(I,J) = 11 ELSE if (IMG(I,J).EQ.30).OR.(IMG(I,J).EQ. 15) THEN IMG(I,J) =10 ENDIF ENDIF 240 CONTINUE 360 CONTINUE C wPJTE(9,*y’FINISHED SEPARATING PORES IN IMAGE “,NINC C C ALL PIXELS STILL SET TO 10 ARE RE-SET TO PORE VALUES, C AS ARE ANY PIXELS WITH VALUES GREATER THAN 30. C THOSE PIXELS SET TOll ARE SET TO SEPR (20) C DO 250 I=3,NX-t-4 DO 250 J=3,NY÷4 IF (IMG(I,J).EQ. 10).OR.(IMG(I,J).GT.30) THEN IMG(I,J) = 255 ELSE if (IMG(I,J).EQ.11) THEN IMG(I,J) = 20 ENDIF ENDIF 250 CONTINUE C C THIS NEXT SECTION CHECKS TO ENSURE EACH PORE THROAT IS AT C LEAST TWO PIXELS IN WIDTH. C DO 270 M=1,2 DO 260 I=4,NX+3 DO 260 J=4,NY+3 if (IMG(I,J).EQ.SEPR) THEN NPOR=0 IF (IMG(I+1,J).EQ.POR) NPOR=NPOR+1 IF (JMG(I-l,J).EQ.POR) NPOR=NPOR+1 IF (NPOR.EQ.2) THEN IMG(I+1,J) = SEPR IMG(I- 1 ,J) = SEPR ENDIF NPOR=0 if (IMG(I,J÷1).EQ.POR) NPOR=NPOR+1 1 42 IF (IMG(I,J-1).EQ.POR) NPOR=NPOR+1 IF (NPOR.EQ.2) THEN IMG(I,J÷1) = SEPR IMG(I,J-1) = SEPR ENDIF ENDIF 260 CONTINUE 270 CONTINUE C C ALL PIXELS STILL SET TO 20 ARE ASSUMED TO BE WITHIN C THE PORE THROATS AND ARE LEFT SET TO THIS VALUE. C C C THE SUBROUTINE PRO IS NOW CALLED TO LOCATE AND RETURN C THE LOCATIONS OF THOSE POINTS iN THE PORES WHICH ARE C FARTHEST FROM THE EDGES OF THE PORES, JE THE CENTERS C OF THE PORES. C 719 SEPR=20 TMP = MDL(NINC) CALL PR(IMG,POR,TMP) WRITE(9,*)”AFTER PR” C C THE SUBROUTINE PRO IS NOW CALLED TO LOCATE AND RETURN C THE LOCATIONS OF THOSE POINTS IN THE PORE THROATS WHICH ARE C FARTHEST FROM THE EDGES OF THE PORE THROATS, IE THE CENTERS C OF THE THROATS. C WRITE(9,*)”BEFORE PR” CALL PR(IMG,SEPR,TMP) WRITE(9,*)”AFTER PR” C OPEN(6,FILE=OUTP(NINC)) C C THE SUBROUTINE PRSZ IS CALLED TO DETERMINE THE SIZE OF THE C PORES AND PORE THROATS PRESENT IN THE IMAGE. C WRITE(9,*)’?BEFORE PRSZ” CALL PRSZ(IMG,PRDIST,THDIST,NPR,NTHR,MXRAD) WRITE(9,*)”AFrER PRSZ” C C THE SUBROUTINE SA IS CALLED TO DETERMINE THE POROSITY, PORE C SURFACE TO VOLUME RATION, PORE THROAT SURFACE TO VOLUME RATIO, C AND THE SURFACE TO VOLUME RATIO OF THE TOTAL VOID SPACE. C CALL SA(IMG,POROS,PSV,TSV,THSV) wRrrE(9*)”AFTER SA” C C THE POROSITY, SURFACE TO VOLUME RATIOS, AND THE PORE AND PORE C THROAT SIZES ARE WRHTEN INTO THE FILE PRSZ.”MDL” C OUTP(NINC) = “PRSZ.”//MDL(NINC) OPEN(7,FILE=OUTP(NINC)) 143 WR1TE(7, 109) 109 FORMAT(”THE CURRENT IMAGE IS “,A20) WRITE(7,111) POROS 111 FORMAT(”POROSITY OF MODEL “,F12.2) WRITE(7,1 12) PSV,THSV,TSV 112 FORMAT(”SA/VOL RATIO FOR PORES ,F8 5/,”SA/VOL RATIO OF & PORE THROATS”,F8 5/,”SA/VOL RATIO OF TOTAL VOID SPACE tt,F8 5/) wRITE(7,*)”RAD PORE MN THROAT WRITE(7,95) 95 FORMAT(’ DIST DIST”,/i) DO 100 I=1,MXRAD WRITE(7,96)I,PRDIST(I),THDIST(I) 96 FORMAT(14,6X,15,4X,I5,4X,14) 100 CONTINUE CLOSE(7) C C THE IMAGE IS WRITTEN TO THE OUTPUT FILE. C DO 370 I=4,NX+3 WRITE(6,9 1 )(IMG(I,J),J=4,NY÷3) 370 CONTINUE C C THE INPUT AND OUTPUT FILES ARE CLOSED. C 13 CLOSE(S) CLOSE(6) C C THE NEXT STATEMENT IS THE END OF THE LOOP WHICH C IS REPEATED FOR EACH IMAGE. C 380 CONTINUE STOP END C C THE SUBROUTINES CALLED IN THE MAIN PROGRAM ARE INCLUDED. C PR.F CONTAINS PRO AND FINDS THE CENTERS OF PORES AND PORE C THROATS, WRITES OUT THE CENTERS AND CHANGES THE PIXEL C CORRESPONDING TO A CENTER TO THE VALUE 200. C SA.F CONTAINS SAO WHICH DETERMINES THE SURFACE TO VOLUME C RATIO OF THE PORE SPACE AND THE POROSITY OF THE IMAGE. C PORSIZ.F CONTAINS PRSZQ AND DETERMINES THE MAX RADIUS OF A C CIRCLE WHICH IS CONTAINED ENTIRELY WITHIN THE PORE OR PORE C THROAT AND CENTERED ON A PIXEL OF VALUE 200. C SPHTHR.F CONTAINS SPHTHR() IS CALLED BY PRSZO TO DO THE C COMPARISON BETWEEN THE IMAGE AND A REFERENCE CIRCLE GIVEN C THE COORDINATES OF THE CENTER OF THE PORE OR PORE THROAT. C 144 INCLUDE PR.F INCLUDE SA.F INCLUDE PORSIZ.F INCLUDE SPHTHR.F 145 PRSZ.F C C MARCH, 1992 ALICE CHAPMAN C ROCK PHYSICS LAB, UBC C C THIS SUBROUTINE ACCEPTS A CONFOCAL IMAGE FILE IN WHICH C THE PORES HAVE BEEN SEPARATED, AND THE PIXELS AT THE C CENTER OF PORES AND PORE THROATS HAVE BEEN RESET TO C THE VALUE 200. AN ARRAY CONTAINING A CIRCLE AND A C REFERENCE ARRAY ARE READ IN FROM THE FILE SPHS. C C THE CIRCLE IS PRODUCED BY SPHERE.F AND IS USED AS A C COMPARISON TO THE IMAGE TO DETERMINE THE SIZE OF THE C PORES OR PORE THROATS IN THE IMAGE GIVEN THE COORDINATES C OF THE POINT WITHIN THE PORES OR PORE THROATS WHICH IS C FARTHEST FROM THE EDGE. C C C C THE VARIABLES ARE DECLARED. C SUBROUTINE PRSZ(IMG,PRDIST,THDIST,NPR,NTHR,MXRAD) C INTEGER*2 IMG(530,780), SPH(200,200) INTEGER SOL,SEPR,RAD,RAD2,CTR,MXRAD INTEGER PRDIST( 101),THDIST(101) INTEGER NTHR,CHKSPH( 101) ,NPR REAL RFIT,MXFIT CHARACTER*20 MDL,INP C C*************************** C C VARIABLE DESCRIPTION: C C IMG(530,780): THE CONFOCAL IMAGE IS PLACED IN THIS ARRAY C WITH THREE ROWS OF SOLiDS SURROUNDING IT C SPH(200,200): AN ARRAY CONTAINING A SET OF CONCENTRIC C CIRCLES WITH THE OUTSIDE ROW OF PIXELS IN EACH C CIRCLE BEING NUMBERED (RAD+0.5) WERE RAD WAS THE C RADIUS OF THE CIRCLE IN QUESTION. IE: THE C CENTERMOST PIXEL IS NUMBERED 1 SINCE RAD = 0.5. C ALSO, THE CENTER IS AT SPH(100.5,100.5) C CHKSPH(I): A COMPARISON ARRAY WHICH RECORDS HOW MANY C PIXELS RE WITHIN A CIRCLE OF RADIUS (1-0.5) C SOL: THE NUMBER WHICH PIXELS WITHIN THE SOLID PORTION C OFTHEIMAGEWILLBESETTO(0) C POR: THE NUMBER WHICH PIXELS WITHIN THE PORES OF THE C IMAGE WILL BE SET TO (255) C SEPR: THE VALUE WHICH PIXELS IN THE PORE THROATS WILL C BE SET TO TO SEPARATE THE PORES (20) C CrR: THE NUMBER TO WHICH THE PIXEL AT THE CENTRE OF A C PORE OR PORE THROAT IS SET C NPR: THE NUMBER OF PORES WHOSE SIZE WAS MEASURED 146 C NTHR: TH NUMBER OF PORE THROATS WHOSE SIZE WAS MEASURED C RAD: THE SIZE OF THE CIRCLE WHICH FITS ENTIRELY WITHIN C THE PORE OR THROAT BEING MEASURED. (SIZE: ACTUAL C RADIUS +0.5) C PRDISTjTHDIST: AFTER THE SIZE OF THE PORE OR PORE THROAT C IS DETERMINED THE SIZE IS RECORDED IN PRDIST (PORES) C AND THDIST (PORE THROATS). PRDIST(RAD) OR THDIST(RAD) C IS INCREMENTED BY ONE AFTER EACH MEASUREMENT. C MXRAD: THE INTEGER REPRESENTING THE SIZE OF THE LARGEST PORE C SIZE MEASURED. MXRAD = RAD1US-i-0.5 C MDL: THE NAME OF THE MODEL, WHICH IS USED TO NAME SOME C OUTPUT FILES C C*************************** C C C THE WIDTH (NX) AND LENGTH (NY) OF THE CONFOCAL IMAGE ARRAY C AND THE WORKING ARRAY (NTX) AND (NTY) ARE SET. C NX = 512 NY = 768 NTX = NX + 6 NTY = NY +6 C C THE ARRAYS IN THE PROGRAM ARE INiTIALIZED TO ZERO. C DO 10 1=1,101 PRDIST(I) =0 THDIST(I) 0 CHKSPH(I) = 0 10 CONTINUE C C THE SOLID, POR, SEPARATION, AND PORETHROAT CENTER VALUES FOR C THE IMAGE FILE ARE INITIALIZED. C SOL=0 SEPR =20 POR = 255 CTR = 200 C C THE COMPARISON SPHERE AND ARRAY (CHKSPH) ARE NOW READ IN FROM C THEFILE”SPHS’ C INP = ‘SPHS’ OPEN(7,FILE=INP) DO 201=1,200 READ(7,91 )(SPH(I,J),J=1 ,200) 20 CONTINUE 91 FORMAT(10014) C DO 301=1,101 READ(7,94)J,CHKSPH(I) IF (I.NE.J) THEN 147 WPJTE(9*) “THERE IS A PROBLEM WITH READING IN THE” WRITE(9,*)”CHKSPH ARRAY FROM THE FILE: SPHS.” WRITE(9,*)”DO YOU WISH TO QUIT? IF SO, ENTER A 1” READ(9,*)K IF (K.EQ.1) THEN CLOSE(7) GOTO 11 ENDIF ENDIF 30 CONTINUE 94 FORMAT(218) C CLOSE(7) C C THE COUNTING VARIABLES USED TO KEEP TRACK OF THE MAXIMUM C RADIUS OF AND THE NUMBER OF PORE AND PORE THROATS MEASURED C ARE INITIALIZED TO ZERO. C MXRAD =0 NPR=0 NTHR =0 C C THE SIZE OF THE PORES AND PORE THROATS IN WHICH THE POINTS C FARTHEST FROM THE EDGES OF THE PORES HAVE BEEN SET TO CTR C ARE NOW MEASURED. C DO 40 I=4,NX+3 DO 40 12=4,NY+3 IF (IMG(I,12).EQ.200) THEN J=I K =12 C C SPHTHR IS CALLED TO DETERMINE THE SIZE OF THE LARGEST CIRCLE C CENTERED AT IMG(J,K) WHICH FITS ENTIRELY WITHIN THE PORE/THROAT. C THE SIZE WHICH IS RETURNED IS CHECKED AND, IF IT IS LARGER THAN C ANY PREVIOUSLY RETURNED, MXRAD IS RE-SET TO EQUAL RAD. C CALL SPHTHR(J,K,IMG,SPH,NREF,RFIT,CHKSPH) IMG(I,12) = 210 RAD = INT(RFIT) IF (RAD.GT.MXRAD) MXRAD = RAD IF (RAD.GT.0) THEN C C IF THE REFERENCE NUMBER RETURNED BY SPHTHR IS EQUAL TO POR, C THEN THE SIZE OF A PORE WAS MEASURED AND NPR AND PRDIST(RAD) ARE C INCREMENTED. IF NREF WAS EQUAL TO SEPR, A PORE THROAT WAS C MEASURED AND NTHR AND THDIST(RAD) ARE INCREMENTED. IF NREF WAS C 0, THOUGH, THIS INDICATES THE PORE OR PORE THROAT WAS IN CONTACT C WITH THE EDGE OF THE IMAGE AND THEREFORE ITS MEASUREMENT WILL C NOT BE RECORDED. C IF (NREF.EQ.POR) THEN NPR=NPR+ 1 PRDIST(RAD) = PRDIST(RAD) + 1 148 ELSE IF (NREF.EQ.SEPR) THEN NTHR = NTHR + 1 THDIST(RAD) = THDIST(RAD) + 1 ENDIF ENDIF END1F ENDIF 40 CONTINUE C C THE IMAGE IS NOW CHECKED FOR PORE/PORE THROATS THAT WERE C MEASURED BUT WHOSE CENTERS WERE LATER INCLUDED IN ANOTHER C PORE AND THEREFORE THEMEASUREMENT MUST BE REMOVED FROM THE C DISTRIBUTION. C DO 50 I=4,NX÷3 DO 50 12=4,NY+3 IF (IMG(I,I2).EQ.220) THEN J=I K=12 C C SPHTHR IS CALLED TO DETERMINE THE SIZE OF THE LARGEST CIRCLE C CENTERED AT IMG(J,K) WHICH FITS ENTIRELY WITHIN THE PORE/THROAT. C CALL SPHTHR(J,K,IMG,SPH,NREF,RFIT,CHKSPH) IMG(I,I2) = NREF IF (NREF.EQ.0) THEN IMG(I,I2) = IMG(I+1,12) ENDIF RAD = INT(RFIT) IF (RAD.GT.0) THEN C C IF THE REFERENCE NUMBER (NREF) RETURNED BY SPHTHR IS EQUAL TO C POR, THEN THE SIZE OF A PORE WAS MEASURED AN]) NPR AND C PRDIST(RAD) ARE DECREASED BY 1. IF NREF WAS EQUAL TO SEPR, A PORE C THROAT WAS MEASURED AND NTHR AN]) THDIST(RAD) ARE DECREASED BY C 1. IF NREF WAS 0, THOUGH, THIS INDICATES THE PORE OR PORE THROAT C WAS IN CONTACT WITH THE EDGE OF THE IMAGE AND THEREFORE ITS C MEASUREMENT WILL NOT BE RECORDED. C IF (NREF.EQ.POR) THEN NPR=NPR- 1 PRDIST(RAD) = PRDIST(RAD) - 1 ELSE IF (NREF.EQSEPR) THEN NTHR = NTHR -1 THDIST(RAD) = THDIST(RAD) - 1 ENDIF ENDIF ENDIF ELSE IF (IMG(I,12).EQ.210) THEN IMG(I,12) = 200 ELSE 149 IF (IMG(I,12).EQ.205) THEN IF (IMG(I-1,12).GT.O) THEN NREF = IMG(I- 1,12) ELSE IF (IMG(I+1,12).GT.O) THEN NREF = IMG(I+1,12) ENDIF ENDIF IMG(I,12) NREF ENDIF ENDIF ENDIF 50 CONTINUE 11 RETURN END 150 PR.F C C MARCH, 1992 ALICE CHAPMAN C ROCK PHYSICS LAB, UBC C C THIS SUBROUTINE TAKES CONFOCAL IMAGES IN WHICH C THE PORES HAVE BEEN SEPARATED AND FINDS THE PIXELS C WiTHIN THE PORES OR PORE THROATS WHICH ARE FARTHEST C FROM THE EDGES. IF THE VALUE OF NREF IS 255, THE C VALUE CORRESPONDING TO PORES N THE IMAGE, THE CENTER C OFTHEPORESAREFOUNDAND,IFITIS2O,THECENTEROF C THE PORE THROATS ARE FOUND. C C THE POINT IN THE CENTERMOST PORTION OF THE PORES AND C PORE THROATS IS FOUND USING THE WATERSHED ALGORITHM. C N THIS ALGORITHM, THE PIXELS ON THE OUTER EDGES OF C THE PORE/THROAT ARE SET TO 50, THEN INNER ONES IN CONTACT C WITH THE 50’S ARE SET TO 51...ETC., UNTIL THE PORES/THROATS C ARE COMPLETELY RE-NUMBERED. THEN THOSE PIXELS GREATER C THAN 50 WHICH DO NOT HAVE A LARGER NUMBER WITHIN A BOX C WHICH EXTENDS 5 PIXELS ON EACH SIDE OF THEM ARE ASSUMED C TO BE NEAR THE CENTERS. THE PIXELS UP TO 20 PIXELS ON C EITHER SIDE OF SUCH A PIXEL ARE CHECKED TO DETERMINE C IF THEY ARE PART OF A CLUSTER OF SIMILARLY VALUED PIXELS C INCLUDiNG THE ORIGINAL PIXEL. IF SO, THE COORDINATES OF C THE PIXELS IN THE CLUSTER ARE AVERAGED TO DETERMINE THE C ACTUAL CENTER OF THE PORE/THROAT. C C C C SUBROUTINE PR(IMG,NREF,MDL) C C THE VARIABLES ARE DECLARED. C INTEGER*2 IMG(530,780) INTEGER POR,SOL,SEPR,CTR( 10000,2),BDRY INTEGER NSPH,MAX,NREF CHARACTER*20 MDL,OUTP C C*************************** C C VARIABLE DESCRIPTION: C C IMG(530,780): THE CONFOCAL IMAGE IS PLACED IN THIS ARRAY C WITH THREE ROWS OF SOLIDS SURROUNDING IT C POR: THE NUMBER WElCH PIXELS WITHIN THE PORES OF THE C IMAGE ARE SET TO (255) C SOL: THE NUMBER WHICH PIXELS WITHIN THE SOLID PORTION OF C THE IMAGE ARE SET TO (0) C SEPR: THE NUMBER WHICH PIXELS IN THE SMALLEST PART OF THE C PORE THROATS ARE SET TO TO SEPARATE THE PORES (20) C NREF:THE NUMBER WHICH INDICATES WHETHER PORE OR PORE THROAT 151 C CENTERS WILL BE DETERMINED. IF NREF = POR, THEN PORE C CENTERS WILL BE FOUND, AND IF NREF = SEPR, PORE THROAT C CENTERS WILL BE FOUND. C CTR(K,I): A WORKING ARRAY WHICH HOLDS THE COORDINATES OF THE C PIXELS WiTHIN THE PORES OR PORE THROATS WHICH ARE C FARTHEST FROM THE EDGES OF THE PORES. C NSPH:THE TOTAL NUMBER OF PORES OR THROATS WHOSE CENTERS ARE C RECORDED IN CTR. C BDRY:THE INfl1AL NUMBER WHICH PORE PIXELS BORDERING THE SOLID C ARE SET TO IN THE WATERSHED ALGORITH1vI. PORE PIXELS WHICH C BORDER THOSE WHICH ARE SET TO BDRY ARE THEN SET TO BDRY-i-1 C ETC, ETC. C MAX: THE LARGEST NUMBER USED IN THE WATERSHED ALGORITHM. (IE: C THE LARGEST NUMBER WHICH A PORE PIXEL IS RE-SET TO. C MDL: THE NAME OF THE MODELS, USED TO NAME OUTPUT FILES C C*************************** C C THE WIDTH (NX) AND LENGTH (NY) OF THE CONFOCAL IMAGE ARRAY C AND THE WORKING ARRAY (NTX) AND (NTY) ARE SET. C NX = 512 NY = 768 NTX = NX + 6 NTY = NY +6 C C THE ARRAYS IN THE PROGRAM ARE INITIALIZED TO ZERO. C DO 201=1,2 DO 20 J=1,10000 CTR(J,I) = 0 20 CONTINUE C C THE PORE/SOLID, SEPARATION, AND BOUNDARY VALUES FOR THE IMAGE C FILE ARE INITIALIZED. C POR=255 SOL=0 SEPR =20 BDRY =50 C C ThE “WATERSHED ALGORiTHM IS NOW APPLIED. THIS SETS THE VALUES C OF THE PIXELS WITHIN THE PORES TO PROGRESSWELY LOWER OR HIGHER C VALUES DEPENDING ON THEIR PROXIMITY TO THE SOLID. C WRITE(9,*)”STARTING WTSHD PORTION” DO 40 I=BDRY,255 NTST=0 ITMP = I DO 50 J = 4,NX+3 DO 50 K= 4,NY+3 IF (I.EQ.BDRY) THEN IF (IMG(J,K).EQ.POR).OR.(IMG(J,K).LT.25) THEN IF (IMG(J,K).NE.NREF) THEN 1 52 NTST= 1 IF (IMG(J+1,K).EQ.NREF) IMG(J+1,K) = ITMP iF (IMG(J-1,K).EQ.NREF) IMG(J-1,K) = ITMP IF (IMG(J,K+1).EQ.NREF) IMG(J,K+1) = ITMP IF (IMG(J,K-1).EQ.NREF) IMG(J,K-1) = ITMP ENDIF ENDIF ELSE IF (IMG(J,K).EQ.(ITMP-1)) THEN NTST= 1 IF (IMG(J+1,K).EQ.NREF) IMG(J+1,K) = ITMP IF (IMG(J-1,K).EQ.NREF) IMG(J-1,K) = ITMP IF (IMG(J,K+1).EQ.NREF) IMG(J,K+1) = ITMP IF (IMG(J,K-1).EQ.NREF) IMG(J,K-1) = ITMP ENDIF ENDIF 50 CONTINUE IF (NTST.EQ.0) THEN MAX = I-i GOTO 111 ENDIF 40 CONTINUE C C THE POSSIBLE CENTERS OF PORES ARE NOW SELECTED, THEIR C COORDINATES STORED IN THE ARRAY CR C 111 NSPH=0 WPJTE(9,*)”PAST WTSHD PORTION” DO 751= MAX,BDRY,-1 DO 70 J=4,NXi-3 DO 70 K=4,NY+3 C C EACH PIXEL IN THE IMAGE IS CHECKED UNTIL ONE WHICH HAS C BEEN ALTERED BY THE WATERSHED ALGORITHM IS FOUND. THE C PARAMETER NPOR IS THEN SET TO THE VALUE OF THIS PIXEL C AND THE TEST VARIABLE, 1TST IS SET TO 0. C IF (IMG(J,K).EQ.I) THEN NPOR = IMG(J,K) ITST=0 C C THIS LOOP CHECKS A BOX IN THE IMAGE WHICH EXTENDS 5 PIXELS C TO EACH SIDE OF THE PIXEL SELECTED. IF THE VARIABLES J2, C K2 GO OUT OF RANGE, THE PIXEL IS NOT CHECKED, OTHERWISE iT C IS CHECKED TO DETERMINE IS IT IS LARGER THAN THE CHOSEN C PIXEL’S VALUE OF NPOR. IF SO, ITST IS INCREMENTED BY 1. C DO 80 J2=J-5,J-t-5 DO 80 K2=K-5,K+5 IF ((J2.LT.l).OR.(J2.GT.NTX)) GOTO 80 IF ((K2.LT.1).OR.(K2.GT.NTY)) GOTO 80 IF (IMG(J2,K2).GT.NPOR) THEN 1 53 ITST = ITST+1 ENDIF 80 CONTINUE C C IF THERE WERE NO PIXELS LARGER THAN THE ONE SELECTED IN C THEAREAAROUNDIT,THENITISASSUMEDTOBEATORNEAR C A PORE/PORE THROAT CENTER AND ThE NUMBER OF CENTERS FOUND C , NSPH, IS INCREMENTED BY ONE. THE POSSIBLE CENTER, OR C ORIGIN OF THE PORE/THROAT IS THEN SET TO THE COORDINATES C OF THE CHOSEN PIXEL. (ORU = J, ORIK = K) NORI IS THE C NUMBER OF PIXEL POSITIONS AVERAGED TO DETERMINE THE ACTUAL C CENTER OF THE PORE THROAT AND IS INiTIALIZED TOP 1. C IF (ITST.EQ.0) THEN NSPH = NSPH + 1 ORIJ=J ORIK=K NORI= 1 C C THE AREA AROUND THE PIXEL IN A BOX EXTENDING UP TO 20 PIXELS C TO EITHER SIDE IS NOW CHECKED AND ANY PIXELS OF THE SAME C VALUE ARE AVERAGED TO DETERMINE THE ACTUAL CENTER OF THE C PORE/THROAT. ONCE THE PIXEL HAS BEEN ALTERED, iT IS RE-SET TO C 240 SO THAT ThE SAME PORE/THROAT CENTER WILL NOT BE “FOUND” FROM C MORE THAN ONE ORIGINATING PIXEL. C FOR EACH LOOP IN WHICH L IS INCREMENTED, THE PORTION OF THE C BOX TO BE CHECKED IS ONE PIXEL FARTHER FROM THE CENTER THAN C BEFORE AND WITHIN EACH LOOP, ITST IS ORIGINALLY SET TOO, AND C IF A PIXEL OF VALUE NPOR IS FOUND, RE-SET TO 1. IF NO PIXELS C OF VALUE NPOR ARE FOUND, THE LOOP IS TERMINATED. C DO 90 L=1,20 N2=L*2 1TST=0 DO 100 J2=J-L,J+L,N2 DO 100 K2=K-L,K+L C C THE NEXT TWO STATEMENTS CHECK IF J2 OR K2 IS OUT OF RANGE C OF THE ARRAY BOUNDARIES AND, IF SO, SKIPS CHECKING THE C VALUES OF THE PIXELS AT THOSE POSITIONS. C IF ((J2.LT.2).OR.(J2.GT.NTX-1)) GOTO 100 IF ((K2.LT.2).OR.(K2.GT.NTY-1)) GOTO 100 IF (IMG(J2,K2).EQ.NPOR) THEN DO 105 J3=J2-1,J2+1 DO 105 K3 = K2-1,K2+1 IF (IMG(J3,K3).EQ.240) THEN ORIJ = (ORIJ*REAL(NORI) + REAL(J2))/(1.O -i-REAL(NORI)) ORIX = (ORIK*REAL(NORI) + REAL(K2))/(1.O +REAL(NORI)) NORI=NORI÷ 1 ITST= 1 GOTO 100 ENDIF 105 CONTINUE 1 54 ENDIF 100 CONTINUE DO 110 K2=K-L,K+L,N2 DO 110 J2=J-L+1 ,J+L- 1 C C THE NEXT TWO STATEMENTS CHECK IF J2 OR K2 IS OUT OF RANGE C OF THE ARRAY BOUNDARIES AND, IF SO, SKIPS CHECKING THE C VALUES OF THE PIXELS AT THOSE POSITIONS. C IF ((J2.LT.2).OR.(J2.GT.NTX-1)) GOTO 110 IF ((K2.LT.2).OR.(K2.GT.NTY-1)) GOTO 110 IF (IMG(J2,K2).EQ.NPOR) THEN DO 115 J3=J2-1,J2+1 DO 115 K3=K2-1,K2-f-1 IF (IMG(J3,K3).EQ.240) THEN ORIJ = (ORIJ*REAL(NORI) + REAL(J2))/(1.0 -i-REAL(NORI)) ORIK = (ORIK*REAL(NORI) + REAL(K2))/(1.0 +REAL(NORI)) NORI=NORI+ 1 ITST = 1 IMG(J2,K2) = 240 GOTO 110 ENDIF 115 CONTINUE ENDIF 110 CONTINUE C C IF ITST IS STILL EQUAL TOO, THEN NO PIXELS OF VALE NPOR HAVE C BEEN LOCATED IN THIS LOOP AND ALL PIXELS SURROUNDING THE PORE C CENTER ARE ASSUMED TO HAVE BEEN FOUND SO THE LOOP IS C TERMINATED. C IF (ITST.EQ.0) GOTO 121 90 CONTINUE C C THE J AND K COORDINATES OF THE CENTER OF THE PORE/THROAT ARE C ROUNDED UP IF THE DECIMAL PORTION IS GREATER THAN 0.5, DOWN C OTHERWISE AND THEN CHANGED TO INTEGERS AND PLACED IN THE CTR C ARRAY. C 121 IF ((ORIJ + O.5).GT.(INT(ORIJ+O.5))) THEN ORU = ORIJ + 0.5 ENDIF IF ((ORIK + O.5).GT.(INT(ORIK+0.5))) THEN ORIK=ORIK+0.5 ENDIF CrR(NSPH,1) = INT(ORIJ) CTR(NSPH,2) = INT(ORIK) IMG(INT(ORIJ),INT(ORIK)) = 200 ENDIF ENDIF 70 CONTINUE 75 CONTINUE WRITE(9,*)“FINISHED WITH CENTERS” C 1 55 C THE ARRAY CTR IS NOW CHECKED TO ENSURE THAT REPEATED CTR C COORDINATES ARE REMOVED AND THAT ANY COORDiNATES WHICH C CORRESPOND TO PORE THROATS ARE ALSO REMOVED. C DO 160 L= 1,NSPH IF (CTR(L,1).GT.0) THEN IF (CTR(L,2).GT.0) THEN M = CTR(L,1) N = CTR(L,2) DO 170 J=L÷1,NSPH IF (CTR(J,1).EQ.M) THEN IF (CTR(J,2).EQ.N) THEN CTR(J,1) =0 CTR(J,2) =0 ENDIF ENDIF 170 CONTiNUE ENDIF ENDIF 160 CONTINUE C C THE ARRAY CTR IS NOW CHECKED TO ENSURE THAT CTR COORDINATES C WHICH ARE NOT VIABLE (IE: LT 0 OR EQUAL TO ANOTHER CUR) ARE C REMOVED. C NCHK=0 DO 180 I=1,NSPH IF (CTR(I,1).GT.0) THEN IF (CTR(I,2).GT.0) THEN NCHK = NCHK + 1 IF (NCHK.LT.I) THEN CTR(NCHK,1) = CTR(I,1) CTR(NCHK,2) = CTR(I,2) CTR(I,1) =0 CTR(I,2) = 0 ENDIF ENDIF ENDIF 180 CONTINUE NSPH = NCHK C C THE PIXELS WHICH HAVE NOT BEEN RE-SET TO THEIR ORIGINAL C VALUE OF NREF ARE NOW RE-SET. C DO 300 I=4,NX+3 DO 300 J=4,NY+3 IF (IMG(I,J).GT.SEPR).AND.(IMG(I,J).NE.200) THEN IF (IMG(I,J).NE.POR) THEN IMG(I,J) = NREF ENDIF ENDIF 300 CONTINUE 1 56 C C THE PIXELS WHICH CORRESPOND TO THE CENThR COORDINATES C DETERMINED WITH THIS PROGRAM ARE NOW SET TO THE VALUE C 200 C DO 310 I=1,NSPH J = CTR(I,1) K = CTR(I,2) IMG(J,K) = 200 310 CONTINUE RETURN END 1 57 SA.F C C MARCH, 1992 ALICE CHAPMAN C ROCK PHYSICS LAB, UBC C C THIS SUBROUTINE ACCEPTS A CONFOCAL IMAGE FILE iN WHICH THE C PORES HAVE BEEN SEPARATED. THE SURFACE AREA, VOLUME, C AND THE SURFACE AREA TO VOLUME RATIO OF THE ENTIRE IMAGE C IS DETERMINED AS WELL AS THE POROSITY. C C C*************************** C SUBROUTINE SA(IMG,POROS,PSV,TSV,THSV) C C THE VARIABLES ARE DECLARED. C INTEGER*2 JMG(530,780) INTEGER SOL,SEPR,TVOL,TSA,PVOL,THVOL INTEGER PSA,THSA REAL POROS ,PSV,TSV,THSV CHARACTER*20 MDL C C*************************** C C VARIABLE DESCRIPTION: C C IMG(530,780): THE CONFOCAL IMAGE IS PLACED IN THIS ARRAY C WITH THREE ROWS OF SOLIDS SURROUNDING IT C SOL: THE NUMBER WHICH PIXELS WITHIN THE SOLD) PORTION OF C THEIMAGEWILLBESETTO(0) C SEPR: THE NUMBER WHICH PIXELS IN THE SMALLEST PART OF THE C PORE THROATS WILL BE SET TO TO SEPARATE THE PORES (20) C POR: THE NUMBER WHICH PIXELS WITHIN THE PORES OF THE IMAGE C WILL BE SET TO (255) C CrR: THE NUMBER TO WHICH THE PIXEL AT THE CENTRE OF A PORE C OR PORE THROAT IS SET C PVOL: THE NUMBER (VOLUME) OF PORE PIXELS IN THE IMAGE C THVOL:THE NUMBER (VOLUME) OF PORE THROAT PIXELS IN THE IMAGE C TVOL: THE TOTAL NUMBER (VOLUME) OF PIXELS IN THE IMAGE C PSA: THE SURFACE AREA OF PORES IN THE IMAGE C THSA: THE SURFACE AREA OF PORE THROATS IN THE IMAGE C TSA: THE TOTAL SURFACE AREA OF PORES AND THROATS IN THE C POROS: THE POROSITY OF THE IMAGE: IE THE RATIO OF THE NUMBER C OF iMAGE PORE AND PORE THROAT PIXELS IN THE IMAGE TO C THE TOTAL NUMBER OF PIXELS IN THE IMAGE. C PSV: THE SURFACE TO VOLUME RATIO (PSA/PVOL) OF THE PORES C THSV: THE SURFACE TO VOLUME RATIO (THSAITHVOL) OF THE PORE C THROATS C TSV: THE OVERALL SURFACE TO VOLUME RATIO OF THE VOID SPACE C IN THE IMAGE: (PSA + THSA)/(PVOL + THVOL) C MDL: THE NAME OF THE MODEL, WHICH IS USED TO NAME SOME C OUTPUT FILES C 1 58 c*************************** C C THE WIDTH (NX) AND LENGTH (NY) OF THE CONFOCAL IMAGE ARRAY C AND THE WORKiNG ARRAY (NTX) AND (NTY) ARE SET. C NX = 512 NY = 768 NTX = NX + 6 NTY = NY +6 C C C THE SOLID, SEPARATION, PORE AND CENTER VALUES FOR THE IMAGE FILE C ARE INITIALIZED. C SOL=0 SEPR =20 POR=255 CTR=200 C C THE COUNTERS FOR THE VOLUMES AND SURFACE AREAS OF PORES AND PORE C ThROATS IN THE MODEL ARE IN111ALIZED TO ZERO. C TVOL=0 PVOL=0 THVOL=0 PSA =0 THSA=0 TSA=0 C C THE FOLLOWING LOOP CHECKS EACH PIXEL IN THE IMAGE AND INCREMENTS C THE VOLUMES AND SURFACE AREA COUNTERS APPROPRIATELY. C DO 10 I=4,NX+3 DO 10 J=4,NY+3 C C THE TOTAL VOLUME COUNTER (TVOL) IS INCREMENTED FOR EACH PIXEL. C TVOL = TVOL + 1 ITST=0 C C THIS SERIES OF IF/ThEN STATEMENTS SETS ITST TO 2 IF THE PIXEL IS C IN A PORE AND TO 3 IF iT IS IN A PORE THROAT. C IF (IMG(I,J).EQ.CTR) THEN DO 20 12=1-1,1+1 DO 20 J2=J- 1 ,J+ 1 IF (IMG(I,J).EQ.POR) THEN 1TST=2 ELSE IF (IMG(I,J).EQ.SEPR) 1TST =3 ENDIF 20 CONTINUE ELSE IF (IMG(I,J).EQ.POR) THEN 1 59 ITST =2 ELSE IF (IMG(I,J).EQ.SEPR) 1TST = 3 ENDIF ENDIF C C IFTHEPIXELWASINAPORETHROAT,ITSTIS3ANDTHVOLIS C INCREMENTED BY ONE. THE FOUR ADJOINING PIXELS ARE THEN CHECKED C AND, IF THE PIXEL CHECKED IS A SOLID, THSA IS INCREMENTED BY 1. C SIMILARLY, IT THE PIXEL WAS IN A PORE, ITST IS 2, AND PVOL IS C INCREMENTED BY ONE. THE FOUR ADJOINING PIXELS ARE THEN CHECKED C AND, IF THE PIXEL CHECKED IS A SOLID, PSA IS INCREMENTED BY 1. C IF (ITST.EQ.3) THEN THVOL = THVOL + 1 IF (IMG(I÷1,J).EQ.SOL) THSA = THSA + 1 IF (IMG(I-1,J).EQ.SOL) THSA = THSA + 1 IF (IMG(I,J+1).EQ.SOL) THSA = THSA + 1 IF (IMG(I,J-1).EQ.SOL) THSA = THSA + 1 ELSE IF (ITST.EQ.2) THEN PVOL = PVOL + 1 IF (IMG(I+1,J).EQ.SOL) PSA = PSA + 1 IF (IMG(I-1,J).EQ.SOL) PSA = PSA + 1 IF (IMG(I,J+1).EQ.SOL) PSA = PSA + 1 IF (IMG(I,J-1).EQ.SOL) PSA = PSA + 1 ENDIF ENDIF 10 CONTINUE C C THE TOTAL SURFACE AREA (TSA) OF THE VOID SPACE IS CALCULATED. C TSA = PSA + THSA C C THE POROSITY AND THE SURFACE AREA TO VOLUME RATIO OF THE PORES C AND PORE THROATS ARE CALCULATED AND PRINTED OUT. C TSV = REAL(TSA)/REAL(PVOL + THVOL) PSV REAL(PSA)/REAL(PVOL) THSV = REAL(THSA)/REAL(THVOL) POROS = REAL(PVOL + ThVOL)/REALfVOL) RETURN END 1 60 1.91. 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 :cj3 3 3 3 3 3 3 3 3 3 3 3 3 3 :N0IJpJ3SGIHVfl1VA 1Hivan flN’DIN‘TfN’HN’(loJ)HJSN”TcI’Nu‘XVNJI9aLNI I1N‘(JoT)HaSNHD‘1SL‘GVIITIDLLNI (ooVoo)Has‘(o8L’oc)ov’Izrnou’..ii (HasiH3’1IiI’iaIN’Has’oNI’IIN’IIN)llUHdSItll1flO1HflS ****************************** V13iDHfl‘fiVISDISAHJ)J3(j{:yJAJ1}{333J’pr S1V011iL 310d103ZIS31113MIELL3(I01A1TlVI’ 41ThdC135flSISITUIfN’TIN SLLVNJcflI0ODTc1UN3DLDNISflS3I0d3HLMIG35013N3KITWrIcIWO3 SIHDIHML313IIDIJiLI03ZIS3H1S3M1WHL3GN0IS’I3AIV’1fl3LLIVd SITUN0LLSflONIS1IXIJ3HL3C[fflDNI(flflOMH3fflA313II3IIU10 co÷sniciviKLoiGNdSIOD5I3mNnN3111•S3oI3ZioaNnoI0N3va VNOSI3HNflNMI(ULLV3IGMI3IVHJSNIS3’T3II3HL1V0IIU 3I0cT/3I0J311101HJSXIaLVN3111NIRLIMU31’IdcIflSS313II3‘EIHL S3IVJNO3‘1V0IHL3I0dI03110dVtIO‘IaLM333111I0SLVMCII0OD AQNVX3H1N3A19N31IA‘TllALflOIHflS/NVIO0JcISIRL (co-i)SflIUVN10fl3ll3VN11UJAkaIV1V0IRLflI0dI0310J 311110S’13X1JANVNAkOHSCfI033NHDIHAkAVfflVNOSIIIVdWO3y:(I)HdSN 9J411‘AVINVEI9VVJI31{LNI 1V0llU3I0J10aI0d3H1dOITLLN3C)1RL10S1VNK0033HJ,:TIN’UK SGMflOHAVRIV3H1NllUI/*ATILLSSIH3IHMITLLN33 31-ILIA10H3DMVISKI1S3H1IVd311105WSIllSflHI‘HdS NI313TID3DN3Ia1aI1S30aVl3H1dOSflICIVI3RLSII+XVN:yi HCTSAV’RIV3111NI S313ll333N3lad3I3111dOIaLM333111dOSLVNIl003:TTT)TJJ 1V0ll1L3l0dVSVM1J’O=daINllGNV 3I0dVSVM.IIN3H1‘cc=d3INdlU3IflSV3NM33HSVH 1V0lHL3I0dI03I0JVI3RL3HMSLV3KIMI3fflAfl43fTJ (co+ snicivi1VfLL3V:3ZIS)UlflSV3N9N13H1V0IHII03l0d 3111MIHIIMA13ULN3SlIdHDIHAkWT3llD3111dO3ZIS3111 (co-i)SflIUVINIV1I33VdO3’1DllDVNIIUIM3H(flflOHS S13XIJANVNA&OHS([I0D3IHDIIIMAVIVN0SNVcWI0DV:(j)j (cool’cool)HaSIVSIIaLN3D3111co=(IVI3DNISI 0113551‘13X1J1SowIaLM3D3HV31NOIIS3flONI313II3 RLdOSfllCIV13H1SIC[VI3l3HhX(co+uva)(J3i3HWflN ONI3H3’I3ll3HDV3NI513X1JdOM0l3GISIflO3111HJJh. S3’IDII3DIULM3DNODdO13SV9NJNIVIMODAVIIVNV IIONIc[MflOIflS501105dOSM0I331RLHJJM AVR1VSITUNI(133VicISI3OVIAII1VDOdNOD3111:(O8L’oEc)9jAJ 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 C*************************** C C THE CENTER COORDINATES OF THE REFERENCE CIRCLE (NK1,NL1), THE C MAXIMUM DIMENSIONS OF THE IMAGE ARRAY (MXX,MXY), AND THE DISTANCE C FROM THE CENTER OF THE REFERENCE CIRCLE TO THE EDGE OF THE ARRAY C SPH(MAX) ARE INiTIALIZED. C NK1 = 100 NL1 = 100 RFIT = 1 MXFIT= 1 MAX =99 MXX = 530 MXY=780 C C THE ARRAY NSPHQ IS INITIALIZED C DO 10 1=1,101 NSPH(I) = 0 10 CONTINUE C C THE AREA IN A BOX EXTENDING UP TO 20 PIXELS TO EiTHER SIDE OF THE C CENTER IF NOW CHECKED AND PIXELS WHICH ARE PORE OR PORE THROAT C PIXELS RESULT IN THE APPROPRIATE ARRAY (NSPH) BEING INCREMENTED C BY ONE. THE ARRAY POSITION WHICH IS INCREMENTED IS THAT WHICH C CORRESPONDS TO THE RADIIJS+1.5 OF THE LARGES CIRCLE CONTAINING C THE PIXEL WHICH WAS CHECKED. C FOR EACH LOOP IN WHICH L IS INCREMENTED, THE PORTION OF THE BOX C TO BE CHECKED IS ONE PIXEL FARTHER FROM THE CENTER THAN IN THE C PREVIOUS LOOP. WITHIN EACH LOOP, NSPH(L+2), THE MAXIMUM SIZE OF C CIRCLE WHICH WOULD FIT WITHIN THE AREA CHECKED IN THAT LOOP, IS C COMPARED TO THE REFERENCE ARRAY CHKSPH(L+2) AND, IF THEY ARE C NOT EQUAL, THE LOOP IS TERMINATED. C DO2OL= 1,MAX N2 = L*2 NSZ=L+2 DO 30 NI = NI1-L,M1-f-L,N2 IF ((NI.LT.1).OR.(NLGT.MXX)) GOTO 30 DK = NI-Mi NK=NK1 +DK DO 30 NJ = NJ1-L,NJ1+L IF ((NJ.LT.1).OR.(NJ.GT.MXY)) GOTO 30 DL=NJ-NJ1 NL=NL1 +DL NIN = SPH(NK,NL) + 1 IF (IMG(NI,NJ).EQ.200) THEN IMG(M,NJ) = 205 ELSE IF (IMG(NI,NJ).EQ.210) THEN IMG(NI,NJ) = 220 ENDIF 1 62 ENDIF IF (IMG(NI,NJ).GT.0) THEN NSPH(NIN) = NSPH(NIN) + 1 ENDIF 30 CONTINUE DO 40 NI = NI1-L+1,NI1+L-1 IF ((NI.LT.1).OR.(NI.GT.MXX)) GOTO 40 DK = NI-Nil NK=NK1 +DK DO 40 NJ = NJ1-L,NJ1+L,N2 IF ((NJ.LT.1).OR.(NJ.GT.MXY)) GOTO 40 DL=NJ-NJ1 NL=NL1 +DL IF (IMG(NI,NJ).EQ.200) THEN IMG(NI,NJ) = 205 ELSE IF (IMG(NI,NJ).EQ.210) THEN IMG(NI,NJ) = 220 ENDIF ENDIF N1N = SPH(NK,NL) + 1 IF (IMG(NI,NJ).GT.0) THEN NSPH(NIN) = NSPH(NIN) + 1 ENDIF 40 CONTINUE IF (NSPH(NSZ).LT.CHKSPH(NSZ)) GOTO 11 20 CONTINUE C C THIS IF BLOCK ENSURES THAT ALL POSSIBLE PORE SIZE VALUES C ARE CHECKED. C 11 IF (L.GT.99) THEN L= 101 ELSE L=L+2 ENDIF C C THIS LOOP DETERMINES THE MAXIMUM SIZE OF CiRCLE WHICH IS C COMPLETELY WITHIN THE PORE OR PORE THROAT. C DO 50 I=2,L IF (NSPH(I).EQ.CHKSPH(I)) THEN RFIT = I-I ENDIF 50 CONTINUE C C THE VALUE OF NREF IS DETERMINED: IF THE CENTER USED IN THIS C SUBROUTINE WAS THE CENTER OF A PORE, NREF IS SET TO 255, AN]) C IFITWASTHECENTEROFAPORETHROAT,ITISSETTO2O. C DO 60 I NI1-1,NIl+1 DO 60 J=NJ1-l,NJ1+l IF (IMG(I,J).EQ.255) THEN NREF = 255 1 63 ELSE IF (IMG(I,J).EQ.20) NREF =20 ENDIF 60 CONTINUE C C THIS FINAL CHECK ENSURES THAT ANY CIRCLES WHICH TERMINATE AGAINST C THE EDGES OF THE IMAGE ARE NOT INCLUDED IN THE PORE SIZE C DISTRIBUTION, AS THESE PORES ARE NOT COMPLETELY INCLUDED WITHiN C THE IMAGE. C IF ((NI1-RFIT).LT.4) THEN NREF =0 ELSE IF ((NI1+RFIT).GT.515) THEN NREF=0 ELSE IF ((NJ1-RFIT).LT.4) THEN NREF=0 ELSE IF ((NJ1+RFIT).GT.771) NREF =0 ENDIF ENDIF ENDIF RETURN END 164 APPENDIX 3 Berea 100 Confocal Pore Size Distributions 1 65 Berea 100 Pore Size Distributions 0.12 0.08 0.04 0.00 106 0.12 0.08 0.04 0.00 Image ABR1 ii’ Image BBR1 II 0 Q E 0 0 0 10 Pore Radius (m) A io5 Pore Radius (m) 1 66 Berea 100 Pore Size Distributions 0.10 0.08 0.06 0.04 0.02 0.00 0.08 106 0.06 0.04 0.02 0.00 - Image CBR1 Image DBR1 10 11 0 0 0 0 C.) E 0 10 Pore Radius (m) —•1 I io5 Pore Radius (m) 1 67 0 I. C) 0 Image FBR1 L 10 Berea 100 Pore Size Distributions 0.10 0.08 0.06 0 I. C.) I 0.04 0.02 10-s Pore Radius (m) Image GBR1 0.15 0.10- 0.05 0.00 -I 106 LL i04 Pore Radius (m) 1 68 Berea 100 Pore Size Distributions 0.08 0.06 0.04 0.02 0.00 10.6 Image HBR1 Image JBR1 I I’ C C) I C •— C) Pore Radius (m) 0.10 0.08 0.06 0.04 0.02 0.00 106 Pore Radius (m) 1 69 Berea 100 Pore Size Distributions 0.08 0.06 0.04 0.02 0.00 -1 106 0.12 - 0.08 Image KBR1 Image LBR1 10 I 0 Q J 0 •— ‘S J Pore Radius (m) 0.04 0.00 10556 I Pore Radius (m) 1 70 0C) E 0 0 .— IS 0 11) 0 Berea 100 Pore Size Distributions Image MBR1 ‘I’ll rh Pore Radius (m) 0.08 0.06 0.04’ 0.02 0.00 106 0.08 — 0.06 0.04 0.02 0.00 106 Image NBR1 Pore Radius (m) 10 171 0 C) a) E 0 C C.) II j Berea 100 Pore Size Distributions Image OBR1 0.10 0.08 0.06 0.04 0.02- 0.00 - — 106 0.12 0.08 0.04 Pore Radius (m) Image PBR1 0.00 106 ‘ii... 10 Pore Radius (m) 1 72 Berea 100 Pore Size Distributions 0.00- 106 0.08 0.06 0.04 0.02 0.00 Image QBR1 Image RBR1 0.12 0.08 0.04 iIlA0.— a) 0 0 .— C) 0 10 Pore Radius (m) 10 I 10 Pore Radius (m) 1 73 APPENDIX 4 S55 Confocal Pore Size Distributions 1 74 Carbonate (S55) Pore Size Distributions Image A78 0.2 — C C.) .) 0.1• 0 0.0 - 106 0.12 0 .— 0 10 Pore Radius (m) 10 Image B78 0.08 0.04’ Pore Radius (m) 1 75 0 .-4 ‘-I C’S E 0 0 C.) C’S J Carbonate (S55) Pore Size Distributions Image C78 0.15 10 Pore Radius (m) 0.10• 0.05 0.00• 0.15 0.10 0.05 0.00• 106 Image D78 Pore Radius (m) 10 1 76 0 . 0 Carbonate (S55) Pore Size Distributions Image E78 0.15 - 0.10- 0.05 - 0.00 4 106 0.3 L Pore Radius (m) Image F78 0 0 0 0.2 0.1 10 Pore Radius (m) 10 177 Carbonate (S55) Pore Size Distributions I C IS C 0.2 0.3 Image G78 Image H78 10 0.1 10 Pore Radius (m) 0.2 0.1• I Pore Radius (m) 1 78 0 .— -‘S 0 ‘-4 E 0 Carbonate (S55) Pore Size Distributions Image J78 0.12 — 0 0.08 0.04 0.00 -. 106 0.2 0.1’ 10 Pore Radius (m) Image K78 0.0• Pore Radius (m) 10’ 1 79 0 C) 1) E 0 I 0.2 Pore Radius (m) L Carbonate (S55) Pore Size Distributions Image L78 0.2 — 0.1 0.0 106 10 Pore Radius (m) Image M78 0.1 1 80 Carbonate (S55) Pore Size Distributions C C I 0.2 Image N78 Image 078 10 0.1- 1 L i05 Pore Radius (m) C C Q C 0.2 0.1’ 0.0 106 10 Pore Radius (m) 181 C I j A Carbonate (S55) Pore Size Distributions Image P78 0.12 C I. C.) j 0.08 0.04 io-5 Pore Radius (m) 10 Image Q78 0.12 0.08 0.04 0.00 106 Pore Radius (m) 1 82 Carbonate (S55) Pore Size Distributions Image R78 0.2 0.1• 0.0• Image S78 10 C 0 I’ 0 I 0.2 1 io4 C 0 0 0) E 0 io5 Pore Radius (m) 0.1• Pore Radius (m) 1 83 0 .— 0 .— C.) I A Carbonate (S55) Pore Size Distributions Image T78 0.10 — 0.08 0.06 0.04 0.02 0.00 106 0.3 0.2 0.1• 0.0 Pore Radius (m) i04 Image U78 106 10 Pore Radius (m) io4 1 84 APPENDIX 5 S77 Confocal Pore Size Distributions 1 85 Carbonate (S77) Pore Size Distributions Image Fl 0.3 - :z 0.0- II1 106 io io io3 Pore Radius (m) Image Gi 0.5 0.4 0.3 J 0.2 106 io4 io3 Pore Radius (m) 1 86 Carbonate (S77) Pore Size Distributions Image Hi I Pore Radius C 0 C.) I C 0 .— C) I 1 • (m) 0.4 0.3 0.2 0.1• 106 0.6 — 0.5 0.4 0.3 0.2 0.1• 0.0 10 Image Ii 106 Pore Radius io (m) 10 1 87 0C.) j 0 .— C.) J Carbonate (S77) Pore Size Distributions Image Ji 0.4 0.3 0.2 0.1 I 0.0 — 106 I 10 Pore Radius (m) Image Ki 0.3 0.2 0.1 0.0 106 10 Pore Radius (m) 1 88 0 .— 0 0 I Image Li 106 Carbonate (S77) Pore Size Distributions 0.2 0.1 0.3 0.2 0.1• 0.0• 106 io-4 io-3 Pore Radius (m) Image Mi IA i05 io4 Pore Radius (m) 10 1 89 0 .,- I. io-5 Pore Radius (m) Carbonate (S77) Pore Size Distributions Image Ni 0.3 - ________ 0.2 - 0.1 0.0 —‘ 106 0.3 0.2 0.1• o•o• r. 106 Pore Radius (m) Image Oi 0 I 11 10 1 90 Carbonate (S77) Pore Size Distributions C () E C Image P1 io5 1o 1o Pore Radius (m) 0 I. I 0.2 0.1• 106 0.2 0.1• 0.0• r 106 Image Q1 10 Pore Radius (m) 191 Carbonate (S77) Pore Size Distributions 0 .— C) 0 Image Ri 0.3 0.2 0.1 0.0• • rrr 106 io i- Pore Radius (m) i03 1 92 APPENDIX 6 Carbonate Samples NMR T1 and Pore Size Distributions 193 SAMPLE S31 Continuous Relaxation Spectrum 0.08 0.06 { 0.02 0.00 - 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 05 0.4 < 0.3 0.2 - 0.1 0.0 - 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 1 94 SAMPLE S31 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 •.. 108 10.6 io Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 C -4 t0.2 A 106 10_s io Pore Radius (m) 1 95 SAMPLE S32 Continuous Relaxation Spectrum 0.08 - 0.06 — 0.04 0.02 0.00 I .....• 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 05 0.4 0.3 0.2 ____,__I_,,FTTDTflI ..........,. I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 1 96 SAMPLE S32 Pore Size Spectra Fast Diffusion Method 0.08 o 0.06 0.04 C 0.02 0.00 1O8 106 io4 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.2 - C I.Q a) 0.1 C 0.0 I i0 106 io5 io Pore Radius (m) 1 97 SAMPLE S33 Continuous Relaxation Spectrum 0.08 EA 0.00 .. I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum - 05 j < 0.3 a) 0.2 .... . 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 1 98 0 .— I. C.) C 0 El SAMPLE S33 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 iO8 iO Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0*.J I i07 ;(.6 Pore Radius (m) 10 1 99 SAMPLE S34 Continuous Relaxation Spectrum 0.08- 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 a) 05 . 0.4W 0.3 0.2 0.1- 0.0- —‘-—= 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 200 0 .— 0 0 0 0 I SAMPLE S34 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 0..8 iO Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0*j\. 10 106 Pore Radius (m) 10 201 SAMPLE S35 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 o.00 I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 a) 05 0.4 0.3w a) 0.2W .= .. ... 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 202 0 C) 0 SAMPLE S35 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 i0..8 10.6 ° Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0 C) J 0.4 0.3 0.2 0.1 0.0 *..I io7 106 Pore Radius (m) 203 SAMPLE S36 Continuous Relaxation Spectrum 0.08 0.06 0.00 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 0• 0.4 0.3 0.2 0.1 0.0 ==c—,—,—,-’-,- 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 204 0 .— c) C 0I. c) 1) E C SAMPLE S36 Pore Size Spectra Fast Diffusion Method 0.08 - 0.06 0.04 0.02 0.00 i08 O Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0*.. 106 io4 Pore Radius (m) 205 SAMPLE S37 Continuous Relaxation Spectrum 0.08 0.06 0.00• 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 ti) 05 0.4 03 0.2 .= 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 206 SAMPLE S37 Pore Size Spectra Fast Diffusion Method 0.08 o 0.06 O.04 0 0.02 0.00- 108 io 106 io5 io Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0 0.2 0.1- . A i07 106 Pore Radius (m) 207 SAMPLE S38 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.00 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) a) a) .— a) 0.6 0.5 0.4 0.3 0.2 0.1 0.0 GWQJ.. oooa 1 10 100 1000 208 SAMPLE S38 Pore Size Spectra Fast Diffusion Method 0.08 i08 io7 i06 io io4 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0 0.2 I 0.1• . A..1 i07 106 io5 io4 Pore Radius (m) 209 SAMPLE S39 Continuous Relaxation Spectrum 0.08 I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 I) 0 0.4 0.3 a) 0.2 0.1• 0.0• DC 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 210 SAMPLE S39 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 0 .— cj 0 0 I. 1)$ 0 108 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0*../ Pore Radius 10 (m) 10 211 SAMPLE S40 Continuous Relaxation Spectrum 0.08 0.06 - 0.04 0.02 0.00 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) I a) a) .. a) 0.8 0.6 0.4 0.2 0.0 .DD 1 10 100 1000 212 0.08 SAMPLE S40 Pore Size Spectra Fast Diffusion Method 0.06 .— C.) c 0.04 , C > 0.02 0.00 10.8 O 10.6 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0 C.) 11 0.3 0.2 0.1 0.0 __ 106 Pore Radius (m) 213 SAMPLE S41 Continuous Relaxation Spectrum O.08 0.06 0.04 0.02 0.00 I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 a 0.4 < 0.3 p. 0.2 0.1• 0.0 ;I.fl.I.•.I.fl.1fl. i:i:i;ixR—————r——r——i—t—i—n 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 214 0.08 SAMPLE S41 Pore Size Spectra Fast Diffusion Method 0.06 0.04 0.02 0.00 0 I. c) 0 0 I’ 0 108 io7 106 io io Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0 io Pore Radius 10 (m) 10 215 SAMPLE S42 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.00• I I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 a) a 0.4 < 0.3 a) 0.2 0.1• 0.0 .IiI000aaDo 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 216 SAMPLE S42 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 C I. C) 1 C C I 0 I 108 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0_/.A iø io Pore Radius (m) 217 SAMPLE S43 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.00- 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 0.4 • 0.3 JL 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 218 SAMPLE S43 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 C C I. c) C C 0 Q C i08 iO Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0 10 A • 1116 Pore Radius (m) 10 219 0.06 - 0.04 0.02 - 0.00 - Spin-Lattice Relaxation Times (ms) SAMPLE S44 Continuous Relaxation Spectrum 0.08 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum I 1) 0.6 0.5 0.4 0.3 0.2 0.1 0.0 Looaa j...Lq.. 1 10 100 1000 220 00 SAMPLE S44 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 - 0.00 10.8 ° Pore Radius (m) Pore Size Spectrum Slow Diffusion Method C 0 I 0.3 0.2 0.1 0.0 *..I 10 1O6 Pore Radius (m) i04 221 SAMPLE S45 Continuous Relaxation Spectrum 0.08 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.5 - 0.4 - . 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 222 0 C.) C 0 C.) z C > SAMPLE S45 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 iO8 0 0 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0*.i . io6 io Pore Radius (m) 223 SAMPLE S46 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.0o I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 0.5 0.4 1 10 100 1000 Spin-Lattice Relaxation Times (ms) 224 SAMPLE S46 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 C .— C C I i0.8 1O6 0 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0 1A! i- Pore Radius (m) 225 SAMPLE S47 Continuous Relaxation Spectrum 0.08 0.06 0.04 - 0.02 - 0.00 - 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) E — ci) 0.8 0.6 0.4 0.2 0.0L 1 10 100 1000 226 SAMPLE S47 Pore Size Spectra Fast Diffusion Method 0.08 o 0.06 0.04 0 0.02 0.00 I 108 iO-7 10 io5 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0 : :: 0.0• A io-7 io-4 Pore Radius (m) 227 SAMPLE S48 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 - 0.00 - ——I 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) I 0.8 0.6 0.4 0.2 0.0 JDODDD 1 10 100 1000 228 0 C) 0 0 . r-) I SAMPLE S48 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 10.8 O. ° Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0 *..\. i07 io Pore Radius (m) 10 229 SAMPLE S49 Continuous Relaxation Spectrum 0.08 0.06 - 0.04 0.02 0.00 1 10 100 1000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum Spin-Lattice Relaxation Times (ms) a 0.8 0.6 0.4 0.2 0.0*..E. 1 10 100 1000 230 0 .— C.) 0 0 C.) j SAMPLE S49 Pore Size Spectra Fast Diffusion Method 0.08 0.06 0.04 0.02 0.00 0.3 0.2 0.1 0.0 1- 108 10 106 io5 io4 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 10.6 Pore Radius (m) i04 231 SAMPLE S54 Continuous Relaxation Spectrum 0.07 0.00 I 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum I) 05 0.4W 0.3W ...... 1 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) 232 SAMPLE S54 Pore Size Spectra Fast Diffusion Method 0.07 0.06 •, 0.05 C) 0.04 000 E 0.03 0.02 0.01- i07 106 IIY Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.8 0.6 C) 0.4 C 0.2 0.0 106 io5 io4 Pore Radius (m) 233 SAMPLE 855 Continuous Relaxation Spectrum 0.07 - 0.00 I I 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.6 0 -e 05- I ____ _____ 0.0 — .....__r_..r_r_r.......rITrTrrlIdIrrrrrrr •.I - 1 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) 234 SAMPLE S55 Pore Size Spectra Fast Diffusion Method 0.06 0.05 0.04 0.03 0.02 0.01 0.00 0 C.) 0 0 C.) C.) z 0 > 102 1111 101 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.3 0.2 0.1 0.0 106 Pore Radius (m) 10 235 SAMPLE S77 Continuous Relaxation Spectrum 0.08 0.06 0.04 0.02 0.00 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) Discrete Relaxation Spectrum 0.8 0.6 0.4 1 10 100 1000 10000 Spin-Lattice Relaxation Times (ms) 236 SAMPLE S77 Pore Size Spectra Fast Diffusion Method O.08 0.06 g 0.04 C 0.02 0.00 i0 10.6 Pore Radius (m) Pore Size Spectrum Slow Diffusion Method 0.6 0.5 0 .— 0.4 ) 0.3 I 0.2 0.1• I Pore Radius (m) 237 APPENDIX 7 NMR Ti distributions and M(O) vs. Sw graphs for Partially Saturated carbonate sampies 238 Sa m pl e S3 1 N or m al iz ed T i Sp ec tr a a t Pa rt ia l Sa tu ra tio ns 0. 08 Sa tu ra tio n 0. 98 4 0. 85 7 0. 77 4 0. 71 9 0. 65 3 0.0 6 0. 04 0. 02 0. 00 — 0 - - 1 10 10 0 Sp in -L at tic e Re la xa tio n Ti m e (m s) 10 00 F’. ) CA ) CD 22 00 20 00 18 00 M (0) 16 00 14 00 12 00 Sa m pl e S3 1 C or re la tio n o f Sa tu ra tio n w ith M (O ) 0. 6 0. 7 0. 8 0. 9 Sa tu ra tio n 1. 0 F’ ) I Sa m pl e S4 1 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns 0. 08 0. 06 0. 04 0. 02 0. 00 Sa tu ra tio n - - 0. 98 3 0. 92 0 0. 83 2 0. 79 6 0. 67 3 0. 56 0 1 10 10 0 Sp in -L at tic e R el ax at io n Ti m e (m s) 10 00 1) 13 00 12 00 11 00 10 00 90 0 80 0 70 0 60 0 Sa m pl e S4 1 C or re la tio n o f Sa tu ra tio n w ith M (O ) M (0 ) 0. 5 0. 6 0. 7 0. 8 0. 9 Sa tu ra tio n 1. 0 F\ ) I\) Sa m pl e S4 2 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns 0. 08 Sa tu ra tio n 0. 98 9 0. 87 2 0. 71 1 0. 61 1 0. 51 2 0. 06 0. 04 0. 02 0. 00 — 0 - - - 1 10 10 0 Sp in -L at tic e R el ax at io n T im e (m s) 10 00 I’) c) Sa m pl e S4 2 C or re la tio n o f Sa tu ra tio n w ith M (O ) M (O ) 12 00 80 0 40 0 0 .4 0. 6 0. 8 Sa tu ra tio n 1. 0 r’3 - a I Sa m pl e S4 4 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns 0. 06 0. 05 0. 04 0. 03 0. 02 0.0 1 0. 00 - 0 - - - - 0 - Sa tu ra tio n 0. 91 4 0. 81 2 0. 71 1 0. 61 4 0. 50 9 0. 32 4 1 10 10 0 Sp in -L at tic e R el ax at io n T im e (m s) 10 00 C, ’ 0) M (O ) 16 00 12 00 80 0 40 0 Sa m pl e S4 4 C or re la tio n o f Sa tu ra tio n w ith M (O ) 0. 3 Sa tu ra tio n 0. 5 0. 7 0. 9 Sa m pl e S5 4 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns Sa tu ra tio ns 0. 99 5 0. 92 0 — G — 0. 88 6 — . — — — 0. 79 4 0. 73 0 — 0 -— — 0. 63 2 0. 55 7 0. 46 4 0. 36 2 - I- 0. 23 7 ro —1 Sp in -L at tic e Re la xa tio n Ti m e (m s) ii) 0. 08 0. 06 0. 04 0. 02 0. 00 10 10 0 10 00 10 00 0 I.’ ) - 70 00 60 00 50 00 M (0) 40 00 30 00 20 00 10 00 Sa tu ra tio n Sa m pl e S5 4 C or re la tio n o f Sa tu ra tio n w ith M (O ) 0. 2 0. 4 0. 6 0. 8 1. 0 F’ ) - (.0 Sp in -L at tic e R el ax at io n Ti m e (m s) — 0 - - . 1. — - - I Sa m pl e S5 5 N or m al iz ed T i Sp ec tra a t Pa rt ia l Sa tu ra tio ns 0. 07 0.0 6 0. 05 0. 04 0. 03 0. 02 0.0 1 0. 00 Sa tu ra tio n 0. 98 4 0. 87 0 0. 77 8 0. 69 7 0. 46 5 0. 40 5 0. 30 8 0. 18 4 0. 03 5 1 10 10 0 10 00 10 00 0 Sa m pl e S5 5 C or re la tio n o f Sa tu ra tio n w ith M (O ) 60 00 50 00 40 00 M (O ) 30 00 20 00 10 00 0 0. 0 0. 2 0. 4 0. 6 0. 8 Sa tu ra tio n 1. 0 I\ ) 01 0

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
United States 4 1
China 3 27
Japan 2 0
Russia 1 0
Poland 1 0
City Views Downloads
Ashburn 3 0
Beijing 3 0
Tokyo 2 0
Unknown 2 0
Redmond 1 1

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}

Share

Share to:

Comment

Related Items