Temporal seismic wave velocity variations at Láscar volcano

We report on the first study using Seismic Wave Interferometry to determine variations of seismic velocities through time, in the vicinity of Láscar volcano in Chile. Seismic Wave Interferometry has been used as a powerful tool to determine spatial and temporal changes of seismic velocities within the Earth. Spatial variations of seismic velocities are related to heterogeneities of material properties, which are expected to occur in a complex structure. However, temporal changes are indicative of dynamic process within the elastic media, and thus, this tool can be used to monitor dynamic processes at volcanic zones. We find consistent variations on three stations close to the volcano, with dv/v of ±0.6%, most likely related to the inflation/deflation process due to fluid movement of magmatic or hydrothermal origin within the volcanic structure. During the observed period of velocity variation, OVDAS reported an increase of volcanic activity evidenced by the increase of the number of long period seismic events, increase of gas emissions and the formation of incandescence above the crater. We suggest that this tool can contribute to the understanding of volcano related dynamic processes, as well as for routine volcano monitoring purposes.


Introduction
Láscar is a stratovolcano located at 23.95° S, 67.88° W, and 5,592 meters of altitude, at approximately 350 km from the trench separating the Nazca and South American plates (Fig. 1).It has been active since at least 50 ka.Its recent activity is characterized by the continuous emission of gases and occasional sub-plininians eruptions, as in April 1993.In 2006, the Observatorio Volcanológico de los Andes del Sur (OVDAS) reported an increase of activity characterized by a plume emitted by the main cone, up to an altitude of 3,000 above the crater.Most recently, between March and April 2013, OVDAS reported an increment of the number of seismic Long period events, together with a degasification process and incandescence near the crater of Láscar volcano.
Due to the availability of seismic data from 5 broad band stations, near the crater, it poses an opportunity to test if with Seismic Wave Interferometry it is possible to detect changes of velocity around the volcano.If it is possible, it would contribute towards the understanding of the associated processes and could become a tool for monitoring purposes.
Seismic Wave Interferometry is based upon the concept that ambient noise, originating everywhere throughout the Earth's surface, acts as a continuous source of waves propagating along every possible direction.When considering a pair of stations, the cross-correlation of their signals produces constructive interference for waves propagating along the great circle path between both stations, and destructive interference among all others paths.Therefore, considering long lapse times, the cross-correlation signal between two stations becomes an empirical Green's function for a virtual source located at one station and recorded at the other.The Green's function includes the main wave and its coda.The coda represents the contribution of waves originating at the virtual station and being scattered in its vicinity before arriving at the other station.Slight waveform variations through time along the coda, or stretching, is due to variations of seismic velocity through time within the volume between both stations.

Temporal seismic wave velociTy variaTions aT láscar volcano
The Seismic Wave Interferometry methodology consist in quantifying the stretching within the coda, giving the amount seismic velocity variation through time (Campillo, 2006;Snieder, 2006).This methodology has been applied to several volcanic centers worldwide.For instance, at the Piton de la Fournaise volcano, in La Reunion Island, a decrease of 0.1% of seismic velocities before eruptions were observed (Brenguier et al., 2011).At the Merapi volcano, Sens-Schönfelder and Wegler, 2006, observed a strong seasonal variation of seismic velocities, and modelled as due to a hydrological effect.Since OVDAS is monitoring several volcanic centers in Chile, we asked if this methodology can contribute towards monitoring velocity changes around these volcanoes, by testing on specific case, namely Láscar.

Data and processing
Láscar volcano is being monitored since 2013, by OVDAS, with 5 broad band stations (Fig. 1), recording at 100 samples per second.Stations LAS and QUE are located at approximately 4.5 km from the crater, LEJ at 5.5 km, PUN at 2.5 and TAL at 20 km.Here we use vertical components between the period of February 12 to July 11, 2013.We remove the instrumental response and resample at 50 Hz.
We first generate time series of 1 day duration, and use a 1-bit normalization filter in amplitude to eliminate the presence of seismic events (Shapiro y Campillo, 2004;Bensen et al., 2007).We then calculate the cross-correlation between all pair of stations (Fig. 2), and autocorrelation for all stations.The cross-correlation can be considered as the Green's function where one source acts as a virtual source (Wapenaar et al., 2010).For the autocorrelation, the virtual source is located at the same station.In this case, the coda generated is due to scattering within a volume surrounding the station.Therefore, velocity variations occurring in a region close to a station, would affect both the autocorrelation and crosscorrelations associated to such station.
The range of frequencies to use depends directly on the geometry of the network, and the object of study.Sens-Schönfelder, 2008, suggests to use low frequencies for cross-correlations, due to interference generated by high frequencies on the far field.In our case we concentrate on the 0.3-6 Hz band for cross-correlations and 4-6 Hz for autocorrelations.
Once the correlation functions are obtained for each day, during the complete period of 6 months, we stack them to obtain the average correlation function for such period.The main feature of the average cross-correlation, is the surface wave propagating at an average velocity of v s =4.2 km/s (Fig. 2).To the daily correlation function we subtract the average correlation function, to identify anomalous features through time.
Slight changes of seismic velocities during the period of observation, would produce slight variations of the coda within such period.If there is a slight decrease/increase of seismic velocity in a region, waves passing through this region will be stretched like an accordion, increasing/decreasing the period of waveforms.Seismic Wave Interferometry concentrates on the analysis of the coda, to identify such waveform variations, calculating the amount of stretching, dt/t, (Snieder, 2006).The amount of stretching is obtained by modeling the waveforms for different stretching factors, and choosing the one with the highest correlation factor R between model and observations.Here we select values of stretching only when the correlation function R>0.8.Knowing the stretching factor dt/t, the velocity change through the region is given by dv/v=-dt/t.
The source of waves within the coda of the correlation function, are waves scattered in the vicinity of the involved stations.The longer the time within the coda, the larger the region involved.
Here we choose windows of 5-10 s and 5-7 s in the correlation matrix (Fig. 3), representing roughly a volume of 10-15 km radius, which is expected to be the region most affected by volcano processes.

Results and Discussion
Our main result is the calculated stretching factor dt/t shown in figures 4 and 5 for cross-and auto-correlations.Blue dots indicate the stretching values with the highest correlation factors R (for R>0.8).There are common features on the observed velocity variation through time, obtained from the cross-and auto-correlation for different stations, as shown in figure 6.If these features are due to a real phenomenon, and not to a processing artifact, then, it would prove to be a powerful tool for studying volcanic processes near volcanoes, and also serve as a monitoring tool.Among these common features one can observe a period between March and April 2013, with an abrupt increase of velocity (~0.6% at QUE station) followed by a slow decrease.These changes could be due to a deflation process followed by an inflation process related to a loading and unloading of a hydrothermal and/or magmatic system beneath the Láscar volcano.
It is difficult at this stage, and without strong independent evidence, to state that these features are a real phenomenon.However, it is reasonable to expect the presence of f luids near the crater of an active volcano.These fluids can move due to complex interactions between magma chambers and/ or hydrothermal sources.The increase or decrease of the amount of fluids in the media, will produce changes in the average seismic velocity, which can be detected by seismic waves passing through these regions.Changes of velocity of the order of less than 1% can not be detected measuring the variations of the arrival time of main seismic phases.However, these changes of velocity are capable to deform the waveform of high frequency seismic waves, in a way that can be modeled by the stretching of waveforms within the high frequency coda.The variation of seismic velocities due to the presence of fluids has been shown by Sens-Schönfelder, 2008, at Merapi volcano, using this same technique.
If our observations are due to a real phenomenon, in the following we speculate if it can be related to other inferred processes.
The deformation of the crust is sensitive to processes within, including the migration of fluids.In the past decade a great progress has been made with geodetic measurements, especially near volcanoes.However, it has been shown using InSAR that the deformation near Láscar volcano is considerable less than other zones within the Andes (Pritchard and Simons, 2004;Henderson y Pritchard, 2013;Fournier et al., 2010).Therefore, it is likely that the source of the inferred velocity change can be due to fluid motion in a very local region, and by small amounts, and can not be detected on a larger scale as observed with InSAR.
On smaller scales, studies using Magnetotelluric (Díaz et al., 2012), suggest the presence of hydrothermal source at a depth of 1 km slightly to the south of the crater.This can clearly affect seismic velocities, but within a small region, therefore it is an interesting possibility.
There are seismic studies related to the Láscar eruption in 1994-1995 (Hellweg, 2000), suggesting the presence of fluids due to complex and time varying frequencies of tremors.However, with these observations it was not possible to infer the location or geometry of the region with fluids.Another interesting analysis related to the eruption of 1994-1995, was done by Paves et al., 2006, using InSAR images, inferring a subsidence of 17 mm near the crater.This subsidence was modeled by Paves et al. 2006, as due to a volume change of 2×10 3 m 3 , of a spherical source located at 180 m below the crater.This is a very local source, and difficult to find a relation with the source of fluid motion 20 years later.

Conclusions
We use Seismic Wave Interferometry to determine variations of seismic velocities through time, in the vicinity of Láscar volcano in Chile.We find consistent variations on three stations close to the volcano, with dv/v of ±0.6%, most likely related to the inflation/ deflation process due to fluid movement of magmatic or hydrothermal origin within the volcanic structure.During the observed period of velocity variation, OVDAS reported an increase of volcanic activity evidenced by the increase of the number of Long period seismic events, increase of gas emissions and the formation of incandescence above the crater.We suggest that this tool can contribute to the understanding of volcano related dynamic processes, as well as for routine volcano monitoring purposes.

FIG. 1 .
FIG. 1. Station distribution around Láscar volcano.The study zone is marked by a yellow square.

FIG. 2 .
FIG. 2. Stacked cross-correlation for each station pair.The red line shows the surface wave arrival.

FIG. 3 .
FIG. 3. Correlation matrix for each station pair.Correlation matrix (xcorr) using the vertical component (Z) for each station pair.

FIG. 6 .
FIG. 6. Up: Relative velocity changes from the cross-correlation between LASZ-PUNZ, LASZ-QUEZ and PUNZ-QUEZ (the last Z means vertical component of the station) and their autocorrelations.Down: Cumulative long period and volcano tectonic events from February to July 2013.