# Deep Self-supervised Learning with Visualisation for Automatic Gesture
Recognition

Fabien ALLEMAND
Department ARTEMIS
Télécom SudParis
<EMAIL_ADDRESS>
Alessio MAZZELLA
Department ARTEMIS
Télécom SudParis
<EMAIL_ADDRESS>
Jun VILLETTE
Department ARTEMIS
Télécom SudParis
<EMAIL_ADDRESS>
Decky ASPANDI LATIF
Department ARTEMIS
Télécom SudParis
<EMAIL_ADDRESS>
Titus ZAHARIA
Department ARTEMIS
Télécom SudParis
<EMAIL_ADDRESS>

###### Abstract

Gesture is an important mean of non-verbal communication, with visual modality
allows human to convey information during interaction, facilitating peoples
and human-machine interactions. However, it is considered difficult to
automatically recognise gestures. In this work, we explore three different
means to recognise hand signs using deep learning: supervised learning based
methods, self-supervised methods and visualisation based techniques applied to
3D moving skeleton data. Self-supervised learning used to train fully
connected, CNN and LSTM method. Then, reconstruction method is applied to
unlabelled data in simulated settings using CNN as a backbone where we use the
learnt features to perform the prediction in the remaining labelled data.
Lastly, Grad-CAM is applied to discover the focus of the models. Our
experiments results show that supervised learning method is capable to
recognise gesture accurately, with self-supervised learning increasing the
accuracy in simulated settings. Finally, Grad-CAM visualisation shows that
indeed the models focus on relevant skeleton joints on the associated gesture.

_K_ eywords Deep learning $\cdot$ Automatic gesture recognition $\cdot$ Self-
supervised learning $\cdot$ Visualisation

##  1 Introduction

Gesture is an important means of non-verbal communication facilitating human
interactions. Certain conditions of disability or other constraints force some
people to communicate using sign languages. However, sign language can also be
used for human to machine interactions. In this context, most smart assistants
use voice as input, but sign language recognition, as shown in Figure 1b,
should be used to enhance accessibility.

With the recent development of wearable, virtual reality and spatial
computing, hand tracking and gesture recognition have become key method in
human-machine communications. Meta offers virtual reality headsets capable of
real-time hand tracking (Figure 1a), and the latest generations of Apple
products (Apple Watch and Apple Vision Pro) are attempting to revolutionise
the way we interact with our electronic devices using hand gestures.

Other applications are to be expected, particularly in the automotive
industry, where vehicles are tending to become increasingly intelligent and
autonomous thanks to powerful on-board computers connected to numerous
sensors. It’s not impossible to imagine a vehicle sound or climate control
system using gesture recognition as shown in Figure 1c. Perhaps one day, there
will be possibility of the utilisation of such system similar to the computers
in futuristic films, controlled solely by voice and gestures.

(a) Real time hand tracking for Virtual Reality (VR)

(b) Sign language recognition for human interactions

(c) Hand gesture recognition for human to machine communication

Figure 1: Several applications of gesture recognition.

In machine learning context, there are several types of automatic learning:
Supervised Learning (SL), Unsupervised Learning (UL), reinforcement learning
and Self-Supervised Learning (SSL). Although the latter method has proved its
worth in other computer vision applications, the work is based on the
observation that the use of self-supervised learning methods in the field of
automatic gesture recognition is still limited. The aim of this work is
therefore to apply a self-supervised learning method using deep neural
networks to gesture recognition on 3D skeleton data and to compare the results
obtained with a conventional supervised learning method.

##  2 Methodology

###  2.1 Deep Learning for Gesture Recognition

The goal of our work is to compare the performances of several neural networks
architectures for sign language recognition in order to find the type of model
that best fit the data containing hands gestures. For this task, we choose
three well known neural network architectures and trained then using a
supervised approach. Here is a brief description of the three types of neural
networks we used:

  * •

Fully Connected (FC) model: Figure 2a shows the architecture of the FC neural
network in blue. There are three layers composed of artificial neurons. Each
neuron from a layer is connected to every neuron in the following layer. This
model is the most simple one.

  * •

Convolutional Neural Network (CNN) model: CNN architecture [3] is known to
perform well in processing two dimensional data such as images as it does not
process each value of the input independently. Instead it looks at patches of
values allowing it to learn information based on position. The convolution
operation used for that learns how to reduce the size of the input data while
preserving the most relevant information, in a way, it encodes the input data.
The light blue layers in Figure 2b depict the convolution layers and the head
of the neural network (dark blue layers) are fully connected layers that
perform the prediction on the encoded data.

  * •

Long Short-Term Memory (LSTM) model: LSTM [4] are part of a larger family of
models called recurrent neural networks. These models are known to work best
for input data that corresponds to sequences. The architecture of an LSTM
neural network is visualised in Figure 2c.

(a) Fully Connected (FC) model.

(b) Convolutional Neural Network (CNN) model.

(c) Long Short-Term Memory (LSTM) model.

Figure 2: Deep learning models used for sign language recognition.

Supervised learning is the most widespread machine learning technique. It
consists of the utilisation of the data with labels in a feedback loop that
allow the neural network to learn how to best predict the correct label given
an input sampled from the data (see Figure 3a).

The learning process may vary a little between two neural network
architectures but it always consists in adjusting the weights of the neural
network in order to minimize the loss: a measure telling how far the
prediction is from being correct. After the training, the model can be used to
make accurate predictions as shown in Figure 3b.

(a) Training step.

(b) Prediction step.

Figure 3: Supervised learning process for class prediction.

###  2.2 Self-supervised Learning

An other aspect of this work is to compare the performances of the
aforementioned neural networks when trained with regular supervised techniques
and with self-supervised techniques, a learning approach that is relevant when
dealing with a lot of unlabelled data (as it could be the case with hand
gesture) but has never been used for gesture recognition.

Self-supervised learning consists in pre-training the supervised model (that
requires expensive labelled data) with easily accessible unlabelled data[5].
This can be done by training a model following an unsupervised method (such as
input reconstruction as explained in Figure 4) requiring unlabelled data then
taking a relevant part of this model and pursue its training with labelled
data. The overall process is summarized in Figure 5. This process consists of
training the model to produce latent features [1, 2] in unsupervised way (i.e.
in this case is reconstruction as pretext task). Then, the trained models used
further on downstream task - supervised gesture recognition.

By doing so, we expect the unsupervised learning to help the model understand
the data it is working with and the supervised learning to be more efficient.
That is to say, with a fixed amount of labelled data, we expect the self-
supervised model to be more accurate than the standard supervised learning
model thanks to the unsupervised pre-training. Expected results are presented
in Figure 6.

(a) Training step.

(b) Reconstruction step.

Figure 4: Unsupervised learning process for input reconstruction. Figure 5:
Overview of self-supervised learning for class prediction. Figure 6:
Supervised versus self-supervised learning, expected results (the values are
for illustration purpose).

Self-supervised learning can also benefit from contrastive learning as
demonstrated in [6]. This technique allows the model to learn a representation
of an input sample in a latent space based on its similarity and dissimilarity
with other samples as presented in Figure 7a. Using contrastive learning
should help the model to learn a representation of the data in the latent
space where similar data are closer together as depicted in Figure 7b. We
expect that this coherence in the model will help it perform better during the
downstream task.

(a) Contrastive learning.

(b) Contrastive learning expected results (visualisation in the latent space).

Figure 7: Contrastive learning overview.

###  2.3 Grad-CAM

visualisation is an essential part of training deep learning models,
especially when working on new types of data. It deepens our understanding of
the way the models learn and can help explain their observed behaviors and
decisions. That is why we applied the Gradient-weighted Class Activation
Mapping (Grad-CAM) algorithm [7] to our CNN model.

Figure 8: Grad-CAM applied to CNN for visualisation.

The main idea of the Grad-CAM method is to visualise which parts of the input
data are more ’important’ for predictions. The process starts with an input
image, which in this case is the 2D array representing the values of the 3D
coordinates of all 79 joints of skeleton for 100 frames. The image is passed
through the layers of the CNN and the model makes a prediction. To understand
why the model made a particular prediction, the following steps are involved :

  * •

Gradient computation: Computing the gradients of the predicted class score (in
this case, one hand Mono or two hands Bi) with respect to the feature maps of
the last convolutional layer.

  * •

Global average pooling: Performing global average pooling on these gradients
to obtain the weights for each feature map. These weights represent the
importance of each feature map for the prediction.

  * •

Weighted sum: Computing a weighted sum of the feature maps using the weights
obtained in the previous step. This weighted sum highlights the important
regions in the feature maps that contribute to the prediction.

  * •

ReLU activation: Applying a ReLU activation function to the weighted sum to
obtain the final Grad-CAM heat-map. The ReLU function ensures that only
positive contributions are considered, highlighting the regions that
positively influence the prediction.

The resulting heat-map is then overlaid on the input image (with interpolation
to fit its size). To take the visualisation further, a max-pooling is done
over the three values of the three dimensional coordinates of each joint on
the heat-map, and from these values (one per joint), the ten most significant
joints for each frame are extracted and highlighted in green in the animation
of the skeleton. Figure 8 is an overview of the application of the Grad-CAM
algorithm to the skeleton data.

##  3 Task setting

###  3.1 Dataset

The skeleton animation data used for this work was provided by MocapLab
company 111https://www.mocaplab.com/fr, a high fidelity motion capture studio
located in Paris. Each skeleton sequence consists of a csv file containing the
three dimensional coordinates of each joint depicted in the skeleton over
time. In other words, each sample can be seen as an array of $3\times n+1$
columns and $t$ rows where $n$ is the number of joints in the skeleton and $t$
is the duration of the movement. (The extra column contains the time stamp for
each position). The first batch of files provided by MocapLab contained
skeletons composed of $n=79$ articulations. Figure 9 shows the first ten time
steps of the recording of a movement for the two first joints.

The dataset contained 111 skeleton files. An additional file containing the
labels for each skeleton file was also provided by the motion capture studio.
A skeleton file could be either labelled as Mono or Bi depending if the
movement of the skeleton was a one handed or a two handed sign.

Figure 9: Example of skeleton file (Avion.csv).

###  3.2 Data Pre-processing

Before being fed to the different machine learning models, the skeleton files
require some pre-processing. As all the models we tried required data of a
fixed size, every file was extended (or padded) to match the length of the
longest sequence. In other words, every skeleton file was transformed into an
array of $3\times n+1$ columns and $t_{max}$ rows by adding some rows full of
zeros, $t_{max}$ being length of the longest skeleton file in the dataset. In
the dataset we worked with, we found $t_{max}=100$. Figure 10 shows the
results of the pre-processing on a skeleton file. This two dimensional data
sample could be fed to the CNN as is but the two other models (FC and LSTM)
required the data to be flattened. Each skeleton file was then transformed
into a vector of size $3\times n\times t_{max}=237000$ (as shown in Figures
2a, 2c).

Figure 10: Examples of a matric of 3D skeleton data which has been transformed
following our data pre-processing scheme.

###  3.3 Experiment Settings

The target of the models is either 0 or 1, representing classes Mono (one hand
movement) and Bi (two hands movement) respectively. The 111 labelled skeleton
files available in the dataset were used as follows:

  * •

Supervised learning (all three models): 66 samples were used for training, 11
for validation and 34 testing.

  * •

Self-supervised learning (applied to the FC and CNN models): 101 samples were
used for the unsupervised training. In this simulated settings, the supervised
model for the downstream task was trained with 5 samples for the train dataset
and 5 samples for the validation dataset. The 101 data used for the
unsupervised part were reused for the final test.

We used accuracy as the metric to evaluate each model. As a way to analyze and
explain the results of the predictions, we applied the Grad-CAM method to our
trained CNN in supervised learning. The source code of the methods (including
the demonstration) can be found in this link
222https://github.com/FABallemand/ProjetCassiopee.

##  4 Results

###  4.1 Supervised Learning

According to our results, the three models we tested are able to learn
properly to do the binary classification on the skeleton data. All of them can
predict with high accuracy if a movement requires one or two hands. We
observed a test accuracy of 100% for the CNN and LSTM models. The FC model
also performed very well with a test accuracy of 97% as shown in Table 1. In
this respect, it should also be noted that the time taken to train each model
was under ten minutes, with the FC being the fastest model to train since it
is the one containing the lowest number of parameters.

| FC | CNN | LSTM
---|---|---|---
Test accuracy | 97% | 100% | 100%
Accuracy (Test+Validation) | 96.1% | 97.4% | 96.1%
Misclassified data |  Venir, Chambre, Couche |  Venir, Petit |  Venir, Chambre, Chef
Table 1: Results comparison for supervised learning on MocapLab skeleton data.

The CNN has overall the best performances with only two misclassified data
over the test and validation datasets, whereas the other models have three
misclassified data (exhibited in Table 1). The loss curves in Figure 11b
suggest that the CNN is also the most stable model.

(a) Result of FC

(b) Result of CNN

(c) Result of LSTM

Figure 11: Training results (accuracy and loss graph for training and
validation data during training, accuracy and confusion matrix on test data).

The data Venir is misclassified by all three models (label 0, predicted as 1).
One explanation could be that the movement uses only one hand but the unused
hand is not as low as in the other movements. Additionally, the moving hand
comes near the head for less than ten frames, meaning it is basically at the
same height level as the unused hand for the most part of the gesture (cf.
Figure 12).

(a) Frame 22

(b) Frame 28

(c) Frame 35

Figure 12: Venir skeleton visualisation.

###  4.2 Self-supervised Learning

Training method | Data utilisation | Test accuracy
---|---|---
Supervised CNN (10% of data) | 5% train, 5% validation, 90% test | 83%
Self supervised CNN | 90% for unsupervised learning, 5% train, 5% validation | 93%
| Note : data used for unsupervised training is also used for testing |
Table 2: Result comparison for SL and SSL for the CNN model on Mocaplab
skeleton data.

The experimentation with our three models has demonstrated the successful
execution of binary gesture classification tasks. On this, we conducted an
experiment using our CNN, training it with 5 data samples in both the training
and validation datasets, while reserving the remaining 101 samples for
testing. The results revealed a test accuracy of 83% (Table 2), which, as
anticipated, was lower than the previous outcomes.

Nonetheless, in the simulated scenario, self-supervised learning method can
effectively capitalize on the scarcity of labelled data by using unlabelled
data. We employed SSL to train our CNN, and achieved a test accuracy of 93%
(Table 2), which indeed increase the accuracy obtained as opposed to sole use
of SL, hence this demonstrates the potential of SSL. Finally, other
exploration experiment conducted on a larger dataset for multi-class object
classification is presented in Appendix A.

###  4.3 Visualisation of the models focus

In Figure 13, we can see two examples of Grad-CAM visualisation, 13a is for a
one hand sign (Mono label) and 13b is for a two hands sign (Bi label). The
example for Mono shows that the CNN focuses on the joints of the moving hand
to make its prediction. Most of the one hand sign data show similar patterns.
The fact that the model performs the prediction by analyzing the joints in the
moving hands (joints) provides some insights where the model is focusing (i.e.
the relevant joints) [8]. However, all two hands signs Grad-CAM visualisation
results are identical as shown on the Figure 13b, meaning the focus of the
model is focusing on the static area (such as legs and torsos). As such,
improving the explainability of the predictions of the CNN for Bi labelled
data can be considered as a relevant task in the future.

(a) visualisation for Mono.

(b) visualisation for Bi.

Figure 13: Grad-CAM based visualisation of models focus with respect to the
joint coordinates (green dots are the point where models focus) for one and
two hands sign.

##  5 Conclusion

In this work, we explore the application of deep self-supervised learning to
the field of automatic gesture recognition, specifically focusing on hand
gestures using 3D moving skeleton data. Our objective is to assess the
performance of various neural network architectures — fully connected network,
convolutional neural networks, and long short-term memory networks — to
recognize gestures by analyzing 3D joint coordinates over time. This approach
provides a comprehensive representation of motion, which is crucial for
accurately capturing and understanding complex gestures.

Our findings indicate that supervised learning techniques yield high accuracy
in gesture classification with this data type, where CNN and LSTM models
achieve perfect test accuracy of 100%, and the FC model follows closely with
97% accuracy. The CNN model demonstrate superior performance, particularly in
handling the temporal and spatial complexities of 3D skeleton data, making it
well-suited for the task.

The application of self-supervised learning methods shows promising potential,
particularly in scenarios with limited labelled data. By pre-training models
on simulated task of using unlabelled 3D skeleton data, we observe enhanced
performance when transitioning to supervised tasks. Although the results for
self-supervised learning are slightly less conclusive, they underscore the
importance of leveraging large amounts of unlabelled data to improve model
training efficiency and accuracy.

Our experiments with Grad-CAM provide valuable insights into model
interpretability and understanding the decision-making process of our CNN
model, particularly for gestures involving single-hand movements. However, the
visualisation for two-handed gestures indicates a need for further refinement
in model interpretability, as the focus tends to misaligned with the intended
regions of interest.

Moving forward, future work could explore more complex classification tasks
involving multiple gesture categories and a larger dataset to improve
generalization and model robustness. Additionally, investigating advanced
self-supervised learning techniques and integrating multi-modal data sources
could further enhance the performance of gesture recognition systems,
broadening their application in areas such as sign language recognition,
virtual reality, healthcare, and human-computer interaction.

As such, this work demonstrates the significant potential of deep learning and
self-supervised learning techniques in advancing automatic gesture
recognition, paving the way for more accessible and intuitive human
interactions and human-computer interfaces.

## Acknowledgments

We would like to thank Zied LAHIANI at Department of ARTEMIS for the technical
support. The authors would also like to thank Remi BRUN and Boris DAURIAC from
Mocaplab studio for providing the captured high quality data which is used for
study case on this work.

## References

  * [1] Decky Aspandi, Adria Mallol-Ragolta, Björn Schuller, and Xavier Binefa.  Latent-based adversarial neural networks for facial affect estimations.  In 2020 15th IEEE International Conference on Automatic Face and Gesture Recognition (FG 2020), pages 606–610. IEEE, 2020.
  * [2] Decky Aspandi, Federico Sukno, Björn W Schuller, and Xavier Binefa.  Audio-visual gated-sequenced neural networks for affect recognition.  IEEE Transactions on Affective Computing, 14(3):2193–2208, 2022\.
  * [3] Ryan Nash Keiron O’Shea.  An introduction to convolutional neural networks.  arXiv:1511.08458, 2015.
  * [4] Eric Rothstein Morris Ralf C. Staudemeyer.  Understanding lstm – a tutorial into long short-term memory recurrent neural networks.  arXiv:1909.09586, 2019.
  * [5] Tim Schneider, Chen Qiu, Marius Kloft, Decky Aspandi Latif, Steffen Staab, Stephan Mandt, and Maja Rudolph.  Detecting anomalies within time series using local neural transformations.  arXiv preprint arXiv:2202.03944, 2022.
  * [6] Ting Chen, Simon Kornblith, Mohammad Norouzi, and Geoffrey Hinton.  A simple framework for contrastive learning of visual representations.  In ICML, pages 1597–1607. PMLR, 2020.
  * [7] Ramprasaath R Selvaraju, Michael Cogswell, Abhishek Das, Ramakrishna Vedantam, Devi Parikh, and Dhruv Batra.  Grad-cam: Visual explanations from deep networks via gradient-based localization.  In ICCV, pages 618–626, 2017.
  * [8] Ciheng Zhang, D. Aspandi, and Steffen Staab.  Predicting eye gaze location on websites.  In Proceedings of the 18th International Joint Conference on Computer Vision, Imaging and Computer Graphics Theory and Applications, VISIGRAPP 2023, Volume 4: VISAPP, Lisbon, Portugal, February 19-21, 2023, pages 121–132, 2023.
  * [9] Umberto Michelucci.  An introduction to autoencoders.  arXiv:2201.03898, 2022.
  * [10] Raja Giryes Dor Bank, Noam Koenigstein.  Autoencoders.  arXiv:2003.05991, 2020.
  * [11] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 770–778, 2016.

##  Appendix A Self-supervised Learning for Object Classification

The experiment from Section 4.2 shows that unlabelled data can be leveraged in
the context of hand gesture recognition thanks to SSL. Having at our disposal
only a few skeleton files labelled for binary classification, we also wanted
to know if SSL could be scaled up with more data and for multi-class
classification. For that matter, we used the RGB-D Object Dataset, a very
large dataset of RGB, depth and segmented images (each of a resolution of
640x480 as shown in Figure 14) with object class and position labels. This
dataset contains more than 124000 data samples distributed across 53 object
classes. We implemented SL for object classification, UL for image
reconstruction (with and without contrastive learning and SSL with object
classification as the downstream task, as described in Section 2.2, using only
RGB images and class labels. We used an auto-encoder [9, 10] for image
reconstruction during the unsupervised task with ResNet-18 [11] as the encoder
part. The same pre-trained ResNet-18 was then fine tuned and used for the
downstream task of image classification.

Figure 14: Data sample from the RGB-D Object dataset after pre-processing
step.

The results presented in Table 3 show that SL achieve great to moderate
accuracy results when trained on the entire or on half of the dataset.
However, the results of SSL are below our expectations. While there is no
noticeable accuracy loss, results of SSL are on par with the results obtained
with SL with an accuracy of around 50%. In this particular case, the effect of
SSL still is not achieved as larger training time and more (unlabelled) data
than SL are required. This thus yields the slightly reduced accuracy obtained.

It should be noted that this behavior may be due to our implementation of SSL
on this dataset. While our model has a lot in common with the SimCLR framework
presented in [6], our implementation is different notably at the level of the
contrastive learning. Our implementation relies on labels in order to select
data samples similar or different than the sample currently processed whereas
SimCLR consider every data sample of a batch to be different from the others.
Our implementation should enhance the separability of workion in the latent
space but it requires labels which makes it unsuitable if the goal is to
leverage large amount of accessible unlabelled data. Due to some technical
difficulties and a lack of time, we did not investigate the issue further.
However, we believe that this task can be further be explored to reveal more
result.

| Training Method | Model | Data | Metric
---|---|---|---|---
Accuracy (%) | F1
1 | Supervised Learning | Pretrained ResNet18 | Split A + B | 93.47 | 0.9320
2 | Supervised Learning | Pretrained ResNet18 | Split A | 49.15 | 0.3540
3 | Unsupervised Learning |  | Autoencoder with
---
pretrained ResNet18
Split B |  |
4 | Unsupervised Learning with Contrastive Loss |  | Autoencoder with
---
pretrained ResNet18
Split B |  |
5 | Self-supervised Learning |  | Pretrained ResNet18
---
from 3
Split A | 48.29 | 0.3412
6 | Self-supervised Learning |  | Pretrained ResNet18
---
from 4
Split A | 48.37 | 0.3446
Table 3: Results comparison for SL and SSL for object classification on RGB
data from the RGB-D Object dataset.


# Visualizing particle networks in granular media by in situ X-ray computed
tomography

Matthias Ruf  Institute of Applied Mechanics (CE), University of Stuttgart,
Stuttgart, Germany  Corresponding author: Matthias Ruf,
<EMAIL_ADDRESS>Kianoosh Taghizadeh  Institute of
Applied Mechanics (CE), University of Stuttgart, Stuttgart, Germany
Multiscale Mechanics (MSM), University of Twente, Enschede, The Netherlands
Holger Steeb  Institute of Applied Mechanics (CE), University of Stuttgart,
Stuttgart, Germany  SC SimTech, University of Stuttgart, Stuttgart, Germany

(December 23, 2021; December 23, 2021)

###### Abstract

In the present contribution, cylindrical samples consisting of monodisperse
soft (rubber) and stiff (glass) particles are pre-stressed under uniaxial
compression. Acoustic P-waves at ultrasound frequencies are superimposed into
prepared samples with different soft-stiff volume fraction. Earlier
investigations showed the importance of particles networks, i.e. force chains,
in controlling the effective mechanical properties of particulate systems,
e.g. elastic moduli. Measured P-wave modulus showed a significant decline
while more soft particles are added, i.e. higher rubber fractions, due to a
change in microstructure. However, for small contents of soft particles, it
could be observed that the P-wave modulus is increasing. For the understanding
of such kinds of effects, detailed insight into the microstructure of the
system is required. To gain this information and link it to the effective
properties, we make here use of high-resolution micro X-ray Computed
Tomography (µXRCT) imaging and combine it with the classical stiffness
characterization. Both performed in situ meaning inside the laboratory-based
XRCT scanner. With µXRCT imaging, the granular microstructure can be
visualized in 3d and characterized subsequently. By post-processing of the
µXRCT data, the individual grains of the particulate systems could be uniquely
identified. Finally, the contact network of the packings with low and high
rubber contents which connects the center of particles was established to
demonstrate the network transition from stiff- to soft-dominated regimes. This
has allowed for unprecedented observations and a renewed understanding of
particulate systems. It has been demonstrated that µXRCT scans of particles
packings can be analyzed and compared in 3d to gain extensive information on
the scale of the single particles. Here, the in situ setup and workflow from
the start of acquiring images in situ till the post-processing of the image
data is explained in detail and demonstrated by selected results.

_Keywords_ granular media $\cdot$ wave propagation $\cdot$ X-ray computed
tomography

##  1 Introduction

Dense granular media are of widespread importance in a number of applications,
ranging across time and length scales from geophysical earthquakes in the San
Andreas fault zone, to plastic powder sintering in the 3d printer on your
desk. Granular mixtures are of interest for a large number of fields,
materials, and applications, including mineral processing, environmental
engineering, geomechanics and geophysics, and have received a lot of attention
in the last decades. A specific example in geotechnical engineering is the
increasing incorporation of recycled materials (e.g. shredded or granulated
rubber, crushed glass) often used into conventional designs and soil
improvement projects [4, 14, 29]. Moreover, sophisticated mixtures of asphalt
and concrete are widely used to construct roads [16, 17, 53, 68, 69].

The mechanical behavior of granular materials is highly nonlinear and involves
plastic deformations also for very small strain due to rearrangements of
particles [27, 32, 34]. On the other hand, the concept of an initial purely
elastic regime at small strains for granular assemblies is an issue still
under debate in the soil mechanics community [6, 10, 67]. On the other hand,
approaches that neglect the effect of elastic stored energy are also
questionable, i.e. all the work done by the internal forces is dissipated.
Features visible in experiments, like the propagation of acoustic waves, can
hardly be described without considering an elastic regime [58, 59]. The small-
strain stiffness provides useful soil information, which is relevant to a wide
range of engineering applications including the design of foundations
subjected to dynamic loading, process monitoring, liquefaction assessment, and
soil improvement control. Within the small strain region, the geomaterial
exhibits linear-elastic behaviour and the moduli are independent of strain
amplitude [9, 21, 27]. The elastic range depends on confining pressure for
non-plastic soils.

Probing a granular media with (ultra) sound waves gives useful information on
the state, the structure and the mechanical properties of the bulk media as
there is a one to one relation between a wave speed and small-strain stiffness
of packings. Estimation of stiffness has traditionally been made in a triaxial
apparatus using precise displacement transducers or resonant column devices
[1, 47]. However, these methods have a disadvantage of destructing samples,
whereas ultrasonic measurements (by propagation of an elastic wave) are widely
accepted for their rapid, non-destructive, and low-cost evaluation methods.
Earlier studies using a wave propagation technique have shown that the
dissipative, elastic and lightweight properties of materials (like soils,
asphalt, etc.) can be enhanced by deliberately adding dissipative, soft, light
inclusions of various types and compositions. Therefore, such a tailored
system with improved material properties leads to new application designs.

There exists a diverse array of experimental imaging techniques which can be
exploited in order to investigate particulate systems. One of them is
attenuation-based micro X-ray Computed Tomography (µXRCT) [7, 8, 22, 46, 56].
µXRCT imaging is in general a non-destructive imaging technique that offers
the possibility to visualize the internal structure of objects. In contrast to
other microscopy imaging methods like optical microscopy, it provides a
3-dimensional (3d) representation of the investigated sample. The method is
based on the physical effect of X-ray attenuation which depends among other
things on the atomic number and consequently on the chemical composition in
each material point. The final, so called reconstructed 3d volume resulting
from an XRCT scan is typically represented as a stack of 2-dimensional (2d)
gray values images. Through subsequent post processing of the 3d raw image
data, different material phases can be identified and subsequently segmented
[18, 45, 49, 63]. Based on the segmented 3d volume, quantifications on
different length scales can be performed. One example to be mentioned on the
macroscopic level is the determination of the volume fractions of the
individual phases as well as their distribution. But it is also possible to
obtain detailed quantities such as the coordination number of individual
particles in the case of granular materials. Consequently, experimental
results and observed phenomena obtained by classical characterization methods
can be better and more comprehensively understood by help of µXRCT. If µXRCT
imaging is combined with classical characterization methods and performed at
the same location we talk about in situ µXRCT imaging [5]. Motivated by the
example of a discovered not fully understood effect shown by tailored granular
media [62], a workflow for in situ investigation is demonstrated in this
contribution.

In the following, we start with a short overview of image-based
characterization using X-ray computed tomography. The combination (in situ) of
this method with the traditional ultrasonic through-transmission measurement
technique is the focus of Section 3. We apply this to investigate monodisperse
particulate systems consisting of weak (rubber) and stiff (glass) particles of
different volume fractions. We show a possible realization including the
required hardware as well as the overall workflow. Based on this, some
selected results of the investigation are provided in Section 4 to give a
motivation of what is in general possible. The acquired experimental data for
one selected mixture is open access published for demonstrative purposes. The
proposed in situ workflow is discussed in Section 5. A summary of the
presented work is given in Section 6.

##  2 Image-based characterization using X-ray computed tomography

Figure 1: Illustration of the principle of attenuation-based X-ray computed
tomography (XRCT) imaging based on the scan of a cylindrical sample of
monodisperse stiff (glass) and soft (rubber) particles.

The capability which micro X-ray Computed Tomography (µXRCT) as a non-
destructive imaging technique offers are meanwhile well known in many research
domains and it is becoming more and more a routine microscopy technique [11,
26, 28, 33, 36, 50, 56, 57]. The method is based on the revolutionary
discovery of X-ray radiation by Wilhelm Röntgen in 1895 [46]. X-rays are high-
energy electromagnetic waves with a wavelength in the range of 0.01 nm to 10
nm. They have a long absorption length and interact with matter through
different mechanisms (photoelectric effect, Compton scattering,
Rayleigh/Thomson scattering, and pair production) [8]. Depending on the
material penetrated, X-rays are attenuated to different degrees. Assuming a
monochromatic X-ray beam, the Beer-Lambert law formulates mathematically the
transmitted intensity $I$ of an X-ray beam penetrating an object along the
straight ray path from 0 to $\tilde{s}$ by

$I(\tilde{s})=I_{0}\mathrm{e}^{-\int_{0}^{\tilde{s}}\mu(s)\mathrm{d}s}$ (1)

where $I_{0}$ is the initial intensity of the X-ray beam and $\mu(s)$ are the
unknown local linear attenuation coefficients along the ray path [7, 8, 56].
If an object is irradiated from one side and the attenuated X-ray beam is
measured on the opposite side by an X-ray detector, a so-called radiogram or
projection image results. The distribution of the X-ray intensity $I$ is
typically represented as a gray value image as displayed in Figure 1(a)
showing one projection image of a cylindrical sample containing monodisperse
soft (rubber) and stiff (glass) particles. According to Equation (1) each
pixel of such an image represents the intensity $I$ of all attenuated X-ray
photons which are captured by the specific sensor pixel during the adjusted
exposure time. The darker the gray value the less the intensity $I$ and vice
versa in the radiogram. Since each pixel contains just the integral
information for one specific X-ray path, it is not possible to back-calculate
the unknown linear attenuation values in the material points from one
radiogram. However, if a sufficiently large set of projections is taken from
different directions, cf. Figure 1(b), enough information is available to
calculate the attenuation coefficient for each material point, cf. Figure
1(c). The 3d distribution of the calculated attenuation coefficient then
represents the inner structure of the investigated object. It is to mention,
due to different reasons, that the absolute values of the 3d images are mostly
only correlating with the theoretical attenuation coefficients and are not
identical.

The general method for the determination of the 3d structure based on numerous
projection images is called reconstruction [7, 22]. To capture projections
from different directions, in µXRCT imaging typically the sample is rotated,
cf. Figure 1(a) and (b). The number $n$ of required equidistant projection
angles $\varphi_{i}$ over one turn of the sample correlates with the number of
detector pixel columns. The size/discretization of the 3d volume is given by
the resolution of the detector. The base area (x- and y-direction) corresponds
with the number of detector pixels in the horizontal direction and the maximum
height (z-direction) with the number of vertical detector pixels, cf. Figure
1(c). Corresponding to the term pixel in 2d, in 3d the term voxel is used
[45]. Assuming an identical voxel edge length, the edge length of one voxel is
given by the ratio of the detector pixel size and the applied geometric
magnification $M_{geo.}$. The adjustment of the magnification in a cone beam
system is performed by geometric magnification. If we assume the source-
detector-distance ($SDD$) is fixed, the geometric magnification $M_{geo.}$ is
adjusted by varying the source-object-distance ($SOD$) and is given by the
ratio $M_{geo.}=SDD/SOD$. The maximum theoretically achievable spatial
resolution correlates with the underlying focal spot size of the used X-ray
source for the adjusted power setting [56]. A micro-focus X-ray source
typically has a focal spot size in the micrometer range and results in a
theoretical system resolution which is also in the micrometer range.

For the reconstruction itself, it is distinguished between two groups of
reconstruction techniques. On the one side, there are analytical
reconstruction techniques and on the other side algebraic ones, e.g. [7, 22].
While the first group of methods tries to solve the inverse problem
analytically, the second group treats it as an optimization problem. Both have
their respective advantages and disadvantages. The shape of the X-ray beam has
a significant influence on the reconstruction, in particular for the
analytical techniques. The ideal case of a parallel X-ray beam is only present
in beamlines of synchrotron radiation facilities [42] and not given in
laboratory XRCT systems. To account for the conical beam shape in laboratory
XRCT systems, the so-called FDK algorithm is mostly applied. The FDK approach
is a practical analytical cone-beam algorithm which goes back to Feldkamp
L.A., Davis L.C., and Kress J.W. [13]. Independent of the used reconstruction
techniques, the final raw data set consists of a stack of slices called
tomograms, cf. Figure 1(c). The voxel values, typically represented as gray
values, represent their influence on X-ray attenuation. The brighter the are
the higher their influence on the X-ray attenuation corresponding with the
material density. The bright gray voxels in Figure 1(c) correspond to the
stiff (glass) particles and the dark gray ones to the (soft) rubber particles
and the nearly black voxels reproduce the air phase. The reconstructed data
set contains a certain amount of noise due to the procedure as well as
potentially other hardly avoidable artifacts [8, 46, 56].

By subsequent image processing, based on the 3d raw volume, different kinds of
information can be extracted from the data set. For this purpose, a so-called
segmentation of the different phases is usually carried out. Segmentation is
the partitioning of a gray-scale image into disjoint regions that are
homogeneous with respect to some characteristics [18, 46, 49, 63].

In the simplest case, this can be achieved by defining a threshold value,
since the individual phases differ in terms of the average gray value. To
determine the related threshold values, a histogram of the underlying gray
value distribution of the raw data set is usually employed. The most well-
known approach based on this scheme is the “Otsu method” [40]. However, due to
the intrinsic present noise in the raw data set, incorrect assignment of
different individual voxels occurs. To account for this, various filtering
techniques are used [18, 49, 63]. To mention is that in the last years more
and more image segmentation using machine learning techniques has been
robustly studied as a new approach due to its well-known benefits [23, 25].
Which phases have to be separated depends on the scientific question.
Exemplary, in Figure 1(d) the separation of the soft (rubber) and stiff
(glass) particles, as well as the pore space (transparent), is showed.

By repeating the scanning process as well as the image processing, the
evolution of physical processes (e.g. deformation, crack initiation, and
growth) over time can be observed. The final result is one stack of tomograms
(3d image) for each time step. What kind of processes can be observed depends
on the required acquisition time for one complete scan. While acquisition
times of less than one second for one complete scan are achievable in
beamlines of synchrotron radiation facilities [5, 15, 42], the acquisition
time in a standard lab-based XRCT system is between some minutes up to several
hours, correlating with the required scan quality. This is one limitation of
physical processes which can be studied. However, processes which can be
stopped in an equilibrium state are also open for investigation in a
laboratory system. For time-resolved scans, the term 4-dimensional (4d)
scanning (space + time) is frequently used.

##  3 In situ experimental testing

The combination of XRCT imaging with classical mechanical characterization
methods is an excellent technique to better understand the behavior of
materials as well as physical processes in general, cf. [5, 54]. In this
context often the term “in situ X-ray computed tomography” is used. “In situ”
is a Latin phrase and translates literally to “on-site" or “on place" and is
the antonym of “ex situ”. In the experimental mechanics’ context, it describes
the way the measurement is taken, highlighting that the measurements are
acquired in the same place the phenomenon is occurring without
removing/installing the sample each time. Thus, in situ XRCT means that the
mechanical characterization is performed inside the scanner. It offers the
possibility to enrich the information on the macroscopic level, e.g. resulting
from acoustic measurements as demonstrated within this contribution, with
information on the microscopic level resulting from 3d imaging. In the
following, we provide an overview of the applied workflow and used hardware to
investigate particulate systems with acoustic wave measurements combined with
in situ XRCT imaging.

###  3.1 Stiffness determination based on wave propagation measurement

Mechanical (sound) waves are disturbances that propagate through space and
time in a medium in which deformation leads to elastic restoring forces. This
produces a transfer of momentum and energy from one material point to another,
usually involving little or no associated mass transport if the amplitude is
small enough. The P-wave, or primary wave, is the fastest and the first wave
detected by seismographs. They are able to move through both solid rock as
well as through liquids. P-waves are compressional or longitudinal waves that
oscillate the ground back and forth along the direction of wave travel, in
much the same way that sound waves (which are also compressional) move air
back and forth as the waves travel from the sound source to a sound receiver.
In a longitudinal wave, the particle displacement is parallel to the direction
of wave propagation [35, 38, 52, 60, 37].

Figure 2: Schematic view of wave propagation measurement setup [62].

Quantitative ultrasonic measurement has been widely-used by different
disciplines ((geo-)physics, soil and geomechanics, materials sciences,
mechanical and civil engineering) to describe the small strain stiffness
behavior of a particulate system, as can be found in literature [19, 20, 64,
55]. Velocity testing by through transmission ultrasound methods has gained
popularity due to its relative ease of obtaining the modulus of a sample.
Figure 2 shows a schematic view of the used low-frequency ultrasound setup
integrated into an oedometer cell including its electrical pieces, cf. Figure
3 and Figure 4. For our highly attenuating particulate system, it consists of
a pair of 100 kHz P-wave broadband piezoelectric transducers (Olympus-
Panametrics Videoscan V1011), an ultrasonic square wave pulser/receiver unit
(Olympus-Panametrics 5077PR) and a digital oscilloscope (PicoScope 5444B).
Piezoelectric transducers are used to determine the small-strain compression
stiffness, $M$, of the granular system by determining the velocity of
mechanical waves (“speed of sound”) through the tested samples. The
transducers are generally used in pairs when one transducer operates as a
transmitter and the other as a receiver. The transmitting transducer is
generally embedded at one end of the particulate soil sample and the receiving
is, aligned with the transmitter, embedded at the other end. This kind of
pair-wise transducer arrangement is also named ultrasonic through transmission
measurement technique in the literature.

This allows to probe the stiffness of the soil or particulate sample along a
given stress path. The transmitting transducer transform the input electric
signal to an acoustic wave and sends it through the medium. The receiving
transducer receives the propagated mechanical sound wave and transforms it
back to an electronical signal which could be acquired by the digitizer,
compare [30, 48] for technical details on the transducers. Transducers are
often incorporated in geomechanical and geophysical laboratory testing
equipment such as in triaxial or oedometers cells.

The P-wave velocity ($V_{P}$) can be calculated from the travel time
($t_{P}$), given the height ($H$) of the sample in the actual configuration,
$V_{P}=H/t_{P}$. Knowing the elastic wave velocity ($V_{P}$) and knowing the
total mass density ($\rho_{B}$) of the particulate system, the longitudinal
P-wave stiffness of the sample ($M$) is determined, $M=\rho_{B}V_{P}^{2}$. In
our investigations, the sample consists of rubber and glass particles
$\rho_{B}=\nu\rho_{r}+(1-\nu)\rho_{g}$, where $\rho_{r}$ and $\rho_{g}$ are
the mass densities of the rubber and glass particles, respectively. The volume
fraction of the rubber particles is denoted with $\nu$.

###  3.2 XRCT system and oedometer cell

The experimental study was performed in an open, modular, and flexible lab-
based micro X-ray Computed Tomography (µXRCT) system with an intrinsic large
space allowing for the installation of a mechanical testing device required
for in situ investigations. A detailed description of the system can be found
in [43]. The system presented in [43] was extended by the integration of a
load frame for XRCT investigations, cf. Figure 3(a).

Figure 3: Experimental in situ setup: (1) X-ray source, (2) X-ray flat panel
detector, (3) Universal testing machine with synchronized rotary tables, (4)
Oedometer cell with integrated P-wave transducers containing monodisperse
glass and rubber particles, (5) Square wave/pulser receiver unit and PC
oscilloscope for ultrasonic measurements.

For this, a Zwick 1445 10 kN Universal Testing Machine (UTM), refurbished from
Doli Elektronik GmbH, Germany, and equipped with a modern EDC 222V controlling
system, is employed. In the UTM, two rotatory tables with high loading
capacities (XHuber 1-Circle Goniometer 411-X3W2) are integrated. The rotatory
tables are pre-stressed with the controlled force $F_{stat.}$ applied to the
sample during the image acquisition process. The X-ray source of the system is
an open micro-focus tube with tungsten transmission target (FineTec FORE
180.01C TT) from FineTec FineFocus Technology GmbH, Germany. For the image
acquisition, a Dexela detector 1512NDT with GOS-based DRZ Standard
scintillator option from PerkinElmer, Inc., Waltham, MA, USA is applied.

To perform ultrasonic through transmission measurements in situ, a low X-ray
absorbing oedometer cell was designed and manufactured in-house. Compared to a
triaxial cell, only the vertical stress in the axial direction is controlled.
In Figure 3 the oedometer cell installed in the XRCT-system and loaded by the
employed universal testing device is presented.

Figure 4: Three-quarter section view of the low X-ray absorbing oedometer cell
with integrated ultrasonic transducers for in situ measurements.

Further, the electrical devices for the ultrasonic testing corresponding to
the schematic Figure 2 can be observed. A detailed view of the oedometer cell
design is provided in Figure 4. The oedometer cell consists of three
components. The top (1) and bottom (3) parts, both with integrated P-wave
transducers, transfer the applied axial compression force $F_{stat.}$ to the
sample and results in the axial stresses denoted as pressure $p_{stat.}$. As
transducers, the pair of 100 kHz P-wave broadband piezoelectric transducers
(Olympus-Panametrics Videoscan V1011) is employed, cf. Subsection 3.1. For the
coupling of the transducers to the Poly (methyl methacrylate) (PMMA) holders
of the oedometer cell, an adequate couplant fluid is used. On both sides the
wave travels through a 10 mm long PMMA distance, before and after the wave
enters and leaves the sample. The contact pressure $p_{cont.}$ in between the
transducers and the PMMA holders can be adjusted by an adjusting wheel. The
sample is held in a rigid confining ring (2) made out of PMMA which prevents
lateral displacement of the investigated mixture. Between the PMMA ring and
the top and bottom part is a small gap in the radial direction to ensure that
the emitted waves propagate through the sample. PMMA is used for the confining
ring, as it has a very low attenuation coefficient. The inner diameter $D$ of
the PMMA ring is 80 mm. The sample height $H=L$ can be varied from 70 mm to 80
mm.

Two capture the whole cell content (diameter 80 mm; height 80 mm), a geometric
magnification $M_{geo.}=1.36$ was set for all particulate systems. This leads
to a field of view of 106.92 mm $\times$ 84.48 mm. Since we are not interested
in small features (bead diameter 4 mm) all projection images were acquired in
$2\times 2$ detector binning mode. This means, that the detector’s full
resolution of 1944 $\times$ 1536 pixels with 74.8 µm pixel pitch and 14-bit
pixel depth is reduced to 972 $\times$ 768 pixels with a pixel size of 149.6
µm. On one side, this significantly improves the signal-to-noise ratio (SNR),
and on the other side, the final tomogram data set size is reduced by factor
8. Both simplify the subsequent image processing significantly without the
loss of information required for the presented study. For all scans, an
acceleration voltage of 110 kV with an acceleration flux of 110 µA for the
X-ray source were set. Combined with a detector exposure time of 1000 ms and
1440 equidistant projection angles. Further, a detector bad pixel compensation
as described in [43] was employed using two different projection positions for
each projection angle. The final tomogram stacks have a size of 972 $\times$
972 $\times$ 768 voxels with the uniform voxel edge length of 110 µm. The 3d
reconstruction of all scans was performed with the software Octopus
Reconstruction (Version 8.9.4-64 bit) [66] using the Filtered Back Projection
(FBP) method [22] in combination with the FDK reconstruction algorithm [13].
To account for common artifacts in XRCT imaging (ring artifacts and beam
hardening), different types of filters were used. Based on the tomograms,
segmentation is performed to distinguish between the rubber and glass beads
and the remaining pore space. The segmentation workflow is described in detail
in Section 3.3.

###  3.3 Image analysis

Figure 5: A top view of a xy-layer of particles for samples prepared with
$\nu=0.1$ (left column) and 0.5 (right column). Top to bottom shows the steps
taken from a reconstructed raw image (first row) to segmented (second row),
and labeled (third row). The 3d raw data set for $\nu=0.5$, here presented
with a single slice (“slice_xy_0322.tif”), can be found in [44]. Figure 6: In
situ workflow summary.

Here, we present steps taken using the Dragonfly software version 2020.2 for
Windows from Object Research Systems (ORS, Montreal, QC, Canada) [39], for the
segmentation and the subsequent quantification of mixture samples. An
important step before interpreting the data is to identify individual
particles. To increase the quality of images and reduce unwanted noise, one
can apply different image filtering techniques before segmenting the phases,
i.e. separating different objects in an image. Here, a 3d-median smoothing
function with a kernel size of seven was applied on image stacks to improve
the resolution of pixels. Next, pixels of images were manually divided into
three different regions of interest (air, rubber and glass) which are used as
boundaries of objects. After segmenting samples into three separate regions
(air, rubber, and glass), watershed transform technique is employed on defined
markers of individual particles to label every particle in the separated
regimes of rubber and glass. The watershed is a classical region-based
algorithm that has its origins in mathematical morphology [3, 12, 51, 65] used
for segmentation. Starting from user-defined markers, the watershed algorithm
treats pixels values as a local topography (elevation). The algorithm floods
basins from the markers until basins attributed to different markers meet on
watershed lines. In many cases, markers are chosen as local minima of the
image, from which basins are flooded.

Figure 5 shows a z-direction view of a layer of particles prepared at two
rubber fractions, $\nu=0.1$ (left) and 0.5 (right), under 800 N load in the
oedometer cell. Light and dark gray pixels represent glass and rubber
particles respectively in the first row. Second row shows the segmented pixels
in the regions of interests (rubber (blue) and glass (red)) for which deep
learning was employed. Thanks to the watershed transform technique, particles
of segmented regions are labeled individually (third row).

###  3.4 Investigated particulate samples and measurement workflow

All investigated samples are made of monodisperse glass and rubber beads with
4 mm diameter. Further information on the material characteristics for both,
glass and rubber beads, can be found in [62]. Particles are poured carefully
into the oedometer cell described in Section 3.2. Samples are prepared at
different rubber volume fractions, $\nu=v_{r}/v_{c}$, from $\nu=0$ to 1 in 0.1
increments, where $v_{r}$ is the volume possessed by the rubber particles and
$v_{c}$ is the total volume of particles in oedometer cell. All mixtures are
tested in the custom-made oedometer cell with the sample height and diameter
equal to 80 mm, cf. Figure 3 and Figure 4 ex-situ without imaging. The
prepared samples are uniaxially compressed in the axial direction via the top
piston of the oedometer cell in subsequent force increments from $F=\mbox{200
N}$ to 1000 N. At each instant step, the system was relaxed for some time at
the defined force to overcome the creep behavior of soft particles. Then, at
each intermediate load step, a high voltage burst signal (P-wave) is excited
from the top cap (the sound source) transducer and the bottom transducer
collects the signal. To remove the influence of network configurations and
user errors on outputs, each experiment was repeated at least 3 times by
mounting and dismounting particles into the cell. For the calculation of the
system stiffness besides the wave propagation time $t_{P}$ the exact height
$L$ of the compressed packing is required. This is determined based on the
loading protocol of the UTM. Since about 2 h are needed for one XRCT scan,
only a subset of the possible combinations of the parameter space (rubber
volume fraction $\nu$ and preload $F_{stat.}$ were scanned in situ. XRCT
imaging was performed for a preload $F_{stat.}=\mbox{400 N}$ and
$F_{stat.}=\mbox{800 N}$ in combination with rubber volume fractions
$\nu=\\{0.1,0.2,0.3,0.4,0.5,0.6\\}$. Otherwise, the procedure is identical to
the previously performed ex-situ measurements. During the image acquisition,
it is switched from force-control to displacement-control to avoid slight
movements of the particles due to relaxation processes. The selection is based
on previous knowledge from comparable experiments in a triaxial cell without
XRCT imaging [62] as well as the before performed ex-situ measurements in the
oedometer cell. A summary of the applied in situ workflow is given in Figure
6.

##  4 Results - insight into granular media

Measurement of ultrasound velocity $V_{P}$ provides complementary information
about material properties and combined use of velocity and topological
structure in seismic analysis provides greater insight into the granular
packings. Here, we first report the results on the bulk stiffness of granular
mixtures with diverse rubber content obtained by the experimental tests.
Particular attention is devoted to the dependence of sound velocity on the
applied load and soft-stiff composition of samples, since this is an important
controllable experimental parameter [24, 31, 62]. After that, X-ray images
taken during wave propagation are analysed to obtain micro-insights of
packings.

###  4.1 Material characterization by ultrasonic measurements

(

(

a) b)

Figure 7: (a) In situ measured P-wave $M$ modulus against rubber fraction
$\nu$. (b) P-wave modulus against uniaxial applied load from $F=200$ N to 1000
N. The underlying measurement data for $\nu=0.5$ is provided in [44].

Figure 7 shows the compressional modulus determined for samples prepared with
rubber fractions from $\nu=$ 0.1 to 0.6 at five different load steps acquired
in situ. The underlying data of the mixture with $\nu=0.5$ together with the
recorded XRCT data set is for illustrative purposes available in [44]. $M$
modulus shows high stiffness for rubber content up to $\nu=0.3$ where the bulk
behavior of samples is highly controlled by the stiff phase. Previous
experimental studies under triaxial stress conditions showed already that a
small amount of soft particles enhances the effective stiffness of the medium
[62]. Thus, the effective stiffness of bidisperse granular mixtures consisting
of stiff and soft particles does not follow a simple mixture rule. The highest
modulus is observed at $\nu\approx$ 0.2. Thus, granular mixtures can be
manipulated to obtain aggregates with even higher stiffness, but lighter and
more dissipative thanks to rubber, when appropriate external conditions are
matched (in this case the pressure) [24]. Increasing the amount of rubber
particles, here 0.3 < $\nu$ < 0.6, reduces the effective stiffness where a
phase transition from stiff to soft phase occurs.

###  4.2 Image characterization

(

(

Figure 8: (a) 3d view of labeled particles for a sample prepared with
$\nu=0.5$. (b) The network of contacts of the sample in which different
particles colors represent a number of contact each particle carries.
Particles carrying zero, one, two, three, four, five, six, and above six
number of contacts are colored by black, magenta, blue, red, green, yellow,
orange, and white, respectively.

(

(

Figure 9: 3d view of (a) glass and (b) rubber particles network for a sample
prepared with rubber fraction $\nu=0.1$, respectively. Particles center are
marked with different colors based on their number of contacts and connected
with black lines. Particles carrying zero, one, two, three, four, five, six,
and above six number of contacts are colored by black, magenta, blue, red,
green, yellow, orange, and white, respectively.

(

(

Figure 10: 3d view of (a) glass and (b) rubber particles network for a sample
prepared with rubber fraction $\nu=0.5$, respectively. Particles center are
marked with different colors based on their number of contacts and connected
with black lines. Particles carrying zero, one, two, three, four, five, six,
and above six number of contacts are colored by black, magenta, blue, red,
green, yellow, orange, and white, respectively.

a) b) a) b) a) b)

Recent progresses in the application of $\mu$XRCT imaging in the field of
geomechanics have allowed all the individual particles of packings in a test
sample to be seen and identified uniquely in 3d. Combining such imaging
capabilities with experiments carried out “in situ” within an imaging setup
(shown in Figure 3) has led to the possibility of directly observing the
topology of packings.

Granular matter has a heterogeneous nature, this is expressed as force chains
through the medium. These force chains spread the forces through the contact
with another particle [60, 61]. This causes that neighboring particles could
have forces of different magnitudes and non-isotropic distributions of stress
throughout the medium [2]. This heterogeneity is manifest of the fact that
granular matter exhibits a strong configuration and history dependence. One of
the interesting outcomes of XRCT imaging of particulate systems is its
potentiality in depicting a particles contact network (e.g. force chains)
which reveals micro-structural information. A 3d view of labeled particles is
shown in Figure 10(a) for a sample prepared with $\nu=0.5$; and the contact
network of the sample is illustrated in Figure 10(b). Center of particles,
marked with different colors, are connected with black lines. The color code
of the number of contacts follows as black (zero), magenta (one), blue (two),
red (three), green (four), yellow (five), orange (six), and white (above six).

Next, to enhance micro-structural investigation, glass and rubber particles
are separated as a result of labelling every particles for samples with volume
fraction of $\nu=0.1$ (Figure 10) and 0.5 (Figure 10). Looking at glass
(Figure 10(a)) and rubber (Figure 10(b)) networks of $\nu=0.1$, it is not
surprising to see the scatter of rubbers among glass beads. These rubbers
which do not carry number of contacts above three are called rattlers.
Rattlers do not transfer force between each others as their number of contact
are not sufficient to form a chain of particles; whereas, moving to a sample
with higher rubber fraction $\nu=0.5$ (Figure 10), we do not see rattles
neither in the glass and nor in the rubber networks. This observation explains
why $M$ modulus remains fairly constant for samples with low rubber volume
fraction ($\nu\leq 0.3$) in Figure 7 since their bulk behavior is controlled
by stiff particles. But in case of intermediate regime, $\nu$ between 0.4 and
0.6, forces are distributed among glass and rubber chains which leads to a
smaller $M$ modulus.

##  5 Discussion

As is almost always the case, no experimental investigation is perfect. The
problems that arise must be critically questioned with regard to their
influence on the scientific question. The motivation for the presented study
was given by investigations of the same kinds of mixtures in a conventional
triaxial cell with well-defined boundary conditions [62]. Since it is quite
more complicated to design and build such a cell that can be used in an XRCT
device, we went for an oedometer cell. As the results in Figure 7 show, we can
observe the same phenomenon. As a consequence, it is an admissible
simplification for the problem under investigation.

One common significant uncertainty and difficulty is associated with
determining the exact wave propagation time by ultrasonic measurements
required to calculate the sample stiffness. Suggested criteria and
recommendations vary depending on the installation, application, and input
signal. The most common methodology is to interpret the received signal in the
time domain. It is typical to consider the first peak at the receiver
transducer as arrival time and the required time difference minus the needed
time in the other parts as the travel time of the signal within the packings
[41] which was consistently done within this study. However, in particular,
for high rubber packings, choosing a peak is not easy since the wave contains
high-frequency tones raised from the cluster of rubber particles. Thus, a low-
pass filter was applied to remove unwanted noises. For low rubber volume
fractions, this is nevertheless a difficult undertaking and not completely
objective. Although the same experimental protocol and measurement were
applied, a user error is inevitable. To account for this, the before performed
ex-situ measurements were repeated three times and are consistent with the in
situ acquired data shown in Figure 7.

To avoid movements of the particles during the image acquisition period taking
up to 2 h, it was switched from force-control to displacement-control. This
was done to avoid slight movements of the particles due to creeping processes
caused by the (viscoelastic) rubber particles. As a result, a corresponding
relaxation can be observed in the measured load, which, however, is less
critical with regard to the movement of individual beads. This can be
confirmed by the final quality of the 3d images as well as by comparing the in
situ ultrasonic measurement results with the before performed ex-situ
measurements. Consequently, the relatively long XRCT scan time has no
influence on the imaging. In the shown concept, the static preload $F_{stat.}$
is applied using a load frame, here a modified traditional universal testing
machine. Consequently, two rotary tables must be used to rotate the pre-
stressed sample. While this concept generally offers extremely great
flexibility (uniaxial tension/compression and torsion) compared to specially
designed load cells with integrated actuators, cf. [5, 54] and the literature
cited therein, it also presents some difficulties. The major problem is the
accurate alignment of the two rotatory tables in the micrometer range. Each
eccentricity error leads either to an undefined stress state and/or to bad
image quality due to potential movement errors in soft samples since the error
is multiplied with the set geometric magnification ($M_{geo.}$). However,
because an extremely low geometric magnification is applied here
($M_{geo.}=1.36$) and, in addition, the individual particles are comparatively
large with a diameter of 4 mm, this influence is to be classified as minor on
the underlying scientific question. Although it is slightly present in the
samples with a very low rubber fraction. To get completely rid of this error
source the alignment of the rotary tables should be further improved which is
technically challenging or alternatively, a load cell with integrated actuator
could be employed. Latter leads however to a significantly more complicated
and expensive cell design. That is the reason why this approach was
intentionally not used. Common artifacts in XRCT imaging (e.g. noise, ring
artifacts, system alignment errors, and beam hardening) could be reduced to a
minimum by the combined usage of hardware precautions and software-based
corrections. In general, it is advantageous to use relatively large particles
as in our approach, as this significantly reduces many potential sources of
error. For instance, alignment errors are less significant and here a $2\times
2$ detector binning was possible which significantly reduces the intrinsic
noise in the data sets. This simplifies the subsequent image processing due to
a significant better signal-to-noise ratio.

Segmentation of XRCT data is doubtlessly the most critical step before
quantitative analyses. Due to the good quality of the raw data sets
(justifying the scanning time of 2 h), and the fact that only three high-
contrast phases (air, glass, and rubber) had to be segmented, this could be
accomplished with a straight forward traditional workflow. The separation of
the single particles was performed using the well-known watershed algorithm.
Since the exact composition of each mixture is known (total particle volume
and number of particles), the segmentation and subsequent separation of the
single particles could be validated quite easily on a global level. The exact
deformation of the individual particles is not of interest for the presented
study. Consequently, the image processing and the quantification can be judged
as very reliable.

It can be concluded that the presented in situ setup in combination with the
applied workflow is adequate to get a more comprehensive understanding of
particulate systems as presented in this contribution. Further, it
demonstrates that an investigation of such kinds of questions are accessible
with a laboratory-based XRCT system if some boundary conditions are fulfilled,
e.g. an intrinsic large space within the system.

##  6 Summary

Understanding the response of granular-based systems in applications requires
a detailed grasp of the connection between the basic ingredients (particles)
and the macroscale properties of the systems considered. These are complex
systems and an understanding of the overall behavior cannot be gained by
studying individual particles. While significant progress has been made during
the last decades on understanding relevant physical mechanisms, there are
still many open questions, starting from the physics of particle interactions
to general features of multiscale models that will bridge the different
spatial and temporal scales of interest. The purpose of this research was to
describe an in situ approach for the combination of traditional experimental
techniques in combination with 3d imaging in order to explore a micro-macro
relation of granular assemblies.

The given contribution started initially with a short overview about
attenuation-based XRCT imaging provided in Section 2. Based on this, a
detailed explanation of how X-ray computed tomography and wave propagation
measuring technique can be combined was demonstrated in Section 3. For this, a
low X-ray absorbing oedometer cell with integrated P-wave transducers was
designed and built. The cell was integrated into a modular µXRCT-system which
was extended by a load frame using a refurbished universal testing machine.
This approach enables applying in a flexible manner different kinds of load
cases needed for advanced mechanical in situ investigations. A smart and
elegant approach to segment different material phases of the 3d tomograms was
explained step-by-step. The sample material chosen for the study is a
composition of monodisperse granular particles, glass, and rubber beads; not
only because of the ubiquity of granular materials but also because it is a
paradigm for complex disordered media and the unique characteristics exhibited
by various physical phenomena associated with mechanical waves in it
(dispersion, scattering attenuation, intrinsic attenuation, diffusion, weak
localization, energy transfer across different frequencies, etc.).

In Section 4, the effect of rubber volume content variation of the studied
glass-rubber mixtures on the compressive elastic modulus (P-wave modulus) by
means of wave propagation at different uniaxial preload levels was shown. It
was demonstrated that the modulus of the mixtures can be increased by adding
low portions of soft particles since it is dominated by stiff particles.
Adding more soft particles into assemblies led to a transition from stiff to
soft dominated regime where bulk behavior of samples is mainly controlled by
soft phase. The macro behavior of such granular samples can hardly be
described without micro information. To gain more insights, XRCT imaging was
performed in situ for two of the five examined preloads. Analyzing images
helped to understand the transition mechanism from stiff to soft controlled
regimes by separating glass and rubber networks. It was found that due to the
isolation of rubber particles, rattlers, the samples were dominated by glass
grains. Whereas, a sample with $50\,\%$ of rubber and glass has distributed
number of contacts, fairly, among themselves.

The combination of material and image characterization experimental techniques
is a complex procedure that was explained carefully in this article. However,
it could be shown what incredible potential it has. Especially since many
problems from the field of mechanics are accessible in laboratory-based XRCT
systems. Numerical simulations, based on the Discrete Element Method (DEM),
have revealed the outmost role of the microstructure in characterizing the
elastic behavior of granular soils. To further advance our understanding, one
could simulate particles networks obtained by image processing using DEM.

## Acknowledgments

The authors acknowledge the scientific discussion with Stefan Luding. Also, we
thank Ralf Plonus for his technical support in modification of devices. M.R.
and H.S. acknowledge funding from the German Science Foundation (DFG) through
Project No. STE 969/13-1. K.T. and H.S. acknowledge funding by the German
Science Foundation (DFG) through the project STE-969/16-1 within the SPP 1897
“Calm, Smooth and Smart”. H.S. thanks the DFG for supporting this work under
Grant No. SFB 1313 (Project No. 327154368).

### Conflict of interest

The authors declare no potential conflict of interests.

## Supporting information

For demonstrative purposes, the µXRCT images (reconstructed data sets,
projection data sets, and metadata) of the investigated glass-rubber mixture
with a rubber volume content $\nu=0.5$ as well as the corresponding ultrasonic
measurement data are openly available in the Data Repository of the University
of Stuttgart (DaRUS) at https://doi.org/10.18419/darus-2208, [44].

## References

  * Atkinson and Sallfors [1991] J. H. Atkinson and G. Sallfors.  Experimental determination of soil properties. General report to session 1.  In _Proc. 10th ECSMFE, Florence_ , volume 3, pages 915–956, 1991\.
  * Bassett et al. [2012] D. S. Bassett, E. T. Owens, K. E. Daniels, and M. A. Porter.  Influence of network topology on sound propagation in granular materials.  _Physical Review E_ , 86(4):041306, 2012.  doi:10.1103/PhysRevE.86.041306.
  * Beucher [1992] S. Beucher.  The watershed transformation applied to image segmentation.  _Scanning Microscopy_ , 1992(6):28, 1992.
  * Bosscher et al. [1997] P. J. Bosscher, T. B. Edil, and S. Kuraoka.  Design of highway embankments using tire chips.  _Journal of Geotechnical and Geoenvironmental Engineering_ , 123(4):295–304, 1997.  doi:/10.1061/(asce)1090-0241(1997)123:4(295).
  * Buffiere et al. [2010] J.-Y. Buffiere, E. Maire, J. Adrien, J.-P. Masse, and E. Boller.  In situ experiments with X ray tomography: an attractive tool for experimental mechanics.  _Experimental Mechanics_ , 50(3):289–305, Jan. 2010.  doi:10.1007/s11340-010-9333-7.
  * Burland [1989] J. B. Burland.  Ninth Laurits Bjerrum memorial lecture: “small is beautiful”—the stiffness of soils at small strains.  _Canadian Geotechnical Journal_ , 26(4):499–516, 1989.  doi:10.1139/t89-064.
  * Buzug [2008] T. M. Buzug.  _Computed Tomography_.  Springer Berlin Heidelberg, 2008.  doi:10.1007/978-3-540-39408-2.
  * Carmignato et al. [2018] S. Carmignato, W. Dewulf, and R. Leach, editors.  _Industrial X-Ray Computed Tomography_.  Springer International Publishing, 2018.  doi:10.1007/978-3-319-59573-3.
  * Chang et al. [1995] C. S. Chang, S. J. Chao, and Y. Chang.  Estimates of elastic moduli for granular material with anisotropic random packing structure.  _International Journal of Solids and Structures_ , 32(14):1989–2008, 1995.  doi:10.1016/0020-7683(94)00225-L.
  * Clayton [2011] C. R. I. Clayton.  Stiffness at small strain: research and practice.  _Géotechnique_ , 61(1):5–37, 2011.  doi:10.1680/geot.2011.61.1.5.
  * Cnudde and Boone [2013] V. Cnudde and M. N. Boone.  High-resolution X-ray computed tomography in geosciences: A review of the current technology and applications.  _Earth-Science Reviews_ , 123:1–17, Aug. 2013.  doi:10.1016/j.earscirev.2013.04.003.
  * Digabel and Lantuéjoul [1978] H. Digabel and C. Lantuéjoul.  Iterative algorithms.  In _Proc. 2nd European Symp. Quantitative Analysis of Microstructures in Material Science, Biology and Medicine_ , volume 19, page 8. Riederer Verlag, 1978.
  * Feldkamp et al. [1984] L. A. Feldkamp, L. C. Davis, and J. W. Kress.  Practical cone-beam algorithm.  _Journal of the Optical Society of America A_ , 1(6):612, Jun. 1984.  doi:10.1364/josaa.1.000612.
  * Garga and O’Shaughnessy [2000] V. K. Garga and V. O’Shaughnessy.  Tire-reinforced earthfill. part 1: Construction of a test fill, performance, and retaining wall design.  _Canadian Geotechnical Journal_ , 37(1):75–96, 2000.  doi:10.1139/t99-084.
  * Hasan et al. [2020] S. Hasan, V. Niasar, N. Karadimitriou, J. Godinho, N. Vo, S. An, A. Rabbani, and H. Steeb.  Direct characterization of solute transport in unsaturated porous media using fast X-ray synchrotron microtomography.  _Proceedings of the National Academy of Sciences_ , 117(38):23443–23449, Sep. 2020.  doi:10.1073/pnas.2011716117.
  * Heimdahl and Drescher [1999] T. C. Heimdahl and A. Drescher.  Elastic anisotropy of tire shreds.  _Journal of Geotechnical and Geoenvironmental Engineering_ , 125(5):383–389, 1999.  doi:10.1061/(asce)1090-0241(1999)125:5(383).
  * Hınıslıoğlu and Ağar [2004] S. Hınıslıoğlu and E. Ağar.  Use of waste high density polyethylene as bitumen modifier in asphalt concrete mix.  _Materials Letters_ , 58(3-4):267–271, 2004.  doi:10.1016/S0167-577X(03)00458-0.
  * Iassonov et al. [2009] P. Iassonov, T. Gebrenegus, and M. Tuller.  Segmentation of X-ray computed tomography images of porous materials: a crucial step for characterization and quantitative analysis of pore structures.  _Water Resources Research_ , 45(9), Sep. 2009.  doi:10.1029/2009wr008087.
  * Jia and Mills [2001] X. Jia and P. Mills.  Sound propagation in dense granular materials.  In _Powders and Grains 2001_ , pages 105–112. Balkema, 2001.
  * Jia et al. [1999] X. Jia, C. Caroli, and B. Velicky.  Ultrasound propagation in externally stressed granular media.  _Physical Review Letters_ , 82(9):1863, 1999.  doi:10.1103/PhysRevLett.82.1863.
  * Jovičić and Coop [1997] V. Jovičić and M. R. Coop.  Stiffness of coarse-grained soils at small strains.  _Géotechnique_ , 47(3):545–561, 1997.  doi:10.1680/geot.1997.47.3.545.
  * Kak and Slaney [1987] A. C. Kak and M. Slaney.  _Principles of Computerized Tomographic Imaging (Classics in Applied Mathematics)_.  Society for Industrial and Applied Mathematics, 1987.  ISBN 9780898714944.
  * Karimpouli et al. [2017] S. Karimpouli, P. Tahmasebi, H. L. Ramandi, P. Mostaghimi, and M. Saadatfar.  Stochastic modeling of coal fracture network by direct use of micro-computed tomography images.  _International Journal of Coal Geology_ , 179:153–163, Jun. 2017.  doi:10.1016/j.coal.2017.06.002.
  * Kim and Santamarina [2008] H.-K. Kim and J. C. Santamarina.  Sand-rubber mixtures (large rubber chips).  _Canadian Geotechnical Journal_ , 45(10):1457–1466, 2008.  doi:10.1139/T08-070.
  * Kodym and Španěl [2018] O. Kodym and M. Španěl.  Semi-automatic CT image segmentation using random forests learned from partial annotations.  In S. Wiebe, H. Gamboa, A. Fred, and S. B. i Badia, editors, _Proceedings of the 11th International Joint Conference on Biomedical Engineering Systems and Technologies - BIOIMAGING_ , pages 124–131. INSTICC, SciTePress, 2018.  doi:10.5220/0006588801240131.
  * Kotwaliwale et al. [2011] N. Kotwaliwale, K. Singh, A. Kalne, S. N. Jha, N. Seth, and A. Kar.  X-ray imaging methods for internal quality evaluation of agricultural produce.  _Journal of Food Science and Technology_ , 51(1):1–15, Aug. 2011.  doi:10.1007/s13197-011-0485-y.
  * Kruyt et al. [2010] N. P. Kruyt, I. Agnolin, S. Luding, and L. Rothenburg.  Micromechanical study of elastic moduli of loose granular materials.  _Journal of the Mechanics and Physics of Solids_ , 58(9):1286–1301, 2010.  doi:10.1016/j.jmps.2010.06.003.
  * Landis and Keane [2010] E. N. Landis and D. T. Keane.  X-ray microtomography.  _Materials Characterization_ , 61(12):1305–1316, Dec. 2010.  doi:10.1016/j.matchar.2010.09.012.
  * Lee et al. [1999] J. H. Lee, R. Salgado, A. Bernal, and C. W. Lovell.  Shredded tires and rubber-sand as lightweight backfill.  _Journal of Geotechnical and Geoenvironmental Engineering_ , 125(2):132–141, 1999.  doi:10.1061/(asce)1090-0241(1999)125:2(132).
  * Lee and Santamarina [2005] J.-S. Lee and J. C. Santamarina.  Bender elements: performance and signal interpretation.  _Journal of Geotechnical and Geoenvironmental Engineering_ , 131(9):1063–1070, 2005.  doi:10.1061/(asce)1090-0241(2005)131:9(1063).
  * Lee et al. [2007] J.-S. Lee, J. Dodds, and J. C. Santamarina.  Behavior of rigid-soft particle mixtures.  _Journal of Materials in Civil Engineering_ , 19(2):179–184, 2007.  doi:10.1061/(asce)0899-1561(2007)19:2(179).
  * Magnanimo and et al [2008] V. Magnanimo and et al.  Characterizing the shear and bulk moduli of an idealized granular material.  _EPL (Europhysics Letters)_ , 81(3):34006, 2008\.  doi:10.1209/0295-5075/81/34006.
  * Maire and Withers [2013] E. Maire and P. J. Withers.  Quantitative X-ray tomography.  _International Materials Reviews_ , 59(1):1–43, Dec. 2013.  doi:10.1179/1743280413y.0000000023.
  * Makse et al. [2004] H. A. Makse, N. Gland, D. L. Johnson, and L. Schwartz.  Granular packings: Nonlinear elasticity, sound propagation, and collective relaxation dynamics.  _Physical Review E_ , 70(6):061302, 2004.  doi:10.1103/PhysRevE.70.061302.
  * Misra and Nejadsadeghi [2019] A. Misra and N. Nejadsadeghi.  Longitudinal and transverse elastic waves in 1d granular materials modeled as micromorphic continua.  _Wave Motion_ , 90:175–195, 2019.  doi:10.1016/j.wavemoti.2019.05.005.
  * Mizutani and Suzuki [2012] R. Mizutani and Y. Suzuki.  X-ray microtomography in biology.  _Micron_ , 43(2-3):104–115, Feb. 2012.  doi:10.1016/j.micron.2011.10.002.
  * Mouraille and Luding [2008] O. Mouraille and S. Luding.  Sound wave propagation in weakly polydisperse granular materials.  _Ultrasonics_ , 48(6-7):498–505, 2008.  doi:10.1016/j.ultras.2008.03.009.
  * Nesterenko [1984] V. F. Nesterenko.  Propagation of nonlinear compression pulses in granular media.  _Journal of Applied Mechanics and Technical Physics_ , 24(5), 1984.  doi:10.1007/BF00905892.
  * Object Research Systems Inc., Montreal, Canada() [ORS] Object Research Systems (ORS) Inc., Montreal, Canada.  Dragonfly 3.1 (computer software).  Available online: http://www.theobjects.com/dragonfly.
  * Otsu [1979] N. Otsu.  A threshold selection method from gray-level histograms.  _IEEE Transactions on Systems, Man, and Cybernetics_ , 9(1):62–66, Jan. 1979.  doi:10.1109/tsmc.1979.4310076.
  * O’Donovan et al. [2015] J. O’Donovan, C. O’sullivan, G. Marketos, and D. M. Wood.  Analysis of bender element test interpretation using the discrete element method.  _Granular Matter_ , 17(2):197–216, 2015.  doi:10.1007/s10035-015-0552-6.
  * Rahimabadi et al. [2020] P. S. Rahimabadi, M. Khodaei, and K. R. Koswattage.  Review on applications of synchrotron-based X-ray techniques in materials characterization.  _X-Ray Spectrom._ , 49(3):348–373, Mar. 2020\.  doi:10.1002/xrs.3141.
  * Ruf and Steeb [2020] M. Ruf and H. Steeb.  An open, modular, and flexible micro X-ray computed tomography system for research.  _Review of Scientific Instruments_ , 91(11):113102, 2020.  doi:10.1063/5.0019541.
  * Ruf et al. [2021] M. Ruf, K. Taghizadeh, and H. Steeb.  micro-XRCT data sets and in situ measured ultrasonic wave propagation of a pre-stressed monodisperse rubber and glass particle mixture with 50 % volume rubber content, 2021.  DaRUS. doi:10.18419/DARUS-2208.
  * Russ and Neal [2016] J. C. Russ and F. B. Neal.  _The Image Processing Handbook_.  CRC Press, 7 edition, 2016.
  * Russo [2018] P. Russo, editor.  _Handbook of X-ray Imaging: Physics and Technology_.  CRC Press, 2018.  ISBN 978-1-4987-4152-1.
  * Santamarina and Aloufi [1999] J. C. Santamarina and M. Aloufi.  Small strain stiffness: a micromechanical experimental study.  _Proceedings of Pre-failure Deformation Characteristics of Geomaterials_ , pages 451–458, 1999.
  * Sawangsuriya [2012] A. Sawangsuriya.  Wave propagation methods for determining stiffness of geomaterials.  In _Wave Processes in Classical and New Solids_. IntechOpen, 2012\.  doi:10.5772/48562.
  * Schlüter et al. [2014] S. Schlüter, A. Sheppard, K. Brown, and D. Wildenschild.  Image processing of multiphase images obtained via X-ray microtomography: a review.  _Water Resources Research_ , 50(4):3615–3639, Apr. 2014.  doi:10.1002/2014wr015256.
  * Schoeman et al. [2016] L. Schoeman, P. Williams, A. du Plessis, and M. Manley.  X-ray micro-computed tomography (µCT) for non-destructive characterisation of food microstructure.  _Trends in Food Science & Technology_, 47:10–24, Jan. 2016.  doi:10.1016/j.tifs.2015.10.016.
  * Serra [1986] J. Serra.  Introduction to mathematical morphology.  _Computer Vision, Graphics, and Image Processing_ , 35(3):283–305, 1986.  doi:10.1016/0734-189X(86)90002-2.
  * Shrivastava and Luding [2017] R. K. Shrivastava and S. Luding.  Effect of disorder on bulk sound wave speed: a multiscale spectral analysis.  _Nonlinear Processes in Geophysics_ , 24(3):435–454, 2017.  doi:0.5194/npg-24-435-2017.
  * Siddique and Naik [2004] R. Siddique and T. R. Naik.  Properties of concrete containing scrap-tire rubber–an overview.  _Waste Management_ , 24(6):563–569, 2004.  doi:10.1016/j.wasman.2004.01.006,.
  * Singh et al. [2014] S. S. Singh, J. J. Williams, P. Hruby, X. Xiao, F. D. Carlo, and N. Chawla.  In situ experimental techniques to study the mechanical behavior of materials using X-ray synchrotron tomography.  _Integrating Materials and Manufacturing Innovation_ , 3(1):109–122, Apr. 2014.  doi:10.1186/2193-9772-3-9.
  * Somfai et al. [2005] E. Somfai, J.-N. Roux, J. H. Snoeijer, M. Van Hecke, and W. Van Saarloos.  Elastic wave propagation in confined granular systems.  _Physical Review E_ , 72(2):021301, 2005.  doi:10.1103/PhysRevE.72.021301.
  * Stock [2008a] S. R. Stock.  _MicroComputed Tomography: Methodology and Applications_.  CRC Press, 2008a.  ISBN 9781420058765.
  * Stock [2008b] S. R. Stock.  Recent advances in X-ray microtomography applied to materials.  _International Materials Reviews_ , 53(3):129–181, May. 2008b.  doi:10.1179/174328008x277803.
  * Taghizadeh et al. [2015] K. Taghizadeh, N. Kumar, V. Magnanimo, and S. Luding.  Understanding the effects of inter-particle contact friction on the elastic moduli of granular materials.  In _IOP Conference Series: Earth and Environmental Science_ , volume 26, page 012008. IOP Publishing, 2015.  doi:10.1088/1755-1315/26/1/012008.
  * Taghizadeh et al. [2017] K. Taghizadeh, S. Luding, and V. Magnanimo.  DEM applied to soil mechanics.  _ALERT Doctoral School 2017 Discrete Element Modeling_ , page 129, 2017.
  * Taghizadeh et al. [2021a] K. Taghizadeh, R. K. Shrivastava, and S. Luding.  Stochastic model for energy propagation in disordered granular chains.  _Materials_ , 14(7):1815, 2021a.  doi:10.3390/ma14071815.
  * Taghizadeh et al. [2021b] K. Taghizadeh, H. Steeb, and S. Luding.  Energy propagation in 1d granular soft-stiff chain.  In _EPJ Web of Conferences_ , volume 249, page 02002. EDP Sciences, 2021b.  doi:10.1051/epjconf/202124902002.
  * Taghizadeh et al. [2021c] K. Taghizadeh, H. Steeb, S. Luding, and V. Magnanimo.  Elastic waves in particulate glass-rubber mixtures.  _Proceedings of the Royal Society A_ , 477(2249):20200834, 2021c.  doi:10.1098/rspa.2020.0834.
  * Tuller et al. [2015] M. Tuller, R. Kulkarni, and W. Fink.  Segmentation of X-ray CT data of porous materials: a review of global and locally adaptive algorithms.  pages 157–182. American Society of Agronomy and Soil Science Society of America, Oct. 2015.  doi:10.2136/sssaspecpub61.c8.
  * Van den Wildenberg et al. [2016] S. Van den Wildenberg, A. Tourin, and X. Jia.  Sound velocity fluctuations in confined granular materials: Coarse-graining lengths and elastic heterogeneities.  _EPL (Europhysics Letters)_ , 115(3):34005, 2016\.  doi:10.1209/0295-5075/115/34005.
  * Vincent and Soille [1991] L. Vincent and P. Soille.  Watersheds in digital spaces: an efficient algorithm based on immersion simulations.  _IEEE Transactions on Pattern Analysis & Machine Intelligence_, 13(06):583–598, 1991.  doi:10.1109/34.87344.
  * Vlassenbroeck et al. [2007] J. Vlassenbroeck, M. Dierick, B. Masschaele, V. Cnudde, L. Van Hoorebeke, and P. Jacobs.  Software tools for quantification of X-ray microtomography at the UGCT.  _Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment_ , 580(1):442–445, Sep. 2007.  doi:10.1016/j.nima.2007.05.073.
  * Walton [1987] K. Walton.  The effective elastic moduli of a random packing of spheres.  _Journal of the Mechanics and Physics of Solids_ , 35(2):213–226, 1987.  doi:10.1016/0022-5096(87)90036-6.
  * Xiao et al. [2007] F. Xiao, S. Amirkhanian, and C. H. Juang.  Rutting resistance of rubberized asphalt concrete pavements containing reclaimed asphalt pavement mixtures.  _Journal of Materials in Civil Engineering_ , 19(6):475–483, 2007.  doi:10.1061/(asce)0899-1561(2007)19:6(475).
  * Zornberg et al. [2004] J. G. Zornberg, A. R. Cabral, and C. Viratjandr.  Behaviour of tire shred sand mixtures.  _Canadian Geotechnical Journal_ , 41(2):227–241, 2004.  doi:10.1139/t03-086.