1import math
  2import os, simpy
  3from random import expovariate, uniform
  4from sim_tools.time_dependent import NSPPThinning
  5import pandas as pd
  6# from random import expovariate
  7from utils import Utils, SeededDistribution
  8from class_patient import Patient
  9# Revised class for HEMS availability
 10from class_hems_availability import HEMSAvailability
 11from class_hems import HEMS
 12from class_ambulance import Ambulance
 13from datetime import timedelta
 14import warnings
 15import numpy as np
 16from math import floor
 17warnings.filterwarnings("ignore", category=RuntimeWarning)
 18# import all distributions
 19import ast
 20from numpy.random import SeedSequence
 21from typing import List, Dict, Tuple
 22from numpy.random import SeedSequence, default_rng
 23import random
 24
 25import logging
 26logging.basicConfig(filename='log.txt', filemode="w", level=logging.DEBUG, format='')
 27
 28class DES_HEMS:
 29    """
 30        # The DES_HEMS class
 31
 32        This class contains everything require to undertake a single DES run for multiple patients
 33        over the specified duration.
 34
 35        NOTE: The unit of sim is minutes
 36
 37    """
 38
 39    def __init__(self,
 40                run_number: int,
 41                sim_duration: int,
 42                warm_up_duration: int,
 43                sim_start_date: str,
 44                amb_data: bool,
 45                random_seed: int,
 46                demand_increase_percent: float,
 47                activity_duration_multiplier: float,
 48                print_debug_messages: bool):
 49
 50        self.run_number = run_number + 1 # Add 1 so we don't have a run_number 0
 51        self.sim_duration = sim_duration
 52        self.warm_up_duration = warm_up_duration
 53        self.sim_start_date = sim_start_date
 54
 55        self.random_seed = random_seed
 56        self.random_seed_sequence = SeedSequence(self.random_seed)
 57
 58        self.print_debug_messages = print_debug_messages
 59
 60        self.utils = Utils(
 61            master_seed=self.random_seed_sequence.spawn(1)[0],
 62            print_debug_messages=self.print_debug_messages
 63            )
 64
 65        self.utils.setup_seeds()
 66
 67        self.demand_increase_percent = demand_increase_percent
 68
 69        # Option for sampling calls per day by season or quarter
 70        self.daily_calls_by_quarter_or_season = 'quarter'
 71
 72        # Option to include/exclude ambulance service cases in addition to HEMS
 73        self.amb_data = amb_data
 74        #self.debug(f"Ambulance data values is {self.amb_data}")
 75
 76        self.all_results_location = self.utils.ALL_RESULTS_CSV
 77        self.run_results_location = self.utils.RUN_RESULTS_CSV
 78
 79        self.env = simpy.Environment()
 80        self.patient_counter = 0
 81        self.calls_today = 0
 82        self.new_day = pd.to_datetime("1900-01-01").date
 83        self.new_hour = -1
 84
 85        self.hems_resources = HEMSAvailability(env=self.env,
 86                                               sim_start_date=sim_start_date,
 87                                               sim_duration=sim_duration,
 88                                               print_debug_messages=self.print_debug_messages,
 89                                               master_seed=self.random_seed_sequence,
 90                                               utility=self.utils
 91                                               )
 92
 93        # Set up empty list to store results prior to conversion to dataframe
 94        self.results_list = []
 95
 96        # Set up data frame to capture time points etc. during the simulation
 97        # We might not need all of these, but simpler to capture them all for now.
 98
 99        # This might be better as a separate 'data collection' class of some sort.
100        # I am thinking that each resource will need its own entry for any given
101        # patient to account for dual response (ambulance service and HEMS)
102        # stand downs etc.
103        self.results_df = None
104
105        # self.inter_arrival_times_df = pd.read_csv('distribution_data/inter_arrival_times.csv')
106
107        self.activity_duration_multiplier = activity_duration_multiplier
108
109        # self.seeded_dists = self.utils.build_seeded_distributions(
110        #     self.utils.activity_time_distr,
111        #     master_seed=self.random_seed
112        #     )
113
114    def debug(self, message: str):
115        if self.print_debug_messages:
116            logging.debug(message)
117            #print(message)
118
119
120    def calls_per_hour(self, quarter: int) -> dict:
121        """
122            Function to return a dictionary of keys representing the hour of day
123            and value representing the number of calls in that hour
124        """
125
126        # self.debug(f"There are going to be {self.calls_today} calls today and the current hour is {current_hour}")
127
128        hourly_activity = self.utils.hourly_arrival_by_qtr_probs_df
129        hourly_activity_for_qtr = hourly_activity[hourly_activity['quarter'] == quarter][['hour','proportion']]
130
131        calls_in_hours = []
132
133        for i in range(0, self.calls_today):
134            hour = pd.Series.sample(hourly_activity_for_qtr['hour'], weights = hourly_activity_for_qtr['proportion'],
135                                random_state=self.utils.rngs["calls_per_hour"]).iloc[0]
136            #self.debug(f"Chosen hour is {hour}")
137            calls_in_hours.append(hour)
138
139        calls_in_hours.sort()
140
141        #self.debug(calls_in_hours)
142
143        d = {}
144
145        hours, counts = np.unique(calls_in_hours, return_counts=True)
146
147        for i in range(len(hours)):
148            d[hours[i]] = counts[i]
149
150        return d
151
152    def predetermine_call_arrival(self, current_hour: int, quarter: int) -> list:
153        """
154            Function to determine the number of calls in
155            24 hours and the inter-arrival rate for calls in that period
156            Returns a list of times that should be used in the yield timeout statement
157            in a patient generator
158
159        """
160
161        hourly_activity = self.utils.hourly_arrival_by_qtr_probs_df
162        hourly_activity_for_qtr = hourly_activity[hourly_activity['quarter'] == quarter][['hour','proportion']]
163
164        d = self.calls_per_hour(quarter)
165
166        ia_time = []
167
168        for index, row in hourly_activity_for_qtr.iterrows():
169            hour = row['hour']
170            if hour >= current_hour and hour in d:
171                count = d[hour]
172                if count > 0:
173                    scale = 60 / count  # mean of exponential = 1 / rate
174                    calc_ia_time = self.utils.rngs["predetermine_call_arrival"].exponential(scale=scale)
175                    tmp_ia_time = ((hour - current_hour) * 60) + calc_ia_time
176                    current_hour += floor(tmp_ia_time / 60)
177                    ia_time.append(tmp_ia_time)
178
179        return ia_time
180
181
182    def generate_calls(self):
183        """
184            **Call generator**
185
186            Generate calls (and patients) until current time equals sim_duration + warm_up_duration.
187            This method calculates number of calls per day and then distributes them according to distributions determined
188            from historic data
189
190        """
191
192        while self.env.now < (self.sim_duration + self.warm_up_duration):
193            # Get current day of week and hour of day
194            [dow, hod, weekday, month, qtr, current_dt] = self.utils.date_time_of_call(self.sim_start_date, self.env.now)
195
196            if(self.new_hour != hod):
197                self.new_hour = hod
198
199                #self.debug("new hour")
200
201            # If it is a new day, need to calculate how many calls
202            # in the next 24 hours
203            if(self.new_day != current_dt.date()):
204                # self.debug("It's a new day")
205                # self.debug(dow)
206                # self.debug(f"{self.new_day} and {current_dt.date}")
207
208                # Now have additional option of determining calls per day by quarter instead of season
209                #self.calls_today = int(self.utils.inc_per_day(qtr, self.daily_calls_by_quarter_or_season) * (self.demand_increase_percent))
210
211                # Try with actual sampled values instead
212                self.calls_today = int(self.utils.inc_per_day_samples(qtr, self.daily_calls_by_quarter_or_season) * (self.demand_increase_percent))
213
214                # self.debug(f"{current_dt.date()} There will be {self.calls_today} calls today")
215                #self.debug(f"{current_dt.date()} There will be {self.calls_today} calls today")
216
217                self.new_day = current_dt.date()
218
219                ia_dict = {}
220                ia_dict = self.calls_per_hour(qtr)
221
222                #self.debug(ia_dict)
223
224                # Also run scripts to check HEMS resources to see whether they are starting/finishing service
225                yield self.env.process(self.hems_resources.daily_servicing_check(current_dt, hod, month))
226
227
228            if self.calls_today > 0:
229                if hod in ia_dict.keys():
230                    minutes_elapsed = 0
231                    for i in range(0, ia_dict[hod]):
232                        if minutes_elapsed < 59:
233                            # Determine remaining time and sample a wait time
234                            remaining_minutes = 59 - minutes_elapsed
235                            # wait_time = random.randint(0, remaining_minutes)
236                            wait_time = int(self.utils.rngs["call_iat"].integers(0, remaining_minutes+1))
237
238                            yield self.env.timeout(wait_time)
239                            minutes_elapsed += wait_time
240
241                            self.env.process(self.generate_patient(dow, hod, weekday, month, qtr, current_dt))
242
243                            [dow, hod, weekday, month, qtr, current_dt] = self.utils.date_time_of_call(self.sim_start_date, self.env.now)
244
245
246                next_hr = current_dt.floor('h') + pd.Timedelta('1h')
247                yield self.env.timeout(math.ceil(pd.to_timedelta(next_hr - current_dt).total_seconds() / 60))
248
249                    #self.debug('Out of loop')
250            else:
251                # Skip to tomorrow
252
253                self.debug('Skip to tomorrow')
254
255                [dow, hod, weekday, month, qtr, current_dt] = self.utils.date_time_of_call(self.sim_start_date, self.env.now)
256                next_day = current_dt.floor('d') + pd.Timedelta(days=1)
257                self.debug("next day is {next_day}")
258                yield self.env.timeout(math.ceil(pd.to_timedelta(next_hr - current_dt).total_seconds() / 60))
259
260
261    def generate_patient(self, dow, hod, weekday, month, qtr, current_dt):
262
263        self.patient_counter += 1
264
265        # Create a new caller/patient
266        pt = Patient(self.patient_counter)
267
268        # Update patient instance with time-based values so the current time is known
269        pt.day = dow
270        pt.hour = hod
271        pt.weekday = weekday
272        pt.month = month
273        pt.qtr = qtr
274        pt.current_dt = current_dt
275
276        #self.debug(f"Patient {pt.id} incident date: {pt.current_dt}")
277
278        # Update patient instance with age, sex, AMPDS card, whether they are a HEMS' patient and if so, the HEMS result,
279        pt.ampds_card = self.utils.ampds_code_selection(pt.hour)
280        #self.debug(f"AMPDS card is {pt.ampds_card}")
281        pt.age = self.utils.age_sampling(pt.ampds_card, 115)
282        pt.sex = self.utils.sex_selection(pt.ampds_card)
283        hems_cc_or_ec = self.utils.care_category_selection(pt.ampds_card)
284        pt.hems_cc_or_ec = hems_cc_or_ec
285        self.debug(f"{pt.current_dt}: {pt.id} Pt allocated to {pt.hems_cc_or_ec} from AMPDS {pt.ampds_card}")
286
287        pt.pt_outcome = self.utils.pt_outcome_selection(pt.hems_cc_or_ec, int(qtr))
288        self.debug(f"{pt.current_dt}: {pt.id} Pt allocated to patient outcome: {pt.pt_outcome}")
289
290        self.add_patient_result_row(pt, "arrival", "arrival_departure")
291
292        if self.amb_data:
293            # TODO: We'll need the logic to decide whether it is an ambulance or HEMS case
294            # if ambulance data is being collected too.
295            self.debug("Ambulance case")
296            pt.hems_case = 1 if self.utils.rngs["hems_case"].uniform(0, 1) <= 0.5 else pt.hems_case == 0
297        else:
298            pt.hems_case = 1
299
300        if pt.hems_case == 1:
301            #self.debug(f"Going to callsign_group_selection with hour {pt.hour} and AMPDS {pt.ampds_card}")
302            # pt.hems_pref_callsign_group = self.utils.callsign_group_selection(pt.ampds_card)
303            # pt.hems_pref_callsign_group = self.utils.callsign_group_selection(pt.hems_cc_or_ec)
304            pt.hems_pref_callsign_group = self.utils.callsign_group_selection()
305            #self.debug(f"Callsign is {pt.hems_pref_callsign_group}")
306
307            # Some % of calls have a helicopter benefit
308            # Default to y for all patients
309            helicopter_benefit = 'n'
310
311            # Update for REG patients based on historically observed patterns
312            with open("distribution_data/proportion_jobs_heli_benefit.txt", "r") as file:
313                expected_prop_heli_benefit_jobs_reg = float(file.read().strip())
314            if pt.hems_cc_or_ec == 'REG':
315                helicopter_benefit = 'y' if self.utils.rngs["helicopter_benefit_from_reg"].uniform(0, 1) <= expected_prop_heli_benefit_jobs_reg else 'n'
316
317            # Update for CC patients based on historically observed patterns
318            with open("distribution_data/proportion_jobs_heli_benefit_cc.txt", "r") as file:
319                expected_prop_heli_benefit_jobs_cc = float(file.read().strip())
320            if pt.hems_cc_or_ec == 'CC':
321                helicopter_benefit = 'y' if self.utils.rngs["helicopter_benefit_from_cc"].uniform(0, 1) <= expected_prop_heli_benefit_jobs_cc else 'n'
322
323            # Update for EC patients based on historically observed patterns
324            with open("distribution_data/proportion_jobs_heli_benefit_ec.txt", "r") as file:
325                expected_prop_heli_benefit_jobs_ec = float(file.read().strip())
326            if pt.hems_cc_or_ec == 'EC':
327                helicopter_benefit = 'y' if self.utils.rngs["helicopter_benefit_from_ec"].uniform(0, 1) <= expected_prop_heli_benefit_jobs_ec else 'n'
328
329            # Following conversation with HT 21/5
330            # Also count as having had a helicopter benefit if the patient is conveyed by HEMS
331            # This is consistent with how things are done in the golden codes paper
332            # https://static-content.springer.com/esm/art%3A10.1186%2Fs13049-023-01094-w/MediaObjects/13049_2023_1094_MOESM1_ESM.pdf
333            # while having not always been coded in the historic dataset
334            # (TODO: SR 30/5 I have commented this out for now and we should revisit this at some point -
335            # as we might be slightly overestimating heli benefit patients as a result? These need to be
336            # more closely/cleverly linked or just covered by the historic dataset instead)
337
338            if pt.hems_result == "Patient Conveyed by HEMS":
339                helicopter_benefit = 'y'
340
341            pt.hems_helicopter_benefit = helicopter_benefit
342            self.add_patient_result_row(pt, pt.hems_cc_or_ec, "patient_care_category")
343            self.add_patient_result_row(pt, pt.hems_helicopter_benefit, "patient_helicopter_benefit")
344
345            # pt.hems_pref_vehicle_type = self.utils.vehicle_type_selection(pt.hems_pref_callsign_group)
346            pt.hems_pref_vehicle_type = self.utils.vehicle_type_selection_qtr(pt.hems_pref_callsign_group, int(qtr))
347            # pt.hems_pref_vehicle_type = 'helicopter'
348            #pt.hems_pref_callsign_group = '70'
349            #pt.hems_helicopter_benefit = 'y'
350
351            self.add_patient_result_row(pt, pt.hems_pref_callsign_group, "resource_preferred_resource_group")
352            self.add_patient_result_row(pt, pt.hems_pref_vehicle_type, "resource_preferred_vehicle_type")
353
354            if (pt.hems_cc_or_ec == 'CC' or pt.hems_cc_or_ec == 'EC') and self.utils.rngs["know_cc_ec_benefit"].uniform(0, 1) <= 0.5:
355                hems_res_list: list[HEMS|None, str, HEMS|None] = yield self.hems_resources.allocate_resource(pt)
356            else:
357                # Separate function to determine HEMS resource based on the preferred callsign group
358                # (which is sampled from historical data), the preferred vehicle type
359                hems_res_list: list[HEMS|None, str, HEMS|None] = yield self.hems_resources.allocate_regular_resource(pt)
360
361
362            self.debug(f"{pt.id} hems_res_list: {hems_res_list}")
363
364            hems_allocation = hems_res_list[0]
365
366            # This will either contain the other resource in a callsign_group and HEMS category (EC/CC) or None
367            hems_group_resource_allocation = hems_res_list[2]
368
369            self.add_patient_result_row(pt, hems_res_list[1], "resource_preferred_outcome")
370
371            if hems_allocation is not None:
372                #self.debug(f"allocated {hems_allocation.callsign}")
373
374                self.env.process(self.patient_journey(hems_allocation, pt, hems_group_resource_allocation))
375            else:
376                #self.debug(f"{pt.current_dt} No HEMS resource available - non-DAAT land crew sent")
377                self.env.process(self.patient_journey(None, pt, None))
378
379
380    def patient_journey(self, hems_res: HEMS|None, patient: Patient, secondary_hems_res: HEMS|None):
381        """
382            Send patient on their journey!
383        """
384
385        #self.debug(f"Patient journey triggered for {patient.id}")
386        #self.debug(f"patient journey Time: {self.env.now}")
387        # self.debug(hems_res.callsign)
388
389        try:
390            patient_enters_sim = self.env.now
391
392            # Note that 'add_patient_result_row' does its own check for whether it's currently within
393            # the warm-up period, so this does not need to be checked manually when adding result
394            # rows in this section
395            # not_in_warm_up_period = False if self.env.now < self.warm_up_duration else True
396
397            # Add variables for quick determination of whether certain parts of the process should
398            # be included per case
399            hems_case = True if patient.hems_case == 1 else False
400            hems_avail = True if hems_res is not None else False
401
402            if not hems_avail:
403                self.debug(f"Patient {patient.id}: No HEMS available")
404                self.add_patient_result_row(patient, "No HEMS available", "queue")
405                self.add_patient_result_row(patient, "No Resource Available", "resource_request_outcome")
406
407            # if hems_avail
408            else:
409
410                # Record selected vehicle type and callsign group in patient object
411                patient.hems_vehicle_type = hems_res.vehicle_type
412                patient.hems_registration = hems_res.registration
413                patient.callsign = hems_res.callsign
414
415                self.add_patient_result_row(patient, patient.hems_vehicle_type, "selected_vehicle_type")
416                self.add_patient_result_row(patient, patient.hems_callsign_group, "selected_callsign_group")
417
418                #patient.hems_result = self.utils.hems_result_by_callsign_group_and_vehicle_type_selection(patient.hems_callsign_group, patient.hems_vehicle_type)
419                #self.debug(f"{patient.hems_cc_or_ec} and {patient.hems_helicopter_benefit}")
420                #patient.hems_result = self.utils.hems_result_by_care_category_and_helicopter_benefit_selection(patient.hems_cc_or_ec, patient.hems_helicopter_benefit)
421
422                # Determine outcome
423                # If we know a resource has been allocated, we can determine the output from historical patterns
424                if hems_res:
425                    patient.hems_result = self.utils.hems_results_by_patient_outcome_and_time_of_day_and_quarter_and_vehicle_type_and_callsign_group_probs(
426                        patient.pt_outcome,
427                        int(patient.qtr),
428                        patient.hems_vehicle_type,
429                        patient.hems_callsign_group,
430                        int(patient.hour)
431                    )
432                else:
433                # Default to what is recorded when no resource sent
434                    patient.hems_result = "Unknown"
435
436                self.debug(f"{patient.current_dt}: PT_ID:{patient.id} Pt allocated to HEMS result: {patient.hems_result}")
437                self.add_patient_result_row(patient, "HEMS Resource Available", "resource_request_outcome")
438                self.add_patient_result_row(patient, hems_res.callsign, "resource_use")
439                self.debug(f"{patient.current_dt} Patient {patient.id} (preferred callsign group {patient.hems_pref_callsign_group}, preferred resource type {patient.hems_pref_vehicle_type}) sent resource {hems_res.callsign}")
440                self.add_patient_result_row(patient, hems_res.callsign, "callsign_group_resource_use")
441
442                # Check if HEMS result indicates that resource stood down before going mobile or en route
443                no_HEMS_at_scene = True if patient.hems_result in ["Stand Down"] else False
444                # Check if HEMS result indicates no leaving scene/at hospital times
445                no_HEMS_hospital = (
446                    True if patient.hems_result
447                    in ["Stand Down",  "Landed but no patient contact", "Patient Treated but not conveyed by HEMS"]
448                    else False
449                    )
450
451            #self.debug('Inside hems_avail')
452            if self.amb_data:
453                ambulance = Ambulance()
454
455            patient.time_in_sim = self.env.now - patient_enters_sim
456
457            if hems_case and hems_avail:
458                #self.debug(f"Adding result row with csg {patient.hems_callsign_group}")
459                self.add_patient_result_row(patient, "HEMS call start", "queue")
460
461            if self.amb_data:
462                self.add_patient_result_row(patient,  "AMB call start", "queue")
463
464            # Allocation time --------------
465            # Allocation will always take place if a resource is found
466
467            if hems_case and hems_avail:
468                # Calculate min and max permitted times.
469                allocation_time = (
470                    self.utils.activity_time(patient.hems_vehicle_type, 'time_allocation')
471                    * self.activity_duration_multiplier
472                    )
473                self.add_patient_result_row(patient, allocation_time, 'time_allocation')
474                self.add_patient_result_row(patient, "HEMS allocated to call", "queue")
475                #self.debug(f"Vehicle type {patient.hems_vehicle_type} and allocation time is {allocation_time}")
476                yield self.env.timeout(allocation_time)
477
478            if self.amb_data:
479                #self.debug('Ambulance allocation time')
480                yield self.env.timeout(180)
481
482            patient.time_in_sim = self.env.now - patient_enters_sim
483
484            # Mobilisation time ---------------
485
486            # Calculate mobile to time at scene (or stood down before)
487            if hems_case and hems_avail:
488                mobile_time = (
489                    self.utils.activity_time(patient.hems_vehicle_type, 'time_mobile')
490                    * self.activity_duration_multiplier
491                    )
492                self.add_patient_result_row(patient, mobile_time, 'time_mobile')
493                self.add_patient_result_row(patient,  "HEMS mobile", "queue")
494                yield self.env.timeout(mobile_time)
495
496            if self.amb_data:
497                # Determine allocation time for ambulance
498                # Yield until done
499                #self.debug('Ambulance time to going mobile')
500                self.debug('Ambulance mobile')
501                yield self.env.timeout(1)
502
503            patient.time_in_sim = self.env.now - patient_enters_sim
504
505            # On scene ---------------
506
507            if hems_case and hems_avail and not no_HEMS_at_scene:
508                tts_time = (
509                    self.utils.activity_time(patient.hems_vehicle_type, 'time_to_scene')
510                    * self.activity_duration_multiplier
511                    )
512                self.add_patient_result_row(patient, tts_time, "time_to_scene")
513                self.add_patient_result_row(patient,  "HEMS on scene", "queue")
514                yield self.env.timeout(tts_time)
515
516            if self.amb_data:
517                    # Determine allocation time for ambulance
518                    # Yield until done
519                    #self.debug('Ambulance time to scene')
520                    yield self.env.timeout(20)
521
522            patient.time_in_sim = self.env.now - patient_enters_sim
523
524            if self.amb_data:
525                self.debug('Ambulance stand down en route')
526
527
528            # Leaving scene ------------
529
530            if hems_case and hems_avail and not no_HEMS_at_scene:
531                tos_time = (
532                    self.utils.activity_time(patient.hems_vehicle_type, 'time_on_scene')
533                    * self.activity_duration_multiplier
534                    )
535                self.add_patient_result_row(patient, tos_time, "time_on_scene")
536                yield self.env.timeout(tos_time)
537
538                if no_HEMS_hospital:
539                    self.add_patient_result_row(patient, f"HEMS {patient.hems_result.lower()}", "queue")
540                else:
541                    self.add_patient_result_row(patient, "HEMS leaving scene", "queue")
542
543            if self.amb_data:
544                #self.debug('Ambulance on scene duration')
545                yield self.env.timeout(120)
546                self.debug('Ambulance leaving scene time')
547
548            patient.time_in_sim = self.env.now - patient_enters_sim
549
550            # Arrived destination time ------------
551
552            if hems_case and hems_avail and not no_HEMS_hospital:
553                travel_time = self.utils.activity_time(patient.hems_vehicle_type, 'time_to_hospital') * self.activity_duration_multiplier
554                self.add_patient_result_row(patient, travel_time, 'time_to_hospital')
555                self.add_patient_result_row(patient, "HEMS arrived destination", "queue")
556                yield self.env.timeout(travel_time)
557
558            if self.amb_data:
559                #self.debug('Ambulance travel time')
560                yield self.env.timeout(30)
561                self.debug('Ambulance at destination time')
562
563            patient.time_in_sim = self.env.now - patient_enters_sim
564
565            # Handover time ---------------
566
567            # Not currently available
568
569            # Clear time ------------
570
571            if hems_case and hems_avail:
572                clear_time = (
573                    self.utils.activity_time(patient.hems_vehicle_type, 'time_to_clear')
574                    * self.activity_duration_multiplier
575                    )
576                self.add_patient_result_row(patient, clear_time, 'time_to_clear')
577                self.add_patient_result_row(patient, "HEMS clear", "queue")
578                yield self.env.timeout(clear_time)
579
580            if self.amb_data:
581                #self.debug('Ambulance clear time')
582                yield self.env.timeout(60)
583
584            patient.time_in_sim = self.env.now - patient_enters_sim
585
586            if self.amb_data:
587                self.debug('Ambulance clear time')
588
589            #self.debug(f"Depart for patient {patient.id} on run {self.run_number}")
590
591            self.add_patient_result_row(patient, "depart", "arrival_departure")
592        finally:
593            # Always return the resource at the end of the patient journey.
594            if hems_res is not None:
595                self.hems_resources.return_resource(hems_res, secondary_hems_res)
596                self.debug(f"Attempting to return {hems_res} and {secondary_hems_res} to resource store")
597                self.add_patient_result_row(patient, hems_res.callsign, "resource_use_end")
598
599    def add_patient_result_row(self,
600                               patient: Patient,
601                               time_type: str,
602                               event_type: str,
603                               **kwargs) -> None :
604        """
605            Convenience function to create a row of data for the results table
606
607        """
608
609        results = {
610            "P_ID"              : patient.id,
611            "run_number"        : self.run_number,
612            "time_type"         : time_type,   # e.g. mobile, at scene, leaving scene etc.
613            "event_type"        : event_type,  # for animation: arrival_departure, queue, resource_use, resource_use_end
614            "timestamp"         : self.env.now,
615            "timestamp_dt"      : self.sim_start_date + timedelta(minutes=self.env.now),
616            "day"               : patient.day,
617            "hour"              : patient.hour,
618            "weekday"           : patient.weekday,
619            "month"             : patient.month,
620            "qtr"               : patient.qtr,
621            "registration"      : patient.hems_registration, # registration of allocated/attending vehicle
622            "callsign"          : patient.callsign, # callsign of allocated/attending vehicle
623            "callsign_group"    : patient.hems_callsign_group, # callsign group of allocated/attending vehicle
624            "vehicle_type"      : patient.hems_vehicle_type, # vehicle type (car/helicopter) of allocated/attending vehicle
625            "hems_res_category" : patient.hems_category,
626            "ampds_card"        : patient.ampds_card,
627            "age"               : patient.age,
628            "sex"               : patient.sex,
629            "care_cat"          : patient.hems_cc_or_ec,
630            "heli_benefit"      : patient.hems_helicopter_benefit,
631            "hems_result"       : patient.hems_result,
632            "outcome"           : patient.pt_outcome,
633            "hems_reg"          : patient.hems_registration
634        }
635
636        #self.debug(results)
637
638        # Add any additional items passed in **kwargs
639        for key, value in kwargs.items():
640             results[key] = value
641
642        if self.env.now >= self.warm_up_duration:
643            self.store_patient_results(results)
644
645
646    def store_patient_results(self, results: dict) -> None:
647        """
648            Adds a row of data to the Class' `result_df` dataframe
649        """
650
651        self.results_list.append(results)
652
653
654    def convert_results_to_df(self, results: dict) -> None:
655        self.results_df = pd.DataFrame(results)
656        try:
657            self.results_df.set_index("P_ID", inplace=True)
658        # TODO - Improve error handling here
659        except KeyError:
660            pass
661
662    def write_all_results(self) -> None:
663        """
664            Writes the content of `result_df` to a csv file
665        """
666        # https://stackoverflow.com/a/30991707/3650230
667
668        # Check if file exists...if it does, append data, otherwise create a new file with headers matching
669        # the column names of `results_df`
670        if not os.path.isfile(self.all_results_location):
671           self.results_df.to_csv(self.all_results_location, header='column_names')
672        else: # else it exists so append without writing the header
673            self.results_df.to_csv(self.all_results_location, mode='a', header=False)
674
675    def write_run_results(self) -> None:
676        """
677            Writes the content of `result_dfs` to csv files that contains only the results from the
678            single run
679
680            Note that this cannot be done in a similar manner to write_all_results due to the impacts of
681            the parallisation approach that is taken with joblib - depending on process timings it can lead to
682            not all
683        """
684
685        self.results_df.to_csv(f"{Utils.RESULTS_FOLDER}/output_run_{self.run_number}.csv", header='column_names')
686
687    def run(self) -> None:
688        """
689            Function to start the simulation.
690
691        """
692        self.debug(f"HEMS class initialised with the following: run {self.run_number}, duration {self.sim_duration}, warm-up {self.warm_up_duration}, start date {self.sim_start_date}, demand increase multiplier {self.demand_increase_percent}, activity duration multiplier {self.activity_duration_multiplier}")
693
694        # Start entity generators
695        self.env.process(self.generate_calls())
696
697        # Run simulation
698        self.env.run(until=(self.sim_duration + self.warm_up_duration))
699
700        # Convert results to a dataframe
701        self.convert_results_to_df(results=self.results_list)
702
703        # Write run results to file
704        self.write_all_results()
705        self.write_run_results()