Microscopic deformations in MDF swelling: a unique 4D-CT characterization

Medium-density fiberboard (MDF), a wood-based material that consists of a tight random network of wood fibers, deforms more than wood when exposed to water. For the first time, the microscopic deformations of MDF were tracked during swelling. A hygroscopic swelling setup imposing the material to deform throughout tomographic acquisition was used coupled to X-ray microtomography. An advanced reconstruction algorithm enabled reconstruction of images free of motion artefacts, and state-of-the-art digital volume correlation was applied to determine the mechanical strain fields at high resolution. Wood fiber bundles were then segmented from single fibers with deep learning using the UNet3D architecture. Combined with the strain fields, this segmentation showed that wood fiber bundles were the drivers of MDF swelling. This contrasts with the hygroscopic behavior of wood, where structured wood swells less than single fibers, which might be caused by a difference in penetration and distribution of the adhesive, in and on the wood fiber cell wall. The unique characterization of MDF’s dynamic behavior can already be used to develop manufacturing strategies to improve water resistance, therefore widening the uses of natural fiber-based materials.

production of about 410 million m 3 per year [2]. However, distinctive fragilities limit further potential application of engineered wood. For instance, wood is well known to swell while it absorbs water [3], a principle that was already used some 5500 years ago by the ancient Egyptians to fracture rocks [4]. But medium-density fiberboard (MDF: $ 100 million m 3 produced/year [2]) is known to swell more than solid wood in the panel thickness direction. This specific disadvantage, which does not erase the other benefits that MDF has over solid wood, makes it so that MDF is still rarely used in humid environments, with the exception of highly-modified and costlier products, e.g. Tricoya.
MDF consists of a mat of separated wood fibers randomly arranged horizontally, [...], pressed with a combination of heat (180-210 C) and pressure (0.5-5 MPa) [5]. The structure is held together by a thermosetting adhesive resin, usually urea-formaldehyde, which is present in relatively low amounts ( $ 10 %wt) and distributed on the loose fibers before compression. In addition to the adhesive, hydrophobic agents (e.g. paraffin) are also mixed with the wood fibers. Nonetheless, MDF remains prone to swelling due to water absorption from both liquid water and water vapor. The irreversible swelling affects the aesthetics and mechanical properties of products made with MDF (typically furniture and flooring), resulting in a reduced service life under moist conditions. That is problematic, since MDF is still poorly recycled [6], although research is being done in that respect.
MDF has already been the focus of studies to improve its resistance to humidity. The structuring element of MDF is the wood fiber, of which the morphometry is expected to primarily condition the panels' mechanical behavior. The fibers' morphometry is affected by the heavy industrial processing involved in the making of MDF [7]. Therefore, fibers need to be studied while already embedded in MDF panels. X-ray micro-tomography (lCT) is a technique that uses the penetrating ability of X-rays to reconstruct 3D images of a sample's internal structure. Since lCT can non-destructively resolve features down to the lm, which is the scale of the wood fibers, it is an excellent tool to characterize MDF. With lCT, the diameter of wood fibers embedded in MDF was found to be distributed around 20 lm [8]. The proportion of residual wood fiber bundles, results of an imperfect processing into single fibers, could be evaluated as well [9].
Together, these studies produced inputs for eventual structural models of MDF, which could in turn be used in simulations. These simulations could then lead to new manufacturing processes that would result in humidity-resistant MDF. However, no study was yet able to embody swelling simulations [10][11][12]. This could be due to computational limitations, as simulating a dense random network of wood fibers (each with various shapes and dimensions) demands a large number of computations. Additionally, it could be due to the lack of thorough characterization of the swelling itself, as to achieve and verify good-quality simulation dynamic characterization of the swelling in the three spatial dimensions is key. While lCT can be used to observe time-dependent phenomenons, to the best of our knowledge such characterizations have never been done for MDF. A reason could be the complexity associated with capturing the relatively fast dynamics of swelling MDF, a material that is already challenging for static observations.
In lCT, higher spatial resolution images require longer exposure times, in the order of a second per projection. This leads to scanning times between tens of minutes and hours to retrieve a good quality 3D image, something that impairs the observation of fast processes. Under specific conditions, faster scanning can be achieved [13], yet these conditions do not apply to MDF. The remaining option, acquiring long exposure high-resolution images of MDF while it is deforming, resulted so far in blurry reconstructed 3D images. However, the latest developments in reconstruction algorithms include motion compensation [14,15], which finally allows to reconstruct sharp, high quality images despite the material deformation.
In this work, we present an unprecedented 4dimensional observation of MDF swelling at high resolution. A custom-developed experimental setup mounted on the rotation stage of the lCT scanner allows the sample to undergo hygroscopic swelling during the acquisition of a high resolution lCT dataset. Sharp reconstructions of the MDF, obtained using motion compensated CT reconstruction to compensate for the global swelling motion, were used as input for regularized digital volume correlation (DVC) and deep learning segmentation. This analysis revealed that these segmented areas correspond to distinct zones of preferential swelling, which constitute a basis for further investigations of MDF.

Medium-density fiberboard samples
Four replicate samples, size 3 Â 3 Â 9 mm 3 , were extracted from a single commercial hydrophobicgrade 9 mm Umidax MDF panel (3 samples for timelapse photography and 1 for 4D lCT). The panel itself was provided by the Belgian MDF manufacturer Unilin, and chosen considering that 9 mm is the median thickness of MDF panels produced in Europe [16]. For representativeness, the samples were uncut in the thickness direction. In the other two directions, the sample dimensions were chosen so that at a lCT voxel size of 2.5 lm, necessary for the study of wood fibers [17], the field of view would encompass the full sample width. Taking advantage of the characteristic symmetry of MDF panels in that direction relative to the center of the panel [8], only half of the sample height was scanned with lCT. Cuboid samples were prepared using a table saw, which allowed us to cut such small samples without apparent damage to the remaining, brittle structure. The samples were kept in the experiment room, at ambient temperature and relative humidity, from at least a night before the start of the experiments.

Experimental setup
A lCT setup was designed to enable swelling during acquisition, as presented in Fig. 1. The sample was set in a climate chamber in which it could rotate close to the X-ray tube. The climate chamber had a moist air inlet, and an outlet for passive exhaust. A commercial relative humidity generator GenRH (Surface Measurement Systems Ltd., UK) provided the moist air input, at 95% relative humidity (RH), and ambient temperature.

Climate chamber
A climate chamber was designed to fit a single sample with size of a few millimeters at most [18]. Provided a moist air input, the climate chamber maintained an atmosphere of set relative humidity and temperature. The air flow around the sample in that climate chamber is laminar, minimizing unwanted sample motion [19]. This is important for lCT, as displacements should be kept below the voxel pitch, which is only a few micrometers. Inside the chamber, the samples' rotation axis is as close as possible to the X-ray source to reach maximal geometrical magnification. Two polyimide windows were installed on either side of the climate chamber, to allow the X-rays to pass through with minimal attenuation. For X-rays, wood is indeed a relatively low absorbing material, and retrieving a high quality signal requires to prevent high X-ray absorption in the field of view. A custom mounting system was designed, which allows the climate chamber to hang from the top of an X-ray tube. As such, the only moving part during the experiments was the MDF sample that rotated.

Relative humidity generator: GenRH
The GenRH creates a moist air input. At the outset of the device's tube, there was a RH and temperature probe (see Fig. 1). This probe (measuring RH 0 ) allows for closed-loop control of the GenRH over the air flow it produced. The closed-loop was set with a tolerance of 1% in RH, while the temperature was left ambient. Additionally, an independent RH and temperature probe (measuring RH 1 ) was added in the climate chamber for validation. This standard DHT22 sensor (Adafruit industries, USA) was placed inside the climate chamber, at the passive air outlet. Readings are logged on a Raspberry Pi 3 at 1 Hz. Finally, the GenRH's RH control was set as a single isothermal step at 95% RH, lasting 5 h.

Photography experiment
A photographic experiment separate from the CT scanner was first performed to characterize the sample's global swelling. An EOS 6D DSLR camera (Canon, Japan, with an EF 100 mm macro focus lens) was used in place of the X-ray detector, and paired with a time-lapse controller (Captur module -timer, Hähnel Industries Ltd., Ireland), and no X-rays were used. One standard photograph was taken every 30 s. This experiment, replicated for 3 similar samples, resulted in 1000 photographs of each sample. It described their evolution for a little over 5 hours, during which the GenRH actively provided air at 95% RH.
These photographs were then processed by a Python code developed in-house using scikit-image [20], illustrated in Fig. 2. After a crude segmentation of the sample's edges, Hough lines [21] were fitted to the edges, filtered for verticality (where the vertical direction is taken with respect to the horizontal plane, or equivalently, to the panel surfaces), and merged under a proximity condition: when the distance between the two lines was less than 30 pixels, and the angle between them less than 30 degrees. This eventually resulted in two (left and right) meta Hough lines per photograph (shown in Fig. 2), of which the length was computed. The maximum of these two values was taken as a measure of the sample's vertical thickness in a given picture. Finally, the median height among the 3 samples was taken at each time point, and divided by the initial sample height to retrieve a global vertical strain e z . The resulting swelling evolution is presented in Fig. 3.

lCT experiment
The photographic experiment was replicated using the lCT system Nanowood [22], built at the Centre for Xray tomography of Ghent University (UGCT, https:// www.ugct.ugent.be). The scans were performed using the Hamamatsu L9181-2 source operated at a tube voltage of 50 kV, and a target power of 5 W. For each scan, 2001 projections were acquired over a 360 degrees rotation. Each projection required 1 s of exposure time. As a result, one single scan required about 36 minutes of scanning time. As presented in Fig. 3, one scan was acquired prior to the swelling experiment. Then, 5 scans were acquired successively, while the GenRH provided an active input to the climate chamber. Finally, a last scan was acquired at 1 h 30 min after scan V. This was after the 4 h 15 min mark overall, i.e. when the sample was expected to be stable enough for standard reconstruction. These scanning parameters enabled the so-called ''scan I'' to capture the extent of the swelling onset (as seen on Fig. 3), and required no compromise on scan quality. The swelling onset corresponds to the first instants (a little more than 20 min) after the humidity reached a Fig. 1 Schematic representation of the experimental setup, in line with the X-ray micro-tomography system (represented by the system X-ray source-sample-X-ray detector or ''camera'').
A picture of one of the replicate samples, seen through the climate chamber's polyimide window is seen on the right Fig. 2 Illustration of the major image processing steps for the photography experiment. At the last step, two meta Hough lines (red and blue) are computed to follow the lateral edges of the samples, and give the value for its thickness peak inside the climate chamber, when the sample swells the fastest. Volumes were initially reconstructed with the Octopus Reconstruction software [23], licensed by Tescan XRE (https://www. XRE.be, part of the Tescan Orsay Holding a.s.), and the standard FDK reconstruction algorithm [24]. The reconstructed voxel size was approximately 2.5 lm 3 , and 16-bit images were obtained. However, due to the sample swelling during tomographic acquisition, motion artifacts were present requiring the use of motion-compensated reconstruction, particularly at the swelling onset.

Motion-compensated reconstruction
Due to the non-instantaneous nature of CT acquisition and its subsequent numerical reconstruction, motion actually leads to intricate image artefacts that do not reflect the sample deformation. Since fibre structures are backprojected strongest when they align with the optical axis of the system, the local orientation of the fibre structures determines which projection time represents them most. The predominantly vertical motion of the sample in this sample, where 360 rotations were used for circular cone beam reconstructions, results in a line doubling as the structures align twice at different heights. This may appear as sharper detail in the reconstruction, although the sample contains no actual corresponding structures. This also disturbs the DVC algorithm, and prevents it from converging to a high-quality solution. By using MCR, motion artifacts are mitigated and the reconstructed position of the structure does not depend on the orientation, allowing for better-quality deformation analysis through DVC. For MCR, the preliminary reconstructions were produced using CTrex (1 iteration of the SART algorithm [25], relaxation factor of 0.6), a CT reconstruction framework developed at UGCT, which also contains methods for motion compensation [14]. The mean vertical slices of all reconstructions were registered using ImageJ's plugin register virtual stack slices (2D?T). This technique registers the dynamic slices by matching corresponding SIFT features [26]. The resulting affine transformation matrices were interpolated for all projection times to yield a steady affine motion [27]. The interpolation assumes that the bottom face of the sample, which is fixed to the sample holder and can thus not be affected by the swelling, is fixed in space. Finally, a second reconstruction set was produced by compensating (warping) for this interpolated motion during the reconstruction process, referencing a volume's state to the initial scan. The warped volumes, thus obtained by MCR, contain much less motion artefacts.

Digital volume correlation
DVC, based on the theory of DIC [28] but developed for 3D datasets, is a technique to estimate a displacement field, from which a strain field can in turn be estimated. The strain is a description of the material  deformation caused by an external event (e.g. an external load, a change in temperature, a chemical reaction, etc.). The study of the strain is crucial in material science, because it describes the intrinsic material response to perturbation. In the context of this study, the computation of the strain allows to investigate the mechanisms underlying MDF swelling at the microscopic scale (e.g. how much it swells, preferential swelling areas/structures, etc.). To enable this investigation, the strain field was computed using an in-house code [29]. A Gaussian pyramid with 4 levels was used. The deformation was modelled by using cubic B-splines transformation with node spacing of 8 voxels at the finest pyramid level (x2 at each higher pyramid level). The cost function consists of two terms: the similarity metric, namely sum of squared differences (SSD), and the total variation strain regularization term.
While volumes corresponding to the full field of view of our system were motion compensated, this DVC step and subsequent analyses were performed on a region of interest (ROI) manually selected in the bulk of the sample due to memory requirements. This ROI, which remains sufficiently large to represent the MDF fibrous structure, was taken in the inside at necessary distance ( $ 1 mm) from the sample's cut edges as well as the panel's surface. This was done to ensure a representative analysis for the bulk of the material. Near this surface, the density is indeed highest and wood fibers are almost completely collapsed, which is not representative for the bulk of the material, which is in the end what determines behavior. Moreover, at 2.5 lm resolution it is hard to discern any structure in this boundary area.

Fiber bundle segmentation
The strain fields were compared with a binary segmentation of the MDF structure, which isolated ''fiber bundles'' from the rest of the ROI (i.e. single wood fibers). These bundles (an example is shown in Fig. 5) are ''array''-like structures of wood fibers, partly retaining the properties of solid wood. During the manufacturing process of MDF, larger wood pieces are thermo-mechanically pulped to retrieve individual wood fibers. However, this processing is not perfect, and some fiber bundles can be found in MDF. Because such bundles are expected to behave differently than single fibers, we worked on algorithms to identify these structures in lCT scans of MDF. An optimal deep learning segmentation method was identified using the UNet3D architecture [30]. All the characterization details and different algorithms can be found in [31]. This UNet3D segmentation was achieved in Dragonfly 2020.1 (Object Research Systems Inc, Montreal, Canada -free of charge for noncommercial use). The network's first layer counted 32 filters, and it was organized in 5 levels. It weighed 3 Á 10 8 nodes, and expected input patches of 128 Â 128 Â 5 voxels 3 , where the smallest dimension was taken along either horizontal direction (in the plane of the MDF panel). The network was trained in 200 epochs, on just 9 manually-annotated slices.

Evaluation of the experimental setup
From Fig. 3, we observed that the precision with which conditions could be imposed in the climate chamber fluctuated, and could have deviated from the target 95% RH by maximum 10%. Since the second sensor DHT22 was at the passive outlet of the climate chamber, the reading RH 1 might have been degraded by the proximity to a non-climatized environment. On the other hand, the air flow from the GenRH (200 sccm -standard cubic centimetres per minute, associated to RH 0 ) might have been too little to prevent air from entering the climate chamber at the passive exhaust. We thus expect the real relative humidity value around the sample to be in between the two sensor readings (i.e. around 90% RH). However, the precision on environmental conditions is not determinant in the present study, as its focus was to acquire and analyze for the first time a 4D CT dataset of MDF swelling rapidly. As can be seen in Fig. 3, the sample was still swelling significantly. In future and for experimental repeatability, the setup could be improved by replacing the GenRH with a moist-air generator that could also input the DHT22 RH readings to its closed-loop control. A machine with the ability to actively cool the air flow produced (the GenRH is only able to heat) would also result in a better accuracy on the set RH. With the setup described here, a comparison between the original and swollen states of a MDF sample is shown in Fig. 4.

Photography experiment: surface deformation
The swelling profile captured in photography is presented in Fig. 3. After an initial phase of rapid swelling lasting about 20 minutes, the sample slowly swelled for the remaining of the experiment. The swelling regime can be modelled as a piecewise linear function of time. During the first 20 minutes, the sample swelled at a rate of 2.7%/h. In the second phase, it swelled ten times slower, at 0.27%/h. The swelling onset is therefore expected to contain the most severe motion artefacts in the tomography experiment. Yet, the swelling onset is also expected to yield the most valuable input for material investigation. Arguably, the first instants of swelling would condition the rest of the swelling experiment, something that could only be verified using 4D CT. Remark that the piecewise linear approximation is indicative, but a logarithmic model -assuming an asymptotic swelling value, or a square root model [3] seem otherwise more appropriate. Because swelling was already slowing down before the 5 h mark, an asymptotic plateau was expected in the hygroscopic swelling profile. 5 h was therefore considered an appropriate length for the experiment, considering that the last lCT scan acquired would by then experience little motion, and therefore exhibit minimal motion artifact.

lCT experiment: 3D internal deformation
The setup presented in Fig. 1 was used to acquire, in combination with the scanning schedule shown in Fig. 3, 14007 projections ( $ 70 GB of data). The first 12006 projections were almost perfectly consecutive. Taken in separate batches of 2001 projections, 6 full scans could be independently reconstructed over 360 degrees rotations. In a lCT image, the grey value represents the linear X-ray attenuation intensity of a given voxel. Here, a voxel's grey value mainly accounts for the photoelectric and Compton effects, i.e. is a function of (principally) atomic number and electron density, respectively [32]. More precisely, given the relatively low difference in elemental composition (mainly C, H, and O in similar proportions throughout the material), the intensity would be primarily affected by density. Therefore, the produced grey scale images, such as those shown in Fig. 5, can be seen as an approximate volumetric map of local density. In stable conditions, such a high-resolution map of density allows to observe material structures in detail. Notice how standard ''filtered-back projection'' (in reality, FDK for circular cone beam geometry) reconstructs scans I to V as blurry due to the swelling during acquisition. Due to blur, microscopic features such as individual wood fibers were not identifiable (Fig. 5). Nevertheless, meso-scale features (such as cracks, or high-density intrusions) remained discernable.
To assess the effect of MCR on a scan that was less impacted by motion artifacts than scan I at the swelling onset, Fig. 6 shows the effect of MCR on scan V. While the microscopic structure of MDF can then be observed before MCR, hinting that the sample is indeed relatively stable then, the algorithm undoubtedly improves image sharpness and thus the quality of microscopic feature observation.

Motion-compensated reconstruction
To retrieve microscopic sharpness in the dynamic scans, vertical slices from all separate volumes were registered. The corresponding (interpolated) vertical strain is shown in Fig. 3. It is in good accordance with the macroscopic strain that had been observed with photography. However, towards the end of the 5 h experiment, the sample swelling measured with tomography was not stabilizing. This was not clear until the data could be processed, i.e. months after the experiment. Sample variability (given only one sample was observed with tomography) could be the cause, but it could also result from the inclination of the sample in photography. As photography can only acquire surface information of the sample, its reading may deviate from a volumetric measurement when the sample tilts about the vertical axis. Because the more the sample swells the more it is likely to tilt, this effect could be more pronounced at the end of the experiment. Therefore, the lCT input is expected to yield the true swelling values.
Using the mesoscopic strain (the MCR curve in Fig. 3), we were able to warp the volume into a time series with improved spatial resolution. The importance of motion compensation can be judged in Fig. 5, where the same vertical slice through scan I is shown before and after MCR. Fiber bundles at the top center or bottom left for instance, were blurred before MCR, but are easily identifiable after. At the top of the sample, cell wall doubling can be seen, which is corrected after MCR. While the improvement in image quality is clear, Fig. 5 cannot explicitly show that DVC algorithms failed to converge before MCR.  However, in this study MCR was absolutely necessary for that reason.

Swelling strain and fiber bundles
After MCR, the vertical strain field was computed inside a representative ROI. The strains corresponding to the largest deformation are first compared, i.e. the strain field computed between scans 0 and VI. In this work, all DVC strain fields were computed with respect to the pre-experiment scan 0. Therefore, the notation e VI z for instance, refers to the strain field between scans 0 and VI. As can be seen in Fig. 7 the high strains seem to coincide with fiber bundles (more structured wood fibers) present in the MDF. These bundles are commonly found in MDF and are residues from the mechanical processing. Several works already endeavored to segment them in lCT scans of MDF [9,12], which led us to pursue these investigations with advanced algorithms. In Fig. 7, the automatic (deep learning) segmentation of these fiber bundles is shown in yellow, which highlights the bundle boundaries. In this ROI, 42 fiber bundles were identified, and accounted for 10% of the ROI's volume.
The strain values distribution, in and out of fiber bundles is also shown in Fig. 7. The e VI z strain values are masked by the automatically-segmented fiber bundle phase and arranged into a histogram. The strains are noticeably higher inside fiber bundles than in the rest of the material, which is a surprise. Indeed, solid wood is expected to swell less than single wood fibers. This is due to the middle lamella (a hydrophobic structure surrounding wood fibers) that constrains the fiber's swelling, and is generally lost during thermo-mechanical pulping [33]. The wettability of wood fibers is primordial for the good penetration of adhesive resin. Thus we hypothesize that in fiber bundles, the middle lamella constitutes a barrier to the penetration of adhesive resin in the limited mixing time during the production of MDF. This leaves unhardened wood fibers in the middle of fiber bundles, which are eventually more prone to swelling than fibers well penetrated by additives, when exposed to water for longer.

Evolution through time
Having acquired an almost-continuous time series of the material deformation, the same DVC algorithm was applied sequentially to all motion-compensated scans. The swelling dynamics are investigated: the strain evolution through time is shown in Fig. 8. The cumulative strain is shown in absolute values, using the dynamic range pertinent to e VI z . Naturally, the strains before scan VI are more homogeneous, and it is harder to discern swelling areas. However, this visualization allows to compare in absolute values, the swelling state during e.g. scan I, and scan VI (the final state). The strain values at the time step i were also normalized by their mean valueê i z . This visualization allows to identify which MDF structures are swelling at an instant i. Peaks in these normalized values are therefore indicative of swelling onset.
Moreover, a three dimensional rendering of the higher strain values (thresholded above 1.035) is shown in Fig. 9, alongside a view of the largest fiber bundles in the ROI.
Finally, the evolution of the strains through time is shown in Fig. 10. It reveals that fiber bundles already swelled more at the swelling onset. Through time, this distinction between the swelling behavior of bundles and of the rest is either maintained, or increased. On Fig. 10, it can be seen that the rate of swelling (i.e. the Fig. 7 A e VI z strain overlaid with the MDF structure (greyscale) on a representative slice. The contours of automaticallysegmented fiber bundles are drawn in yellow. B e VI z strain values grouped inside and outside of the fiber bundles, showing the difference in relative swelling slope of the mean strain curve) is always higher in the fiber bundles. As such, the fiber bundles are understood to truly initiate swelling locally. The swelling can then ''propagate'' to the other fiber structures inside MDF, which is remarkable.

Conclusion
Wood-based panels such as MDF are part of a growing industry, important for generalizing the use and manufacturing of eco-based building materials. To improve these materials, dynamic characterizations are required. In the present study, we extracted millimetric samples of commercial MDF panels to analyze with X-ray micro-tomography and a custommade climate chamber. The outcome of these experiments is a first-ever 4D dataset that encodes the evolution of the samples while swelling. Exploiting this dataset required the use of state-of-the-art algorithms in motion-compensated reconstruction, digital volume correlation and deep learning segmentation. Yet this unique analysis revealed that contrary to expectation, the swelling of wood fiber bundles present in MDF precede the swelling of single fibers. Due to the importance of the swelling onset to the macroscopic swelling of MDF panels, this result is of immediate interest to fiberboard manufacturers the cumulative strain was normalized by the mean strain valuê e z , to give a qualitative insight through time. However, the variability of this second measure decreases at higher mean strain values, so its information is less sensitive for the later scans. The absolute timings of each scan can be seen on   wishing to mitigate the water-induced swelling of MDF. Moreover, the research can now be continued to further investigate the core differences between wood fibers and wood fiber bundles embedded in MDF.
Our approach presents both novelty and relevance to the study of wood-based materials. We hope that such work can also be continued with the characterization of more green materials, in environmental conditions even more challenging (e.g. involving liquid water absorption, to conform with industrial norms such as EN 317). Furthermore, the application of MCR to one of the most complex case studies may convince the community about the feasibility of more experiments for the benefit of material science.