# Comparison of various internal fibre orientation determination methods for predicting the mechanical properties of one Douglas-fir timber board

Helene Penvern, Guillaume Pot, Jeffrey Viguier

<figure><img src="image_1.png" /></figure>

## To cite this version:

Helene Penvern, Guillaume Pot, Jeffrey Viguier. Comparison of various internal fibre orientation determination methods for predicting the mechanical properties of one Douglas-fir timber board. 2025. hal-05265237

HAL Id: hal-05265237

https://hal.science/hal-05265237v1

Preprint submitted on 23 Sep 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Comparison of various internal fibre orientation
determination methods for predicting the
mechanical properties of one Douglas-fir timber
board

Helene Penvern¹, Guillaume Pot¹, Jeffrey Viguier¹

¹Arts et Metiers Sciences et Technologies, LaBoMaP, 13 rue Porte de Paris, Cluny, F-71250, France.

Contributing authors: helene.penvern@ensam.eu;

Abstract

Technologies utilising laser dot scanning to measure surface fibre orientation have evolved into some of the most advanced system-ready solutions for machine strength grading. These methods have significantly enhanced predictive models of timber strength; however, they typically rely on surface-based interpolations that may not fully reflect internal fibre architecture.

This study proposes an adaptation of the method by Penvern et al (2023), employing gradient structure tensor (GST) analysis to determine the normal direction of growth layers with improved computational efficiency. By comparing this approach with existing internal fibre orientation determination (IFOD) techniques—including those integrating laser dot measurements and GST-based analysis—their capacity to predict the mechanical behaviour of sawn timber is assessed.

The evaluation is grounded in local bending stiffness profiles derived from digital image correlation (DIC) measurements and finite element (FE) simulations. Our findings show that precise characterization of internal fibre orientation, especially when informed by laser-based measurements, leads to genuinely accurate predictions of timber mechanical performance. On another note, GST analysis proves to be an efficient tool for determining wood's local radial direction and offers a promising approach for timber strength grading based on tomography.

Ultimately, this study aims not only to extract meaningful data on internal fibre architecture for modelling purposes, but also to demonstrate that fibre orientation is the critical determinant of mechanical behaviour in timber containing knots. The proposed method thus lays the groundwork for generating
reliable internal fibre orientation data and for developing efficient, data-driven
models that infer internal orientation from accessible surface or tomographic mea-
surements—thereby enabling more accurate predictions of timber strength and
stiffness.

**Keywords:** grain direction, tomography, GST, growth rings, stiffness
# List of Symbols

($\vec{x}, \vec{y}, \vec{z}$) Coordinate system of the wood specimen: $\vec{x}$ denotes the board's longitudinal direction, $\vec{y}$ is the width-wise direction, and $\vec{z}$ is the thickness-wise (planing) direction.

$L, R, T$ longitudinal, radial and tangential wood local directions

$\alpha$ in-plane fibre angle, defined as the angle between the major axis of the elliptical laser dot scatter and $\vec{x}$ longitudinal direction of the specimen

$\gamma$ out-of-plane "diving" angle, i.e. the fibre angle relative to the $xy$ planing plane

$\vec{u}$ $xy$ projected fibre direction in relation to $\alpha$

$\vec{v}$ 3D fibre direction

$\vec{n}$ normal direction of the growth-layers

$\epsilon_{xx,DIC}$ longitudinal strain observed on the board surface using DIC

$\epsilon_{xx,s,DIC}$ same longitudinal strain smoothed by a convolutions filter of dimension $20 \times 20 \text{ mm}^2$

$\epsilon_{xx,FEM}$ longitudinal strain obtained in FEM on corresponding external surface

$\kappa$ local $y$-displacement curvature

$E_{b,bench}$ localised bending modulus of elasticity computed from board inertia and bending bench local deflection

$E_{b,\kappa,w_x,DIC}$ localised bending modulus of elasticity computed from curvature over a longitudinal window of $w_x$ measured on the board surface using DIC

$E_{b,e,w_x,DIC}$ localised bending modulus of elasticity computed from strain measured on the board surface using DIC over a longitudinal window of $w_x$

$E_{b,\kappa,w_x,FEM}$ localised bending modulus of elasticity computed from curvature simulated using FEM on the corresponding surface over a longitudinal window of $w_x$

$E_{b,e,w_x,FEM}$ localised bending modulus of elasticity computed from strain simulated using FEM on the corresponding surface over a longitudinal window of $w_x$

# 1 Introduction

## 1.1 Background

The accurate prediction of the mechanical behaviour of sawn timber is a critical aspect of wood engineering, with significant implications for both structural applications and material science. Wood as a building material has seen increasing demand in recent years, driven by its aesthetic appeal and its status as a carbon-negative material. Structural wood strength-to-weight ratio can rival that of steel. However, unlike steel, wood mechanical properties are very variable and difficult to assess, leading to conservative safety standards that result in costly designs. A deeper understanding of sawn timber properties could enable a better non-destructive grading and more efficient use of this material, ultimately reducing construction costs while maintaining safety.

Wood is a highly orthotropic material; due to its anatomical structure, it is significantly stronger and stiffer along the grain than in the transverse directions. For defect-free, straight-grained specimens, a relationship can be established between wood
strength or stiffness and the angle between the fibre direction and the applied load, which can be derived from material stiffness matrix rotation (Olsson et al, 2013; Pot et al, 2024). This relationship is usually well described by the well known empirical Hankinson formula (Hankinson, 1921). Detailed values of anisotropic material parameters for several wood species can be found in (Kretschmann, 2010) or (Guitard, 1987), along with their variation as a function of influencing factors such as relative humidity, density, and moisture content.

However, the presence of knots complicates this relationship. Knots are considered the most critical defects in structural timber (Bodig and Jayne, 1982; Guindos and Guaita, 2012), significantly reducing its strength. This impact is largely attributed to the severe local fibre deviation around knots. Understanding local fibre orientation and integrating this information into predictive models for the mechanical properties of individual boards is therefore at the forefront of current research.

Recent strength grading machines are equipped with laser dots which allow to measure fibre on board surfaces. When spotting a wood surface with a laser dot, the lighted spot gets an elliptical shape due to light dispersion. As long as the ratio of the major to minor axes of the laser dots' dispersion ellipses is large enough, the ellipses' major axes give a reliable indication of the fibre direction projected onto the visible surface (Besseau et al, 2020). This phenomenon is called the tracheid effect as it was first observed on softwood (Nyström, 2003; Simonaho et al, 2004). This technology has led to the development of more precise strength grading methods that utilise laser dot measurements. These methods rely on rather direct regression and correlation analysis (Viguier et al, 2015, 2017; Rais et al, 2021; Jehl et al, 2011), or rely on beam theory to compute stiffness profiles, enabling more accurate strength predictions compared to previous approaches (Olsson et al, 2013; Viguier et al, 2017; Olsson et al, 2018; Briggert et al, 2020; Olsson et al, 2022). Still, those rely on interpolations to estimate local material orientations within the board from surface fibre directions.

Subsequent studies have aimed to refine these interpolation models by improving the definition of internal fibre orientation. One approach proposed by Briggert et al (2018) involved using the ellipse ratio to determine a diving angle ; however, this approach revealed to be ineffective because the laser dot ellipse ratio is highly influenced by other factors, such as wood moisture content (Purba et al, 2020) and the density ratio between earlywood and latewood (Besseau et al, 2020). Hu et al (2018b) have integrated the position of the pith into interpolation models to better capture fibre trajectories, showing slight improvements in the prediction of the bending strength (Olsson et al, 2022). Some studies have leveraged Foley's internal fibre orientation model to estimate fibre deviation around knots (Foley, 2003; Lukacevic et al, 2019; Habite and Olsson, 2022), however raising the difficulty of adjusting the model parameters.

Recent advances in X-ray tomography techniques for sawn timber analysis have enabled the direct acquisition of knot positions and growth layer structures, providing valuable input for fibre orientation modelling (Jacquin et al, 2017; Duchateau et al, 2013; Longuetaud et al, 2023, 2022; Huber et al, 2022). A persistent challenge lies in the limited resolution of X-ray tomography, which prevents simultaneous visualization of both the full geometry of a knot and the detailed structure of surrounding fibres. To
overcome this, Hu et al (2022) or Müller et al (2018) conducted multiscale tomography
studies of branch-stem junctions to better capture the actual fibre paths near knots.
Multiscale observations of the branch-stem junction have provided crucial insights into
fibre orientation patterns near knots. In particular, Hu et al (2022) revealed abrupt
fibre orientation changes at the upper part of the junction, supporting Müller et al
(2018)'s concept of sacrificial tissue formation there.

Although anatomically relevant, these studies are limited to local observations and
lack a global view of fibre orientation in sawn timber, which is essential for accurate
mechanical behaviour prediction at the board scale.

On the other hand, Huber et al (2022) proposed a model based on X-ray tomography, in which internal fibre orientation around knots is derived from angular interpolation and the computation of a flow field relative to the distance to the segmented knots. This approach thus partially follows the widely used flow-grain analogy (Bodig and Jayne, 1982; Lang and Kaliske, 2013; Foley, 2003). In addition, this study introduced the use of the gradient structure tensor (GST) to estimate the local radial and tangential directions in tomographic longitudinal slices—specifically in regions assumed to be undisturbed by knots—by identifying the direction of maximum gradient as orthogonal to the growth layers. Based on the reconstructed geometry of the timber board—including pith location, knots, and local fibre orientations—Huber et al (2022) developed 3D and 1D finite element models (FEM) to simulate four-point bending tests. The resulting mechanical behaviour was evaluated by comparing simulated and experimental values of bending strength, global modulus of elasticity, and local modulus of elasticity.

However, these global estimators offer only limited insight into the local accuracy
of the reconstructed fibre orientation fields. Complementary validation strategies that
provide spatially resolved comparisons—such as matching simulated finite element
models (FEM) stiffness profiles with those derived from digital image correlation (DIC)
measurements—may therefore provide a more reliable basis for assessing internal fibre
orientation models, as demonstrated in (Hu et al, 2018b) and (Lukacevic et al, 2019).
Furthermore, all (Huber et al, 2022; Hu et al, 2018b; Lukacevic et al, 2019), express
the need of accurate fibre reconstruction for fibre orientation calibration.

To address this need, several destructive methods for determining internal fibre ori-
entation have been developed, which rely directly on serial sectioning. For instance, Hu
et al (2018a) used in-plane laser dot measurements from two serially sectioned orthog-
onal planes of a supposed symmetrical knot, thus enabling a 3D reconstruction of fibre
orientation from both in-plane measurements. Later, Penvern et al (2024) proposed a
simplified method based on successive planing in one direction and reconstructing fibre
orientation from visible growth layers, under the assumption that fibres are tangent to
these layers. This method was further extrapolated by Hu et al (2025), who used X-
ray tomography to extract growth layers and reconstruct fibre directions for an entire
spruce board, enabling the generation of full-volume internal fibre orientation data.

Serial sectioning, as commonly applied in materials science, offers an effective
alternative when non-destructive methods are limited in the specific type of inter-
nal information they can provide (Zankel et al, 2014). In such contexts, successive
2D observations—acquired destructively—can be recombined into a rich 3D dataset.
Here, this principle is applied by adding in-plane fibre orientation measurements to the greyscale images used to visualise growth layers, to allow genuine internal fibre orientation computation.

More recently, Huber and Olofsson (2023) suggested extending the GST approach to estimate fibre direction as the direction of minimum gradient in the tomographic volume, under the assumption that fibres follow paths of minimal curvature within growth layers. However, this idea has not yet been evaluated or validated.

## 1.2 Aim and Scope

This study pursues two main objectives:

### Methodological development

The study proposes a new adaptation of the method by Penvern et al (2024) to accelerate internal fibre orientation determination in sawn timber boards. Instead of relying on growth layer segmentation and reconstruction, the new method uses GST analysis to estimate the local radial direction—defined as the normal to growth layers—thereby significantly reducing computational cost. This use of GST was previously demonstrated by Huber et al (2022) in longitudinal 2D tomographic slices in fibre-straight regions.

### Mechanical validation

The study evaluates the effectiveness of this new fibre orientation determination (IFOD) method against other methods—namely that of Penvern et al (2024) and the GST-based full-volume approach suggested by Huber and Olofsson (2023)—by comparing resulting local mechanical behaviour. This is done through analysis of bending stiffness profiles obtained via DIC and FEM simulations incorporating the various orientation datasets.

The different IFOD methods considered in this study, based on either destructive serial sectioning (DSS) or computed tomography (CT) data, are summarised in Table 1, with bibliographical references. The method introduced in the present work is highlighted in gray, and a short identifier is assigned to each method for reference in the following sections.

# 2 Materials and experimental procedure

One 2775 mm long board of Douglas fir grown in France (*Pseudotsuga menziesii* Mirb. Franco) was selected from a larger batch for its prominent central knot and significant global fibre deviation, as visible in Figure 1a. The board was dried and squared before being scanned in an industrial scanner. The board was then weighed with a *Ohaus I-20W* scale and its dimensions measured with a measuring tape and its density accordingly computed. The moisture content was also assessed with a *testo 606-2*. Results are given in Table 2. A reference number was written on the tree crown side of the board. This numbering makes it possible to identify the board within the larger batch for further study, and used to easily identify the board orientation, as described in Figure 1b.
<table><caption>**Table 1:** Internal fibre orientation determination methods summary</caption><thead><tr><th>denomination</th><th>IFOD methods</th><th>reference</th></tr></thead><tbody><tr><td>sym,DSS</td><td>destructive serial sectioning method with laser dot measurements of the projected angles of the fibres of a supposedly symmetrical knot in orthogonal planes</td><td>(Hu et al, 2018a)</td></tr><tr><td>seg,DSS</td><td>destructive serial sectioning method with laser dot measurements of in-plane fibre angle and normal to growth layer defined by principal component analysis (PCA) after reconstruction of growth layers by planed surfaces colour scans segmentation</td><td>(Penvern et al, 2024)</td></tr><tr><td>seg,CT/DSS</td><td>destructive serial sectioning method with laser dot measurements of in-plane fibre angle and normal to growth layer defined by PCA after reconstruction of growth layers by board X-ray tomography segmentation</td><td>(Hu et al, 2025)</td></tr><tr><td>GST,DSS</td><td>destructive serial sectioning method with laser dot measurements of in-plane fibre angle and normal to growth-layer computed by planed surfaces gray scans GST analysis extrapolated by flow analogy extrapolation and without optic flow extrapolation</td><td>present study</td></tr><tr><td>GST,DSS,init</td><td>non-destructive method fully based on board X-ray tomography GST analysis, here applied on gray scans resulting from destructive serial sectioning</td><td>(Huber and Olofsson, 2023)</td></tr><tr><td>GST,CT*</td><td>non-destructive method fully based on board X-ray tomography GST analysis, here applied on gray scans resulting from destructive serial sectioning</td><td>(Huber and Olofsson, 2023)</td></tr></tbody></table>

The board was monitored by image correlation during non-destructive elastic bending tests as described in paragraph 2.1, to measure its stiffness. Finally, the board was entirely planed and scanned, as described in paragraph 2.2, to yield a complete greyscale DSS image stack, along with a full in-plane fibre angle field mapped over the entire board.

<table><caption>**Table 2:** Dimensions, weight and moisture content of selected board</caption><thead><tr><th>dimensions (mm)</th><th>weight (kg)</th><th>density (kg m<sup>-3</sup>)</th><th>moisture content (%)</th></tr></thead><tbody><tr><td>2775 x 146 x 41</td><td>8.50</td><td>511.7</td><td>9.8</td></tr></tbody></table>

## 2.1 Bending Test and Digital Image Correlation monitoring

To measure the actual stiffness of the beam locally, it was decided to use stereo digital image correlation to monitor the board displacements during 4-point bending tests. The loading span was 625 mm, with a support span of 1875 mm, under symmetrical loading conditions, as presented in Figure 2. After painting a speckle pattern on the board, it was placed on a 3R RP100 ELM bending bench, equipped with a LaVision StrainMaster Portable system, composed of a laptop computer to run, a device control
<figure><img src="image_3.png" /><figcaption>(a) Board enhanced visualization using linearly interpolated laser dot angle measurements (color-mapped from pink to green) and the corresponding in-plane fibre orientation field represented as streamlines, displayed side by side with the original greyscale scan image of the board.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 1:** (a) Visualization of investigated board (#39) and (b) associated reference frame.</figcaption></figure>

unit with integrated AD converter, two 9 Mpx cameras, two continuous LED illumination units and a tripod. The *LaVision StrainMaster Portable* system was used to collect and treat, via the *LaVision* software *DaVis 10.2* the correlated images.

The orientation and spacing of the cameras (36 cm) were chosen as a compromise between stereo correlation accuracy and the size of the monitoring window within the actuator frame, in the constant moment zone. As a result, only approximately 20 cm of the board could be observed at a time. The related Cameras' resolution after calibration was 12.184 px/mm. For the digital correlation displacement and strain calculation, a typical subset size of 31 px and a step size of 10 px were chosen allowing a resolution of 1.096 px/mm. Areas affected by edge effects, corresponding to an approximate 5 mm border of the board, were removed from the strain and displacement maps.
<figure><img src="image_5.png" /><figcaption>Fig. 2: Bending bench initial configuration</figcaption></figure>

Because of this limited window size, it was chosen to reproduce multiple bending tests, at regular intervals, along the beam. Initially, the numbered board's end was placed at the limit of the bending bench support, with the reference number towards the cameras. The board was bent with a loading speed of 5 mm/min up to a maximum stress of ~15 MPa, within the elastic range, before being released and moved to be bent and monitored again every further ~14 cm. The initial set-up is presented in Figure 2, also showing in red the proportion of the area covered by DIC monitoring through the successive bending tests.

Seven bending tests were performed, covering a total length of 1195 mm. For the bending monitoring, pictures of the board were taken every ~3 MPa increase of the bending stress, starting from the unloaded board as new reference frame, and the corresponding force recorded.

The translation monitoring was also achieved by image correlation, pictures of the board being taken every ~3 cm displacement step. This information was used to reposition the bending test results in relation to the initial configuration of the board. It was ensured that the DIC coordinate system shared the same *x* – *axis* and *y* – *axis* as the board thanks to the initial DIC calibration. Thus, the data was repositioned into the board's global reference frame by precisely determining the *x* – *axis* and *y* – *axis* offsets between the same pattern (usually in the vicinity of a knot or crack) in the board's greyscale image—aligned with the global reference frame—and in the successive bending DIC reference images, which have coordinates relative to the board's initial configuration.

## 2.2 All through board planing and scanning operations: destructive serial sectioning

In order to reconstruct the board's internal fibre orientation field as in (Penvern et al, 2024), an all through process of successive planings and scannings was applied. First, a 10 mm thick plywood board and a 16 mm thick MDF board of the same exact width and length of the considered timber board were glued onto its top face. This created a plywood-MDF-board collage. MDF was chosen as the interface with the board to prevent laser diffusion in transparency through the considered wood material during the final scanning operations, while the plywood provided greater rigidity during scanning
and planing. The planings, in 1 mm increments, were carried out using an SCM S520 planer.

Prior and between each planing, the board was scanned in an industrial timber grading scanner. It is equipped with a near infrared (830 nm) laser dots and line device (cf. Figure 1a). The gray images, obtained from the laser line lighting and the laser dots ellipses were both acquired by the same matrix cameras. This ensures that the destructive serial sectioning data, i.e. the board stacked gray images and the fibre in-plane angle measurements, share the same coordinate system. Only the bottom face matrix camera was used, as it ensured to have same camera focus for each successive scanning of the planed face of the destructive serial section. Because of unevenness of the laser line lighting, the gray images needed lighting bias correction. This was achieved by removing the related mean crosswise intensity variation along the board.

The detected laser dots elliptical contours were fitted by an ellipse using Python OpenCV library cv2.fitEllipse function to retrieve each ellipse's centre, angle to the board's axis, and minor and major axis lengths. Fitted laser dots ellipses mean dimensions were 14.81 px x 9.52 px. Given a recorded bottom matrix camera resolution of 4.21 px/mm, this gives an average ellipse size of 3.5 mm x 2.3 mm with a laser dot projected every 4.7 mm approximately.

The longitudinal resolution depends on a roller encoder and was therefore the same for all the data, fixed at 1 px/mm. A roller presser ensured that there was no slipping. In this respect, every image was built during the board's longitudinal displacement through the scanner at a speed of 62 m/min, with a line image being taken every 1 mm. Thus, the greyscale image resolution was 4.21 px/mm × 1 px/mm, while the laser dot resolution was approximately 0.18 px/mm × 1 px/mm.

To determine the fibre orientation from these two datasets, it was therefore decided to work at an intermediate resolution of 1 mm × 1 mm, resulting in a referenced regular mesh of 1 mm × 1 mm × 1 mm for the board after the successive planing/scanning operations.

The regularised in-plane fibre angle field was obtained by linear interpolation using the scipy.griddata function on the laser dot measurement scattered data, while the regularised greyscale destructive serial sectioning of the board was obtained by linear interpolation of the corrected greyscale images of the surfaces revealed at each 1 mm planing step, using the cv2_resize function with the cv2 INTER_AREA interpolation scheme.

For more details on this scanner and possible data processing, the reader can refer to (Longuetaud et al, 2022).

# 3 Method

## 3.1 Generic Internal Fibre Orientation Determination Method

To get the internal fibre orientation inside the board, the idea developed in (Penvern et al, 2024) is to make use of the fact that fibres follow growth layers. Thus, when having the projected fibre angle $α$ obtained by laser dot measurements, and after determining the local normal direction to the growth layers either by PCA analysis of the manually segmented and reconstructed growth-layer 3.1.1 or by GST analysis
<figure><img src="image_6.png" /><figcaption>Fig. 3: Principle of defining the local fibre orientation from the projected angle measured by laser points and the normal direction to the growth layers</figcaption></figure>

of the gray images stack 3.1.2, it is possible to compute the local fibre direction as illustrated in Figure 3.

Let $\vec{v}$ be the 3D fibre direction in the specimen, this direction $\vec{v}$ can be decomposed in an out-of-plane direction $v_z \vec{z}$ and an in-plane direction $\vec{u}$, computed from the in-plane fibre angles measurements as :

$$ \vec{u} = \begin{pmatrix} \sin \alpha \\ \cos \alpha \\ 0 \end{pmatrix}_{(\vec{x}, \vec{y}, \vec{z})} \qquad (1) $$

If $\vec{n} = \begin{pmatrix} n_x \\ n_y \\ n_z \end{pmatrix}_{(\vec{x}, \vec{y}, \vec{z})}$ is the growth layer normal direction, then $\vec{v}$, parallel to the growth-layer, is orthogonal to $\vec{n}$, i.e. $\vec{v} \cdot \vec{n} = 0$, or $(\vec{u} + v_z \vec{z}) \cdot \vec{n} = 0$, hence :

$$ \vec{v} = \begin{pmatrix} v_x \\ v_y \\ v_z \end{pmatrix}_{(\vec{x}, \vec{y}, \vec{z})} = \vec{u} + v_z \vec{z} = \begin{pmatrix} \sin \alpha \\ \cos \alpha \\ -\frac{\vec{u} \cdot \vec{n}}{n_z} \end{pmatrix}_{(\vec{x}, \vec{y}, \vec{z})} \qquad (2) $$

The following sections thus outline the different internal fibre orientation methods based on how the normal direction to the growth layer, $\vec{n}$ is determined.

### 3.1.1 seg, DSS IFOD by destructive serial sectioning growth-layer segmentation and in-plane fibre angle laser dot measurements

The method developed in (Penvern et al, 2024) based on destructive serial sectioning and growth layer reconstruction, as illustrated in Figure 4, was applied. In this
<figure><img src="image_7.png" /><figcaption>**Fig. 4:** Illustration of IFOD method by Penvern et al. (2023) applied to the board of present study with (a) initial gray scan of the board, (b) semi-manual growth layer segmentation, (c) growth layer reconstruction, and (d) streamlines computed from fibre orientation field overlaid with knot masks from destructive 3D-scan intensity threshold.</figcaption></figure>

method, the growth layers are reconstructed by first applying a watershed segmentation to the corrected DSS greyscale image stack. The segments are then manually adjusted, merged and identified according to the growth layer they belong to (Figure 4b). Finally, point clouds of the growth layers boundary surfaces (Figure 4c) are generated, so that normal directions to those growth layers may be computed using PCA. The obtained normal directions, computed at each point of the growth layers' point clouds, forming a scattered dataset, are then linearly interpolated using `scipy.griddata` function to the shared 1 mm by 1 mm by 1 mm regular grid of the board's volume.

Combining this data to the laser dot fibre angle using Equation 2, a high-resolution coordinated dataset of fibre direction of the full board is computed. Illustration of the fibre orientation result, showing fibre streamlines using *PyVista* Python library on the fibre orientation field is visible in Figure 4d.

Although the method proves effective, the semi-manual segmentation of a board's growth layers requires approximately two weeks, limiting its applicability to large-scale studies.
### 3.1.2 GST,DSS,init IFOD by GST analysis of the DSS image stack and laser dot in-plane fibre angle measurements

The Gradient Structure Tensor (GST) is a powerful mathematical tool used to extract the principal directions of a material from tomographic images. It is particularly useful in analysing anisotropic structures, such as wood, composite materials, and biological tissues. It provides local orientation information by capturing the dominant direction of intensity variations.

In the current study, the Gradient Structure Tensor was computed on board's 1 mm × 1 mm × 1 mm regularised DSS image stack using a $\sigma = 3$ px Gaussian filter by means of the OrientationPy Python module developed by the EPFL Center for Imaging (2022).

$$GST = \begin{bmatrix} (J_x(p))^2 & J_x(p)J_y(p) & J_x(p)J_z(p) \\ J_x(p)J_y(p) & (J_y(p))^2 & J_y(p)J_z(p) \\ J_x(p)J_z(p) & J_y(p)J_z(p) & (J_z(p))^2 \end{bmatrix} \quad (3)$$

where $J_x$, $J_y$ and $J_z$ are the partial Gaussian derivatives of the image intensity $I$:

$$J_x = w_{\partial x, \sigma} * I ; \quad J_y = w_{\partial y, \sigma} * I ; \quad J_z = w_{\partial z, \sigma} * I \quad (4)$$

Eigenvalue decomposition of the GST computed on the greyscale DSS image stack yields eigenvectors that represent the principal directionalities of the overall spatial organisation within the greyscale field. Thus, OrientationPy allows to compute 2D or 3D image directionalities, either in a "fibre mode" - this outputs the GST eigenvector field of lowest eigenvalue, i.e. the direction of the lowest gradient, as along a fibre; or in a "membrane mode" - this outputs the GST eigenvector field of highest gradient eigenvalue, i.e. the direction of highest gradient, as in the direction normal to a membrane. This "membrane mode" is therefore the mode used in this method to determine the direction normal to growth layers, i.e. the wood local radial direction. The computation of this direction in the whole board with the current 1 mm resolution grid requires less than two minutes with GST, underlining the method's high efficiency and suitability for large-scale applications.

### 3.1.3 GST,DSS IFOD by GST analysis of extrapolated DSS image stack and in-plane fibre angle laser dot measurements

To investigate the possibility of reducing GST analysis edge effects, extrapolating the DSS image stack in the $z$-direction using the cv2.calcopticalFlowFarneback function from the OpenCV Python library was explored. The cv2.calcopticalFlowFarneback() is a dense optical flow algorithm that estimates pixel displacement between two consecutive frames in an image sequence. The function works by analysing patterns in the intensity variations of the images to compute a per-pixel flow field. By calculating the flow between consecutive stack images, it is possible to get an estimation of how the growth layers would extend in additional frames. Two additional frames in both positive and negative $z$-directions were thus above and beneath the board. In this way, the DSS image stack as input array of the GST is extrapolated by fluidic analogy where the Gaussian
<figure><img src="image_8.png" /><figcaption>Fig. 5: Illustration of optic flow extrapolated DSS image stack</figcaption></figure>

derivative filter overlaps the board broad external surface. Where the GST is applied on the initial gray images stack, the method will be refered to as *GST,DSS,init*, while it will be refered to *GST,DSS* when applied on the extrapolated stack. Stack extrapolation is achieved in a few seconds on a laptop computer.

### 3.1.4 GST,CT\* internal fibre orientation determination by sole GST analysis of the DSS image stack

When no in-plane fibre angle data is available and the wood direction is to be determined solely from X-ray computed tomography, rather than assuming that fibres are merely parallel to the growth layers, they may, as suggested by Huber and Olofsson (2023), be considered as following the path of least curvature within them. The fibre orientation may then be directly computed from the GST analysis of the X-ray CT scans of the board 3.1.4.

This method was applied to the optic flow-extrapolated DSS image stack obtained; however, the same approach could be applied non-destructively using X-ray computed tomography (CT), as demonstrated in (Huber and Olofsson, 2023), hence the designation *GST,CT\**. The normal to the growth layer direction, i.e. the wood local radial direction may be defined by the GST eigenvector field of highest eigenvalues, while the fibre direction, i.e. the wood local longitudinal direction may in turn be defined by the GST eigenvector field of lowest eigenvalues. The GST being symmetric, the eigenvectors are orthogonal, ensuring that radial and longitudinal directions of the board are orthogonal. Here, the local radial direction field was obtained as previous, using the "membrane mode" of the *OrientationPy* module with $\sigma = 3$ px, while the
local fibre direction field was obtained using the "fibre mode" of the image stack directionalities computation of the above mentioned OrientationPy module with similar structure tensor computation parameters.

## 3.2 Localised Bending MoE Profiles from Beam Theory

Beam theory provides the following relationships for a homogeneous isotropic beam subjected to pure bending in the *xy*-plane with a constant moment *M*(*x*) of *z*-axis :

$$ \kappa(x) = \frac{M(x)}{E_b(x)I_z} \qquad (5) $$

$$ \sigma_{xx}(x, y) = \frac{M(x) \bar{y}}{I_z} \qquad (6) $$

$$ \varepsilon_{xx}(x, y) = \frac{M(x)}{E_b(x)I_z} \bar{y} \qquad (7) $$

where:

* $\kappa(x)$ is the curvature of the beam at position *x* along its length,

* $M(x)$ is the internal bending moment at location *x*,

* $E_b(x)$ is the localised bending MoE (along the beam axis),

* $I_z$ is the second moment of area of the cross-section with respect to the z-axis (perpendicular to the bending plane),

* $\varepsilon_{xx}(x, y)$ is the corresponding normal strain,

* $\bar{y}$ is the vertical coordinate within the cross-section, measured from the neutral axis (which corresponds to zero strain).

### 3.2.1 Digital Image Correlation Analysis

Thus, according to the beam theory, the strain slope coefficient $\frac{\varepsilon_{xx}(y)}{\bar{y}}$ and the curvature $\kappa$ should be equal (to $\frac{M}{EI_z}$) when considering a homogeneous material. In our case, because wood is a heterogeneous and orthotropic material, these two variables differed significantly. This discrepancy is evident when plotting these variables — normalised by the force exerted by the bending test bench and averaged over a longitudinal window of $w_x = 90$ mm — after compiling the results of successive bending tests, as visible in Figure 6a and 6c.

Nonetheless, both variables can be used to compute a local bending modulus of elasticity (MoE) to provide a stiffness profile of the board, either through the DIC longitudinal strain or transversal displacement results. Thus, in the DIC monitored *xy*-plane of the board, in the area of pure bending of moment $M(F)$, two localised MoEs were computed :

$$ E_{b,\varepsilon,w_x,DIC}(x) = \frac{M(F)}{I_z \cdot \left\langle \frac{\varepsilon_{xx,s,DIC}(x,y)}{\bar{y}} \right\rangle_{w_x}(x)} \qquad (8) $$

$$ E_{b,\kappa,w_x,DIC}(x) = \frac{M(F)}{I_z \cdot \langle \kappa_{DIC}(x,y) \rangle_{w_x}(x)} \qquad (9) $$
<figure><img src="image_9.png" /><figcaption>**Fig. 6:** Displacement curvature (a) and strain slope coefficient (c) from successive bending tests with examples of linear (b) and circle (d) fits obtained with the 3rd bending DIC measurements of the $w_x$ = 90 mm wide window centred at $x$ = 1245 mm</figcaption></figure>

The averaged strain slope coefficients $\left\langle \frac{\epsilon_{xx,s,DIC}(x,y)}{\bar{y}} \right\rangle_{w_x}$ were computed by linear regressions performed using the ` numpy.polyfit ` function on the DIC smoothed longitudinal strain $\epsilon_{xx,s,DIC}(x,y)$ over a window of length $w_x$ = 90 mm. The smoothed strain was obtained by applying a convolution uniform filter (`scipy.ndimage uniform filter`) of 19 px × 19 px on the raw DIC strain measurements $\epsilon_{xx,DIC}$. The DIC resolution being 1.096 px/mm, this corresponds to a uniform filter of size 20.8 mm × 20.8 mm. This is similar to the computation made in (Hu et al, 2018b) where the regression was run over the average strain value of the surrounding area of about 20 mm × 20 mm, and then the localised MoE averaged over a longitudinal resolution window of $L_r$ = 50 mm, for the comparison of DIC measured and FE simulated stiffness profiles of two boards of Norway spruce (*Picea abies* (L.)
H. Karst). This averaging window length has also been extended to 90 mm by Olsson et al (2022) in the definition of their indicative properties for the prediction of local bending strength, with application on Douglas fir sawn timber. An example of such regression is given in Figure 6b. Figure 6a illustrates the computed DIC averaged strain slope coefficients normalised by a constant moment corresponding to a 1 kN exerted force, $\langle \frac{\epsilon_{xx,s,DIC}(x,y)}{y} \rangle_{F=1kN, w_x=90mm}$, of the successive bending tests in relation to the board $x$-axis.

The local curvatures $\langle \kappa_{DIC}(x, y) \rangle$ were computed on the basis of the DIC displacement measurements directly, using G. Taubin algebraic circle fit (circle-fit.TaubinSVD) supplied by Klear (2023) over the $y$-wise mean transversal $y$-displacement values over the same $w_x = 90$ mm longitudinal window as for the strain slope coefficients. An example of circle fit is given in Figure 6d. Figure 6c illustrates the computed DIC averaged curvature normalised by a constant moment corresponding to 1 kN exerted force, noted $\langle \kappa_{DIC}(x, y) \rangle_{F=1kN, w_x=90mm}$, of the successive bending tests of board #39 again, in relation to the board $x$-axis.

In parallel with the processing of the DIC results, the recorded local deflection $w_l$ as a function of the force exerted in the various bending tests (see Figure 2) was measured on the bending bench, and processed in order to obtain reference values of local stiffness $E_{b,bench}$.

$$E_{b,bench} = \frac{l_1^2}{16I_z} \left\langle \frac{\Delta w_l}{\Delta M} \right\rangle \qquad (10)$$

where $l_1 = 5/6 \cdot \text{loading span} = 521$ mm is the base length for determining the local deflection and $\langle \frac{\Delta w_l}{\Delta M} \rangle$ is the mean deflection to moment slope computed by `numpy.polyfit` over the subsequent increments of the successive bending tests. This made it possible to verify the actual correlation between the stiffness profiles obtained from the DIC longitudinal strain and transverse displacement results (plotted in Figure 7) and the reference local MoE values from the bending bench, indicated by cyan dots.

It is interesting to note that, although all these quantities are related to localised elastic moduli, the two quantities $E_{b,\epsilon,w_x,DIC}(x)$ and $E_{b,\kappa,w_x,DIC}(x)$ (Figure 7), just as $\langle \frac{\epsilon_{xx,s,DIC}(x,y)}{y} \rangle_{F=1kN, w_x=90mm}$ and $\langle \kappa_{DIC}(x, y) \rangle_{F=1kN, w_x=90mm}$ (Figure 6 produce significantly different outcomes. In particular, the MoE profile computed from the curvature exhibits infinite local values, caused by local inversion of curvature in the knot vicinity. This phenomenon was previously observed and discussed in the numerical study by Pot et al (2024). Nevertheless, since these localised moduli are derived from two different datasets and discussed by different authors in the literature, Hu et al (2018b); Pot et al (2024); Huber et al (2022), they can be used to add details that will prove useful in justifying the reliability of the investigated fibre orientation reconstruction methods for predicting the mechanical properties of sawn timber.

### 3.2.2 Finite Element Simulation Analysis

To get the bending profiles of the investigated board, finite element simulations were carried out using the *Abaqus* software.
<figure><img src="image_10.png" /><figcaption>**Fig. 7:** Overlaid Flexural Moduli, computed from DIC strains ($E_{b,ε,DIC}$ in dashed lines), or from the DIC displacement curvatures ($E_{b,κ,DIC}$ in continuous lines), over the successive DIC monitored bending phases (color shade)</figcaption></figure>

A three-dimensional finite element model of the board was created in *Abaqus* based on its measured dimensions $l_x, l_y, l_z$. The material is defined as orthotropic and the material parameters defined as proposed by Guitard (1987), taking into account the global density of the board and its moisture content (MC). The estimated shear and elastic moduli are given in Table 3.

<table><caption>**Table 3:** Board material parameters according to Guitard (1987)</caption><thead><tr><th>density (g/cm³)</th><th>MC (%)</th><th>E<sub>L</sub></th><th>E<sub>R</sub></th><th>E<sub>T</sub></th><th>G<sub>LR</sub></th><th>G<sub>LT</sub></th><th>G<sub>RT</sub></th><th>ν<sub>LR</sub></th><th>ν<sub>LT</sub></th><th>ν<sub>RT</sub></th></tr></thead><tbody><tr><td>511.7</td><td>9.8</td><td>16190</td><td>1222</td><td>804</td><td>1055</td><td>859</td><td>104</td><td>0.39</td><td>0.44</td><td>0.51</td></tr></tbody></table>

A structured mesh of hexahedral elements with an element size of 2 mm was generated, chosen as a compromise between resolution and computational resources, along a regular rectangular parallelepiped defined by the board's dimensions.

Given the orthotropic nature of wood, a discrete field was created to define local material axes for each finite element. A KD-tree search using `scipy.spatial.KDTree.query` was employed to match each element to its nearest material orientation reference point. The element's principal material direction was aligned with the mean fibre direction within the element volume, while the secondary direction was oriented along the mean normal direction to the growth layers, and the third direction was set orthogonal to the first two.

Default longitudinal and radial directions were set to $\vec{x}$ and $\vec{z}$, respectively, consistent with the plain-sawn configuration of the board. Where in-plane angle values were missing, at the board's edges, due to laser dot spacing, the fibre direction was computed as $\vec{n} \wedge (\vec{x} \wedge \vec{z})$.
<figure><img src="image_11.png" /><figcaption>Fig. 8: Finite element model loading set-up</figcaption></figure>

To simulate a pure bending scenario, the board was subjected to the following constraints:

* a fully fixed constraint at the "left" edge midpoint ($x = 0$), i.e. $(0, l_y/2, l_z/2)$, where all translational degrees of freedom were set to zero;

* a partially fixed constraint on the "right" edge ($x = l_x$), where only transverse displacements were constrained to prevent rigid body motion while allowing bending deformation;

* an additional constraint at the origin $(0, 0, 0)$ to prevent out-of-plane rotation by fixing displacement in the z-direction.

A distributed pressure load was applied to both the left and right faces of the board. The loading followed an analytical expression designed to induce a constant bending moment equivalent to the maximum moment applied in the bending bench for a 1 kN applied force, ensuring consistency with experimental conditions.

These loading conditions are illustrated in Figure 8.

A static bending step was defined in *Abaqus* to simulate the board's response under the applied loading. Upon completion of the simulation, nodal data—including displacement and strain fields—were extracted and compared to the experimental DIC results. Local bending moduli of elasticity, $E_{b,ε,FEM}$ and $E_{b,κ,FEM}$, were computed in the same way as their DIC counterparts, $E_{b,ε,DIC}$ and $E_{b,κ,DIC}$, based on simulation-derived local curvature $\langle \kappa_{FEM}(x, y, z = l_z) \rangle w_x(x)$ and the slope of longitudinal strain $\langle \frac{\varepsilon_{xx,FEM}(x,y,z=l_z)}{y} \rangle w_x(x)$ at surface nodes corresponding to the DIC-monitored board face.

$$E_{b,ε,w_x,FEM}(x) = \frac{M}{I_z \cdot \left\langle \frac{\varepsilon_{xx,FEM}(x,y,z=l_z)}{y} \right\rangle_{w_x}(x)} \quad (11)$$

$$E_{b,κ,w_x,FEM}(x) = \frac{M}{I_z \cdot \langle \kappa_{FEM}(x, y, z = l_z) \rangle_{w_x}(x)} \quad (12)$$

The sliding averaging window size $w_x$ was chosen equal to 90 mm as previous.
# 4 Results and discussion

## 4.1 Fibre orientation analysis

Fibre orientation results of the various IFOD methods are illustrated in Figure 9.

<figure><img src="image_12.png" /><figcaption>**Fig. 9:** **a.** colour scan and fibre orientation fields coded in RGB values as in **c.** colour circle, at *z* = 10 mm, with **b.** a zoom in the vicinity of the central knot ($x \in [1050, 1250]$ mm, $y \in [0, 100]$ mm) according to *seg,DSS* (red frame), *GST,DSS* with (green frame) and without (orange frame) image stack extrapolation, and *GST,CT*\* (blue frame) IFOD methods; **c.** overlay of *GST,CT*\* fibre orientation field and colour scan zooms</figcaption></figure>

It is noteworthy that the methods *seg,DSS* and *GST,DSS*, whether with or without DSS image stack extrapolation, yield overall very similar results (Figure 9 **a.**). However, significant discrepancies are observed in the knot vicinity, particularly within its boundaries (Figure 9 **b.**). A closer inspection of the knot region reveals that, while the general background distribution of fibre orientation remains consistent across methods, the *seg,DSS* approach produces a much more fluctuating orientation pattern.

This behaviour may reflect genuine local variability in fibre orientation. Nevertheless, it is more plausibly attributed to limitations in the segmentation step, particularly in accurately delineating the growth layer boundaries near complex geometries such as knots. Consequently, incorporating GST-based analysis to define the direction normal to the growth layers — and thus the radial direction of the wood — appears to be a valuable improvement to the method approaching the knot.
<figure><img src="image_13.png" /><figcaption>**Fig. 10:** growth layer normal orientation fields zoom in the vicinity of the central knot ($x \in [1050, 1250]$mm, $y \in [0, 100]$mm), coded in RGB values as in colour circle at $z = 10$ mm, according to *seg,DSS* (red frame), *GST,DSS* with (green frame) and without (orange frame) serial sectioning extrapolation, and *GST,CT*\* (blue frame) IFOD methods</figcaption></figure>

Regarding the differences observed within the knot (i.e., the branch area), it is important to note that the growth layer segmentation was ineffective in this region, as the boundaries of the growth layers were no longer distinguishable. In the *seg,DSS* method, the directions normal to the growth layer are then interpolated across the knot region, leading to artificially smoothed fibre orientations. In contrast, the *GST,DSS* method reveals abrupt directional changes within the knot. These discontinuities are only artefacts resulting from the presence of cracks. These cracks introduce strong directional signals that do not correspond to any actual fibre orientation or growth layer, but to fracture patterns, thereby complicating the interpretation. Hence, both *seg,DSS* and *GST,DSS* give erroneous fibre direction within the knot area due to false growth layer normal field definition (see Figure 10).

However, the impact of the crack becomes critical when estimating fibre orientation based solely on GST analysis, as in the *GST,CT*\* method. When examining the results in terms of fibre orientation, the influence of the crack is evident, with abrupt changes in direction along its path (Figure 9 d.). Yet, independent of this effect, a distinct line can be observed in *GST,CT*\* fibre orientation result along the middle of the board, corresponding to the top of the growth layer surfaces. When the growth layer aligns with the scanning plan, determining the direction of least intensity variation from sole GST analysis becomes more challenging, as visible in Figure 11. The greyscale variability within the growth layers of the scanner image—particularly between consecutive sections—reduces the contrast between growth layers (see red rectangle) can thus mislead the GST analysis.

Nevertheless, this figure also illustrates the good relative positioning of the DSS image stack. If some misalignments could be observed (see blue rectangles), these are not uniform across the board, and therefore difficult to assess. These misalignments in the image stacking may have led to small inconsistencies in the reconstructed fibre orientation field, for all IFOD methods, whether for the GST analysis or for the growth layer segmentation. These two artefacts-contrast variability and gray images
<figure><img src="image_14.png" /><figcaption>**Fig. 11:** extrapolated greyscale DSS image stack orthogonal slices with cross lines colour corresponding to orthogonal views frame colour **a. in the knot vicinity**, at *x* = 1146 mm and at *y* = 48 mm, with some misalignments highlighted by blue squares ; **b. in the clear wood**, at *x* = 1959 mm and at *y* = 72 mm, with middle blur highlighted by red square.</figcaption></figure>

stack misalignments-could be mitigated by higher resolution X-ray tomography. Such tomographies provide clearer differentiation between growth layers and noise, yet with the underlying difficulty of repositioning the scanner laser dot data in relation to the tomographic data as in the GST,CT/SS IFOD method investigated by Hu et al (2025).

The fibre orientation in the vicinity of the central knot is illustrated in Figure 12 by the streamlines obtained using *PyVista* by integrating the *seg,DSS* fibre orientation field from seeds placed within the central knot volume in both directions. This field was chosen because it does not account for the presence of the crack, thus providing a clearer visualisation of the underlying fibre orientation. Above the knot, the fibres are flat, $v_z = 0$, but converging, with a sharp junction point. While below the knot, fibres have an out-of-plane angle, fibre merging more gently with the branch. This corresponds to the branch-stem junction described by Hu et al (2022) or Müller et al (2018) using multi-scale tomography, then that of Hu et al (2025) using Penvern et al (2024) fibre orientation determination method. The fibres from below gradually fuse with the branch, while the more lateral fibres extend around it, wrapping the branch knot until they reach a junction point where fibre orientation diverges significantly. At this point, some fibres turn sharply—up to 90°—to integrate with the branch, while others continue their upward trajectory. It is important to recall that no interpretation of the fibre orientation within the knot should be made, due to the unreliability of the computed orientation field in this area.

## 4.2 DIC and FEM bending MoE profiles analysis

Directly comparing fibre orientation remains challenging, as there is no tangible reference for validation. In contrast, assessing investigated IFOD methods based on their ability to reproduce the mechanical behaviour of the board provides a reliable basis for evaluation.

Figures 13a and 13c show the computed localised bending MoE profiles computed on FEM strain and curvature simulation results respectively with *seg,DSS*,
<figure><img src="image_15.png" /><figcaption>**Fig. 12:** fibre orientation field streamlines in the central knot volume ($x \in [1050, 1250]$mm, $y \in [0, 100]$mm, $z \in [0, 41]$mm)</figcaption></figure>

**GST,DSS, GST,DSS, init** and **GST, CT*** IFOD methods. The modelled curvature is also given in Figure 13b as more legible than the curvature computed MoE. The curves obtained from DIC display large black patches, which correspond to the reconnection or merging zones between successive bending test measures.

### 4.2.1 GST,DSS and seg,DSS results

What is the most striking at first is the good correlation between the DIC results and the FE modelled results, especially with regards with the curvature or curvature computed MoE (Figures 13a and 13b). Determination coefficients for the modelled curvature profiles, visible Table 4, are equal to 0.754, 0.799 and 0.902 for this study **GST,DSS**, **GST,DSS,init** and **seg,DSS** IFOD methods respectively.

Likewise, the modelled bending MoEs computed from FEM strain results showed increased correlations (0.826 > 0.754, 0.859 > 0.799), except for the **seg,DSS** IFOD method, which showed a slight decrease in correlation, with a determination coefficient of $R^2 = 0.882$. In all cases, the coefficients of determination correspond to a root mean
<figure><img src="image_16.png" /><figcaption>(a) localised bending MoE profiles computed on curvature.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>(b) Curvature profiles for a moment of 312.5 N m ($F = 1$ kN).</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>(c) localised bending MoE profiles computed on strain.</figcaption></figure>

<caption>**Fig. 13:** Comparison of FEM and DIC profiles: (a) Apparent modulus from curvature, (b) Curvature profiles, (c) localised modulus from strain.</caption>
square error (RMSE) ranging from 1.0 to 1.4 GPa and a mean absolute percentage error (MAPE) below 10%. These errors can be primarily attributed to the intrinsic sparsity and variability of the DIC measurements. On the one hand, this indicates that the internal fibre orientation is relatively well captured by all three methods; on the other hand, it confirms that fibre orientation is both necessary and suitable to predict the local mechanical behaviour of the board.

Notably, similar RMSE and MAPE values for the local MoE based on strain were reported by Lukacevic et al (2019) for two Norway spruce (Picea abies (L.) H. Karst) boards, following FEM-DIC comparison after modelling the internal fibre orientation using the mathematical fibre deviation model proposed by Foley (2003). The modelling parameters in that study were optimised to fit the experimental data, and the resulting errors ranged from an RMSE of 1.25 GPa to 2.59 GPa, and a MAPE between 8.96% and 13.30%.

These results suggest that Foley's model is an interesting approach for modelling internal fibre orientation. However, the higher RMSE observed in the latter case could potentially be reduced by integrating more accurate fibre deviation data—such as that provided by the methods assessed in the present study—rather than relying solely on parameter optimisation.

Finally, no determination coefficients are reported for the bending modulus of elasticity computed from curvature, as this modulus tends toward infinity and thus renders such coefficients meaningless. Similarly, MAPE and RMSE based on curvature are omitted due to invalid computations resulting from division by zero.

The primary limitation of these IFOD approach lies in their reliance on destructive serial sectioning. Since the method requires successive planing of the board, mechanical tests must be restricted to the elastic domain. Inducing fracture during mechanical testing would compromise the subsequent tomographic reconstruction of fibre orientation, making it incompatible with the method. Consequently, the bending tests performed in this study had to remain within low strain ranges.

The main significant difference between the **GST,DSS** and **seg,DSS** FEM computed profiles occurs, as could be expected, within the knot area. An explanation for this discrepancy can be found in previous paragraph, which illustrates the fibre orientation in the knot vicinity. In the case of GST determined fibre orientation, the bulk knot region is indeed somewhat broader, with a sharp contrast to the surrounding wood, and an inside meander of fibre direction balancing between +z̄ and -z̄. This stems from the fact that the GST analysis is influenced by the presence of drying cracks at the core of the knot. In contrast, the method proposed by Penvern et al (2024) relies on a simple linear interpolation of the fibre direction within knots, as the segmentation of normal growth layers is interrupted near these areas due to the impossibility of delimiting the growth rings in these zones. Nonetheless, it is interesting to note that the sharp variation in fibre direction determined by GST within the knot appeared to lessen the displacement and deformation values observed in the FEM simulations, and led to a reduced correlation between the DIC and FEM profiles. In this context, the interpolation strategy used by Penvern et al (2024), which avoids over-interpretation of the knot core by smoothing the fibre field across disrupted zones, proves advantageous.
The impact of the DSS image stack extrapolation by optic flow extrapolation appears to be negligible, showing neither improvement nor degradation of the results. However, it would be of interest to further refine this extrapolation process, as the generated objects, while visually meaningful (as seen in Figure 5), currently miss an effective extrapolation of the board; as visible in Figure 10, this extrapolation does not correctly extend the growth layer dark circles.

### 4.2.2 GST,CT\* results

Interestingly, simulations based on fibre orientation derived from pure GST analysis in **GST,CT** yield significant correlations, with $R^2 = 0.720$ for MoE profiles (strain-based) and $R^2 = 0.422$ for curvature profiles when comparing FEM and DIC. These results are particularly noteworthy given the known limitations of the method, including defects such as cracks, intensity gradient variability within growth layers, and stack misalignments. In addition, the GST analysis parameters employed in this study were not fully optimised.

This method has the great advantage of being non-destructive, and could be conducted on X-ray computed tomography which would partly resolve current intern growth layer intensity variation and misalignments of presently used destructive serial sectioning. Factors such as the smoothing scale, gradient magnitude thresholds, and filtering strategies were chosen based on preliminary tests and heuristic adjustments. Iterative optimisation of these parameters, possibly incorporating adaptive multi-scale analyses, could also enhance the fidelity of fibre orientation estimation and further improve the correlation between simulation and experiment.

Nonetheless, it should be noted that computing the elastic modulus from strain measurements inherently yields flatter and smoother profiles, which could partly account for the observed high correlation. This observation is consistent with the findings of Pot et al (2024), who showed that the computation of localised bending MoE from curvature leads to a very different behaviour compared to localised bending MoEs derived from strain, as in (Olsson et al, 2018).

The significant influence of knots on the bending stiffness of the board is in all profiles clearly evidenced, and accurately reproduced by finite element simulations when using fibre orientation fields obtained through the **GST,DSS** and **seg,DSS** methods. For instance, starting from a local elastic modulus of approximately 14 500 MPa, values drop to around 7000 MPa when considering averaged strain profiles, and even further to 5000 MPa when evaluating the curvature-based bending stiffness of the board. It is worth emphasizing that this initial value of 14 500 MPa is notably lower than the 16 190 MPa predicted according to Guitard (1987) relations for this specific board wood species, density and moisture content. This discrepancy arises from the fact that the board exhibits spiral grain, visible in Figures 1a or 4, deviating from the straight fibre assumption. This observation is further supported when analysing stiffness profiles derived from the **GST,CT\*** method, which relies solely on growth layer gray tomographic imaging to define fibre orientation. In this case, the modulus starts at 16 190 MPa in knot-free regions and remains globally overestimated along the
board. Thus, a strong limitation to the GST,CT* method consists in its incapability to capture the influence of spiral grain, leading to a systematic overestimation of stiffness.

<table><caption>**Table 4:** R² determination coefficients comparing FEM and DIC results for modelled curvature κ and average bending modulus E<sub>b,avg</sub> for the different IFOD methods, with additional RMSE and MAPE for E<sub>b,avg</sub>.</caption><thead><tr><th></th><th>R<sup>2</sup>(κ)</th><th>R<sup>2</sup>(E<sub>b,avg</sub>)</th><th>RMSE<sub>E<sub>b,avg</sub></sub> [GPa]</th><th>MAPE<sub>E<sub>b,avg</sub></sub> [%]</th></tr></thead><tbody><tr><td>seg,DSS</td><td>0.902</td><td>0.882</td><td>1.288</td><td>9.17</td></tr><tr><td>GST,DSS, init</td><td>0.799</td><td>0.859</td><td>1.069</td><td>6.72</td></tr><tr><td>GST,DSS</td><td>0.754</td><td>0.826</td><td>1.394</td><td>9.81</td></tr><tr><td>GST, CT*</td><td>0.422</td><td>0.720</td><td>2.648</td><td>28.39</td></tr></tbody></table>

# 5 Conclusion

This study aimed to accelerate and validate the internal fibre orientation determination (IFOD) in sawn timber boards by adapting and extending existing methods. A novel approach based on gradient structure tensor (GST) analysis of destructive serial section (DSS) greyscale images was proposed. This method estimates the local radial direction—i.e., the normal to the growth layers—without requiring explicit segmentation, thereby significantly reducing processing time from several weeks to just a few minutes. Despite this simplification, the derived fibre orientation fields yield mechanical predictions comparable to more labour-intensive methods.

The mechanical relevance of these fibre orientation fields was evaluated both qualitatively, by comparison with literature-reported fibre deviation patterns around knots, and quantitatively, by comparing FEM-predicted and DIC-measured bending stiffness. High determination coefficients (R² ≈ 0.8–0.9) were obtained using DSS-based methods, which benefit from the precise fibre orientation measurement enabled by laser dot tracking on planed sections. These strong correlations were achieved despite the absence of explicit knot modelling and without optimisation of local mechanical properties, underlining the crucial role of internal fibre orientation modelling for accurate stiffness prediction in sawn timber.

In addition, the direct estimation of fibre orientation as the direction of minimal intensity gradient in GST analysis was investigated on DSS data, with the goal of enabling its transfer to non-destructive CT imaging. This approach yielded promising correlations with mechanical tests (R² ≈ 0.5), demonstrating its potential for practical application in a non-destructive context.

These findings highlight that proper characterisation of fibre orientation is a key factor for predicting timber strength and stiffness, and for improving mechanical grading approaches. They also point to several areas for further improvement to increase predictive accuracy:

* The explicit modelling of knots and their effects on fibre flow and local stress;

* The representation of cracks as discontinuities in FEM simulations;
* The refinement of elastic property assignment methods, potentially using vibration-based characterisation or high-resolution density maps;

* Using a CT scan rather than the DSS image stack for GST analysis which allows better contrast and no misalignments, but with the difficulty then of relocating the laser dot angle measurement data obtained from DSS.

Nonetheless, some limitations remain:

* The inability to reliably determine fibre orientation within knot volumes;

* The constraints inherent to destructive serial sectioning, which require testing to remain in the elastic regime to preserve boards for scanning;

* The current limitation of the **GST,CT*** method, which cannot account for spiral grain and systematically overestimates stiffness when fibre orientation is inferred only from growth layer intensity gradients.

The proposed DSS-based method thus lays the groundwork for generating reliable internal fibre orientation data, which can in turn support the development of efficient, data-driven models that infer internal orientation from accessible surface or tomographic measurements. This may include the adaptation of existing mathematical fibre deviation models (e.g., Foley-type) or the development of new, wood-specific modelling approaches. These efforts offer promising pathways toward improving timber strength prediction and enhancing mechanical grading procedures in industrial contexts.

In the longer term, reconstructing fibre orientation across entire logs could enable the development of growth models linking trunk geometry (branch size, position, orientation) to fibre architecture. Such models would allow the prediction of future board quality from a given tree, offering valuable perspectives for both forestry management and wood product industries.

## Declarations

This research was funded by the French National Research Agency (EffiQuAss project ANR-21-CE10-0002).

The authors have no relevant financial or non-financial interests to disclose.

The authors are deeply grateful to **Thomas Bultel** and to **Ameur Ben Chaabane** for their help in the investigations.

The authors contributed to this article as follows:

* Conceptualisation: GP, JV, HP;

* Investigation: HP, JV, GP;

* Data Curation: HP;

* Software: HP, GP, JV;

* Formal Analysis, writing – original draft and visualisation: HP;

* Funding Acquisition and project administration: GP;

* Methodology: HP;

* Supervision and writing – review & editing: GP, JV;
## References

Besseau B, Pot G, Collet R, et al (2020) Influence of wood anatomy on fiber orientation measurement obtained by laser scanning on five European species. Journal of Wood Science 66(1):74. https://doi.org/10.1186/s10086-020-01922-y, URL https://jwoodscience.springeropen.com/articles/10.1186/s10086-020-01922-y

Bodig J, Jayne BA (1982) Mechanics of wood and wood composites. Van Nostrand Reinhold

Briggert A, Hu M, Olsson A, et al (2018) Tracheid effect scanning and evaluation of in-plane and out-of-plane fiber direction in norway spruce timber. Wood and Fiber Science 50(04):411-429. URL https://wfs.swst.org/index.php/wfs/article/view/2728

Briggert A, Olsson A, Oscarsson J (2020) Prediction of tensile strength of sawn timber: definitions and performance of indicating properties based on surface laser scanning and dynamic excitation. Materials and Structures 53:1-20

Duchateau E, Longuetaud F, Mothe F, et al (2013) Modelling knot morphology as a function of external tree and branch attributes. Canadian Journal of Forest Research 43(3):266-277

Foley C (2003) Modeling the effects of knots in structural timber. PhD thesis

Guindos P, Guaita M (2012) A three-dimensional wood material model to simulate the behavior of wood with any type of knot at the macro-scale. Wood Science and Technology 47. https://doi.org/10.1007/s00226-012-0517-4, URL https://doi.org/10.1007/s00226-012-0517-4

Guitard D (1987) Mécanique du matériau bois et composites. Cépaduès

Habite T, Olsson A (2022) Modelling of knots and 3d fibre orientation within timber boards based on information obtained from optical scanning. SSRN Electronic Journal https://doi.org/10.2139/ssrn.4110671

Hankinson (1921) Investigation of crushing strength of spruce at varying angles of grain. Air service information circular 3(259):16

Hu M, Briggert A, Olsson A, et al (2018a) Growth layer and fibre orientation around knots in norway spruce: a laboratory investigation. Wood Science and Technology 52(1):7-27. https://doi.org/10.1007/s00226-017-0952-3, URL https://doi.org/10.1007/s00226-017-0952-3

Hu M, Olsson A, Johansson M, et al (2018b) Modelling local bending stiffness based on fibre orientation in sawn timber. European Journal of Wood and Wood Products 76(6):1605-1621. https://doi.org/10.1007/s00107-018-1348-2, URL https://doi.org/10.1007/s00107-018-1349-1
//doi.org/10.1007/s00107-018-1348-2

Hu M, Olsson A, Hall S, et al (2022) Fibre directions at a branch-stem junction in nor-way spruce: a microscale investigation using x-ray computed tomography. Wood Sci-ence and Technology 56(1):147-169. https://doi.org/10.1007/s00226-021-01353-y, URL https://doi.org/10.1007/s00226-021-01353-y

Hu M, Olsson A, Abdeljaber O, et al (2025) Determining growth surfaces and fiber orientation in norway spruce sawn timber using x-ray computed tomogra-phy and optical scanning. Construction and Building Materials 482:141734. https://doi.org/https://doi.org/10.1016/j.conbuildmat.2025.141734, URL https://www.sciencedirect.com/science/article/pii/S0950061825018859

Huber JA, Broman O, Ekevad M, et al (2022) A method for generating finite element models of wood boards from x-ray computed tomography scans. Computers & Structures 260:106702. https://doi.org/https://doi.org/10.1016/j.compstruc.2021.106702, URL https://www.sciencedirect.com/science/article/pii/S0045794921002248

Huber JAJ, Olofsson L (2023) Evaluation of knots and fibre orientation by gradi-ent analysis in x-ray computed tomography images of wood. In: 3rd ECCOMAS Thematic Conference on Computational Methods in Wood Mechanics (CompWood 2023), September 5-8, 2023, Dresden, Germany, International Center for Numerical Methods in Engineering (CIMNE), pp 143-144

EPFL Center for Imaging (2022) Orientationpy. URL https://epfl-center-for-imaging.gitlab.io/orientationpy/introduction.html

Jacquin P, Longuetaud F, Leban JM, et al (2017) X-ray microdensitometry of wood: A review of existing principles and devices. Dendrochronologia 42:42-50. https://doi.org/https://doi.org/10.1016/j.dendro.2017.01.004, URL https://www.sciencedirect.com/science/article/pii/S1125786516300820

Jehl A, Bleron L, Meriaudeau F, et al (2011) Apport de la mesure de pente de fil au classement mécanique des bois de structures. In: CFM 2011-20ème Congrès Français de Mécanique, AFM, Maison de la Mécanique, 39/41 rue Louis Blanc-92400 Courbevoie

Klear M (2023) circle-fit: A python package for circular fitting. URL https://pypi.org/project/circle-fit/, licence: MIT License (Copyright (c) 2022 Michael Klear Permission

Kretschmann DE (2010) Mechanical properties of wood. Environments 5:34

Lang R, Kaliske M (2013) Description of inhomogeneities in wooden structures: mod-elling of branches. Wood Science and Technology 47(5):1051-1070. https://doi.org/10.1007/s00226-013-0557-4, URL https://doi.org/10.1007/s00226-013-0557-4
Longuetaud F, Pot G, Mothe F, et al (2022) Traceability and quality assessment of Douglas fir (Pseudotsuga menziesii (Mirb.) Franco) logs: the Tree-TraceDouglas database. Annals of Forest Science 79(46). https://doi.org/10.1186/s13595-022-01163-7, URL https://hal.science/hal-03658479

Longuetaud F, Schraml R, Mothe F, et al (2023) Treetrace spruce. https://doi.org/10.57745/WKLTJI, URL https://doi.org/10.57745/WKLTJI, dataset

Lukacevic M, Kandler G, Hu M, et al (2019) A 3d model for knots and related fiber deviations in sawn timber for prediction of mechanical properties of boards. Materials and Design 166:107617. https://doi.org/10.1016/j.matdes.2019.107617, URL https://www.sciencedirect.com/science/article/pii/S0264127519300541

Müller U, Gindl-Altmutter W, Keckes J (2018) How softwood tree branches are attached to stems: hierarchical extension of shigo's stem-branch model. Trees 32(4):1113-1121. https://doi.org/10.1007/s00468-018-1699-7, URL https://doi.org/10.1007/s00468-018-1699-7

Nyström J (2003) Automatic measurement of fiber orientation in softwoods by using the tracheid effect. Computers and Electronics in Agriculture 41(1-3):91-99. https://doi.org/10.1016/S0168-1699(03)00045-0, URL https://linkinghub.elsevier.com/retrieve/pii/S0168169903000450

Olsson A, Oscarsson J, Serrano E, et al (2013) Prediction of timber bending strength and in-member cross-sectional stiffness variation on the basis of local wood fibre orientation. European Journal of Wood and Wood Products 71(3):319-333. https://doi.org/10.1007/s00107-013-0684-5, URL https://doi.org/10.1007/s00107-013-0684-5

Olsson A, Pot G, Viguier J, et al (2018) Performance of strength grading methods based on fibre orientation and axial resonance frequency applied to norway spruce (picea abies l.), douglas fir (pseudotsuga menziesii (mirb.) franco) and european oak (quercus petraea (matt.) Liebl./Quercus robur l.). Annals of Forest Science 75(4):102

Olsson A, Pot G, Viguier J, et al (2022) Performance of timber board models for prediction of local bending stiffness and strength—with application on douglas fir sawn timber. Wood and Fiber Science 54(4):226-245

Penvern H, Demoulin L, Roux B, et al (2023) 3d fibre orientation reconstruction around a knot in douglas fir. In: 3rd ECCOMAS Thematic Conference on Computational Methods in Wood Mechanics (CompWood 2023), September 5-8, 2023, Dresden, Germany, International Center for Numerical Methods in Engineering (CIMNE), pp 137-138
Penvern H, Demoulin L, Pot G, et al (2024) A laboratory method to determine 3d fibre orientation around knots in sawn timber: case study on a douglas fir specimen. Wood Science and Technology 58(5):1735-1760. https://doi.org/10.1007/s00226-024-01583-w, URL https://doi.org/10.1007/s00226-024-01583-w

Pot G, Duriot R, Girardon S, et al (2024) Comparison of classical beam theory and finite element modelling of timber from fibre orientation data according to knot position and loading type. European Journal of Wood and Wood Products 82(3):597-617

Purba CYC, Viguier J, Denaud L, et al (2020) Contactless moisture content measurement on green veneer based on laser light scattering patterns. Wood Science and Technology https://doi.org/10.1007/s00226-020-01187-0, URL https://doi.org/10.1007/s00226-020-01187-0, aCL

Rais A, Bacher M, Khaloian-Sarnaghi A, et al (2021) Local 3d fibre orientation for tensile strength prediction of european beech timber. Construction and Building Materials 279:122527. https://doi.org/10.1016/j.conbuildmat.2021.122527, URL https://www.sciencedirect.com/science/article/pii/S0950061821002877

Simonaho SP, Palviainen J, Tolonen Y, et al (2004) Determination of wood grain direction from laser light scattering pattern. Optics and Lasers in Engineering 41(1):95-103. https://doi.org/10.1016/S0143-8166(02)00144-6, URL https://linkinghub.elsevier.com/retrieve/pii/S0143816602001446

Viguier J, Jehl A, Collet R, et al (2015) Improving strength grading of timber by grain angle measurement and mechanical modeling. Wood Material Science & Engineering 10(1):145-156. https://doi.org/10.1080/17480272.2014.951071, URL http://dx.doi.org/10.1080/17480272.2014.951071, aCL

Viguier J, Bourreau D, Bocquet JF, et al (2017) Modelling mechanical properties of spruce and Douglas fir timber by means of X-ray and grain angle measurements for strength grading purpose. European Journal of Wood and Wood Products 75(4):527-541. https://doi.org/10.1007/s00107-016-1149-4, URL http://link.springer.com/article/10.1007/s00107-016-1149-4, aCL

Zankel A, Wagner J, Poelt P (2014) Serial sectioning methods for 3d investigations in materials science. Micron 62:66-78. https://doi.org/10.1016/j.micron.2014.03.002, URL https://www.sciencedirect.com/science/article/pii/S096843281400050X