1 Introduction

Post-earthquake reconnaissance and damage survey activities have demonstrated to be valuable sources of information to assess local seismicity and structural vulnerability. In the last decades, several studies have focused on the effects of seismic events on reinforced concrete and masonry buildings, to deepen the knowledge of their structural behaviour and to develop sound risk assessment methodologies. Concerning a particularly relevant but vulnerable typology, namely masonry churches, the observational data collected in the aftermath of significant Italian earthquakes, starting from the 1976 Friuli Earthquake (Doglioni et al. 1994), allowed the definition and continuous improvement of solid procedures for rapid damage survey at territorial scale. This is evidenced, for instance, in the applications after the 1997 Umbria–Marche Earthquake (Lagomarsino and Podestà 2004a, b), the 2002 Molise Earthquake (Lagomarsino and Podestà 2004c), the 2003 Valle Scrivia Earthquake (Ruggieri et al. 2022), the 2009 L’Aquila Earthquake (Da Porto et al. 2012; Lagomarsino 2012; De Matteis et al. 2016, 2019a, b), the 2012 Emilia Earthquake (Sorrentino et al. 2014), the 2016–2017 Central Italy seismic sequence (Hofer et al. 2018; De Matteis and Zizi 2019; Penna et al. 2019; Cescatti et al. 2020; Canuti et al. 2021; Ceroni et al. 2022; Sisti et al. 2023) and the 2017 Ischia Earthquake (Salzano et al. 2022). The survey form for masonry churches provided by the Italian Civil Protection supports the damage assessment by investigating the local response of the main structural components, namely dividing the structure into several distinct macro-elements, i.e. façade, nave, lateral naves, transept, dome, apse, roof, lateral chapels and bell tower (Fig. 1), that tend to exhibit independent behaviours when subjected to earthquake loading (MiBACT 2015). Indeed, past earthquakes have demonstrated that monumental masonry buildings exhibit complex structural behaviour far away from the so-called box behaviour (Giuffrè 1991), also addressed as integral or global behaviour. Therefore, the main goal of the existing works has been to identify the different damage patterns and mechanisms observed in each macro-element, with the corresponding damage level, and to estimate the global damage level through a weighted average of each church within the investigated building stock (DPCM 2011).

Fig. 1
figure 1

Subdivision of the church into structural macro-elements

Full size image

In this context, the gathered data are generally exploited to perform statistical analysis to provide damage distributions, to assess the predictive performance of existing vulnerability models as well as to develop novel methodologies suitable to peculiar building typologies and specific datasets in hand, supporting the application of rapid assessment and mitigation strategies. To provide the frequency of damage occurrence, the most common and simple tools for data interpretation are the Damage Probability Matrices (DPMs) which cluster the investigated buildings or components with likely homogenous behaviour according to increasing damage extents. Commonly, damage is classified according to the European Macroseismic Scale, EMS-98 (Grünthal 1998). This provides a qualitative evaluation method according to a 6-grade system of the severity of damage to structural and non-structural elements, namely the damage level ranges from 0 (No damage) up to 5 (Destruction or collapse). Numerous studies have reported DPMs for different types of churches under specific seismic events, for instance, one-nave churches (De Matteis and Zizi 2019; Ceroni et al. 2022; Ruggieri et al. 2022) or three-nave churches (De Matteis et al. 2016, 2019b). Besides the DPMs, data are generally correlated with a ground motion parameter characterising the seismic input. In this regard, a well-established approach consists of using empirical functions that correlate the expected mean damage with the macroseismic intensity, taking into account the vulnerability and ductility indexes related to the structure. The original formulation of this vulnerability curve was provided by (Lagomarsino and Podestà 2004b; Lagomarsino and Giovinazzi 2006). Distinct values of vulnerability and ductility indexes for different masonry building typologies (e.g. church, palace, tower, monastery/convent, mosque or castle type) are suggested in the literature (Lagomarsino et al. 2004; Despotaki et al. 2018). Although the use of macroseismic intensity for predictive vulnerability models is widespread, recently, several methodologies are exploring the use of Peak Ground Acceleration (PGA) as a ground motion parameter and PGA-based approaches are becoming very popular (De Matteis and Zizi 2019; Ceroni et al. 2022; Sisti et al. 2023). Among these methodologies, fragility curves play a key role in forecasting future seismic damage scenarios and are widely proposed for masonry churches (Hofer et al. 2018; De Matteis et al. 2019b; Cescatti et al. 2020; Ceroni et al. 2022; Sisti et al. 2023).

Even though a very large number of studies can be found in the literature for masonry churches as a whole, the aforementioned analyses are rarely conducted to investigate, independently, the behaviour of single macro-elements. Nonetheless, the macro-elements influence in governing the seismic response of such buildings is notorious, as recently highlighted by the analyses reported in (Ceroni et al. 2022; Sisti et al. 2023). Among the macro-elements, historic masonry bell towers have demonstrated to be very susceptible to structural damage, thus likely to be highly vulnerable. For instance, data regarding the 2009 L’Aquila earthquake and the 2016 Central Italy sequence show the activation of damage mechanisms in 60% and 70% of the towers, respectively (De Matteis et al. 2016; De Matteis and Zizi 2019). Moreover, in the case of Central Italy, 16% of towers within the investigated dataset (11 out of 68), as reported in (De Matteis and Zizi 2019), exhibits partial and total failure, presenting higher damage levels when compared with other structural macro-elements. Furthermore, bell towers are generally characterised by high vibration periods, differently from other structural components, suggesting the possibility of a dynamic response completely different from the rest of the building.

DPMs for distinct macro-elements, including masonry towers, are provided by (De Matteis et al. 2016; De Matteis and Zizi 2019), while (Hofer et al. 2018; Ruggieri et al. 2022) also investigate distinct damage mechanisms. A noteworthy contribution on the seismic response of masonry bell towers is the work of (Curti et al. 2008), who exploited the collection of a wide dataset of historic masonry towers affected by different earthquake events. This was used to calibrate the vulnerability and ductility indexes according to distinct damage mechanisms by fitting the vulnerability curve proposed by (Lagomarsino and Podestà 2004b; Lagomarsino and Giovinazzi 2006) to the observed data. After this seminal application, similar studies have recently focused on the validation of the existing formulations against datasets of specific events, as well as on the calibration of new formulation coefficients (Canuti et al. 2021; Ceroni et al. 2022).

Based on these promises, the present work presents a statistical analysis of data gathered from an extensive literature review on studies addressing post-earthquake surveys of historic masonry bell towers in Italy. The collected database, described in Sect. 2, is exploited, in Sect. 3, to provide the global vulnerability of the investigated towers through DPMs for different damage mechanisms and for homogenous clusters of towers. A comparison with results of similar studies is also included. Subsequently, vulnerability curves are obtained, in Sect. 4, using the existing formulations suggested for historic towers found in the literature (Lagomarsino et al. 2004; Curti et al. 2008), and then compared against the observational data. A regression analysis is performed by fitting the samples through a function with the same form of the existing equation, and new coefficients are proposed for the specific dataset in hand. In Sect. 5, fragility curves are developed by correlating the cumulative probability function with the PGA; while the former is obtained by applying the binomial probability density function, the latter is retrieved using an empirical correlation between PGA and macroseismic intensity (PGA-I) proposed by (Faenza and Michelini 2010). Finally, in Sect. 6, the main conclusions of this work are drawn.

2 Database of the investigated building typology

The database collected in the present study encompasses 129 historic masonry bell towers, hit by six moderate to strong earthquakes for which post events reconnaissance reports and scientific publications have been found, upon an extensive literature review, ensuring a reliable identification of the occurred damage mechanisms and extents. In particular, the classification of the collapse mechanisms follows the macro-element approach presented by the Italian Guidelines for Cultural Heritage (DPCM 2011). The search focused on Italian events due to the combination of a large exposed heritage and the significant seismic hazard. This has allowed the identification of a relevant number of cases analysed and reported. The database, openly available as an annex to the present work, gathers together and expands similar existing collections that have mainly focused on a single seismic event, by simultaneously considering different earthquake scenarios. This not only provides numerous samples for future studies but also ensures the comparison between different characteristics of the seismic input and the affected towers.

Details on such a large database in terms of characteristics of the events, the towers and the emerged damage scenarios are provided in the following sections.

2.1 Earthquake events

Table 1 summarises the characteristics of the selected earthquakes, namely location of the affected area, date and time, epicentral coordinates, hypocentral depth and moment magnitude, obtained from the Italian Earthquake Catalogue (Rovida et al. 2020), while Fig. 2 shows the map of the epicentres. The earthquake events were selected upon a review of technical and scientific publications, based on the availability of data regarding structural damage observed in historic masonry towers.

Table 1 Earthquake events considered in the present study
Full size table
Fig. 2
figure 2

Distribution of the analysed earthquakes events

Full size image

A synthetic measure that combines the seismic hazard with exposure and vulnerability of the built environment is the macroseismic intensity. Different scales of intensities have been proposed over the years, e.g., Mercalli-Cancani-Sieberg (MCS), Medveded-Sponheuer-Karnik (MSK), Modified Mercalli (MM), and European Macroseismic (EMS-98). Among them, the MCS intensity scale (\(I_{MCS}\)), with a maximum of 12, is here adopted to characterise the seismic events at the specific location of each analysed tower. Generally, damage to buildings and structures is expected for intensities greater than 5. Table 2 provides a summary of the maximum intensity registered for each earthquake together with the corresponding source according to the Italian Macroseismic Database (DBMI15) (Locati et al. 2022).

Table 2 Maximum MCS intensity associated to each earthquake event
Full size table

2.2 Historic masonry towers

A total of 129 historic masonry towers, affected by the aforementioned events, is collected, grouping the towers by the corresponding earthquake. The groups of the analysed towers are distributed over the Italian territory in Fig. 3. In particular, 31 towers were struck by the 1976 Friuli seismic event, 8 towers by the 1997 Umbria–Marche, 18 towers by the 2002 Molise, 6 towers by the 2009 L’Aquila, 55 towers by the 2012 Emilia, and the remaining 11 towers by the 2016 Central Italy sequence. The distribution of the analysed towers considering the different MCS intensities recorded in the location is shown in Fig. 4 upon a classification in five distinct groups for increasing levels of MCS intensity: Group I includes 40 towers with \(5 \le I_{MCS} < 6\), Group II comprises 29 towers with \(6 \le I_{MCS} < 7\), Group III is composed of 25 towers characterised by \(7 \le I_{MCS} < 8\), Group IV includes 22 towers with 8 \(\le I_{MCS} < 9\), and finally Group V consists of 13 towers with \(I_{MCS}\) equal to or greater than \(9\). It should be noted that whenever the MCS intensity associated with the tower location was not available in the sources listed in Table 2, the intensity associated with the Municipality where the tower is located was assigned.

Fig. 3
figure 3

Geographical locations of towers across the Italian territory, grouped by earthquake event

Full size image
Fig. 4
figure 4

Map highlighting the tower distribution grouped by MCS intensity

Full size image

2.3 Observed damage mechanisms and extents

The damage mechanisms and extents for each historic masonry tower included in the present database are collected from the results of post-earthquake survey activities conducted by field engineers and researchers focusing on masonry churches using the macro-element approach, according to the Italian Guidelines for Cultural Heritage (DPCM 2011). Each macro-element can be activated during the ground shaking leading to different failure mechanisms. According to the EMS-98 (Grünthal 1998), the damage level ranges from 0 up to 5. Therefore, during the post-event survey, a grade for the damage is assigned to each failure mechanism related to the macro-elements. It is worth noting that the number of collapse mechanisms expected for churches, thus considered in the survey form, has increased over time, highlighting the lessons learned from past seismic experiences. Just as an example, 18 failure mechanisms have been considered for the fast damage assessment of masonry churches after the 1997 Umbria–Marche Earthquake, while 28 mechanisms have been taken into account after the 2012 Emilia Earthquake. Historic masonry towers represent one individual macro-element with two expected failure mechanisms, related to their two main structural components, the main body of the tower and the belfry, as shown in Fig. 5. These two distinct failure mechanisms, commonly highlighted as Mec27 and Mec28 in the survey form provided by the Italian Civil Protection, are in this work renamed MecA and MecB, respectively. However, damage to non-structural elements, which may be present in the tower, such as pinnacles or other emerging elements, are not included in the data collection.

Fig. 5
figure 5

Failure mechanisms associated with the tower macro-element: a tower mechanism MecA and b belfry mechanism MecB

Full size image

Although historic towers are apparently regular structures typically with a heavy weight that may be beneficial for structural stability, experience shows that, when subjected to ground shaking, they exhibit damage patterns that can vary widely due to specific structural features (e.g. the slenderness, the presence of openings, the configuration of the belfries, the interaction with adjacent buildings and soil). The aforementioned MecA and MecB, indeed, collect a series of typical failure mechanisms that occurred in historic masonry towers (Fig. 6).

Fig. 6
figure 6

Typical damage mechanisms observed in historic masonry towers

Full size image

In particular, several combinations of damage patterns can affect the main body of the tower. In-plan shear cracking may lead to a rocking mechanism, thus affecting the entire superstructure. The crack patterns and the subsequent formation of a hinge depend on the orientation of the shear plans as they can form a horizontal axis of rotation along one side, converge towards a corner or fragment the panels upon the formation of double-diagonal (X-type) cracks. Very rare is the case of failure due to sliding. Nonetheless, some examples of towers that suffered this type of damage exist, showing remarkable crack patterns but without causing total failure. Cracking in the belfry is likely the most common damage observed in bell towers. In most cases, the damage is due to rocking of the upper part as a consequence of the formation of horizontal cracks in the masonry pillars of the belfry, often involving some detachment of masonry corners in the level below. In some cases, the tower can undergo a shear failure right below the belfry combining the two mechanisms. Other possible combinations of mechanisms depend on the distribution of the openings. Indeed, vertical cracks may occur in the presence of aligned openings and lead to a splitting of the tower.

It is worth noting that the identification of the two aforementioned mechanisms (i.e. MecA and MecB) for masonry bell towers dates back to the seminal studies on the macro-element behaviour of masonry churches, in the aftermath of the 1976 Friuli Earthquake (Doglioni et al. 1994). The masonry towers affected by such severe event and included in the database have been investigated by (Curti et al. 2008), according to this distinction in two mechanisms: the tower and the belfry. For the 1997 Umbria–Marche Earthquake and the 2002 Molise Earthquake, information about the damage to the towers derives from technical reports, (Spence et al. 1998) and (Cifani et al. 2005), respectively. For the following events included in the database, the damage records are extracted from scientific publications. In particular, for the 2009 L’Aquila Earthquake, the main sources are (Augenti and Parisi 2010; Brandonisio et al. 2013; Criber et al. 2015). For the towers affected by the 2012 Emilia Earthquake, relevant information is provided by (Paupério et al. 2012; Sorrentino et al. 2014; Valente et al. 2017; Ferrari 2020). Finally, the damage observed in the towers after the 2016 Central Italy Earthquake is obtained from the works of (De Matteis and Zizi 2019; Giordano et al. 2019; Clementi et al. 2020; Jain et al. 2020; Acito et al. 2021). For around 65% of the identified case studies, the aforementioned sources provide the damage levels, e.g. (Curti et al. 2008; Sorrentino et al. 2014; Criber et al. 2015; Ferrari 2020), while in the remaining percentage (around 35%) the damage level is assigned based on to the visual inspection of the gathered photographic documentation, e.g. (Cifani et al. 2005; Augenti and Parisi 2010; Brandonisio et al. 2013), according to the EMS-98 (Grünthal 1998). These cases are highlighted in the database with the letter c (Annex).

Finally, it should be noted that for a few towers included in the database no damage level is reported associated with one or both of the failure mechanisms (Appendix). This solution is adopted to highlight the cases for which no clear information is obtained from the sources regarding the activation of a mechanism and its severity. These towers belong to the 1976 Friuli and 2012 Emilia earthquakes affecting the entire range of MCS intensities, specifically between 5 and 9. The final number of towers is reported in Table 3, including the subgroups associated with each earthquake, each MCS intensity group and each failure mechanism (MecA, MecB).

Table 3 Number of towers included in the database to which a damage score has been assigned, grouped by earthquake event, damage mechanism and macroseismic intensity
Full size table

3 Damage probability matrices

Following an approach inspired by the Italian Guidelines on Cultural Heritage (DPCM 2011), besides the assessment of the individual mechanisms, a global assessment is conducted. Thus, after the identification of each occurring mechanism and the evaluation of its damage level according to the EMS-98 grading system \(\left( {d_{i} } \right)\), a global damage index for the whole tower is computed \(\left( {i_{d} } \right)\) with a maximum equal to the unit value, based on the following equation:

$$i_{d} = \frac{1}{5} \frac{{\mathop \sum \nolimits_{i = 1}^{2} \rho_{i} d_{i} }}{{\mathop \sum \nolimits_{i = 1}^{2} \rho_{i} }}$$
(1)

where \(\rho_{i}\) denotes the weight to be assigned to each damage mechanism, based on its importance for the structural behaviour of the whole building. In this study, an importance weight equal to the unit value is considered for both damage mechanisms, as commonly done for the global assessment of churches, for instance in (De Matteis and Zizi 2019). An alternative, not pursued here, would be to give more relevance to the tower mechanism, MecA, as the consequences of failure are likely to be greater. Therefore, the expression becomes:

$$i_{d} = \frac{1}{10} \mathop \sum \limits_{i = 1}^{2} d_{i}$$
(2)

The global damage index \(\left( {i_{d} } \right)\) is, thus, the average of the damage levels observed in the activated damage mechanisms. The global damage index provides the final global damage level \(\left( {D_{k} } \right)\) to assign to each tower, according to the correlations provided in (Lagomarsino and Podestà 2004b), as reported in Table 4.

Table 4 Correlation between the global damage index and the global damage level
Full size table

The individual (from mechanisms MecA and MecB) and global damage (for the tower) levels can be sorted to provide Damage Probability Matrices (DPMs), which are common tools to graphically observe the damage distribution within a certain dataset of investigated structures. In particular, when the DPMs are developed for the individual damage mechanisms, the mean damage \(\mu_{d}\) is given by:

$$\mu_{d} = \frac{{\mathop \sum \nolimits_{i = 1}^{n} d_{k,i} }}{n}\quad k = \left[ {0, \ldots , 5} \right]$$
(3)

where \(d_{k}\) is the damage level attributed to either the lower part of the tower or the belfry, depending on the damage mechanism under study, and \(n\) denotes the number of towers investigated. When the DPMs are developed for the entire macro-element, the mean damage \(\mu_{D}\) is given by:

$$\mu_{D} = \frac{{\mathop \sum \nolimits_{i = 1}^{n} D_{k,i} }}{n}\quad k = \left[ {0, \ldots , 5} \right]$$
(4)

where \(D_{k}\) is the global damage level obtained for each tower, and \(n\) denotes the number of towers analysed.

It is worth mentioning that the global damage level is computed only for those towers whose damage levels (ranging from 0 to 5) are available for both damage mechanisms (i.e. 106 towers out of the total 129). When the damage level related to one of the two distinct failure mechanisms is not known (i.e. empty cell in the database), the global damage level is not computed. It is stressed that an empty cell does not indicate “no damage” but indicates lack of information. The DPMs are compared with the estimation obtained using the Binomial Density Probability Function (BDPF) that allows to predict the probability of damage at a given level, \(p_{k}\), based on the knowledge of the mean damage \(\mu\) (or the binomial coefficient \(\mu /5\)) only, as follows:

$$p_{k} = \frac{5!}{{k! \left( {5 - k} \right)!}} \cdot \left( { \frac{\mu }{5} } \right)^{k} \cdot \left( { 1 - \frac{\mu }{5} } \right)^{5 - k} \quad k = \left[ {0, \ldots , 5} \right]$$
(5)

Very similar values are obtained for the mean damage associated to the two damage mechanisms and to the whole tower macro-element: 2.69 and 2.53 for MecA and MecB, respectively, and 2.77 for MecA + MecB. The comparisons between the observed DPMs and the BDPFs are reported in Fig. 7, considering the entire set (i.e. \(5 \le I_{MCS} < 11\)).

Fig. 7
figure 7

DPMs for distinct damage mechanisms (a and b), and for the entire macro-element (c)

Full size image

Although similar values of the observed mean damage are obtained, the distributions of the collected levels of damage are quite different. In particular, for the MecA, a unimodal almost symmetrical distribution of the damage with a peak at level 3 emerges, while MecB exhibits a growing trend in the distribution with a large number of towers presenting a damage score equal to 5. The last case, MecA + MecB, reveals an almost uniform distribution of damage with slightly larger occurrence of D3 and slightly more cases with low damage than severe damage. The displayed damage histograms for MecA and MecA + MecB are relatively well fitted by the binomial distribution, whereas the distribution of MecB seems not to match with the BDPFs.

It is worth noting that the value obtained for the global mean damage is larger than the individual mechanisms mean damage. To better investigate the unexpected trend, the analysis is repeated by taking into account the same number of towers (i.e. 106) instead of the 110 and 112 towers originally considered for MecA and MecB, respectively. Thus, all the towers for which no information is available for one of the two mechanisms are excluded from the calculation. The new DPMs and the BDPFs are illustrated in Fig. 8.

Fig. 8
figure 8

DPMs for distinct damage mechanisms (a and b), and for the entire macro-element (c) considering 106 towers only

Full size image

No significant variation in the mean damage values and corresponding damage distributions emerged and the mean global damage continued to be larger than the mean damage obtained for the single individual mechanisms. This is likely ascribed to the equation used to combine the damage levels into the global damage index (Eq. 2) and to the related correlation with the global damage level (Table 4). The overestimation of the global damage with respect to the average of the individual mechanisms’ values provided by this approach, when adopted for the case of towers only, is likely influenced by the weight assigned to each mechanism and by the correlation used to transform the continuous global damage index into the discrete global damage level. In the present study, both weights and correlation are maintained as in the original approach by (Lagomarsino and Podestà 2004a, b) and in similar works as in (De Matteis and Zizi 2019).

To further increase the reliability of the statistical analysis, the DPMs are also developed taking into account the five aforementioned MCS intensity groups, as shown in Fig. 9a–q. The results show an unexpected increase in the mean damage for the lowest MCS intensity (Group I, \(5 \le I_{MCS} < 6\)), which is higher than the mean damage obtained for Group II (\(6 \le I_{MCS} < 7\)). This result could be due to two main non-mutually exclusive reasons, on one hand, it is possible that some towers, due to their characteristics and conditions, have shown significant damage in an area that has been then characterised by low macroseismic intensity, on the other hand, the way the database has been collected, by relying on published surveys of damaged towers, may have led to a bias in availability of results. In fact, many cases of undamaged towers, which would have reduced the mean damage for a given intensity, may have been overlooked. However, an increase in the observed mean damage with increasing intensity appears from Group II, through Group III (\(7 \le I_{MCS} < 8\)), to Group IV (8 \(\le I_{MCS} < 9\)), which is expected and consistent with similar results obtained from the literature. Group V (\(I_{MCS} \ge 9\)) presents a lower value of mean damage when compared with the mean damage of Group IV. This is likely due to the very small, non-representative, number of towers included in Group V.

Fig. 9
figure 9

DPMs for different MCS intensities

Full size image

The DPMs here developed based on the observational data of earthquake-damaged bell towers can be compared with similar ones found in the literature, as reported in Table 5 (De Matteis et al. 2016; De Matteis and Zizi 2019; Cianchino et al. 2021). Recent studies (Ceroni et al. 2022; Sisti et al. 2023) have reported the DPMs specifically for bell towers hit by the Central Italy seismic event as function of the PGA instead of the MCS intensity, thus, they cannot be directly compared with the other datasets and have not been presented hereafter. Among the investigated studies, (De Matteis et al. 2016) analysed 64 three-nave masonry churches affected by the 2009 L’Aquila earthquake to provide the DPMs for the entire set and for each structural macro-element including bell towers. Their results obtained for the towers as entire macro-elements indicate that 60% of the towers suffered damage. Similarly, (De Matteis and Zizi 2019) developed DPMs based on the damage data of 68 one-nave masonry churches collected during the inspection activities after the 2016–2017 Central Italy seismic sequence. According to their results obtained for the macro-element, structural damage was observed in 70% of the towers.

Table 5 DPMs for towers as entire macro-elements (MecA + MecB)
Full size table

The DPMs of the present database can be also compared with the results of (Hofer et al. 2018; Ruggieri et al. 2022), who specifically developed DPMs for the distinct damage mechanisms. In particular, (Hofer et al. 2018) reported the damage distribution for the belfry mechanism by analysing a database of 196 masonry churches after the 2016 Central Italy Earthquake, whereas (Ruggieri et al. 2022) provided the damage distribution for both damage mechanisms included in a database of 90 one-nave masonry churches inspected after the 2003 Valle Scrivia Earthquake. This was a moderate event that occurred on April 11, 2003, in the Piemonte region, with epicentral coordinates of 44.68 (Lat.) and 8.82 (Long.), hypocentral depth of 26 km and magnitude of 4.9, which caused significant damage to buildings and churches. It is worth noting that such towers and in general this event was not included in the database, because the available information was insufficient. The results are summarised in Tables 6 and 7, including a comparison with the present database, taking into account the damage mechanisms separately.

Table 6 DPMs for tower damage mechanism MecA
Full size table
Table 7 DPMs for belfry damage mechanism MecB
Full size table

Even though the comparison should be carried out according to the distinct macroseismic intensity groups for a better characterisation of the samples, this information is not readily available in the reviewed papers. However, this preliminary comparative analysis is mainly conducted to further investigate the reasons of the aforementioned unexpected trends in the damage levels. From the comparison, it emerges that the present database includes a larger number of towers. In this regard, the existing studies analysed masonry churches, and the exact number of towers for which the DPMs have been developed is not always explicitly specified, although not all the churches present a tower. A few exceptions are worth to be mentioned: in the work of (Hofer et al. 2018) the number of activated belfry mechanisms is reported (i.e. in 17 out of 196 churches), whereas in the work of (Ruggieri et al. 2022) a total percentage of towers within the dataset (i.e. about 70%) is provided without the distinction of the activated damage mechanisms.

Moreover, from the comparison, it is also noteworthy that the present database presents a significant percentage of cases with severe damage (from D3 to D5), and at the same time, a rather low number of towers that are undamaged (D0). This, on one hand, highlights the importance of providing a database that encompasses several events, since a single event may present a narrower range of macroseismic intensity and fewer cases with higher damage extent. On the other hand, the results further stress the aforementioned risk of availability bias. In this light, data collection and analysis must address not only the towers characterised by the occurrence of damage but also those with no damage, since the latter are paramount for an unbiased assessment.

As an example, (Cianchino et al. 2021), analysing the damage data of 87 churches hit by the 2016 Central Italy earthquake, reported a low percentage of mechanism activation, with about 74% of the towers presenting no damage. This result is quite different from similar DPMs provided for the same seismic event. The comparison of distinct works on the effect of this same earthquake demonstrates the need of an adequate definition of representative samples for the analyses. It is also worth noting that the existing studies for the same seismic event may refer to inspections carried out at different times during the earthquake sequence, thus leading to dissimilar results. In this regard, a further source of discrepancy can be linked to the typologies of towers included in the analysis. Indeed, some of the works (De Matteis et al. 2016; De Matteis and Zizi 2019; Cianchino et al. 2021; Ceroni et al. 2022) seem to develop DPMs including various types of towers (e.g. tower, bell cell, and bell gable), differently from the present database which encompasses bell towers only. The DPMs for the two damage mechanisms (MecA and MecB) are instead compared with those reported in the work of (Canuti et al. 2021), who developed the DPMs accounting for distinct groups of macroseismic intensity (see Tables 8 and 9). The comparative evaluation confirmed the discussed data availability bias especially for the lowest macroseismic intensity. It is worth noting that, in some cases (Hofer et al. 2018; Canuti et al. 2021; Ruggieri et al. 2022), the percentages are approximated values, being extrapolated from the graphs.

Table 8 DPMs for the tower mechanism (MecA), grouped by macroseismic MCS intensity
Full size table
Table 9 DPMs for the belfry mechanism (MecB), grouped by macroseismic MCS intensity
Full size table

4 Vulnerability curves

Alternatively to the DPMs, vulnerability curves are a powerful tool to rapidly estimate the expected mean damage, \(\mu_{d}\), for a given level of macroseismic intensity, \(I_{MCS}\). To build a vulnerability model, the function originally proposed by (Lagomarsino and Podestà 2004b; Lagomarsino and Giovinazzi 2006) is usually assumed, as follows:

$$\mu = 2.5 \left[ {1 + tanh\left( {\frac{{I_{MCS} + 6.25V_{I} - 13.1}}{Q}} \right)} \right]$$
(6)

where the mean damage \(\mu\) is correlated with increasing levels of macroseismic intensity \(I_{MCS}\) through the definition of two physical parameters, namely, the ductility index, Q, and the vulnerability index, \(V_{I}\). Q controls the rate of increase of the damage with the intensity, typically defined according to the typological building classification. \(V_{I}\) is generally an average value for a specific building typology that can be further modified by taking into account a series of parameters that influence the structural response of the building under investigation, thus requiring a detailed knowledge about the structure for its complete definition. The abovementioned formulation for the vulnerability curve has been adapted for different building typologies by calibrating the ductility and vulnerability indexes based on real observed data collected from the past occurrence of seismic events. In particular, values for vulnerability and ductility indexes for historic masonry bell towers are recommended by (Lagomarsino et al. 2004; Curti et al. 2008) and have been here used as a reference (Table 10). (Lagomarsino et al. 2004) proposed a vulnerability function for the whole tower macro-element, also recently reported by (Despotaki et al. 2018), while (Curti et al. 2008) developed vulnerability functions specifically addressing both the damage mechanisms of the towers.

Table 10 Vulnerability and ductility indexes for historic masonry towers
Full size table

Figure 10a, b and c show the vulnerability curves obtained by applying these existing equations against the present dataset, namely the curves are compared with the values of mean damage observed in the investigated towers for different groups of MCS intensity, further grouped by seismic events, through a colour key.

Fig. 10
figure 10

Comparison between the existing vulnerability functions and the observed mean damage for distinct mechanisms (a and b) and for the whole macro-element (c)

Full size image

It should be noted that the towers characterised by maximum MCS intensity have been distinguished into two subgroups: 9 \(\le I_{MCS} < 10\) and 10 \(\le I_{MCS} < 11\). The same plots also report the values of the mean damage considering the entire set of earthquakes for each group of MCS intensity as the Global Mean. This value is deemed important to overcome a possible bias given by the different number of towers observed in distinct earthquakes. Nonetheless, it is worth noting that this approach does not completely solve the unexpected high mean damage in the range of 5 \(\le I_{MCS} < 6\) and the extremely low mean damage in the range of 10 \(\le I_{MCS} < 11\). Samples in these ranges may be not fully representative of the expected real population of damaged towers; therefore, for the subsequent calibration of new vulnerability curves, they are discarded.

The existing predictive functions for the individual damage mechanisms seem to fit well with the observed mean damage, especially in the range 6 \(\le I_{MCS} < 10\). A certain scatter between the predicted and the observed damage levels emerges but the same is reduced when the Global Mean values are considered. The predictive function for the entire macro-element, instead, appears to be non-conservative, producing an underestimation of the mean damage.

In light of these considerations, further analyses are conducted to fit the observed mean damage using the well-known formulation of the vulnerability curve but recalibrating its parameters. Initially, to explore a better fit against the mean damage observed in the investigated towers, the three existing vulnerability functions are adjusted by varying the parameters (α, β, Q), as shown in Eq. (11), one by one, while keeping the vulnerability index. This procedure is similar to (De Matteis et al. 2019b).

$$\mu = 2.5 \left[ {1 + tanh\left( {\frac{{I_{MCS} + \alpha V_{I} - \beta }}{Q}} \right)} \right]$$
(11)

The results obtained for the vulnerability curves with a range of parameters are illustrated in Fig. 11a–i. The lower and upper bound for the parameters (α, β, Q) of the three existing formulations and the intervals in their variation, considered for the definition of the curves are reported in Table 11 (here ∆ is the step adopted for the different curves shown). It can be observed that while the variations of the parameters α and β produce a shift of the vulnerability curves, the variation of the ductility index Q generates a change in the slope of the curves. Among the calibration coefficients, α and β play a major role in fitting the curves to the present dataset.

Fig. 11
figure 11

Parametric analyses of the existing vulnerability functions for distinct damage mechanisms, (a–f), and for the whole macro-element, (g–i)

Full size image
Table 11 Constraints adopted in the analysis
Full size table

The parametric analyses performed on the vulnerability equations provide a graphic tool useful to identify the best curve according to specific data in hand. Moreover, a regression analysis is conducted to develop new curves fitting the observed mean damage, by choosing α as a variable to calibrate. It should be noted that the observed mean damage used for data fitting is referred to towers grouped by MCS intensity and earthquake, and the definition of the parameter to calibrate is suggested upon the parametric analyses as shown in Fig. 11, where the effect on the curve of increasing α corresponds to the effect of reducing β and vice-versa. The results related to the individual damage mechanisms (MecA and MecB), Eqs. (7)–(8), show similar formulations to the proposed by (Curti et al. 2008), confirming their capability for damage prediction. When fitting the dataset for the entire tower macro-element (MecA + MecB), instead, the obtained Eq. (10) produces a significant improvement with respect to the original formulation proposed by (Lagomarsino et al. 2004), for the specific dataset in hand. Figure 12 shows the performance plots, namely predicted mean damage versus observed mean damage values, obtained by applying both the existing and proposed formulations, together with the Residual Sum of the Square (RSS) metric. The latter is expressed as follows:

$$RSS = \mathop \sum \limits_{i = 1}^{n} \left( {y_{i} - \widehat{{y_{i} }}} \right)^{2}$$
(12)

where n is the number of observations; \(y_{i}\) denotes the actual observed damage extent, and \(\widehat{{y_{i} }}\) denotes the predicted value estimated via the vulnerability function. This metric shows a significant improvement of the predictive capability upon recalibration of the vulnerability function coefficients, considering both the Global Mean values (plotted as diamonds in Fig. 12), very well predicted, and the mean values for the samples grouped by earthquake (plotted as circles in Fig. 12) that, despite a larger dispersion, are better predicted as well.

Fig. 12
figure 12

Performance plots before (a) and after (b) recalibration. RSS estimated with the Global Mean values

Full size image

To highlight the statistical relevance of the present database and to further validate the existing and proposed vulnerability functions, multiple statistical analyses are performed by randomly generating 10,000 distinct subsets of historic masonry towers. In each subset, 10 towers are randomly removed from the initial whole dataset, thus resulting in smaller subsets of 100, 102, and 96 towers, for MecA, MecB and MecA + MecB, respectively. The Global Mean damage collected from all the subsets are firstly computed and then clustered accounting for each macroseismic intensity as shown in Fig. 13a, c and e. The Global Mean damage of the subsets are characterised by a very narrow interquartile range and, in general, by a rather narrow distribution. Nonetheless, a few significant outliers emerge. It should be also noted that the distribution of the Global Mean damage for high macroseismic intensity features the largest spread. To statistically measure the performance of the existing and novel formulations, the RSS metric is evaluated over all the random subsets, and the obtained distributions are illustrated in Fig. 13b, d, and f. The RSS distributions present lower mean values for MecA than for MecB and MecA + MecB. Finally, the RSS comparative results confirm a significant improvement upon recalibration of the curve for MecA + MecB instance, as their mean decreases and their distribution becomes narrower. A recalibration of the curve coefficients for MecB data has been attempted as well but has led to a negligible improvement of the performance.

Fig. 13
figure 13

Validation of the existing (a–d) and the proposed (e–f) vulnerability functions against 10′000 subsets generated by randomly removing 10 towers each time

Full size image

To further test the robustness of the formulations, randomly generated samples of different sizes are investigated. For the sake of brevity, the results hereafter discussed and reported in Fig. 14 refer to the tower macro-element only, i.e. MecA + MecB, for which the calibration has been proposed. In particular, 3 and 20 towers each time are randomly removed from the initial whole dataset, thus resulting in smaller subsets of 103 and 86 towers, respectively. For the case of 3 towers excluded, all the possible combinations of 3 samples out of the original database are considered, summing up to 192,920 different analysed subsets (Fig. 14a). For 20 towers excluded, the number of all the possible combinations led to an unmanageable computational burden, therefore, 10,000 distinct random combinations are generated (Fig. 14b). The results confirm the rather narrow distribution of the Global Mean damage that tends to spread more with more outliers as the number of removed samples increases and, in general, for higher intensity, likely due to the cardinality of the subsets for each macroseismic intensity level that significantly reduces when more samples are removed.

Fig. 14
figure 14

Validation of the proposed vulnerability function using randomly generated subsets of different sizes: a 3 towers excluded; and b 20 towers excluded

Full size image

These analyses suggest that the database used for the training ensures a sufficiently robust vulnerability curve, in the case of the recalibrated one, and allows to generally conclude that the three investigated curves, for MecA, MecB and MecA + MecB, present an acceptable performance. At the same time, the presence of outliers in the distribution of the Global Mean damage values and the RSS values, whose extent increases when more samples are removed, call for future extensions of the database, fostered by its distribution in open access, and a recalibration of the curves if needed.

5 Fragility curves

To define the fragility curves, the two existing vulnerability functions for the mechanisms MecA and MecB, Eqs. (7)–(8), and the novel proposed formulation for the global behaviour of the macro-element, Eq. (10), can be exploited for the estimation of the expected mean damage for each intensity (\(5 \le I_{MCS} \le 12\)). The expected mean damage values are introduced in the BPDFs, Eq. (5), for the definition of the probability of exceeding a certain level of damage for increasing values of intensity, as in (Lagomarsino and Podestà 2004b). Afterwards, the cumulative density function can be estimated as follows:

$$P_{k} = \mathop \sum \limits_{k = 1}^{5} p_{k}$$
(13)

The fragility curves correlate the probability of increasing damage levels with a pre-selected ground motion parameter, typically the Peak Ground Acceleration (PGA). Although several studies have presented fragility curves considering different parameters including the macroseismic intensity (De Matteis et al. 2019b), in the present work, the PGA is preferred and this is derived based on the existing PGA-I empirical formulation proposed by (Faenza and Michelini 2010), expressed as:

$$I_{MCS} = 2.58 \cdot log PGA + 1.68 $$
(14)
$$log PGA = - 0.65 + 0.39 \cdot I_{MCS}$$
(15)

It should be noted that empirical methods are affected by the quality of input data for model development and, to some extent, by the amount of data. The use of PGA-I empirical correlation demonstrated to be of great potential to extrapolate the fragility curves with a certain level of confidence when PGA is not available. For alternative similar laws the interested reader is referred to (Faccioli and Cauzzi 2006; Gomez-Capera et al. 2020; Masi et al. 2020). The adopted approach is preferred at this stage of the work, given the limited availability of data. Nonetheless, other strategies are recommended for future improvement of the fragility curves, such as the use of more robust methods to estimate the PGA, for instance, based on input parameters such as the magnitude and epicentral distance (Sabetta and Pugliese 1987), or the use of instrumentally recorded PGA nearby the towers.

Figure 15a, and b show the fragility curves developed for the case of damage level associated with distinct damage mechanisms, and Fig. 15c illustrates the fragility curves built for the case of global damage level associated with the tower as an entire structural macro-element. It can be observed that the probability of exceeding a certain level of damage for MecA is generally lower than the probability obtained for MecB, namely the frequency of damage occurrence is higher for the belfry mechanism, confirming the expected significant vulnerability of such elements. Regarding the fragility curves produced for the entire macro-element, they exhibit an intermediate behaviour.

Fig. 15
figure 15

Fragility curves for distinct damage mechanisms (a and b) and for the entire macro-element (c). Here, solid lines are the curves from the Authors’ database and the dashed lines are from (Marotta et al. 2021). The observed damage points are included for a validation purpose

Full size image

After having obtained the fragility curves, the observed damage points are then reported in Fig. 15, upon the correlation of their macroseismic intensity with the PGA computed using the same empirical equation, to discuss the validity of the adopted approach and to further elaborate possible biases included in the collected data. It is worth noting that the estimation of PGA using the empirical equation for the database at hand presents a gap in the range of 0.3–0.4 g. However, the comparisons indicate that the fragility curves for the distinct damage mechanisms obtained from the BPDFs overestimate the probability of the mean damage for D1 and underestimate the probability of the mean damage for D4–D5 in the range of PGA between 0.1 and 0.2 g. In the remaining bin of PGA (0.5 g), the comparison for the tower mechanism MecA presents acceptable prediction, while the fragility curves for the belfry mechanism MecB overestimate the probability of the mean damage along all the possible damage scenario extents (D1–D5). Regarding the belfry mechanism, fragility curves have been already developed by (Marotta et al. 2021) using a fitting procedure based on the lognormal cumulative distribution method over their damage data. These proposed fragility curves, presented in Fig. 15b, underestimate the probability of the mean damage of the present database, especially for D3–D5, when compared with the fragility curves developed in this study. Finally, when comparing the fragility curves with the observed mean damage for the entire macro-element (Fig. 15c), the results show very good agreement.

The development of fragility curves based on BPDFs may explain the poorer accuracy against MecB observed damage, considering the emerged difference between these curves and the DPMs for the case of MecB (see Sect. 3). The application of suitable alternative approaches for developing the fragility curves, directly relying on the observed data fitting, as the lognormal cumulative distribution function (Baker 2015), have been considered but excluded at this stage due to the characteristics of the database. Indeed, a preliminary attempt of using this approach to the present database indicated a need for future expansion of the dataset, once more data is available. This will make it more representative of the real population of towers in seismic-prone areas and will allow to produce more robust predictive models.

6 Conclusions

The present paper collected, upon an extensive literature review, a database of 129 historic masonry bell towers struck by six past seismic events that occurred in Italy, reporting their geographic location, the Mercalli-Cancani-Sieberg (MCS) intensity \(I_{MCS}\) characterising the seismic input, and the observed damage using a macro-element approach. Damage Probability Matrices (DPMs) have been developed for the entire macro-element and for two distinct damage mechanisms that can be activated in historic masonry towers (i.e. lower part of the tower and belfry) and have been compared with the results obtained from previous studies on single seismic events. In general, an expected trend (i.e. increasing observed mean damage with increasing intensity) has been found between the range 6 \(\le I_{MCS} < 10\), whereas a likely availability bias has been detected in the MCS intervals 5 \(\le I_{MCS} < 6\) and 10 \(\le I_{MCS} < 11\).

The database has been used as a testbed for simplified models existing in the literature, providing a valuable rapid vulnerability assessment tool for historic masonry towers at a territorial scale. The results revealed that, while the formulations of vulnerability curves proposed for individual damage mechanisms performed well against the collected data, the formulation proposed for the entire macro-element had a comparatively lower predictive capability. A wide range of curves has been also reported to support future users in fitting specific data in hand. Still, a new calibration of the vulnerability curve coefficients has been carried out to propose an improved formulation for the entire macro-element. The validation of the existing and newly proposed vulnerability curves has been extensively carried out against a large set of samples generated from the initial database by randomly removing 3, 10 and 20 towers. This allowed to explore the statistical relevance of the dataset at hand, and at the same time to better interpret the obtained results.

Finally, fragility curves have been derived by correlating the cumulative probability function, in the form of the Binomial Density Probability Function (BDPF) based on the expected mean damage, with the PGA, retrieved using an existing empirical correlation between PGA and macroseismic intensity (PGA-I). The proposed fragility curves indicated the higher susceptibility to damage of the belfry when considering the distinction between the two failure mechanisms.

The comparison of the DPMs with others existing in literature as well as the comparison of the fragility curves with the observed damage call for a future careful analysis of the representativeness of the datasets used to develop essential rapid predictive models such as the vulnerability and the fragility curves. To this end, an extension of the database once more data is available, especially to include undamaged samples and to increase the number of cases for lower and higher macroseismic intensities, is an essential future scope of the work to ensure that the samples are fully representative of the entire data population. Additionally, the development of more effective predictive models may benefit from the implementation of more robust strategies for the estimation or instrumental measurement of the PGA, combined with the availability of a larger database which will allow the development of fragility curves based on observed data fitting.