capbtabboxtable[][]

# High angular sensitivity X-ray phase-contrast microtomography of soft-tissue
through a two-directional beam-tracking synchrotron set-up

Carlos Navarrete-León  Department of Medical Physics and Biomedical
Engineering. University College London, Malet Place, Gower Street, London,
WC1E 6BT, United Kingdom  P. Stephen Patrick  Centre for Advanced Biomedical
Imaging, Division of Medicine, University College London, Paul O’Gorman
Building, 72 Huntley Street, WC1E 6DD, United Kingdom  Adam Doherty
Department of Medical Physics and Biomedical Engineering. University College
London, Malet Place, Gower Street, London, WC1E 6BT, United Kingdom  Harry
Allan  Department of Medical Physics and Biomedical Engineering. University
College London, Malet Place, Gower Street, London, WC1E 6BT, United Kingdom
Silvia Cipiccia  Department of Medical Physics and Biomedical Engineering.
University College London, Malet Place, Gower Street, London, WC1E 6BT, United
Kingdom  Shashidhara Marathe  Diamond Light Source, Harwell Science and
Innovation Campus, Fermi Avenue, Didcot, OX11 0DE, United Kingdom  Kaz Wanelik
Diamond Light Source, Harwell Science and Innovation Campus, Fermi Avenue,
Didcot, OX11 0DE, United Kingdom  Michela Esposito  Department of Medical
Physics and Biomedical Engineering. University College London, Malet Place,
Gower Street, London, WC1E 6BT, United Kingdom  Charlotte K. Hagen  Department
of Medical Physics and Biomedical Engineering. University College London,
Malet Place, Gower Street, London, WC1E 6BT, United Kingdom  Alessandro Olivo
Department of Medical Physics and Biomedical Engineering. University College
London, Malet Place, Gower Street, London, WC1E 6BT, United Kingdom  Marco
Endrizzi  Department of Medical Physics and Biomedical Engineering. University
College London, Malet Place, Gower Street, London, WC1E 6BT, United Kingdom

###### Abstract

Two-directional beam-tracking (2DBT) is a method for phase-contrast imaging
and tomography that uses an intensity modulator to structure the X-ray beam
into an array of independent circular beamlets that are resolved by a high-
resolution detector. It features isotropic spatial resolution, provides two-
dimensional phase sensitivity, and enables the three-dimensional
reconstructions of the refractive index decrement, $\delta$, and the
attenuation coefficient, $\mu$. In this work, we report on the angular
sensitivity and the spatial resolution of 2DBT images in a synchrotron-based
implementation. In its best configuration, we obtained angular sensitivities
of $\sim$20 $n$rad and spatial resolution of at least 6.25 $\mu$m in phase-
contrast images. We also demonstrate exemplar application to the three-
dimensional imaging of soft tissue samples, including a mouse liver and a
decellularised porcine dermis.

## Introduction

X-ray phase-contrast tomography (XPCT) is a non-destructive imaging technique
that enables the 3D visualization of materials and tissues composed of low-Z
elements. Through generating contrast also from the phase shift induced in the
x-ray wavefront, this technique allows for the visualization of details
otherwise undetectable by using a conventional, attenuation-based approach to
generating image contrast [1]. Notable examples have been demonstrated in the
biomedical field using synchrotron radiation, including brain [2, 3], lung [4,
5], kidney [6, 7], breast [8, 9, 10], and oesophagus [11] imaging, amongst
others.

Several methods have been developed at synchrotron radiation facilities for
XPCT including crystal-based interferometric methods, propagation-based
imaging methods, analyzer-based imaging methods, grating-based interferometric
methods, speckle-based imaging methods, and non-interferometric mask-based
methods [12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22].

Here we focus on two-directional beam-tracking (2DBT), which belongs to the
category of non-interferometric mask-based methods. The method was first
proposed by a patent in the mid-90s [23] and it shares similarities with the
Shack-Hartmann wave-front sensor. It requires the use of a single optical
element (modulator), which structures the beam into an array of independent
beamlets that are then resolved by a high-resolution detector. The modulator
is placed before the sample, contrary to other implementations [24, 25],
meaning that all photons reaching the sample contribute to image formation,
limiting the absorbed dose.

In this approach, the radiation field intensity modulation, coupled with a
dedicated data analysis methodology, provides a way to measure how the sample
affects the intensity and position of each X-ray beamlet. A shift in their
position is interpreted as a refraction effect, which is related to the first
derivative of the phase shift imposed by the sample. A change in intensity is
interpreted as attenuation of the X-ray beam.

This approach provides aperture-driven [26, 27, 28] and isotropic [29] spatial
resolution; and with dedicated mask designs alongside efficient acquisition
schemes, scanning time can be improved for optimal acquisition [30, 31], and
fly-scan data acquisition schemes. We note that the approach allows also for
X-ray dark-field imaging [32]. The method was recently demonstrated in a
compact laboratory set-up [33], within a small ($<$1 m) footprint and by using
a low power (10 W) source. Here we report on an implementation that made use
of synchrotron radiation, and that we found suitable for high-sensitivity
measurement of phase gradients, providing excellent contrast for the
visualisation of the morphology in soft tissue samples. We report on the
angular sensitivity of this approach, characterised as a function of exposure
time and system geometry, as well as its spatial resolution, estimated with
Fourier ring correlation on tomographic reconstructions. We demonstrate
exemplary application to the three-dimensional imaging of soft tissue samples,
including both a formalin-fixed sample of mouse liver, as well as a
decellularised, iodine-stained porcine dermis. For both biological samples,
phase-contrast imaging enhanced the visibility of key physiological features
above attenuation-based images acquired with equivalent x-ray exposure.

Figure 1: Schematic diagram of the two-directional beam-tracking experimental
set-up.

## Materials and methods

### 2DBT x-ray set-up

The XPCT set-up is presented in Figure 1. The experiments were carried out at
Diamond Light Source Beamline I13-2. The angular sensitivity measurements were
done with a mean energy of 16 keV from a filtered pink beam with a Silicon
mirror and filters of 1.34 mm Pyrolytic graphite, 1.4 mm Aluminium, and 0.042
mm Niobium. For the biological specimens, the mean energy was increased to 27
keV to reduce sample damage by changing to a Platinum mirror and filters of
1.34 mm Pyrolytic graphite and 3.2 mm Aluminium.

The sample was placed roughly 221 m from the source, and the modulator was
placed 15 cm upstream of the sample. The modulator is fabricated with laser-
ablation from a 100 $\mu$m thick tungsten foil (Goodfellow), and has a period
of 50 $\mu$m. The apertures have a conical shape with diameters of 15 $\mu$m
in the front and 30 $\mu$m in the back. The detector is a pco.edge 5.5 camera
coupled to a scintillator-objective combination with an effective pixel size
of $2.6\times 2.6$ $\mu$m2 and a field of view of $6.6\times 5.6$ cm2.

### Angular sensitivity measurements

The angular sensitivity was assessed with a custom-built phantom composed of
soda-lime glass microspheres of 50 $\mu$m diameter (Fischer Scientific,
monodisperse) embedded in wax and polyethylene foam. To study the sensitivity
as a function of the object-to-detector distance ($z_{od}$), the sample was
imaged at the following distances: $z_{od}$={2.5, 7.5, 17.5, 37.5, 77.5} cm,
by moving the detector. The sample was moved in a $10\times 10$ grid with
steps of 5 $\mu m$ in the $xy$ plane. At each sample position, $10\times 0.1$
s frames were acquired, and flat and dark images were taken before and after
each scan. The frames were used to assess sensitivity as a function of the
exposure time. The assessment of the angular sensitivity during imaging was
carried out by calculating the mean and standard error of the standard
deviation of the measured refraction angles in an area without the sample, for
which eight different windows of $5\times 40$ pixels were used.

### Tomography of unstained and stained ex-vivo tissues

Two biological samples were imaged: a mouse liver and a hernia mesh,
consisting of decellularised porcine dermis. The liver was fixed in 4% para-
formaldehyde for 24 hours upon dissection from a 2 month old C57BL/6 mouse
(Charles River Laboratories), then stored in 0.9% saline. The decellularised
porcine dermis (XenmatrixTM, Bard), was stained in 3% Lugol’s iodine solution
(Scientific Laboratory Supplies) in phosphate buffered saline (Gibco) for 24
hours before storage in 0.9% saline. Both liver and decellularised dermis were
prepared for imaging by embedding in 1% agar (Thermo Scientific Chemicals).
The samples were between 3.1 and 3.8 mm wide and they were scanned by
acquiring 1200 projections while rotating over 180° in a fly-scan fashion with
an exposure time of 0.15 s per projection. This was repeated at different
modulator sub-pitch displacements to increase sampling. The modulator was
raster-scanned in $8\times 8$ positions, by using 6.25 $\mu$m displacements
both in $x$ and $y$. This led to a total exposure time of $1200\times 8\times
8\times 0.15$ s = 3.2 h for each sample. Flat and dark images were acquired at
each modulator position, before and after rotating the sample. Note that this
is different from radiography, in which it was the sample that was moved. The
detector was placed 128 cm away from the sample to further increase the
angular sensitivity, which was also assessed for this configuration with eight
different windows of $8\times 8$ pixels.

### Data Analysis

The transmission, refraction in $x$, and refraction in $y$ images were
obtained by selecting a window of $20\times 20$ pixels around each beamlet and
comparing the intensities with ($I_{s}(x,y)$) and without ($I_{0}(x,y)$) the
sample in the beamlet. The transmission was calculated by dividing the sum of
the intensities in the windows: $t=\sum_{xy}I_{s}(x,y)/\sum_{xy}I_{0}(x,y)$,
and the two refraction images by measuring the displacements ($\Delta_{x}$,
$\Delta_{y}$) between the beamlets with a subpixel cross-correlation algorithm
[34]. This was performed for all images acquired at each sample or modulator
position, which were then stitched to obtain an image with higher sampling
[33].

Assuming small refraction angles and under a geometrical optics approximation,
the refraction angle, $\alpha_{xy}$, is related to the displacements
$\Delta_{x}$ and $\Delta_{y}$ and the orthogonal gradients of the phase shift,
$\Delta\Phi_{x,y}$, by:

$\alpha_{x,y}=\frac{\Delta_{x,y}}{z_{od}}=\frac{\Delta\Phi_{x,y}}{k},$ (1)

where $z_{od}$ is the object-to-detector distance and $k$ is the wavenumber.
This allows us to obtain the phase shift $\Delta\Phi$ by integrating the two
gradients through a Fourier space method [35]).

The retrieved quantities $t$ and $\Delta\Phi$ are linked to integrals along
the photon path of the linear attenuation coefficient ($\mu$) and the real
part of the refractive index ($\delta$) in the following way:

$-\ln t(x,y)=\int_{o}\mu(x^{\prime},y^{\prime},z^{\prime})dz$ (2)

$-\frac{\Delta\Phi(x,y)}{k}=\int_{o}\delta(x^{\prime},y^{\prime},z^{\prime})dz.$
(3)

Therefore, for the biological specimens, volumes of $\mu$ and $\delta$ were
obtained from the projections taken at different viewing angles using the
filtered back projection (FBP) implementation of the Astra toolbox [36].

For both the $\mu$ and $\delta$ volumes, the spatial resolution of slices in
the three orthogonal planes was estimated using an implementation of Fourier
ring correlation (FRC) [37] provided as part of the BIOP ImageJ plugin [38].
Independent inputs were provided to the algorithm by reconstruction of two
volumes, each using half of the available projections. To reduce the noise of
the FRC estimate, the resultant curves of 5 adjacent representative slices
from the middle of the volume were averaged. Resolutions are stated using the
3-$\sigma$ criterion, expressing the spatial frequency at which the FRC curve
exceeds by 3 standard deviations the expected correlations within the random
background noise [39].

Figure 2: Angular sensitivity of the method as a function of (a) object-to-
detector distance and (b) exposure time for different system configurations.
(c) Phase images of the phantom (polyethylene foam, and microspheres and air
bubbles embedded in wax) are shown for increasing object-to-detector
distances. The improvement in angular sensitivity reveals interfaces in the
foam and small bubbles in the wax substrate, as pointed out by the arrows. d)
An inset in the phase image is shown, along with the two refraction images,
for further demonstration of thin wax deposits on the substrate being unveiled
with increasing sensitivity.

## Results and discussion

### Angular sensitivity

The results from the angular sensitivity measurements are reported in Figure
2. The sensitivity is shown for different object-to-detector distances (Fig.
2a) and for an increasing number of integrated frames for both energy
configurations (Fig 2b). We observe that the smallest resolvable angle
decreases proportionally with increasing propagation distance ($\propto
1/z_{od}$) within the range of propagation distances explored, as expected
from the geometrical optics approximation used in Equation 1. We note a small
deviation from this trend at 77.5 cm of propagation distance and we observe it
is associated with a small decrease in visibility, from 85% to 82%. Longer
propagation distances offer a further increase in the angular sensitivity,
however, the angular sensitivity is expected to increase at a relatively
slower rate beyond this point. In addition, we also observed that when
integrating a few frames, the sensitivity goes inversely with the square root
of the number of counts ($\sqrt{N}$). This is the expected trend from Poisson
statistics, indicating that this is the dominant noise source up to
accumulations of 400 ms; beyond this point, additional noise sources become
significant in limiting the smallest measurable refraction angle. The smallest
angle with a mean energy of 16 keV was measured with a combination of
$z_{od}=77.5$ cm and $10\times 0.1$s frames, for which we report an angular
sensitivity of $21.6\pm 0.2$ $n$rad. In the conditions for tomographic imaging
at 27 keV, the angular sensitivity benefits from the increased propagation
distance of $z_{od}=128$ cm, and we measured $35\pm 2$ $n$rad with only
$1\times 0.15$s frame.

The effect of increasing angular sensitivity on image quality can be observed
in Figures 2c,d, where the integrated phase images are presented along with
both refraction images of the smaller, highlighted region of interest. The
increasing angular sensitivity unveils various interfaces in the foam and
small bubbles in the wax substrate, as pointed out by the arrows in Fig. 2c.
This is further demonstrated with the insets in Fig. 2d, where the two
refraction images show increasingly lower noise levels as propagation distance
is increased, which in this case reveals thinner deposits of wax on the
substrate as angular sensitivity increases.

Figure 3: X-ray phase-contrast and attenuation tomography of biological soft-
tissues shown as 3D-rendered volumes and in representative perpendicular
planes. (a, b) phase-contrast, and (d) attenuation tomography of mouse liver.
(c) Line profile of signal through indicated liver cross sections of b and d
showing improved contrast to noise ratio in the phase contrast. (e) Fourier
ring correlation curve calculated from two independent reconstructions of (b)
and (d) showing an increased resolution for phase-contrast. (f) 3D render of
phase contrast tomography of decellularised porcine dermis, showing
perpendicular cross-sections of (g) phase-contrast and (h) attenuation.

### Tomography of biological soft tissues

The attenuation and phase contrast tomographic reconstructions of the liver
tissue and decellularised dermis samples are presented in Figure 3. The
significant increase in contrast-to-noise ratio achieved by means of phase
contrast is evident across all slices and for both stained and unstained
samples. While the image quality enhancement from phase contrast was
particularly evident in the unstained (i.e. poorly absorbing) liver tissue,
even tissue optimised for attenuation-based imaging using an iodine stain
still showed notable improvement with phase contrast. In terms of
physiological features identifiable in the liver sample, phase-contrast gave
clear visualisation of the hepatic portal vein, arteries, and bile ducts (Fig
3a,b), while these were not discernible on the attenuation-based images (Fig.
3c). For the stained sample of the decellularised dermis, while the lumen of
the hair follicles was discernible in both phase (Fig. 3f,g) and attenuation-
based images (Fig. 3h), phase-contrast improved the clarity of the structure
of the hair follicles including the layers of the root sheath, sebaceous
gland, and fibre alignment around the dermic sheath (Fig 3f-h).

Figure 3e displays the FRC curves obtained from the unstained liver tissue
axial slices, as illustrated in Figure 3b and Figure 3d. The attenuation curve
intersects the threshold at 0.4 px-1, indicating a spatial resolution of 16
$\mu$m. Meanwhile, the phase-contrast curve fails to intersect with the
threshold, suggesting that in this case the spatial resolution is sampling
limited and is at least equal to the sampling pixel size of 6.25 $\mu$m. We
note that the FRC curves in the orthogonal cross-sections showed comparable
trends. We interpret the disparity between the FRC curves obtained through the
attenuation- and phase-contrast tomography as a consequence of the noise and
contrast dependence inherent in the FRC resolution metric. High-spatial
frequency image features are unable to surpass the noise threshold in the
noisier, low-contrast attenuation volume, whereas the much higher signal-to-
noise ratio achieved with phase contrast allows for separating even the
smallest features from the background. We also note that by splitting the
projection dataset to obtain independent volumes, the angular tomographic
sampling has also been halved to 600 projections. Although the full dataset
was largely oversampled in terms of viewing angles, the halved dataset is
undersampled with respect to the Nyquist sampling theorem for samples between
500 and 600 pixels of width. As such, the FRC result should still be
considered a conservative estimate of the achievable resolution.

## Conclusion

We have here studied the angular sensitivity and spatial resolution in images
obtained through a 2DBT synchrotron set-up and shown the potential of the
method for volumetric imaging of soft tissues with poor attenuation contrast.
We report angular sensitivities of $\sim$20 $n$rad at 1s exposure time, 77.5
cm of propagation distance, and 16 keV mean energy; and $\sim$35 $n$rad at 150
ms, 128cm, and 27 keV mean energy. Our results indicated that the geometrical-
optics approximation used by the phase-retrieval algorithms is well satisfied
within 1 m of propagation distance. We have also shown sub-aperture spatial
resolution in phase-contrast tomography, which was observed to be limited by
sampling to a factor 2.4x better than the apertures in the modulator. We note
that this indicates a spatial resolution improvement with respect to what was
previously modelled [26]. These results provide a basis for future
experimental designs with the 2DBT method, especially for identifying optimal
trade-offs between angular sensitivity, spatial sampling, and acquisition
time. They also provide a basis for comparison with similar imaging methods.

## Acknowledgements

We gratefully acknowledge Diamond Light Source for time on Beamline I-13 under
Proposal MG30748. This work was supported by the EPSRC EP/T005408/1; Wellcome
Trust 221367/Z/20/Z; the National Research Facility for Lab x-ray CT (NXCT)
through EPSRC grants EP/T02593X/1 and EP/V035932/1. PSP acknowledges funding
by MRC grant MR/R026416/1, EPSRC Knowledge Exchange and Innovation Fund award
EP/R511638/1, and Wellcome Trust/UCL Devices & Diagnostics Therapeutic
Innovation Network Pilot Data fund (562448 – Linked to Lead 553191). CKH was
supported by the Royal Academy of Engineering under the Research Fellowship
scheme. AO is supported by the Royal Academy of Engineering under their Chair
in Emerging Technologies scheme (CiET1819/2/78).

## References

  * [1] D. M. Paganin and D. Pelliccia, “Chapter two - x-ray phase-contrast imaging: a broad overview of some fundamentals,” in Advances in Imaging and Electron Physics (M. Hÿtch and P. W. Hawkes, eds.), vol. 218 of Advances in Imaging and Electron Physics, pp. 63–158, Elsevier, 2021.
  * [2] B. Pinzer, M. Cacquevel, P. Modregger, S. McDonald, J. Bensadoun, T. Thuering, P. Aebischer, and M. Stampanoni, “Imaging brain amyloid deposition using grating-based differential phase contrast tomography,” NeuroImage, vol. 61, no. 4, pp. 1336–1346, 2012.
  * [3] M. Töpperwien, F. van der Meer, C. Stadelmann, and T. Salditt, “Correlative x-ray phase-contrast tomography and histology of human brain tissue affected by alzheimer’s disease,” NeuroImage, vol. 210, p. 116523, 2020.
  * [4] R. P. Xian, C. L. Walsh, S. E. Verleden, W. L. Wagner, A. Bellier, S. Marussi, M. Ackermann, D. D. Jonigk, J. Jacob, P. D. Lee, and P. Tafforeau, “A multiscale x-ray phase-contrast tomography dataset of a whole human left lung,” Scientific Data, vol. 9, June 2022.
  * [5] J. Reichmann, S. E. Verleden, M. Kühnel, J. C. Kamp, C. Werlein, L. Neubert, J.-H. Müller, T. Q. Bui, M. Ackermann, D. Jonigk, and T. Salditt, “Human lung virtual histology by multi-scale x-ray phase-contrast computed tomography,” Physics in Medicine and Biology, vol. 68, p. 115014, may 2023.
  * [6] M.-C. Zdora, P. Thibault, W. Kuo, V. Fernandez, H. Deyhle, J. Vila-Comamala, M. P. Olbinado, A. Rack, P. M. Lackie, O. L. Katsamenis, M. J. Lawson, V. Kurtcuoglu, C. Rau, F. Pfeiffer, and I. Zanette, “X-ray phase tomography with near-field speckles for three-dimensional virtual histology,” Optica, vol. 7, pp. 1221–1227, Sep 2020.
  * [7] C. L. Walsh, P. Tafforeau, W. L. Wagner, D. J. Jafree, A. Bellier, C. Werlein, M. P. Kühnel, E. Boller, S. Walker-Samuel, J. L. Robertus, D. A. Long, J. Jacob, S. Marussi, E. Brown, N. Holroyd, D. D. Jonigk, M. Ackermann, and P. D. Lee, “Imaging intact human organs with local resolution of cellular structures using hierarchical phase-contrast tomography,” Nature Methods, vol. 18, pp. 1532–1541, Nov. 2021.
  * [8] Y. Zhao, E. Brun, P. Coan, Z. Huang, A. Sztrókay, P. C. Diemoz, S. Liebhardt, A. Mittone, S. Gasilov, J. Miao, and A. Bravin, “High-resolution, low-dose phase contrast x-ray tomography for 3d diagnosis of human breast cancers,” Proceedings of the National Academy of Sciences, vol. 109, pp. 18290–18294, Sept. 2012.
  * [9] P. Baran, S. Mayo, M. McCormack, S. Pacilè, G. Tromba, C. Dullin, F. Zanconati, F. Arfelli, D. Dreossi, J. Fox, Z. Prodanovic, M. Cholewa, H. Quiney, M. Dimmock, Y. Nesterets, D. Thompson, P. Brennan, and T. Gureyev, “High-resolution x-ray phase-contrast 3-d imaging of breast tissue specimens as a possible adjunct to histopathology,” IEEE Transactions on Medical Imaging, vol. 37, no. 12, pp. 2642–2650, 2018.
  * [10] L. A. Peña, S. Donato, D. Bonazza, L. Brombal, F. Martellani, F. Arfelli, G. Tromba, and R. Longo, “Multiscale x-ray phase-contrast tomography: From breast CT to micro-CT for virtual histology,” Physica Medica, vol. 112, p. 102640, Aug. 2023.
  * [11] C. K. Hagen, P. Maghsoudlou, G. Totonelli, P. C. Diemoz, M. Endrizzi, L. Rigon, R.-H. Menk, F. Arfelli, D. Dreossi, E. Brun, P. Coan, A. Bravin, P. D. Coppi, and A. Olivo, “High contrast microstructural visualization of natural acellular matrices by means of phase-based x-ray tomography,” Scientific Reports, vol. 5, Dec. 2015.
  * [12] U. Bonse and M. Hart, “An x-ray interferometer,” Applied Physics Letters, vol. 6, no. 8, pp. 155–156, 1965.
  * [13] A. Momose, T. Takeda, Y. Itai, and K. Hirano, “Phase–contrast x–ray computed tomography for observing biological soft tissues,” Nature Medicine, vol. 2, pp. 473–475, Apr. 1996.
  * [14] A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, and I. Schelokov, “On the possibilities of x-ray phase contrast microimaging by coherent high-energy synchrotron radiation,” Review of Scientific Instruments, vol. 66, pp. 5486–5492, Dec. 1995.
  * [15] K. A. Nugent, T. E. Gureyev, D. F. Cookson, D. Paganin, and Z. Barnea, “Quantitative phase imaging using hard x rays,” Physical Review Letters, vol. 77, pp. 2961–2964, Sept. 1996.
  * [16] D. Chapman, W. Thomlinson, R. E. Johnston, D. Washburn, E. Pisano, N. Gmür, Z. Zhong, R. Menk, F. Arfelli, and D. Sayers, “Diffraction enhanced x-ray imaging,” Physics in Medicine and Biology, vol. 42, pp. 2015–2025, Nov. 1997.
  * [17] A. Momose, “Phase-sensitive imaging and phase tomography using x-ray interferometers,” Optics Express, vol. 11, p. 2303, Sept. 2003.
  * [18] K. S. Morgan, D. M. Paganin, and K. K. W. Siu, “X-ray phase imaging with a paper analyzer,” Applied Physics Letters, vol. 100, p. 124102, 03 2012.
  * [19] S. Berujon, H. Wang, and K. Sawhney, “X-ray multimodal imaging using a random-phase object,” Physical Review A, vol. 86, Dec. 2012.
  * [20] A. Olivo, F. Arfelli, G. Cantatore, R. Longo, R. H. Menk, S. Pani, M. Prest, P. Poropat, L. Rigon, G. Tromba, E. Vallazza, and E. Castelli, “An innovative digital imaging set-up allowing a low-dose approach to phase contrast applications in the medical field,” Medical Physics, vol. 28, pp. 1610–1619, Aug. 2001.
  * [21] K. S. Morgan, D. M. Paganin, and K. K. W. Siu, “Quantitative single-exposure x-ray phase contrast imaging using a single attenuation grid,” Optics Express, vol. 19, p. 19781, Sept. 2011.
  * [22] F. A. Vittoria, M. Endrizzi, P. C. Diemoz, A. Zamir, U. H. Wagner, C. Rau, I. K. Robinson, and A. Olivo, “X–ray absorption, phase and dark–field tomography through a beam tracking approach,” Scientific Reports, vol. 5, Nov. 2015.
  * [23] S. W. Wilkins, “Improved x-ray optics, especially for phase contrast imaging,” 1995.
  * [24] G. B. Provinciali, A. Cedola, O. de La Rochefoucauld, and P. Zeitoun, “Phase-contrast tomography with x-ray hartmann wavefront sensor,” in International Conference on X-Ray Lasers 2020 (D. Bleiner, ed.), SPIE, July 2021.
  * [25] G. B. Provinciali, M. Piponnier, L. Oudjedi, X. Levecq, F. Harms, A. Cedola, O. de La Rochefoucauld, and P. Zeitoun, “High-sensitivity x-ray phase imaging system based on a hartmann wavefront sensor,” Condensed Matter, vol. 7, p. 3, Dec. 2021.
  * [26] P. C. Diemoz, F. A. Vittoria, and A. Olivo, “Spatial resolution of edge illumination x-ray phase-contrast imaging,” Opt. Express, vol. 22, pp. 15514–15529, Jun 2014.
  * [27] C. K. Hagen, F. A. Vittoria, M. Endrizzi, and A. Olivo, “Theoretical framework for spatial resolution in edge-illumination x-ray tomography,” Phys. Rev. Appl., vol. 10, p. 054050, Nov 2018.
  * [28] M. Esposito, L. Massimi, I. Buchanan, J. D. Ferrara, M. Endrizzi, and A. Olivo, “A laboratory-based, low-energy, multi-modal x-ray microscope with user-defined resolution,” Applied Physics Letters, vol. 120, p. 234101, 06 2022.
  * [29] C. Navarrete-León, A. Doherty, A. Astolfo, C. K. Hagen, P. Munro, A. Olivo, and M. Endrizzi, “Two-directional beam-tracking for phase-sensitive x-ray tomography with laboratory sources,” AIP Conference Proceedings, vol. 2990, p. 050001, 09 2023.
  * [30] G. Lioliou, O. R. i Morgó, S. Marathe, K. Wanelik, S. Cipiccia, A. Olivo, and C. K. Hagen, “Cycloidal-spiral sampling for three-modal x-ray CT flyscans with two-dimensional phase sensitivity,” Scientific Reports, vol. 12, Dec. 2022.
  * [31] G. Lioliou, C. Navarrete-León, A. Astolfo, S. Savvidis, D. Bate, M. Endrizzi, C. K. Hagen, and A. Olivo, “A laboratory-based beam tracking x-ray imaging method achieving two-dimensional phase sensitivity and isotropic resolution with unidirectional undersampling,” Scientific Reports, vol. 13, May 2023.
  * [32] E. S. Dreier, C. Silvestre, J. Kehres, D. Turecek, M. Khalil, J. H. Hemmingsen, O. Hansen, J. Jakubek, R. Feidenhans’l, and U. L. Olsen, “Single-shot, omni-directional x-ray scattering imaging with a laboratory source and single-photon localization,” Optics Letters, vol. 45, p. 1021, Feb. 2020.
  * [33] C. Navarrete-León, A. Doherty, S. Savvidis, M. F. M. Gerli, G. Piredda, A. Astolfo, D. Bate, S. Cipiccia, C. K. Hagen, A. Olivo, and M. Endrizzi, “X-ray phase-contrast microtomography of soft tissues using a compact laboratory system with two-directional sensitivity,” Optica, vol. 10, pp. 880–887, Jul 2023.
  * [34] M. Guizar-Sicairos, S. T. Thurman, and J. R. Fienup, “Efficient subpixel image registration algorithms,” Optics Letters, vol. 33, p. 156, Jan. 2008.
  * [35] C. Kottler, C. David, F. Pfeiffer, and O. Bunk, “A two-directional approach for grating based differential phase contrast imaging using hard x-rays,” Opt. Express, vol. 15, pp. 1175–1181, Feb 2007.
  * [36] W. van Aarle, W. J. Palenstijn, J. De Beenhouwer, T. Altantzis, S. Bals, K. J. Batenburg, and J. Sijbers, “The astra toolbox: A platform for advanced algorithm development in electron tomography,” Ultramicroscopy, vol. 157, pp. 35–47, 2015.
  * [37] R. P. Nieuwenhuizen, K. A. Lidke, M. Bates, D. L. Puig, D. Grünwald, S. Stallinga, and B. Rieger, “Measuring image resolution in optical nanoscopy,” Nature methods, vol. 10, no. 6, pp. 557–562, 2013.
  * [38] A. Herbert and O. Burri, “Fourier ring correlation imagej plugin.” https://github.com/BIOP/ijp-frc, 2016.
  * [39] M. Van Heel and M. Schatz, “Fourier shell correlation threshold criteria,” Journal of structural biology, vol. 151, no. 3, pp. 250–262, 2005.


# The TerraByte Client: providing access to terabytes of plant data

Michael A. Beck1, Christopher P. Bidinosti, Christopher J. Henry, Manisha
Ajmani
University of Winnipeg, Winnipeg, MB Canada
Email<EMAIL_ADDRESS>
###### Abstract

In this paper we demonstrate the TerraByte Client, a software to download
user-defined plant datasets from a data portal hosted at Compute Canada. To
that end the client offers two key functionalities: (1) It allows the user to
get an overview on what data is available and a quick way to visually check
samples of that data. For this the client receives the results of queries to a
database and displays the number of images that fulfill the search criteria.
Furthermore, a sample can be downloaded within seconds to confirm that the
data suits the user’s needs. (2) The user can then download the specified data
to their own drive. This data is prepared into chunks server-side and sent to
the user’s end-system, where it is automatically extracted into individual
files. The first chunks of data are available for inspection after a brief
waiting period of a minute or less depending on available bandwidth and type
of data. The TerraByte Client has a full graphical user interface for easy
usage and uses end-to-end encryption. The user interface is built on top of a
low-level client. This architecture in combination of offering the client
program open-source makes it possible for the user to develop their own user
interface or use the client’s functionality directly. An example for direct
usage could be to download specific data on demand within a larger
application, such as training machine learning models.

## Availability

The _TerraByte Client_ is available at the TerraByte homepage
https://terrabyte.acs.uwinnipeg.ca/resources.html.

##  I Introduction

Digital agriculture is an upcoming field with the goal of using modern
technologies such as robotics, phenotyping, and machine learning to advance
and transform the way we grow food [1, 2, 3, 4, 5, 6, 7]. The vision of
autonomous agents making plant-by-plant decisions within the field motivates
this research. However, we often encounter, especially in the area of machine
learning, a bottleneck of _labelled data_ on which models can be trained and
improved upon [8, 9]. Labelled data, at a minimum, is the data itself (for
example an image) plus information about what the data represents (for example
the label “Soybean” attached to an image). Often, the information provided –
metadata – goes beyond simple labelling (for example also tracking the age of
the Soybean plant, time and day of imaging, geolocation, etc.). In machine
learning, for example with convolutional neural networks (CNN), we require
thousands if not hundreds of thousands of labelled images to perform well.
Labelled data is often created by manually labelling existing datasets and
scaled to millions of images by crowdsourcing the efforts [10, 11, 12].
However, in the case of plant identification, the knowledge required for
accurate labelling is more complex and uncommon and thus this approach of
crowdsourcing is not available.

Still, when looking at similar applications we can expect strong results, once
the challenge of having labelled image-data is resolved. CNNs have been
applied very successfully in image-classification tasks with results that
eventually outperformed human-level accuracy (see for example [13, 14, 15]).
What these models have in common is that they were trained upon a vast amount
of image data. The most prominent of these datasets being ImageNet [16, 17],
which consists of more than 1.2 million labelled images.

Figure 1: The architecture that provides plant-data to researchers, consisting
of the TerraByte Client, the server back-end, and an S3-compatible object
storage. In this image the buffer containing part 12 of a dataset is being
downloaded by the client, while in parallel the next part of the dataset is
being fetched from the object storage by the server into the second buffer.

At the TerraByte project [18] we have the goal to provide labelled plant-data
at a similar scale, thus, fuelling the development of better models. For this
purpose, we are deploying and developing several systems to collect plant-data
in the lab or in the field. Since 2019 we have collected over half a million
images in the field and over half a million images in the lab (as of the end
of 2021) [19]. The images collected in the lab further contain multiple plants
per image, which had been automatically cropped out, resulting in more than
1.5 million images of individual plants. All lab-images are fully labelled,
providing information about each plant’s species, planting date, imaging date
and more. This data is not only critical for machine learning purposes, but
also allows to filter the data by these metadata fields. The species of
plants, the age of the plant, and the imaging date, are just three example
fields we want to mention here; a full list of the metadata collected is given
in Table 1.

Generating and labelling data alone would not be sufficient, however, to
sustainably support machine learning research in digital agriculture. As
outlined, for example, by Lobet [20], data-collections must have long-term
access, storage, and curation. To this end we work together with _Compute
Canada_ 111Compute Canada’s services are actually merged into those of the
_Digital Research Alliance of Canada_. As the day of writing, this merge has
not happened, yet, thus we stick to _Compute Canada_ in this paper. The merger
has no effect on where or how to access our data. to host our data in a
permanent and unchanging location.

To provide this data to other researchers we are performing several steps: (1)
our data is regularly backed up on an object storage system hosted by _Compute
Canada_(CC). (2) The respective metadata is ingested on a database that runs
on a server-instance on CC. This instance also hosts the back-end application
through which data is delivered from the object storage to the end-user. (3)
We have developed a client program, the _TerraByte Client_ (TC), that has been
successfully tested on Windows, MacOS, and Unix. This client communicates with
the back-end application and can download data to the researcher’s system. The
TC can also be used to query the database without having to fully download the
dataset in its entirety. Thus, the end-user can quickly check how many images
are available under a specified set of filters (for example, how many images
of 2-3 week old soybeans are available). Furthermore, samples of such
specified datasets can be downloaded within seconds to confirm that these are
indeed images that suit the end-user’s goals. By providing these functions, we
think a dataset of this size is easier to navigate and data relevant to a
specific research question can be found faster. To reduce friction further, we
developed a graphical user interface that enables all functionality of the TC;
we also deliver the program as an installation-free standalone executable for
windows-based machines.

The rest of this paper is structured as follows: In Section II we give an
overview on how the TC is operated to obtain data. It follows a short
discussion on the client’s performance in Section III. The architecture of the
TC and the server’s API is provided in Section IV. Section V offers an outlook
on how the TC and the available datasets will evolve.

##  II Operation

Table 1:

Field |  Description
---|---
camera & lens |  Name of camera used and lens/model used
camera_pose |  Coordinates and orientation (pan and tilt) of the camera, when the image was taken
date & time |  Date and time when the images was taken
institute & room |  Location of the imaging system (more than one imaging system will be in place)
horizontal and vertical dimensions |  Dimensions of the image taken
tags |  For future use

Metadata fields specific for original (uncropped) images.

Field |  Description
---|---
plant_id |  Unique identifier for each individual plant
label |  Common name label for species of plant, e.g. “Soybean”
scientific_name |  E.g. “Fallopia convolvulus”
planting_date |  Date at which the individual plant was planted
age |  Age of the plant when image was taken in days
position_id |  An identifier that corresponds to a specific location inside the imagable volume at which the plant was located when the image was taken
x_min, x_max, y_min, y_max |  Coordinates of the subimage within the original image

Metadata fields specfic for cropped out single-plant images

(see Figure 3).

The interface of the TC is shown in Figure 2; it consists of several tabs each
relating to one class of data and one additional tab providing version-info
and useful links. In Figure 2 the active tab is set to _EAGL-I Data_ , which
corresponds to images that have been taken in the lab by our robotic imager
EAGL-I (see [19]). The tab itself is structured into three sections; via the
top two sections the user can define the data they are interested in, whereas
the third section contains controls for configuring the tool itself or for
interacting with the server-backend.

The filtering options of the first section relate directly to the metadata
information that is saved along the image data. Table 1 gives a detailed
description of all metadata tracked and we have highlighted fields in boldface
that are currently being used by the TC. All filters default to the least
restrictive setting and combine in an AND-relation. For example, setting a
maximal age of 10 days and a specific plant ID will only give images of that
plant ID in which the plant was less than 10 days old. Or stated differently,
activating additional filters will generally result in less data fitting the
filters. The only exception to this is the multi-select box in which different
plant-species can be chosen. Choosing no plant at all in this box deactivates
the filter completely, whereas multiple selections will result in every plant
that fits into any of the selected species being returned. Thus, within the
species-filter we have an OR-relation, while the species-filter itself is
combined via AND with the other filters.

Figure 2: The TerraByte Client user interface currently showing the _EAGL-I
Data_ tab, relating to images taken in the lab.

Through the second section in the interface the user can define what types of
data they are interested in. This relates to the way the robotic system EAGL-I
takes images and we give a brief explanation of that here: EAGL-I moves a
camera inside an imageable volume to different positions and from each
position one image is taken. Each of those images can contain multiple plants.
In a post-processing step each plant is cropped out from the original image
separately. This usually leads to more than one cropped out plant image per
original image (see Figure 3). We call the cropped out images _single plant
images_ and the original images _multiple plant images_. Further, a copy of
the original multiple plant image is created and annotated with the boundaries
of the cropped out single plant image. The four selections that can be made in
the second section of the interface relate to these three image types and the
respective metadata, which is given in a JSON file. Selecting several
filetypes will result in all checked filetypes to be returned from the server.
Finally, the second section contains a drop-down list that contains pre-
compiled datasets. These datasets are permanently stored on the server-backend
(in contrast to user-defined datasets that are being pulled from the object
storage) and thus can be delivered immediately by the server.

Once the user has defined a dataset through the first two sections, they can
query the server in three different ways through the third section of the
interface.

  * •

First, the user can check how many files fulfill the specified filter, without
committing to any downloads, by using the “Check Query” button. The
information returned by the server is displayed in the TC as the number of
files and in how many parts these files would be downloaded.

  * •

By using the “Get Sample” button the TC will download a small randomly
selected sample per selected filetype. For example, if “Single Plant Images”
and “Multiple Plant Images” is selected the TC will download 20 random images
of each filetype.

  * •

Finally, the user can download all files resulting from a query by using the
“Download” button. This will download all files that fit the specified filters
and filetypes.

Three additional buttons are available:

  * •

“Set Credentials”: This opens a dialog in which the user can enter or change
their user credentials consisting of username and password. The credentials
are saved locally and are reused at future launches of the TC.

  * •

“Set Sample Path”: This allows the user to select in which folder samples
should be saved to.

  * •

“Set Download Path”: This allows the user to select in which folder downloads
through the “Download”-button will be saved to.

The functionality of the TC when using the “Fielddata”-tab is similar to the
one under the “EAGL-I Data”-tab. Future datasets that will be available
through the TC client will be accessed through additional tabs, one per
dataset, each containing a similar interface.

Figure 3: An original image (A) overlayed with the areas that are being
cropped out as single-plant images (B-D)

##  III Performance

The TC’s performance depends on several factors that lie outside of its
control:

  * •

The end-systems download-speed and write-speed to drive

  * •

The current load on Compute Canada’s internal network

  * •

The number of clients connected and using the server back-end

To investigate the impact of parallel usage, we tested the system with up to 6
users downloading files at the same time. Only 2 users (user 2 and user 3)
shared the same access network. In this test we asked the users to download a
predefined dataset of 6193 (542 MB) and measured the time until the data was
available on their system. This was compared to the time it took user 1 to
download the same files, while no other user is active, which took 130
seconds. We saw that the impact of parallel usage was measurable, yet small,
compared to the impact the end-to-end bandwidth between the client’s host and
the server has (see Table 2). For example, users 2, 3, and 6 encountered a
similar download time independent of whether they downloaded pre-compiled or
not pre-compiled datasets. This tells us that the bottleneck in their download
speed was not the process of fetching data from the object-storage. Indeed,
user-defined (not pre-compiled) datasets impose additional workload on the
server that consists of retrieving the data from the object storage and tar-
balling the files before sending them to the TC (see the next Section for
details). We also saw that downloading pre-compiled datasets is significantly
faster if the end-user has enough download bandwidth.

In our testing we further noticed that the performance of the TC and the
server-backend fluctuates even if the available bandwidth between client and
server remains stable and no other clients are active. We have observed a
decrease of download speed up to a factor of 2 in rare instances. When
investigating these temporary drops in performances we noticed that the
transfer-speed between object storage and server-backend lead to these longer
delays. We thus assume that the overall load on Compute Canada’s internal
network can impact the system’s overall performance.

Table 2:

| Download Time | Download Speed |  Download Time(precompiled)
---|---|---|---
User 1 | 316 seconds | 1.7 MB/s |  36 seconds
User 2 | 957 seconds | 0.56 MB/s |  840 seconds
User 3 | 1,006 seconds | 0.53 MB/s |  847 seconds
User 4 | 1096 seconds | 0.49 MB/s |  49 seconds
User 5 | 300 seconds | 1.8 MB/s |  58 seconds
User 6 | 972 seconds | 0.55 MB/s |  965 seconds

Performance Test of 6 users downloading 6193 files (542 MB) in parallel

We also performed a long running speed-test for downloading 95K single plant
images (9.7 GB) for a single active user, where the data was not precompiled.
We repeated that download at several different days to rule out fluctuations
stemming from Compute Canada’s internal traffic. The results had been
consistent within 1 minute. Similarly we downloaded 20K multiple plant images
(34.3 GB) and 15K field data images (105 GB) to compare the download speeds
with respect to different filetypes (single plant images are smaller files).
The results can be found in Table 3. We can observe that the end-to-end
download speed for single plant images lies far below the download speed for
the other two categories. The reason for that is the overhead of copying many
indiviual small files from the object storage to the server-backend within
Compute Canada’s network. For each file we can observe a significant overhead
spent on initializing the transfer, resulting in poorer performance.

Table 3:

| Download Time | Download Speed
---|---|---
Single plant images | 80 minutes | 2 MB/s
Multiple plant images | 34 minutes | 18 MB/s
Field data images | 95 minutes | 19 MB/s

Performance test for data from different categories.

##  IV Architecture

In this section we give a high-level description of the client and server
architecture, starting with the client.

###  IV-A Client Architecture

We focus the discussion on how downloading a user-defined dataset is
implemented on the client side. Downloading samples instead or querying the
database without starting a download at all are just variations and
simplifications of a full download.

The main functionality is defined in $\mathtt{simple\\_client.py}$, which
provides two key methods: $\mathtt{send\\_req}$ and $\mathtt{get\\_files}$.
The former is a wrapper to send HTTPS POST- and GET-requests to the server’s
API. It is used to communicate requests to the database or request a file to
be downloaded. In the former case the return of this method is a JSON-
formatted object containing the response of the server. When the user wishes
to download a dataset $\mathtt{send\\_req}$ is used to inform the server about
the query parameters. In response the client obtains a unique randomized job-
id together with the number of parts the download will consist of. Then the TC
main-loop starts the $\mathtt{get\\_files}$ method as a subprocess for this
job-number and the list of parts.

The $\mathtt{get\\_files}$ method repeatedly asks the server if a specific
part of a specific job-id is ready to download. To not overload the server
with these queries a dynamic backoff period is used, which increases between
successive requests. This backoff period resets after a request yields a
positive answer from the server. If the server answers that the part is ready
for download, $\mathtt{get\\_files}$ continues by a direct request to the tar-
file via $\mathtt{send\\_req}$. The such received data is then unpacked into
the download folder on the client side (specified by the user (see Section
II). If the server repeatedly answers that the data is not ready, the
$\mathtt{get\\_files}$ method will eventually reach a maximal number of tries
and abandon the download of that part. Once a part is downloaded (or the
maximal number of tries exceeded) the main-loop ends that subprocess and moves
on to download the next part in a new subprocess. This is repeated until all
parts are processed.

Since the main functionality of the client is contained within
$\mathtt{simple\\_client.py}$ an operation outside the provided graphical user
interface is possible. Indeed, the graphical user interface is just a
convenient way to construct well-formed HTTPS requests to the server. It is
thus possible to skip its usage entirely for the purpose of accessing the
server-side API directly within custom code. This allows, for example, to
fetch data from the server on-demand from within other methods, such as the
iterative training of machine learning models.

###  IV-B Server Architecture

The server is hosted on Compute Canada on an end-system separate from the data
object storage. For communication with the clients the server runs a stack of
Flask, MongoDB, and Celery. Flask provides an API over HTTPS, which is also
used to deliver the data to the clients. Besides requests that are related to
user-management the server handles two types of requests: database queries and
download requests.

Database queries from the client are a collection of query parameters which
are sanitized and transformed into MongoDB-queries server-side. The such
formed MongoDB-query is used to obtain the list of files the user wishes to
download in the future. The server also creates a job-id and allocates
resources for the upcoming download process. A summary of the file-list is
sent back to the client immediately, together with the job-id. Then, using
Celery, the server starts a subprocess that downloads the first part of the
requested data from the object storage. Once the first part is fully fetched
from the object storage those files are tar-balled and are now available for
download by the client. Without waiting for the client, the server immediately
starts the download of the second part following the same process. However,
the third part will only be fetched from the object storage after the first
part has been downloaded by the client and the tar-archive has been deleted at
the server-instance. The server continues this double-buffered approach,
downloading new parts from the object-storage while serving the previously
fetched part to the client (see Figure 1). This architecture has two
advantages: First, it serves the first parts after a short time to the client
thus being a responsive design. Second, we do not have to download the
requested dataset in its entirety to the server-instance. In fact, this would
be impossible for datasets that exceed the server-instance’s storage (the
object storage can host hundreds of TB of data, whereas the server instance
has a limited disk space of 5 TB).

##  V Future Work

We presented the TerraByte Client a tool for downloading user-defined datasets
from the TerraByte data-server. We continue our goal to provide labelled
plant-data to researchers and industry with the purpose of progressing digital
agriculture. To that end, we will extend the client as well as the data-portal
along several dimensions:

  * •

More variety in field-data: Future iterations of the data-portal will contain
drone data of field images, including test-plots of research facilities and
fields of local farmers. Further, we are developing a data-rover that will
collect image data in test-plots. More plant-variety in field-data is planned
as well.

  * •

Beyond RGB-data: We have several systems that go beyond simple RGB-imaging
that are being prepared to ingest data on a large into our database. Amongst
these systems are a 3D-scanner for phenotyping purposes, a hyperspectral
camera, and a photogrammetry station.

  * •

Additional imagers and plant variety: Next to continuing the operation of
EAGL-I, we are setting up additional imaging systems that will image a wider
variety of crops and weeds inside our laboratories.

We will further work on improving the overall performance of the client and
the server. Specifically, we will increase the end-to-end download speed for
user-defined datasets. One options is to have datasets prepared and saved
server-side, so that they can be downloaded later at a faster speed. For this
purpose, each user will have a limited allocation of disk-space on the server,
they can still continue to use the client in the way outlined above as well
without space-limitations. Another options is to create archives on the object
storage each containing a subset of the data, such that if one file inside a
subset is queried it is very likely that other files from that archive are
also queried (e.g., one specific day of imaging). Since the data querried
usually is highly correlated along such a clustering, transferring an archive
from the object storage to the server-backend can yield many relevant files
while only paying the transfer overhead once. Indeed, we observed in tests a
tranfer speed of more than 150 Mb/s from object storage to the server
instance. This would shift the bottleneck of end-to-end transfer to the
connectin between the TC and the server instance (for which we observed upload
speeds of 20 Mb/s). Further improvements can be achieved by implementing a
cache on the level of downloaded archives and a cache on the level of
extracted (or copied) files to speed-up follow-up queries.

## References

  * [1] Ilker Ünal and Mehmet Topakci.  Design of a remote-controlled and GPS-guided autonomous robot for precision farming.  International Journal of Advanced Robotic Systems, 12:194, 2015\.
  * [2] Avital Bechar and ClÃ©ment Vigneault.  Agricultural robots for field operations. part 2: Operations and systems.  Biosystems Engineering, 153:110–128, 2017.
  * [3] Avital Bechar and ClÃ©ment Vigneault.  Agricultural robots for field operations. part 2: Operations and systems.  Biosystems Engineering, 153:110–128, 2017.
  * [4] Redmond Ramin Shamshiri, C. Weltzien, I. A. Hameed, I. J. Yule, T. E. Grift, S. K. Balasundram, L. Pitonakova, Desa Ahmad, and Girish Chowdhary.  Research and development in agricultural robotics: a perspective of digital farming.  International Journal of Agricultural and Biological Engineering, 11:1–14, 2018.
  * [5] Tom Duckett, Simon Pearson, Simon Blackmore, and Bruce Grieve.  Agricultural robotics: The future of robotic agriculture.  CoRR, abs/1806.06762, 2018.
  * [6] J.E. Relf-Eckstein, Anna T. Ballantyne, and Peter W.B. Phillips.  Farming reimagined: a case study of autonomous farm equipment and creating an innovation opportunity space for broadacre smart farming.  NJAS - Wageningen Journal of Life Sciences, 90-91:100307, 2019.
  * [7] Manlio Bacco, Paolo Barsocchi, Erina Ferro, Alberto Gotta, and Massimiliano Ruggeri.  The digitisation of agriculture: a survey of research activities on smart farming.  Array, 3-4:100009, 2019.
  * [8] Konstantinos G. Liakos, Patrizia Busato, Dimitrios Moshou, Simon Pearson, and Dionysis Bochtis.  Machine learning in agriculture: A review.  Sensors, 18(8), 2018.
  * [9] Kirtan Jha, Aalap Doshi, Poojan Patel, and Manan Shah.  A comprehensive review on automation in agriculture using artificial intelligence.  Artificial Intelligence in Agriculture, 2:1–12, 2019.
  * [10] Bryan C Russell, Antonio Torralba, Kevin P Murphy, and William T Freeman.  LabelMe: a database and web-based tool for image annotation.  International Journal of Computer Vision, 77:157–173, 2008.
  * [11] Michael Buhrmester, Tracy Kwang, and Samuel D Gosling.  Amazon’s mechanical turk: A new source of inexpensive, yet high-quality data?  Perspectives on Psychological Science, 6:3–5, 2016.
  * [12] Christopher J Rapson, Boon-Chong Seet, M Asif Naeem, Jeong Eun Lee, Mahmoud Al-Sarayreh, and Reinhard Klette.  Reducing the pain: A novel tool for efficient ground-truth labelling in images.  In Proc. IEEE IVCNZ, pages 1–9, 2018.
  * [13] Karen Simonyan and Andrew Zisserman.  Very deep convolutional networks for large-scale image recognition, 2014\.
  * [14] Yaniv Taigman, Ming Yang, Marc’Aurelio Ranzato, and Lior Wolf.  DeepFace: Closing the gap to human-level performance in face verification.  In Proc. IEEE CVPR, 2014.
  * [15] Christopher J. Henry, Christopher D. Storie, Muthu Palaniappan, Victor Alhassan, Mallikarjun Swamy, Damilola Aleshinloye, Andrew Curtis, and Daeyoun Kim.  Automated LULC map production using deep neural networks.  International Journal of Remote Sensing, 40:4416–4440, 2019.
  * [16] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Delving deep into rectifiers: Surpassing human-level performance on imagenet classification.  In Proceedings of the IEEE International Conference on Computer Vision (ICCV), December 2015.
  * [17] Olga Russakovsky, Jia Deng, Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bernstein, Alexander C. Berg, and Li Fei-Fei.  ImageNet large scale visual recognition challenge.  International Journal of Computer Vision, 115:211–252, 2015.
  * [18] Michael Beck, Christopher Bidinosti, and Christopher Henry.  Terrabyte homepage, May 2021c.
  * [19] Michael A. Beck, Chen-Yi Liu, Christopher P. Bidinosti, Christopher J. Henry, Cara M. Godee, and Manisha Ajmani.  An embedded system for the automated generation of labeled plant images to enable machine learning applications in agriculture.  PLOS ONE, 15(12):1–23, 12 2020.
  * [20] Guillaume Lobet.  Image analysis in plant sciences: Publish then perish.  Trends in Plant Science, 22(7):559–566, 2017.