Using cellular automata to simulate wildfire propagation and to assist in fire prevention and fighting

Cellular Automata have been successfully applied to simulate the propagation of wildfires with the aim of assisting fire managers in defining fire suppression tactics and in designing fire risk management policies. We present a Cellular Automata designed to simulate a severe wildfire episode that took place in Algarve (southern Portugal) in July 2012. During the episode almost 25 thousand hectares burned and there was an explosive stage between 25 and 33h after the onset. Results obtained show that the explosive stage is adequately modeled when introducing a non-local propagation rule where fire is allowed 5 to spread to the nearest and next nearest cells depending on wind speed. When the rule is introduced deviations in modeled time of burning from estimated time based on hotspots detected from satellite have a root mean square difference of 8.7 hours for a simulation period of 46h (less than 20%). The simulated pattern of probabilities of burning as estimated from an ensemble of 100 simulations show a marked decrease out of the limits of the observed scar, indicating that the model represents an added value for fire fighting in what respects to the choice of locations to allocate resources for fire combat. 10

steep terrain. It was the largest wildfire in Portugal in 2012, contributing to more than 22% of the total amount of 110, 232 ha of burned area (ICNF, 2012) in that year.
The terrain is prevailingly steep, with slopes of 20% located in the higher altitude region in the northern part of the Tavira municipality with hilltops reaching up to 541 meters. The altitude and slope decrease towards the southeast area of Tavira and southwest area of the São Brás de Alportel municipalities, having slopes between 0 and 20% and lower altitudes reaching sea 5 level at several locations (Viegas et al., 2012).
Since 2012 was a year of extreme drought, the meteorological background conditions were very prone to the occurrence of large fire events (Trigo et al., 2013). The region of Tavira is characterized by Mediterranean climate, the maximum monthly temperature in August ranging from 25 o C to more than 30 o C, with maximum absolute temperatures around 39 o C, and mean relative humidity below 65% (ANPC, 2012). In 2012, the precipitation in Tavira was 45% below the normal record and the 10 study area had a soil water content value below 10% at the time of the fire (Viegas et al., 2012). The wildfire propensity was further aggravated by above average precipitation in 2010 and 2011, that favored vegetation growth and fuel build up.
The fire propagated in two distinct phases. In the first stage, from 13:00 UTC on July 18 to 17:00 UTC on July 19, the fire burned about 5, 000 ha, representing one fifth of the total burned area. In this phase, the fire advanced through rugged terrain, and wind direction was highly variable, causing frequent shifts in the direction of maximum spread, which was mainly 5 towards south/southeast until it reached the Leiteijo stream, where it gained speed under the influence of topography. Then, in a transition stage, around 16:30 UTC on July 19, the fire started spreading through steep slopes along the Odeleite stream.
Spotting occurred up to hundreds of meters due to low fuel moisture, and multiple spot fires were recorded. Fire suppression was difficult due to the steepness of terrain and frequent wind direction changes, and operations were focused on life and property salvation.

10
In the second stage from 17:00 UTC to 24:00 UTC on July 19 the fire turned into a major conflagration, greatly increasing its propagation speed and burning about 20, 000 ha in 7 hours.
When the fire reached the Odeleite stream it became orographically channeled, as an increase in wind speed led to fast and intense fire growth towards south, where copious amounts of fuel loading were present. The fire split into two advanced sections heading west and east to the São Brás de Alportel and the Tavira municipalities, with a 10 km wide fire front. In 15 addition, spotting now created new fires up to two kilometers ahead of the fire front. All these factors allowed rapid propagation of the fire front while turning suppression extremely difficult.

Input data
A study area of 30 km x 30 km was defined centered on the burned area ( Figure 1) and fine-scaled raster data from various sources were collected and pre-processed in a common format suitable as input for the wildfire simulations. Data include the 20 ignition points, the start and end times of the fire event, the fire perimeters, the burned areas, the surface wind speed and direction, the topography, and information about the landcover (vegetation type, vegetation density, areas burnt in previous wildfires, waterlines and roads).
Patch-slope information was derived from elevation data as obtained from the Digital Elevation Model provided by the Shuttle Radar Topography Mission (Farr et al., 2007).

25
Hourly wind data were obtained from a regional weather simulation performed with the Weather Research and Forecast model (WRF), version 3.1.1 (Skamarock et al., 2008). The quality of the simulation was previously assessed for wind (Cardoso et al., 2012;Soares et al., 2014). WindNinja (version 2.1.3) (Forthofer, 2007) was then used to spatially model the hourly wind input data taking into account the interaction with topography. The temporal behavior of the wind field was then validated against the information contained in the report by the Portuguese Authority for Nature and Forest (ICNF) (ICNF, 2012).  without vegetation, areas with sparse, normal and dense vegetation ( Figure 2, bottom left panel). As described in section 3.1, for the different categories of vegetation type and density, values of the associated probability factors, respectively p veg and p den , were empirically assigned or taken from literature (Alexandridis et al., 2008). Assigned values are listed in Table 1.
Roads and waterlines inside the simulation area ( Figure 2, top right panels) were also included in the model by assigning low values to both probability factors p veg and p dens , with p veg = p dens . Primary, secondary and tertiary roads were assigned the  Active fire data as identified from satellite were used for the quality assessment of the CA model simulations by evaluating temporal and spatial discrepancies between active fire observations and simulated fire growth. For this purpose, we used the MODIS (MODerate Resolution Imaging Spectroradiometer) active fire product that provides hotspots detected at 1km × 1km pixel resolution, at the time of the satellite overpass. The MODIS sensor on the Terra and Aqua satellites supply daytime and nighttime observations at four nominal acquisition times, thus providing information about the geographical location, date, and 5 time of the detected active fires (Giglio et al., 2003).
For each satellite overpass totally or partially covering the total burned area by the Tavira fire, we used the centroids of the active fire footprints (  with a probability p burn which is a function of the variables that affect fire spread.

10
In its basic formulation, probability p burn of a given cell depends: 1) on a constant reference probability that a cell in the neighborhood of a burning cell (containing a given type of vegetation and density) starts burning at the next time step under no wind and flat terrain, p 0 ; 2) on the vegetation type, p veg , and vegetation density, p den ; 3) on topography, p s ; and 4) on wind fields, p w , as follows:

15
As described in section 2.2, in order to account for the effect of vegetation, both type and density were stratified into discrete classes, and for each class a constant probability factor was assigned as specified in Table 1.
The effect of the wind is modeled as: where c 1 and c 2 are adjustable coefficients, V is the wind speed, and θ is the angle between the wind direction and the fire 20 propagation direction. As expected, p w increases when wind and fire directions are aligned.
The probability factor that models the effect of the terrain elevation is given by: where θ s is the slope angle of the terrain and a s is a coefficient that can be adjusted from experimental data. Slope angle θ s was derived from elevation data, E, according to where D is equal to the size L of the square cell when the two neighboring cells are adjacent or to √ 2L when the two cells are diagonal. As expected this topography effect is higher when the fire spreads uphill. The model with the new propagation rule N 2 will be hereafter referred to as the modified model.

Simulations
The landscape was discretized into square cells with size of 100 m and the model free parameters were set according to the model parameters along the boundary of the final observed scar. However, this setting along the scar boundary is not an artificial device since it reflects the known a posteriori fact that the shape of the scar resulted from effective fire combat in these locations.
Two different ensembles of 100 simulations were generated, one with the baseline fire spread model and the other with the modified model. Results obtained at four selected stages of the fire are displayed in Figure 6. When using the baseline rule 5 ( Figure 6, left panels), and excepting for the slot at 25h (after ignition) where there is a fair agreement between the simulated burned area and the front lines of the fire as indicated by the hotspots identified by satellite, the simulated burning is well behind the fire front, an indication that the modeled propagation of the fire is too slow. A strong contrast is observed when using the modified model ( Figure 6, right panels). In this case, the modeled burned areas spread much closer to the fire front as defined by the hotspots. The exception is the slot at 25h (after ignition), where the modeled propagation of fire is faster than 10 the one suggested by the location of the hotspots. On the other hand, it is worth emphasizing that the explosive behavior of fire between slots at 25h and at 33h is very well simulated when using the new wind propagation rule.

Quality assessment
Burned area in each one of the two ensembles was identified by assuming that a given pixel is a burned one when the modeled probability that it burned is larger than 0.2. Each pixel identified as burned was assigned the respective time step as an indi-   . Each triplet of columns corresponds to the burned cells identified, respectively, in the intervals [0,6[, [6,12[, [12,18[, [18,24[, [24,30[, [30,36[ and [36,42[ hours. restricting the burned area to cells with a burning probability larger than 80% (Figure 9, right panel). Since fire containment was mainly due to actions by firemen along the perimeter, results indicate that unconstrained simulations represent a very useful tool to assist decision makers during a fire event, by providing indications about locations of low burning probability to be selected as appropriate to allocate resources for fire combat.