From e4807b6b2cfc8ae901be0f7889c974c9481a3cdf Mon Sep 17 00:00:00 2001
From: Luis Aleixo <luis.aleixo@cern.ch>
Date: Mon, 28 Feb 2022 11:32:47 +0100
Subject: [PATCH] Added method to generate the BLO curve

---
 cara/short_range_plots/results.py | 37 +++++++++++++++++--------------
 cara/short_range_plots/scripts.py | 18 +++++++++++++++
 2 files changed, 38 insertions(+), 17 deletions(-)

diff --git a/cara/short_range_plots/results.py b/cara/short_range_plots/results.py
index ce33297b..56c4b37e 100644
--- a/cara/short_range_plots/results.py
+++ b/cara/short_range_plots/results.py
@@ -39,25 +39,25 @@ from cara.monte_carlo.data import symptomatic_vl_frequencies
 #                             thickness = [2, 2])
 
 # print('\n<<<<<<<<<<< Peak viral concentration with short range interactions for baseline scenarios >>>>>>>>>>>')
-concentration_curve(models=[exposure_module_with_short_range(
-                                    activity='Light activity',
-                                    expiration={"Speaking": 1, "Breathing": 2},
-                                    mask='No mask',
-                                    sr_presence=[(10.5, 11.0)],
-                                    sr_activities=['Speaking']),
-                            exposure_module_without_short_range(
-                                activity='Light activity',
-                                expiration={"Speaking": 1, "Breathing": 2},
-                                mask='No mask',)
-                            ],
-                            labels = ['Concentration with short range interactions', 'Background (long-range) concentration'],
-                            labelsDose = ['Dose (full)', 'Dose (long-range)'],
-                            colors = ['salmon', 'royalblue'],
-                            linestyles = ['-', '--'],
-                            thickness = [2, 2])
+# concentration_curve(models=[exposure_module_with_short_range(
+#                                     activity='Light activity',
+#                                     expiration={"Speaking": 1, "Breathing": 2},
+#                                     mask='No mask',
+#                                     sr_presence=[(10.5, 11.0)],
+#                                     sr_activities=['Speaking']),
+#                             exposure_module_without_short_range(
+#                                 activity='Light activity',
+#                                 expiration={"Speaking": 1, "Breathing": 2},
+#                                 mask='No mask',)
+#                             ],
+#                             labels = ['Concentration with short range interactions', 'Background (long-range) concentration'],
+#                             labelsDose = ['Dose (full)', 'Dose (long-range)'],
+#                             colors = ['salmon', 'royalblue'],
+#                             linestyles = ['-', '--'],
+#                             thickness = [2, 2])
 
 
-print('\n<<<<<<<<<<< Dose vs SR exposure time >>>>>>>>>>>')
+# print('\n<<<<<<<<<<< Dose vs SR exposure time >>>>>>>>>>>')
 #Always assume 1h for the short range interactions.
 #Always assume that in each model there is only ONE short range interaction.
 # plot_vD_vs_exposure_time(exp_models = [
@@ -80,3 +80,6 @@ print('\n<<<<<<<<<<< Dose vs SR exposure time >>>>>>>>>>>')
 #                          time_in_minutes=True,
 #                          normalize_y_axis=True)
 
+print('\n<<<<<<<<<<< BLO curve >>>>>>>>>>>')
+generate_BLO_curve(activity='Speaking')
+    
\ No newline at end of file
diff --git a/cara/short_range_plots/scripts.py b/cara/short_range_plots/scripts.py
index 64550008..1b384c2c 100644
--- a/cara/short_range_plots/scripts.py
+++ b/cara/short_range_plots/scripts.py
@@ -223,3 +223,21 @@ def plot_vD_vs_exposure_time(exp_models: typing.List[mc.ExposureModel], labels,
     plt.ylabel('Mean cumulative dose\n(infectious virus)', fontsize=12)
     plt.legend()
     plt.show()
+
+def generate_BLO_curve(activity):
+    diameters = short_range_expiration_distributions[activity].diameter
+    BLOcurve = BLOmodel(expiration_BLO_factors[activity]).distribution(diameters
+            )/BLOmodel(expiration_BLO_factors[activity]).integrate(0.1,100)
+    # you have to divide by the integral to get the normalized distribution,
+    # so that you can compare with the normalized histogram
+
+    plt.xscale('log')
+    plt.yscale('log')
+
+    plt.plot(diameters,BLOcurve,'.',ms=5)
+    plt.hist(diameters,bins=100,range=(0,100),density=True)
+    plt.title(activity)
+
+    plt.xlabel("diameters (microns)")
+    plt.ylabel("Normalized distribution and histogram")
+    plt.show()
\ No newline at end of file