Volcanic ballistic projectile deposition from a continuously erupting volcano: Yasur Volcano, Vanuatu

Volcanic ballistics are the main hazard to life and infrastructure from Strombolian eruptions, which are a tourist drawcard, exposing people to this hazard. Most research to date has been to understand this style of eruption and how ballistics form and travel. However, little focus has been placed on how ballistics are distributed within ballistic fields or the inclusion of this data into hazard and risk assessments. In this study we used a UAV to image the ballistic field, and cameras to record eruptions at Yasur Volcano, Vanuatu from 28 July – 2 August and 17 – 19 October 2016. We present the mapped distributions from the two trips, how the field changes with distance and direction from the vent, and how eruption dynamics influence these changes. Our evidence for directionality results in considerable variation in summit ballistic hazard and is an important consideration for ballistic hazard and risk assessments.

Risk reduction and mitigation measures (such as land use planning, exclusion zones, protective-wear, shelters, communication and education products) are Volcanic ballistic projectile deposition from Yasur Volcano Fitzgerald et al., 2020 key to reducing risk to visitors, though to be most effective they must be supported by robust hazard and risk assessment . VBP hazard assessments determine the likelihood of VBP-producing eruptions and the likely impacted areas [Alatorre-Ibargüengoitia et al. 2012; Thouret et al. 2000]. VBP risk assessments go one step further and determine the likelihood of specific consequences occurring (e.g. fatalities, casualties, damage) due to the exposure to VBPs [Blong 1996]. Increasing tourism and other activities on active volcanic summit areas leads to greater exposure to the hazard, and there are increasing regulatory requirements and societal expectations to inform users about the risk to which they may be exposed [Deligne et al. 2018;Jolly et al. 2014]. These have driven the requirement to more accurately assess risk so that it may be evaluated and treated within a risk management framework. Yet, our understanding of VBP hazard is relatively limited. A lack of detailed mapping of how VBPs are distributed in space (VBP fields) during explosive eruptions is notable with only a few fully mapped VBP fields reported [e.g. Fitzgerald et al. 2014;Gurioli et al. 2013;Swanson et al. 2012]. The main reasons for this are fieldwork time constraints (mapping every individual VBP on foot is time-consuming), the high risk involved in mapping proximal vent areas especially during unrest or eruption phases, the resolution limitations of readily available remote sensing imagery, and the limited geological preservation of VBP fields due to cover from ash and lapilli and subsequent erosion. Instead, published VBP distribution maps often only show the maximum travel distance or extent of the field [Minakami 1942;Nairn and Self 1978;Robertson et al. 1998;Yamagishi and Feebrey 1994]. This does not give a complete understanding of the hazard and can lead to inaccurate estimations of risk, particularly without the knowledge of spatial density. Sufficiently high resolution aerial photos help get a wider understanding of field characteristics, making mapping the field possible and allowing for targeted field investigations to supplement and check the accuracy of photo mapped data, maintain detail and reduce uncertainty in identifying VBPs. As such, aerial photos have been used to map VBP fields Kaneko et al. 2016]. Previously this has been accomplished using planes and helicopters. However, the use of these aircraft is expensive and can be a limitation for scientists with limited resources.
Drones, or remotely piloted aircraft, are a tool that can be used to assess geological hazards [Gomez and Purdie 2016] while reducing the risk to scientists. They offer particular capabilities that lend themselves to volcano observation and monitoring. This includes: (1) access to dangerous or hard to reach areas while keeping scientists at a safer distance from the hazard; (2) portability; (3) ability to produce higher resolution outputs compared with those taken by larger aircraft such as helicopters, by being able to get closer to fea-tures; and (4) relatively low cost. For example, drones have been used in volcanic research to thermally map a geothermal valley [Harvey et al. 2016], observe crater activities in active volcanoes [Jordan 2015;Patterson et al. 2005], map active lava flows and model future flow paths [Turner et al. 2017], determine the likely location of magma and dykes using aeromagnetic surveys [Kaneko et al. 2011], measure the gas composition [Shinohara 2013] and flux from a volcanic plume [McGonigle et al. 2008], 3D image a volcanic plume [Gomez and Kennedy 2018], and survey land changes after volcanic eruptions [Nakano et al. 2014]. Two MSc theses have assessed the methodology of using drones to map VBP fields [Gates 2018;Pitchika 2017], though to the best of our knowledge, no-one has utilised a drone to map VBP fields and assess VBP hazard.
In this paper we present the results of two field campaigns from 28 July to the 2 August and 17 to the 19 October 2016 to map the VBP field and understand the VBP hazard at Yasur. The study objectives were to: (1) use drones and photogrammetry to determine the distribution and varying spatial intensity of the hazard (both the number and size of VBPs) constrained between the two flights; and (2) evaluate explosion frequency and dynamics from video footage to understand how this influences the VBP hazard distribution at Yasur. The results presented here will add to the limited data available on VBP distributions from Strombolian eruptions, improve our knowledge of how VBP fields vary with distance and direction from vent and over time, and improve understanding on the causes of heterogeneity in hazard around a vent as well as within a field. The paper is structured into four main sections. We describe our methodology, followed by results which are subdivided into Spatial distribution containing mapping results and Eruption dynamics containing video and observation derived results. Then we discuss the relationship between VBP size and distance from vent, how spatial density changes with direction around the vents, differences between the two data sets, and finally, we describe the hazard/risk implications of our results.
(at the summit; 169°26' 43.06" E 19°31' 36.74" S), situated within the >20 ka Siwi caldera at the north-western edge of the Yenkahe resurgent dome [Battaglia et al. 2016;Brothelande et al. 2016;Gaudin et al. 2014;Merle et al. 2013;Métrich et al. 2011]. At its base, the cone is 1500 m in diameter, while the summit is irregular in height-highest in the west and south and decreasing in elevation towards the north and east. Present-day Yasur is composed of a 660 m diameter crater subdivided into two sub-craters with three main vent areas (A and B in the southern

Page 184
Volcanica 3(2): 183 -204. doi: 1 .3 9 9/vol. 3. 2.1832 4 sub-crater, and C in the northern sub-crater) [Battaglia et al. 2016;Kremers et al. 2012;Meier et al. 2016;Nabyl et al. 1997;Oppenheimer et al. 2006;Spina et al. 2015, Figure 1]. It is evident that the positions of the vents migrate over time and the two southern vent areas (A and B) are alternately named by different authors. Drone footage from this study revealed a change in the number and position daily of vents in both the North and South craters. Thus, we have decided to adopt the Jolly et al. [2017] naming convention of simply North and South Crater. Firth et al. [2014] identified 3 stages of eruption from tephra stratigraphy and literature review: pre-700 AD, a variable frequency and intensity eruptive episode; from~700 AD to 1270 AD a higher magnitude lower frequency episode; and from 1270 AD to the present a persistently active low magnitude high frequency episode. This most recent eruptive stage has been predominantly Strombolian with intermittent Vulcanian activity [Bani and Lardy 2007;Firth et al. 2014;Nabyl et al. 1997]. Explosion frequency has been variable over the current high frequency stage. Meier et al. [2016] reported explosions several times a minute over a 15day period in 2008, while Marchetti et al. [2013] and Battaglia et al. [2016] reported a frequency of at least one explosion per minute (over 4 days in July -August 2008 and the entire year of 2008, respectively). Kremers et al. [2013] report similar results of explosion recurrence under one minute (from Vent Area A) over two weeks in August-September 2008. However, Vent Area B frequency ranged from minutes to days and Vent Area C had strong explosions every 10 minutes. Vent Area B explosion recurrence was reported to be on the minutes end of the range by Oppenheimer et al. [2006] with 13 explosions occurring over a 15-minute period in January 2005. Over a two-hour period in September 2011, Bani et al. [2013] recorded 200 explosions from Vent Area B using a thermal infrared thermometer. Infrasound recordings reported by Pichon [2005] showed 20-1300 explosions per day.

Methodology
VBP fields can cover large areas (e.g. 6 km 2 in area recorded by Fitzgerald et al. [2014]), making them difficult to map in detail due to the time needed to record all pertinent information from each VBP (e.g. dimensions, lithology, crater formation, impact angle) and the sheer number of VBPs that are contained in the field . Aerial photography and photogrammetry allow for mapping to be done remotely, in detail over the whole area and in a much shorter time frame. In our study, a DJI Phantom 3 drone was flown over the crater and flanks of Yasur Volcano, capturing 3,863 images and covering an area of 0.82 km 2 from NE to S ( Figure 2) at heights above the ground up to 60 m with a 12-megapixel camera. Time and equipment re-strictions meant limiting flying to a 135°swath of the volcano, with the NE-S area chosen due to it encompassing the viewing areas, track, and carpark. Eight orthophotographs and accompanying DEMs were created from these images using Structure from Motion (SfM) software (PhotoScan) in high enough resolution to map individual VBPs. Flights were flown on autopilot (using the Map Pilot application) with flight paths created on an iPad in the field. Using Map Pilot ensured that photos were taken with sufficient overlap for SfM to be used and that no areas were missed.
Geospatial data was collected using a Trimble R8 Real-Time Kinematic Global Navigation Satellite System (RTK GNSS), with a base station set up at the edge of the car park ( Figure 2) and linked to the drone imagery through ground control points taken at ground targets visible in the aerial images.
Although RTK GNSS produces a typical location error in X, Y < 2-5 cm and < 10 cm in Z compared to single GNSS and GPS units, the SfM photogrammetric method relies on statistical approximation of the camera characteristics [Westoby et al. 2012] and on the resolution of the imagery [Gomez 2014]. Therefore, even with GNSS RTK tie points, error can spread in between the ground control points and precisely aligning the produced orthophotographs and DEMs from the two trips can be challenging once imported into the ArcGIS environment. To make visual comparisons between the two orthophoto datasets easier and rectify the misalignment, the orthophotos and DEM sets were stacked in the GIS software ENVI and exported as combined TIFFs. The stacks were then imported into ArcGIS where they were split into smaller areas (between 4 and 20 depending on the size of the area the original orthophoto covered) and georeferenced to each other. Splitting the layers resulted in more accurate georeferencing as it reduced the warping produced when the larger datasets were used.
As the area captured in the orthophotos was too large to map all VBPs in the available timeframe, twentythree 20ˆ20 m squares (an area of 400 m 2 for each square) were visually mapped from the orthophotos from the July-August trip and thirty from the October trip (as the latter covered a larger area). The squares were positioned at 100 m increments and along transects 22.5°apart radiating from a central point between the two craters (hereafter called the crater mid-point) to capture as much of the orthophoto area and variation in VBP distribution around the volcano as possible ( Figure 2). A central point between the two craters was chosen as it was not possible to distinguish which VBP came from which vent. Spatial density contours were then applied to the mapped distributions using the spline tool in ArcGIS and best fit estimates. In addition to the location being recorded, VBP dimension (taken from the largest axis) and whether an impact crater has been produced and its dimensions was also noted. Identification of VBPs was based on mul- tiple factors though not all were required for identification. Factors included (1) morphology (a small topographic increase caused by VBP sitting on the surface, distinct VBP rounding or angularity), (2) colour differences (darker or lighter clasts compared to the ground surface), (3) round depressions around or near a VBP indicating impact with the ground, and (4) disturbed material around a VBP indicating an impact. Figure 3 illustrates this process, showing mapping squares preand post-VBP identification. Comparison of the July-August and October orthophotos helped to identify new VBPs that had landed post July-August data collection and also confirmed if something was a VBP (i.e. seeing it in a slightly different light or angle in a different orthophoto could help to confirm it was a VBP and not a lighter patch of ground or autobreccia from the old lava flows on the flanks). Unfortunately due to time constraints in the field, we could not ground-truth the mapping squares to assess for resolution and mapping error.

Presses universitaires de �rasbourg
Along with the drone flights, three GoPro Hero 3c ameras with orthorectified lenses were set up around the crater rim to record eruption dynamics and frequency for 10 hours over 30 July to 1 August 2016 (Figure 2). Each camera could see into the crater though not to the vents as crater walls blocked these from view. When possible, we assessed the size and directionality of the explosions from the GoPro video footage. Multiple locations around the crater and concurrent filming from the cameras meant that directionality could be checked at another angle, reducing stereological effects. Directionality was classified (by sight) on the angle of the eruption jet in relation to landmarks around the crater (i.e. the jet was directed towards the viewing areas to the east and south). Explosions were sized based on the height that VBPs or a sustained tephra plume reached (small: a third of the way up the crater; moderate: between a third of the crater height and the crater rim; large: above the crater rim). An anemometer (Kestrel 5500 L) was also installed at the southern viewing area ( Figure 1C) to record wind speed and direction every 1 to 5 minutes over 29 July to 1 August. Additionally, systematic observational logs were collected for one hour each day from the seating/viewing area in the SSE from 5 September to 16 October 2016 in between the two imaging trips by B. Simons. This record provided the source crater that any bombs were coming from, which direction they were landing and if they landed outside the crater rim.
As part of concurrent studies [Iezzi et al. 2019;Jolly et al. 2017;Matoza et al. 2017], a temporary infrasound and seismic network was set up around the cone. Two stations from this network (YIB22 and YS01), lo-cated~700 m from the vent area, are used in this pa- per to support changes in eruptive activity observed in video footage from the same time period. These stations are also used to place the explosion events selected for analysis in this study in context with the continuous ongoing activity. As a proxy for event size, we use peak event amplitudes from unfiltered infrasound waveforms from YIB22 in a time-window froḿ 5 to`10 s around an automatic network-coincident STA/LTA trigger. The STA/LTA detection was performed with 0.1-50 Hz filtered data using an STA length of 0.5 s, an LTA length of 40 s, and a coincident STA/LTA threshold of 3 on at least 2 stations of a 3element array. Real-Time Infrasonic Amplitude (RIAM) and Real-Time Seismic Amplitude (RSAM) were also estimated as the 10-minute mean of absolute amplitude of 0.1-50 Hz filtered data [Endo and Murray 1991].

Results and analysis
This section is divided into two subsections: Spatial distribution, which reports results from the aerial mapping; and Eruption dynamics which conveys the results from the GoPro footage and visual observations.

Spatial distribution
We present two VBP distributions from Yasur volcano. The first is of VBPs deposited between 28 July and 19 October 2016 (hereafter referred to as the two-month distribution), by mapping VBPs that appear only in the October images and that are not present in the July-August images. The second was the distribution of all observable VBPs in October, which includes the longerterm distribution of VBPs deposited potentially over months to years.
From the two-month distribution, 1,550 new VBPs were identified from 23 mapping squares (each square 400 m 2 in area). VBP diameters ranged from 5 cm (the minimum size that could be distinguished) to 304 cm with a mean of 43 cm. The number of VBPs per m 2 drops rapidly as distance from the crater increases, from 182ˆ10´2 m´2 at 200 m from the crater midpoint (defined in section 2) to an average of 6.691 0´2 m´2 at a distance of 400 m ( Figure 4A, 5D). Greater than 500 m from the crater mid-point the number of VBPs decreases to between 0 and 0.5ˆ10´2 VBPs per m 2 in most azimuths. The distribution of VBPs also shows directionality in the S and SSE azimuths where a greater number of VBPs were deposited between 300 and 500 m from the crater mid-point than in other azimuths (VBP densities ranging 124-3.41 0´2 m´2 compared to 23 -0ˆ10´2 m´2, respectively) ( Figure 4A). The E and ESE azimuth also have more VBPs at 400 m from the crater mid-point than the SE. Conversely, Figure 5A shows the SE having the largest median VBP diameter at 300 m that decreases away towards the ESE and SSE, creating a lobe that might indicate directionality. Median diameter was used as there are large VBP outliers that skew the mean value. Analysis of the size distribution also revealed an increase in median diameter at 500 m following a trend of decreasing size in the SSE azimuth.
The total October distribution shows similar patterns as the two-month distribution. From 30 mapping squares (400 m 2 areas) 5481 VBPs were mapped. VBPs ranged from 3 to 353 cm in diameter with an average of 34 cm. VBP spatial density decreases with distance from the crater from an average of 138.5ˆ10´2 m´2 at 200 m from the crater mid-point to 4ˆ10´2 m´2 at 700 m, though not as rapidly as seen in the twomonth distribution ( Figure 4B and 5D). When analysed in azimuths we see a general decrease in density along all azimuths, though at some point along all but the E and ENE azimuths, an increase in density is observed. For example, in the S azimuth, density decreases from 72.25ˆ10´2 m´2 at 400 m from the crater mid-point to 16.5ˆ10´2 m´2 at 500 m distance. Density then increases to 32.75ˆ10´2 m´2 at 600 m from mid-point and then decreases again with distance. Similar to the two-month distribution, VBP deposition is greater towards the S-SSE and E-ESE in contrast to the SE at most distances except at 400 m from mid-point where Presses universitaires de �rasbourg it is at similar densities to the S and SSE ( Figure 4B).
Analysing the October size distribution at varying distances and azimuths from the crater shows higher median VBP diameters to the SE at 300 and 500 m from crater mid-point that decrease away towards the ESE and SSE ( Figure 5B). At 400 m from the crater mid-point, the ENE has the largest median diameter that decreases towards the E and NE, potentially indicating two lobes towards the SE and ENE ( Figure 5B). The October distribution shows a general decrease of median VBP size with distance from crater mid-point until 600 m where VBP size begins to increase (Figure 5C). Increasing diameters following a general decreasing trend is also seen in the S, SSE, and SE azimuths ( Figure 5B).

Eruption dynamics
Over 10 hours (between 30 July and 1 August) of GoPro footage was collected, capturing 758 explosions. On average, 68 explosions occurred per hour (hr´1), with 42 hr´1 from South Crater and 27 hr´1 from North Crater. Each explosion was classified small, moderate or large, the vent origin noted (Table 1), and what style of eruption it was. As described earlier, the larger the explosion, the higher the VBPs are ejected vertically. Subsequently, a larger explosion has the ability to eject VBPs further from the vent area. Eruption styles included ballistic, ballistic with a tephra plume, tephra plume, and gas. To compare eruption frequency over the three days we normalised the number of small, moderate and large explosions each day to a 13-minute time period as this was the longest period of time explosions were Presses universitaires de �rasbourg recorded on 1 August. Therefore, the normalised eruption frequency provides the average number of explosions that occur over 13 minutes on that day. A marked increase in the frequency of large events is noted on 1 August (on average 13 large explosions per 13-minute period) compared to the previous days (on average 1 large explosion per 13-minute period on 30 July and 5 large explosions per 13-minute period on 31 July). Additionally, an increase in the number of explosions from the south crater vs the northern crater was noted. The frequency increased from 5 explosions per 13 minutes (30 July) and 14 explosions per 13 minutes (31 July) to 20 explosions per 13 minutes on 1 August. Seismic and infrasound data over this time also show an increase in explosivity ( Figure 6). This can be more clearly seen in the infrasound record ( Figure 6A and 6B) where on day 5 (1 August) there is a noticeable increase in peak event and total (RIAM) waveform amplitude. Highlighted in blue in Figure 6A are the video observed explosions from Table 1. On average small explosions produce a peak pressure of 6.71 Pa, increasing to 13.6 Pa for moderate explosions and 27.7 Pa for large explosions. Eruptive activity therefore is not steady and can fluctuate over hours.
From our analysis of the GoPro footage over 30 July to 1 August period, most explosions were vertically directed (40 %). However, when the jet was angled, a SE directionality was the favoured direction from South Crater and S from North Crater (Figure 7). The SE directionality from South Crater is also evident in the acoustic radiation recordings over the same time pe-riod, presented in Jolly et al. [2017]. No relationship between explosivity and directionality was observed with between one third and one half of all explosion sizes (that directionality could be ascertained) having vertical jets and the rest directed towards an azimuth. Between the two trips from 5 September to the 16 October, observations were made daily for 1 hour a day of where VBPs were landing (rather than orientation of the jet as for the GoPro footage), indicating directionality favouring the south. This period of time is when most of the VBP deposition occurred that is recorded in the two-month mapped distribution and thus it was important to get an idea of the eruption directionality within this time. Figure 7C shows the azimuths where VBPs ejected from South Crater landed, with 39 % landing to the south. Only two observations were made of VBP directionality from North Crater. On both occasions they were directed to the south. During the second trip from 16-20 October, general visual observations of eruption dynamics showed directionality to the south and east from South Crater and to the east from North Crater [Gomez and Kennedy 2018].
We also noted a changing directionality throughout individual eruption events on some occasions. For example, VBPs would be ejected toward the north at the start of the explosion and subsequently move towards the south in a spraying motion. This was also noted by Gaudin et al. [2014] in Strombolian explosions, where the mean ejection angle shifted up to 40°in a single explosion. They theorise that this is due to the changing location of the bubble rupture point as the rupture area Presses universitaires de �rasbourg continues to open and increase in size. These instances have not been included in the directionality data presented here for simplicity, however, it is important to recognise that directionality is not fixed and can change even during one eruption event.
Not all explosions at Yasur have a VBP component, with 18 % of the explosions recorded on the GoPro videos having no observable VBP component at all. This is an important consideration when using eruption frequency to assess VBP hazard.
Summarising our key results, we find that spatial density of VBPs decreases with distance from the crater, with this occurring more sharply in the two-month distribution than the October distribution, and that median VBP diameter decreases with distance from crater.
Directionality is evident in the spatial density mapping (S-SSE with a minor increase in density to the E), VBP size distribution analysis (SE and minor ENE) and video analysis (SE from South Crater and S from North Crater) and the daily visual observations, though they do not all agree. Video footage showed an increase in the frequency of large explosions from the 30 July to the 1 August, with large explosions capable of ejecting VBPs above the crater rim and creating a hazard for visitors. This is also seen in the infrasound and seismic record where amplitude increased on 1 August. Table 1 -Number and normalised frequency of small, moderate and large events, and the number and normalised frequency of explosions from each crater recorded on GoPros between 30 July and 1 August. Explosion frequency has been normalised to a 13 minute time period (the length of time explosions were recorded for on 1 August) as each day the length of time explosions were recorded for differed. Peak pressure recorded on infrasound station YIB22 was found for each video recorded explosion and the average and standard deviation of all explosions categorised into small, moderate or large was calculated. : Infrasound and seismic amplitude metrics from a temporary co-located infrasound station (YIB22) and broadband seismometer (YS01) approximately 700 m from the summit vents. A shows peak event amplitudes for unfiltered infrasound data at YIB22. Blue symbols in A show infrasound amplitudes at YIB22 corresponding to the explosions analysed in video data in this study (see Table 1). Vertical dashed lines represent times of network completion for coincident triggers. B shows 10-minute Real-Time Infrasonic Amplitude (RIAM) and C shows Real-Time Seismic Amplitude (RSAM) values at YIB22 and YS01, respectively.

Discussion
While VBPs are primarily transported through the air on ballistic trajectories, two other factors can cause further transport after impact. The first is bouncing or rolling. This has been observed on Ngauruhoe [Allen 1948], Tungurahua [Bernard 2018], Kīlauea [Swanson et al. 2012] and Stromboli [Taddeucci et al. 2017] and can result in the VBP coming to a complete stop tens to hundreds of metres from its original landing position [Bernard 2018;Pistolesi et al. 2008;Rosi et al. 2006]. This typically occurs on the slopes of a volcano, observed on slopes >15-20°at Alaid Volcano [Steinberg and Lorenz 1983] and between 31°and 38°o n Etna ]. Secondly, VBPs can fragment on impact with the ground and the ejecta from this can travel metres away from the original landing/fragmentation site (>5 m by Taddeucci et al. [2017]; tens of metres at Stromboli by Rosi et al. [2006]; 28 m at Shinmoedake by Maeno et al. [2013]; 100 m at Ngauruhoe by Nairn and Self [1978]). Therefore, mapped VBP distributions will often reflect all transport mechanisms as it is often difficult to determine if these additional transport mechanisms have occurred, especially in areas of high VBP spatial densities.
The VBP distributions mapped from the aerial pho-  Table 1 as the plots only show directed explosions, not those vertically directed or where direction could not be ascertained.

Presses universitaires de �rasbourg
[C] shows the directionality from South Crater seen from the observational logs from 5 September to 16 October 2016. The rose diagram bin thicknesses and gaps between bins are not representative of the resolution of the angle of explosion direction.
tos are the products of both the North and South Craters. The frequency of explosions from both craters, the lack of temporal sampling and the varying explosion directionalities make it difficult to assign specific vents to a distribution. However, it is likely that the higher densities in the S-SSE come from South Crater and those in the E from the North Crater. This deduction is based on the dominant explosion directionalities observed from video records and the observational logs and the greater distance VBPs from North Crater need to travel when directed toward the S-SSE.

Decreasing median VBP size with distance from the crater
The relationship between VBP size and distance travelled from source has been discussed in the literature with both increasing and decreasing size with distance observed. Taddeucci et al. [2017] summarise this literature to say that an increase in VBP size with distance is observed in many ash and block rich Vulcanian eruptions [Druitt et al. 2002;Fitzgerald et al. 2014;Minakami et al. 1969;Self et al. 1980;Steinberg 1974;Steinberg 1976] and a decrease in size with distance is usually observed in phreatomagmatic eruptions [Lorenz 1970;Novellis and Luongo 2006;Self et al. 1980;Waitt et al. 1995 However, VBP size-distance relationships have not been discussed in relation to Strombolian eruptions. As Strombolian explosions tend to be frequent (up to 9 hr´1 at Stromboli [Harris and Ripepe 2007]) and typically deposit bombs closer to the vent, the risk of being impacted while conducting fieldwork is high and is likely the main reason the size-distance relationship has not been investigated. A literature search resulted in two references where a relationship was found. Gurioli et al. [2013] note that in the 21 January 2010 eruption at Stromboli, thermal video shows 'a leading spray of smaller bombs, quickly followed by an emission of larger bombs that attained lower heights and fell closer to the vent than those of the first spray'. Andronico and Pistolesi [2010] report pumiceous clast sizes between 7 and 20 cm in the summit area of Stromboli (~790 m.a.s.l.) that reduce to an average long axis of the five largest clasts of 11 cm at 650 m.a.s.l. from the 24 th November 2009 eruption. Both papers show a reverse distribution. Bombrun et al. [2015] also recorded lower initial ejection velocities of larger projectiles compared to smaller projectiles at Stromboli volcano. Strombolian eruptions are driven by gas slugs, providing the Presses universitaires de �rasbourg gas stream needed to reduce drag on smaller projectiles [Harris et al. 2012;Taddeucci et al. 2015]. This is likely the reason for the reverse distribution described in this paper and seen in the aforementioned publications.
Rolling and bouncing of VBPs on slopes post-impact may be a factor influencing the VBP size distribution. The base of the slope may have a greater proportion of larger VBPs than on the slope itself due to the momentum of rolling and bouncing VBPs. This is observed in rockfalls and talus slopes [Evans and Hungr 1993]. However, we observe a more complex relationship. Both increases and decreases in median VBP diameter are observed at the base of the cone compared to those on the cone in Figure 5A and 5B. Increases are noted in the SSE and ENE azimuths in the October distribution and in the SSE and SE in the two-month distribution. Decreases in median diameter are noted in all other azimuths. Whilst post-impact rolling and bouncing might influence the size distribution of VBPs on slopes causing larger VBPs to accumulate at the base of slopes, it is not clearly established within this data. Though increases in the median VBP diameter at the base of slopes were observed in some azimuths, this does not change the general decreasing size with distance relationship observed in the mapped distribution over the whole mapped area (from the vent area to 700 m distance). Nevertheless, larger VBPs rolling downhill on the slope could bias the average diameter towards higher values than that expected if only considering the point of first impact.
There is potential for VBP fragmentation to occur on impact with the ground which may affect the recorded size distribution or what is thought to be the size distribution of original VBPs. However, due to the high spatial density of VBPs on the cone it is near impossible to determine whether they have fragmented or not and therefore how fragmentation contributes to the size distribution. This may be a future research avenue in which the contribution of fragmentation to the size distribution is quantified.

Changing spatial density with direction around the crater
We see a change in spatial density from the south to the east, with higher densities in the S-SSE and ESE-E than the SE. Four main factors might influence where VBPs cluster or appear to cluster: (1) topography [Kilgour et al. 2019], (2) preservation of the deposit (3) VBP fragmentation and (4) eruption directionality [Gurioli et al. 2013;Kilgour et al. 2010]. A topographic high or steep sided crater wall may block VBP deposition as VBPs may not be able to reach heights capable of landing on the other side of the topography or angled VBPs may get blocked by steep sided crater walls. There does not appear to be any topographic shadowing affecting VBP deposition in the SE from South Crater or North Crater. The height of the crater rim is the same in the SE as the ESE (~320 m), though is 10 m and 35 m higher in the SSE and S respectively when measured from the crater mid-point. Therefore, if the SE were shadowed we should see similar densities to the SE in the ESE and lower densities in the SSE and S. Analysing the ability of VBPs to be ejected from the vents uninhibited by crater walls reveals that VBPs ejected at angles ě34°f rom horizontal could escape the crater from the ESE to the SSE (Figure 8). The angles were measured from the middle of the North and South Craters respectively, rather than the all-encompassing crater mid-point to be as accurate as possible. In Strombolian explosions ejection angles are typically close to vertical (within 5°from vertical at Stromboli ; around 72°a t Mt Etna [Gouhier and Donnadieu 2011]; 70-85°at Stromboli [Pistolesi et al. 2008]) with 80-90 % of particles ejected within 20°of the ejection axis Gouhier and Donnadieu 2010]. If we use the lower end of the reported average ejection angles (72°) and add on 20°for the dispersion cone, all VBPs would successfully clear the ESE-SSE crater walls.
Topography may also impact the preservation of the original spatial density (landing positions). Steeper slopes may cause the occurrence of bouncing or rolling downhill post-impact Steinberg and Lorenz 1983] which may present as greater spatial densities at the base of slopes than on the slopes themselves. Analysis of this relationship shows mixed results. In the October distribution ( Figure 9B), we do not see higher densities at the base of the cone in the S and SE at 500 m distance (16.5 and 3.75/m 2 respectively) compared to on the flank at 400 m (72.25 and 69.5/m 2 respectively) where there are greater slope angles. It is also not seen in the E at 250 m at the top of the flank (94.75/m 2 ) compared to the base of the slope partway down the cone at 300 m (79.25/m 2 ). However, we do see this relationship in the E and ESE (97.75 and 117.75/m 2 respectively) at 400 m at the base of the cone compared to 79.25 and 63.5/m 2 at 300 m respectively. Higher spatial densities on the flanks compared to the base of the cone are also seen in the two-month distribution ( Figure 9A). In the ESE we see a spatial density of 12.75/m 2 at the top of the flank (at 300 m from the crater mid-point) compared to 4/m 2 at the base (at 400 m). The SE also has this pattern at 400 m and 500 m. It is possible that higher densities at the base of the cone compared to the flank occur in the E, SSE, and S though we do not have flank densities to compare this to. Additionally, the azimuths that have higher spatial densities at the base of the cone are not the same ones that have higher median diameters at the base, leading to further uncertainty over the contribution of rolling and bouncing to larger scale spatial densities and size distribution.
We also see a lower spatial density in the SE at 300 m from the crater mid-point compared to the same distance in the ESE and SSE ( Figure 4B). This area is on the side of the crater rim dipping in toward the vents Presses universitaires de �rasbourg ( Figure 10). Landing on the vent-dipping side might explain lower densities in the SE at 300 m as the VBPs might land and then either roll back into the crater or stop on lower part of the inner crater rim. None of the other mapping squares at a 300 m distance are on a slope dipping towards the crater-the others are either at the base of the inner rim or are on the crater rim where it is flatter. We see three examples of empty impact craters with similar sized VBPs downhill (in toward the inner rim/crater) in the mapping square, as well as multiple cases of VBPs being on the edge of their impact craters in the direction of maximum slope towards the base of the inner rim ( Figure 10). It is evident that bouncing and rolling occurs and affects spatial distribution patterns, but it is unclear to what extent the overall spatial distribution is affected.
The deposition of tephra fall has the ability to impede preservation of the VBP field and can affect the estimations of hazard and risk. Yasur not only erupts VBPs but also smaller sized tephra fall (ash and lapilli). We see a complete burial of the VBPs in the July-August images and deposition of new VBPs by the twomonth distribution in the SE 200 m from the crater mid-point. Tephra fall deposition could also contribute to the lower spatial densities seen in the SE where windblown tephra may have been deposited. To assess this possibility, we looked at the number of VBPs that were mapped in the July-August images but were no longer visible in the October images ( Figure 11) and the wind directionality over the mapping period. Of the VBPs mapped in the July-August orthophotos, 50 % were not visible in the October orthophotos. When divided into distance and azimuth from vent area, it is apparent that the SE does not have an increased percentage of missing VBPs compared to other azimuths and thus it is unlikely that SE ash and lapilli tephra fall directionality or preservation issues are a main factor in the VBP azimuthal spatial distributions observed. Additionally, the prevailing wind direction for Tanna is from the E-SE [Vanuatu Meteorology and Geohazards Department 2019a] which deposits windblown tephra fall in the opposite direction to the W to NW. The anemometer deployed on the south rim of the crater recorded N-SE winds throughout the July-August trip, though predominantly E-SE [see Jolly et al. 2017], that again would not be responsible for depositing ash and lapilli tephra fall preferentially to the SE. However, deposition of tephra fall in proximal areas of the cone is still likely, especially when wind speeds are low and this along with VBP deposition would explain the 50 % of VBPs that could not be identified in the second survey in October.
As mentioned in the previous sections, VBPs can fragment on impact with the ground producing multiple smaller clasts that cluster around the area of impact, increasing the local spatial density. VBPs can also fragment in flight from collisions with other VBPs or due to ductile deformation [Taddeucci et al. 2017]. To investigate whether the mapped spatial densities were influenced by fragmentation, we divided each mapping square from the October distribution into four 10ˆ10 m squares and calculated the standard deviation for each mapping square. High standard devia-tions could indicate large amounts of smaller particles created by fragmentation, if supported by visual evidence of fragmentation (clustering of small particles near an impact site). We found a range of standard deviations from 0.43 to 27.68 in each of the mapping squares. However, the high standard deviations were not found to be linked with visual evidence of fragmentation, and therefore we cannot conclude that fragmentation is a dominant processes controlling our size distribution. Video footage of the flank and the entire ballistic trajectory would help to corroborate the fragmentation and high standard deviation theory.
Directed eruptions can produce asymmetric VBP fields Gurioli et al. 2013;Houghton et al. 2011;Kilgour et al. 2010]. From the observational logs taken over the two-month period ( Figure 7C), more VBPs were reported landing to the south of the vents in South Crater (i.e. both on the southern crater wall and outside the crater rim to the south). These eruption directionality observations match the increased spatial density in the S and SSE in the mapped distributions, indicating that eruption directionality is the likely cause of the observed VBP spatial distributions. The E-ESE increase in spatial density seen in the mapped distribution can also be explained Presses universitaires de �rasbourg by explosion directionality, with 27 % of the directed explosions recorded on video between 30 July and 1 August from North Crater going to the east ( Figure 7A). Higher spatial densities in the SSE and E are also noted in maps produced following increased activity in 1994 [Global Volcanism Program 1995]. However, we observed from the 3-day video footage captured during the July-August fieldwork a predominantly SE direc-tionality (of those eruptions that were directed) from South Crater ( Figure 7B) which is not supported in the mapped distributions (we would see an increased spatial density in the SE). This could be due to the short timeframe of recorded video (3 days) not representing the longer-term directionality. Salvatore et al. [2018] recorded variations in jet directionality over short timescales of hours to days at Stromboli from a Presses universitaires de �rasbourg four-year record and that this was the result of changes to vent size and morphology. Over longer timescales (months to years) vents were observed to migrate and merge. Both of these factors could cause the eruption directionality and asymmetry seen in our mapping.
In summary, VBP spatial density changes with direction from the crater. This is likely caused by explosion directionality and minimally by slope changes and VBP fragmentation on impact. Topographic shadowing was not observed from the NE to S (though may affect the unmapped western side which is higher than the east and south), and ash and lapilli tephra fall and VBP burial of the VBP field does not explain the observed changes in spatial density with direction but does explain the lack of preservation of many VBPs. As our mapping of the field is limited to the NE to S (clockwise) we do not capture the VBP field around the entire volcano. However, we do capture the areas in which tourists and guides utilise and the predominant eruption directionalities, recorded during our field campaigns. Further work is needed to map the VBP distribution from the S to the NE (clockwise) to ascertain the full picture of the VBP hazard at Yasur. 4.3 Contrasting size, density, and directionality data A larger mean VBP diameter is observed at 300 m distance in the SE in both distributions and at 500 m in the October distribution compared to other azimuths. In the SE azimuth, we also observe the lowest VBP spatial density and highest proportion of eruption directionality from the GoPro videos. The contradicting Go-Pro eruption directionality and spatial density, as explained earlier, is likely due to the GoPro footage representing a small timeframe and not the longer-term preferential directionality. However, this does not explain the large clast size anomaly in the SE. We are unsure as to what is causing this. It could be the result of larger less frequent eruptions ejecting larger VBPs further and wider, depositing large clasts in the SE as well as in all other directions, in comparison to the more common smaller eruptions that are directed in other directions, depositing smaller VBPs in these areas. Unfortunately, our GoPro footage did not record any such sufficiently large events where large VBPs reached these distances.

Differences in the density with distance trend between the two-month and October distributions
A decrease in the VBP spatial density with increasing distance from the vent area is evident in both the twomonth and October distributions. There are low VBP spatial densities in distal areas for the most part because fewer VBPs are ejected to this distance in the small frequent explosions. However, these densities are also affected by the presence of vegetation and the location of the ash plain at these distances. Though both distributions show this decreasing density with distance relationship, they do not follow the same trend. A steeper decrease is seen in the two-month distribution. This could be due to two factors. Firstly, the October distribution covers a longer time period where larger eruptions likely occurred over months to years prior to the survey date, ejecting a greater number of VBPs further than in the two-month period. Secondly, proximal ash and lapilli fall could be covering older VBPs in the two-month distribution (though the VBPs were still deposited within the two months).

Hazard and risk implications
Our results may have implications for risk management. Understanding which areas have higher densities of VBP impact and how the spatial density decreases with distance could be used to determine where transient visitors may be safer to visit or where infrastructure is sited. The spatial density of VBPs reported here should be considered conservative due to burial of proximal deposits by subsequent tephra fall, masking the true deposition; hence the hazard to visitors may be underestimated. Preservation issues with VBP deposits (e.g. burial of deposits as witnessed here, but also erosion of craters by weather and other factors) show how time sensitive VBP mapping is to capture the true distribu-

Presses universitaires de �rasbourg
Page 197 Volcanic ballistic projectile deposition from Yasur Volcano Fitzgerald et al., 2020 tion from an eruption.
This VBP hazard research does not encompass the hazard from larger eruptions which would likely produce a larger VBP field that may extend into the jungle at Yasur. Reports from a more explosive eruption period in June-July 1999 describe ejected VBPs reaching as far as 600 m from the nearest edge of the crater [Global Volcanism Program 1999]. This would require fieldwork at greater distances from the crater and modelling to understand the potential distributions. Further assessment would also be beneficial over a longer time period (e.g. 12-24 months), like that performed by Salvatore et al. [2018], to capture any temporal variability in eruption directionality and potentially eruption size, which may reduce the spatial variability seen in this study. We note the Vanuatu National Disaster Management Office, the government agency responsible for disaster risk reduction, have a permanent exclusion zone that includes the crater and the top of the flanks extending around from the viewing points. Additional 'zones' are restricted during periods of elevated seismic activity (monitored by Vanuatu Meteorology and Geohazards Department [VMGD]) associated with larger or more frequent explosions (Table 1 and Figure 6 show the correlation between elevated seismoacoustic activity and an observed increase in the proportion of larger explosions). Danger Zone A includes the cone and volcano viewing is only allowed from the carpark. Danger Zone B is an~500 m radius from the edge of the crater rim that includes the car park and parts of the ash plain [Vanuatu Meteorology and Geohazards Department 2019b]. As shown in Figure 6 and Table 1, acoustic sensors are highly complementary to the seismic and visual observations with regards to energy dynamics of eruptions and VBP hazard. Acoustic monitoring at the sole seismic station monitored by VMGD may add value to monitoring and risk decision making, though further investigation of its limitations is needed which is outside the scope of this study.
As well-known with volcanoes, past behaviour does not necessarily dictate future behaviour. This needs to be considered when using eruption directionality to assess potential VBP hazard areas (i.e. areas exposed to volcanic ballistic projectile deposition), which could then be used to inform where access on the volcano may be tolerable within risk management decision making. Our data shows that the VBP hazard is highly spatially variable, in part due to explosion directionality, but that directionality may vary significantly temporally. Therefore directionality results from this study alone should not be used as the basis for hazard and risk decisions. Longitudinal studies on directionality patterns are needed to determine the stability, usefulness and representativeness of directionality in VBP hazard and risk decision-making. In the interim, it may be more conservative and pragmatic to use spatial density and how this changes with distance from the vent to describe the hazard intensity and for a radius of a cer-tain distance around the vent based on spatial density to be used to define 'hazard areas' so that any future change in directionality is accounted for. In the past hazard areas have been defined by the maximum travel distance of a specific diameter VBP. However, as discussed earlier, eruptions can produce normal and reverse VBP distributions and thus setting a high or low diameter threshold could cause over or underestimation of travel distance depending on the eruption style. Combining these two approaches where spatial density of VBP deposition above a dangerous diameter is used as a measure to define hazard areas could also be considered in the future. Further investigation is needed to determine what measure is best to apply to characterise VBP hazard.
The VBP hazard mapping and analysis presented in this paper is intended to improve our understanding of the extent and distribution of VBP hazard from Strombolian eruptions. This study should not be used in isolation to guide or inform risk management decisions at Yasur or any other volcano. Rather this study provides data and analysis for one step (of many) in the process of assessing risk to visitors and infrastructure, which can contribute to informing decisions as to how the risk should be managed.

Conclusions
Understanding where VBPs land, their spatial density and size, and how this changes with distance and direction from the vent can inform more effective risk management. This is especially important on volcanoes that are frequently visited such as Yasur Volcano, Vanuatu. We used a drone to image part of the VBP field around Yasur and GoPro cameras and daily observations to record explosions. Mapping and analysis of the images and video revealed: 1) The spatial density of the VBPs decreased with increasing distance away from the crater 2) The median VBP diameter decreased with distance from the crater (evident up to 600 m distance), similar to that seen from Strombolian eruption deposits at Stromboli and also from phreatomagmatic eruptions around the world. This is likely due to the presence of a gas jet decreasing drag on particles, allowing smaller ones to travel higher and further than larger ones.
3) Higher spatial densities of VBPs are seen in the S-SSE than in other directions. We attribute this to directionality in the explosions and not to topographic shadowing.
4) A different directionality was observed in the Go-Pro videos than in the daily observations and VBPs mapped from aerial photos. The mapping represents a longer time frame than the observations or Presses universitaires de �rasbourg videos and likely represents the longer-term directionality trend, showing that directionality is time variable.
5) Field preservation issues are apparent, with the burial of VBPs by further deposition of tephra fall. Therefore, VBP mapping is time sensitive after an eruption and delays in mapping will likely be accompanied by decrease in field preservation and result in an underestimation of VBP hazard.
When creating VBP hazard and risk maps or making risk management decisions at continuously erupting volcanoes it is important to remember that the hazard can vary over relatively small areas and over different time frames. Assessing the hazard over as much of the volcano as possible and over as long a time frame as possible will produce more effective results to subsequently base successful risk management decisions from. assistance, and Sandrine Cevuard at VMGD for monitoring support. Seismic and infrasonic data analysis included the use of ObsPy [Beyreuther et al. 2010]. We are grateful to the editor Jamie Farquharson and two reviewers Benjamin Bernard and Fabian Wadsworth for their constructive comments and suggestions.

Author contributions
R.H. Fitzgerald analysed GoPro footage and catalogued explosions, mapped the VBP distributions and analysed the results, and wrote the manuscript. C. Gomez processed drone imagery and created the orthophotos and DEMs. R.S. Matoza provided the seismic and infrasound data and created Figure 6. B. Simons provided observation logs between two trips. E. Garaebiti coordinated fieldwork and provided local hazard management context. All authors assisted in discussions and writing of the manuscript.