Quantifying Porosity through Automated Image Collection and Batch Image Processing: Case Study of Three Carbonates and an Aragonite Cemented Sandstone

: Modern scanning electron microscopes often include software that allows for the possibility of obtaining large format high-resolution image montages over areas of several square centimeters. Such montages are typically automatically acquired and stitched, comprising many thousand individual tiled images. Images, collected over a regular grid pattern, are a rich source of information on factors such as variability in porosity and distribution of mineral phases, but can be hard to visually interpret. Additional quantitative data can be accessed through the application of image analysis. We use backscattered electron (BSE) images, collected from polished thin sections of two limestone samples from the Cretaceous of Brazil, a Carboniferous limestone from Scotland, and a carbonate cemented sandstone from Northern Ireland, with up to 25,000 tiles per image, collecting numerical quantitative data on the distribution of porosity. Images were automatically collected using the FEI software Maps, batch processed by image analysis (through ImageJ), with results plotted on 2D contour plots with MATLAB. These plots numerically and visually clearly express the collected porosity data in an easily accessible form, and have application for the display of other data such as pore size, shape, grain size/shape


Introduction
Automated image acquisition and image analysis through scanning electron microscopy (SEM) has a long history, with early pioneers [1] presenting an example of automated image acquisition, image analysis and a grid-based graphical representation of clay particle orientation collected through SEM analysis.In more recent years, with improvements in computer power and storage capabilities, it has become possible to use SEM to acquire high-resolution images over increasingly larger areas [2][3][4][5].Pore image analysis from SEM images for a variety of porous media have been performed by a number of authors [6][7][8], and are routinely used in the analysis of petrographic thin sections examined by optical microscopy [9,10].Image analysis has an important role in supporting microscopic observations, being used to improve pore characterization of sedimentary rocks (pore shape, size, orientation, distribution, etc.).Image analysis provides reproducible and representative numerical data for the petrographic study of carbonate and other geological materials in thin sections [11,12].
Porous media (natural and manmade) and particularly geoscience-relevant porous materials are ubiquitous, with relevant examples being sandstone, carbonate, mudrock, concrete, cement and mortar.Knowledge of the porosity, permeability and flow properties of such materials can be critical in understanding their structural stability, effectiveness as aquifers and barriers to fluids.Any method that can potentially provide quantitative data on porosity at the micrometer to millimeter scale is therefore welcome and can provide useful information suitable for modeling flow through such media.
For a review of many of the methods available for the study of porosity and pore structure, see [13].The latter authors recognized three broad routes to obtaining data on porosity, namely: (i) direct measurement; (ii) imaging (using a variety of techniques) and; (iii) down hole porosity logs (wireline logging techniques).Direct measurement and down hole porosity logs typically provide information on porosity, with no additional information on the morphology or texture of pores and pore networks, which can be extracted from imaging techniques.Of the imaging techniques, X-ray tomography (XRT) and scanning electron microscopy (SEM) are perhaps arguably the most significant; with XRT capable of producing micron-scale realistic 3D porosity models of samples, and SEM 2D porosity reconstructions.SEM images, from polished thin sections, can be used with software such as Pore Architecture Model (PAM) and Pore Analysis Tool (PAT), to produce porosity/permeability model data from stochastically produced reconstructions [14].Other potential imaging techniques include 3 H and 14 C PMMA autoradiography [15,16] and confocal microscopy [17].
The new technique presented herein typically utilizes many thousand images collected across a fixed grid network, which is a relatively cost-and time-effective technique.This takes advantage of modern SEM's capability for automated image collection of high-resolution images over largescale areas, and digital storage, combined with readily available freeware for batch image analysis.The production of reproducible porosity data that are realistic, maintain original spatial relationships and are representative of the real world are therefore easily achievable.Although no data on factors such as pore shape/size and connectivity are presented here, such could be obtained relatively simply through further manipulation using image analysis.The technique provides information on vertical and lateral changes in porosity, which are useful for the modeling of fluid flow through porous media, as well as helping to explain observed flow patterns, potential sources of blockage and porosity pathways.

Materials and Methods
Samples examined include a laminated carbonate (laminite) and a coquina from onshore analogues of the Pre-salt Cretaceous of Brazil [18][19][20].With supplementary material including a Carboniferous limestone from Spireslack, Ayrshire, Scotland [21], and a Quaternary aragonite cemented sandstone from County Antrim, Northern Ireland.Polished thin sections were prepared and imaged within an FEI Quanta FEG 650 SEM (FEI, Hillsboro, USA), and imaged with a quad/or CBS (concentric) backscattered electron (BSE) detector in low vacuum mode at 0.83 Torr, with an operating voltage of 20 kV and spot size of 4.5.Whole slides were imaged using the FEI Maps software (version 2.1, 64 bit, FEI, Hillsboro, OR, USA), following a simple workflow (Figure 1), collecting hundreds to thousands of images (tiles) with a horizontal field width (HFW) of 50 to 2000 microns, over a regular grid (Table 1).Map collection times varied from approximately 20 min to 35 h (Table 1).Examples illustrated herein have pixel resolutions ranging from 6.51 nm to 2.6 µm (Table 1).Higher pixel resolutions can be achieved through smaller tile sizes and/or use of higher tile pixel dimensions, therefore the technique is capable of recording microporosity and macroporosity.Individual tiles were batch processed with ImageJ (version 1.47 for Windows, 64 bit, free software, National Institutes of Health, Bethesda, MD, USA) following the protocol illustrated in Figure 2; being thresholded for porosity, binarized (Figure 3) and batch processed returning a percentage porosity value for each tile.Thresholding was carried out visually using the default setting and manually adjusting to incorporate the maximum degree of observed porosity.This was deemed sufficient for the current investigation, as although slight differences in the placement of thresholding would affect the actual percentages measured, herein the absolute value is of less importance than the variability between adjacent tiles.The latter is not affected by this method of thresholding, as the selected threshold is automatically applied to all tiles during batch processing.Percentage porosity values for each tile were saved as a single data column, then copied and pasted as a column in a Microsoft Excel workbook (version 12.0.0.0,Microsoft, Redmond, WA, USA), starting at cell position A1.In position B1 place the formula: where # is replaced with the number of tiles within each row of the collected Maps montage.The formula is then copied from B1 and pasted (starting at B1) in an area equivalent to the number of columns and rows recorded in Maps.Data is now organized in an array of the same dimensions as the grid of image tiles recorded within Maps, although collected rows are transposed into columns within the spreadsheet (Figure 4a,b).The array is then saved in the Excel file and transferred to MathWorks MATLAB R2014b (version 8.4.0.150421, 64 bit, MathWorks Inc., Natick, MA, USA).In MATLAB the plotting function "Contourf" is used to produce colored contoured representations of the variation in recorded porosity.Plotted data require reorientation and resizing to reflect the original disposition of the montaged BSE parent image (Figure 4c).
In addition, some data were also collected on the distribution of silicates within the coquina thin section, siderite within the Spireslack sample and mean grayscale across all four of the test samples.These were processed in the same fashion, with the exception that mean grayscale did not require binarizing, but was instead batch processed directly from the original BSE tiles.Mineral phases were identified through energy dispersive X-ray analysis (EDX), using Oxford Instruments AZtec software (version 3.3, 64 bit, Oxford Instruments, Abingdon, UK), and the characteristic gray level of the particular mineral phase.
As well as contour plotting of the data collected, some additional parameters (mineral phase coverage) were scatter plotted against porosity, using standard Excel graphing functions.

Results
The application of the techniques outlined above were successfully applied, with the construction of color contoured maps for the distribution of porosity, as well as a number of other parameters.Utilization of percentage data collected from individually scanned tiles was also used to data mine additional information and in some cases plot relationships between porosity and other parameters.

Porosity Percentage
Porosity percentage maps were successfully calculated and plotted for all four samples (Figures 5-7), which all showed different distributions in the occurrence of pores.In the case of the laminite, higher levels of pores were found to be restricted to the more thinly laminated part of the thin section (Figure 5b,c).Porosity mainly occurs between 0-30% with an average of 18% and mode of 10%

Results
The application of the techniques outlined above were successfully applied, with the construction of color contoured maps for the distribution of porosity, as well as a number of other parameters.Utilization of percentage data collected from individually scanned tiles was also used to data mine additional information and in some cases plot relationships between porosity and other parameters.

Porosity Percentage
Porosity percentage maps were successfully calculated and plotted for all four samples (Figures 5-7), which all showed different distributions in the occurrence of pores.In the case of the laminite, higher levels of pores were found to be restricted to the more thinly laminated part of the thin section (Figure 5b,c).Porosity mainly occurs between 0-30% with an average of 18% and mode of 10% (Table 2).Smaller, but significant areas with up to 83% porosity occur (Tables 2 and 3).These are associated with thin laminae in the order of 500 µm thick and appear to be laterally discontinuous at the sub-millimeter scale.
The percentage porosity map for the coquina (Figure 6b) exhibits a low bulk porosity (Tables 2  and 3), with isolated apparently disconnected areas with relatively higher porosity (10-20%), although an area that appears more connected (higher percentage of pores) occurs towards the top of the thin section the majority of the sample has a porosity of less than 10% (Table 3).
A porosity map for the Spireslack limestone (Figure 6d) is unusual, as the distribution of pores appears more heterogeneous than homogeneous, with a degree of apparent vertical to steeply inclined fabric that occurs nearly perpendicular to sedimentary bedding.The average and modal porosity values are both low (Table 2), with over 98% of the area having less than 10% porosity.
The aragonite cemented sandstone is particularly interesting, as percentage porosity maps display at least two forms of porosity.Large voids where tiles comprise 100% porosity, and zones or bands with 20-50% porosity (or locally higher), occurring within bivalve shells or as bands parallel to bedding (Figure 7b-d).Average porosity for this sample is under 25% (Table 2), representing the typical background matrix porosity, while the significant occurrence of larger vugs is indicated by the increased frequency of porosity values in the range of 91-100% and higher maximum porosity values (Tables 2 and 3).These reflect different degrees of cementation with an acicular aragonite cement phase (Figure 3b).Porosity maps for the aragonite cemented sandstone were taken at three different tile resolutions, 2 mm, 1 mm and 518 µm (Figure 7b-d).In all three cases, variations in porosity were clearly picked out and show the same frequency trends (Figure 8a).

Mineral Distribution
In the Spireslack limestone sample, the contouring technique clearly displays variable concentrations in the occurrence of diagenetic siderite crystals (Figure 9b).The latter are not apparent from the large field BSE map (Figure 6c), or from inspection of the BSE images at higher magnifications (Figure 9c).Casual inspection of the contoured siderite percentage map (Figure 9b) suggests that on the right-hand side of the slide there is a relatively siderite depleted zone which roughly correlates with areas of higher porosity (Figures 6c and 9b).This is, however, not the case across the whole of the slide, where patches of high siderite and porosity concentration appear to coincide.A scatter plot of percentage porosity against that of siderite abundance per tile (Figure 10a) exhibits a well-developed clustering of points with an apparent positive correlation, with good correlation significance (Spearmans rank = 0.466, n = 25,080, p < 0.001).It was also noted that the technique could be applied to the coquina, with the distribution of silicate rich phases (quartz and clay) in reference to that of porosity.Comparison of the two appear to exhibit a level of positive correlation (Figures 6b and 9a), suggesting that the two may well be linked.A scatter plot of percentage porosity against percentage silicate coverage for individual tiles also shows a positive significant relationship (Figure 10b).

Mean Grayscale
The use of mean grayscale plots also in many cases illustrates information concerning the structure of the thin sections, amplifying and supporting results observed from porosity maps and

Mineral Distribution
In the Spireslack limestone sample, the contouring technique clearly displays variable concentrations in the occurrence of diagenetic siderite crystals (Figure 9b).The latter are not apparent from the large field BSE map (Figure 6c), or from inspection of the BSE images at higher magnifications (Figure 9c).Casual inspection of the contoured siderite percentage map (Figure 9b) suggests that on the right-hand side of the slide there is a relatively siderite depleted zone which roughly correlates with areas of higher porosity (Figures 6c and 9b).This is, however, not the case across the whole of the slide, where patches of high siderite and porosity concentration appear to coincide.A scatter plot of percentage porosity against that of siderite abundance per tile (Figure 10a) exhibits a well-developed clustering of points with an apparent positive correlation, with good correlation significance (Spearmans rank = 0.466, n = 25,080, p < 0.001).It was also noted that the technique could be applied to the coquina, with the distribution of silicate rich phases (quartz and clay) in reference to that of porosity.Comparison of the two appear to exhibit a level of positive correlation (Figures 6b and 9a), suggesting that the two may well be linked.A scatter plot of percentage porosity against percentage silicate coverage for individual tiles also shows a positive significant relationship (Figure 10b).

Mean Grayscale
The use of mean grayscale plots also in many cases illustrates information concerning the structure of the thin sections, amplifying and supporting results observed from porosity maps and distribution maps for elemental and mineral phases.Examples include the location of porous laminae (Figure 11a), the location of shelly versus silica rich layers (Figure 11b), and vertical zones that cross-cut bedding (Figure 11c).The latter in particular shows a strong difference in grayscale, which ties in with areas that are relatively depleted in siderite (Figures 9b and 11c).The grayscale map of the aragonite cemented sandstone from Co. Antrim (Figure 11d) is particularly interesting as it provides a similar degree of information on the distribution of porosity, compared to the corresponding map for percentage porosity based on the same tile size (Figure 7c).In addition, areas interpreted as more highly aragonite cemented, or areas of bivalve shell, are picked out in bright green (Figure 11d).In (A) blue areas correspond to high porosity regions within the laminite, while in (B) blue areas relate to localities with either high porosity and/or high silica contents.The green and orange areas are carbonate-rich.In (C) darker/lower values for average grayscale indicate a vertically orientated feature (right-hand side), which appears to correspond to an area that is relatively depleted in siderite.The darker mean grayscale in (D) is aligned with areas that are high in porosity, while isolated patches of brighter green match areas that are particularly well cemented with aragonite (AC), or represent areas of bivalve shell (S).Direction of scanning as in associated porosity maps.In (A) blue areas correspond to high porosity regions within the laminite, while in (B) blue areas relate to localities with either high porosity and/or high silica contents.The green and orange areas are carbonate-rich.In (C) darker/lower values for average grayscale indicate a vertically orientated feature (right-hand side), which appears to correspond to an area that is relatively depleted in siderite.The darker mean grayscale in (D) is aligned with areas that are high in porosity, while isolated patches of brighter green match areas that are particularly well cemented with aragonite (AC), or represent areas of bivalve shell (S).Direction of scanning as in associated porosity maps.

Beam Stability Issues
Montages collected using Maps software can take many tens of hours (Table 1), and it is therefore essential that the electron beam is stable, exhibiting low levels of drift.Instability will result in areas of identical composition producing a variance in BSE intensity, which will affect image analysis, compromising the ability to threshold BSE grayscale images.Field emission SEM is best suited for such tasks, being characterized by highly stable electron sources, with little inherent drift.Conversely, W-hairpin filaments are likely to be less suitable for scans over more than several hours, with a tendency for beam drift, shorter lifetime and requirement for frequent replacement; making SEM's fitted with such sources less suitable for large montages.

Preparation Considerations
Well-polished thin sections are best suited for such studies, as they exhibit high levels of planarity, low levels of surface roughness and therefore require minimal effort in focusing.Polished resin impregnated blocks can also be suitable, although polished surfaces are not necessarily parallel to the sample base, necessitating use of an interpolated focus setting (based on three points of focus), which increases the time required for scanning.Large-area ion beam-milled surfaces can also be imaged using the technique described herein, although these may require more complex autofocus solutions owing to the shape of the milled surface (concave dish due to rotation during milling), considerably increasing the time required to image large surfaces.Some problem areas were noted during the collection of these data.Namely the need for cleanliness, as fluff, hairs and grease can easily occur.Where these are present, they may locally increase the count for porosity.However, in the current cases, although present to a degree, they appear to have had no major substantial effect (see Figure 3e,f).

Potential Problem Areas
The technique does require operator input into thresholding, to select features of interest, which can lead to discrepancies between users.It is hard to quantify to what degree this will affect results, however this does not significantly impact on the usability of the technique, as the main advantage is visualization of variability in parameters such as porosity, which is not greatly affected by minor differences in the positioning of threshold limits.Where precision is more of an issue, various automated thresholding techniques may be applied (e.g., Huang, Li, Mean, Otsu, Yen thresholding, etc.), which may be more applicable.In cases where it is desirable to constrain any degree of error, such as operator, counting, or specimen error, further steps may be taken [22][23][24].In terms of porosity, many tiles return high values (80-100%), which is much higher than porosity traditionally recorded from whole rock analysis.This is a function of the small size of tiles, which are used to calculate porosity from tile to tile, with porosity values referring to micron-to millimeter-sized areas within each recorded section.Averaging of porosity across whole slides (Table 2) returns values that are representative of the whole rock.Figures obtained for percentage porosity, from individual tiles, should not be taken at face value (in terms of traditionally quoted values for whole rock porosity), but can be used to demonstrate how porosity changes across examined thin sections.

Relevance of Technique
Frequency data of percentage pore coverage for all four samples (Figure 8b and Table 3) clearly illustrates the dominance of areas with lower porosity (0-20%), which is also confirmed by average and modal porosity analysis of individual tiles (Table 2).All samples are characterized by a steady (but low) abundance of percentage porosity values in the 21-90% range, within the Antrim sample a secondary increase in frequency around a porosity of 91-100%.The individual profiles do differ in detail, potentially reflecting differences in the density/distribution of pores.Although further image analysis would be required to fully interpret these differences, it is clear that the technique outlined in the current paper represents a useful source of information that can be used to compare samples and observe variation within single samples.It does not, however, imply pore connectivity, which would require further detailed analysis of the data.

Mineral Distribution
Although the examples here are mainly used to illustrate variation in porosity, the technique can be easily applied, in conjunction with spot EDX analysis, to graphically illustrate variations in the distribution for mineral phases (Figure 9a,b).Color contoured maps can be used directly to determine and illustrate variability in mineral phase distributions, or individual tiles can be successfully used to produce scatter plots of variation in mineral phase concentration plotted against other criteria such as porosity (or other measured criteria) to statistically determine relationships between sets of parameters.Such relationships are valuable, given the high statistical significance that can be attached to data sets containing large numbers of measurements (Table 1).Both methods therefore have high value for supplementing information obtained from traditional petrographic analysis.

Mean Grayscale
Maps based on mean grayscale (Figure 11) have the advantage of not requiring thresholding or binarizing.They are therefore marginally quicker to construct and are not affected by operator decision processes, thus being more reproducible.Although some of the grayscale maps appear to show close similarity to maps for porosity (Figure 11a,b,d), information that can be extracted is less targeted, potentially reflecting a combination of multiple parameters including mineral composition, as well as porosity.In the case of the Antrim aragonite cemented sandstone, the mean grayscale map (Figure 11d) exhibits components that match porosity, some that reflect areas with higher levels of cementation or calcite shelly material.Mean grayscale contoured maps cannot be taken as direct analogues for porosity, and may prove harder to take simple interpretations from.Factors that may have to be taken into consideration include the degree of BSE contrast between mineral phases present, and the level of contrast between pores and matrix.

Other Parameters
It is worthwhile noting that any parameter that can be measured through image analysis can be plotted and therefore visualized using this process.Maps can therefore be easily constructed, with little additional effort, for other parameters such as pore size, shape, perimeter and orientation, as well as the corresponding measurements related to particles.

Effect of Individual Tile Size
The importance in the careful selection of the tile size used to acquire images for analysis is clearly illustrated in the case of the aragonite cemented sandstone (Figure 7).Here as tile size was decreased from 2 mm to 1 mm and finally 518 µm (Figure 7b-d) the increase in the quality of the recorded porosity data is clearly seen.Scans taken with an individual tile field of view of 2 mm clearly depict the larger pore systems, although it takes a tile size of 518 µm to clearly visualize and record the occurrence of more subtle variations in porosity inside bivalve shells and within bedding across the sample.This is perhaps not unexpected, and illustrates that time spent collecting higher resolution data is fruitful, but at the expense of both time and cost.Comparison of the minimum, maximum, average and modal percentage porosity values for the three tile size montages of the aragonite cemented sandstone (Table 2), however, exhibits very little variation between the montages.This indicates that general data extracted on percentage porosity coverage is not affected by the size/resolution of individual tiles.For the Antrim sample (Figure 8a), relative increase in frequency, at either end of the percentage porosity spectrum, in comparison to values between 21-90% porosity, between the three trends (2 mm, 1 mm, 518 µm) is to be expected, as smaller tile sizes would naturally favor the recording of areas with low porosity, as well as large cavities with high porosity.

Future Potential
The technique is also equally applicable to other materials, including sandstone and mudrock, as well as, among others, but not limited to, other porous media such as cement, concrete and mortar.Any material that requires characterization at the sub-micron to centimeter scale can potentially be investigated using this technique.The range of pore scales that can be examined is a function of the size of the individual scanned tiles, as well as the pixel dimensions of these tiles, which determines pixel resolution.By varying these parameters, it is possible to examine a range of pore sizes from nanopores, through micropores to macropores, although the size of the pores targeted will determine the time taken and the area that it is possible to image (see Table 1).Useful data can be collected on porosity distribution from relatively large areas (up to several centimeters squared).Such data is valuable for modeling the flow of fluids both across and parallel to bedding and in evaluating the likely connectedness of pore spaces within carbonates and other porous media.This technique is useful as it allows for the collection of quantifiable measurements of porosity at a sub-micron resolution over millimeter to centimeter squared areas, while still preserving relative spatial relationships, shown particularly well in Figure 5c.Porosity maps indicate preferential porosity alignments that might explain location of permeable pathways which in the rocks described may be very localized elements, which could subsequently be targeted for even higher resolution studies.Resolutions achieved over large areas should be as good, if not better, than that achieved through X-ray tomography.3D information on porosity distribution is required to determine the level of actual pore connectivity, and the 2D contour maps cannot directly infer flow pathways, only showing 2D variability in the distribution of pore percentage coverage.In theory, 3D reconstructions of porosity and permeability, using the Pore Architecture Modeling (PAM) and Pore Analysis Tool (PAT) methods [14] based on the collected binary images, can be used to further refine the analysis, through more detailed modeling in areas of variable porosity as identified through the current method.Indeed, with further development and full automation, it would be possible to apply PAM/PAT models to each tile analyzed, leading to the potential of producing many thousand such stochastic poro-perm models per square centimeter.

Conclusions
This technique generates realistic numerical data on porosity at a sub-micron resolution, over areas of a centimeter or so in scale, it can therefore be used to characterize the nature of pore systems from the nanopore to macropore scale.Such data can be used as the basis in modeling flow through natural porous media.As well as porosity, the technique has potential for mapping other features such as pore shape, size, grain shape, size, orientation, variations in mineral phase distributions and the overall distribution and texture of pore systems.Although other methods are available, such as X-ray tomography, and 3D pore modeling from 2D SEM images, the current method is relatively simple to perform and to an extent less costly.Samples can be automatically run overnight, and data manipulation and plotting takes a minimum of time thanks to batch processing and ease of use.

Supplementary Materials:
The following are available online at www.mdpi.com/2076-3263/7/3/70/s1,Excel data on porosity values for tiles, for each of the examples illustrated.

Figure 1 .Figure 1 .
Figure 1.Screenshots illustrating workflow associated with acquisition of backscattered electron (BSE) tile sets used as data source for porosity analysis.Antrim polished slide used as illustration.(A) Optical view of polished slide within chamber, taken with "Color Nav-Cam", for the purpose of locating sample.(B) BSE overview image with whole slide selected (red arrow), by drawing box over area of interest and selecting add tiles; generating a 15 × 34 grid of tiles approximately 2 mm wide (yellow arrow), obtained with a scan rate of 10 µs.Blue numbered arrows indicate direction and order of tile collection.(C) Selected area for higher resolution scan (red arrow); 19 × 47 grid of tiles approximately 1 mm wide (yellow arrow), with a scan rate of 5 µs.Blue arrows as in (B).(D) Selected area for higher resolution scan (red arrow); 39 × 48 grid of tiles approximately 518 µm wide (yellow arrow), with a scan rate of 5 µs.Blue arrows omitted due to position of yellow arrow, however tile collection details as in (C).Tile sets collected in (B)-(D) by selecting area in Maps, with tile parameters read from microscope control system, and pressing acquire within Maps control area.

Figure 2 .Figure 2 .
Figure 2. Flow diagram showing methodology for batch processing of image tiles within ImageJ, to extract numerical porosity data for contouring within MATLAB.Red boxes highlight screenshots of workflow within ImageJ.

Figure 4 .
Figure 4. Example of how grid network organization changes during various stages of analysis.(A) Layout of grid, used to construct BSE montage, taken with Maps software.Images collected in rows (red arrow), starting top left, working towards the right, then moving down one row (sequence 1, 2, 3, 4, 5, etc.).(B) Data transferred into an Excel matrix, arranged in columns (red arrow), starting top left and working downwards, then progressively moving to columns to the right.Intermediate stage, in ImageJ, where data saved as a single column not illustrated.(C) Organization of displayed data in MATLAB contoured map.Note that data is still displayed in columns, but has been inverted in comparison to the Excel data sheet.Final stage involves rotating contoured map 90° clockwise and resizing to the same dimensions as the original BSE montage collected in (A).

Figure 4 .
Figure 4. Example of how grid network organization changes during various stages of analysis.(A) Layout of grid, used to construct BSE montage, taken with Maps software.Images collected in rows (red arrow), starting top left, working towards the right, then moving down one row (sequence 1, 2, 3, 4, 5, etc.).(B) Data transferred into an Excel matrix, arranged in columns (red arrow), starting top left and working downwards, then progressively moving to columns to the right.Intermediate stage, in ImageJ, where data saved as a single column not illustrated.(C) Organization of displayed data in MATLAB contoured map.Note that data is still displayed in columns, but has been inverted in comparison to the Excel data sheet.Final stage involves rotating contoured map 90 • clockwise and resizing to the same dimensions as the original BSE montage collected in (A).

Figure 5 .
Figure 5. (A,B) Backscattered electron (BSE) images of vertical section through Cretaceous laminite carbonate, (A) overview of whole thin section, showing location of (B).Arrow indicating starting point and direction of tile collection (image rotated 90° clockwise from original scanned orientation).(C) Percentage porosity for corresponding area in (B), comprising data from 6798 tiles.Note low porosity in lower section, with higher porosity laminae within the upper section.

Figure 5 .
Figure 5. (A,B) Backscattered electron (BSE) images of vertical section through Cretaceous laminite carbonate, (A) overview of whole thin section, showing location of (B).Arrow indicating starting point and direction of tile collection (image rotated 90 • clockwise from original scanned orientation).(C) Percentage porosity for corresponding area in (B), comprising data from 6798 tiles.Note low porosity in lower section, with higher porosity laminae within the upper section.

Figure 6 .
Figure 6.(A) Backscattered electron (BSE) images of vertical section through Cretaceous coquina carbonate.(B) Corresponding percentage porosity map, showing isolated high porosity areas within a low porosity matrix.Pores may be connected in 3D, although no evidence is seen for this in the 2D plot.(C) Backscattered electron (BSE) images of vertical section through Lower Carboniferous limestone, from Spireslack, Scotland.(D) Percentage porosity map across area illustrated in (C), showing heterogeneous porosity distribution and a tendency towards elongation in the vertical direction.Arrows indicate starting point and direction of tile collection (images rotated 90° clockwise from original scanned orientation).Note images and data maps (A,B) comprise 5,829 tiles, (C,D) 25,080 tiles.

Figure 6 .
Figure 6.(A) Backscattered electron (BSE) images of vertical section through Cretaceous coquina carbonate.(B) Corresponding percentage porosity map, showing isolated high porosity areas within a low porosity matrix.Pores may be connected in 3D, although no evidence is seen for this in the 2D plot.(C) Backscattered electron (BSE) images of vertical section through Lower Carboniferous limestone, from Spireslack, Scotland.(D) Percentage porosity map across area illustrated in (C), showing heterogeneous porosity distribution and a tendency towards elongation in the vertical direction.Arrows indicate starting point and direction of tile collection (images rotated 90 • clockwise from original scanned orientation).Note images and data maps (A,B) comprise 5,829 tiles, (C,D) 25,080 tiles.

Figure 7 .
Figure 7. (A) Backscattered electron (BSE) image of small vertical section through Quaternary aragonite cemented sandstone, containing bivalves and highly variable in porosity (P = large pores).Arrow indicating starting point and direction of tile collection in (B) and (C).(B) Corresponding percentage porosity map (individual tiles 2 mm).(C) Porosity map, as in (B) but individual tiles 1 mm.(D) Porosity map, selected area view marked by red dashed box in (A), red arrow indicating starting point and direction of tile collection, individual tiles 518 µm.Note progression from (B) to (D) showing change in clarity for which different sized pores are observed, and hence the importance of tile size chosen.

Figure 7 .
Figure 7. (A) Backscattered electron (BSE) image of small vertical section through Quaternary aragonite cemented sandstone, containing bivalves and highly variable in porosity (P = large pores).Arrow indicating starting point and direction of tile collection in (B,C); (B) Corresponding percentage porosity map (individual tiles 2 mm); (C) Porosity map, as in (B) but individual tiles 1 mm; (D) Porosity map, selected area view marked by red dashed box in (A), red arrow indicating starting point and direction of tile collection, individual tiles 518 µm.Note progression from (B) to (D) showing change in clarity for which different sized pores are observed, and hence the importance of tile size chosen.

Geosciences 2017, 7 , 70 13 of 20 Figure 9 .
Figure 9. Contour plots for various mineral distributions.(A) Coquina map showing percentage coverage by silica rich particles (quartz and clays).Note co-associations between silica particles and porosity (Figure 6b).(B) Percentage coverage of siderite from Spireslack limestone, across same area as in Figure 6d.Note vertically orientated concentrations in occurrence of siderite, corresponding to poorly visible brighter areas in Figure 6c.(C) Exploded area from Figure 6c, illustrating the occurrence and shape of siderite crystals within the limestone.Mineral phases identified through spot EDX analysis and grayscale.

Figure 9 .
Figure 9. Contour plots for various mineral distributions.(A) Coquina map showing percentage coverage by silica rich particles (quartz and clays).Note co-associations between silica particles and porosity (Figure 6b).(B) Percentage coverage of siderite from Spireslack limestone, across same area as in Figure 6d.Note vertically orientated concentrations in occurrence of siderite, corresponding to poorly visible brighter areas in Figure 6c.(C) Exploded area from Figure 6c, illustrating the occurrence and shape of siderite crystals within the limestone.Mineral phases identified through spot EDX analysis and grayscale.

Figure 10 .
Figure 10.(A) Plot of percentage porosity versus percentage siderite, from Spireslack limestone, showing strong clustering with a positive correlation.(B) Plot of percentage porosity versus percentage silicates (quartz and clays), from coquina sample.Data appears visually to be clustered (with some lateral "wings") and a positive correlation.Both show statistically significant positive relationship with Pearson r and Spearman rank correlation.Also note the high sample number (n) particularly in (A).

Figure 10 .
Figure 10.(A) Plot of percentage porosity versus percentage siderite, from Spireslack limestone, showing strong clustering with a positive correlation.(B) Plot of percentage porosity versus percentage silicates (quartz and clays), from coquina sample.Data appears visually to be clustered (with some lateral "wings") and a positive correlation.Both show statistically significant positive relationship with Pearson r and Spearman rank correlation.Also note the high sample number (n) particularly in (A).

Figure 11 .
Figure 11.Additional data maps, representing the mean grayscale across various thin sections.(A) Laminite, (B) coquina, (C) Spireslack limestone and (D) aragonite cemented sandstone from Antrim.In (A) blue areas correspond to high porosity regions within the laminite, while in (B) blue areas relate to localities with either high porosity and/or high silica contents.The green and orange areas are carbonate-rich.In (C) darker/lower values for average grayscale indicate a vertically orientated feature (right-hand side), which appears to correspond to an area that is relatively depleted in siderite.The darker mean grayscale in (D) is aligned with areas that are high in porosity, while isolated patches of brighter green match areas that are particularly well cemented with aragonite (AC), or represent areas of bivalve shell (S).Direction of scanning as in associated porosity maps.

Figure 11 .
Figure 11.Additional data maps, representing the mean grayscale across various thin sections.(A) Laminite, (B) coquina, (C) Spireslack limestone and (D) aragonite cemented sandstone from Antrim.In (A) blue areas correspond to high porosity regions within the laminite, while in (B) blue areas relate to localities with either high porosity and/or high silica contents.The green and orange areas are carbonate-rich.In (C) darker/lower values for average grayscale indicate a vertically orientated feature (right-hand side), which appears to correspond to an area that is relatively depleted in siderite.The darker mean grayscale in (D) is aligned with areas that are high in porosity, while isolated patches of brighter green match areas that are particularly well cemented with aragonite (AC), or represent areas of bivalve shell (S).Direction of scanning as in associated porosity maps.

Table 1 .
Image acquisition parameters used within Maps during the current work.

Table 1 .
Image acquisition parameters used within Maps during the current work.

Table 2 .
Percentage porosity (minimum, maximum, average and modal), based upon data from individual tiles from each sample set examined.

Table 3 .
Porosity frequency of tiles within specific percentage ranges (at 10% intervals) for all samples examined.Tile size indicated in brackets.