Echoes of ECOs: gravitational-wave signatures of exotic compact objects and of quantum corrections at the horizon scale

Gravitational waves from binary coalescences provide one of the cleanest signatures of the nature of compact objects. It has been recently argued that the post-merger ringdown waveform of exotic ultracompact objects is initially identical to that of a black-hole, and that putative corrections at the horizon scale will appear as secondary pulses after the main burst of radiation. Here we extend this analysis in three important directions: (i) we show that this result applies to a large class of exotic compact objects with a photon sphere for generic orbits in the test-particle limit; (ii) we investigate the late-time ringdown in more detail, showing that it is universally characterized by a modulated and distorted train of"echoes"of the modes of vibration associated with the photon sphere; (iii) we study for the first time equal-mass, head-on collisions of two ultracompact boson stars and compare their gravitational-wave signal to that produced by a pair of black-holes. If the initial objects are compact enough as to mimic a binary black-hole collision up to the merger, the final object exceeds the maximum mass for boson stars and collapses to a black-hole. This suggests that - in some configurations - the coalescence of compact boson stars might be almost indistinguishable from that of black-holes. On the other hand, generic configurations display peculiar signatures that can be searched for in gravitational-wave data as smoking guns of exotic compact objects.


I. INTRODUCTION
The relativistic collision of two compact objects is the Rosetta Stone of the strong-gravity regime. The gravitational-wave (GW) signal emitted during the process contains a wealth of information on the nature of the colliding bodies. Following the recent LIGO detections [1][2][3], in the next years GW astronomy will deepen our understanding of the gravitational interaction and of astrophysics in extreme-gravity conditions to unprecedented level, playing a role similar to that of atomic spectroscopy in advancing quantum theory during the past century.
The comparison to atomic spectroscopy seems particularly apt at least in two respects: (i) the post-merger ringdown phase is governed by a series of damped oscillatory modes [4][5][6] that can be computed very precisely in perturbation theory, and are akin to the energy levels of the hydrogen spectrum; (ii) the precise modelling of the gravitational waveform allows us to search for smokinggun anomalies due to new physics, similarly to the celebrated Lamb shift in atomic spectroscopy.
GW spectroscopy will play an increasingly important role as more and more events at large signal-to-noise ratio are detected. These observations provide novel ways to test strong gravity [7][8][9][10], black-hole (BH) no-hair results [11], the existence of event horizons [12], possible quantum effects at the horizon scale [12,13], dark matter and environmental effects [13,14], and also exotic compact objects (ECOs) which might reveal themselves for the first time in the GW band [12,15,16]. All these tests require a precise modeling of the gravitational waveform in strong-gravity processes.
It has been recently argued [12] that the post-merger ringdown phase of an ECO in the high-compactness limit is initially almost identical to that of a BH, and that any correction at the horizon scale due to a surface [17] or to quantum effects [18][19][20][21] will reveal itself in secondary pulses that appear in the late-time ringdown waveform. This result was obtained by studying the radial plunge of a test particle into a thin-shell wormhole. If the wormhole throat is located at r 0 ∼ 2M , the initial ringdown signal is due to the vibration modes of the photon sphere (PS), and is the same as those of BHs, even though their quasinormal mode (QNM) spectrum -defined as the poles of the relevant Green's function [5,6] -differ dramatically. BHs have QNMs which can be identified with the PS. There being no other scale in the problem and with ingoing conditions at the horizon, the PS modes are identical to the QNMs and no other mode is excited. For ECOs, on the other hand, the PS modes still exist and they ring in the same way as BHs, but are not QNMs, as they do not belong to the spectrum of the relevant operator anymore. Instead, the spectrum contains a series of trapped modes, which describe the vibration of the inner stable PS, which is absent in BH spacetimes.
With the exception of BSs -which we know to form naturally as a consequence of collapse of scalar fieldsthere is no known formation mechanisms for ECOs. In addition, there are some indications that horizonless, ultracompact objects (what we have termed as ECOs), are linearly or nonlinearly unstable [22] although the timescales involved are model-dependent or unknown. Nevertheless, the time has come to expect the unexpected, and a good understanding of alternatives allows us to search for new physics in GW data. The understanding of ECOs is also important to quantify our confidence in the existence of BHs and event horizons.
The purpose of this work is to extend the analysis of Ref. [12] in different independent directions. On the one hand, we consider different models of ECOs, generic trajectories of test particles, as well as the scattering of Gaussian wavepackets off these objects, showing that a generic feature of microscopic-scale corrections near the horizon is the presence of a modulated series of "echoes"of the PS vibration modes. On the other hand we study, for the first time, the head-on collision of two solitonic boson stars (BSs) with a self-interacting potential including terms up to sixth order in the scalar field [23,24]. The latter are chosen because they can reach a compactness comparable to that of the PS, are relatively easy to evolve numerically, and can naturally form in dynamical scenarios [25]. To the best of our knowledge, no other model of ECO is known to enjoy all these properties. Thus, solitonic BSs stand out as the most natural model of ECOs and an important question is whether they can mimic the GW signal of a BH-BH coalescence. Through this work we use c = G = 1 units.

II. ECHOES OF ECOS
In this section, we investigate the (ringdown) response of several models of ECOs in different scattering processes. Most of our results are derived for the same wormhole model studied in Ref. [12], for thin-shell gravastars [17,26], and for a simple toy model of an empty, spherical thin shell of matter (model II of matter-bumpy BH in Ref. [13] with M = 0). The list is not meant to be exaustive, but merely to show that, qualitatively, the response of ultracompact objects is universal and simple, regardless of the specifics of the object. All these models are characterized by some exotic form of matter that prevents gravitational collapse and, most importantly, have a radius r 0 that can be arbitrarily close to the would-be Schwarzschild radius. We focus on models in which with M , which can qualitatively describe putative microscopic corrections at the horizon scale [cf., e.g., Refs. [17][18][19][20]27] for some proposals].
In the spherically symmetric case, the line element for these models can be collectively written as where F and B depend on the model [12,13,22]. Matter is localized only in the region r ≤ r 0 whereas, in the region r > r 0 , Birkhoff's theorem guarantees that spherically symmetric solutions are described by the Schwarzschild metric, F (r) = B(r) = 1 − 2M/r. Details on each model are given in Appendix A.

A. Scattering of wavepackets
The most relevant signatures of these models are already evident in the simplest scattering process, namely a test scalar wavepacket being scattered off the gravitational potential of the ECO. The scattering is governed by the (Klein-Gordon) master equation [6] 2Φ = 0. Using angular variables (θ, φ) on the sphere, and expanding the scalar in spherical harmonics as Φ = where dr/dr * = √ F B. In the exterior region, r > r 0 , the tortoise coordinate r * and the potential V l coincide with their Schwarzschild values, whereas their expressions in the interior region r < r 0 are model-dependent [cf. Appendix A]. The potential for some representative cases is shown in Fig. 1.
We solve Eq. (3) with initial conditions The waveform obtained by solving Eq. (3) with r g = 10M, σ = 6M is shown in Fig. 2 for various models and compared to the BH case.
As discussed in Ref. [12], the initial ringdown signal is basically identical to that of a BH. This part of the signal corresponds to unstable modes at the outermost PS, which is present when r 0 < 3M , and which are associated with the maximum of the potential V l near r ∼ 3M . The outermost PS is typically unstable on short timescales, explaining the rapid damping of this ringdown stage.
On the other hand, when the Schwarzschild horizon is replaced by a surface (as, e.g., in the gravastar case) or by a throat (as in the wormhole case), the potential also develops a minimum (i.e, an innermost stable PS) which can trap low-frequency modes [12,15,[28][29][30] (cf. Fig. 1). This inner PS can also be thought of as being caused by the centrifugal barrier, and it may become nonlinearly unstable [12]. These modes make their way to the waveforms in Fig. 2 in the form of "echoes" of the initial PS modes after they leak through the potential barrier: the radiation pulse generated at the potential barrier peak (the PS modes) is then trapped in a semi-permeable cavity bounded between the two PSs. Indeed, the time delay between two consecutive echoes is roughly the time that light takes for a round trip between the potential barrier. In general, this delay time reads  [12] and of star-like ECOs with a regular center [22]. The precise location of the center of the star is modeldependent and was chosen for visual clarity. The maximum and minimum of the potential corresponds approximately to the location of the unstable and stable PS, and the correspondence is exact in the eikonal limit of large angular number l. In the wormhole case, modes can be trapped between the PSs in the two "universes". In the star-like case, modes are trapped between the PS and the centrifugal barrier near the center of the star [28][29][30]. In all cases the potential is of finite height, and the modes leak away, with higher-frequency modes leaking on shorter timescales.
where r min is the location of the minimum of the potential shown in Fig. 1. If we consider a microscopic correction at the horizon scale ( M ), then the main contribution to the time delay comes near the radius of the star and therefore, where n is a factor of order unity that takes into account the structure of the objects. For wormholes, n = 8 to account for the fact that the signal is reflected by the two maxima in Fig. 1, whereas for our thin-shell gravastar model and the empty-shell model it is easy to check that n = 6 and n = 4, respectively. The results shown in Fig. 2 for = 10 −6 M are perfectly consistent with this picture, with the wormhole case displaying longer echo delays than the other cases with the same compactness.
Our results show that the dependence on is indeed logarithmically for all the ECOs we studied. As argued in Ref. [12], the logarithmic dependence displayed in Eq. (6) implies that even Planckian corrections ( ≈ L P = 2 × 10 −33 cm) appear relatively soon after the main burst of radiation, so they might leave an ob-servable imprint in the GW signal at late times. From Eq. (6), a typical time delay reads where M 30 := M/(30M ). The picture of GW signal scattered off the potential barrier is also supported by two further features shown in Fig. 2, namely the modulation and the distortion of the echo signal. In general, modulation is due to the slow leaking of the echo modes, which contain less energy than the initial one. In the wormhole case, this effect is stronger due to the fact that modes can also leak to the "other universe" through tunneling at the second peak of the potential. While the amplitude of the echoes is model-dependent, for a given model it depends only mildly on . Distortion is also due to the potential barrier, which acts as a low-pass filter and reflects only the low-frequency, quasibound echo modes. This implies that each echo is a low-frequency filtered version of the previous one and the original shape of the mode gets quickly washed out after a few echoes 1 .

B. Waves generated by infalling or scattered particles
The features above are observed in a simple scattering process, but are also evident in the GW signal produced by head-on collisions or close encounters, in the testparticle limit. The latter differ from the radial plunge studied in Ref. [12] in that their pericenter r min > 3M , i.e. the particle does not cross the radius of the PS (in fact, scattered particles in the Schwarzschild geometry can never get inside the r = 4M surface). In order to compute the GW signal, we use the Regge-Wheeler-Zerilli decomposition reviewed in Appendix B (cf. Ref. [31] for details).
We have studied the GW emitted during collisions or scatters between point particles and ECOs; again the general qualitative features are the same as those discussed in Section II A and independent of the nature of the ECO. To be specific, we show in Fig. 3 the Zerilli wavefunction for a point particle plunging into (left panel) or scattering off a wormhole with = 10 −6 M , with initial Lorentz boost E = 1.5. The coordinate system we use is such that the particles are moving along the equator, and it differs -by a π/2 rotation -from the coordinate axis used in Ref. [12]. As such, the l = 2 Zerilli-Moncrief wavefunction, for example, has contributions from azimuthal numbers m = 0, ±2. Note also that it is easy to The waveform for the radial infall of a particle with specific energy E = 1.5 into a wormhole with = 10 −6 M , compared to the BH case. The BH ringdown, caused by oscillations of the outer PS as the particle crosses through, are also present in the wormhole waveform. A part of this pulse travels inwards and is absorbed by the event horizon (for BHs) or then bounces off the inner (centrifugal or PS) barrier for ECOs, giving rise to echoes of the initial pulse. This is a low-pass cavity which cleans the pulse of high-frequency components. At late times, only a lower frequency, long-lived signal is present, well described by the QNMs of the ECO. Right panel: the same for a scattering trajectory, with pericenter rmin = 4.3M , off a wormhole with = 10 −6 M . The main pulse is generated now through the bremsstrahlung radiation emitted as the particle approaches the pericenter. The remaining main features are as before. We show only the real part of the waveform, the imaginary part displays the same qualitative behavior.
express these results in a rotated frame [32,33], and we checked that the waveforms agree up to numerical errors with our previous study [12] 2 .
The left panel of Fig. 3 shows the l = m = 2 GW waveform generated by point particle plunging radially into a BH and a wormhole with = 10 −6 M . The main features were already shown before [12], but they are much clearer here: as the particle crosses the outer PS, a pulse of radiation is emitted. This pulse, identical for all ultracompact ECOs, has both an outgoing and ingoing component. When the central object is a BH, the ingoing pulse disappears towards the event horizon and the outgoing pulse is all that an outside observer receives. When the central object is an ECO, the ingoing pulse is now trapped between the outer and inner PSs and bounces back and forth, "echoing" through the cavity. Because the pulse is of relatively high frequency, there is leakage at each bouncing, and outside observers are able to detect many echoes. At late times the echo acquires a smaller frequency component, identical to the QNMs of the ECO.
This same process is triggered by more generic orbits, even by particles which do not cross the outer PS, but that approach it sufficiently closely. This time however, it is the bremsstrahlung radiation that excites the outer PS modes. An example is shown in the right panel of Fig. 3 for a pericenter at r min = 4.3M . The remaining steps are the same as before: echoes are observed, delayed by an amount that depends exclusively on , or in other words, on how big the cavity is. The amplitude of the echoes does not seem to be sensitive to but only to the details of the process exciting the main pulse.
In other words, the angular momentum of the infalling particle does not influence significantly the echo structure, and we therefore expect it to be a generic feature of ECOs: even particles on the last stages of merger will excite echoes as the they plunge through the PS.
bottom right panel of Fig.4 in Ref. [12] refers to a Lorentz factor E = 1.01 and not to E = 1.5 as reported in the paper. This has since been corrected in an Erratum. FIG. 4. GW signal produced by a test particle falling radially into a wormhole with E = 1.01. We consider the same setup as in Ref. [12] but for a wormhole without a PS. Without outer (and inner) PS, the ringdown signal is, clearly, different from that of a BH. Because there is no longer a good resonating cavity, echoes do not appear to be excited.
Finally, it is clear that a crucial ingredient for the appearance of echoes in the GW signal is the presence of a PS in the spacetime as well as a sufficiently large "cavity". As shown in Fig. 4, ECOs without a PS display a ringdown different from that of BHs even at early times. Furthermore, because of the absence of trapped states in the spectrum [12,22] (which, in turn, is due to the absence of a potential well, or an inner PS), the late-time ringdown is simply characterized by a damped sinusoid, without the echo structure.

III. HEAD-ON COLLISIONS OF COMPACT BOSON STARS
In this section we go beyond the point-particle limit previously considered to study the head-collision of two equal-mass ultracompact objects which are initially at rest. First numerical studies in this scenario [34,35] already showed the solitonic behavior of boson stars. Later studies have investigated the outcome of highly relativistic collisions between BSs [36] in the context of the socalled Hoop Conjecture (crucial for the trans-planckian collision problem). The orbiting case was considered within the conformally flat approximation, which neglects gravitational waves, in [37]. Studies concerning gravitational waveforms, directly motivated by GW science, include Refs. [32,38]; nevertheless, studies of collisions of BSs aimed at understanding how well their signal can mimic BHs, are missing. As discussed in the introduction, solitonic BSs are a natural candidate for this purpose.
BSs are equilibrium, self-gravitating solutions of the Einstein-Klein-Gordon theory with a minimally-coupled, complex scalar field (cf. Ref. [25] for a review), The Einstein equations read G ab = 8πT Φ ab , with whereas the Klein-Gordon equation is Φ = dV d|Φ| 2 Φ, together with its complex conjugate. We consider solitonic BSs supported by the self-interacting potential [23,24] The mass m S of the scalar is related to the mass parameter µ above through µ = m S / . Here, σ 0 is a constant, generically assumed to be of the same order as µ. This is the simplest potential that can support, in the absence of gravity, nontopological solitonic solutions [23,24]. Solitonic BSs can be very compact, with the minimum radius of stable, spherically-symmetric configurations comparable to the radius of the PS 3 , i.e. R ≈ 3M [15]. Their maximum mass reads [23,39] 3 Strictly speaking, BSs extend all the way to infinity; however, the scalar energy density decreases exponentially at large distances and it is standard practice -which we follow -to define its radius as the point at which 99% of the BS mass is contained.
where the scaling of M max with µ is exact, while the scaling with σ 2 0 is only approximate and valid when σ 0 ∼ µ m P (with m P the Planck mass). The field equations for the static, spherically symmetric case are given in Appendix C, while the numerical setup of the simulations is briefly described in Appendix D. A representative mass-radius relation for solitonic BSs is shown in Fig. 5.
The field equations for the solitonic potential (10) are stiff, and the scalar field has a very steep profile across a surface layer of thickness ∼ µ −1 . This stiffness makes the numerical integration particularly challenging. Here, we use the method presented in Ref. [15] to prepare initial data for the spherically symmetric case. The initial state for the binary head-on collision is simply constructed by superposition of two solitonic BSs. The simulations presented below refer to BSs with the same mass and radius, initially at rest and separated by a distance ≈ 2.7R; we have also tried different configurations finding qualitatively similar results.
The theory (8) is symmetric under Φ → −Φ and under a U (1) transformation. By using these symmetries, we can straightforwardly change the sign of the Noether charge and/or the phase of the scalar field in either of the two BSs in the initial data. We will present results for three different configurations covering the most extreme interaction dynamics between BSs: (1) a binary of two identical BSs in phase (bs-bs), (2) a binary with a BS and an anti-BS with opposite Noether charge (bs-abs), and (3) a binary with two BSs in phase opposition (bsbsop) corresponding to a shift of π in the phase of one of the stars. Notice that these cases have also been studied in the context of binary mini-BSs [32,38].
The gravitational waveforms of these three configurations are presented in Fig. 6 for two representative values of the total mass. In the left panel, we consider the evolution of two relatively light BSs with M/R ≈ 0.118 (red marker in Fig. 5), whose head-on collision could form in principle a BS with nearly maximum mass and with M/R ≈ 1/3 (blue marker in Fig. 5). In all the cases considered here the waveforms display a qualitative behavior different from that of the head-on collisions of two (nonspinning) BHs with the same masses. We argue that this is mainly due to three reasons. First, in order to form a BS which does not exceed the maximum mass in Fig. 5, the compactness of the initial BSs has to be relatively low, so they start merging much before the BH case. Secondly, the radius of the final BS (if it eventually forms after a transient stage) roughly coincides with the PS, and does not meet the requirements to mimic well a BH ringdown at early times [12] (as also shown in Fig. 4). Finally, the scalar fields composing the BSs interact in a non-trivial way during the merger which depends on its specific configuration (i.e., phase shift and Noether charge sign). This behavior manifests more clearly in the maximum of the scalar field norm, which is displayed in Fig. 7. The only scenario with a clear merger and relaxation to another BS configuration is the (bs-bs) case. The opposite Noether charges of the (bs-abs) case annihilate during the merger, dispersing and radiating all the scalar field. The (bs-bsop) case is probably the most exotic scenario, since the scalar field interaction induces a repulsive force. Therefore, the stars suffer several inelastic collisions, bouncing back and forth, before losing all their kinetic energy. At late times the two stars are at rest, next to each other, without merging. In the right panel of Fig. 6 we show the case of two BSs with moderately high compactness, M/R ≈ 0.184 (green marker in Fig. 5). For all the configurations considered, the final product of the merger is a Schwarzschild BH. However, the relative phase of the scalar field and the BS charges have again a dramatic impact in the waveform. When the two BSs have initially the same phase, their waveforms are qualitatively different from that of two colliding BHs, regardless the relative sign of the BS charges. In particular, high frequency oscillations, resulting from the interaction of the BSs scalar fields, appear soon after the merger. The resulting star promptly collapses to a BH, producing the well-known ringdown signal after the merger. However, when the initial BSs have opposite phases the interaction between the scalar fields does not produce any interference pattern and the signal is much more similar to that of a BH. In this latter case, we expect that the small differences arise from the compactness of the initial BSs (which is anyway smaller than the BH case roughly by a factor 2.7) and by the related fact that the mass of the final BH is not exactly the same in the two cases.
Finally, notice that further increasing the total mass of the system will result in a faster collapse of the final object into a BH, reducing the anomalous signatures of the scalar fields interference on the emitted GWs.

IV. DISCUSSION AND CONCLUSION
The recent direct detection of GWs [1,2] has opened two intriguing opportunities related to ECOs, namely constraining these objects as alternatives to BHs and using them as a proxy to probe quantum corrections at the horizon scale. In this paper, we have explored some GW signatures that emerge in both scenarios, finding a number of interesting results.

A. Quantum corrections
First of all, we have extended and clarified the picture proposed in Ref. [12]. The ringdown signal of an ECO whose compactness is parametrically close to that of a BH displays some universal features. The main burst of radiation during the early-time ringdown phase is only associated with the vibration modes of the PS and does not depend on whether or not the spacetime has a horizon. The pulse of radiation at the PS travels unimpeded towards the event horizon, if the object is a BH. Thus the late stages in the dynamics of BHs are simple. By contrast, ECOs have an outer and inner PS: they represent a cavity for perturbations, able to trap them and leak them away on very large timescales. An outside observer will see characteristic imprints of the absence of horizons under the form of distorted echoes of the original pulse, which can live orders of magnitude longer than the timescales usually associated to BHs. The scattering The bs-abs configuration has opposite Noether charges that annihilate soon after they merge, destroying the stars and dispersing/radiating their scalar fields. The scalar field interaction in the bs-bsop is repulsive and larger than the gravitational attraction, so the system undergoes several inelastic collisions -which compresses the star and leads to small bumps on the scalar field norm that can be observed at t ≈ {60, 160, 240}-before relaxing to a binary at rest with the surfaces barely touching.
of wavepackets or of point particles confirms this picture and shows that the modulation of the ringdown signal is associated with the reflection of the main burst of radiation off the potential barrier at the PS, producing a characteristic train of "echoes" at interval ∆t. For a given model, the amplitude of the echoes depends only mildly on the compactness in the → 0 limit. Reflections of the potential barrier give also rise to a distortion of each echo mode, since high-frequency components are filtered out. This is an unfortunate feature for GW searches, because it implies that the signal is quickly washed out after a few echoes. In this context, it will be interesting to find an analytical template for the late-time ringdown waveform, in order to search for these echoes in actual GW data through matched filters. It is also interesting to study the echo structure in the presence of rapidly rotating objects. This will presumably give rise to rich structure in the echoes, since the lifetime of the main burst generated at the PS can be much longer for spinning BHs.

B. Hairy black-holes
The discussion above also sheds some light on the issue of looking for hairy BHs with GWs [11]. Some -if not all -of these solutions are associated with an extra scale in the problem. For example, a minimally coupled massive scalar field theory gives rise to non-trivial spinning hairy BHs, which describe a scalar "cloud" outside the horizon [40], the spatial extent of which is linked to the mass scale of the field. What our results teach us is that the GW response of such geometries may be identicalat early times -to those of Kerr, if their near-horizon geometry is sufficiently close to it: the early time response depends mostly on the properties of the PS. Only at late times will the effect of a different geometry or environment become noticeable (see also the overview [13]). The lesson is therefore that more sensitive detectors are needed to probe the late-time behavior of the dynamical response of hairy BHs, the bonus being that GR will be also tested during the process.

C. Boson star strawmen
Having investigated the ringdown signatures of microscopic deviations at the horizon scale, the complementary problem of ECOs as BH mimickers is also interesting. BHs enjoy two remarkable features. Not only does their compactness exceed that of neutron stars, quark stars and BSs, but they also have an arbitrary mass. The former feature implies that two BHs typically merge when they are just a few Schwarzschild radii apart, whereas the latter feature implies that the merger product is a (stable) BH whose mass is roughly the total mass of the binary (modulo GW energy loss).
In light of our results, this simple consideration shows that the limitation of ECOs are more theoretical than phenomenological. Contrived models like wormholes and gravastars can be as massive and compact as BHs and can therefore mimic the inspiral phase up to the merger. However, their dynamics in a comparable-mass, twobody collision and their formation in dynamical processes are difficult to study. The formation and dynamics of BSs, on the other hand, are relatively simple and well established, but these models are limited by their relatively low maximum compactness and maximum mass. Static BSs seem to be considerably less compact than a Schwarzschild BH, and only very fine-tuned models (marginally) possess a PS [15]. This seems to be associated with the finite compressibility of a scalar field, although it would be interesting to find a general bound "à la Buchdahl" [41] for generic BSs. Thus, while BSs are viable ECOs, they are of less interest to mimic quantum corrections at the horizon scale and to test the findings of Ref. [12] in the comparable-mass regime.
Nonetheless, we find that -in some configurationsthe collision of solitonic BSs can mimic that of two BHs. In particular, two BSs with opposite phase and relatively large compactness produce an initial waveform which is similar to that of two BHs with the same mass, but they merge to form a BH. On the other hand, if the final object is a BS, the initial compactness needs to be sufficiently small and this produces qualitative differences in the pre-merger waveforms. Furthermore, for generic configurations, the GW signal is markedly different from that of a pair of BHs and even the outcome of the merger does not need to be either a BH or a BS.
A natural extension of our work is to study the quasicircular coalescence of two solitonic BSs and to check whether the picture that emerges from the head-on collisions remains valid, especially for what concerns the role of the initial phases of the BSs and the existence of configurations that can mimic the entire BH-BH coalescence waveform. In this context, an important discriminator is provided by the tidal deformability that enters the inspiral waveform at fifth post-Newtonian order (cf., e.g, Ref. [42] and Ref. [43] for a review). Work on the tidal deformability of BSs is underway (see also the recent Ref. [44]). For the case of gravastars, the tidal Love numbers vanish in the BH limit [45,46] which suggests that some -but not all -models of ECOs can be constrained by a GW measurement of the tidal deformability.
Finally, although the effects discussed here might be rather exotic, we advocate a proactive view: not only has GW astronomy the potential to constrain these models, but there is the exciting prospect for novel, unexpected detections.

ACKNOWLEDGMENTS
We benefited from discussions with all the participants of the Unifying tests of General Relativity workshop at Caltech. We are grateful to Niayesh Afshordi for useful correspondence and suggestions. V.C. and S.H. acknowledge financial support provided under the European Union's H2020 ERC Consolidator Grant "Matter and strong-field gravity: New frontiers in Einstein's theory" grant agreement no. MaGRaTh-646597. Research at Perimeter Institute is supported by the Government of Canada through Industry Canada and by the In this appendix we provide some details on the models of ECOs used in the main text. All models are spherically symmetric and described by the line element (2). Some models are discontinuous across the surface and require a thin shell of matter at r = r 0 . In this case, Israel's junction conditions [47] relate the discontinuities of the extrinsic curvature on the surface with the stressenergy tensor of the thin layer. From these conditions, the surface energy Σ and surface pressure p of the shell read [26] denotes the discontinuity of a generic function A(r) across the shell.

A toy model: empty thin shell
The simplest model that displays ringdown echoes is an empty thin shell of matter located at r = r 0 . The line element is Eq. (2) with Equations (A1) imply Σ > 0 for any r 0 > 2M , whereas the dominant energy condition on the shell, |p| ≤ σ, implies r 0 ≥ 25M/12 ≈ 2.08M . Note that this solution is unstable against radial perturbations when r 0 2.37M [48].
In order to avoid a discontinuity of the metric at the shell, here we consider a slightly different model in which the metric is smooth everywhere [13]. We take the ansatz (2) with F (r) = B(r) = 1 − 2m(r)/r and where erf is the error function. This metric describes matter fields localized at r 0 , which we take to be r 0 = 2M with a spatial extent L. The metric above corresponds to a particular case of model II of matter-bumby BHs studied in Ref. [13].

Thin-shell, traversable wormholes
We consider the same model of a traversable wormhole [49,50] used in Ref. [12], which is obtained by identifying two Schwarzschild metrics with the same mass M at the throat r = r 0 > 2M . In Schwarzschild coordinates, the two metrics are identical and described by Eq. (2) with F = B = 1 − 2M/r. Because Schwarzschild's coordinates do not extend to r < 2M , we use the tortoise coordinate dr/dr * = ±F , where the upper and lower signs refer to the two different universes connected at the throat at r 0 = 0. The surgery at the throat requires a thin shell of matter with surface density and surface pressure [50] respectively. The weak energy condition is violated (Σ < 0), whereas the strong and null energy conditions are satisfied when the throat is within the PS, r 0 < 3M .

Thin-shell gravastars
We consider the thin-shell gravastar model [26] studied, e.g., in Ref. [51], which is described by the line element (2) with when Λ = 6M/r 3 0 , both F and B are continuous across the shell (more generic, thin-shell gravastar models have been recently studied in Ref. [46]). Note that the empty shell model discussed above is a particular case of this gravastar model when Λ = 0. Although F and B are continuous, their derivatives are not and this requires a thin shell with vanishing energy density and negative surface pressure [45]. 4. Boundary conditions at the shell Some of the models presented above require thin shells of matter across which the metric functions are discontinuous. Gravitational perturbations of discontinuous geometries can be studied by using the thin-shell formalism developed in Ref. [51]. In the main text we considered the scattering of scalar wavepackets, for which the junction conditions are easier to obtain.
In the frequency domain (FD), the Klein-Gordon equation on the background (2) can be written as with dr/dr * = √ F B and V l (r) = F l(l + 1) Here, primes stand for derivative with respect to r. Therefore, the potential is generically discontinuous if F , B, or their derivatives have a jump across the shell. Given that Eq. (A6) is homogeneous, we can assume that X lmω is continuous without loss of generality. In this case, even if V l has some finite jump, the integral of V l X lmω across the shell is continuous. Therefore, by integrating all terms in Eq. (A6) across the shell we obtain the junction condition i.e., the first derivative (with respect to the tortoise coordinate) of the scalar wavefunction is continuous. However, note that in Schwarzschild coordinates the first derivative might be discontinuous due to the definition of r(r * ). For example, for the wormhole considered above dr/dr * = ±F and therefore the first derivative (with respect to the Schwarzschild coordinate) changes sign across the shell. In the main text, we imposed the junction condition (A8) when solving Eq. (3). A similar procedure was done for gravitational perturbations, where now assumptions on the stiffness of the matter at the shell have to be dealt with. In practice, we implemented a condition akin to (A8) for the Zerilli-Moncrief wavefunction [12].

Appendix B: Regge-Wheeler-Zerilli formalism
Our numerical perturbation theory results shown in Fig. 3 are found by solving the first-order field equations in Regge-Wheeler gauge using a FD code. In what follows one would typically use Schwarzschild r to express spatial dependence. However, given that in this work we consider the wormhole model in addition to a Schwarzschild background, we make use the tortoise coordinate, as explained in Sec. A 2. In the both the wormhole and BH cases r * → ∞ as r → ∞ in the primary universe. On the other hand, in the BH case r * → −∞ as r → 2M , while for the wormhole model r * → −∞ corresponds r → ∞ in the other universe. In both cases, for a given radiative lm mode the field equations reduce to a single 1+1 wave equation, Here V l (r * ) is either the Zerilli potential (l + m even) or the Regge-Wheeler potential (l + m odd). The particle is assumed to be confined to the equator θ = π/2 and then the source contains terms proportional to the Dirac delta function and its first derivative where r * p = r * p (t) is the particle's radial location measured in r * . The time dependent functions G lm (t) and F lm (t) result from the tensor spherical harmonic decomposition of the stress-energy tensor of the point mass and subsequently evaluating r → r p (t) and ϕ → ϕ p (t). Their specifics will be discussed in further detail below.
We move to the FD with a Fourier transform and the inverse relations Given these, the master equation (B1) takes on the following FD form We assume retarded boundary conditions and asymptotically unit-amplitude homogeneous solutions, X ± lmω (r * → ±∞) = e ±iωr * .
The solution to Eq. (B5) follows from the method of variation of parameters, and W lmω is the (constant-in-r * ) Wronskian. Extending the integrals in Eq. (B8) over all space provides the normalization coefficients Finally, inserting the specific form of the source from Eqs. (B3) and (B2) this becomes an integral over time [31]. While a full solution to Eq. (B5) requires the functions c ± lmω (r * ), the constants C ± lmω are all that are required to compute the total radiated energy and angular momentum, as well as the waveform at infinity.
From a practical standpoint, repeated evaluations of the integral (B9) for a range of l, m, and ω make up the brunt of our calculation. The integral converges if the source coefficients G lm (t) and F lm (t) die off as t → ±∞, or equivalently, as r * → ±∞. These source coefficients are unique to the specific master function used, and while all master function sources decay rapidly at the horizon (for the BH case), at infinity they have differing behaviors. For bound motion it is best to use the Zerilli-Moncrief [52] (ZM) function (for l + m even) and the Cunningham-Price-Moncrief [53] (CPM) function (for l + m odd) because they allow for simple time domain reconstruction of the metric perturbation amplitudes. However, for unbound motion these variables are less than ideal. The ZM source tends to a constant at infinity while the CPM source falls off slowly as r −1 p . Therefore, in the even-parity sector it is better to use Zerilli's original variable [54], which has a source that decays as r −1 p . Meanwhile, in the odd-parity sector Regge and Wheeler's [55] original variable is preferable, since its source decays as r −3 p . These original variables are essentially the time derivatives of the ZM and CPM variables, which accounts for the more-rapid source fall-off. It is possible to take further time derivatives and define new master functions with source terms that decay even faster. In particular, the time derivative of the Zerilli function has a source which falls off as r −3 p at large distance, making it far more efficient than the Zerilli variable. In an upcoming work [56] the relation between all these master functions will be explored and the specific details of their sources will be given.
When using the ZM and CPM variables, the total energy radiated for a given lm mode to infinity is The gauge invariant waveform can also be computed from our code As r * → ∞ the FD particular solutions go to X + lmω (r * → ∞) = C + lmω e iωr * (B11) Thus, in order to evaluate these at retarded time u = t − r * and r * → ∞ we compute Ψ lm (u, r * → ∞) = 1 2π ∞ −∞ C + lmω e −iωu dω.
As shown in Ref. [57], these can be summed over l and m to form the transverse-traceless metric perturbation. In practice, the integrals (B10) and (B12) must be discretized. In order to obtain our results we sampled ∆ω as small as 2.5 × 10 −5 /M in a frequency range as wide as −1/M ≤ ω ≤ 1/M , skipping only the zero frequency mode which contributes no radiation and provides only a constant offset to the waveform.