Robust numerical analysis of fibrous composites from X-ray computed tomography image data enabling low resolutions

X-ray computed tomography scans can provide detailed information about the state of the material after manufacture and in service. X-ray computed tomography aided engineering (XAE) was recently introduced as an automated process to transfer 3D image data to finite element models. The implementation of a structure tensor code for material orientation analysis in combination with a newly developed integration point-wise fibre orientation mapping allows an easy applicable, computationally cheap, fast, and accurate model set-up. The robustness of the proposed approach is demonstrated on a non-crimp fabric glass fibre reinforced composite for a low resolution case with a voxel size of 64 μ m corresponding to more than three times the fibre diameter. Even though 99.8% of the original image data is removed, the simulated elastic modulus of the considered non-crimp fabric composite is only underestimated by 4.7% compared to the simulation result based on the original high resolution scan.


Introduction
Image analysis of X-ray computed tomography images can become a standard part of the development process for anisotropic material design [1]. However, this requires faster approaches for setting up computational material models and a holistic concept including all process steps from image acquisition and reconstruction to segmentation and modelling. Consequently, to allow for efficient computed tomography based finite element model analyses of fibre reinforced polymer composites a number of challenges must be overcome. Among these, accurate mechanical performance predictions from large field of view image data stands out. There is an immediate need for methods that can provide reliable predictions of the response of complex and large composites structures to in-service loads across applications [2][3][4]. For this purpose, we propose a highly automated image based numerical analysis method and assess its sensitivity to image resolution.

X-ray computed tomography
X-ray computed tomography can improve the accuracy of finite element models by adding information about the fibre and bundle orientations in fibrous composites after manufacture. In contrast to classical material characterisation with imaging, e.g., scanning electron microscopy, the destruction of the sample can be avoided. The advantages of X-ray computed tomography for analysing key material parameters for complex microstructures have been emphasised in many studies [5][6][7]. It can be combined with tension, bending, and compression experiments, as well as impact and fatigue testing [8].

Fibre orientation estimation
There are two fundamentally different approaches to estimate fibre orientations from X-ray computed tomography images. The first approach is based on tracking the individual fibres, e.g., using the method by Emerson et al. [9] or commercial software such as Avizo™ [10]. Alternatively, local structural orientations can be estimated without tracking individual fibres, e.g., using structure tensor analysis [10][11][12] or, as recently proposed, digital image correlation [13]. These models can provide orientation distributions similar to those estimated using fibre tracking methods, while being both more robust when fibres are difficult to detect [10] and computationally much cheaper [11].

Modelling based on tomography data
Image based finite element models pose new challenging tasks such as 3D boundary segmentation and/or 3D image segmentation (voxelwise labelling) -tasks that are often demanding regarding software, hardware, and time. There are two common approaches for creating models based on tomography data [14]. One is a voxel-based method where every finite element represents a single voxel from the image data [15]. Also Straumit et al. [16] presented a voxel-based mesh, which models a quasi unidirectional composite containing flax fibres. Similar to this study, it determines the fibre volume fraction based on the X-ray computed tomography data, maps fibre orientations, and uses mostly data sheet values. Due to uncertainty of the natural fibres properties, the longitudinal fibre modulus is calculated by the experimentally determined modulus of the quasi unidirectional material. Voxel-based approaches often allow for a fast process but can easily lead to large model sizes and lack of an accurate boundary representation between objects in the scanned structure. The other option is based on a segmented physical boundary between objects [17,18].
Since the fibre orientations in composites are important to the material properties, information about local fibre orientation is critical to model accuracy. The orientations can be modelled with a sine function approximation [19,20], or by extracting a centroid line of the fibre bundle [21]. Ewert et al. [22] compared a segmented boundary-based mesh with a voxel-based mesh. Both approaches deliver similar stiffness results. It was found that a rough surface can lead to a non-physical centroid line. Hence, a sine function approximation or a centroid line based on a surface segmentation may be valid if no other information is available. However, X-ray computed tomography-based models should aim for more accurate solutions. Auenhammer et al. [23] have introduced the X-ray computed tomography aided engineering (XAE) framework to address image-based modelling as full process including image acquisition and reconstruction as well as segmentation. Despite the achieved accuracy, the proposed approach requires extensive computational resources and can only process limited model sizes.

Objective
Traditional 3D image based modelling faces challenges such as limited field of view, complex model creation processes, and high software and hardware cost which in turn lead to accuracy and efficiency issues. However, there is an extensive demand in research and industry to exploit the advantages X-ray computed tomography based modelling has to offer. Here, we aim to pave the way for easily applicable and accurate image based modelling of large and complex composite structures. Therefore, we present a novel open-access, automated, fast, precise, and robust approach that can benefit from the recent progress in Xray computed tomography and serve as a standard for advanced material development. The method is assessed in regard to sensitivity in elastic modulus and stress concentration predictions to image resolution for a textile reinforced composite.

Material
The non-crimp fabric reinforced glass fibre composite for this study was manufactured with E-glass fibres SE2020 from 3B-fibreglass infused by the epoxy resin RIMH035c/037 from Hexion GmbH [24]. The fibre diameter is 17 μm. The sample consists mainly of unidirectional bundles.
To provide stability before infusion they are held together by backing bundles and stitching threads. On the outside, a layer with ±45 • is added with the resulting layup [±45 • , b 80 • /0 • , b 80 • /0 • ] s . Here, b 80 • /0 • indicates a layer of backing bundles orientated at ±80 • , supporting the unidirectional bundles. The total thickness of the sample is approximately 4.5 mm with a unidirectional bundle thickness of approximately 0.7 mm. The linear density of the unidirectional bundles is 2400 tex while the backing and biaxial bundles have 68 tex and 600 tex, respectively. Further details about the material system can be found in Ref. [11]. This material lay-up was milled into tensile test specimens, specified in Ref. [25].

XAE -X-ray computed tomography aided engineering
XAE promotes a holistic view on all aspects of image-based modelling. Poor quality scans or reconstructions will drastically limit the quality of the final model. Further, the meshing, mapping, and modelling approaches are directly effecting the options for images acquisition.

Image acquisition
The 3D image data of the sample was acquired with a Zeiss Xradia Versa 520 X-ray computed tomography scanner. The main settings were an exposure time of 24 s, a number of projections of 3601, an acceleration voltage of 50 kV, an optical magnification of 0.4 × , and a voxel size of 8 μm. Three single scans were taken to increase the field of view resulting in a total scan time of 72 h.
Raw data from a 3D X-ray computed tomography scan are 2D projections consisting of line integrals. Therefore, they need to be reconstructed. Reconstruction is essential to X-ray computed tomography analysis, as the weakest link of all work steps defines the quality of the final numerical results. In the presented work, the Zeiss Xradia Versa 520 built-in reconstruction software using the Feldkamp algorithm for cone beams [27] is used to obtain the 3D images. In Fig. 1 a beam-cone effect from the scanning itself at the left and right ends of the scan is visible. A different reconstruction algorithm may suppress this type of artefact [28][29][30]. Alternatively, this part can easily be cropped away from the analysed region. The effects are described in the results chapter.

Automated segmentation
The simplest way to segment the image data is using thresholding to separate the 3D image into fibre bundles and matrix. This approach was previously used in Refs. [23,31], where it was combined with local orientation mapping. However, this type of segmentation does not allow for separated individual fibre bundle structures.
Automatic segmentation of individual bundles (instance segmentation) is challenging as bundles are often not clearly separated and individual fibres may move between bundles. Thus, even the human eye is sometimes not able to accurately separate the bundles. To segment the individual fibre bundles as non-overlapping meshes, this study uses an image segmentation approach based on sparse layered graphs [32]. In the following, both the thresholding and graph-based bundle segmentation methods are investigated.

Thresholding segmentation
With thresholding segmentation, objects in an image are separated based on their grey-scale value. In the investigated case, glass fibres and resin show distinct grey-scale values, and thereby can be easily segmented. The thresholding segmentation in this study was carried out with Avizo™, but there are many other options that offer a thresholding operator. Avizo™ offers the advantage that a geometrical surface in a mesh form can be generated.

Individual bundle segmentation
The graph-based method for segmenting individual bundles works by radially sampling the image intensities along each centre line orthogonal to the line, which are used to create a new "unfolded" volume for each fibre bundle. Then, the individual bundle surfaces can be segmented with binary optimisation using the approach by Li et al. [33], which uses graph cuts. This approach was extended by Jeppesen et al. [32] to ensure that the surfaces do not overlap.
A number of parameters were used to ensure that the bundle surfaces were accurately segmented. The set of parameters is related to the radial sampling and includes the number of sampling angles, the number of steps per angle, and the radii at which to start and stop the sampling. These parameters determine the sample resolution and are chosen based on the desired resolution of the surface meshes and available system resources. For the segmentation in this study, 120 angles were sampled with 100 steps per angle at 500 steps along each line, resulting in a quad surface mesh with 60,000 vertices per bundle. To ensure accurate surface segmentation, a surface smoothness was applied following the approach of Li et al. [33] and a minimum exclusion margin following the approach of Jeppesen et al. [32] avoiding bundle overlap. The slgbuilder Python package was used to construct the energy function and obtain the solution.

Meshing
Two different models are prepared. One is set up with a mesh, based on the thresholding segmentation. It is used for stiffness and stress investigations, while the mesh, based on the individual bundle segmentation, is only used for stress investigation as, for simplicity, the biaxial and backing bundles are not included.

Thresholding based mesh
The surface, based on thresholding segmentation created with Avizo™ was imported as mesh into the finite element pre-processor ANSA™. There, the surface mesh quality was improved and the average element size was set to 200 μm. Afterwards, a solid mesh for the fibre bundles, as well as for the surrounding resin rich areas, was created with tetrahedral elements. The generated finite element model comprises approximately 2.9 million elements. Compared to the previous study [23], the model volume size was increased by approximately 50%, while the element number was kept nearly constant. An automated meshing process is regarded as too difficult at this stage, given the complexity of the fibre bundle structure. Therefore, the meshing process including manual operations takes approximately 2 h for an experienced operator.

Individual bundle based mesh
In the alternative model, the quad surface mesh provided by the individual bundle segmentation approach was further smoothed to make it more suitable for use with Abaqus™. A Python script was developed that converts the quad mesh to a mesh file with an Abaqus™ syntax. The script also includes an adjustable element number reduction. The mesh was imported in ANSA™ and further manipulated. It remains for future work to remove the mesh generation in a finite element pre-processor and fully automate the 3D meshing in a Python script. A tetrahedral mesh could be generated very easily in ANSA™. However, for this study it was decided to take advantage of the segmented individual bundles and create a hexahedral mesh. However, creating a hexahedral mesh is challenging in this case as the structure is still complex. Two different meshes with an average length of 100 μm and 200 μm, respectively were created resulting in approximately 2.3 million elements and 0.5 million elements, respectively. An experienced operator needs approximately 2 h for the meshing of one model (displayed in Fig. 2).

Structure tensor analysis
Structure tensor analysis 1 provides a fast and accurate way to estimate local fibre orientations in fibre-reinforced composites [10,11]. In the context of 3D images (volumes), a structure tensor is a 3-by-3 matrix that captures the local orientation and potentially the anisotropy around a point in space. The structure tensor, S, can be computed for a point in a volume, V, as T is the gradient of V and the summation/ integration is limited to a certain neighbourhood around the point. This study uses the algorithm by Jeppesen et al. [11,12], specifically the structure-tensor Python package, where Gaussian kernels are used both for calculating the gradients and integration. Once S has been calculated for every voxel in the volume, the dominant orientation of fibre-like structures can be found analytically by eigend ecomposition of S. This yields three eigenvalues and three mutually orthogonal eigenvectors for each point at which S is calculated and decomposed. The eigenvector v, corresponding to the smallest eigenvalue, points towards the direction of least change in pixel intensities. Thus, for any points where fibre-like structures are visible in the image, v describes the fibre orientation at that point. In other words, the local fibre orientations can be estimated by calculating v for all voxels that contain fibres.
This structure tensor analysis approach takes two parameters, σ ST and ρ ST [11]. The σ ST parameter is known as the noise scale and represents the standard deviation of the first-order derivative of Gaussian The values are calculated with the fibre radius r fibre given in voxels as: The structure tensor code supports parallel computation on graphics processing units which decreases the computation time further (in many cases by over an order of magnitude), allowing very large volumes to be analysed quickly.

Mapping
In the mapping step, the estimated local fibre orientations are assigned to the finite element model. For that, the first eigenvectors of the structure tensor computations are taken. The mapping code itself is combined with the structure tensor orientation analysis, as the structure tensor analysis is independent on the mesh generation.
Two different mapping schemes have been established. Both can map the material orientations on continuum (2D, 3D) and structural (beam, shell …) elements for first and second-order formulation. The two schemes differ in their mapping destination. One maps the material orientation element-wise, while the other one allows for an integration point-wise mapping.
The Python code in a Jupyter notebook for the element-wise scheme is designed for an Abaqus™ mesh file syntax. The first step is to calculate the centre of gravity for each element. Next, the estimated local orientation of the voxel position closest to the centre of gravity is extracted (orientations were estimated in form of the eigenvector of the structure tensor with the lowest eigenvalue for each fibre voxel position). Then, the extracted local orientation is assigned as the local element orientation. The choice of representing the material orientation of an element with an average length of 200 μm by the orientation of a voxel with the size 8 μm is only admissible since the structure tensor analysis summarises the orientation of a region around the given position (determined by ρ ST ). The entire mapping process, including the import of the mesh file and export of the new mesh file with the mapped orientation, takes approximately 3 min for the specific model size.
The integration point-wise mapping works in a similar fashion; only the mapping target is different. Therefore, the integration point coordinates need to be supplied. A dummy simulation is used to output the integration point coordinates as an Abaqus™ DAT-file. The mapped orientations are then output in a Fortran 77 data file designed for an efficient compilation in the Abaqus™ subroutine ORIENT. In the subroutine, the integration point-wise orientations are assigned. An example of the mapping is visualised in Fig. 3. The time required for the integration point-wise mapping is approximately the same as for the element-wise mapping. The Jupyter notebook that combines the structure tensor analysis and orientation mapping, along with further details, is available in Ref. [36].

Material modelling
The simple load case, a uniaxial tensile test, allows for a linear elastic material model for both the fibre bundles and the resin rich areas. The material orientation of the fibre bundles is reflected by an orthotropic material model. The required nine stiffness parameters were calculated with the micro-mechanics model by Chamis [37]. The determined material parameters are listed in Table 1. The unidirectional and the biaxial bundles are made from the same fabric while the backing bundles consist of different fibres.
The required fibre volume fraction inside the fibre bundles is calculated according to Equation (3), where W FAW represents the fabric areal weight, n the number of bundle layers across the structure thickness, ρ fibre the glass fibre density, χ bundle the volume of the fibre bundles, and A surf the area of the scanned surface parallel to the main fibre direction. W FAW and ρ fibre can be extracted from the data sheet supplied by the material supplier, while A surf and χ bundle are determined based on the segmentation outcome. The fibre volume fraction inside the unidirectional fibre bundles for the investigated sample was determined to be 55%.

Results
The model, with its stress state results shown in Fig. 4a  models are analysed for an uni-axial tensile strain of 0.25% with constrained displacement at one end. The models are simulated with the commercial finite element solver Abaqus™.
The presented model is based solely on computed tomography image analysis and data sheet values. Consequently, the simulated elastic tensile modulus is sensitive to the thresholded bundle volume and the fibre orientation analysis. There are other options to predict the elastic tensile modulus which are less time consuming and possibly more accurate. The strength of XAE lies in the qualified stress prediction of real samples. Hence, the primary goal for the elastic tensile modulus analysis is not to present a new tensile modulus prediction method but to use it as a measure to compare different input image quality and fibre orientation analysis parameters. However, the result of the simulated elastic tensile modulus is in reasonable agreement with experimentally determined tensile moduli of twelve specimens from the same material batch with a mean value of 35.2 ± 0.3 GPa. The physical specimens were tested according to Ref. [39]. The simulation results predict a tensile modulus 3.3% lower than the experimentally measured moduli. It has to be kept in mind that an XAE result reflects the situation in one single part and is therefore limited to predict a mean tensile modulus value for several specimens.
In the following, different parameters that influence the stress distribution and tensile modulus of the numerical results are discussed in detail. Table 2 provides an overview of the different settings. The outcome of the simulation results is quantified with the tensile modulus of the full model and the mean stress in the local coordinate system, relating to local fibre orientation only for the unidirectional bundles. Additionally, the median of the absolute angle between the fibres in the unidirectional bundles and the z-axis (nominal fibre orientation) is used as a measure of the fibre misalignment [12]. This angle is usually noted as θ [40].
For an accurate and fast XAE process, the image acquisition, image analysis, and modelling need to be aligned. This directly influences the choices for resolution, field of view, mesh size, element formulation, and structure tensor analysis parameters. The entire parameter study was deliberately set up with the intention to identify and push the boundaries of the XAE process rather than focus on the most accurate solution.

Structure tensor parameters
In this analysis, the values for ρ ST and σ ST were set 25% smaller and larger than for the reference case. A comparison of the stress state in a cross-section along the specimen length direction (i.e., parallel to the unidirectional bundles) is shown in Fig. 4 (top row). The depicted crosssection lies at a boundary layer between backing bundles and unidirectional bundles, where stress concentrations become visible. To further quantify the stress state, histograms with a normalised integration point count were created (Fig. 4 bottom row) for the stress distribution in the unidirectional bundles. Both histograms and cross-section stress plots show little difference in the stress results for variations of σ ST and ρ ST . Although the values of ρ ST and σ ST were significantly alternated, neither major differences in the stress state nor tensile modulus (Table 2) were found. This demonstrates the advantage of using this structure tensor analysis in the XAE process.

Mapping scheme
So far, only results for second-order tetrahedral elements with an integration point-wise mapping have been presented. As described in the previous chapter, also an element-wise mapping scheme has been developed. Here, the standard segmentation and simulation results are compared with first-order tetrahedral elements equipped with an  Table 1 Calculated material properties of the fibre bundles by the model by Chamis [37], where the coordinate 1 is the local fibre direction while 2 and 3 is the two transverse directions. The fibres in the unidirectional and biaxial bundles are given with E f = 81 GPa, while the fibres in the backing bundles are listed with E f = 72 GPa. The other constituent properties are E m = 3. 16 GPa, ν f = 0.22, and ν m = 0.366 [24]. element-wise mapping scheme. The material orientation analysis was the same in both cases. As for the integration point-wise mapping scheme, a second-order element formulation is also possible for the element-wise mapping scheme, but not considered here due to the focus on fast computation. The results from the integration point-wise and element-wise schemes are very similar, as shown in Fig. 5. The stress distribution of the element-wise scheme with second-order elements appears smoother. The elastic tensile modulus computed with the element-wise mapping scheme differs by only 0.6%. Moreover, the stress distributions in the histograms are very similar.

Resolution
To investigate the influence of variations in resolution, new 3D images with reduced resolution were artificially created. For these, a binning of 2, 4, and 8 was applied. A binning of 2 means that 2 voxels in all three spatial direction are combined, corresponding an averaging of the intensity in 8 voxels, which leads to a data reduction by the factor 8. Consequently, a binning of 4 reduces the size by 64 and a binning of 8 even by a factor of 512. The binning of 2 can be also expressed as a voxel size of 16 μm, whereas a binning of 4 and 8 translates to voxel sizes of 32 μm and 64 μm, respectively. With a voxel size of 32 μm and 64 μm, it is impossible for the human eye to detect single fibres with a diameter of 17 μm. This is contradictory to other attempts in research where high resolution is desired to ensure a correct fibre orientation analysis. This study aims to go in the other direction and find the minimal required resolution to save computation time and to find a possible measure for a field of view enlargement.
In Table 3 the file size reduction as well as the information loss due to increased voxel size are specified. Moreover, the computation time on a computer with 32 GB RAM and an Intel i7-8750 central processing unit with six cores for the material orientation analysis with the structure tensor code is listed. Size reduction is found to drastically lower the computational cost even though it does not directly scale with file size.
The reconstruction artefact stemming from the beam cone effect (Fig. 1) was enhanced at lower resolutions. Therefore, it was decided to remove the top and bottom domain of the model. The model length for this investigation was then 25 mm instead of 34 mm. The obtained changes for tensile modulus and mean stress inside the unidirectional bundles for the reference model are insignificant, which indicates that the beam cone effect does not have an influence for the 8 μm voxel size model.
With larger voxels the orientation distribution over a larger domain is represented. This means that, for instance, elements close to the boundaries of a unidirectional bundle can be assigned an off-axis fibre orientation. In the case of 64 μm voxel size this led to unacceptable results following Equation (2) for the structure tensor parameters ρ ST and σ ST . The situation was improved by deviating the structure tensor parameter (specified in Table 2). However, the effect of the smeared orientation information persists for all three lower-resolution cases and has two consequences.
Firstly, the tensile modulus is reduced ( Table 2) from 34.0 GPa to 32.8 GPa, and the mean stress in the unidirectional bundles decreases from 112 MPa to 109 MPa for the standard run compared to the case of 64 μm voxel size. In Fig. 6 (top row) this can be observed by stress plots that go more towards green colours. Secondly, the smeared orientations (due to the lower resolutions) lead to a less distinct stress distribution. The peaks in the stress histograms in Fig. 6 (bottom row) flatten as the voxel size increases. Consequently, the reduced resolution results in a lower modulus and mean stress as well as smoother stress distribution. Nevertheless, the concentrated high and low stresses in the unidirectional bundles as well as in the backing bundles, remain visible for all four simulation results (Fig. 6 top row). Even though the backing bundles have an approximate thickness of only 100 μm, the structure tensor analysis and the mapping algorithm prove their ability to capture orientations with the large voxel size of 64 μm. In terms of stiffness, the simulation run based on a voxel size of 64 μm underestimates the tensile modulus by 4.7% and the mean stress in local fibre direction by 2.9% compared to the 8 μm simulation run. Such errors can be acceptable in many cases, in particular considering the advantages that come with a low resolution in terms of scanning and computation times or increased field of view. A voxel size enlargement to 16 μm is nearly unrestrictedly admissible while for 32 μm and 64 μm a more detailed assessment must be performed.

Mesh size and type
In the model built by first-order hexahedral elements (depicted in Fig. 2) for the bundles, neither backing nor biaxial bundles are included. However, the stress concentrations in the unidirectional bundles near the backing bundles are still clearly visible (Fig. 7). This supports the conclusion from Ref. [31] that the backing bundles themselves are not responsible for the stress concentrations, but the bundle waviness that they induce. Furthermore, the stress values closely follow the results from the full model with second-order tetrahedral elements including the backing bundles. The mean value of the stress in the local fibre direction in the unidirectional bundles differs at most by only 0.4 MPa compared to the standard model computed with second-order tetrahedral elements. No major differences between the results from the two models with an average element length of 100 μm and 200 μm, respectively, can be identified. This is despite the fact that the former consists of more than 4 times more elements.

Discussion
From an image analysis point of view, segmenting non-crimp fabric reinforced fibre composites with fibre bundles in contact over large areas is a challenging task. Due to the parallel orientation of adjacent bundles, segmentation of individual parallel fibre bundles is more difficult for non-crimp fabric than woven composites, where neighbouring bundles are oriented differently [41,42]. Even so, the sparse layered graph-based segmentation method was successfully applied to non-crimp fabric reinforced composites for individual bundle segmentation. The segmented individual bundles allow for a hexahedral mesh, albeit the meshing task is challenging. The investigated load case, however, does not justify the time-consuming hexahedral meshing process. Nevertheless, it is believed that the individual bundle modelling approach can be advantageous for other load cases, e.g., compression load.
Two newly developed material orientation mapping schemes have been presented. The integration point-wise scheme is specifically designed for Abaqus™, while the element-wise scheme can be also used for other solvers, e.g., LS-Dyna™. Even though the element-wise mapping scheme can also process higher order elements, it is recommended to use the integration point-wise scheme to fully take advantage of the Table 2 Overview of the simulation runs where all models were built on the thresholding segmentation method and meshed with tetrahedral elements with an average element length of 200 μm. The values of the mean σ 11 in local fibre direction only refer to the stresses in the unidirectional bundles, while the fibre misalignment represents a direct consequence of the image analysis parameters and the chosen voxel size. The top row specifies in which figure the single models are shown. In the second column the standard model is listed for reference. Marked in darker grey are the parameters that were actively changed compared to the standard run (listed in the second column). Marked in lighter grey are the parameters calculated from Equation (2).

Fig. 5.
Comparison of the influence of the mapping scheme on the outcome of the stress distribution in a y-z cross-section at the interface between backing bundles and unidirectional bundles (top row) and the stress distribution in the unidirectional bundles depicted in histograms (bottom row). The element-wise mapping simulation was run with first-order elements, while the integration point-wise simulation was run with second-order elements. method in case higher order elements are desired. Both schemes can be used for continuum as well as structural elements and are open access.
The implementation of the structure tensor method in the XAE process provides a fast and accurate material orientation analysis. Together with the novel material orientation mapping schemes, it allows for a precise prediction of the tensile modulus and stress distribution inside the fibre bundles. It is open access and easy to apply as it only requires two parameter settings. Furthermore, the insensitivity to the choice of these parameters on the stiffness and stress levels has been shown. Moreover, the method allows for a reduction of the required resolution.
The most extreme case with a voxel size of 64 μm may still produce quality results compared to the original image with a voxel size of 8 μm. This is inline with results by Karamov et al. [10] achieved with their structure tensor code, applied on images of glass fibres with a voxel size of 16 μm. Aiming for lower resolutions provides new opportunities.
Firstly, since 99.8% of the data size of the original image data is removed, the computation time is dramatically reduced. Compared to the originally used fibre tracking in Avizo™ in Ref. [23] with approximately 2 h on a special machine with a large graphics processing unit, this analysis runs in a few seconds on a laptop central processing unit.
Even the original image data with a voxel size of 8 μm can be run on a laptop central processing unit in approximately 20 min. Secondly, another possible direction is to exploit the fibre orientation analysis and mapping capability to increase the field of view of the scan. With a voxel Table 3 Effect of artificial down-scaling of the original X-ray computed tomography image scanned with a voxel size of 8 μmm on the file size and the computation time of the structure tensor orientation code run on a computer with 32 GB RAM and an Intel i7-8750 central processing unit with six cores.  Fig. 6. Visualisation of the 3D image resolution influence on the predicted stress distribution in a y-z cross-section at the interface between backing bundles and unidirectional bundles (top row) and the stress distribution in the unidirectional bundles depicted in a histogram (middle row). Examples for the image resolution are shown in the bottom row. Note: Due to the amplification of the beam cone artefact with lower resolutions the model length was reduced from 34 mm to 25 mm. size of 64 μm the field of view can theoretically be increased by a factor of 512. A shift by two orders of magnitude in the possible size of the scanned volume can revolutionise X-ray computed tomography-based modelling. The presented models consist of approximately 2.9 million elements. An increase by two orders of magnitude would require an unrealistic number of 290 million elements. This leads to a situation where the X-ray computed tomography scanning is not the limiting factor. In other words, X-ray computed tomography analysis can provide information about fibre orientation on a scale and field of view that cannot be achieved with common meshing and modelling approaches. Consequently, a demand is created for new methods to benefit from the possibilities provided by the advances in computed tomography and subsequent image analysis. As non-crimp fabric reinforced composites are considered as one of the most challenging fibrous composite layup for image analysis, it is believed that the presented method can be used for other types of composites, such as short fibre injection moulded and woven composites. Moreover, the method was deliberately designed to be widely applicable and serve as basis for more complex load cases, as the code is published open access and all numerical analysis of fibrous composites can benefit from accurate fibre orientation prediction.

Conclusion
The main goal of this study is to present a robust method for numerical analysis of fibrous composites from X-ray computed tomography that allows for low resolution 3D image data. The implementation of the fast, computationally cheap, and robust structure tensor method for material orientation as basis for the newly presented mapping schemes has proven insensitive to parameter settings and image resolution regarding simulated tensile stiffness and stress distributions. The investigation for low resolution 3D images has even moved the boundaries for modelling based on X-ray computed tomography data. In the demonstrated case, the field of view can be increased by more than two orders of magnitude, with practically maintained accuracy in numerical prediction of the elastic modulus of the composite structure.