Soil-structure interaction assessment combining deconvolution of building and field recordings with polarization analysis: application to the Matera (Italy) experiment

In this study, the wavefield radiated from a building to its surroundings is identified and extracted from M4.6 earthquake recordings collected by sensors installed in a building and on the nearby athletic field in Matera (Italy) using a new approach for soil-structure interaction assessment. The proposed approach for earthquake data analysis combines in an innovative way two methods already used in seismology and engineering seismology: deconvolution and polarization analysis. The approach enables the identification, reconstruction, and characterization of the wavefield radiated from a vibrating building into its surroundings, and the estimation of the amount of energy associated with it. The approach consists of four steps: (1) estimation of the resonant frequencies of the building, (2) deconvolution of the earthquake recordings from a building and its surroundings, (3) identification of the seismic phases, reconstruction of the signal transmitted from the building to its surroundings, and estimation of its energy, and (4) polarization analysis. Analysis of recordings of the M4.6 event highlighted that the motion related to the wavefield radiated from the building to the ground was mostly linearly polarized in the radial and transverse planes, while a clear ellipticity was observed only in the horizontal plane. The wavefield radiated from the building might be dominated by unconventionally polarized surface waves, i.e., quasi-Rayleigh waves or a combination of quasi-Rayleigh and quasi-Love waves. The results indicated that the energy transmitted from the analyzed vibrating building to its surroundings was significant and decreased the ground motion shaking due to the out-of-phase motion.


Introduction
The destruction in Mexico City caused by the 1985 Michoacan earthquake (M8.0) at an epicentral distance of about 350 km highlighted the impact of the interaction between soil and built structures on damage in the case of strong ground motion (Wirgin and Extended author information available on the last page of the article Bard 1996).In the following years, soil-structure (SSI), site-city (SCI), and structuresoil-structure interaction (SSSI) were extensively studied, from both engineering and seismological perspectives.SSI, SCI, and SSSI have been investigated in engineering seismology using analytical models (e.g., Paolucci 1993;Guéguen et al. 2002), numerical simulations (e.g., Kham et al. 2006;Liang et al. 2013;Isbiliroglu et al. 2013;Kumar and Narayan 2017) and laboratory experiments (e.g., Pitilakis et al. 2008;Chandra and Guéguen 2017).However, there are presently few published examples of SSI, SCI, and SSSI studies that involve empirical data necessary to validate and supplement the aforementioned types of studies (e.g., Jennings 1970;Wirgin and Bard 1996;Guéguen and Bard 2005;Parolai et al. 2005;Petrovic and Parolai 2016;Petrovic et al. 2018).
By analyzing earthquake and explosion data from two experiments at the Volvi test site in northern Greece, Guéguen and Bard (2005) showed that the wavefield radiated from a building into the ground via the foundations is significant and contaminates the free-field motion.Large variations in the amplification and resonant frequency were also observed by Chávez-García and Cárdenas-Soto (2002) who analyzed the horizontal-tovertical spectral ratio of microtremors in the so-called "free-field".They attributed these phenomena to the presence of neighboring buildings which significantly affected the free-field ground motion.In some recent studies, Petrovic and Parolai (2016) and Petrovic et al. (2018) showed, by analyzing recordings of weak and strong motions that a significant amount of energy is radiated as seismic waves from the structure to the ground during an earthquake.
Analysis of the response of a building and wave propagation analysis within the soilstructure system is one way to analyze SSI.Seismic interferometry by deconvolution (e.g., Snieder and Şafak 2006) has become a widely used approach to evaluate the response of a building to an input impulse.The approach has been applied to both noise (e.g., Nakata and Snieder 2013;Bindi et al. 2015;Mordret et al. 2017) and earthquake recordings (e.g., Nakata et al. 2013;Pianese et al. 2018;Guéguen et al. 2019;Skłodowska et al. 2021) collected in buildings.The advantage of the deconvolution approach is that it is independent of the source characteristics of the analyzed earthquake (Vasconcelos and Snieder 2008).In some recent SSI studies, the deconvolution approach has been applied to earthquake recordings from sensors installed in buildings and nearby boreholes (Petrovic and Parolai 2016;Petrovic et al. 2018), thus extending previous applications to both building (e.g., Snieder and Şafak 2006) and borehole data (e.g., Parolai et al. 2009).Additionally, the first attempt to study wave propagation from a building to the soil by applying a deconvolution approach to seismic noise recordings of sensors installed both in buildings and in the surroundings was made by Petrovic et al. (2019).
Although SSI, SCI and SSSI studies have been performed for years, analysis of the polarization of the wavefield radiated from a vibrating structure into its surroundings has gained limited attention so far.Even though the importance of polarization in earthquake data analysis has been underlined in the past (e.g., Vidale 1986), there are few studies focusing on SSI or SCI assessment.Recently, Kumar and Narayan (2019) modeled the polarization of incident S-waves in SCI effects and observed considerable SCI effects on the free-field motion with different S-wave polarization.However, there are still very few studies on the polarization of the wavefield radiated back from built structures to the ground (e.g., Chávez-García and Cárdenas-Soto 2002).In one such study, Cardenas et al. (2000) concluded that in the case of the Jalapa building in Mexico City, the polarization analysis of the wavefield radiated from the structure to the ground did not provide clear results due to the influence of the incident wavefield on the recorded motion.However, the energy radiated from a vibrating structure might have an impact on the neighboring buildings and therefore, knowledge of the wavefield polarization needs to be explored further to better understand the physics behind this phenomenon.
In this study, we introduce a new approach for analyzing soil-structure interaction to evaluate the wavefield radiated from a structure into its surroundings.The approach is based on combined analysis of earthquake recordings in a building and its surroundings.In contrast to previously proposed approaches, the focus of this study is on polarization analysis of the wavefield radiated from the building.The approach enables the identification of wave types (e.g., P-, S-, Rayleigh, and Love waves) and the estimation of the energy of the radiated wavefield.
The proposed approach is an innovative combination of the deconvolution of earthquake data collected in a building and its surroundings, the application of constrained deconvolution (Bindi et al. 2010) to separate the wavefield components, and polarization analysis for wavefield characterization.The approach consists of four main steps: (1) estimation of the resonant frequencies of the structure to define the frequency band in which SSI might occur; (2) deconvolution of the recordings from sensors installed both in the building and its surroundings; (3) identification of the part of the seismic wavefield transmitted from the building to the ground and estimation of the associated energy; (4) polarization analysis of the reconstructed signal related to the wavefield radiated from the building.We apply the approach to recordings of an earthquake that occurred during an SSI experiment in Matera (Italy) in 2019, where a large number of sensors were installed both in a building and on a nearby athletic field.

Site description
In October 2019, an SSI experiment was carried out in Matera, southern Italy (Fig. 1a).Matera is a site of particular interest in terms of local seismic response because the first resonant frequency of most buildings in the city is relatively close to that of the soil (Gallipoli et al. 2020).Stratigraphic sequences in most parts of the Matera city show strong impedance contrasts due to the presence of terrestrial deposits above calcarenites and Altamura limestone (e.g., see Fig. 1A in Gallipoli et al. 2020).The average shear wave velocities for the subsoil layers of the test site were estimated by Gallipoli and Lupo (2012) through down-hole tests (Table 1).Based on extensive measurements carried out as part of the CLARA project, the soil resonant frequency of the test site is estimated to be around 1.6 Hz (CLARA project 2020).
The test site chosen for the experiment consisted of a 7-floor reinforced concrete residential building and a nearby athletic field (Fig. 1b).The instrumented building was composed of three very similar rectangular units (two external and one in the middle) connected with joints.The length of each unit was 23 m, making in total 69 m in the longitudinal direction (Y).The width of each unit was 12 m in the transverse direction (X).The average floor height was 3 m.

Array description
During the experiment, a dense seismometer array was installed in both the building and on the athletic field (Figs.1b, 2).The experimental setup consisted of 16 three-component digitizer-sensor combinations (from here on simply "sensors"): three CUBE and eight Reftek digitizers, both connected to LE-3Dlite 1 s seismometers with sampling rates of 200 Hz and 250 Hz, respectively, and five Lunitek Sentinel Geos with sampling rate 250 Hz.
For the analysis of wave propagation through the building, the vertical arrays within the building always included at least sensors placed on the highest possible floor and in the basement.Since the building has a pitched roof, it was impossible to install sensors on the roof, so the highest sensor was installed in the stairwell between the 7th floor and the roof.Due to the difference in elevation between the athletic field and the bottom of the analyzed building, the sensors on the athletic field were installed about 7.5 m above the lowest sensors in the building.
For the analysis, all sensors are divided into four arrays: Array 1, 1A, 2, and 2A.Array 1 (Fig. 2, blue triangles) is composed of three sensors installed in the stairwell of Unit 46 (Sensor 101, 7th floor; Sensor 102, between the 3rd and 4th floors, Sensor 103, in the basement) and two sensors installed on the athletic field at 20 m (Sensor 104) and 50 m from the building (Sensor 105), perpendicular to the building's vertical array and the longer side of the building.Array 2 (Fig. 2, green triangles) consisted of four building sensors installed

Data description
Seismic noise recordings were carried out on five consecutive days (22.10.2019 − 26.10.2019).During this period, on 25.10.2019 at 4:31 am (UTC), a M4.6 earthquake with an epicentral distance approximately 145 km southwest of the Matera test site was recorded by the sensors (Fig. 1a).
For the analysis, the data recorded at a sampling rate of 250 Hz is resampled to 200 Hz to make the sampling rate uniform between the sensors.Since the data come from different instruments, the instrumental response is deconvolved from the recordings after removing the mean and trend.An example of 2-hour noise recordings and the earthquake time histories with their corresponding Fourier amplitude spectra is presented in Fig. 3.The Fourier amplitude spectra of the noise recordings are calculated by averaging spectra from 20 s-long time windows (360 windows in total).

Methods
A fundamental assumption for the proposed approach is that the earthquake recordings originate from an event located far enough away that the input can be assumed to consist of only vertically propagating waves.The method requires at least two sensors located in the structure (top and bottom) and one located on the surface in the surroundings.

Resonant frequencies of the building
We suggest using the spectral ratio method to estimate the resonant frequencies of the building, due to its simplicity and the small number of required sensors (minimum of two, at the top and bottom of the structure).The method is based on the division of the Fourier spectra of the same components of simultaneous recordings from two sensors installed on the top and the bottom of the building (e.g.Guéguen and Bard 2005;Parolai et al. 2005;Skłodowska et al. 2021).Apart from the peaks related to the resonant frequencies, the method gives the frequency band in which most of the vibrational energy of the building is concentrated.

Deconvolution of the recordings from the building and its surroundings
To estimate the band-limited Green's functions describing wave propagation from the building to its surroundings, we filter the earthquake data with a causal filter (e.g. a Butterworth filter) applied over the frequency range determined in the first step.In the frequency domain, the deconvolution in a 2D profile (i.e. of vertically and horizontally distributed sensors) is defined as where u 1 (x 1 , z 1 , f ) and u 2 (x 2 , z 2 , f ) are the Fourier transforms of the recordings at locations (x 1 , z 1 ) and (x 2 , z 2 ) , respectively, and x and z indicate the horizontal and vertical positions of the considered sensor.The inverse Fourier transform of S(f) is the deconvolved wavefield.
(1) To avoid instabilities resulting from division by very small values, we perform the deconvolution using a regularization technique (based on Clayton and Wiggins (1976)): where the asterisk indicates the complex conjugate and is the regularization parameter defined as a percentage of the average spectral power.For the deconvolution analysis in this study, we deconvolve the earthquake recordings in the building and on the athletic field using the signal recorded at the top of the building as a reference (i.e.u 2 (x 2 , z 2 , f ) in Eq. 1).
A similar method for deconvolution-based wave propagation analysis proposed by Petrovic and Parolai (2016) uses recordings of a vertical array of sensors installed in a building and nearby borehole.This study, instead, enables the combined analysis of earthquake recordings in two dimensions: vertical (in the instrumented building) and horizontal (in the surroundings of the building).

Phase identification, wavefield reconstruction, and energy estimation
For further data analysis identification of seismic phases of the deconvolved wavefield in the surroundings of a building that relate to energy transmitted from the building to the ground is required.We use an analytical model to describe the wave propagation from the top of the studied structure to its surroundings in order to simplify the seismic phase identification and interpretation of the results.The analytical model defines the case-specific transfer function between the field sensors recording the ground motion x(t), and the reference station located at the top of the building recording the motion y(t) (Appendix Fig. 11).The model is based on a simplified geometry and the assumption that the input wavefield, x 0 (t) , originates from an earthquake at a large epicentral distance that can be approximated as vertically propagating plane waves only.
The constrained deconvolution based on the approach proposed by Bindi et al. ( 2010) is used to reconstruct the wavefield radiated from the building to its surroundings.By applying the positivity constraint around the identified seismic phase of the deconvolved wavefield, the seismic phase of interest is separated.In this study, instead of applying the Landweber method as in Bindi et al. (2010), we convolve the spectra of the identified phases with the spectra of the recording from the top of the building.This reconstructs the ground motion in the field related to a particular seismic phase (i.e.related only to the energy transmitted from a vibrating building to its surroundings).
Finally, the spectral energy of the signal is estimated by calculating the integral of the power spectral density function within the analyzed frequency band, which allows the quantification of the impact of the vibrating building on its surroundings.

Polarization analysis
The polarization analysis is based on the analysis of particle motion of the reconstructed signal that relates to energy radiated from the building into the ground.For this purpose, we analyze the three-component reconstructed field motion (estimated in the previous step) decomposed into three planes of particle motion: radial (X-Z), transverse (Y-Z), and horizontal (X-Y).
( 2) The proposed analysis is similar to the surface waves analysis of Yoshida and Sasatani (2008).In the selected time windows (e.g. the beginning of the signal, the most energetic part, and the coda) we decompose and analyze the particle motion of the reconstructed wavefield transmitted from a vibrating building in each of the three planes.This provides information about the polarization of the signal and allows for the evaluation of different wave types (e.g.P-, S-, and surface waves) that compose the reconstructed wavefield radiated from the building into its surroundings.

Resonant frequencies of the building
To estimate the resonant frequencies of the building we calculate the spectral ratios using both seismic noise and earthquake recordings without smoothing the data (Fig. 4).We calculate the spectral ratios of the seismic noise as averaged spectral ratios of 4320 moving windows (20 s long) from 24 h of recordings.For the spectral ratio calculation, after the removal of the mean, trend, and instrument response, we filter both the seismic noise and earthquake data with a fourth-order Butterworth bandpass filter with corner frequencies 0.5 Hz and 40 Hz.
The fundamental frequencies of the building identified from noise spectral ratios are 2.4 Hz in the X direction and 2.5 Hz in the Y direction.The second peak of both horizontal components is at a similar frequency (2.7 Hz), but the peak in the X direction disappears towards the center of the building (Arrays 2 and 2A).We observe the third peak at 3.8 Hz in all spectral ratios in the X direction.In the Y component spectral ratios, the 3.8 Hz peak disappears towards the middle of the structure.The fourth peak at 6.3 Hz is present in all the spectral ratios for both horizontal components.In the Y component spectral ratio, there is also a peak at around 7.8 Hz, which was not present in the X component spectral ratio.Compared to the noise spectral ratios, the earthquake spectral ratios of the X component do not show a clear peak at 2.4 Hz in Arrays 1 and 1A.Towards the center of the building (Arrays 2 and 2A), however, the peak is more pronounced.In the Y direction, the noise and earthquake spectral ratios show similar peaks present at 2.5 Hz and 2.7 Hz. Around 3.8 Hz, the earthquake spectral ratios are amplified in the X component of Arrays 2 and 2A and in the Y component of Arrays 1 and 1A.At higher frequencies, up to 10 Hz, the earthquake spectral ratios have similar shapes to those calculated from the noise recordings.The vertical component shows no significant amplification in the spectral ratios for both seismic noise and earthquake data.
We observe that the most significant differences in amplitude between the seismic noise and earthquake spectral ratios are in the X component.The spectral ratios of both noise and earthquake recordings show that most of the vibration energy for both horizontal components was between 2 Hz and 10 Hz for both inputs, therefore these frequencies are used to bandpass filter the data before the next step.

Deconvolution of the field and building recordings
Prior to deconvolution, all data are bandpass filtered between 2 Hz and 10 Hz with a fourthorder Butterworth filter (see previous Section).The regularization parameter (Eq.2) is set to 1% of the average spectral power of the reference signal u 2 to eliminate numerical insta- bilities resulting from division by small amplitudes.The chosen value should have little effect on the frequency band of interest.An example of the deconvolved wavefields of the simultaneous earthquake recordings from the instrumented building and the nearby athletic field (Array 1) using the sensor at the top of the building as the reference is presented in Fig. 5.The deconvolved wavefields from Array 2 are very similar, and therefore, for brevity, only results from Array 1 are shown here.
The deconvolved wavefields of recordings from sensors installed in the building show symmetrical peaks in the acausal and causal parts for both horizontal components (Fig. 5a,  b).In the deconvolved wavefields of the sensors installed on the athletic field, the symmetry is lost and they contain numerous phases which are not observed in the building.The deconvolved wavefields of Sensors 104 and 105 have dominant peaks that arrive simultaneously in the acausal part, regardless of the distance from the building.Additionally, the dominant peaks of the field sensors arrive closer to time 0 than those of Sensor 103, installed at bottom of the building (dashed vertical lines in the horizontal components, Fig. 5a, b).
In the vertical component (Z), the separation of acausal and causal phases is not visible.Some peak separation is visible in the deconvolved wavefields at the greatest distance from the building (Sensors 105 and 205), however, the resolution of the deconvolved wavefields is not sufficient to properly identify peak separation.

Seismic phase identification
To identify the seismic phases related to the energy radiated from the vibrating building, we develop an analytical model defining the transfer function between the top of the building and the field sensors (see Appendix 1 for a full derivation).The fitting of the model with the transfer function obtained in the previous step is done in the frequency domain using a grid search approach (e.g., Parolai et al. 2010) by comparing the absolute Fourier amplitude spectra within the considered frequency band.The search is based on finding three model parameters 1 , 3 , and 3 (Fig. 6). 1 is the time delay corresponding to wave propagation with velocity V 1 in the building.2 is the time delay of a wave propagating with velocity V 2 from the bottom of the building to the field sensor.3 is related to the wave travel time in the soil considering the vertical elevation of the athletic field with respect to the base of the building and velocity V 3 .The parameters that provide the best-fit model are chosen based on the minimum rootmean-squared log error (RMSLE) between the absolute spectra (Appendix Eq. 20).For the models, the densities of the building and the soil are approximated as 1 = 250 kg∕m 3 and 2 = 1800 kg∕m 3 (approximate density of wet sand, Lanzo and Silvestri (2016)).The den- sity values only affect the impedance contrast (i.e. the reflection coefficient) and the shear modulus, and not the time delay, which is the primary information sought in this approach.We do not consider attenuation either, since the focus is on identifying the seismic phases related to the energy transmitted from the building to the ground, rather than the amplitude changes.Therefore since we invert for velocities only, the use of approximate and fixed density values is sufficient for the aim of the analysis.
The estimated velocities of the wave propagation in the building and the soil are presented in Table 2.The model is much less sensitive to the velocity of the vertical propagation of the input wave in the soil, V 3 , than to the horizontal propagation velocity, V 2 (i.e., similar RMSLE values for different 3 values considered in the grid search-see Appendix Fig. 12).V 3 has, therefore, a relatively low impact on RMSLE, and the convergence to the minimum value is worse than for the V 2 search.This explains the greater variation of the estimated V 3 values of each component and array compared to V 1 and V 2 (Table 2).
The obtained velocities in the building ( V 1 ) vary between 253 ± 8 m/s to 271 ± 8 m/s, which is a range similar to those available in the literature for reinforced concrete structures (e.g., 322 m/s in Snieder and Şafak (2006), 214 m/s in Nakata et al. (2013), or 200-276 m/s in Bindi et al. (2015)).
The output of the analytical models is simpler than the real data, both in the frequency and time domains (Fig. 7).In the frequency domain, for both horizontal components, the model fits only some of the troughs in the real data results.In the X direction, for the sensors with the largest distance from the building (i.e.Sensors 105 and 205-Fig.7e, i), the model fits the first trough at 2.5 Hz well, which is not the case for Sensor 104, which is the closest sensor to the building (Fig. 7a).The second trough defined by the analytical transfer function for the X component is at around 7.5 Hz for all sensors.In the transfer functions of the real data, however, it is observed only in Array 1 (Sensors 104 and 105) and not in Array 2 (Sensor 205).
In the Y component, similar to the X component, the first trough at 2.5 Hz is captured by the model for Sensors 105 and 205 only (Fig. 7g, k).The second trough at 7.6 Hz is the only trough in the considered frequency range shown by the model for Sensor 104 (Fig. 7c).The models for Sensors 105 and 205 have a double trough at around 7.5 Hz and 8.6 Hz in the Y direction.They do not, however, accurately describe the troughs in the transfer functions of the real data.In general, the transfer functions of the real data are more complex than the best-fit analytical transfer functions.
In the time domain, the models identify three peaks of the deconvolved wavefields.The first two relate to the wave propagation within the building, in both the acausal (Appendix Eq. 16), and the causal part (Appendix Eq. 17).The third relates to the part of the wavefield transmitted from the building to the ground (Appendix Eq. 18).
For the two horizontal components the analytical deconvolved wavefields fit the real data deconvolved wavefields for Sensor 205 (Fig. 7j, l) more accurately than those of Sensor 105 (Fig. 7f, h), where we observe a bigger phase shift of the peaks between the real data and the best-fit model.The analytical solution of Sensor 104 (Fig. 7b, d) accurately delimits phases in the time domain for both components, even though there is a small phase shift of the first peak in the X component of the deconvolved wavefield.We observe a shift of the first peak (Array 1) or the two first peaks (Array 2) in the analytical solution with respect to the real data deconvolved wavefields for the X component of all the sensors.This is in contrast to the Y component, where we observe a good fit.
By comparing the real data with the analytical deconvolved wavefields, the third peak identified by the model relates to energy transmitted from the building to the athletic field.The time delays of the identified phases are around 0.27 s (X component) and 0.28 s (Y) for Sensor 104; 0.54 s (X) and 0.50 s (Y) for Sensor 105; and 0.55 s (X) and 0.52 s (Y) delay for Sensor 205 (Fig. 7).

Reconstructed wavefield
We selected the time intervals for the positivity constraint for the constrained deconvolution method containing the peaks related to the energy radiated from the building by comparing the deconvolved wavefields with the analytical models.We select time intervals from 0.42 s to 0.60 s for Sensor 105 (Array 1), and from 0.39 s to 0.62 s for Sensor 205 Array 2) for all components (Fig. 8a).For Sensor 104, the peaks describing the seismic phases are too close to the dominant peak in the causal part of the deconvolved wavefield and therefore peak separation for the constrained deconvolution is not possible and is not considered in this study.Since the results from field sensors far from the building are very similar, we present only results from sensor 105 in Array 1 (Fig. 8).We show the comparison between the Fourier amplitude spectra of the deconvolved and the constrained wavefield in Fig. 8b-d.We obtain the spectra of the reconstructed wavefield radiated from the building to the athletic field (green lines in Fig. 8e-g) by convolving the Fourier amplitude spectra of the constrained wavefield (dashed blue lines in Fig. 8b-d) with the spectra of the recording from the top of the building ( magenta lines Fig. 8e-g).
For both horizontal components, the troughs in the transfer function around the resonant frequencies of the building (2-3 Hz) have lower amplitude compared to the constrained spectra (Fig. 8b, c).These troughs have a visible effect on the spectra of the reconstructed horizontal components.The spectra of the reconstructed signal (i.e.related only to the energy transmitted from the building to its surroundings-green lines in Fig. 8e-g), therefore, have amplitudes above the recorded field signal (orange lines in Fig. 8e-g) around the resonant frequencies of the structure.
The three components of the reconstructed motion of Sensor 105 from Array 1 are presented in the time domain with the recorded ground motion for comparison in Fig. 9a.

The energy of the radiated wavefield
To evaluate the influence of the building on ground motion in the nearby athletic field, we estimate the energy of the wavefield radiated from the studied building into its surroundings.The energy is calculated as the integral of the power spectral density of the reconstructed wavefield within the frequency band of interest (2-10 Hz).
Within the considered frequency band, the energy radiated from the building at Sensor 105 is around 59% of the energy of the signal recorded by the same sensor in the X component.The Y component of the signal is less affected by the wavefield radiated from the building and is around 17% of the recorded signal.The vertical component of the building motion has the lowest impact on the recording in the athletic field (around 3%).In Array 2, the energy of the reconstructed signal transmitted from the building to the ground is approximately 30% of the recorded signal in the X direction and less, around 24%, in the Y direction.The smallest percentage of the calculated recorded signal energy, around 1%, is for the Z component.
Comparing the energy of the signal recorded by the sensors installed on the athletic field with the "real" input at the surface (considering free surface amplification)-i.e.without the impact of the existing building (obtained by subtracting the spectra of the reconstructed signal related to the wavefield radiated from the building from the spectra of the field recording), shows that in Array 1 the input is around 153% of the recorded signal in the X component, 122% in the Y component and 99% in the vertical component.In Array 2, the energy is smaller for the X component at 135%, and similar for the other components: 121% for the Y component and 99% for the vertical.These values indicate that the input motion is locally diminished, resulting in lower horizontal component ground shaking.

Polarization analysis
To analyze the particle motion of the reconstructed signal related to the wavefield radiated from the building, we select three time intervals of duration 5 s starting at 20 s (P1), 35 s (P2), and 55 s (P3) (Fig. 9a).The intervals are chosen to cover the parts of the signal where most of the recorded energy is located (P1), where most of the reconstructed signal energy is located (P2), and the coda (P3).Figure 9b-j shows the particle motion decomposed into the three planes of motion radial X-Z, transverse Y-Z, and horizontal X-Y.
For the Z component, the reconstructed signal has significantly lower amplitude compared to the recorded signals.This is clearly visible in the decomposed trajectories.In the radial and transverse planes, the polarization of the motion is almost completely linear with a very small contribution from the vertical component, while we only observe the clear elliptical motion in the horizontal plane.The incremental color change highlights the development of elliptical motion in the horizontal plane in the P2 interval.The gray scatters in Fig. 9a-j represent the particle motion of the recorded signal for the corresponding planes and time intervals.In contrast to the trajectories of the reconstructed signal, there is no clear polarization of the recorded signals.Polarization analysis results from Sensor 205 from Array 2 are not presented in this article since they are similar to Sensor 105.

Discussion
Compared to some other case studies where earthquake recordings in buildings are analyzed and where significant co-seismic resonant frequency deviation is observed (e.g., Guéguen and Colombi 2016;Astorga et al. 2018;Pianese et al. 2018;Skłodowska et al. 2021), the variation of the resonant frequencies of the Matera building during the recorded earthquake is small.The shape of the spectral ratios differs between earthquake and noise inputs, however, which could be explained by the excitation of different modes by different inputs (Ditommaso et al. 2012).
The deconvolved wavefields of all the building sensor arrays obtained by deconvolving by the earthquake signal at the top of the building show symmetry between the causal and acausal parts in both horizontal components.These peaks could be associated with the up-and down-going waves propagating within the structure (Snieder and Şafak 2006).The lack of peak separation in the deconvolved wavefield of the vertical component is most likely caused by higher wave propagation velocity in this direction, combined with the low sampling rate of the instruments for such a configuration.For the vertical component, Parolai et al. (2009) suggested that the wave propagation velocity should be close to the P-wave velocity.For the Matera experiment, this leads to arrivals that are separated in time by less than the sampling interval (Bindi et al. 2015) and in space by less than the vertical sensor spacing.
The symmetry between the causal and acausal parts of the deconvolved wavefields of the building sensors is lost in the deconvolved wavefields of the athletic field sensors.The athletic field sensors were elevated above the height of the second floor of the building, shown in the deconvolved wavefields as delays in the arrival time compared to the sensor at the bottom of the building (Fig. 5a, b).This delay is most likely due to the near-vertical incidence of the seismic waves, arriving first at the bottom of the building and only later reaching the sensors on the athletic field.The numerous other wiggles in the deconvolved wavefield might be a result of the influence of local secondary sources (not fully removed due to the limited amount of earthquake data), the complicated geometry of the test site, and the low energy of the recorded earthquake.
The analytical model developed here helps to identify the seismic phases in the complex deconvolved wavefield.In both time and frequency domains, the analytical transfer functions and deconvolved wavefields are much simpler than for the real data.For sensors installed on the athletic field, the real data transfer functions show several features (e.g., many troughs) that are not reproduced by the simple analytical model.For the two sensors furthest away from the building (Sensors 105 and 205), the third peaks of the analytical deconvolved wavefield, related to the transmitted energy, fit the peaks of the real data deconvolved wavefield better in the X than in the Y direction.This better fit could be explained by the orientation of the analyzed arrays, which are aligned in the X direction of the building.
The apparent better fit of the time domain analytical deconvolved wavefield compared to the frequency domain analytical transfer function is caused by the model assumptions.In the time domain, only three peaks of the model (two associated with the up-and down-going wavefields in the building, and a third related to energy transmitted from the building to its surroundings) are considered in the RMSLE.This means that a good fit is obtained only for those three phases, and all other phases of the deconvolved wavefield are neglected.In the frequency domain, however, the grid search of the parameters that describe the best-fit analytical transfer function is performed for the whole time range.Additionally, the differences between the analytical model and real data for both horizontal components are most likely due to simplified model assumptions (e.g.no coupling between the building and the soil, and no influence of surrounding buildings), and the simplified model geometry of a complex experimental setup.
The clear peaks of the Fourier amplitude spectra in the reconstructed signals imply that the dominant frequencies of the horizontal components are strongly affected by the resonant frequency of the building (Fig. 8e, f).A similar phenomenon is observed by Petrovic and Parolai (2016) and Petrovic et al. (2018), where earthquake recordings from sensors installed in buildings and nearby boreholes are jointly analyzed.
The reconstructed radiated wavefield comprises a significant part of the signal recorded in the surroundings of the building, and its energy is up to 59% of the field signal in the direction of the wave propagation from the building (X).Comparing the energy radiated from the building and the input signal (after removing the influence of the building) shows that in the limited frequency band, the impact of the building on the ground motion is significant for the horizontal components, and it actually decreases the ground motion.This observation could be explained as out-of-phase motion of the ground and the energy transmitted from the vibrating building in Matera at the analyzed distances from the structure.
The presented case study is based on the analysis of a recorded weak motion.In case of strong shaking, the response of the soil and buildings can be characterized by nonlinear changes in the wave propagation velocities and resonant frequencies (e.g., Bonilla et al. 2011;Chandra and Guéguen 2017;Astorga et al. 2018;Skłodowska et al. 2021).Chandra and Guéguen (2017) suggested that also the presence of buildings changes the site response as well as its nonlinear response.Especially in the case of rigid buildings, their laboratory experiment demonstrated a reduction in nonlinear behavior, due to the rocking of the structure.Additionally, soil nonlinearity was observed to be one of the parameters affecting the nonlinear response of the system.Those studies suggest that the ratio of the radiated energy versus the so-called "free-field" energy might vary in a nonlinear way in case of strong shaking.Additionally, depending on the type of the analyzed building, the response of the whole system will vary.Therefore, in order to estimate the energy ratio in the strong motion case, more analyses of data with similar sensor configuration but for recordings with larger amplitude are needed.
In this study polarization analysis of the reconstructed wavefield reveals that particle motion is mostly linear in the radial and transverse planes with a small contribution to the vertical component.We only observe the elliptical motion in the horizontal plane.According to Richart et al. (1970), energy transmitted from a circular footing undergoing vertical oscillations is transmitted into the ground as a combination of P-, S-, and Rayleigh waves.Considering that the particle motion of Rayleigh waves is a combination of radial and vertical components (Richart et al. 1970;Stein and Wysession 2003), the reconstructed motion at the Matera test site is not dominated by classically polarized Rayleigh waves.
The estimated wave propagation velocities of the wavefield radiated from the building (Table 2) are smaller than the S-wave velocity ( V S ) of the soil layers of the Matera test site (Table 1), even after considering the significant uncertainties.Therefore we can exclude the strong influence of compressional P-waves in the radiated wavefield, as they are faster than S-waves.The wave propagation velocities are, however, close to the S-wave velocity estimated for the top sand layer in the test site, V S = 180-200 m/s estimated using SH-wave seismic reflection (Lorenzo Petronio, personal communication).
The observed polarization of the radiated wavefield could be caused by the anisotropy of the test site (e.g., complex topography, different media properties along the wave propagation path) which may generate non-classically polarized surface waves (Tanimoto 2004).Yanovskaya and Savina (2004) define quasi-Rayleigh waves by three planes, where the radial and transverse planes have a phase shift, and particle motion in the horizontal plane is elliptical.Another explanation for the observed polarization could be a combination of quasi-Rayleigh and quasi-Love waves, which according to Tanimoto (2004) do not have strictly transverse motion.A further explanation for the observed elliptical particle motion in the horizontal plane could be non-synchronous oscillation of the building in the X and Y directions.To fully understand this behavior, further studies should be performed on both the anisotropy of the test site and the behavior of the building.
The approach developed in this study is complementary to methods proposed in previous SSI studies.For example, by combining this approach with a joint deconvolution approach for analysis of earthquake recordings from an instrumented building and a nearby borehole (e.g.Petrovic and Parolai 2016), it would be possible to analyze the seismic wavefield radiated from a building to its surroundings both in the horizontal (this study) and vertical (Petrovic and Parolai 2016) dimensions.

Fig. 10
Schematic drawing of the changes to the radiated wavefield caused by buildings.a Amplitude of the earthquake input motion (red solid curved arrows) at the surface in the free field (red dashed arrows).b Amplitude of motion of the wavefield radiated from a vibrating building to its surroundings at the ground surface with maxima, minima, and zeros of the propagating wave indicated.c Surface interaction between the wavefield radiated from a building and the ground motion caused by an earthquake at a certain time t.Blue lines indicate the wavefield radiated by a vibrating building.The black dashed line indicates the modified wavefield with possible areas of amplification and deamplification of the signal caused by constructive and destructive interference.In all plots, X-Y defines the surface plane.The top view shows the ground motion at a certain time t at the surface.Apart from input ground motion indicated with red curved arrows in (a) and (c), solid blue, red, and black dashed lines refer to the surface motion Additionally, if the polarization of the surface ground motion produced by an earthquake (Fig. 10a) interacts with the wavefield transmitted from a vibrating building (Fig. 10b), it may modify the ground motion in the surroundings of the building (Fig. 10c).A radiated wavefield polarized in the horizontal plane could interact with the horizontal motion of the ground resulting in an overall amplification or deamplification of the wavefield in the surroundings of the building (Fig. 10c).Understanding the polarization and quantifying the energy of the transmitted waves could provide information about buildings that act as secondary sources of vibration during an earthquake.Knowledge of the dominant wavelength of the radiated wavefield, combined with polarization analysis, could provide information about locations of positive and negative wavefield interference for a given earthquake motion, which could help to identify locations where a future building should (or should not) be constructed.
Furthermore, after validation of the observations presented in this study (using e.g., simulations) the approach could be extended from single-building to multi-building configuration.The polarization analysis of the multi-building case scenario could potentially be used in research on metamaterials applied to geophysical scale (Roux et al. 2018;Lott et al. 2019;Guéguen et al. 2019), to study if e.g., the observed frequency band gap is related only to surface waves or also to body waves.

Conclusions
In this study, we propose an approach for SSI assessment by combining for the first time, two methods that are already widely used for seismology and engineering seismology studies: deconvolution and polarization analysis.By analyzing earthquake recordings using this approach, it is possible to identify and reconstruct the wavefield radiated from a vibrating building to its surroundings, determine the amount of energy associated with it and estimate its polarization.
The characterization of the radiated wavefield (e.g., knowledge of the dominant wavelength, location of the zeros, maxima, and minima, and polarization) could help to identify locations where a future building should not be constructed or areas where the building properties should be adjusted (stiffening or softening slightly the structure).Additionally, regions likely to be subject to increased shaking could be defined as "no development area" through a district plan.By analyzing soil-structure interaction using the approach proposed in this study, it could be possible to avoid the positive interference of earthquake motion with the radiated wavefield.Moreover, with increased knowledge of the radiated wavefield, new buildings could be designed to generate energy that would cancel out the energy radiated from other buildings.
In the present study, we apply the approach to recordings of an M4.6 earthquake from an SSI experiment in Matera, Italy.In contrast to the common assumption that during an earthquake energy is transmitted from a vibrating building to the ground mainly as surface waves, the presented results indicate that in the case of the Matera test site, the energy is transmitted from the building to the surroundings mostly as linearly polarized waves in radial and transverse directions.Since the analysis is based only on the recordings of a single earthquake, however, more data needs to be analyzed to confirm these observations.For validation of the proposed approach and the obtained results, further analysis with different data sets will be carried out in the future.
where 2 is the time delay of a wave propagating with velocity V 2 from the bottom of the building to the field sensor approximated as √ h 2 2 + h 2 3 .3 is related to the wave travel-time in the soil considering the vertical elevation of the athletic field with respect to the base of the building, h 3 , and the velocity V 3 .
In Fourier domain, u(t), d(t), y(t), and x(t) are expressed by By combining Eqs. 8, 9 and 10, an equation describing the recording at the top of the building in the Fourier domain is obtained which can be rewritten as By substituting X 0 (f ) in Eq. 11 with Eq. 13, X(f) can be expressed as By dividing Eq. 14 by Y(f), the analytical transfer function for a sensor in the field, x(t), with a virtual source at the top of the building, y(t), can be therefore described in the Fourier domain by where and The wave propagation travel time in the building, 1 , and therefore wave propagation velocity, V 1 , was estimated using a simple model similar to the one proposed by Snieder and ( 8)

Fig. 1 a
Fig. 1 a Location of Matera, Italy (square) and the epicenter of the M4.6 earthquake (star) that occurred on 25.10.2019 at 04:31:38 (UTC) and was recorded during the experiment.b Satellite view of the test site.The three-component velocimeters (104, 105, and 205) installed on the field are indicated with white arrows

Fig. 2
Fig. 2 Schematic drawing of the positions of the three-component sensors in the building and on the athletic field: a side view, b top view.Blue triangles represent sensors from Array 1, green -sensors from Array 2, yellow -Array 1A, and red -Array 2A.Names of the sensors from Array 1 and Array 2 are next to the corresponding triangles

Fig. 3
Fig. 3 Example of the X component velocity recordings used in this study.a Example of the M4.6 earthquake recordings from sensors in Array 1. b Corresponding earthquake Fourier amplitude spectra.c Example of noise recordings for the same array.d Corresponding moving window Fourier amplitude spectra.The data from sensors installed in the building are marked with black lines and those from sensors installed on the athletic field with gray lines.The sensor names are above the corresponding time histories.The y-axis shows the vertical (building) or horizontal (athletic field) offset of the sensor with respect to the basement sensor.A negative sensor offset indicates that the sensor is offset horizontally (athletic field) with respect to the lowest sensor in the building

Fig. 4
Fig. 4 Comparison of noise and earthquake spectral ratios of top and bottom recordings for each component for Arrays 1, 1A, 2, and 2A.Red lines indicate the noise spectral ratio.Red shading indicates the ± standard deviation of the average noise spectral ratio.Black lines show the spectral ratio of the earthquake recordings

Fig. 5 Fig. 6
Fig. 5 Deconvolved wavefields, obtained using the top of the building (Sensor 101) as the reference recording, for sensors installed in the building and on the nearby field, filtered between 2 Hz and 10 Hz for Array 1 for the X (a), Y (b), and Z (c) components.Positive sensor offsets indicate that the sensors were located in the building (i.e.vertical distance from the bottom sensor of the building).Negative sensor offsets indicate that the sensors were installed in the field (i.e.horizontal distance from the bottom sensor in the building).The gray background indicates the sensors installed in the building.The red dashed line (a and b) indicates the time delay of the acausal peak, related to the up-going wave of the sensor at the bottom of the building.Note that in (c) the wave propagation velocity is too high to observe peak separation

Fig. 7
Fig. 7 Comparison of the transfer functions and the corresponding deconvolved wavefields of the real data (black lines) and analytical solution (blue lines) for sensors in the field: a-d Sensor 104, e-h Sensor 105 and i-l Sensor 205.The first and second columns present the model fit of the X component in the frequency domain (transfer functions) and time domain (deconvolved wavefields), respectively.The same for the third and fourth columns of the Y component.Red vertical lines indicate the identified seismic phases of the real data deconvolved wavefield related to the wavefield transmitted from the building to the athletic field.The time delays of the selected phase peaks are written on the plots

Fig. 8 Fig. 9
Fig. 8 Example of the constrained deconvolution analysis for Sensor 105 from Array 1. a Deconvolved wavefields of all three components, shown by a black solid line.The time interval for the constrained phase selection is marked by gray shading.The blue dashed line indicates the constrained deconvolution in time.b-d Corresponding Fourier amplitude spectra of the deconvolved wavefields.e-g Comparison of the Fourier amplitude spectra of the signal recorded by the sensor at the top of the building (magenta line), the signal recorded by the sensor in the field (orange line), and the reconstructed motion related to the energy radiated from the building into the soil (green line)