3D-characterization of three-phase systems using X-ray tomography: tracking the microstructural evolution in ice cream

Abstract

The microstructure of food is key to its sensorial perception, and methods to characterize the microstructure are of crucial importance in food engineering. Ice cream is a special example whose microstructure changes dramatically in response to temperature variations. Since ice cream is a multiphase material, the complex interactions among the phases and the physical mechanisms that drive the evolution of microstructure are not yet well understood. This is mostly due to the fact that observing the microstructure with traditional microscopic methods is destructive and does not allow the study of undisturbed samples. With X-ray micro-tomography, it is possible to overcome these limitations and carry out time lapse studies of the evolution of the microstructure of ice cream. Using iodine as a contrast agent, we measured the three-dimensional distribution of the three main phases (air, unfrozen sugar solution, and ice crystals) with a voxel size of 6 μm. An automated routine was developed that allows for the segmentation of the three phases. Based on the three-dimensional data we calculated the temporal evolution of air bubble sizes and ice crystal sizes during cyclic variations of temperature. Under the given temperature variations we find strong hints that for ice crystal coarsening a melt refreeze mechanism and for air microstructure coarsening coalescence are the dominating underlying mechanisms. This method—which can be applied to a plethora of soft multiphase materials—provides new insights into the coarsening mechanisms of multiphase materials and could contribute to a better understanding of complex materials.

---

1 Introduction

The microstructure of food is one of the key factors that determine its sensorial perception. Therefore, controlled creation and stability of food microstructure is of great interest to the food industry. One method to create and stabilize a specific food microstructure is freezing and subsequent frozen storage. Due to the large number of food ingredients, their mutual interactions, and the non-equilibrium nature of food materials, the mechanisms which govern the behavior of food microstructure are not yet fully understood. This is particularly true for frozen food both during the freezing process and during the long term evolution in the frozen state.

Frozen food microstructure is often very complex and can consist of a large number of distinct phases; ice crystals, air, unfrozen freeze concentrated solution (further on termed unfrozen matrix) and fat particles. The unfrozen matrix itself is usually a complex mixture containing small sugars like mono- and oligo-saccharides, water soluble polymers like proteins and hydrocolloids and small molecular weight surfactants. The different phases can be either dispersed or continuous. In both cases they exhibit characteristic sizes in one or more dimensions which typically are on the order of 1 to 100 μm. Ice cream, the main subject of this study, contains four phases: (1) the unfrozen matrix, (2) air bubbles, or more precisely air microstructure with typical sizes between 20 and 150 μm,¹ (3) ice crystals with sizes between 10 and 75 μm, and finally (4) fat particles with sizes between 0.4 and 4 μm as well as fat particle aggregates with a size of around 10 μm.^2,3 It is assumed that the unfrozen matrix is continuous. For ice crystals and air microstructure it is not evident whether they constitute dispersed or continuous phases, and in addition it depends on the exact composition and the environmental conditions (temperature, storage time).

During freezing and thawing this microstructure is changed dramatically.⁴ The reason is that the unfrozen matrix at a typical storage temperature of T = −18 °C is well above the glass transition temperature, which is typically around −35 °C. Above the glass transition temperature quite rapid, irreversible coarsening of the ice crystals and the air microstructure occurs, which is particularly pronounced in cyclic warming and cooling conditions as they are commonly found in commercial freezers. This yields deteriorations of the ice cream microstructure like shrinkage (loss of air content) and a coarse and icy texture caused by ice crystals larger than around 50 μm.^4–7 Therefore it is important for the food industry to achieve a better understanding of the structural changes and the structural stability during these processes.

To a large extent this gap in understanding the mechanisms which govern the behavior of frozen food microstructure is due to the lack of methods which are capable of visualizing the three dimensional (3D) microstructure in a direct and non-destructive way in a realistic product environment. Up to now most of the techniques that are used to investigate the microstructure of frozen food are either two-dimensional or indirect or destructive or rely on preparation techniques that introduce artifacts. A number of different microscopic techniques have been used to characterize the microstructure of frozen food. Among them are cold stage light microscopy coupled with a dispersion technique,⁸ direct light microscopy,⁹ micro slicing techniques to reconstruct 3D microstructure from 2D images¹⁰ and several approaches using electron microscopy.^1–3 The reader can find more detailed information on methods to investigate frozen food microstructure in several reviews.^{2, 3, 8, 11–13}

Direct information about the 3D structure of ice crystals during freezing of sucrose solutions has been reported with help of magnetic resonance imaging (MRI),¹⁴ but the technique is limited due to the low resolution of 78 μm which is larger than most structures of interest. Regardless of the above mentioned limitations all these techniques have yielded a tremendous amount of insight about frozen food microstructure.

As already stated, the typical size of frozen food microstructures is often larger than 10 μm which renders them accessible to the resolution of X-ray micro-tomography with a voxel size in the range of 3 to 6 μm. This resolution limit can in principle be pushed to the sub-micron range in large synchrotron facilities. X-ray micro-tomography is already widely used for the investigation of food structures,^15,16 although its use has been restricted to ambient temperatures. In addition most of the studies have been limited to the investigation of air microstructure, owing to the fact that the contrast between air and condensed materials is high. The ice microstructure for example has only been investigated in an indirect and destructive manner using freeze drying of frozen food materials.^16,17

In this paper we describe the application of X-ray micro-tomography to investigate the microstructure of ice cream and its evolution with time, including ice crystals. For this purpose we apply a cyclic temperature variation further on called heat shock protocol to an ice cream sample. By means of iodine as a contrast agent, we were able to distinguish between three different phases in a model ice cream, namely the ice crystal phase, the unfrozen matrix phase and the air phase. We developed a novel image processing protocol for the segmentation and quantitative evaluation of such a three phase microstructure.

The methods developed in this work are not limited to ice cream and other frozen foods, but will be beneficial for the general study of random three-phase materials, where the same problems in image analysis occur.

2 Materials and methods

As introduced in Section 1, the four main phases in standard ice cream are air bubbles, H₂O ice crystals, unfrozen matrix and fat particles. Although fat particles play a significant role in stabilizing the ice cream microstructure the typical sizes of the fat particles is below 10 μm and therefore below the resolution of the used X-ray micro-tomograph. Due to this limitation we will not consider the fat particles in this paper.

2.1 Model ice cream

2.1.1 Recipe. Table 1 contains the recipe for our model ice cream which is based on typical industrial ice cream recipes used in Europe. All quantities are given in weight percent. The iodine salt is added as a contrast agent and is replacing some of the small molecular weight sugars of industrial ice cream recipes. The total molality of small molecules (below or equal to the size of disaccharides) for our recipe is roughly 1.4 mol kg⁻¹. The melting point of the recipe has been determined with differential scanning calorimetry to be T_m = −2.94 ± 0.3 °C (see Section 2.1.4).

Table 1 Recipe for model ice cream with iodine as contrast agent

Ingredient	Weight percent
Water	65.0
Dry glucose syrup (DE40)	10.1
Whey protein replacer (16% protein)	9.5
Coconut fat	7.8
Sodium iodine	3.2
Skim milk powder	2.3
Dextrose	1.6
Emulsifier/stabilizer	0.5

---

2.1.2 Ice cream production process. All dry ingredients and the fat were mixed with deionized water at T = 65 °C. After a minimal hydration time of 1 h the ice cream mix was homogenized using a two stage high pressure homogenizer (APV, type APV-mix) at 150 bar and 50 bar respectively and pasteurized at T = 86 °C for 30 s. The mix was then rapidly cooled to 4 °C and stored at that temperature over night for maturation. The ice cream mix was then aerated and frozen in a continuous surface scrape heat exchanger (Hoyer, Technohoy MF 50). An open dasher rotating at 750 rpm was used. The back pressure has been set to 1 bar, the outlet temperature to T = −7 °C and the overrun was set to 100% which corresponds to a volume based air fraction of 0.5. The ice cream was then filled into plastic pots of 130 ml and then quiescently cooled to T = −50 °C for storage.

2.1.3 Absorption properties. The linear absorption coefficient μ for the ice crystals and liquid solution are very similar for photon energies between 10 and 75 keV, characteristic for the tomography equipment used in the present study (cf. Section 2.2). For example, at 30 keV, a 40% sucrose solution at 0 °C has a 15% higher absorption coefficient compared to water ice. Therefore we added a small amount of sodium-iodine salt (NaI) to the solution. Iodine is well known for its large X-ray absorption coefficient, and has a very low solubility in solid H₂O ice. Iodine is expected to be freeze-concentrated in the unfrozen matrix, and thus to produce a better X-ray absorption contrast between the pure ice crystals and the liquid solution.

However, care must be taken to avoid high concentrations, to prevent notable modifications of the physical properties of the ice cream. The linear absorption coefficient of iodine is about two orders of magnitude larger compared to the other constituents like sucrose and water. With a density five times larger than that of water, the iodine absorption per mass unit is thus about 20 times larger. To achieve a (desirable) two-fold increase of the mean absorption of an iodine-water solution compared to pure water, a mass fraction of 1/19 or roughly 5%-wt of iodine should be added. Experimental comparison of different concentrations between 0.5 and 5%-wt both in frozen sucrose solutions and in model ice creams led to a good compromise of 3.2%-wt NaI concentration. It is important to note that the actual NaI concentration varies with temperature due to the freeze concentration in the unfrozen matrix.

Fig. 1 Reconstructed X-ray tomography slice of model ice cream with 3.2%-wt NaI added to the recipe. The dark regions correspond to air bubbles, the brightest regions denote unfrozen matrix. The overview image shows boundary effects at the interface, resulting from the extraction of the sample from a larger block with a saw at -60 °C. The subvolume for evaluations was taken from the undisturbed center. The scale bar on the right is 100 μm.

In the following, ice cream with 3.2%-wt NaI was used for all experiments (see Fig. 1 for an example of a raw tomographic slice).

2.1.4 Freezing curve and ice fraction. For part of the image analysis the ice fraction as a function of temperature was used as an input parameter. The ice fraction was calculated from a differential scanning calorimetry (DSC) curve of the ice cream mix. A small aliquot (8.9mg) of the ice cream mix was weighted into a DSC pan. The sample was then cooled from room temperature down to T = −60 °C at −5 °C min⁻¹. In order to enlarge the total amount of frozen water the sample was then partially melted during a first heating ramp with a rate of 1 °C min⁻¹ to T = −5 °C and subsequently cooled again to T = −60 °C, at a rate of −5 °C min⁻¹. During the following second heating ramp at 1 °C min⁻¹ up to T = 10 °C the heat flow was measured. Starting at T = −34 °C around the glass transition temperature of the maximally freeze concentrated solution the heat flow was integrated in order to calculate the enthalpy H of the ice cream mix as a function of temperature. An absolute value of H = 0 was arbitrarily assigned to T = −34 °C. From the position of the maximum heat flow during the heating ramp the melting temperature of the ice cream mix was deduced to be T_m = −2.94 ± 0.3 °C. Using an empiric model proposed by Riedel,¹⁸ the ice fraction was calculated as a function of temperature. The parameters of the model were set to α = 0.39, β = 0, a = 5.99, b = 0.62, c = 1.40, d = 0. The ice fraction (in volume percent) as a function of temperature is shown in Fig. 2. Throughout the paper, ice, air and unfrozen matrix fractions will be given in volume percent. As we will explain in the following sections the ice fraction is needed as an input parameter for the segmentation routine in the image analysis process.

Fig. 2 Ice fraction of the ice cream mix in volume percent as a function of temperature, with reference to the example segmentations shown in Fig. 7.

2.2 Tomography and temperature control

A Scanco μCT40 desktop device inside a cold room was used for the tomographic imaging. The tomograph was equipped with a microfocus X-ray source, operated at 75 kV acceleration voltage. The typical integration time to capture one projection with a nominal pixel size of 6 μm was 900 ms. For tomography, 1000 projections per 180 degrees were recorded.

The temperature inside the sample chamber was indirectly controlled by the lab temperature. We logged the temperature inside the μCT40 every 30s using a UTL-3 data logger (Geotest). The thermocouple of the data logger was measuring the air temperature in the sample chamber, at a distance approximately 10 cm from the ice cream sample. Since the response time scale for temperature changes was on the order of hours, we assume that this temperature reading represents the core temperature of the ice cream within an uncertainty of a few tenth of a degree.

A cylindrical sample holder (10 mm diameter) was filled with a piece of ice cream and placed inside the CT scanner (Fig. 1). The cold lab was programmed to follow a temperature step function, oscillating between −20 °C and −8 °C. The effective temperature in the CT sample chamber is shown in Fig. 3. At the same time, the CT scanner was programmed to start a scan every 4 h. Note that the CT scanning periods are also indicated in Fig. 3 as gray bars.

Fig. 3 Measured temperature response inside the sample chamber of the μCT. The gray bars indicate the periodic μCT scans that were started every 4 h. Two dark blue bars emphasize the time frames for one cold temperature and one high temperature which are referred to in Fig. 7 and 8.

2.3 Image processing

The following evaluation was based on a subvolume of 400³ voxel (or (2.4 mm)³), taken from the undisturbed central part of the sample.

The characteristic structural dimensions of ice cream were close to the resolution limit of the X-ray tomographic device (6 μm isotropic voxel size). To improve the signal quality an edge-preserving smoothing filter was applied before segmentation.¹⁹ Namely we used a three-dimensional anisotropic diffusion filter as implemented in the open source library Insight Tool Kit (ITK).²⁰ The filter requires the tuning of three parameters, which we introduce in a quick summary of the principle of the algorithm.

In anisotropic diffusion filtering, an initial image f(x, y, z, t = 0) evolves over time according to the partial differential equation describing diffusion,

with a spatially varying diffusion constant c = c(x, y, z). The key difference to isotropic diffusion (c = const.), which is equivalent to a convolution with a Gaussian kernel filter in the case of periodic boundary conditions, is that c is usually decreased at the edges in the image, such that no intensity ‘flows’ across the edges. In practice, c is thus defined as a monotonically decreasing function of the gradient magnitude, e.g.where k determines the sensitivity of the diffusion constant to the gradient magnitude: the larger the value of k, the slower c decreases with the gradient magnitude. The normalization with the root mean square of the gradient magnitude, , renders k dimensionless.

The advantage of anisotropic diffusion filters is their edge-preserving property whilst smoothing homogeneous regions, but they require correct tuning: the number of iterations N, the time step Δt, and the diffusion sensitivity k have to be chosen correctly. However, these three parameters are not independent. In the limit k → ∞, the diffusion becomes isotropic and solving the PDE at time NΔt is equivalent to a convolution with a Gaussian with . This means—given Δt is chosen small enough to limit numerical errors—that the filter depends only on two parameters, namely k and the total diffusion time N × Δt.

To optimize the choice of parameters, we considered as a measure for the image quality the contrast-to-noise ratio between two materials m₁ and m₂,

where S represents the mean signal of the object and σ is the standard deviation. Here, due to the small contrast, the interesting materials are ice crystals and liquid solution. In principle, calculating the CNR requires prior knowledge about the sample, i.e. one has to know the precise spatial distribution of the objects. In practice, we used a reference segmentation that was visually checked, and determined the CNR for different sets of parameters based on this reference segmentation. Fig. 4 shows the CNR as a function of , for various choices of k. If k is very small diffusion is suppressed at any point and the filtering has only little effect. For k ≥ 3, all curves coincide. Increasing k further does not improve the CNR. For k ≥ 1 the CNR has a maximum at the total time , further temporal evolution over-smoothes the image and reduces the CNR. A set of parameters that is compatible with above observations is

N = 11, Δt = 0.03, k = 3, (4)

which was used for all following evaluations.

Fig. 4 Contrast to noise ratio as a function of the square root of the total diffusion time . The curves collapse for k ≥ 3, and have a maximum for .

The next step towards quantitative evaluation is the classification of image voxels into the correct materials. A straightforward segmentation with two thresholds necessarily introduces artifacts when applied to a random three phase material, as illustrated in Fig. 5. Consider a boundary between highly absorbing (here: unfrozen matrix) and barely absorbing material (here: air): the continuity of the measured signal leads to intermediate gray values where the absorption changes from high to low; these intermediate values would be classified as the material with absorption coefficient between the two others (here: ice crystals) by simple threshold segmentation. As a result, in the case of ice cream, all air bubbles would be surrounded by a thin ice shell even if in reality they were embedded in unfrozen matrix (Fig. 6).

Fig. 5 Measurement artifacts and how to identify them: the true absorption signal (top) is blurred by the measurement process, which can be described as a convolution with the point spread function of the measurement device (middle panel). A thin spurious layer of material II is introduced between material I and material III when applying a simple threshold-segmentation (arrow). However, these ‘spurious voxels’ of material II can be distinguished from the real voxels of material II by taking into account the gray value gradient (bottom panel): those voxels of material II where the gradient exceeds a certain threshold (horizontal line) have to be assigned to the neighboring materials.

Fig. 6 Principle of the iterative segmentation algorithm. (left) A direct segmentation by simple thresholding falsely assigns the transition voxels between the unfrozen matrix and the air to the ice phase, as depicted in Fig. 5. This artifact results in a thin layer of ice apparently covering the air bubbles. (right) 3D-gradient magnitude within the ice phase. The bright rings around the air bubbles reveal the “spurious voxels”.

To eliminate these artifacts, a two-step, recursive segmentation algorithm was applied. In a first step, the air phase was identified directly using the clearly distinguishable dip in the gray value histogram (see fig.7b). To proceed, the expected ice fraction for the recorded temperature was calculated (see Section 2.1.4) and a tentative second threshold was chosen such that the fraction of voxels that would be in the “ice crystal” class matched the expected ice fraction. As we have argued above this class comprises many spurious, non-ice voxels because they happen to be at the interface between unfrozen matrix and air. What distinguishes these spurious ice voxels from true ice voxels is the fact that the gradient of the gray value is very high compared to normal ice voxels. Therefore, the three-dimensional gradient magnitude can be used to discard the spurious ice voxels, resulting in a decreased ice fraction that does not match the calculated ice fraction. Knowing the number of spurious ice voxels, we can decrease the tentative threshold to increase the set of ice voxels to be tested. This process is iterated until the tentative threshold does not change any more and the segmented ice fraction matches the calculated one.

Fig. 7 Result of the segmentation procedure, illustrated by two example slices, the left one from a scan taken at time t = 52.5 h, and the right one taken at t = 104.5 h. Both time points are indicated by a dark blue bar in Fig. 3. The input parameter for segmentation was the expected ice fraction at the measured temperature, i.e. 48.7% for a mean temperature of −16.0 °C and 21.9% at a mean temperature of −4.3 °C.

2.4 Quantitative analysis

To extract structural parameters, we used the scriptable Image Processing Language (IPL) from Scanco Medical. Since the CNR of ice crystals and sucrose solution was close to one, salt-and-pepper noise was present after the segmentation, which would affect any measurement of size distributions. Therefore, a component labeling procedure was performed and small, non-connected clusters (<4 voxels) were not included in the size calculation, both for the ice phase and the air phase. In addition, small inclusions within both phases that were not connected to the respective background were eliminated.

The volume fraction of each phase was determined by counting the voxels. For the size distribution an algorithm based on the distance transformation was applied.²¹ The algorithm can be visualized by inscribing maximal spheres into the structure, assigning to each voxel the diameter of the maximal sphere it belongs to. The reported mean values correspond to the volume-averaged diameter.

For a correct interpretation of the results, it is important to note that this definition of thickness is sensitive to the smallest extension of a three-dimensional object. There is no unique and straightforward definition for structure size in three dimensions in the literature. Other definitions could take into account the maximum extension in any dimension, for example. In the definition that we used, the deformation of objects in one dimension will be reflected also in the value of thickness, which is reasonable for capturing changes in the structure.

3 Results and discussion

3.1 Overview of structural evolution

30 tomography scans, separated by 4 h, were acquired while the sample was subjected to the heat shock protocol depicted in Fig. 3. The logged temperature showed that the rising flank occupied more time than the cooling, since there was no active heating applied to the cold room; instead, the warming had to take place through heat transfer from the surrounding laboratory, resulting in an undefined temperature during several hours. Those five scans which happened to fall inside a rising flank (t ∈ {0.5, 4.5, 28.5, 64.5, 100.5}) were discarded from the analysis, leaving a set of 25 tomograms that were evaluated according to the protocol described in Section 2.3.

Visual inspection of the tomographic data reveals dramatic effects of the heat shock protocol on the samples, and marked differences between warm and cold periods. To exemplify this, Fig. 7a) shows the raw data for two reconstructed slices, taken at the two time points marked in Fig. 3 (t = {52.5, 104.5}). The gray value corresponds to the absorption coefficient measured within that voxel, leading to dark air bubbles and bright liquid solution. Both slices belong to the exact same analyzed sub-volume, and give a first hint about the structure evolution over time. The left column in Fig. 7 corresponds to a sample temperature of −16.0 °C (cold room set to −20 °C), and the right column to a sample temperature of −4.3 °C (cold room set to −8 °C). Obviously, the air bubbles have undergone a significant coarsening during the heat shock cycle. The images reveal that the shape of the air bubbles is qualitatively different at cold temperatures compared to warm temperatures, a point that will be discussed further in the following paragraphs. Note the different ice fractions (dark grey areas) that can be distinguished visually from the unfrozen matrix (light grey areas) in the gray value images.

The segmentation of the air bubbles was straightforward using the clear minimum in the gray level histogram separating the air phase from the condensed phase as shown in Fig. 7b). The measured air content of the initial ice cream at t = 8.5 h was 27.5%. During the time-lapse study we observed a rise in the air content by a significant amount, up to a value of 42.5% at the end of the experiment. The air content stayed constant during the cold periods, but was rising with time when the temperature was high. There is no physical reason and it has never been observed that the air fraction in ice cream would increase with time. During heat shock treatment actually the opposite—a loss of air—is frequently observed which leads to shrinkage of the ice cream. Therefore the here observed increase of the air fraction must be an artifact of the measurement method. There are two possible reasons for this, which both have to do with the limited resolution of the used microtomograph. First, a significant part of the air fraction can be present in spaces with a dimension smaller than our resolution limit. The true resolution limit to identify structures in tomography is usually assumed to be between two and three voxel dimensions, which for our case results in a spatial resolution of around 15 μm. Second, again due to the resolution limit, a large number of voxels represent a mix of two phases e.g. air and unfrozen matrix. The minimum in the gray level histogram does not necessarily present a voxel containing 50% air and 50% unfrozen matrix. Both effects change over time due to the coarsening of the air miccrostructure. The amount of air in spaces smaller than the resolution limit as well as the total surface area (mixed voxels) between air and the other phases is reduced due to the coarsening process. This is well in line with the observation that the measured air fraction changes significantly only during the high temperature periods. We are aware that this is a limitation of the presented experiment and that this has to be taken into account for the interpretation of the air bubble sizes. However, we would like to remark that other techniques investigating air microstructure of ice cream do not give access to the air fraction at all or only in an indirect way. During ice cream production the air fraction was measured to be 50 ± 5%. Given that ice cream is not a homogeneous material on the investigated subvolume used here the air fraction towards the end of the measurements agrees well with the expected air fraction from the ice cream production.

In order to segment the ice phase we applied the iterative segmentation algorithm described in Section 2.3 to all measurements. The three resulting phases for the examples are shown in Fig. 7c), taking into account the expected ice fractions from the DSC measurements of 48.7% and 21.9%, respectively. The segmentation can be confirmed visually by comparing the gray scale images of the raw data (see Fig. 7a) to the segmentation results (see Fig. 7c). Overall, there is good correspondence, although the ice-fraction might be intuitively perceived too low in the segmented images, at least for the high temperature measurement. The ice fraction is an input parameter for the segmentation, and was determined experimentally as a function of temperature (Section 2.1.4, Fig. 2). One problem is the large uncertainty in the calculation of the ice fraction, particularly at temperatures close to the melting temperature, which is due to uncertainties in the DSC measurement, uncertainties in the used empiric model which for example neglects mixing enthalpy and enthalpic contributions from fat crystallization, and uncertainties in the measurement of the sample temperature. While the former two possibly introduce a systematic error in the function, for the present paper the accurate value for the ice fraction is not important as long as it is determined in a consistent way. Our conclusions drawn from a comparison of the microstructure of the ice phase for different time points do not depend on the exact ice fraction chosen for the two different temperatures. On the other hand, uncertainties in temperature measurements are of more concern. They affect in particular the higher temperatures, where the ice fraction changes rapidly with temperature. An error of ±0.5 °C results in an error in the ice fraction of around ±6% at T = −4.3 °C, but only an error of around ±1% at T = −16.0 °C. This emphasizes the importance of precise temperature measurements.

Based on the segmentation, three-dimensional renderings of small sub-sheets of the two example measurements are shown in Fig. 8, together with the subsequent scan to visualize short term changes. Within the cold period (Fig. 8a), both phases are very stable, and thus there is no discernible difference between two subsequent scans. While the air phase at the higher temperature is characterized by well rounded shapes, determined by minimization of the interface-energy between air and the unfrozen matrix, the shape of the air bubbles in the low-temperature image exhibits rugged, strongly deformed interfaces. This can be interpreted as mechanical deformation during ice crystal growth when cooling down the sample. Partially, this may be caused by ice crystals growing in preferred directions and partially by the density difference between water and ice. The increase of the ice fraction from 21.9% to 48.7% causes an expansion of the non air phase of around 3.5% which can only be accommodated by the compressible air phase. The rugged shape of the air bubble interfaces also does not relax during the time at cold temperature. This demonstrates that the relaxation time of the air bubble interface driven by the surface tension is much longer than our observation time. From this we conclude that Ostwald ripening, a coarsening process also driven by surface tension, does not play a significant role during cold temperature periods. In contrast, at high temperature as shown in Fig. 8b, one can observe rapid changes in both phases. For the air bubbles we observe structural relaxation and in addition coalescence (see arrow in Fig. 8b). Note that this is the first time that a bubble coalescence process has been directly observed in ice cream under realistic storage conditions. Considering the ice phase, at higher temperatures a significant amount of the ice structure melts and therefore the volume fraction of the unfrozen matrix increases, and the remaining ice objects are small and dispersed.

Fig. 8 Three-dimensional rendering of the ice phase (top) and the air phase (bottom), for subsequent time steps. Sub-figures a) and b) refer to the two selected points in the temperature cycle, which are indicated in Fig. 3. The first two measurements (a) were made at low temperature (T = −16 °C), and the last two measurements (b) were taken at high temperature (T = −4 °C), as illustrated by the insets. The temporal separation between the two pairs of scans was four hours. A site where two air bubbles coalesce in (b) was emphasized with a black arrow. For an animation of the complete time series please see the Supplementary Material online. (Dimensions of the volumes: 2.4 × 2.4 × 0.18mm³).

To quantify the changes, the difference images between two subsequent scans can be considered. Fig. 9 shows the disappearing (red) and appearing (green) voxels in the air phase between subsequent scans separated by 4 h. Fig. 9a corresponds to the cold temperature scans shown in Fig. 8a, where only a very small fraction of voxels (+4%) at the boundaries changes. For the image in Fig. 9b showing changes between two subsequent scans at high temperature, a significant fraction of 18% of the air voxels were relocated, for example at the point of coalescence.

Fig. 9 Changes in the air phase between subsequent scans, derived from the images shown in Fig. 8. Voxels marked green represent appearing air phase and voxels marked red represent disappearing air phase. (a) For low temperatures, the bubbles remain mostly stable, which is reflected by the little changes occur during the 4 h between the scans. (b) At high temperature, there are significant changes, the growth/ablation pattern results from both moving bubbles and deforming bubbles.

For a systematic study of the coarsening effects we calculated the distribution of structure thickness (see Section 2.4) for both the air and the ice. The temporal evolution of ice crystal and air bubble sizes are shown in Fig. 10. Over the 120 h of the time lapse study the ice and the air microstructure exhibited significant coarsening which can be seen in the overall increase of mean sizes as well as in the broadening of the size distributions. But there are large differences between the high and low temperature phases as well as between ice and air phase.

Fig. 10 Evolution of typical structural sizes in relation to temperature changes. In (a), the evolution of the mean thickness of ice crystals is shown, while (b) depicts the evolution of air bubble thickness. For both phases in the left column mean sizes (dark blue dots) are plotted together with temperatures (red line) in order to facilitate the distinction between warm and cold periods within the temperature cycle. The histograms behind some selected dots (marked by black lines) are shown on the right hand side. For interpretation of these graphs it is important to note that frequencies mean the fraction of voxels that belong to a certain size class as calculated by the distance transform.

3.2 Evolution of ice crystal thickness

Regarding the ice crystals the most striking feature is the difference of the mean size between low and high temperature situations. This is due to the more than two-fold increase of the ice fraction when going from a warm to a cold period in the temperature cycle (e.g. from 21.9% to 48.7% for the two representative example scans in Fig. 7 and 8). If a structure grows homogeneously in all dimensions, a two-fold volume increase would result in a size increase by a factor of 2^1/3 ≈ 1.3. In contrast to this we observed an increase of the mean ice crystal size by roughly a factor of four. Therefore, the largest part of this increase has to be explained by a connection of previously unconnected ice crystals. The opposite is the case when the temperature rises which means that the ice crystals not only shrink but also break up into a number of fragments during the partial melting. This breakup requires that the crystals have an elongated, non-spherical shape which is also what we observe particularly at low temperatures (see Fig. 7). This itself is a very important finding of the present study. So far it has been generally assumed that the number of ice crystals is a monotonously decreasing function in time. Here we have demonstrated for the first time that the partial melting phase during a heat shock can induce breakup of ice crystals in ice cream with the consequence that the number of ice crystals is increasing during this phase.

During the periods with constant (either high or low) temperature there are no strong changes of the mean size of the ice crystals within the 20 h observation time even though there seems to be a trend, e.g. in the frame from t = 68 h to t = 86 h. This does not exclude that coarsening could be observed during longer periods of constant temperature and more precise measurements. In contrast, over several cycles of temperature variation we do observe a significant coarsening of the ice crystals both comparing subsequent cold periods and subsequent warm periods which results in a net growth of ice crystals over the 120 h. Looking at the warm periods and cold periods separately we observe marked differences in the effective growth rates as can be seen from the two thin blue lines in Fig. 10a. We think that this is a startup effect and that the shapes of the ice crystal size distributions do not reach a steady state during our measurements. As significant coarsening of the mean ice crystal size is observed only over an entire temperature cycle and not during the constant temperature periods we conclude that the main coarsening effect is due to a partial melting refreezing mechanism. Since the ice fraction at a given temperature stays constant this implies that the number of ice crystals over an entire temperature cycle must decrease. This is not in contradiction with our previous observation that the number of ice crystals is fluctuating during the cycle. We would like to remark that this does not prove that a partial melting refreezing mechanism is the dominant coarsening mechanism for all situations as for example ice cream stored in a household freezer is typically subject to faster but smaller temperature fluctuations compared to our heat shock protocol.

The thickness distributions of the ice crystals detail the discussed behavior. While the distributions for t = 9 h and t = 117 h (both from a warm period) have very similar shape and width, the two distributions at t = 53 h and t = 93 h show that a significant broadening occurred during the relatively short time between the second and third cold period. This strengthens our previous idea that the shape of the ice crystal size distributions have not reached a steady state during the observation time of 120 h.

3.3 Evolution of air bubble thickness

In contrast to the behavior of ice crystals during temperature fluctuations, the air bubbles are expected to grow monotonously with time, which is indeed the case for the first 3 cycles (see Fig. 10b). Interestingly, the mean bubble size at the beginning of the third cold period at t = 88 h decreases slightly, which is probably due to deformations of the bubbles during cooling—and therefore deviations from the spherical shape which yields the maximal thickness for a given volume within the frame of the analysis using a distance transform algorithm.

Confirming the visual insight as demonstrated in Fig. 8, the mean bubble size stays almost constant during low temperature periods. A significant increase of the mean bubble size can only be observed during the high temperature periods. Regarding the high temperature periods an interesting feature of air bubble coarsening is that we observe an almost constant and linear mean size increase of 3.8 ± 0.8, 4.0 ± 0.1, and 4.0 ± 0.4 μm h⁻¹ for the three warm periods, respectively. In other words: the coarsening rate at high temperature stays constant over the entire time lapse study.

Looking at the evolution of the air bubble size distribution shown in Fig. 10b, we observe a strong broadening of the distribution in addition to the shift of the distribution maximum to larger sizes. This yields a hint to the dominant mechanism for the coarsening of the air bubbles during our heat shock protocol. An Ostwald ripening process usually leads to a narrow self-similar size distribution where the mean bubble size evolves with time proportional to t^1/3 and one would not expect significant broadening of the distribution over time.²² Both the almost linear coarsening rate at high temperature and the strong broadening of the size distribution together lead us to the hypothesis that for this experiment Ostwald ripening is not a dominant coarsening mechanism and that coalescence represents the main coarsening mechanism for the air microstructure. On the other hand there are at least two arguments against this hypothesis. First, our bubble size distributions might not represent a steady state distribution with respect to Ostwald ripening as our time lapse study only covers a limited time interval. Second, as noted in the beginning of the results section we might miss a significant amount of bubbles smaller than our resolution limit. However we want to remark that the latter argument leads to an even faster coarsening of the air microstructure as is shown in Fig. 10b. This is evident from the increase of the observed air fraction with time which implies that less air is present in very small bubbles at later times which would result in a faster increase of the mean bubble size.

We mentioned above that the air fraction increased during the experiment, from roughly 30% to 40%, and hypothesized this to be due to unresolved air bubbles at early stages. As more and more air is resolvable, this has also an influence on the size distribution. Since the very small bubbles can coalesce with bubbles of any size, they contribute to the coarsening at any scale. For the mean value, this means that we probably over-estimate the bubble size at earlier times, since a fraction of very small bubbles remains unresolved. Given the narrow distribution at early times, this effect is not expected to be significant. The general observations remain valid, and experiments with higher spatial resolution will reveal the impact of unresolved bubbles.

The data we have presented show that the evolution of the morphology can be tracked nondestructively with X-ray tomography. It yields interesting and novel insight into the complex interplay between ice crystals and air bubbles. A systematic study of the structural evolution will provide the basis for understanding the physics behind the coarsening of the different phases.

4 Conclusions

We developed a method to visualize nondestructively the three main phases of ice cream: air, water ice, and unfrozen matrix. The advantage of this method is the inherent three-dimensional character, which gives a deep insight into the structural relationships within the individual phases as well as amongst those phases. Compared to the widely used and established techniques, it is clear that the minimal sample preparation results in a reduction of artifacts. On the other hand, the limited resolution of desktop microtomographs introduces systematic errors as we have shown with the evolution of the total air fraction. This will render it necessary to complement the information e.g. with data from microscopic techniques such as confocal microscopy or fluorescence microscopy which will give a higher resolution and a better idea on the distribution of specific ingredients. This is, however, not a principal limitation of the method, since the resolution of X-ray tomography can be pushed down to sub-micrometer scales by using synchrotron sources.

We solved a general problem that appears in the threshold-based segmentation of three-phase materials, namely the ‘shell’ of falsely assigned voxels of intermediate gray values that appear at the interfaces between high and low gray values. The procedure can readily be applied to similar problems in a variety of applications.

First results show the evolution of the microstructure of an ice cream subjected to a heat shock protocol. The non-destructive measurements allow for a quantification of structural parameters of one and the same subregion of the sample, which gives new insight into the coarsening process of ice cream in spite the discussed difficulties due to the resolution limit. We observed that particularly during cold periods elongated ice crystals can form. From the large difference between the ice crystal mean size at low and high temperatures we could conclude that during partial freezing ice crystals grow together and during partial melting large elongated ice crystals again split up into smaller ice crystals. This implies that the number of ice crystals is fluctuating during a temperature cycle. Since significant ice crystal coarsening occurs only over an entire temperature cycle we could identify a partial melting refreezing mechanism as the dominant coarsening mechanism for the investigated storage condition.

In contrast to the coarsening of the ice crystals, where the amplitude of the temperature fluctuation is the most important parameter, for the coarsening of the air microstructure the time spent at high temperature periods is the most important factor. We observe significant coarsening only during the high temperature periods. The observed coarsening behavior, both the coarsening rate and the evolution of the size distribution, hint on coalescence as the main coarsening mechanism at high temperatures. At low temperatures air microstructure coarsening is strongly suppressed due to the high viscosity of the bulk phase between the air.

For future studies of microstructure of ice cream and other frozen foods we believe that a better resolution and a better contrast between the different phases will overcome major limitations of the present study and will be able to deliver a concise picture of the microstructure evolution. Both can be achieved by using X-ray tomography at a synchrotron source. A less strong heat shock protocol and reduced scan time will in addition allow to follow individual ice crystals and air bubbles and will help to reveal the mechanisms behind the coarsening. Understanding the coarsening process will ultimately open up the possibility to specifically modify production or distribution processes to optimize product quality.

References

1. Y. Chang and R. Hartel, Int. Dairy J., 2002, 12, 463–472.
2. K. G. Berger, B. K. Bullimore, G. W. White and W. B. Wright, Dairy Ind., 1972, 37, 419.
3. K. G. Berger, B. K. Bullimore, G. W. White and W. B. Wright, dairy Ind., 1972, 37, 493.
4. D. Reid, in Frozen Food Technology, ed. C. P. Mallet, Chapman & Hall, London, 1993, ch. Basic physical phenomena in the freezing and thawing of plant and animal tissues, pp. 1–19.
5. W. S. Arbuckle, Ice Cream, AVI Publishing Company Inc., Westport, CT, 1986.
6. C. Clarke, Food Science and Technology, 2006, 20, 31–33.
7. K. Cook and R. Hartel, Compr. Rev. Food Sci. Food Saf., 2010, 9, 213–222.
8. D. P. Donhowe, R. W. Hartel and R. L. B. Jr, J. Dairy Sci., 1991, 74, 3334–3344.
9. A. Caillet, C. Cogn, J. Andrieu, P. Laurent and A. Rivoire, LWT–Food Sci. Technol., 2003, 36, 743–749.
10. S. Ueno, G.-S. Do, Y. Sagara, K.-I. Kudoh and T. Higuchi, Int. J. Refrig., 2004, 27, 302–308.
11. W. S. Arbuckle, Ice Cream Trade Journal, 1960, 56, 62–68–168–172.
12. A. Smith, H. Goff and A. Regand, in Ice Cream II, ed. H. D. Goff and B. W. Tharp, International dairy Federation, Brussels, 2004, ch. Advanced microscopy techniques for ice cream research, pp. 197–208.
13. G. Petzold and J. Aguilera, Food Biophys., 2009, 4, 378–396.
14. R. Mahdjoub, P. Chouvenc, M. Seurin, J. Andrieu and A. Briguet, Carbohydr. Res., 2006, 341, 492–498.
15. P. Falcone, A. Baiano, F. Zanini, L. Mancini, G. Tromba, D. Dreossi, F. Montanari, N. Scuor and M. Del Nobile, J. Food Sci., 2005, 70, E265–E272.
16. R. Mousavi, T. Miri, P. Cox and P. Fryer, Int. J. Food Sci. Technol., 2007, 42, 714–727.
17. R. Mousavi, T. Miri, P. Cox and P. Fryer, J. Food Sci., 2005, 70, E437–E442.
18. L. Riedel, Chem. Mikrobiol. Technol. Lebensm., 1978, 5, 129–133.
19. A. Kaestner, E. Lehmann and M. Stampanoni, Adv. Water Resour., 2008, 31, 1174–1187.
20. T. S. Yoo, Insight into Images: Principles and Practice for Segmentation, Registration, and Image Analysis, AK Peters, 1st edn, 2004.
21. T. Hildebrand and P. Ruegsegger, J. Microsc., 1997, 185, 67–75.
22. I. Lifshitz and V. Slyozov, J. Phys. Chem. Solids, 1961, 19, 35–50.

---

Footnotes

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c2sm00034b

‡ Present address: Swiss Light Source, Paul Scherrer Institut, Villigen, Switzerland; E-mail: bernd.pinzer@psi.ch

---