Data assimilation of structural deformation using finite element simulations and single-perspective projection data based on the example of X-ray imaging

Abstract

Gaining insight in the complex deformation processes inside a vehicle during the process of crashing is to date only possible via numerical simulations. To validate these simulations experimentally, high speed, high-energy X-ray imaging is currently developed. So far, X-ray images have only been compared qualitatively to simulations as quantitative analysis proves challenging: Three dimensional (3D) information has to be reconstructed from two dimensional (2D) X-ray images showing overlapping parts which undergo complex deformation. We propose a novel analysis method which adapts ideas from data assimilation to reconstruct three dimensional motion and deformation from two dimensional measurements using high speed X-ray video data and finite element (FE) simulations. Furthermore, we analyze the accuracy of the approach based on an observing system simulation experiment.

Introduction

Predictive simulations play an increasingly important role in research and development as well as in acceptance testing, especially in automotive applications. Not only structural deformation and material collapse is of interest for recent research, e.g. in the context of weight minimization (e.g., [19, 25]) or crashworthiness optimization (e.g., [22, 26]) of certain vehicle components. Also, the simulation of human body dynamics is of huge importance to understand injury risks during car accidents concerning e.g. internal organs [30] or the head [17]. Finite element (FE) car crash simulations are commonly validated by quantitative and qualitative comparison to experimental crash tests. Many simulation parameters, especially material parameters and models are determined in separate material characterization tests. Other parameters are set as boundary conditions to reproduce a certain experiment. Quantitative analysis focuses on comparing acceleration and impact load curves at the location of test instrumentation in the crash test (e.g., [20]), whereas qualitative evaluation is based on optical comparison of the global deformation in the impact zone [21]. Usually, deformation processes visible on the outside of the vehicle are analysed using high speed imaging systems and methods from photogrammetry (e.g., [27]). However, to date it is not possible to quantitatively validate the deformation processes hidden from optical imaging taking place inside the car during the crash. X-ray technology allows gaining an insight to these hidden internal structures. Using high-energy X-ray computed tomography, it is possible to statically image entire vehicles in three dimensions [29]. However, this process takes several hours and is therefore not feasible to image dynamic deformation processes during crash happening on a timescale of milliseconds. High speed X-ray imaging has proven useful in the context of car crash experiments, however, due to experimental constraints often only a single imaging perspective is available. [1] have successfully acquired and used high speed X-ray videos to investigate engineering processes. [28] apply high speed X-ray video analysis to investigate cervical spine injuries in the context of rollover crashes. Changes to the shape of cervical spine are determined qualitatively by tracking motion of the neck on the X-ray images. The X-ray data generated in this work is used by [33] to validate their human FE model. By qualitatively comparing corresponding frames of the human model simulation with the experimentally acquired X-ray images, they can confirm the representativeness of the overall kinematics around the cervical spine. Our proposed method enables quantitative analysis of structural deformation processes visible on a sequence of 2D single-perspective (or multi-perspective) X-ray images which can be used for validation of predictive crash simulations. To do so, quantitative, time-resolved evaluation of three dimensional motion from two-dimensional X-ray images is necessary. To fill the gap opened by the missing dimension of the experimentally acquired data, including a-priori knowledge of the observed processes is necessary. This approach has been used in a wide variety of research fields: The project XROMM reconstructs skeletal motion of animals from biplanar in vivo X-ray videos [5]. For each video frame, three-dimensional bone models are registered with X-ray images from two perspectives resulting in an accurate animation of the skeletal motion. The animal bones in XROMM are assumed to remain rigid during movement. This assumption acts as a-priori information filling the information gap and enabling the reconstruction in three dimensions. Vehicle components, however, can deform drastically during crash. Thus, for analysis of deformation processes, this approach is not viable. [23] acquire high speed X-ray images of fragments after a high speed impact incident. Using flash X-ray technology few images are acquired at different points in time and from different perspectives. Using assumptions on fragment motion, the 3D flight path of fragments is computed. In the field of medical imaging, particularly for image-guided radiotherapy the ability to track moving objects in real time on X-ray images is crucial. [7] track tumour motion by registering X-ray images to a planning 4D-CT using the structural similarity measure on digitally rendered radiographs from the 4D-CT and the X-ray images acquired during radiotherapy. From a set of volumetric training data [16] compute deformation vector fields between the separate respiration phases using deformable image registration. Using principle component analysis (PCA) a set of eigenvectors and PCA coefficients is determined to model respiratory motion. To gain 3D information from a single X-ray projection during therapy, the parameters of this model are optimized to fit the X-ray image. In this way, the approach of model-based image registration uses a-priori information inherent to an underlying dynamic model to reconstruct three dimensional processes from 2D images. In the context of structural deformations and material collapse, especially in the field of automotive crash testing, FE simulations are commonly used to model the system dynamics. To be able to draw inferences from two-dimensional measurement data on three-dimensional deformation processes, the a-priori information inherent to the dynamic model has to be used, as the two-dimensional X-ray image yields information on the deformation of an object only within the imaging plane. However, deformation processes within the projection axis remain hidden to the observer. In regular X-ray analysis of rigid bodies, depth information in the viewing direction can be derived using the known size of the objects and the fixed magnification factor of the cone beam imaging system. However, this approach fails when the tracked features or markers are too small for a local magnification factor evaluation and the deformations invalidate the rigid body approach. Therefore, this information gap has to be filled using known dynamic correlations between processes happening within the imaging plane and parallel to the projection axis. Generally, these correlations are not inherent to the geometric results of numerical simulations: Geometrically morphing the simulated object to fit the measured X-ray image does not generally yield reliable results in the third dimension, where only limited information is available from the image. In fact, the correlation between the visible and hidden deformation is inherent to the physics-based model underlying the numerical computations. Therefore, the proposed method adapts parameters of the FE simulation to match the simulated geometry with its X-ray projection whilst taking into account the physical correlations in the third, non-measurable dimension. Note that this publication focuses on the methodology of making the physics of the FE simulation accessible as a-priori information for this matching and adaption process. We investigate the method using a simplified test case, where the parameters chosen for adaptation represent generic materials for demonstration purposes. For more complex and realistic simulations, the choice of parameters for modification and the respective confidence intervals remain task of a simulation expert in order not to lose the predictability and reliability of already carefully calibrated simulations. Combining measurement data from X-ray images and information from FE simulations has been investigated in a medical imaging context: Here, FE simulations of the compression process of the female breast during mammography are used to register 3D magnetic resonance scans and 2D X-ray images: [15] optimize unknown input parameters of the FE-simulation considering patient examination setup, whereas [10] optimize patient specific parameters of the FE model of the patient’s breast. These works aim at adapting a three-dimensional FE model in order to match the deformed state visible on acquired X-ray images. However, in the context of mammography, only the final deformed state of the breast is relevant. In contrast to this, many applications require insight in a time-dependent deformation process. In the specific case of car crash testing, the final state is well known through disassembly of the crashed car. However, from this, it is not possible to draw inferences on the high speed dynamic processes that yielded this final state. The task of optimizing highly complex dynamic models based on spatially as well as temporarily discrete measurements is also addressed in the field of data assimilation. This conception originates from the research area of weather forecasting. Here, spatially and temporally incomplete and error-prone data from weather stations is fed into a running simulation of the current weather state, aiming at improving the precision and realism of the forecast (for an overview on different assimilation techniques and their applications refer to e.g., [24, 31]). The potential of transferring the methods of data assimilation from numerical weather forecasting has been highlighted for different fields of application: The great potential of applying the methods of data assimilation for prediction of tumour growth using biophysical models and imaging measurements is emphasised by [32]. [13] propose an ensemble approach to assimilate magnetic resonance images of the human brain with mathematical models describing brain tumour growth. In the wider context of engineering, the methods of data assimilation are also used to adapt general physical models based on measurement data in the context of the concept of ’gray box modelling’ [6]. More specifically, data assimilation has also been successfully applied on FE simulations as dynamic model. [18] apply the Kalman filtering technique to achieve multimodal medical image registration using a finite elemente formulation of a dynamic linear elastic continuum mechanics model for deformation of soft tissues. The method is found to be promising for image-guided procedures where tissues undergo deformation. Apart from the medical sciences, applying data assimilation techniques to adapt FE simulation has also proven useful in completely different fields: The 3D-Var methodology is used to recover ideal FE parameters of a simulation of an electrical machine [3]. Ebrahimian et al. [8] demonstrate a successful application of Kalman filtering to find unknown parameters of FE models simulating the response of a structural system to earthquake ground motions. We propose transferring the aforementioned methods to the field of automotive crash testing: Currently, a fast high-energy radiography system is established, which enables X-ray imaging of relevant vehicle parts during crash experiments at a frame rate of 1 kHz [14]. Using the methods of data assimilation we augment the two dimensional observations from X-ray measurements using the information inherent to the numerical FE simulation. To do so, we adapt chosen parameters of the FE simulation using a Kalman filtering technique, yielding an optimized simulation matching the experimental measurements. In this proof-of-concept simulation study paper, we present first developments of this analysis technique on the basis of a simple structural test case. The paper is organized as follows: the simulated test object and the developed analysis technique are explained in Section 2. Results of the simulation study are displayed and discussed in Section 3. Section 4 concludes the paper and gives an outlook on possible future developments.

Methods

FE simulation of test case

Generally, an FE simulation aims at solving a system of partial differential equations, consisting of the equations of motion, the equations of state, the constitutive equations as well as boundary and initial conditions for a certain object geometry. To be able to set up the equations of state as well as the constitutive equations, material properties of the object need to be considered (for a more detailed explanation on the finite element simulation technique refer to e.g., [4]). Note that FE simulations are always approximations relying on the validity of the underlying equations as well as of the input quantities regarding object geometry, material properties, initial and boundary conditions [9]. Assuming the underlying equations describe the physical behaviour of the object sufficiently well, the FE simulation can be calibrated to match the corresponding real process by adjusting these input quantities.The proposed parameter estimation method aims at adjusting a preselected set of input quantities in order to calibrate an FE simulation of a deformation process using projection data. In order to investigate and evaluate the applicability of the proposed parameter estimation method we chose to carry out a simulation study, such that the results can be compared to a well-defined and known ground truth as opposed to a real-world experiment where deviations between simulated and experimentally observed deformation behaviour are less controllable. Since the proof of concept itself is not dependent on the complexity of the simulation, a simple test case was chosen consisting of a square rod which is bent and rotated simultaneously (see Fig. 1). The physical dimensions of the rod in the simulation amount to $100\text{cm}$ x $10\text{cm}$ x $10\text{cm}$. Note that these dimensions have been chosen to resemble the dimensions of a generic component considered in a car crash experiment. The simulations were performed using the LS-Dyna numerical simulation framework (common in the field of automotive crash testing) over a time frame of $1.0\text{ms}$ with a step size of $0.025\text{ms}$ (explicit analysis). Note that these time scales are generic values not intended to represent a realistic test case. Provided the deformation process of interest can be observed for a sufficiently high number of time steps, the presented analysis method works on any time scale. The rod was modelled using 10000 cube-shaped solid elements (1 cm x 1 cm x 1 cm) and an elasto-plastic material formulation (LS-Dyna *MAT_PIECEWISE_LINEAR_PLASTICITY). A bending force of 1000 N is applied to a node set at both ends of the rod, a node set in the middle of the rod is fixed in space acting as central bearing.

Fig. 1 [Images not available. See PDF.]

FE simulation (LS-Dyna) of a square rod, which is bent and rotated simultaneously. Nodes used for the analysis are marked with their corresponding node IDs

Two simulations with well-defined differences in three chosen parameters (density $\rho$ and Young’s modulus E of the rod, angular velocity $\omega$) are generated. Due to these differences, the resulting spatial trajectories of the nodes on the rod differ in all three dimensions. The suggested method is evaluated in an observing system simulation experiment: The first FE simulation is used as “simulation”, whereas the second FE simulation represents a fictional experiment or “truth run”. Based on this truth run, synthetic measurement data is generated by computing projections of the 3D data on a virtual X-ray detector. To evaluate the influence of measurement uncertainties, Gaussian noise is added to the projected coordinates. These virtual measurements are equivalent to measuring the trajectory of certain existing features on a part (e.g. corners) or X-ray markers introduced to featureless areas (e.g. flat surfaces) by using feature tracking methods on experimentally acquired X-ray images. Figure 2 visualizes the rod geometry at $t_0=0\text{ms}$ and $t_{\mathit{end}}=1.0\text{ms}$ as well as respective simulated X-ray projections. Note that the X-ray projections were added for visualization purposes, the virtual measurement data is computed directly from the respective simulation.

Fig. 2 [Images not available. See PDF.]

3D visualization (top) and respective simulated X-ray projection (bottom) at $t_0=0\text{ms}$ (left) and $t_{\mathit{end}}=1.0\text{ms}$ (right)

State-space model for parameter estimation

Figuratively speaking, the aim of the proposed method is to reconstruct the true 3D position of designated simulation nodes based on their projected position on the 2D X-ray image and the position predicted by the simulation, taking into account measurement and model uncertainties. In order for the reconstructed node position to match the corresponding measured projection, its location is constrained only in two dimensions (detector plane) with a remaining degree of freedom along the connecting line between X-ray source and measured projection (see Fig. 3).

Fig. 3 [Images not available. See PDF.]

Schematic visualizing the proposed approach: Considering a certain 3D node coordinate of the simulation and comparing its projection to the X-ray measurement on the detector reveals a deviation between simulation and reality. By taking into account uncertainties of measurement and FE model, an assimilated coordinate is computed using the Unscented Kalman filter

One approach to locate the node along this line would be to directly infer its position using Bayesian state estimation on the considered node. In doing so, the simulated object is geometrically (and unphysically) morphed to fit its projection without including any information regarding the correlation between the visible deformation of the object (parallel to the detector plane) and the hidden deformation processes (along the line of sight). However, assuming the object’s initial geometry is modelled accurately, deviations between real and simulated 3D node position have to be caused by differences in the deformation processes yielding these node positions. Therefore, the physics-based model underlying the numerical simulation of the deformation process has to be considered. By adapting the dynamic model of the deformation process itself, such that it yields a 3D object geometry that matches the measured projection, we can remove the remaining degree of freedom from the location of the considered node whilst accounting for the deformation mechanics causing differences between model and experiment. For this reason, we propose a parameter estimation approach to adapt the numerical simulation of the deformation process based on the measured projection data. The FE simulation of a complex deformation process generally takes a high number of model parameters and initial conditions as input quantities. For the parameter estimation, a subset of these quantities which is presumed to correlate with the deformation of the considered component is chosen (the choice of parameters is further elaborated in Section 3). The 3D geometry of simulation nodes at estimation step k, $x_k$, can be written as a function of the object’s initial geometry $x_0$, the simulation input parameters ${\theta }_k$ and p and the simulation time t. Here, ${\theta }_k$ denotes the set of estimated parameters at estimation step k and p denotes the set of all other input quantities of the FE model F:

1

$$\begin{array}{c}\hfill x_k=F(x_0,{\theta }_k,p,t).\end{array}$$Simulation time t and estimation step k are related via $t=k\ast s^{-1}$, with s the measurement sampling rate. Therefore, to determine the 3D node geometry resulting from a certain set of time-invariant parameters, the simulation is computed from $t=0$ to $t=k\ast s^{-1}$. The resulting 3D node geometry is related to the observed projection data $y_k$ via the central projection $h_{\mathit{proj}}$:

2

$$\begin{array}{c}\hfill y_k=h_{\mathit{proj}}(x_k,g),\end{array}$$where g denotes the X-ray setup geometry (i.e. source position, detector position and orientation). The state-space model for the parameter estimation is therefore given by:

3

$$\begin{array}{cc}& {\theta }_{k+1}={\theta }_k+w_k\hfill \end{array}$$

4

$$\begin{array}{cc}& y_{k+1}=h_{\mathit{proj}}\circ F(x_0,{\theta }_{k+1},p,t)+{ϵ}_{k+1},\hfill \end{array}$$where the process noise $w_k$ and the observation noise ${ϵ}_{k+1}$ are assumed to be statistically uncorrelated Gaussian noise with zero-mean. Note that Eq. 3 is valid for the estimation of time-invariant model parameters, like those considered in this simulation study, but can be generalised to time-varying parameters.

Unscented Kalman filter

In order to estimate the correct values of the chosen simulation parameters, we use a data assimilation method: The Kalman filter [12] enables estimation of non-measurable system quantities of dynamic systems. It is constructed as mean squared error minimiser taking into account both measurement error and accuracy of the model describing the dynamic system. Normal distribution of noise in the state dynamics and measurements as well as linearity of state dynamics and the observation process are assumptions underlying the Kalman filter. In this case, the observation operator (Eq. 4) is a combination of running a (non-linear) FE simulation and computing a virtual X-ray projection from a focal spot to the detector screen, which mathematically corresponds to a central projection. Due to non-linearity of the central projection and the FE simulation we make use of the Additive Unscented Kalman filter (UKF). This variant of the Kalman filter proposed by [11] enables the application of non-linear observation processes and state dynamics by parametrising state mean and covariance with a set of discretely sampled points. This method cannot only be used for state estimation, but also in the context of parameter estimation (e.g., [2]). The mathematical formulation of the UKF from [11] is put in context of the state-space model explained above as follows: Given an initial set of parameters ${\theta }_0\in {\mathbb{R}}^L$ and initial state covariance ${\widehat{P}}_0^{\theta \theta }\in {\mathbb{R}}^{L\times L}$ (where L is the number of parameters to be estimated), a set of sample points, the so-called sigma points ${\chi }_{k-1}^{(i)}$ with corresponding weights $W_i$ is computed according to:

5

$$\begin{array}{cc}\hfill {\chi }_{k-1}^{(0)}=& {\theta }_{k-1},\hfill \end{array}$$

6

$$\begin{array}{cc}\hfill {\chi }_{k-1}^{(i)}=& {\theta }_{k-1}+{\sqrt{(L+\lambda ){\widehat{P}}_{k-1}^{\theta \theta }}}_i,i=1,...,L,\hfill \end{array}$$

7

$$\begin{array}{cc}\hfill {\chi }_{k-1}^{(i)}=& {\theta }_{k-1}-{\sqrt{(L+\lambda ){\widehat{P}}_{k-1}^{\theta \theta }}}_i,i=L+1,...,2L,\hfill \end{array}$$

8

$$\begin{array}{cc}\hfill W_m^{(0)}=& \lambda /(L+\lambda ),\hfill \end{array}$$

9

$$\begin{array}{cc}\hfill W_c^{(0)}=& \lambda /(L+\lambda )+(1-{\alpha }^2+\beta ),\hfill \end{array}$$

10

$$\begin{array}{cc}\hfill W_m^{(i)}=& W_c^{(i)}=1/(2(L+\lambda )),i=1,...,2L,\hfill \end{array}$$where k is the assimilation time step (starting at $k=1$), $\lambda ={\alpha }^2(L+\kappa )-L$, $\alpha$, $\beta$ are parameters preset in the used implementation (pykalman library, https://pykalman.github.io/; accessed May 11, 2020) and the term ${\sqrt{(L+\lambda ){\widehat{P}}_{k-1}^{\theta \theta }}}_i$ denotes the ith row of the matrix square root. Upon computation of the sigma points, follows the “Time Update” or “Prediction Step”: At first the sigma points are propagated to the next time step using the dynamic state model or transition operator f. This corresponds to updating the time-invariant parameter values according to Eq. 3, therefore:

11

$$\begin{array}{c}\hfill {\chi }_k^{(i)}=f({\chi }_{k-1}^{(i)})=\mathbb{I}{\chi }_{k-1}^{(i)}.\end{array}$$Then, the updated sigma points are used to determine the predicted parameter set ${\theta }_{k|k-1}$ and covariance ${\widehat{P}}_{k|k-1}^{\theta \theta }$, where the subscript $k|k-1$ denotes the prediction at step k using the data from step $k-1$:

12

$$\begin{array}{cc}& {\theta }_{k|k-1}=\sum _i{W}_m^{(i)}{\chi }_k^{(i)},\hfill \end{array}$$

13

$$\begin{array}{cc}& {\widehat{P}}_{k|k-1}^{\theta \theta }=\sum _i{W}_c^{(i)}[{\chi }_k^{(i)}-{\theta }_{k|k-1}]{[{\chi }_k^{(i)}-{\theta }_{k|k-1}]}^T+Q_k,\hfill \end{array}$$where $Q_k$ is the process covariance associated with the errors $w_k$ of the transition function (refer to Eq. 3). Then, the updated sigma points are projected to observation space using the observation operator $h_k({\chi }_k^{(i)})=h_{\mathit{proj}}\circ F(x_0,{\chi }_k^{(i)},p,t)$ (Eq. 4), and subsequently used to compute the observation prediction ${\widehat{y}}_{k|k-1}$ and the corresponding predicted observation covariance ${\widehat{P}}_{k|k-1}^{\mathit{yy}}$:

14

$$\begin{array}{cc}& {\mathcal{Y}}_k^{(i)}=h_k({\chi }_k^{(i)}),\hfill \end{array}$$

15

$$\begin{array}{cc}& {\widehat{y}}_{k|k-1}=\sum _i{W}_m^{(i)}{\mathcal{Y}}_k^{(i)},\hfill \end{array}$$

16

$$\begin{array}{cc}& {\widehat{P}}_{k|k-1}^{\mathit{yy}}=\sum _i{W}_c^{(i)}[{\mathcal{Y}}_k^{(i)}-{\widehat{y}}_{k|k-1}]{[{\mathcal{Y}}_k^{(i)}-{\widehat{y}}_{k|k-1}]}^T+R_k,\hfill \end{array}$$where $R_k$ is the measurement covariance associated with the measurement errors ${ϵ}_k$ (Eq. 4). After this follows the ”Measurement Update” or ”Correction Step”: First, the cross-covariance matrix ${\widehat{P}}_{k|k-1}^{\theta y}$ is computed as follows:

17

$$\begin{array}{c}\hfill {\widehat{P}}_{k|k-1}^{\theta y}=\sum _i{W}_c^{(i)}[{\chi }_k^{(i)}-{\theta }_{k|k-1}]{[{\mathcal{Y}}_k^{(i)}-{\widehat{y}}_{k|k-1}]}^T.\end{array}$$From this, the Kalman gain is determined:

18

$$\begin{array}{c}\hfill K_k={\widehat{P}}_{k|k-1}^{\theta y}{({\widehat{P}}_{k|k-1}^{\mathit{yy}})}^{-1}.\end{array}$$Subsequently, the measurement data $y_k$ (dimension d, not necessarily equal to L) is incorporated to find the updated parameter set ${\theta }_{k|k}$ as well as the updated parameter covariance ${\widehat{P}}_{k|k}^{\theta \theta }$:

19

$$\begin{array}{c}\hfill {\theta }_{k|k}={\theta }_{k|k-1}+K_k(y_k-{\widehat{y}}_{k|k-1}),\end{array}$$

20

$$\begin{array}{c}\hfill {\widehat{P}}_{k|k}^{\theta \theta }={\widehat{P}}_{k|k-1}^{\theta \theta }-K_k{\widehat{P}}_{k|k-1}^{\mathit{yy}}{K}_k^T.\end{array}$$To proceed to the next time step, a new set of sigma points is computed based on the updated parameters and their covariance, followed by applying the prediction and correction steps. In this way, discrete measurements over time are used to iteratively correct a non-linear dynamic model. The Kalman filtering process described above and the role of FE simulation, X-ray projection as well as observation and transition operator within this process are visualized in Fig. 4.

Fig. 4 [Images not available. See PDF.]

Process diagram visualizing the Kalman filtering approach and the role of FE simulation, X-ray projection as well as observation and transition operator within this process

Implementation of observing system simulation experiment

As mentioned above, for the parameter estimation in this simulation study, we chose two material parameters, namely density $\rho$ and Young’s modulus E of the square rod, as well as one initial condition, namely the angular velocity $\omega$ (therefore $L=3$, ${\theta }_k=[{\rho }_k,E_k,{\omega }_k]$). The parameters of the simulation were chosen to be $\rho =2.7\text{g}·{\text{cm}}^{-3},E=70\text{GPa},\omega =0.100{\text{ms}}^{-1}$, whereas the parameters of the fictional experiment or ”truth run” were set to $\rho =2.2\text{g}·{\text{cm}}^{-3},E=60\text{GPa},\omega =0.055{\text{ms}}^{-1}$. Note that these parameters describe generic material properties used as an exemplary test case. The first part of the analysis was performed using node 12221 situated on a corner of the rod (see Fig. 1). The assimilation was performed for 40 assimilation steps with a step size of $0.025\text{ms}$ (corresponding to $s^{-1}$). The interaction with the FE simulation was realised by writing a wrapper around the LS-Dyna numerical simulation framework, modifying input files, running the simulation and analysing output files for every time step. The virtual X-ray geometry was set up with a source detector distance of 600 cm, with the rod positioned at a distance from the source of 300 cm and a virtual detector size of 250 x 250 cm (see Fig. 5). These values are not meant to represent a realistic X-ray setup and are chosen for demonstration purposes only.

Fig. 5 [Images not available. See PDF.]

Virtual X-ray setup used to compute projected coordinates of the considered nodes

To generate synthetic measurement data from the truth run, the spatial coordinates of the considered node were projected onto the detector plane for each time step. Then, measurement noise was introduced to the data by adding normal distributed numbers with a standard deviation of $0.03\text{cm}$ in both dimensions to the projected coordinates (dimension of the observation $d=2$). Due to the generic nature of the test case, this value does not represent a physical measurement error of a standard X-ray detector, but it was chosen in relation to the spatial extent of the node’s trajectory on the virtual X-ray detector. The following covariances were set to initialize the Kalman filter:

21

$$\begin{array}{cc}\hfill R_k=& \left(\begin{array}{ccc}0.{03}^2& 0& \\ 0& 0.{03}^2\end{array}\right)\text{cm},\hfill \end{array}$$

22

$$\begin{array}{cc}\hfill Q_k=& \left(\begin{array}{ccc}0.01& 0& 0\\ 0& 0.01& 0\\ 0& 0& 0.0005\end{array}\right)\left(\begin{array}{c}\text{g}·{\text{cm}}^{-3}\\ \text{GPa}\\ {\text{ms}}^{-1}\end{array}\right),\hfill \end{array}$$

23

$$\begin{array}{cc}\hfill {\widehat{P}}_0^{\theta \theta }=& \left(\begin{array}{ccc}0.5^2& 0& 0\\ 0& {10}^2& 0\\ 0& 0& 0.{045}^2\end{array}\right)\left(\begin{array}{c}\text{g}·{\text{cm}}^{-3}\\ \text{GPa}\\ {\text{ms}}^{-1}\end{array}\right).\hfill \end{array}$$As the parameters remain unchanged when propagated over time, the process covariance $Q_k$ was set to small numbers with respect to the absolute parameter values. The initial state covariance ${\widehat{P}}_0^{\theta \theta }$ was set to a diagonal matrix, neglecting possible correlations between the parameters as a first-order approximation. The variances of the respective parameters (diagonal entries of the covariance matrix) were chosen, such that the deviating parameters underlying the truth run lie just within the interval spanned by the corresponding standard deviation. To evaluate the influence of the number of observations on the accuracy of the adaptation of the FE simulation performed by the Kalman filter, the analysis was repeated with ten nodes measured on the projections (nodes marked in Fig. 1). Here, all input parameters and covariance matrices were kept the same for comparability of the results. Note, that using ten two-dimensional observations increases the dimension of the observation vector $y_k$ to $d=20$, such that the size of $R_k$ increases to 20x20, where the observation error of each node was kept the same:

24

$$\begin{array}{c}\hfill R_k=\left(\begin{array}{cccc}0.{03}^2& 0& \cdots & 0\\ 0& \ddots & \ddots & ⋮\\ ⋮& \ddots & \ddots & 0\\ 0& \cdots & 0& 0.{03}^2\end{array}\right)\text{cm}.\end{array}$$

Results and discussion

In the first part of the analysis, data assimilation was performed on the aforementioned simulation parameters using the projected position of one node on the square rod as observation. For an adapted set of parameters found in assimilation step k, the FE simulation was started again from time $t=0$ up to time $t=k\ast s^{-1}$ using those parameters. The resulting 3D trajectory of the considered node is saved for each assimilation step. These trajectories can then be compared to the trajectory yielded by the initial simulation as well as the respective trajectory from the truth run. Figure 6 shows the trajectories in three dimensions (top) as well as their projections on the virtual X-ray detector (bottom).

Fig. 6 [Images not available. See PDF.]

Three-dimensional visualization of the results of assimilating the trajectory of node 12221 (top) for some of the assimilation steps. Projection on the virtual X-ray detector (bottom). Blue corresponds to the trajectory of this node from the initial simulation, red corresponds to the trajectory from the truth run and green colours correspond to the assimilated trajectories at different assimilation steps. On the detector, the red trajectory corresponds to the noisy measurements, whereas the projected data from the truth run is shown in black

It can be observed that the final assimilated trajectory (bright green in Fig. 6) coincides well with the trajectory of the experiment on the detector. Note that fluctuations due to noise in the experimental trajectory are not expected to be tracked by the assimilated trajectory. Instead, the assimilated trajectory is expected to track the true, noise-free trajectory underlying the measurements. Looking at the projection of the noise-free truth run trajectory, the Kalman filter accounts for the noise correctly and closely traces the noise-free truth run (unknown to the analysis process). This is also validated in the 3D visualization: Again, the final assimilated trajectory is in close accordance to the truth run. This indicates that by using the information inherent to the physical laws underlying the FE simulation, it is possible to correctly augment the two-dimensional projections to three dimensions. Comparing the trajectories yielded by the different assimilation steps, the process of approaching the true parameters of the simulation becomes visible: Whilst the trajectories from the early assimilation steps do not match the ground truth yet, the trajectories from the later assimilation steps approach the ground truth more closely, such that not only the current position is estimated more accurately, but also the positions previous to the current assimilation step. Note that the spatial extent of the considered trajectory in this example only amounts to approximately 0.05 cm x 1.26 cm x 5.00 cm, which represents a rather small deformation process considering the size of the rod. However, by ensuring that the measurement uncertainty on the detector scales reasonably with the spatial extent of the trajectory, the applicability of the proposed method is independent of spatial dimensions. As a second step, the same assimilation method was applied measuring the projected trajectories of ten nodes on the rod simultaneously. For comparison, the resulting trajectory for the same node (12221) is depicted in three dimensions and on the virtual X-ray detector in Fig. 7. Again, the resulting assimilated trajectory closely follows the ground truth. By visual comparison to Fig. 6 it appears that the assimilation using ten nodes results in more accurate results compared to the analysis using just one node.

Fig. 7 [Images not available. See PDF.]

Three-dimensional visualization of the assimilated trajectories of the same node 12221 (top) for some of the assimilation steps. Projection on the virtual X-ray detector (bottom). For this analysis, the trajectories of ten nodes were assimilated simultaneously, e.g. the positions of ten nodes on the detector were measured for each time step. Blue corresponds to the trajectory of this node from the initial simulation, red corresponds to the truth run and green colours correspond to the assimilated trajectories for at different assimilation steps. On the detector, the red trajectory corresponds to the noisy measurements, whereas the projected data from the truth run is shown in black

Figure 8 shows the spatial deviation between initial simulation and truth run compared to the spatial deviation between the final assimilated trajectory and truth run over time and in all three spatial directions for both analyses. Note that for the analysis using ten nodes, mean deviations are computed over all ten nodes for each time step, where the error bars depict the corresponding standard deviation. Obviously, the deviation in x- and y-direction (which correspond to the detector plane) between assimilation and truth run is considerably smaller as compared to the initial simulation. This could also be achieved by simply shifting or deforming the object unphysically to make the projections match. The benefit of the proposed method becomes manifest when considering the third dimension: the deviation between assimilation and truth run is also reduced in z-direction (viewing direction), where no information is available directly from the measurement. This information gap is filled by including the underlying physics of the FE simulation via adaptation of chosen physical parameters.

Fig. 8 [Images not available. See PDF.]

Deviation between assimilated trajectory and truth run (green) as compared to deviation between initial simulation and truth run (blue) in x- (top), y-(middle) and z-direction (bottom) for both analyses (dark blue/dark green: one node, light blue/bright green: ten nodes). For the analysis using ten nodes the deviations were averaged over all ten nodes for each time step, where the error bars depict the corresponding standard deviation

Comparing between both of the analyses (using one vs. ten nodes) confirms the obvious expectation that the mean deviation between assimilated trajectories and ground truth trajectories can be reduced further when using more nodes for observations. Getting a number of the order of magnitude of ten measurement points on the X-ray image for each investigated component seems feasible: In case the component does not offer enough remarkable features, it is common in X-ray testing to attach a number of X-ray markers, e.g. spherical lead markers at predefined positions which can be matched to a node in the FE simulation. Note that including more observations in the assimilation step does not increase the computation time substantially. The parameter values found in the final assimilation step of the analysis using one node amount to $\rho =2.00\pm 0.17\text{g}·{\text{cm}}^{-3}$, $E=65.9820\pm 5.3224\text{GPa}$, $\omega =0.053\pm 0.003{\text{ms}}^{-1}$. For the analysis using ten nodes, the final parameters are found to be $\rho =2.18\pm 0.03\text{g}·{\text{cm}}^{-3}$, $E=59.8447\pm 0.7675\text{GPa}$, $\omega =0.054\pm 0.0004{\text{ms}}^{-1}$. The corresponding uncertainties were determined using the final parameter covariance of the assimilation process: From the diagonal entries corresponding to the variance of the respective parameters, the standard deviation can be computed. The parameters density and angular velocity are in good agreement with those underlying the truth run, whereas the Young’s module could only be approached for the assimilation using ten nodes. Note that by testing the method with just a single measured node trajectory the merit of incorporating the physics intrinsic to the FE simulation becomes obvious: although measuring just a single displacement curve, we can correctly infer density and angular velocity of the rod by making use of the known physical correlation between the simulation parameters and the measured node position on the X-ray detector.

The process of adaptation of the chosen parameters density $\rho$, Young’s modulus E and angular velocity $\omega$ over assimilation time is visualized in Fig. 9. The corresponding uncertainty is determined using the current parameter covariance for each time step of the assimilation process.

Fig. 9 [Images not available. See PDF.]

Visualization of parameter adaptation (top: density $\rho$, middle: Young’s module E, bottom: angular velocity $\omega$) over the assimilation time steps for both analyses (dark green squares: one node, bright green diamonds: ten nodes). Error bars depict the computed uncertainty for each parameter yielded by the current state covariance for each time step in the Kalman filter. Parameter values of the initial simulation are depicted in blue, red corresponds to parameter values of the truth run, green corresponds to the assimilated parameter values at the respective time step

For the assimilation using ten nodes, it can be observed that the parameter value for the angular velocity approaches the value underlying the truth run within few assimilation time steps. The associated uncertainty is reduced considerably during the assimilation, confirming that the parameter could be matched accurately. The parameter value for the density could also be matched accurately to the density value underlying the ground truth. For this parameter the true value was approached slightly more slowly and less accurately than for the angular velocity: Whilst for the angular velocity the true value was found within approximately five assimilation steps and the uncertainty was reduced to approximately 1.0% of the initial uncertainty, approximately seven time steps were needed for the density and the uncertainty was reduced to 5.2% of the initial uncertainty. The parameter value for the Young’s modulus only approaches the true value after approximately 20 assimilation steps, where the uncertainty still amounts to 7.6% of the initial uncertainty. By comparing the uncertainty of the assimilated parameter to its initial uncertainty the quality of the assimilation results can be assessed: If the uncertainty associated to the assimilated parameter still amounts to a high percentage of the initial uncertainty or did not decrease at all, this indicates that the assimilation results for this parameter are not reliable. However, if the uncertainty is reduced to a small fraction of the initial uncertainty, it can be concluded, that the assimilated parameter likely matches the underlying truth. It is worth noting that although the parameters do not match the ground truth yet in the early assimilation steps, the method still yields improvement for these time steps as the reconstructed trajectories are retrieved from a simulation run with the final estimated parameter set. Comparing the assimilation using one and ten nodes again demonstrates the supremacy of using a higher number of observations: The true values of the angular velocity and the density were not only approached more quickly, they could also be matched more accurately as the remaining parameter uncertainty is considerably smaller. Using only one node as observation yielded remaining uncertainties in the assimilated parameter values of 34.2% (density), 53.2% (Young’s modulus) and 7.5% (angular velocity) of the initial uncertainty. For the Young’s modulus, using only one node did not suffice to recover the true parameter value. Note that the results of the Kalman filtering process are strongly dependent on the covariance matrices chosen to initialise the assimilation process. Here, a diagonal matrix was chosen, neglecting possible correlations between the parameters. Whilst this is a valid assumption for the parameter of the angular velocity, the correlation between density and Young’s modulus could be strong enough to have an impact on the results. Approximating the covariance between chosen parameters for assimilation more accurately and assessing the impact on assimilation results could be a topic for further research. Considering the varying performance of the assimilation process regarding the respective parameters, it is also presumed that changes in both the density as well the Young’s modulus have a smaller impact on the trajectory of the observed nodes than the angular velocity: As the parameters can only be adapted based on their correlation with the measurement, parameters with smaller influence on the measurement cannot be recovered as easily. The sensitivity of the measured trajectory to the respective simulation parameters is assessed in Fig. 10. Here, the trajectories of the considered node for the original simulation as well as for a set of adapted simulations are depicted, where only one of the three parameter was changed, whilst the other two parameters were kept equal to the original simulation. The following parameter changes were considered: ${\rho }_{\pm }={\rho }_0\pm \Delta \rho =2.7\pm 0.5\text{g}·{\text{cm}}^{-3}$, $E_{\pm }=E_0\pm \Delta E=70\pm 10\text{GPa}$, ${\omega }_{\pm }={\omega }_0\pm \Delta \omega =0.1\pm 0.045{\text{ms}}^{-1}$, where ${\rho }_0$, $E_0$ and ${\omega }_0$ are the parameter values of the original simulation. In order to quantify the impact of each parameter on the difference between simulated and measured trajectories, the respective parameters were varied by an amount matching the difference between parameter deviation between original simulation and ground truth, such that $\Delta \rho =|{\rho }_0-{\rho }_{groundtruth}|$, $\mathrm{\Delta }E=|E_0-E_{groundtruth}|$ and $\Delta \omega =|{\omega }_0-{\omega }_{groundtruth}|$. Compared to the original trajectory (orange), changing the angular velocity has the strongest impact (dots), while changes in the density (cross marks) and the Young’s modulus (diamonds) result in trajectories which are considerably closer to the original trajectory. This confirms the aforementioned presumption: The assimilation performs best on parameters which correlate strongly with the measurement on the detector, whereas parameters which are less correlated are approached less accurately.

This prompts the question which parameters to choose for assimilation. From the analysis above, it becomes clear that in order to be able to accurately infer the true parameters, these have to be sufficiently correlated to the measurement on the X-ray detector. Note that the degree of correlation between simulation parameters and observation can be influenced to a certain extent by optimizing the experimental X-ray setup.

Fig. 10 [Images not available. See PDF.]

Trajectory of node 12221 on the detector for different simulations: the trajectory from the original simulation is depicted in orange. To analyse the sensitivity of the measurable trajectory to parameter changes, the following parameter adaptations to the original simulation have been performed and plotted: $\rho =2.7\pm 0.5\text{g}·{\text{cm}}^{-3}$ (diamonds), $E=70\pm 10\text{GPa}$ (cross marks), $\omega =0.1\pm 0.045{\text{ms}}^{-1}$ (dots). Negative adaptations are depicted in purple, positive adaptations are depicted in green

At the same time, a potential difficulty when applying the method to real experimental data, where no ground truth is available, is the choice of parameters: in this case, it is not obvious which parameters in the highly complex simulation model are causing differences between observations and numerical computations. To date simulations have to be adapted entirely based on experience of car crash simulation engineers: By comparing experimental data and FE simulation, the expert adjusts the simulation based on a sophisticated guess. Therefore, given the expert advice on the choice of parameters for adaptation, the proposed method constitutes a significant benefit as the chosen parameters are inferred by quantitatively incorporating the available measurement data. It is worth noting that the choice of parameters to be estimated is not restricted to material properties or initial/boundary conditions. In fact, given continuity of the quantity to be inferred, any variable or parameter of the FE simulation can be chosen, e.g. geometric variables like the shell thickness of a certain part. The problem of automatically choosing the right parameters for adaptation could be addressed by future research. It is conceivable that approaches could be transferred from the field of FEM parameter sensitivity analysis. Another challenge comes with the computational cost of the suggested approach as it involves a significant number of repeated simulation runs. Compared to other non-linear parameter estimation methods (e.g. Monte Carlo sampling), the UKF is fairly computationally efficient. However, with complex tasks, this approach is still very time-consuming. Possible solutions to be investigated in future research could be applying the method to sub-domains of the simulation (e.g. when examining a part in the front of the car, the compuational effort could be reduced to a domain surrounding this part).

Conclusions

X-ray technology is a promising tool to make accessible information concerning dynamics and hidden deformation of internal parts of a car during crash. However, from X-ray images only 2D information in the detector plane is available. In order to be able to draw conclusions on dynamics and deformations in three dimensions we augment the two-dimensional observations from X-ray measurements using the information inherent to a numerical FE simulation of the deformation process. To do so, we adapt chosen parameters of the FE simulation transferring methods from data assimilation, yielding an optimized simulation matching the experimental measurements. We perform two simulations of the same simplified test case (a rotating and bending square rod) differing in three chosen parameters to investigate the proposed approach. Artificial measurements are generated by computing two-dimensional projections of certain nodes of the second simulation and adding normal distributed noise. This truth-run simulation can then be used to validate and assess the accuracy of the computed results. The proposed method proves suitable for adapting a three-dimensional numerical simulation of a deformation process based on two-dimensional X-ray images. The 3D trajectories of the considered nodes from the adapted simulation are in closer agreement to the experiment (truth-run) than the original simulation: The mean deviation between assimilated trajectory and experiment is significantly smaller than the deviation between the original simulation and the experiment in all three spatial directions. Apart from yielding an adapted simulation matching the acquired measurements, uncertainties associated to the adapted simulation parameters result directly from the assimilation process and can be used to assess the quality of the computed results. It has to be noted, that only those parameters of the FE simulation that correlate with the observations can be adapted with reasonable accuracy; from the measurements, no information on uncorrelated quantities can be retrieved. The method presented in this paper does not only prove feasible to augment 2D X-ray measurements with a-priori information from FE simulations. It could also provide a pathway towards optimization of FE simulations based on measurements in order to improve their predictive power, which could constitute an important step towards digitalisation of automotive crash testing. If viable, the method can be used to enhance the investigation of dynamic processes in a wide range of applications beyond the scope of car safety development, assuming that a simulation or analytical model of the dynamic process is available.

Acknowledgements

We would like to offer our special thanks to Dr. Benjamin Lang (Fraunhofer EMI) for designing and preparing the FE simulation of the chosen test case. The authors also wish to thank Prof. Leonhard Scheck (Deutscher Wetterdienst), MSc. Thomas Soot (Fraunhofer EMI) and MSc. Michael Schmitt (Fraunhofer EZRT) for fruitful scientific discussions and valuable input.

Open Access

This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Author’s contributions

IB originated the presented method, derived and implemented the algorithms, performed data analysis and interpreted the results and drafted the manuscript. SH, SM and SN conceived of the scientific concept and substantively revised the manuscript. SM and SN supervised the work and were in charge of project administration.

Funding

Open Access funding enabled and organized by Projekt DEAL.

Declarations

Funding

This work was financed in the scope of the internal Fraunhofer research program MAVO fastXcrash and the Fraunhofer EMI financed research project XCC.

Conflict of interest

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

Availability of data and materials

The datasets used and analysed within the presented study are available from the corresponding author on reasonable request.

Code availability

The code is available from the corresponding author on reasonable request.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

References

1. Abt F, Boley M, Weber R, Graf T, Popko G, Nau S (2011) Novel x-ray system for in-situ diagnostics of laser based processes—first experimental results. Physics procedia, lasers in Manufacturing 2011—proceedings of the sixth international WLT conference on lasers in manufacturing, 12:761–770

2. Astroza R, Barrientos N (2020) On the performance of unscented Kalman filters in parameter estimation of nonlinear finite element models. Exp Vib Anal Civil Struct :399–407, Doi: https://doi.org/10.1201/9781003090564-39

3. Bacchus, A; Tounzi, A; Argaud, JP; Bouriquet, B; Biet, M; MacAire, L; Le Menach, Y. Estimation of FEM model parameters using data assimilation and its application to an electrical machine. IEEE Trans Mag; 2016; [DOI: https://dx.doi.org/10.1109/TMAG.2015.2495265]

4. Benson, DJ. Computational methods in Lagrangian and Eulerian hydrocodes. Comput Methods Appl Mech Eng; 1992;MathSciNet ID: 1186728[DOI: https://dx.doi.org/10.1016/0045-7825(92)90042-I]zbMath ID: 0763.73052

5. Brainerd EL, Baier DB, Gatesy SM, Hedrick TL, Metzger KA, Gilbert SL, Crisco JJ (2010) X-ray reconstruction of moving morphology (XROMM): precision, accuracy and applications in comparative biomechanics research. J Exp Zool Part A Ecol Genet Physiol 9999A(May 2009):n/a–n/a. Doi: https://doi.org/10.1002/jez.589

6. Chen, B; Wang, Y; Wang, R; Zhu, Z; Ma, L; Qiu, X; Dai, W. The gray-box based modeling approach integrating both mechanism-model and data-model: The case of atmospheric contaminant dispersion. Symmetry; 2020; [DOI: https://dx.doi.org/10.3390/sym12020254]

7. Choi, JH; Lee, S. Real-time tumor motion tracking in 3D using planning 4D CT images during image-guided radiation therapy. Algorithms; 2018; 11, 10 155. [DOI: https://dx.doi.org/10.3390/a11100155]

8. Ebrahimian, H; Astroza, R; Conte, JP. Extended Kalman filter for material parameter estimation in nonlinear structural finite element models using direct differentiation method. Earthquake Eng Struct Dyn; 2015; 44, pp. 1495-1522. [DOI: https://dx.doi.org/10.1002/eqe.2532]

9. Hiermaier, S. Pred Model Dyn Process Trib Prof Klaus Thoma; 2009; [DOI: https://dx.doi.org/10.1007/978-1-4419-0727-1]

10. Hopp, T; Dietzel, M; Baltzer, PA; Kreisel, P; Kaiser, WA; Gemmeke, H; Ruiter, NV. Automatic multimodal 2D/3D breast image registration using biomechanical FEM models and intensity-based optimization. Med Image Anal; 2013; [DOI: https://dx.doi.org/10.1016/j.media.2012.10.003]

11. Julier, SJ; Uhlmann, JK. New extension of the Kalman filter to nonlinear systems. Signal Process Sensor Fusion Target Recog; 1997; 3068, pp. 182-193.

12. Kalman, RE. A new approach to linear filtering and prediction problems. J Fluids Eng Trans ASME; 1960; 10, 1115/1 3662552.MathSciNet ID: 121033

13. Kostelich, EJ; Kuang, Y; McDaniel, JM; Moore, NZ; Martirosyan, NL; Preul, MC. Accurate state estimation from uncertain data and models: an application of data assimilation to mathematical models of human brain tumors. Biol Direct; 2011; [DOI: https://dx.doi.org/10.1186/1745-6150-6-64]

14. Kurfiß M, Bösl P, Popko G, Moser S, Nau S, Fritsch J (2020) Fast high-energy radiography at 1khz. In: paper presented at the 20th world conference on non-destructive testing

15. Lee, AW; Rajagopal, V; Babarenda Gamage, TP; Doyle, AJ; Nielsen, PM; Nash, MP. Breast lesion co-localisation between X-ray and MR images using finite element modelling. Med Image Anal; 2013; [DOI: https://dx.doi.org/10.1016/j.media.2013.05.011]

16. Li, R; Jia, X; Lewis, J; Gu, X; Folkerts, M; Men, C; Jiang, S. WE-C-204B-05: Real-Time Volumetric Image Reconstruction and 3D Tumor Localization Based on a Single X-Ray Projection Image during Lung Cancer Radiotherapy. Med Phys; 2010; 10, 1118/1 3469377.

17. Mao, H; Zhang, L; Jiang, B; Genthikatti, VV; Jin, X; Zhu, F; Makwana, R; Gill, A; Jandir, G; Singh, A; Yang, KH. Development of a finite element human head model partially validated with thirty five experimental cases. J Biomech Eng; 2013; 10, 1115/1 4025101.

18. Marami, B; Sirouspour, S; Capson, DW. Non-rigid registration of medical images based on estimation of deformation states. Phys Med Biol; 2014; [DOI: https://dx.doi.org/10.1088/0031-9155/59/22/6891]

19. Marklund, PO; Nilsson, L. Optimization of a car body component subjected to side impact. Struct Multidiscip Optim; 2001; 21, 5 pp. 383-392. [DOI: https://dx.doi.org/10.1007/s001580100117]

20. Marzougui D, Kan C, Bedewi N (1996) Development and validation of an NCAP simulation using LS-DYNA3D. In: 4th International LS-DYNA3D conference

21. Marzougui D, Samaha R, Nix L, Kan C (2013) Extended validation of the finite element model for the 2010 Toyota Yaris passenger sedan (MASH 1100kg vehicle). In: Transportation research board 92nd annual meeting

22. Mohammadiha O, Beheshti H, Aboutalebi FH (2015) Multi-objective optimisation of functionally graded honeycomb filled crash boxes under oblique impact loading. Int J Crashworthiness. https://doi.org/10.1080/13588265.2014.970398

23. Moser S, Wickert M, Nau S (2019) High-speed x-ray imaging and 3D analysis of impact-formed fragments. Proceedings of SPIE 11051. In: 32nd International congres on high-speed imaging and photonics (January), 25. Doi: https://doi.org/10.1117/12.2525447

24. Navon, IM. Data assimilation for numerical weather prediction: a review. Data Assim Atmos Ocean Hydrol Appl; 2009; [DOI: https://dx.doi.org/10.1007/978-3-540-71056-1_2]

25. Obradovic, J; Boria, S; Belingardi, G. Lightweight design and crash analysis of composite frontal impact energy absorbing structures. Comp Struct; 2012; 94, 2 pp. 423-430. [DOI: https://dx.doi.org/10.1016/j.compstruct.2011.08.005]

26. Parrish, A; Rais-Rohani, M; Najafi, A. Crashworthiness optimisation of vehicle structures with magnesium alloy parts. Int J Crashworthiness; 2012; [DOI: https://dx.doi.org/10.1080/13588265.2011.648518]

27. Rentschler W, Uffenkamp V (1999) Digital photogrammetry in analysis of crash tests. In: SAE technical papers, 10.4271/1999-01-0081

28. Roberts C, Kerrigan JR (2015) Injuries and kinematics: response of the cervical spine in inverted impacts. In: 24th International technical conference on the enhanced safety of vehicles (ESV) (15-0432):1–27

29. Salamon M, Boehnel M, Reims N, Ermann G, Schmitt M, Uhlmann N, Hanke R, Iis F, Ezrt XrT (2013) Applications and methods with high energy CT systems. In: 5th International symposium on NDT in aerospace (November 2013), pp 13–15

30. Shigeta K, Kitagawa Y, Yasuki T (2009) Development of next generation human FE-model capable of organ injury prediction. In: Proceedings of the 21st enhanced safety vehicle (ESV) conference

31. Wikle, CK. Atmospheric modeling, data assimilation, and predictability. Technometrics; 2005; [DOI: https://dx.doi.org/10.1198/tech.2005.s326]

32. Yankeelov, TE; Quaranta, V; Evans, KJ; Rericha, EC. Toward a science of tumor forecasting for clinical oncology. Cancer Res; 2015; [DOI: https://dx.doi.org/10.1158/0008-5472.CAN-14-2233]

33. Yasuhiro D, Shigeo T, Kazuki O (2017) An evaluation method of cervical spinal injury for car passenger in dynamic rollover using a human FE model. In: 25th International technical conference on the enhanced safety of vehicles (ESV)