Focal mechanism and depth of the 1956 Amorgos twin earthquakes from waveform matching of analogue seismograms

Historic analogue seismograms of the large 1956 Amorgos twin earthquakes which occurred in the volcanic arc of the Hellenic subduction zone (HSZ) were collected, digitized and reanalyzed to obtain refined estimates of their depth and focal mechanism. In total, 80 records of the events from 29 European stations were collected and, if possible, digitized. In addition, bulletins were searched for instrument parameters required to calculate transfer functions for instrument correction. A grid search based on matching the digitized historic waveforms to complete synthetic seismograms was then carried out to infer optimal estimates for depth and focal mechanism. Owing to incomplete or unreliable information on instrument parameters and frequently occurring technical problems during recording, such as writing needles jumping off mechanical recording systems, much less seismograms than collected proved suitable for waveform matching. For the first earthquake, only seven seismograms from three different stations at Stuttgart (STU), Göttingen (GTT) and Copenhagen (COP) could be used. Nevertheless, the waveform matching grid search yields two stable misfit minima for source depths of 25 and 50 km. Compatible fault plane solutions are either of normal faulting or thrusting type. A separate analysis of 42 impulsive first-motion polarities taken from the International Seismological Summary (ISS bulletin) excludes the thrusting mechanism and clearly favors a normal faulting solution with at least one of the potential fault planes striking in SW–NE direction. This finding is consistent with the local structure and microseismic activity of the Santorini–Amorgos graben. Since crustal thickness in the Amorgos area is generally less than 30 km, a source depth of 25 km appears to be more realistic. The second earthquake exhibits a conspicuously high ratio of body wave to surface wave amplitudes suggesting an intermediate-depth event located in the Hellenic Wadati– Benioff zone. This hypothesis is supported by a focal mechanism analysis based on first-motion polarities, which indicates a mechanism very different from that of the first event. A waveform matching grid search done to support the intermediate-depth hypothesis proved not to be fruitful because the body wave phases are overlain by strong surface wave coda of the first event inhibiting a waveform match. However, body to surface wave amplitude ratios of a modern intermediate-depth event with an epicenter close to the island of Milos observed at stations of the German Regional Seismic Network (GRSN) exhibit a pattern similar to the one observed for the second event with high values in a frequency band between 0.05 Hz and 0.3 Hz. In contrast, a shallow event with an epicenter in western Crete and nearly identical source mechanism and magnitude, shows very low ratios of body and surface wave amplitude up to 0.17 Hz and higher ratios only beyond that frequency. Based on this comparison with a modern event, we estimate the source depth of the second event to be greater than 100 km. The proximity in time and space of the two events suggests a triggering of the second, potentially deep event by the shallow first one. Published by Copernicus Publications on behalf of the European Geosciences Union. 1028 A. Brüstle et al.: Focal mechanism and depth of the 1956 Amorgos twin earthquakes

the Santorini-Amorgos graben.Since crustal thickness in the Amorgos area is generally less than 30 km, a source depth of 25 km appears to be more realistic.
The second earthquake exhibits a conspicuously high ratio of body wave to surface wave amplitudes suggesting an intermediate-depth event located in the Hellenic Wadati-Benioff zone.This hypothesis is supported by a focal mechanism analysis based on first-motion polarities, which indicates a mechanism very different from that of the first event.A waveform matching grid search done to support the intermediate-depth hypothesis proved not to be fruitful because the body wave phases are overlain by strong surface wave coda of the first event inhibiting a waveform match.However, body to surface wave amplitude ratios of a modern intermediate-depth event with an epicenter close to the island of Milos observed at stations of the German Regional Seismic Network (GRSN) exhibit a pattern similar to the one observed for the second event with high values in a frequency band between 0.05 Hz and 0.3 Hz.In contrast, a shallow event with an epicenter in western Crete and nearly identical source mechanism and magnitude, shows very low ratios of body and surface wave amplitude up to 0.17 Hz and higher ratios only beyond that frequency.Based on this comparison with a modern event, we estimate the source depth of the second event to be greater than 100 km.The proximity in time and space of the two events suggests a triggering of the second, potentially deep event by the shallow first one.

Introduction
On 9 July 1956 (03:11 GMT), one of the strongest earthquakes in the southern Aegean in the 20th century occurred between the islands of Amorgos and Santorini with a magnitude of M S = 7.4 (Makropoulos et al., 1989).A second strong earthquake (03:24 GMT) with a magnitude of M S = 7.2 occurred just 13 min later in the same region.Both earthquakes were located in the Santorini-Amorgos graben within the central Hellenic volcanic arc (HVA) of the Hellenic subduction zone (HSZ; Fig. 1).The first earthquake was located at the northern flank of the graben southwest of Amorgos (Comninakis and Papazachos, 1986;Makropoulos et al., 1989;Okal et al., 2009).Epicenter locations of the second earthquake vary from the northern graben shoulder southwest of Amorgos (Comninakis and Papazachos, 1986) and the northeastern Santorini-Coloumbo volcanic complex (Makropoulos et al., 1989) to the southern graben shoulder near Anafi (Okal et al., 2009).
The Santorini-Coloumbo volcanic complex is located in the southern extension of the western Santorini-Amorgos graben shoulder.During the last 60 yr the Santorini volcano exhibited only weak background seismicity (Dimitriadis et al., 2005(Dimitriadis et al., , 2009;;Bohnhoff et al., 2006;Hensch, 2009), except for a period of unusually high seismicity lasting from January 2011 to January 2012, with local magnitudes up to 3.2 accompanied by surface deformations within the Caldera as a result of magma upwelling (Newman et al., 2012;Chouliaras et al., 2012;Vallianatos et al., 2013).In contrast, the submarine Coloumbo volcano is characterized Figure 1.Well-located microseismic events observed during the CYCNET (CYClades NETwork) and EGELADOS (Exploring the GEodynamics of subducted Lithosphere using an Amphibian Deployment Of Seismographs) experiments in the Santorini-Amorgos region (Brüstle, 2012).Included are locations of major fault structures in the area according to Perissoratis and Papadopoulos (1999) and Stiros et al. (1994).SAG: southern Amorgos graben; NAG: northern Amorgos graben; CoV: Coloumbo volcano; SaV: Santorini Volcano; AnG: Anidros graben.
by continuous high seismic activity and volcanic unrest.Local seismicity occurs at depths of 3-9 km within a 3 km wide vertical column beneath the northeastern flank of the submarine volcano as a result of magmatic intrusions (Dimitriadis et al., 2005;Bohnhoff et al., 2006;Kolaitis et al., 2007;Dimitriadis et al., 2009;Hensch, 2009).Focal mechanisms of these earthquakes also tend to show NE-SW striking normal faulting mechanisms implying a strong connection to the Kameni-Coloumbo fracture, the southern extension of the western Santorini-Amorgos fault zone.Intermediate-depth seismicity is observed at 120-160 km depth indicating the Wadati-Benioff Zone (WBZ) of the subducting Hellenic slab (Knapmeyer, 1999;Papazachos et al., 2000;Meier et al., 2007;Brüstle, 2012).
The island Amorgos is part of the northern footwall of the graben structure (Fig. 1).Neotectonic investigations on Amorgos infer a tensional NW-SE trending stress pattern (Papadopoulos and Pavlides, 1992) and backtilting (subsidence of the north coast, uplift of the south coast) of the island (Stiros et al., 1994).The southern coastline of the island has been uplifted by about 30 cm as a result of the seismic activity in July 1956.Seismic reflection profiles of the Santorini-Amorgos graben point to a large recent sediment slumping event between Santorini and Amorgos, which is also associated with the seismic activity in July 1956 (Perissoratis and Papadopoulos, 1999).
The shallow origin of at least one of the earthquakes seems to be evident from the reported tsunami waves, which caused damage within a region of about 100 km (Ambraseys, 1960;Stiros et al., 1994).Detailed studies were undertaken  (Ambraseys, 1960;Galanopoulos, 1960;Okal et al., 2009) and studies of sedimentary deposits (Dominey-Howes et al., 2000).
to determine wave heights and wave run-ups on the islands in the southern Aegean.Highest sea wave amplitudes (Fig. 1) were observed on the south coast of Amorgos (20-30 m) and the north coast of Astypalea (10-20 m) by eye witnesses and are seen from marine flood deposits (Ambraseys, 1960;Galanopoulos, 1960;Dominey-Howes, 1996;Dominey-Howes et al., 2000;Okal et al., 2009), whereas on the coastlines of the two islands which faced away from the Santorini-Amorgos graben structure as well as on the northern coastline of Crete, on the Dodecanese islands and the Turkish west coast wave amplitudes of less than 4 m were observed.The only exception is the island Folegandros, where wave amplitudes of up to 14.6 m were observed locally on the western coast.
Both earthquakes caused severe damage on the surrounding islands, especially on Santorini and Astypalea.In total, 53 people were killed, 100 people were injured and 529 buildings were destroyed.More than 3200 buildings were severely damaged.Most of these were on Santorini, where the second event apparently effected the final collapse of houses already damaged by the first event (Ambraseys, 1960;Galanopoulos, 1982).
While the two earthquakes have been unanimously located in the Santorini-Amorgos graben, focal parameters and focal depth are controversial (Table 1) or affected by large uncertainties (Galanopoulos, 1982;Comninakis and Papazachos, 1986;Makropoulos et al., 1989).For the first event, Shirokova (1972) proposed a fault plane solution of a ENE-WSW striking normal fault with a slight dip to the SSE, calculated from first-motion polarities of P phases (Fig. 3); this has been confirmed by stress pattern determination from field observations.Other fault plane solutions indicating a strike-slip solution (Papazachos and Delibasis, 1969;Ritsema, 1974) were judged of poor quality by the authors themselves.The focal mechanism of the second event is essentially unknown.Since the data were recorded on analogue devices and mostly written on paper, the only way to reassess these important earthquakes using modern seismological techniques is by collecting, digitizing and reprocessing historic analogue recordings.Recently, Okal et al. (2009) digitized teleseismic recordings of the first event and inverted them for the spectral amplitudes of mantle Rayleigh and Love waves and obtained a SW-NE striking normal faulting at a focal depth of about 45 km.In the study reported here, analogue seismograms and bulletins of regional stations were collected to refine estimates of focal parameters and focal depth of the two earthquakes.The scanned analogue seismograms were digitized and analyzed for their frequency content using a multiple filter technique according to Dziewonski et al. (1969).The MFT (multiple filter technique) allowed us to identify P , S and surface waves of the first event and to distinguish body waves of the second earthquake from surface wave coda of the preceding one.A grid search was applied to the first event based on cross-correlation and waveform matching of the digitized analogue seismograms with synthetic waveforms to determine the best fitting focal mechanism and hypocentral depth.In addition, a classical focal mechanism determination using the HArdebeckSHearer method (Hardebeck and Shearer, 2002) was carried out based on first-motion polarities taken from the International Seismological Summary (ISS) bulletin (Villaseñor et al., 1997) to clarify nonuniquenesses of the grid search result.
Body wave phases of the second event were masked by the surface wave coda of the preceding earthquake, making it impossible to perform waveform matching.Clues on the focal depth of the second earthquake are gained through comparison of MFT analyses of the historic seismograms with those of a modern event.Constraints on the focal mechanism of the second event were obtained from first-motion polarities published in the ISS bulletin (Villaseñor et al., 1997).

Analogue seismograms and instrument parameters
Based on the ISS bulletin (Villaseñor et al., 1997) of 1956 the two major Amorgos earthquakes were recorded at 317 (first event) and 60 (second event) stations worldwide.In total, more than 80 analogue seismograms from 29 European stations were collected (Fig. 4) to reinvestigate the first and, if possible, the second earthquake of 9 July 1956.Seismograms were collected by the first author and mostly scanned as part of the SISMOgrammi Storici (SISMOS) project (Michelini et al., 2005) or obtained already scanned from the SISMOS database.The collection of analogue seismograms covers dominantly western and northwestern Europe with an azimuthal coverage of 120 • and regional distances of 247 to 2337 km.The digitization of the seismograms was done with a plug-in for the GNU Image Manipulation Program (GIMP) called TESEO (Pintore et al., 2005, Turn the Eldest Seismograms into the Electronic Original ones ).The software allows to apply geometrical corrections of the digitized seismogram (e.g., curvature of the writing arm, skewness of the recording) and to resample the digitized waveform to equally spaced data points.
In addition, bulletins and reports were searched or taken from the SISMOS database to gain detailed information on instrument dimensions and parameters for further processing.This is very important, as seismometer instruments were not standardized at that time.Most observatories had their unique seismograph with individual orientation of the instrument and polarization of the components.To obtain sufficient amplification, the damping of the instruments was either very weak in comparison to today or instruments were recording without any damping (e.g., vertical pendulum at Ebro (EBR) in Spain, Milne seismograph at Malta (MLT)).Also, each observatory adjusted instrument parameters quite frequently to obtain the best resolution on the recording paper (Fig. 5a) At that time, the stations were equipped with mechanical as well as electrodynamic or electromagnetic instruments.The recordings were carried out mechanically or galvanometrically on smoked or photographic paper.Because of the absence of broadband seismometers many observatories were equipped with more than one seismograph, each covering a different frequency band.A very impressing example was the station STU at the Stuttgart observatory (Fig. 5b).During the studied time period about seven seismograph systems were recording simultaneously.
Due to the large magnitude of the first earthquake, high amplitudes, especially of the surface waves, caused clipping of the waveforms or recording needles jumping off the mechanical recording systems, making a considerable number of recordings unusable for further investigations.Moreover, waveforms are often disturbed by time marks, variable paper speed or, in many cases of photo-optical paper recordings, by insufficient exposure time.Finally, missing information on instrument parameters (to calculate the instrument response), on orientation and polarization of the components, or just missing station and instrument name annotations on the seismogram paper are reasons why even good recordings proved unusable.Figure 6 shows a typical recording of the two earthquakes on the vertical component of the 1300 kg Wiechert seismograph at the station in Göttingen (GTT, Germany).Body waves (P 1, S1) and surface waves (sf1) of the first event are clearly visible, while the corresponding phases (P 2, S2, sf2) of the second event are overlain by the surface waves coda of the first event.
For both earthquakes only recordings of the fully mechanical systems of Wiechert or Mainka seismographs at regional distances of about 17 • to 21 • proved suitable for further investigations (Table 2).Table 3 provides a detailed overview of the reasons why other digitized records had to be excluded from the waveform matching procedure.Instrument responses for the usable stations were noted on the seismogram paper or documented in the available bulletins (Table 4).The component polarity of station GTT (marked on the seismogram recording) and COP (Copenhagen, Denmark; Charlier and Van Gils, 1953) was inverse, while the component polarity of station STU (R. Schick, personal communication, 2006) and Jena (JEN, Germany;Charlier and Van Gils, 1953) can be assumed to be correct.Figure 7 shows the digitized waveform sections (black) for the two earthquakes which were suitable for waveform matching.For the first earthquake, they comprise the following phases: P and S phases on all three components of station GTT, on the N component of station COP and the E component of station STU; the P phase only on the E component of station COP and the N component of station STU.For the second earthquake, three of the four digitized waveforms (GTT-Z, GTT-E, JEN-Z) cover the entire time window of body and surface waves, whereas only body waves are available on the N component of station COP.

Instrument responses
In a first step transfer functions were calculated from available instrument parameters.The standard transfer function of a fully mechanical seismograph system is defined by the following expression:  (Gerecke and Güth, 1961).(b) Extract from the bulletin "Seismischer Bericht des Württembergischen Erdbebendienstes" of 1951.At the station STU seven seismograph systems were installed recording on smoked and photo-optical paper (Hiller, 1951).
where d s is the damping constant, ω 0 is the angular eigenfrequency of the instrument and A a frequency-independent calibration factor (Schick and Schneider, 1973).In general, the instrument parameter eigenperiod (T 0 ), damping ratio ( ) and calibration factor (A) are reported in the annual bulletins of the observatories.The damping constant is related to the damping ratio as follows (Galizin, 1914): Figure 8 shows the frequency-dependent normalized transfer functions of two mechanical seismometers.The transfer function of the Z component of the 1300 kg Wiechert seismograph at the station GTT shows a nearly optimal damping and a flat response for frequencies above 0.2 Hz and thus ensures a constant amplification of the recorded frequencies.However, not all instruments were optimally damped back then.An example is the vertical pendulum at the station EBR, which was completely undamped with a resonance at 0.4 Hz.  (3), (4), (5) (4), ( 5

Multiple filter analysis of the analogue seismograms
To identify different phases in the digitized seismograms and also distinguish phase arrivals of the second event from surface wave coda of the first event, we calculated timefrequency spectra with a multiple filter technique (MFT, Dziewonski et al., 1969).The digitized data are fast Fourier transformed (FFT) into the frequency domain and filtered at a series of equidistant center frequencies using narrow Gaussian bandpasses: where ω n is the center frequency of the filter, ω l,n and ω u,n are lower and upper bandpass limits and α controls the resolution of the filter.Note that α is linearly dependent on frequency in order to optimize the time-frequency resolution (Meier et al., 2004).After inverse FFT, the envelope of the filtered waveform is calculated and displayed as a color-coded function of time with the maximum of the envelope normalized to 1. High amplitudes of the envelope are indicated by regions of reddish colors whereas blue colors signify low amplitudes.The time-frequency spectrum of the unfiltered Z component of the 1300 kg Wiechert seismograph at station GTT (Fig. 9a) clearly shows the non-dispersed P and S phases at 4.5 min and 8 min and the dispersion curve of the surface waves at 8-17 min after event time of the first earthquake.The offset of the dispersion curve at frequencies less than 0.02 Hz corresponds to a vertical displacement (clipped recording) of the recorded seismogram at about 12 min.Phase arrivals of the second event are masked in the seismogram by the coda of the surface waves of the first event.However, the signals of P and S phases are identified at 16.5 and 20 min in the time-frequency spectra, while surface waves can not be recognized.The second example shows the unfiltered waveforms of the first and second earthquake recorded by an undamped vertical pendulum at station EBR (Fig. 9b).The amplitude spectrum of this record shows high amplitudes in the range of 0.35-0.45Hz corresponding to the resonance frequency of the calculated transfer function.Although the undamped signal of EBR is not as clear as the one recorded at station GTT, P (5-7 min) and S (8-10 min) phases as well as the dispersion curve (8-16 min) of the surface waves can be identified for the first earthquake.Though the second earthquake exhibits a weaker signal, P (17-19 min) and S (20-22 min) phases can still be identified.However, surface waves are missing as for station GTT, suggesting a hypocenter at greater depths probably located in the subducting Hellenic slab.(1950). 4Written on the seismograms. 5Jensen (1957). 6Charlier and Van Gils (1953).

Focal mechanism and hypocenter depth determination from waveform matching
A grid search based on a waveform matching procedure was attempted in order to obtain best fitting estimates of the focal mechanism and hypocentral depth.The basic idea is to go through all possible hypocentral depths and focal mechanisms characterized by the strike and dip of the fault plane and the slip direction on the fault plane to determine the corresponding moment tensor with unit seismic moment and to calculate synthetic seismograms for this moment tensor and depth.After scaling the observed waveforms to the maximum rms value of all data traces and the synthetic waveforms to the rms value of the corresponding synthetic trace, a misfit over all traces is calculated whose minimum points to the optimal selection of source depth and focal mechanism.
Epicentral coordinates, origin time and earth model are kept fixed during the search.The synthetic seismograms were calculated for the standard 1D-earth model AK135 (Kennett et al., 1995) using the GEMINI (Green function of the Earth by MInor Integration) code (Friederich and Dalkolmo, 1995).The program provides complete seismograms for moment tensor point sources in a spherically symmetric earth model.The source time function is a delta function leading to velocity seismograms for a step function moment tensor source.The synthetics were then filtered in the same way as the data and convolved with the instrument response of the corresponding seismograph calculated from the parameters listed in Table 4. Finally, they were cut to the same time window as the corresponding data traces and interpolated to identical sampling times.Further processing encompassed the removal of the average and of a possible linear trend due to a drift of the instrument as well as bandpass filtering.Since the synthetics are calculated as velocity waveforms, the data were differentiated in case they were recorded as displacements.
A direct matching of the instrument-corrected synthetics with the digitized data turned out to be unsuccessful because of phase mismatches due to either timing errors in the data or deviations of the true earth structure from the standard earth model used.Timing errors are likely for recordings of the 1950s, where, in general, time was noted manually on the seismogram paper from a pendulum clock and time adjustments of the pendulum clock were applied from the telephone or radio time signal once a day (Wielandt, 1996).To compensate for these phase mismatches, we first calculate a normalized cross-correlation function of the processed data and synthetics for lag times ranging from −60 s to +60 s, Then, we identify the time lag associated with the local maximum of the cross-correlation closest to zero lag and time shift the synthetics by this amount.To accommodate amplitude differences between data and synthetics, the data are scaled to the maximum rms value of all data traces (d max ) and the synthetics to the rms value of the corresponding synthetic trace (s max ).Only after these steps, a misfit between data and synthetic trace is calculated as follows: where X k is the misfit for the kth trace, s k (t) is the kth synthetic trace and d k (t) the k − th data trace.The scaling is trace-independent and therefore preserves relative amplitude differences between stations and components.The total misfit is evaluated as the sum of the individual trace misfits.The additional weighting of the trace misfits by the inverse maximum cross-correlation drives the grid search to solutions associated with a high similarity of data and synthetic traces.Once a best fitting focal mechanism and hypocentral depth are found from the grid search, the seismic moment is determined.For this purpose, the energy of all traces is determined by calculating the sum of squares of all waveform samples.The seismic moment is then calculated as the square root of the ratio of the average energy of data and synthetic traces (E d , E s ): The moment magnitude is calculated according to Kanamori (1977) by Rupture area A R and the average displacement D are estimated following Wells and Coppersmith (1994), who proposed empirical relationships between magnitude and rupture parameters for shallow continental intraplate and interplate earthquakes of magnitudes larger than 4.5: lgA R = −2.87(±0.5)+ 0.82(±0.08)M w and lgD = −4.45(±1.59)+ 0.63(±0.24)M w . (8) 6 Results and discussion

First event
The grid search for focal parameters was applied in 30 steps of 12 • for strike and rake, in 30 steps of 3 • for dip and for depths ranging from the surface to 200 km depth.Various low-frequency bandpass filters between 0.02 and 0.1 Hz were tested to reduce the influence of local crustal structures.The best waveform fits were achieved for a bandpass between 0.03 and 0.07 Hz.When experimenting with the grid search, matching both body and surface wave phases lead to unsatisfactory waveform fits.Either a good fit of the body waveforms or a good fit of the surface waveforms could be achieved.A satisfactory fit of both would at least have required a 2-D earth model that reflects the complex crustal and lithospheric structure along the wave paths which cross the Hellenides, Dinarides and Eastern Alps.Since surface waves are much more affected by such lithospheric complexity than body waves, we restricted ourselves to matching P and S wave phases only.Finally, since three slightly differing locations of the earthquake have been proposed in the literature (Comninakis and Papazachos, 1986;Makropoulos et al., 1989;Okal et al., 2009), we performed the grid search for all of them.The smallest waveform misfit achievable for each trial source depth (Fig. 10) varies significantly with depth and exhibits two local minima at 25 and 50 km depth, respectively.2009), (red) Comninakis and Papazachos (1986) and (green) Makropoulos et al. (1989).
Greater source depths can be excluded due to the rapidly growing waveform misfit.As the Moho of the Aegean plate in the central HVA lies at 25-27 km depth (e.g., Papazachos and Nolet, 1997;Bohnhoff et al., 2001), a value of 25 km for source depth appears to be more appropriate (Konstantinou, 2010).An example of a waveform fit for the Comninakis-Papazachos location and a source depth of 25 km is shown in Fig. 11.Owing to the narrow azimuthal coverage and the fact that low-pass filtered data traces are matched, even for the optimal depths of 25 and 50 km, the grid search still yields a multitude of fault plane solutions with acceptably small misfits.To provide a visualization of the variability of these solutions, we have binned the focal sphere and, for each bin, we use color shading to plot the smallest waveform misfit achievable with any fault plane solution whose P axis or T axis pierces the focal sphere within this bin (Fig. 12).P as well as T axes associated with small waveform misfit cluster in the central western part of the focal sphere indicating either normal or thrust faulting mechanisms compatible with the historic waveforms.Solutions with nearly horizontal P (T ) axes spreading around SE belong to solutions with T (P ) axes located in the central western area.Apparently, the historic waveforms are unable to distinguish between normal and thrust faulting, but they require solutions with steeply dipping (> 45 • ) P or T axes roughly pointing west.This behavior is observed for all three locations the grid search was carried out for.
To resolve the question of whether normal faulting or thrust faulting is more appropriate, we collected 42 firstmotion polarities marked as "impulsive" from the ISS bulletin (Villaseñor et al., 1997), and carried out a focal mechanism search using the HASH method as described by  Comninakis and Papazachos (1986) location and a source depth of 25 km.Black traces are data and red traces are synthetics.Amplitude ratio of data and synthetics is correctly reproduced.Hardebeck and Shearer (2002).This method searches the set of all fault plane solutions that explain the polarity data given errors in location, seismic velocity model and polarity.Since the effect of location errors and velocity model on fault plane solutions is hard to estimate, we adopted a heuristic approach by assuming that the location error mostly influences the azimuth of the ray leaving the focal sphere while the seismic velocity model affects the take-off angle at the source.The azimuthal error is assumed to decrease with distance while the error of the take-off angle is assumed to increase with epicentral distance.Maximal assumed errors are two degrees for azimuth and five degrees for take-off angle.HASH finds 545 compatible fault plane solutions which are depicted in Fig. 13 by plotting the piercing points of their P and T axes on the focal sphere.Although the piercing points scatter widely over the focal sphere, the result demonstrates that only strikeslip and normal faulting solutions are compatible with the data.Thrust faulting can be excluded based on the polarity data.Remarkably, the distribution of the P axes on the focal sphere looks similar to the one obtained from the waveform matching grid search.
For the sake of completeness, we also show for the three locations and a depth of 25 km the best fitting fault plane solutions resulting from the waveform grid search (Fig. 14).All of them represent a SW-NE striking normal fault with either a small dextral motion on a SE dipping fault plane The seismic moment M 0 of the these solution ranges between 4.22 and 7.02 × 10 19 Nm.The mean moment magnitude from M w is 7.1.According to Wells and Coppersmith (1994) this corresponds to a rupture  Comninakis and Papazachos (1986) was used for the waveform grid search.Focal sphere is mapped on the horizontal plane using a Lambert-Schmidt projection as commonly done when plotting beach balls.North is directed upwards and east is directed to the right.Color scale is truncated at a waveform misfit of 6.  area A R of about 900 km 2 and an average displacement D of about 1.0 m.Our estimate of moment magnitude is smaller than the magnitude values given by other authors (Table 1).One possible reason could be a directivity effect caused by the very small azimuthal coverage provided by the few usable historic traces.
Among the two possible fault planes, we prefer the dextral SE dipping normal fault mechanism, as all three hypocenter determinations located the earthquake along the northern shoulder of the Santorini-Amorgos graben.This is in agreement with seismic reflection profiles showing a SE dipping normal fault along the northern graben shoulder as well as with neotectonic investigations identifying a NW-SE oriented tensional stress field and a coseismic uplift of the south coast of Amorgos of 0.30 m (Papadopoulos and Pavlides, 1992;Stiros et al., 1994).The fault plane solution shown is also in good agreement with previous focal mechanism determinations of Shirokova (1972) and Okal et al. (2009), and  the local stress field obtained from studies of recent seismicity (Bohnhoff et al., 2006;Friederich et al., 2014).Okal et al. (2009) found a similar focal mechanism, but at hypocentral depths of 46 km.This depth seems to be unlikely, as the Moho has been identified as being not deeper than 20-27 km in the central volcanic arc of the HSZ (e.g., Makris, 1978;Bohnhoff et al., 2001;Sodoudi et al., 2006Sodoudi et al., , 2013;;Tirel et al., 2004;Drakatos et al., 1997;Papazachos and Nolet, 1997;Karagianni et al., 2005;Karagianni and Papazachos, 2007;Endrun et al., 2008).The focal depth of the preferred solution of the present study at about 25 km seems to be more realistic and corresponds to the observed recent seismic activity of the graben region (Bohnhoff et al., 2006;Brüstle, 2012).

Second event
The analogue seismograms as well as their time-frequency spectra exhibit a strong interference of surface wave coda from the first event and body waves from the second event in a time window of 17 min to 22 min after origin time.For this reason, all attempts at matching synthetic waveforms to the data failed.We are therefore also unable to provide an estimate of the seismic moment.Still, the time-frequency spectra (Fig. 9) rather suggest a deep event because surface waves are very weak compared to P and S waves.We want to corroborate this hypothesis in two ways: (1) a comparison of time-frequency spectra of modern shallow and intermediate-depth events in the Aegean recorded by station Moxa (MOX) of the German Regional Seismic Network which is very close to the historic station GTT, and (2) a HASH fault plane solution search based on first-motion polarities available from the ISS bulletin (Villaseñor et al., 1997).
For the comparison of the time-frequency spectra, we selected a shallow event near Western Crete with magnitude 6.1 and a N-S striking sinistral strike-slip mechanism which occurred at a depth of 24 km on 17 March 2004, and an intermediate-depth earthquake which occurred close to the island of Milos on 21 May 2002 at a depth of 84 km with moment magnitude 6.0 and nearly identical sinistral N-S striking strike-slip mechanism.Values of magnitude, depth and focal mechanism were taken from the Swiss Earthquake Service Moment Tensor Catalogue (Braunmiller et al., 2002).Figure 15 shows time-frequency spectra calculated from vertical component data of station MOX for the two events.For the shallow event, the MFT exhibits very low amplitudes of body waves (relative to the surface wave amplitudes) in the frequency band between 0.05 and 0.17 Hz.Beyond 0.17 Hz the spectral amplitudes of body waves become comparable to the surface wave amplitudes.In contrast, for the intermediate-depth Milos event, ratios of body and surface wave amplitude are much higher in the same frequency band between 0.05 and 0.17 Hz. Figure 16 shows a comparison of the MFT spectrum for station MOX of the modern Milos event and that of the second Amorgos event calculated from the record at station GTT.In both spectra, body waves can be well recognized beyond 0.05 Hz, while, similar to the 2002 Milos event, surface wave signals are very weak in comparison.Based on the very weak surface wave amplitudes in the MFT of the GTT record (compared to the Milos event), we suggest that the second Amorgos event was even deeper than the 2002 Milos event, probably located in the Wadati-Benioff Zone of the subducting Hellenic slab, which is located at 120-160 km depth beneath the Santorini-Amorgos graben (Knapmeyer, 1999;Papazachos et al., 2000;Meier et al., 2007;Brüstle, 2012).
This hypothesis is corroborated by a HASH focal mechanism analysis (Hardebeck and Shearer, 2002), as already described for the first event, based on 16 first-motion polarities taken from the ISS bulletin (Villaseñor et al., 1997).Although P and T axes scatter over the focal sphere, the compatible solutions (Fig. 17) indicate a thrusting mechanism rather than a normal faulting or strike-slip mechanism which are characteristic of the shallow seismicity of the Amorgos region.A thrusting mechanism is characteristic for intermediate-depth earthquakes occurring beneath the Cyclades, which preferably show thrust faulting mechanisms with steeply dipping T axes (Friederich et al., 2014).

Conclusions
Refined estimates of the hypocentral depth and focal mechanism of the large Amorgos twin earthquakes of 9 July 1956 (03:11 and 03:24 GMT) were obtained from historic seismograms by a grid search based on the matching of historic, digitized waveforms to corresponding synthetic seismograms.In total, 80 historic records of the two earthquakes were collected either at European observatories or taken from the SISMOS database (Michelini et al., 2005).After scanning the analogue recordings, the TESEO software (Pintore et al., 2005) was used to digitize the seismograms and to correct them for distortions.To obtain instrument responses, information bulletins were searched for relevant parameters.

A. Brüstle et al.: Focal mechanism and depth of the 1956 Amorgos twin earthquakes
Owing to variable recording quality and incomplete documentation of instrument parameters, seven recordings proved suitable for subsequent analysis of the first earthquake while four records could be used for analyzing the second event.High surface wave amplitudes point to a shallow origin of the first earthquake while missing surface waves suggest a deep hypocenter for the second event.These findings are confirmed by a waveform matching grid search for depth and focal mechanism resulting in a SW-NE striking normal faulting mechanism for the first event at a depth of about 25 ± 5 km.The focal mechanism is consistent with both the location of the epicenter on the northern flank of the Santorini-Amorgos graben, and a crustal thickness of less than 27 km of the Aegean crust in the central HVA.The seismic moment of the first event is estimated to be 4.2-7.02×10 19Nm corresponding to a mean moment magnitude of 7.1, a mean slip on the fault of about 1 m and a mean rupture area of about 900 km 2 .Although the Santorini-Amorgos graben is currently characterized by spatially and temporarily varying clusters of microseismic activity confined to the upper crust that are very likely related to fluid flow and the volcanic activity along the fault zone (Dimitriadis et al., 2005;Bohnhoff et al., 2006;Kolaitis et al., 2007;Dimitriadis et al., 2009;Hensch, 2009;Brüstle, 2012;Newman et al., 2012;Chouliaras et al., 2012;Vallianatos et al., 2013), the first of the Amorgos twin earthquakes shows the disquieting evidence for the extreme seismic hazard in this area.
Determination of a meaningful focal mechanism for the second earthquake was only possible using first-motion polarities.Waveform matching failed because of the strong interference of surface wave coda from the first event with body waves from the second one.Nevertheless, amplitude ratios of body and surface waves from the second event suggest a source at intermediate-depths of at least 100 km.This depth range is in agreement with the observed depths of the WBZ of the subducting Hellenic slab beneath the HVA (Knapmeyer, 1999;Papazachos et al., 2000;Meier et al., 2007;Brüstle, 2012).The proximity in time and space of the two events suggests a triggering of the second, potentially deep event by the shallow first one.
Acknowledgements.The study was funded by the German Research Foundation (DFG) within the Collaborative Research Centre 526 "Rheology of the Earth -From the Upper Crust to the Subduction Zone".The majority of figures was generated with the Generic Mapping Tools (GMT) by Wessel and Smith (1995) and the background topography was generated from the ETOPO1 global relief model of Amante and Eakins (2009).We would like to thank Alberto Michelini, Graziano Ferrari and the SISMOS Working group for scanning the seismogram collection at the INGV, for providing an introduction to the TESEO software and providing additionally scanned seismograms and bulletins.We would also like to thank the following persons for providing analogue seismograms, bulletins and the information which was not documented in any bulletin: Ivo Allegretti (University of Zagreb), Manfred Baer (former member of the Swiss Seismological Service SED), Ioannis Baskoutas (National Observatory Athens NOA), Josep Batllo (University of Lisbon, former member of the Ebre Observatory, Spain), Wolfgang Brüstle (Landeserdbebendienst Baden-Württemberg), Wolfgang Brunk (Wiechert'sche Erdbebenwarte Göttingen e.V.), Andreas Tilgner (Georg August University Göttingen), Petra Buchholz (University of Leipzig), Torsten Dahm (German Research Centre for Geosciences GFZ, former member of University of Hamburg), Martin Hensch (University of Iceland, former member of University of Hamburg), Pauline Galea (University of Malta), Klaus-G.Hinzen (University of Cologne), Manfred Joswig (University of Stuttgart), Wolfgang Lenhardt (Zentralanstalt für Meteorologie und Geodynamik, Austria), Nurcan Meral Özel (Kandilli Observatory), Wilfried Steinhoff (Göttingen), Tuncay Taymaz (Istanbul Technical University), Roland Verbeiren (former member of the Royal Observatory of Belgium), Pawel Wiejacz (University of Warzawa) and Rolf Schick (former member of University of Stuttgart).Special Issue: "Rheology of the Earth -observations, laboratory experiments and numerical modelling from the micro-to the macro-scale" Edited by: H. Thybo and W. Friederich

Figure 2 .
Figure 2. Map of locations of the two earthquakes on 9 July 1956 at 03:11 and 03.24 UTC investigated in this study.Also shown are isoseismal lines and tsunami heights (in meters) taken from eye witnesses(Ambraseys, 1960;Galanopoulos, 1960;Okal et al., 2009) and studies of sedimentary deposits(Dominey-Howes et al., 2000).

Figure 3 .
Figure 3. Fault plane solutions for the first earthquake from previous studies.

Figure 4 .
Figure 4. Station map for the 80 analogue seismograms collected from 29 European stations.Seven components of the stations STU, GTT and COP for the first and four components of the stations JEN, GTT and COP for the second earthquake (yellow stars; Comninakis and Papazachos, 1986) were suitable for the present study.

Figure 5 .
Figure 5. (a) Extract of the bulletin of the stations Potsdam (POT), Halle, Plauen and Sonneberg.In 1957, instrument parameters were changed twice at the station POT(Gerecke and Güth, 1961).(b) Extract from the bulletin "Seismischer Bericht des Württembergischen Erdbebendienstes" of 1951.At the station STU seven seismograph systems were installed recording on smoked and photo-optical paper(Hiller, 1951).
, (5) * Both horizontal components on one recording paper.(1) Incomplete waveform recording because of too short an exposure time or bad smoked paper.(2) Waveform clipped because maximum deflection of recording needle reached.(3) Incomplete waveform recording because of end of recording paper.(4) Mixed waveform lines.(5) Incomplete waveform recording because of time mark.(6) Recording stopped because needle jumped off the recording system.(7) Incomplete seismogram paper.(8) Amplification too low.(9) Unclear instrument parameters.(10) No instrument parameters available.(11) Undamped recording.(12) Unknown or unclear component orientation.(13) Component orientation not in N or E. (14) Unknown or unclear seismometer type.(15) Unknown or unclear recording paper type.

Figure 6 .Figure 7 .
Figure 6.Scanned analogue seismogram from the original smoked paper of the first earthquake on 9 July 1956.The Z component of the 1300 kg Wiechert seismograph at the station GTT is at about 2073 km distance from the epicenter.Body waves (P 1, S1) and surface waves (sf1) of the first earthquake are recorded in high quality.The recording is clipped at the largest amplitude of the recording and following signals are recorded with an offset of 5-10 mm.Theoretically calculated phase arrivals (P2, S2, sf2?) show that the waveform of the second earthquake is hidden by the surface waves of the first earthquake.Component polarity is inverse to the real ground motion.("Up" on the seismogram recording represents "down" in real ground motion.)

Figure 8 .
Figure 8. Frequency-dependent normalized transfer functions of a nearly optimally damped (station GTT, red) and an undamped (station EBR, blue) mechanical seismometer.The transfer function of the Z component of the 1300 kg Wiechert seismograph shows a constant plateau for frequencies greater than 0.2 Hz.In contrast, a resonance from 0.18 to 1.26 Hz of the undamped vertical pendulum disturbs the recordings significantly.

Figure 9 .
Figure 9. Digitized waveform (upper part), time-frequency spectra (middle part) and frequency-amplitude spectra (lower part) of digitized unfiltered waveforms of (a) the Z component of a 1300 kg Wiechert seismograph at station GTT and (b) an undamped vertical pendulum at station EBR.Both digitized seismograms clearly show P and S phases and dispersed surface waves from the first earthquake, whereas the body phase arrivals from the second earthquake can only be identified in the time-frequency spectrogram.The absence of surface waves in both recordings from the second earthquake indicates a large hypocentral depth.

Figure 11 .
Figure 11.Example of a waveform fit for the Comninakis and Papazachos (1986) location and a source depth of 25 km.Black traces are data and red traces are synthetics.Amplitude ratio of data and synthetics is correctly reproduced.

Figure 12 .
Figure 12.Smallest misfit achievable for fault plane solutions with (left) P axes or (right) T axes piercing the focal sphere at a given bin.Source depth is 25 km.Location ofComninakis and Papazachos (1986) was used for the waveform grid search.Focal sphere is mapped on the horizontal plane using a Lambert-Schmidt projection as commonly done when plotting beach balls.North is directed upwards and east is directed to the right.Color scale is truncated at a waveform misfit of 6.

Figure 13 .
Figure 13.Compatible fault plane solutions obtained with the HASH method for first-motion polarities of the first Amorgos event taken from the ISS bulletin.Blue and turquoise dots with station names indicate positive and negative polarity, respectively.Red squares mark the piercing points of the P axes and green triangles indicate piercing points of the T axes.The nodal lines represent the average fault plane solution calculated by HASH.

Figure 15 .
Figure 15.Multiple filter technique time-frequency spectra computed for vertical component data of the shallow 2004 western Crete (left) and intermediate-depth 2002 Milos event (right) recorded at station MOX of the German Regional Seismic Network (GRSN).In each subfigure, the top panel shows the data trace, the central panel displays the time-frequency spectrum and the bottom panel shows the amplitude spectrum of the trace.

Figure 16 .
Figure 16.Left: multiple filter technique (MFT) time-frequency spectrum for the second Amorgos event computed from vertical component data of historic station GTT.Right: for comparison, MFT time-frequency spectrum for a modern intermediate-depth event near the island of Milos computed from vertical component data of station MOX of the German Regional Seismic Network (GRSN).In each subfigure, the top panel shows the data trace, the central panel displays the time-frequency spectrum and the bottom panel shows the amplitude spectrum of the trace.

Figure 17 .
Figure 17.Compatible fault plane solutions obtained with the HASH method for first-motion polarities of the second Amorgos event taken from the ISS bulletin.Blue and turquoise dots with station names indicate positive and negative polarity, respectively.Red squares mark the piercing points of the P axes and green triangles indicate piercing points of the T axes.The nodal lines represent the average fault plane solution calculated by HASH.

Table 1 .
Hypocenter parameters of the two earthquakes on 9 July 1956 from previous studies indicate a large uncertainty of hypocentral depth.

Table 2 .
Station Comninakis and Papazachos (1986)ogue seismograms used in the present study.Distances are given relative to the hypocenter location of the first earthquake according toComninakis and Papazachos (1986).

Table 3 .
Recording device, paper style, available component and reasons why records of the twin Amorgos earthquakes were of limited use.

Table 4 .
Instrument parameters that seismograph systems of the recorded analogue seismograms used in the present study.Component polarity in writing direction: left is up (+), right is down (−).