diff --git a/cara/model_scenarios_paper.py b/cara/model_scenarios_paper.py
index 7b038ae4..bb3b460d 100644
--- a/cara/model_scenarios_paper.py
+++ b/cara/model_scenarios_paper.py
@@ -1,5 +1,5 @@
 from cara import models
-from cara.monte_carlo.data import activity_distributions
+from cara.monte_carlo.data import activity_distributions, virus_distributions
 import cara.monte_carlo as mc
 import numpy as np
 
@@ -28,6 +28,60 @@ yann_er = [8396.78166, 45324.55964, 400054.0827]
 
 ######### Standard exposure models ###########
 
+######### Breathing model ###########
+def breathing_exposure():
+    exposure_mc = mc.ExposureModel(
+        concentration_model=mc.ConcentrationModel(
+            room=models.Room(volume=100, humidity=0.5),
+            ventilation=models.AirChange(
+                active=models.SpecificInterval(((0, 24),)),
+                air_exch=0.25,
+            ),
+            infected=mc.InfectedPopulation(
+                number=1,
+                virus=virus_distributions['SARS_CoV_2'],
+                presence=mc.SpecificInterval(((0, 2),)),
+                mask=models.Mask.types["No mask"],
+                activity=activity_distributions['Seated'],
+                expiration=models.Expiration.types['Breathing'],
+            ),
+        ),
+        exposed=mc.Population(
+            number=14,
+            presence=mc.SpecificInterval(((0, 2),)),
+            activity=models.Activity.types['Seated'],
+            mask=models.Mask.types["No mask"],
+        ),
+    )
+    return exposure_mc
+
+######### Talking model ###########
+def talking_exposure():
+    exposure_mc = mc.ExposureModel(
+        concentration_model=mc.ConcentrationModel(
+            room=models.Room(volume=100, humidity=0.5),
+            ventilation=models.AirChange(
+                active=models.SpecificInterval(((0, 24),)),
+                air_exch=0.25,
+            ),
+            infected=mc.InfectedPopulation(
+                number=1,
+                virus=virus_distributions['SARS_CoV_2'],
+                presence=mc.SpecificInterval(((0, 2),)),
+                mask=models.Mask.types["No mask"],
+                activity=activity_distributions['Seated'],
+                expiration=models.Expiration.types['Talking'],
+            ),
+        ),
+        exposed=mc.Population(
+            number=14,
+            presence=mc.SpecificInterval(((0, 2),)),
+            activity=models.Activity.types['Seated'],
+            mask=models.Mask.types["No mask"],
+        ),
+    )
+    return exposure_mc
+
 ######### Breathing model for specific viral load ###########
 def breathing_exposure_vl(vl):
     exposure_mc = mc.ExposureModel(
diff --git a/cara/plot_output.py b/cara/plot_output.py
index 4af9af8b..8febdcfc 100644
--- a/cara/plot_output.py
+++ b/cara/plot_output.py
@@ -4,7 +4,8 @@ Date: <date>
 Code version: <code version>
 Availability: <where it's located> """
 
-from cara.montecarlo import *
+from cara import model_scenarios_paper
+from cara.results_paper import *
 from cara.test_plots import *
 
 # Exhaled virions while talking, seated #
@@ -21,7 +22,23 @@ print('\n<<<<<<<<<<< Vlout for Talking, seated with chosen Cn,L >>>>>>>>>>>')
 
 # Exhaled virions while breathing according to BLO model, seated #
 print('\n<<<<<<<<<<< Vlout for Breathing, seated with chosen Cn,B >>>>>>>>>>>')
-exposure_model_from_vl_breathing_cn()
+#exposure_model_from_vl_breathing_cn()
+print('\n')
+############ Statistical Data ############
+
+############ Breathing model ############
+exposure_mc = model_scenarios_paper.breathing_exposure()
+exposure_model = exposure_mc.build_model(size=SAMPLE_SIZE)
+emission_rate = exposure_model.concentration_model.infected.emission_rate_when_present(cn_B=0.06, cn_L=0.2)
+print('\n<<<<<<<<<<< Breathing model statistics >>>>>>>>>>>')
+print_er_info(emission_rate)
+
+############ Talking model ############
+exposure_mc = model_scenarios_paper.talking_exposure()
+exposure_model = exposure_mc.build_model(size=SAMPLE_SIZE)
+emission_rate = exposure_model.concentration_model.infected.emission_rate_when_present(cn_B=0.06, cn_L=0.2)
+print('\n<<<<<<<<<<< Talking model statistics >>>>>>>>>>>')
+print_er_info(emission_rate)
 
 ############ Used for testing ############
 #exposure_model_from_vl_talking_new_points()
diff --git a/cara/montecarlo.py b/cara/results_paper.py
similarity index 100%
rename from cara/montecarlo.py
rename to cara/results_paper.py