3D Quantitative Mineral Characterization of Particles Using X-ray Computed Tomography

A new method to measure and quantify the 3D mineralogical composition of particulate materials using X-ray computed micro-tomography (CT) is presented. The new method is part of a workflow designed to standardize the analysis of particles based on their microstructures without the need to segment the individual classes or grains. Classification follows a decision tree with criteria derived from particle histogram parameters that are specific to each microstructure, which in turn can be identified by 2D-based automated quantitative mineralogy. The quantification of mineral abundances is implemented at the particle level according to the complexity of the particle by taking into consideration the partial volume effect at interphases. The new method was tested on two samples with different particle size distributions from a carbonate rock containing various microstructures and phases. The method allowed differentiation and quantification of more mineral classes than traditional 3D image segmentation that uses only the grey-scale for mineral classification. Nevertheless, due to lower spatial resolution and lack of chemical information, not all phases identified in 2D could be distinguished. However, quantification of the mineral classes that could be distinguished was more representative than their 2D quantification, especially for coarser particle sizes and for minor phases. Therefore, the new 3D method shows great potential as a complement to 2D-based methods and as an alternative to traditional phase segmentation analysis of 3D images. Particle-based quantification of mineralogical and 3D geometrical properties of particles opens new applications in the raw materials and particle processing industries.


INTRODUCTION
The use of X-ray computed tomography (CT) for 3D visualization of microstructures in geological materials has a variety of applications in the geosciences (Cnudde & Boone 2013), for example, to measure porosity (Zhang et al., 2019) and pore-scale processes (Hasan Sharul et al., 2020;Anduix-Canto et al., 2021), to determine the distribution of and association among minerals (Reyes et al., 2017;Guntoro et al., 2021;Da Ferraz Costa et al., 2022), which are useful properties in ore geology, and to determine particle properties (Wang et al., 2017;Miller & Lin 2018), which are important in the mineral processing industry (Videla et al., 2007).
CT is based on the principle that each component of a sample has a specific X-ray attenuation that is proportional to its electron density and is representative of its composition (Withers et al., 2021). Thereafter, the component should be represented in a 3D image of the sample by voxels with a grey-value that is proportional to its attenuation coefficient (Dhaene et al., 2015). However, in practice, because of imaging artefacts and the lack of chemical information, a direct mineral classification based on the grey-scale of a 3D image is difficult (Boas and Fleischmann 2012;Godinho et al., 2021b). Therefore, it is not surprising that the wider applications of CT in Earth sciences are related to imaging of pores because air is the least X-ray attenuating component of a sample, and thus it is easy to classify. Contrarily, the application of CT in the fields of ore geology and minerals engineering requires identification of the different phases in often-complex samples (Bam et al., 2016;Warlo et al., 2021). To broaden the applicability of CT to such fields, it is necessary to develop standardized and automated methods that allow faster, yet reliable, classification and quantification of mineral phases from grey-scale 3D images.
According to the literature, different solutions have been demonstrated to help phase classification from grey-scale 3D images, for example, by direct correlation of the same slice in 3D with 2D classified images (Furat et al., 2018), by measuring chemical information from transmitted X-rays (Godinho et al., 2021b;Sittner et al., 2021), or by calculation of attenuation coefficients from the grey-scale (Pankhurst et al., 2018;Bam et al., 2020). Additionally, image processing methods are available to reduce quantification deviations caused by imaging artefacts (Lin et al., 2015;Godinho et al., 2019;Ketcham and Mote 2019;Voigt et al., 2020). Nevertheless, the complexity and variety of microstructures in natural materials make it difficult and laborious to standardize the available methods as universal solutions to every material.
In this research work, we present an alternative solution to quantify the mineral composition of particles based on a workflow that has the potential to be standardized and semi-automated. The new method consists of breaking the complexity of a sample into smaller sub-volumes (individual particles) that are simpler to analyse. Based on prelimi-nary mineralogical and microstructural information, the individual particles can be classified using particle specific criteria (an approach also used in 2D automated mineralogy) defined from particle histograms. The results are compared to 2D automated quantitative mineralogy, specifically using mineral liberation analysis (MLA), which is also used as a source of phase and microstructural information. The advantages and limitations of the new method relative to other 3D image segmentation methods and to 2D-based characterization are discussed in detail.

Sample Preparation
The workflow applied in this study is schematized in Figure 1. A rock sample (519 g), from the Snowbird Mine, Montana, USA / Canada (Samson et al., 2004), measuring approximately 7 9 7 9 8 cm was scanned using CT. The rock was crushed to particle sizes < 5 mm using a Retsch GmbH jaw crusher (BB200 Mangan, Germany) and then milled to < 1 mm. The milled particles were homogenized and split using a rotary splitter. To evaluate the effect of particle size in the results, two sub-samples were analysed: (1) the size fraction 710-1000 lm (named '' > 710 lm''), as an example of a sample with narrow size distribution; and (2) the size fraction 63-1000 lm (named ''AllSizes''), as an example of a sample with wider size distribution. The particle size distribution based on equivalent spherical diameter measured from the 3D images is given in Supplementary Figure SI1. For each sub-sample, one 30-mm-diameter grain mount was prepared for MLA (named ''GM > 710'' or ''GM-AllSizes''). GM > 710 was also scanned by CT to assess the representativity of the surface analysed by MLA.
The particle dispersions were prepared using a standardized procedure (Godinho et al., 2021a) that uses sugar particles as spacer in the ratio of 7 g of sugar and 1 g of sample particles. The resin used was a fast-curing acrylic polymer, Paladur (Kulzer, Mitsui Chemical Group). The particles of the material together with the sugar were mixed with methylmethacrylate-copolymer powder in the mass ratio 1:1. The methylmethacrylate liquid resin was added to the solid mixture in the ratio 3 ml to 10 g forming a paste saturated with particles. The final paste was poured into a tube with 12-mm diameter and then left to dry for one hour. The grain mounts were prepared by mixing 1 g of the sample with the same amount of graphite and epoxy resin. To prevent bias caused by settling, the so-called B-sections were used for surface analysis (Schulz et al., 2020).
The rock was scanned at a maximum X-ray energy of 175 kV using 1-mm copper and 0.8-mm-thick steel filters mounted at the source, with voxel size of 40 lm and beam current of 40 W. The two particle dispersions were scanned at a maximum X-ray energy of 160 kV, using a 0.5-mm copper filter, 15 W, 10-lm voxel size, and the pixels of the detector were binned in a 2 9 2 grid. The conditions of all scans were optimized so that the resolution was limited by geometrical factors and not by the spot size of the beam. Reconstructed images were processed as 16bit image processing and visualization using Avizo software (v9.3.0, Thermo Fisher Scientific, Waltham, Figure 1. Overview of the work done in this paper. Red panel: 3D images generated by CT of the sample before and after milling. Green panel: 2D SEM-based mineral analysis as source of preliminary mineralogical and microstructural information. Blue panel: preliminary information is converted into classification criteria following a decision three. Yellow panel: individual particle histograms are quantified according to the decision three. MA, USA) and Dragonfly software (v2021.1, Objects Research Systems, Montreal, Quebec, Canada). Labelling of particles was performed using a deep learning model with 3D sensor architecture as implemented in Dragonfly that was trained using 245,610 voxels, of which 41,770 represented particles.
The SEM used was a FEI Quanta 650F field emission SEM (FE-SEM) equipped with two Bruker Quantax X-Flash 5030 energy dispersive X-ray detectors. Quantitative mineralogical and textural data were extracted from particle maps, which were used to derive the modal mineralogy (Schulz et al., 2020). The grain-based X-ray mapping (GXMAP) mode was used with a resolution of 6 lm / pixel. Backscattered electron (BSE) images were collected with a resolution of 1 lm / pixel. Data were processed using the software package MLA Suite 3.1.4.686.

Workflow Background
The grey-scale of voxelized images obtained by CT are typically assumed proportional to the electron density of sample components . However, the correlation between grey-scale and phase composition can diverge due to image reconstruction artefacts and the spatial distribution of phases in the sample (Pankhurst et al., 2018). For example, in denser sections of a sample, beam hardening causes an artificial increase in the attenuation of a phase relative to lighter sections due to the preferential attenuation of lower energy X-rays (Godinho et al., 2021b). In practice, this broadens the range of grey-values that represent one phase, and it makes the grey-scale dependent on the location of a voxel within a sample. Another example is the boundary between two phases, which in the physical world is sharp but in CT 3D images is represented by a grey-scale gradient several voxels across the interphase due to partial volume and cone beam artefacts (Boas and Fleischmann 2012). In practice, interphases appear to be blurred, which means that it is impossible to define the real sharp boundary that defines 100% of a phase on either side of a boundary. The inaccuracy and bias of traditional segmentation methods that do impose such boundaries to interphases increase as the size of microstructures is close to the voxel size of an image (Godinho et al., 2019;Ketcham and Mote 2019). For instance, the smaller the grain is, the smaller is the ratio of voxels representing the bulk of the grain to voxels representing the interphase, down to a grain size limit where the interphase is dominant relative to the bulk. In this case, the grey-scale becomes dependent not only on the composition, but also on the grain size (Godinho et al., 2019(Godinho et al., , 2021b. A workflow named mounted single-particle characterization for mineralogical analyses (MSPaCMAn) has been shown to reduce the CT image artefacts described above (Godinho et al., 2021a). First, the particles are dispersed throughout a sample using a spacer to prevent particles from touching each other and to reduce the effect of beam hardening because all sections of a sample have similar densities. Additionally, because the particles are separated, the interphases between neighbour particles are reduced, thus segmenting the particles individually is facilitated. The second step is the labelling of the individual particles, which allows processing of their sub-volumes independently. The reduced complexity of each particle relative to the whole sample simplifies the interpretation of the grey-values, and it allows adding the geometrical properties of a particle as classification criteria. Importantly, the segmentation is done at the particle level and not at the grain level (like traditional segmentation methods (Wang et al., 2016;Guntoro et al., 2019;Da Wang et al., 2020)), which reduces the biased decision of defining boundaries at the smaller grain level.
To summarize, the MSPaCMAn workflow consists of (i) preparation of a particle dispersion sample, (ii) segmentation of individual particles, and (iii) phase classification using the histogram of individual particles and preliminary knowledge about the minerals present. In the next section, the fundamental principles of phase quantification at the particle level are described and appended as step (iv) of the MSPaCMAn workflow.

Quantitative Particle Analysis
Class quantification is proposed to be done directly from individual particle histograms that result from step (iii) of MSPaCMAn (as described in Sect. 3.1). From the differences and similarities of the histograms, the particles could be grouped based on their complexity. By correlation of the different groups that can be distinguished, with the prelimi-nary information about the sample, the different classes can be identified. To differentiate classes, various histogram parameters that are indicative of different microstructures in a particle can be derived. G max is the maximum grey-value of all voxels in a particle that can be used as indicative of the densest phase. G mean is the mean grey-value of all voxels, and it typically is close to the grey-value of the most frequent grey-value if the particle is composed of only one phase. G mean /G max informs on whether the particle is composed of one or more phases, i.e. values closer to 1 imply only one phase. G SD is the standard deviation of a histogram, and it gives an idea of whether a particle contains sharp peaks or if the peak frequencies are more distributed through different grey-values, which indicate fine microstructures without a dominant phase. G peak is the most frequent grey-value within a specified search range. For example, using the preliminary knowledge about which phases are present in a sample, their theoretical grey-value can be determined as the centre of the search range. It should be noted that only some examples of histogram derived parameters are presented; thus, other parameters could be more adequate to other materials.
The second step is to use the histogram parameters to create a sequence of logical conditions that allow distinguishing each class. This can be done following a decision tree, whereby each branch contains a set of logical conditions based on measurable properties, which could also include the particle geometrical properties. To help organizing the decision tree, five types of particles with different complexities are explained: type (1) coarse particles composed of one phase; type (2) particles composed of one phase with geometrical properties that affect the grey-scale; type (3) particles composed of two phases and large grains; type (4) particles composed of at least two phases distributed into small grains; and type (5) particles with higher complexity that do not fit the criteria of particle types (1)-(4).
The new quantification method distinguishes voxels representing interphases (affected by partial volume) from voxels representing pure classes. Thereafter, once a particle type and classes are determined, all voxels in a particle histogram are divided into two types: (1) Main peaks correspond to voxels with 100% of a class that is identified using the decision tree. Note that not all peaks necessarily correspond to a class simply because it appears within its expected range of grey-values. (2) Minor peaks, which are all voxels that are not included into the main peaks identified by the decision tree, are assumed to be composed of two classes in a ratio given by: This equation assumes a linear regression between the grey-values of the two classes (G peak1 and G peak2 ) so that the volume fraction of class 1 (V 1 fraction ) in a voxel with any grey-value (G i ) can be calculated (Godinho et al., 2019;Ketcham & Mote 2019). The k is a constant for each pair G peak1 and G peak2 , and it can be measured assuming that V 1 fraction is 100% when G i is equal to G peak1 . In practice, minor peaks represent voxels with a mixed composition at interphases. The underlying hypothesis of this method is that the quantification of the partial volume at interphases, without defining the boundaries between two classes, reduces the quantification bias and reduces uncertainty.
Particle type 1: The simplest case refers to particles with high ratio particle size to voxel size (ParVox), and they are composed of only one class. The histogram of such particles exhibits a single main peak that represents the central region of the particle, whereby G peak appears within a narrow range of grey-values that are specific to each class (Fig. 2). The interphase between particle and background is always present, and it is proportional to the surface area of the particle, e.g. note the red/blue rims around particle 2 and the corresponding segment (same colour) in the histogram. The volume of the interphase is represented by minor peaks on the left side of the main peak, i.e. towards the grey-value that is characteristic of the background. It is also observed that the ratio of peak width to height increases for more attenuating classes, e.g. the peak of particle 2 is wider and shorter than the peak of particle 1. This is expected due to intra-particle beam hardening that is higher for denser particles. To summarize, the quantification follows three steps: (1) the class is identified based on the grey-value of G peak ; (2) the volume of the main peak corresponds to 100% of the class; and (3) voxels at the interphase are composed of some per cent of the phase in the particle and some per cent of background, which are calculated using Eq. 1 and the voxelÕs grey-values.
Particle type 2: This refers to particles composed of only one class and with geometric properties that affect the characteristic G peak of the class. Two specific properties are exemplified in Figure 3, namely shape (particles 3 and 4) and size (particles 5 and 6). This type of particles is characterized by a high ratio of voxels at the interphase particle-background relative to the number of voxels in the bulk. For example, one can note the reduction in the yellow core relative to the blue/red rims in particles 5 and 6 ( Fig. 3) relative to the coarse particle 2 (Fig. 2). Because the main peak in the histogram corresponds to voxels in the bulk, the height of the peak decreases with particle size. Thereafter, below some particle size, the main peak is no longer distinguishable (yellow arrow, Fig. 3). Further decrease in particle size causes a progressive decrease in G max (orange arrow, Fig. 3). Similar to the effect of size on G peak , the thickness of a particle may also affect the ratio of the number of voxels at the interphase relative to the bulk. If both size and shape are below the thresholds affecting G peak , the ratio of surface area to volume (Surf/Vol) can be used as a geometric parameter to define specific classification criteria.
The partial volume in type 2 particles could also be calculated using Eq. 1 as long as the right class has been identified. This might be misleading due to the expected shift in histogram parameters, e.g. G peak and G max . In this work, this type of particles was excluded from the quantification. Nevertheless, it is hypothesized that the shift of histogram parameters may be correlated to geometrical properties affecting a specific class. It is worth noting that by using traditional image segmentation methods, type 2 particles would be misclassified due to the shift in grey-values, and the quantification would be strongly biased because the interphase with the background is larger than the bulk of the particle.
Particle type 3: This refers to particles composed of two classes and large grains. For example, in Figure 4, class B (yellow in 3D and light grey in 2D) is denser than class A (white in 3D and dark grey in 2D). In this type of particles, besides the interphase with the background, another interphase is present in the particle (red line), between classes A and B. Consequently, the right side of the main peak of class A broadens instead of the sharp decline noted in Figure 2 for type 1 particles. The grey-scale range at the interphase between class B and background overlaps with the grey-scale range at the interphase between classes A and B. Thus, in the histogram of all voxels in a particle, it is impossible to know which voxels correspond to each interphase, but this is necessary in order to apply Eq. 1 with the right G peak values (Fig. 5).
One solution is to simplify the particle histogram by dividing it into two: the histogram of surface voxels and the histogram of bulk voxels (Fig. 4). The bulk histogram can be retrieved easily by eroding the voxels on the surface of the particle that represent the gradient towards the background. The bulk histogram in Figure 4 shows two main peaks (100% classes A or B) separated by small peaks that correspond to the interphase between classes A and B. The surface histogram can be calculated by subtracting the histogram of all voxels from the histogram of the bulk. The result contains information from the interphase between both classes and the background. Consequently, two types of gradients overlap to some degree, i.e. the interphase A-background and B-background, which if not corrected can lead to an overestimation of A relative to B. In this work, as a first approach to correct for that overlap, the total amount of surface Figure 2. General example of 3D images of two type 1 particles composed of only one class each but with different densities. Each particleÕs bulk volume is represented by one peak in the histogram (white or yellow for particle 1 or 2, respectively), as well as a gradient towards the background (only noticeable for particle 2 by red/blue rims and the corresponding segment colours in the histogram).
at the interphase A-background and B-background is assumed to be the same as the volume ratio A/B measured in the bulk. Thus, splitting between bulk and surface voxels allows to identify to which interphase a voxel belongs to, which allows using Eq. 1 to calculate the partial volume of each class. In practice, the quantification is based on the segmentation of surface and bulk of the particle instead of the traditional CT methods that would segment classes A and B, i.e. by defining a boundary at some position along the red line.
Particle type 4: This refers to particles composed of at least two classes associated as fine grains with sizes close to or smaller than the voxel size. In this case, G peak can appear at any grey-value between the expected grey-scale ranges of the two classes, depending on which class is dominant and the size of the microstructures. In general, this type of particles Figure 3. Examples of type 2 particles with geometrical properties that affect their histogram parameters: particles 3 and 4 have thin shapes; particles 5 and 6 have small diameters. Particle 2 in the histogram corresponds to particle 2 in Figure 2 that has the same class as particles 5 and 6. . Example of a type 3 particle composed of two grains that are sufficiently large to produce two identifiable peaks in the histogram. Quantification requires differentiating the voxels at the surface from the voxels in the bulk.
has low G mean /G max and high G SD . Classifying these types of particles requires preliminary knowledge, e.g. from classified 2D images, about the typical microstructures in the sample in order to derive suitable histogram parameters that are specific to the microstructure. Once the classes composing the microstructure are identified, their theoretical greyvalues can be estimated and used as input to Eq. 1. Note that separating surface and bulk histograms is also necessary similar to the method followed for type 3 particles. Particle type 5: This refers to particles that do not fit the criteria of any of the other particle types.
For example, particles with more than two classes or particles that have two classes and also have a geometrical property that affects their grey-scale. The partial volume cannot be quantified in this type of particles because Eq. 1 requires attributing one or two classes to the composition of every voxel at an interphase. Alternatively, in these cases the particle can still be analysed using traditional segmentation methods without the partial volume correction. Possible future developments of the method to improve quantification of these type of particles are discussed later in the paper.

Results and Discussion
The scan of the rock showed the spatial distribution of some classes but it did not allow to identify those classes (Fig. 6). The rock was crushed, and a list of mineral phases present in it was obtained using MLA (Table 1, Fig. 7). The composition of the particulate material determined from MLA and CT was compared for two particle size distributions (Fig. 8) in order to evaluate which classes and microstructures can be distinguished and quantified using the new quantitative MSPaCMAn workflow. The preliminary information collected by MLA was used to define the type of particles and the classification criteria used for the analysis of individual particles in the 3D image (Figs. 9, 10). Table 1. List of phases identified by MLA, their abundance in mass per cent, theoretical X-ray attenuation coefficient for the ideal mineral formula ) and expected grey-value for an average mineral composition at the effective energy of the scan (E eff = 73 kV)

Rock Analysis
The rock sample has two eye-distinguishable regions, a white region composed of calcite that has some sub-millimetre veins, and a brown region that contains dark brown crystals (Fig. 6a). With a voxel size of 40 lm (the best resolution possible in the used scanner that included the entire sample in the field of view), three main classes of microstructures based on the grey-values could be distinguished: (1) The lighter main sample matrix (set transparent in Fig. 6b). Note that the white and brown regions visible in Figure 6a are not distinguishable. (2) Denser crystals dispersed throughout the matrix (pink). (3) Regions with grey-values between classes 1 and 2 that are composed of more than one phase (Fig. 6c, slightly different grey-values on top right of Figure 6c and purple-green-orange-red in Figure 6b). The denser phases show irregular distributions of grey-values, which indicate imaging artefacts (highlighted in Supplementary Figure SI2). For example, grains are brighter at their edge, which indicates beam hardening, and different grains have grey-values that depend on their position in the sample and the size of the grain. Thus, classification directly from the grey-scale would misclassify grains  . (a) Heat map of 11,000 individual particle histograms organized from left to right based on ascending equivalent diameter corresponding to sample AllSizes. Warmer colours correspond to higher frequency of voxels, i.e. analogous to a peak in a histogram. (b) Variation of maximum grey-value of each particle as a function of equivalent diameter. Note the increasing trend below the red line. The pink band is likely related to particles containing parisite, and the grey band is likely related to particles containing carbonate matrix.
with different composition as the same class, and grains with the same composition as different classes. To further characterize the mineralogy of the observed 3D microstructures, the sample was crushed and milled < 1 mm, and the particles were analysed by MLA and CT.

Mineralogy and 2D Microstructures
MLA identified at least 11 minerals in the particulate material (Table 1), which appeared in the most common microstructures shown in Figure 7.

Overview of Particle Mounts
The 3D images of the particle dispersions showed an approximately homogeneous distribution of particles throughout the sample (Figs. 8a, b), which reduced the imaging artefacts relative to the rock sample and simplified the segmentation of individual particles. Nevertheless, for sample All-Sizes, it was inevitable that some of the smaller particles touched the larger ones (e.g. red circle, Fig. 8d). Consequently, some segmented regions may have actually contained more than one particle, which may have added to the complexity of the histogram. It should be noted that in most cases, the difference in size of the touching particles was large, thus the particle histogram was dominated by the classes composing the larger particle, and the contribution of the finer particle was imperceptible. For the sample > 710 lm, which had a narrower size fraction, real individual particle segmentation was achieved. In general, the peaks of the histogram were better defined for sample > 710 lm (Fig. 8e), possibly because larger microstructures were preserved, thus the particles had a higher ratio of bulk voxels relative to surface voxels, i.e. less voxels affected by partial volume between particle and background.
The attenuation coefficient of the 11 phases identified by MLA were calculated assuming an effective X-ray energy of 73 kV as measured using a spectral detector (Godinho et al., 2021b;Sittner et al., 2021). Because the grey-scale of the 16-bit image (between 0 and 65,535) is directly propor-tional to the attenuation coefficient, the approximate grey-value corresponding to each phase was calculated (Table 1) and compared to the histogram of the particle dispersions (colour regions, Figure 8e). Note that the attenuation coefficients of the nickel-iron phases and the REE phases were approximate values for the theoretical mineral formula, which might differ from the actual formula of the mineral found in the deposit and can also be different between grains. Out of the 11 phases identified by MLA, only six peaks were distinguishable in the whole sample histograms. Several phases did not seem to be represented by any of the six peaks based on their expected peak positions (e.g. ankerite or gaspeite), other phases seemed to have the peak overlapping with each other (e.g. magnetite, liebenbergite and pyrite), and some phases appeared to have the peak position shifted (e.g. millerite and synchysite). The peak of parisite was very broad, which suggests that different grains may have variable element composition. The peaks of millerite and synchysite appeared to be shifted possibly due to a different elemental composition of the phases relative to their theoretical value calculated.
In conclusion, traditional grey-scale-based mineral classification strategies that were applied to the entire 3D images would only result in a very rough quantitative mineral analysis. From the 11 phases identified in the MLA, only millerite, synchysite and parisite were represented by a distinct range of grey-values. The other phases would require grouping according to the distinguishable peaks in the whole sample histogram (Table 2). Additional uncertainty would be expected because voxels at interphases had grey-values that may be out of the range of the corresponding class. Next, the quantitative MSPaCMAn method described above was applied to samples > 710 lm and AllSizes to improve the accuracy of quantification and to increase the number of phases that can be distinguished.

Particle Histograms vs 3D Microstructures
The first step of quantitative MSPaCMAn is the assessment of similarities and differences between individual particle histograms (Figs. 9, 10a). Ideally, particles with similar compositions would have similar histograms. However, type 2 particles, by definition, have histograms that are dependent not only on their composition and microstructure, but also on their geometrical properties. These particles should be identifiable, and so they can be processed according to specific criteria (Fig. 9). The second step is the grouping of similar histograms and their corresponding microstructures using the 3D image of the particles to confirm the similarity of microstructures (Fig. 10). Third, the 3D microstructures of each group are compared to the known microstructures obtained by MLA (Figs. 7,  10c). Fourth, for each group, the characteristic histogram parameters are found, and the end member values are used to define the classification criteria that are specific of each class (Fig. 11).
Comparing many particle histograms together can be achieved, for example, using heat maps (Figs. 9a, 10a). Each vertical line corresponds to the histogram of a particle, whereby brighter points correspond to higher frequency of a given grey-value in a particle, i.e. higher peaks in traditional histogram representations. By ordering the particles from left to right based on particle size (or other geometrical property), it is possible to also assess the presence of type 2 particles. For example, the heat map of sample AllSizes showed a gradual shift of the histogram peaks towards higher grey-values for coarser particles (white line in Figure 9a, pink band in Figure 9b). This was the consequence of the behaviour described in Figure 3 for type 2 particles. Thereafter, G max increased as a function of the equivalent diameter of the particle up to approximately 110 lm (Fig. 9b). Interestingly, this threshold diameter appeared to be higher for more attenuating phases (note inclination of the red dashed line). Particles with equivalent diameters larger than 110 lm have similar G max , which depends only on the densest phase in the particle, e.g. grey and pink regions in Figure 9b, each corresponding to one class. In conclusion, particles < 110 lm have his-togram parameters that depend not only on the particle composition, but also on its size. In the current work, these particles (about 20% by volume, Supplementary Figure S1) were excluded from the quantification described in the next section.
Because the particles in sample > 710 lm were coarser than 110 lm (Supplementary Figure S1), the histograms were independent of particle size as shown by the lack of gradual variations of the histograms from left to right in Figure 10a. Thereafter, this sample was used to define the classification criteria. The majority of particle histograms had a single peak between 2300 and 5100 (white arrows, Fig. 10b), which was within the range expected for calcite and dolomite, the most abundant phases in the sample. Note the slight variability of the peak position and width ( Fig. 10a and zoom of P1 and P2 in Figure 10b) was likely due to variable ratios of ankerite, calcite and dolomite in microstructures similar to particle 7 in Figure 7. The fine microstructures between ankerite and dolomite implied that the grey-values representing these microstructures were between the values of the pure phases, which coincided with the range expected for calcite. Therefore, the three carbonate phases cannot be distinguished at the scan resolution and was grouped further in a class named carbonate matrix.
Several particle histograms had a peak between 13,000 and 14,000, which was expected for pyrite (e.g. purple arrows, Fig. 10). These peaks can be narrow for type 1 particles (e.g. particle 3, Fig. 7; and P3, Fig. 10) or wider towards lower grey-values for particles containing other phases. For instance, if the veins containing magnetite or gaspeite are not well-resolved (e.g. particle 2, Fig. 7), then a single peak at a position between magnetite and pyrite is observed. A second peak associated with particles containing pyrite appeared typically at the range expected for liebenbergite (e.g. P4, Fig. 10). How- Table 2. Best possible classification of particle mounts using either only the grey-scale (traditional segmentation), or quantitative MSPaCMAn. Phases identified by MLA are coloured as in Figure 7, and the classes defined by MSPaCMAn are coloured as in Figure 10.
The symbol marks the phases that can be identified without grouping ever, it is known from MLA that this phase is mostly associated with the carbonate matrix or liberated and rarely associated with pyrite. Therefore, the second peak was interpreted as corresponding to the typical assemblage of magnetite and gaspeite (e.g. particle 1, Fig. 7), which are always associated with pyrite and manifests as a peak between the greyvalues of magnetite and gaspeite. Because both phases had a similar effect on particle histograms, they were grouped as one class named Mag-Gasp. A repeating histogram feature was identified by orange / yellow arrows. These histograms showed a wide range of grey-values between about 7,700 and 20,000, and the dominant peak appeared at variable grey-values. Such features are typical of type 4 particles. Comparing the microstructures in these parti-cles with MLA imaged an obvious similarity between P5 can be seen in Figure 10c and particles 10 and 11 in Figure 7. The microstructures were composed of the typical jamborite-millerite association. Because the ratio of both phases varied in different particles (compare particle 10 and 11), it was not surprising that the dominant peak can appear at any grey-value between the expected range of jamborite and millerite depending on which phase was dominant. The smaller the grain sizes, the less likely it is that a dominant peak can even be distinguished. This would be the case if all voxels corresponded to interphases between grains inside the particle.
To summarize, by comparing the microstructural information from MLA and CT with the particle histograms, it is possible to distinguish six Figure 11. Scheme of the decision tree used to distinguish the different classes. The boxes are coloured as in Figure 10 according to the class identified in each branch.
individual phases and eight classes following the MSPaCMAn workflow (different colours in Table 2). This is particularly relevant considering that some phases had similar attenuation coefficients or were present as microstructure assemblies represented by grey-values that overlapped with other phases, e.g. jamborite/millerite and pyrite. Relative to traditional classification based only on the greyscale of the whole image, quantitative MSPaCMAn doubled the number of phases that can be quantified individually in this particular mineral system (Table 2). Next, the classification criteria that are specific for the histograms of each class are defined.

Histogram Analysis
The first condition at the top of the decision tree (Fig. 11) is to divide the dataset into particles with equivalent diameters finer than 110 lm and coarser than 110 lm. To prevent phase misclassification when using particle histogram parameters, as G max , G peak and G mean depend on particle size, the finer fraction was not analysed in the current work. The coarser particle size fraction was passed to the next branch of the decision tree. In sample AllSizes, only 3377 particles (about 4% of the total number of particles) obeyed to the size criteria, which corresponded to about 81% of the total volume of particles segmented (Supplementary Figure S1). That means any difference in composition between coarser and finer size classes would be diluted about five times in the uncertainty of the composition of the sample if calculated only from the coarser particle size fraction.
The second layer of the decision tree consists of screening the G peak values. If only one peak is found (type 1), then G peak should appear within the ex-pected grey-scale range of a class (Table 2, Fig. 10). In general, type 1 particles of the classes at the lowest and highest attenuation could be classified easily because they did not overlap with types 3-5 particles. In this particular mineral system, this corresponded to the carbonate matrix and to the parisite classes that were mostly liberated and represented > 94% of the samples masses (Table 3). In principle, other ranges of grey-values, if unique to a class of liberated particles, could be defined, although that was not the case in the studied material. For instance, even though the grey-scale range of the synchysite class was unique, the class did not occur liberated; and the ranges of the other classes were overlapping and required more selective criteria.
The third layer of the decision tree aims at discriminating particles according to the complexity of their microstructures. As defined above, the ratio G mean /G max can be used as an indicator of the presence of more than one phase, i.e. type 1 particles have G mean /G max closer to 1. In this case, 0.8 was set as the threshold of G mean /G max , above which particles were observed to be of type 1. Most particles obeying this criterion were liebenbergite or pyrite, which appeared at close, yet distinguishable, greyranges that were used to discriminate them (Fig. 11). If G mean /G max < 0.8, then particles have at least two phases. The millerite-jamborite microstructures can be distinguished using a G max threshold between 14,500 and 20,000, which indicates the presence of millerite.
The particles that remain unclassified require additional criteria that add specificity to the microstructures expected from MLA analysis. The typical liebenbergite-carbonate matrix association (P2, Fig. 10c) was characterized by the presence of the G peak of the carbonate matrix and a G max < 12,000 that corresponded to small grains of lieben- bergite without a well-defined G peak2 . The typical pyrite-magnetite-gaspeite association (P4, Fig. 10c) was characterized by a broader peak or two narrow peaks at grey-scale ranges lower than those expected for pyrite. Other combinations of two classes (type 3) in a particle are possible and accounted for as long as two peaks are distinguishable and appear within a grey-scale range close to the range expected for a specific class, e.g. large grains of carbonate matrix and parisite in particle P8, Figure 10. It is also possible that no G peak is found if the peaks are not well-defined, which can happen for types 2 and 4 particles if the grain sizes were too small. The selected conditions allowed to classify 1562 of the 1574 particles in sample > 710 lm, and 2668 of the 3377 particles that passed the first size criterion for the sample AllSizes. The particles that did not obey the given classifiers can be either disregarded or segmented using traditional thresholding methods. In sample > 710 lm, the 12 particles that could not be classified were of type 5 (e.g. P6, Fig. 10) and were composed of more than two phases. The grains in these particles were segmented manually by greyscale thresholding similar to traditional segmentation of 3D images; thus, it is associated with the same limitations already discussed above. Nevertheless, because this was only applied to a smaller dataset, the limitations have minor impact.
In sample AllSizes, a higher percentage of particles was not classified because, as observed in Figure 8d, fine particles were often touching the coarse ones. Consequently, if the fine touching particles were composed of a class different from the coarse particle, it added new degrees of complexity that were not visible in the MLA microstructures. For instance, new associations between phases could result in histogram parameters outside the ranges defined in the decision tree. Nevertheless, in most cases, the difference in size of the touching particles was so large that the histogram was dominated by the classes composing the large particle, which in most cases were liberated carbonate matrix. Thus, the application of traditional threshold segmentation to this case was reasonable because the contribution of the fine particles was small.

Particle-based Phase Quantification
Once the classes in each particle were identified using the decision tree, the quantification was done according to the particle type, following the steps described above. The mass content of the different classes quantified from CT and MLA were compared in Table 3 for the two samples. In general, sample AllSizes showed better agreement between the composition measured with CT and MLA (less than 30% difference per class) than sample > 710 lm. To understand the differences, the effect of particle size, resolution and the representability of the data were analysed for both techniques.
To assess the representativity of the surface of the grain mount GM > 710 analysed by MLA, the sample was also analysed by CT (Fig. 12a). One subvolume, two voxels thick (20 lm) that contained the surface of the grain mount (V1), and two other subvolumes parallel to V1 but at different depths (V2 and V3) were compared (Figs. 12b, c). Because the grain mounts were prepared without spacer, the particles were touching each other; thus, quantitative MSPaCMAn could not be used. Therefore, only three classes could be segmented by applying greyscale thresholds. Note that this rough quantification  Figure SI3). Note that only few particles of each mineral were present in the grain mount and not all were present in each cross-section. of three classes was sufficient to assess the variability within the grain mount. Figure 12 shows asymmetries in the distribution of phases throughout the grain mount, which resulted in significant variation of composition at each cross section of the grain mount (Table 4). The lack of representability of each section was caused by the low number of particles carrying some of the phases (Blannin et al., 2021). Based on the CT scan, it was calculated that within a sub-volume 20-lm thick (assuming an electron penetration beneath the surface of up to 20 lm), only about 4.8 mg of material was actually considered in the MLA. This contrasts with a more representative 52 mg of material measured in the particle dispersion of sample > 710 lm with CT. In conclusion, some of the discrepancies between MLA and CT in Table 3 could be due to the lack of statistical representativity of minor phases on the surface analysed by MLA.
In sample GM_AllSizes, the number of particles on the cross section was much higher; thus, the results from MLA were more representative than those from sample GM > 710. On the other hand, the CT scan only measured 8.1 mg of material, about 6.4 time less than sample > 710 lm. This was due to the difficulty to disperse particles with a wider particle size distribution and because the size fraction < 110 lm was not accounted for. Therefore, the new method is more suitable to analyse samples with narrower particle size distributions and with sizes that are significantly larger than the voxel size of the scan. Note that due to less complex microstructure of fine particles, the sources of uncertainty in the CT scan were expected to be reduced, which could also explain the better agreement between CT and MLA for the sample AllSizes.

DISCUSSION
The new quantitative MSPaCMAn workflow was tested to quantify the mineral composition of two particulate samples from the same source rock, but with different particle size distributions, and the results were compared to MLA. First, the novelty and the advantages of the method relative to traditional CT quantitative analysis are discussed below. Second, the limitations of the method are analysed and compared to 2D mineralogy. Third, further developments of the MSPaCMAn workflow are proposed.

Advantages Relative to Traditional CT Quantification
The quantitative extension of the MSPaCMAn workflow has two key differences relative to traditional quantitative 3D image analysis. (1) The classification is undertaken at the particle level, which simplifies the interpretation of simpler particle histograms, as opposed to assigning a uniform interpretation to the more complex whole sample histogram representing all microstructures. (2) The classes are quantified without the need to segment individual grains, which allows calculating the partial volume of voxels at different interphases, thus avoiding biased user input of hard boundaries between classes.
Individual particle classification allows to identify specific microstructures that have a grey-scale signature that overlaps with the signature of other classes or microstructures in the sample, e.g. different classes represented by overlapping grey-ranges that appear as one peak in the histogram of the whole image. The refined classification is achieved by identifying the unique set of histogram-based criteria that are specific of a class or microstructure. Because several parameters can be used, the criteria are more selective than using simple threshold ranges used in traditional segmentation. This improvement is shown, in the analysed samples, to be advantageous for the cases that (a) even though the peaks overlap, the G peak is slightly different (e.g. liberated liebenbergite and pyrite), (b) the typical association with other phases causes a characteristic shift of G peak (e.g. Mag-Gasp association with pyrite), and (c) specific microstructures that cannot be resolved at the scan resolution, yet are represented by a distinctive histogram pattern (e.g. milleritejamborite microstructures). In conclusion, the particle-based classification of this material allowed quantification of jamborite, liebenbergite, pyrite, millerite and the magnetite-gaspeite class, which would not have been distinguishable nor reliably quantified using traditional 3D image processing methods.
It is important to consider that even for samples prepared using the same material, the microstructures may be dependent on the particle size classes, i.e. coarser particles are most likely to contain more complex microstructures, while finer particles are more likely to be liberated. Therefore, it is proposed that narrower particle size distributions should be used in order to narrow the type of microstructures possible, thus improving the likelihood of a class to be recognized by the classification criteria. Additionally, coarser particles are in principle easier to segment. However, above some particle size, the complexity of microstructures may increase to the point where quantitative MSPaCMAn no longer presents an advantage relative to traditional image processing methods. In conclusion, the applicability of the method should be assessed based on particle size distribution, the grain sizes and the resolution of the CT scan.
Another advantage of quantitative MSPaC-MAn is that the volume of each class is quantified directly from particle histograms without grain segmentation. This allows calculating the partial volume of voxels at interphases using Eq. 1, which is automatically applied to every particle according to the classes present. In a different manner, traditional 3D quantitative analysis requires user-biased input to define the boundary between grains, each containing 100 per cent of a class, e.g. grey-scale thresholds or ground truth mask used to train AIbased segmentation models. The quantification bias is also reduced because the segmentation of particles is more accurate, and it generates less interphases than segmentation of individual grains within a particle. Importantly, the boundary of all particles has a common neighbour phase (the sample matrix), which has a constant G peak that is independent of particle composition. This allows treating voxels at the boundary separately from voxels inside the sample, which reduces further the quantitative uncertainty at boundaries. Additionally, avoiding the segmentation of grains also enables the quantification of phases that appear as small grains in type 4 particles, which cannot be segmented accurately by traditional methods.

Comparison between 2 and 3D Mineralogy
Deriving classification criteria from particle histograms is shown to improve the ability to quantify additional phases that would have not been distinguishable from the grey-scale alone. However, this approach also carries risks that must be evaluated per material. To summarize, the accuracy of MSPaCMAn to quantify mineral compositions is affected by three general factors: (1) the contrast between phases and the spatial resolution in relation to the size of particles and grains; (2) the representativity of the preliminary mineralogical information in relation to the material and to the volume actually analysed by the CT; and (3) the complexity of microstructures that may reduce the universality of the classification criteria. Factors 1 and 2 are typically linked because, similar to other imaging techniques, increasing sample representativity implies larger scanning fields of view, worse resolutions and longer scanning and analysis times.
Even though MSPaCMAn improves the capability to distinguish phases and increases the accuracy of quantification relative to other 3D image analysis methods, it is clear that it lags behind 2D MLA when it comes to the ability to identify mineral phases. Depending on the complexity of the material, it is reasonable to expect that techniques capable of measuring direct chemical information or have higher spatial resolution would be able to distinguish more phases than CT-based methods. Nevertheless, this cannot be described as a limitation of MSPaCMAn but rather contextualized in the scope of each material. For example, the need to group dolomite, calcite and ankerite is not due to a physical limitation of the method but rather due to the small grain size and association of dolomite and ankerite. Note that the MLA uses a step size of 1 lm to resolve grains. Arguably, if the CT scan was done with a higher resolution, it could have been possible to discriminate the three carbonate phases because their attenuation coefficients are sufficiently different. However, it would also be expected that due to the smaller fields of view required to increase the resolution to 1 lm, the amount of material analysed may have been less representative, similar to the observations in Figure 12 for MLA.
The classification criteria are only as reliable as the preliminary information used to interpret the particle histograms. For instance, a wrong classification of the MLA would be transferred to a misclassification in the CT. This could have been the cases of synchysite and of millerite, for which the observed G peak did not match the value expected based on the theoretical attenuation coefficient. It is possible that those deviations were due to variable element composition within the same phase that may be sufficient to affect the attenuation coefficient and perhaps the density. Consequently, the value of density used to convert volume into mass may also contribute to some of the discrepancies observed between CT and MLA. In this sample, this was not expected to be significant for synchysite and millerite due to their relatively small amounts. Nevertheless, because some phases must be grouped into classes for which an overall density must be assumed, it is natural to have some differences in results from CT and MLA. For example, due to the high weight of the carbonate matrix in the mineralogy, even small variations of the ratio dolomite, calcite and ankerite that were detectable by the MLA but not by the CT could affect the overall mass balance of all classes.
In conclusion, similar to 2D MLA, quantitative MSPaCMAn is affected by the same sources of uncertainty, e.g. spatial resolution versus field of view and sample representativity, and classificationbiased input. However, the particle-based analysis and the imposition of classification criteria are expected to minimize those sources of uncertainty. The improved accuracy of quantitative 3D mineralogical analysis is a necessary step to the acceptance of CT as a suitable complementary technique to standard 2D automated mineralogy methods. This suitability is of course subject to the adequacy of material and the resolution of the CT scanner to the scientific question to be answered.

Further Developments
This paper lays the ground to broaden the applicability of quantitative 3D automated mineralogy. Further understanding of how the geometric properties of particles affect their histogram parameters may reduce the particle size threshold that is excluded from the quantification. For example, by establishing the function defining the pink region in Figure 9b, along which the equivalent diameter is proportional to G max , parisite could be quantified for particles finer than 110 lm. Naturally, there is a limit where the functions from the different phases overlap, which signifies that the classes cannot be distinguished. Nevertheless, determining the size threshold below which the quantification is expected to be wrong is crucial to any trustable analytical technique, and it represents an additional advantage relative to traditional CT quantification.
This method has the potential to perform better than 2D image methods to analyse coarse particle sizes, for which minor phases may not be representative. However, coarse particles of materials with high complexity are more likely to contain more than two classes per particle. Therefore, extending the method to the quantification of multiple classes per particle would boost its applicability to characterize ore samples. The identification of various G peak per particle is already possible. However, recognizing voxels with partial volume corresponding to each interphase is not yet possible without segmenting the individual grains. The development of spectra deconvolution methods specific to CT particle histograms may solve this challenge.
The application of this method to increasingly complex particulate materials may require decision trees with more branches and possibly more robust criteria to narrow the classification of very specific microstructures. One possibility is to use statistical analysis methods to automatically cluster particle types based on a combination of both geometrical properties and histogram derived parameters. This would also contribute to define the end members of a class, which would substitute the need to define the ranges of a parameter for each class based on visual observation. Therefore, this step is crucial to make the method more automated and less input bias, which are two key necessities on the path of MSPaCMAn to become a standardized workflow that is acceptable and broadly applicable in the raw materials community.

CONCLUSION
The new quantitative extension of the MSPaCMAn workflow consists of breaking the complexity of a sample into smaller sub-volumes (individual particles). Because microstructures in a particle are simpler than those in the whole sample, the classification can be refined using histogrambased particle properties instead of only grey-value thresholds. Additionally, voxels that are affected by partial volume can be distinguished in the particles that contain a maximum of two classes without the need to segment individual grains. Therefore, the partial volume of each class can be quantified with less uncertainty and without the biased definition of boundaries between classes associated with 3D image segmentation. This was demonstrated to be possible using quantitative MSPaCMAn, even for classes with overlapping grey-ranges and for grains close to or below the voxel size of the scan.
Despite the improved ability to distinguish minerals and the increased accuracy of volume quantification relative to traditional quantitative 3D image analysis, MSPaCMAn still does not bring CT as a standalone mineral characterization technique because preliminary 2D microstructural information is necessary. Yet, it is shown that more statistically relevant information can be obtained from a 3D image of particle dispersion than from a cross section of a grain mount analysed by MLA, especially if particles are large and if some phases are present in small amounts.
The overall balance between advantages and limitations of the workflow must be judged case by case based on the particle size, type of microstructures present, phase composition and the scientific question to be answered. Nevertheless, the new method shows potential as an alternative to traditional 3D analysis for particulate materials, and it shows complementary advantages relative to wellestablished 2D-based mineral quantification techniques. In particular, considering that the 3D geometrical properties and the 3D surface liberation of the individual particles can also be measured in addition to the mineralogy, quantitative MSPaC-MAn may become an impactful tool in the raw materials industry and particle technology.

Conflict of Interest
The authors have no competing interests to declare, and no financial or non-financial interests to disclose.

OPEN ACCESS
This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecom mons.org/licenses/by/4.0/.