diff --git a/cara/monte_carlo/data.py b/cara/monte_carlo/data.py
new file mode 100644
index 00000000..79e25517
--- /dev/null
+++ b/cara/monte_carlo/data.py
@@ -0,0 +1,56 @@
+import numpy as np
+
+import cara.monte_carlo as mc
+from cara.monte_carlo.sampleable import Normal,LogNormal,LogCustomKernel
+
+
+# From CERN-OPEN-2021-04 and refererences therein
+activity_distributions = {
+    'Seated': mc.Activity(LogNormal(-0.6872121723362303, 0.10498338229297108),
+                          LogNormal(-0.6872121723362303, 0.10498338229297108)),
+
+    'Standing': mc.Activity(LogNormal(-0.5742377578494785, 0.09373162411398223),
+                            LogNormal(-0.5742377578494785, 0.09373162411398223)),
+
+    'Light activity': mc.Activity(LogNormal(0.21380242785625422,0.09435378091059601),
+                                  LogNormal(0.21380242785625422,0.09435378091059601)),
+
+    'Moderate activity': mc.Activity(LogNormal(0.551771330362601, 0.1894616357138137),
+                                     LogNormal(0.551771330362601, 0.1894616357138137)),
+
+    'Heavy exercise': mc.Activity(LogNormal(1.1644665696723049, 0.21744554768657565),
+                                  LogNormal(1.1644665696723049, 0.21744554768657565)),
+}
+
+
+# From CERN-OPEN-2021-04 and refererences therein
+symptomatic_vl_frequencies = LogCustomKernel(
+    np.array((2.46032, 2.67431, 2.85434, 3.06155, 3.25856, 3.47256, 3.66957, 3.85979, 4.09927, 4.27081,
+     4.47631, 4.66653, 4.87204, 5.10302, 5.27456, 5.46478, 5.6533, 5.88428, 6.07281, 6.30549,
+     6.48552, 6.64856, 6.85407, 7.10373, 7.30075, 7.47229, 7.66081, 7.85782, 8.05653, 8.27053,
+     8.48453, 8.65607, 8.90573, 9.06878, 9.27429, 9.473, 9.66152, 9.87552)),
+    np.array((0.001206885, 0.007851618, 0.008078144, 0.01502491, 0.013258014, 0.018528495, 0.020053765,
+     0.021896167, 0.022047184, 0.018604005, 0.01547796, 0.018075445, 0.021503523, 0.022349217,
+     0.025097721, 0.032875078, 0.030594727, 0.032573045, 0.034717482, 0.034792991,
+     0.033267721, 0.042887485, 0.036846816, 0.03876473, 0.045016819, 0.040063473, 0.04883754,
+     0.043944602, 0.048142864, 0.041588741, 0.048762031, 0.027921732, 0.033871788,
+     0.022122693, 0.016927718, 0.008833228, 0.00478598, 0.002807662)),
+    kernel_bandwidth=0.1
+)
+
+
+# From CERN-OPEN-2021-04 and refererences therein
+virus_distributions = {
+    'SARS_CoV_2': mc.SARSCoV2(
+                viral_load_in_sputum=symptomatic_vl_frequencies,
+                quantum_infectious_dose=100,
+                ),
+    'SARS_CoV_2_B117': mc.SARSCoV2(
+                viral_load_in_sputum=symptomatic_vl_frequencies,
+                quantum_infectious_dose=60,
+                ),
+    'SARS_CoV_2_P1': mc.SARSCoV2(
+                viral_load_in_sputum=symptomatic_vl_frequencies,
+                quantum_infectious_dose=100/2.25,
+                ),
+}
diff --git a/cara/tests/test_predefined_distributions.py b/cara/tests/test_predefined_distributions.py
new file mode 100644
index 00000000..1e066432
--- /dev/null
+++ b/cara/tests/test_predefined_distributions.py
@@ -0,0 +1,47 @@
+import numpy as np
+import numpy.testing as npt
+import pytest
+
+from cara.monte_carlo.data import activity_distributions, virus_distributions
+
+# TODO: seed better the random number generators
+np.random.seed(0)
+
+
+# mean & std deviations from CERN-OPEN-2021-04 (Table 4)
+# NOTE: a mistake was corrected for the std deviation of the
+# "Moderate exercise" case (0.37 in the report, but should be 0.34)
+@pytest.mark.parametrize(
+    "distribution, mean, std",[
+        ['Seated', 0.51, 0.053],
+
+        ['Standing', 0.57, 0.053],
+
+        ['Light activity', 1.24, 0.12],
+
+        ['Moderate activity', 1.77, 0.34],
+
+        ['Heavy exercise', 3.28, 0.72,],
+    ]
+)
+def test_activity_distributions(distribution, mean, std):
+    activity = activity_distributions[distribution].build_model(size=1000000)
+    npt.assert_allclose(activity.inhalation_rate.mean(), mean, atol=0.01)
+    npt.assert_allclose(activity.inhalation_rate.std(), std, atol=0.01)
+
+
+# mean & std deviations from CERN-OPEN-2021-04 (Table 4) - with a refined
+# precision on the values
+@pytest.mark.parametrize(
+    "distribution, mean, std",[
+        ['SARS_CoV_2', 6.59, 1.74],
+
+        ['SARS_CoV_2_B117', 6.59, 1.74],
+
+        ['SARS_CoV_2_P1', 6.59, 1.74],
+    ]
+)
+def test_viral_load_logdistribution(distribution, mean, std):
+    virus = virus_distributions[distribution].build_model(size=1000000)
+    npt.assert_allclose(np.log10(virus.viral_load_in_sputum).mean(), mean, atol=0.01)
+    npt.assert_allclose(np.log10(virus.viral_load_in_sputum).std(), std, atol=0.01)
diff --git a/cara/tests/test_sampleable_distribution.py b/cara/tests/test_sampleable_distribution.py
index e3801c7e..585b4dde 100644
--- a/cara/tests/test_sampleable_distribution.py
+++ b/cara/tests/test_sampleable_distribution.py
@@ -3,9 +3,11 @@ import numpy.testing as npt
 import pytest
 
 from cara.monte_carlo import sampleable
+
 # TODO: seed better the random number generators
 np.random.seed(2000)
 
+
 @pytest.mark.parametrize(
     "mean, std",[
         [1., 0.5],