{'fraction_non_alphanumeric': 0.0492590759556319, 'fraction_numerical': 0.013362834434772482, 'mean_word_length': 4.288683602771362, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 6, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Characterization of echoes: A Dyson-series representation of individual pulses

The ability to detect and scrutinize gravitational waves from the merger and coalescence of compact binaries opens up the possibility to perform tests of fundamental physics. One such test concerns the dark, nature of compact objects: are they really black holes? It was recently pointed out that the absence of horizons -- while keeping the external geometry very close to that of General Relativity -- would manifest itself in a series of echoes in gravitational wave signals. The observation of echoes by LIGO/Virgo or upcoming facilities would likely inform us on quantum gravity effects or unseen types of matter. Detection of such signals is in principle feasible with relatively simple tools, but would benefit enormously from accurate templates. Here we analytically individualize each echo waveform and show that it can be written as a Dyson series, for arbitrary effective potential and boundary conditions. We further apply the formalism to explicitly determine the echoes of a simple toy model: the Dirac delta potential. Our results allow to read off a few known features of echoes and may find application in the modelling for data analysis.

Black holes (BHs) are one of the most intriguing solutions of Einstein field equations. They play a key role in the understanding of the mathematical content of the field equations, but are also a fundamental player in astrophysics and -within the gauge-gravity duality context -in high energy physics. In General Relativity (GR), BHs are very simple and characterized by a one-way membrane called the event horizon. The event horizon and associated boundary conditions is the chief responsible for a number of properties of BHs, including their simplicity [1]. In addition, the "one-way" property of horizons causally disconnects the BH interior from its exterior, shielding outside observers from unknown quantum effects generated close to singularities or those presumably going on at the Cauchy horizon of spinning or charged BHs. Nevertheless, such extraordinary properties naturally raise the question of whether horizons do indeed always form, and where is this dynamical mechanism hiding in the field equations [2,3].
Other questions relate to the nature of compact objects and to tests of gravity in regions of where spacetime becomes so warped as to potentially form horizons and hide singularities, regions with huge intrinsic redshift. There are numerous tests of gravity in the "weak field" regime, and GR seems to successfully describe these regions up to the currently accessible precision. What about strongly warped spacetime regions? What evidence do we have that they exist, how is it quantified, and how can we quantify the evidence for horizons?
These questions all meet at horizons. Fortunately, the access to such fundamental questions has recently been granted through the historical detection of gravitational waves (GWs) by aLIGO [4]. Compact binaries are the preferred sources for GW detectors (and in fact were the source of the first detected events [4,5]). The GW signal from compact binaries is naturally divided in three stages, corresponding to the different cycles in the evolution driven by GW emission [6][7][8]: the inspiral stage, corresponding to large separations and well approximated by post-Newtonian theory; the merger phase when the two objects coalesce and which can only be described accurately through numerical simulations; and finally, the ringdown phase when the merger end-product relaxes to a stationary, equilibrium solution of the field equations [8][9][10].
A. Quantifying the evidence for horizons with precision GW-physics: Multipoles, tidal heating and tidal Love numbers All three stages provide independent, unique tests of gravity and of compact GW sources. Due to their very nature, it is impossible to prove the existence of event horizons [2,11]. However, the evidence for horizons and BHs may be significantly tightened using GW observations: During the inspiral, three different features will play a role, all of which can be used to test the BH-nature of the objects. Firstly, horizonless objects are bound to possess a different multipolar structure from Kerr BHs, impacting the GW phase, and allowing for constraints on possible deviations from the Kerr geometry [12]. Secondly, the presence of horizons translates into BHs absorbing radiation in a very characteristic way. Horizons absorb incoming high frequency radiation, and serve as sinks or amplifiers for low-frequency radiation able to tunnel in. However, horizonless objects are expected not to have significant tidal heating. Thus, a "null-hypothesis" test consists on using the phase of GWs to measure absorption or amplification at the surface of the objects [13]. LISA-type GW detectors [14] will place stringent tests on this property, potentially reaching Planck scales near the horizon and beyond [13]. Finally, non-or slow-spinning BHs were shown to have zero tidal Love numbers. In a binary, the gravitational pull of one object deforms its companion, inducing a quadrupole moment proportional to the tidal field. The tidal deformability is encoded in the Love numbers, and the consequent modification of the dynamics can be directly translated into GW phase evolution at higher-order in the post-Newtonian expansion [15]. It turns out that the tidal Love numbers of a BH are zero [16][17][18][19][20][21], allowing again for null tests. By devising these tests, existing and upcoming detectors can rule out or strongly constraint boson stars [13,[22][23][24] or even generic ultracompact horizonless objects [13,23]. The GW signal overall depends on the entire spacetime structure. The late-time behavior depends, mostly, on the strong-field region and on the boundary conditions set there. Generically, a BH-type ringdown is excited when the spacetime is disturbed at the photosphere. We call this the photosphere ringdown, which has a timescale of order ∼ 5GM/c 3 [2] (we will use geometric units, G = c = 1, from now on). If light or GWs hit the surface of the putative horizonless objects on timescales of this order or smaller, then the ringdown of the spacetime is significantly different from that of a BH. Thus, tests of GR or of the nature of the object can be performed using the standard ringdown tools [1,25,26].
These techniques and tools require very sensitive detectors and hinge on the ability to do precision physics. It turns out that, should horizons not exist, there should be clear signatures of this absence in GW signals. If the outside spacetime is vacuum GR all the way down to a region very close to the horizon, in a way that GWs take a time 5M to hit the surface, then we call these Clean Photosphere Objects (ClePhOs). Because ClePhOs have photospheres which are vacuum in their vicinity, the photosphere ringdown stage is therefore indistinguishable from that of a BH [2,27]. However, on longer timescales the GWs excited at the photosphere have had time to hit the surface and bounce back to the observer, with a fraction of it trapped between the object and the photosphere: the late-time response of a ClePhO is a series of GW echoes [2,[27][28][29][30][31][32]. From this discussion it is clear that echoes should also appear whenever there is structure in the near horizon region, either due to the compact object or due to propagation effects [33].

C. The morphology of echoes and purpose of this work
The non-detection of echoes can place interesting bounds on the surface of ClePhOs and rule out certain models; most importantly, it would quantify in an independent way the evidence for BHs in our universe. Given the tremendous significance of a possible detection of echoes in a GW signal, it is worth devising strategies to detect them in actual stretches of data. Early studies, based on simple templates, have claimed an important statistical evidence for the presence of echoes in the first detections [34,35]. At this point the evidence is not of enough statistical significance to warrant claims of a detection [36,37]. To improve the significance and to understand the underlying physics, a better description of GW echoes is necessary.
An important step in this direction was given recently by Mark and collaborators [38] (see also Refs. [39,40]). These authors have shown how the response of a ClePho can be expressed in terms of the response of a BH, convoluted with an appropriate function which takes into account the different boundary conditions. Doing this, they were able to identify the contribution from isolated echoes and proposed templates for non-spinning BHs which had a very good match with numerical waveforms. Subsequent work has investigated phenomenolog-ical waveforms for detection of echoes [41]. Clearly, the result of these incursions is that a lot more remains to understand. Some of the open issues include, • The very late time behavior of the GW response of ClePhOs is described by the quasinormal modes (QNMs) of the system. The QNMs arise as poles of the relevant Green's function. It is a puzzle that the photosphere modes have no special spectral role, and yet they dominate the response at early times.
• In addition, although the late-time response of Cle-PhOs is presumably governed by its QNMs, no numerical evidence exists of this fact.
• The overall amplitude of sucessive echoes decreases, at least if one is looking for consecutive echoes generated shortly after merger. What type of decay is this, is it polynomial, exponential? Can we characterize the evolution of echoes in a more precise manner?
• The delay between different echoes is a key quantity in any detection strategy. Is the delay really constant or does it evolve in time, and how [42]?
• In a related vein, a generic widening of the pulses, in the time-domain, was observed as time goes by. This is physically intuitive: the pulses are semitrapped within a cavity that let's high frequency waves pass. At late times only low-frequency, resonant modes remain. Hence the pulse is becoming more monochromatic. This is a generic but not yet quantified result.
• Consecutive echoes may be in phase or out of phase, depending on the particular boundary conditions imposed on us by the physical model.
The purpose of this work is to start the discussion on some of these issues. Here we provide an exact expression for each individual echo amplitude, in the form of a Dyson series, that accounts for most of the above points. This formula also works retrospectively: if one echo is measured, then the surface reflectivity becomes fully specified (assuming the potential is known).
The first section consists of the perturbative approach to wave scattering in which: i) the most general appropriate boundary conditions are considered; ii) the Lippmann-Schwinger equation is obtained; iii) the resulting Dyson is resummed into a sum of echoes; iv) the time-dependent solution is obtained through Laplace inversion.
In the second section, we apply this apparatus to a Dirac delta type of potential. The Dyson series does not need to be truncated since the explicit solution is completely attainable. Here we compute the echo frequency amplitude for arbitrary surface reflectivity, and the echo waveform for a set of common boundary conditions.

II. WAVES IN OPEN SYSTEMS
A linear perturbation Ψ(t, x) in a system with potential V (x), supported on R, obeys the wave equation (1)

A. Boundary Conditions
We're interested in open systems, in which waves are only partially bounded and can escape to infinity, in at least of one of the sides. We'll choose +∞ to be the open end. Then, The other boundary condition (BC) can correspond to a reflection at some point x = −L. When the potential vanishes it can be written as where the first term on the rhs corresponds to a wave travelling to the left, out of the system, and the second term is nothing but the reflected wave, travelling to the right, thus, we can identifty R(ω) as the reflectivity associated with the BC at x = −L.
If we do not wish for external influence on the system, the reflectivity should obey otherwise, the reflected wave has a larger amplitude than the outgoing wave, i.e. there's an external input at the left. The condition above is however violated for some well-known systems, such as Kerr BHs, which are known to display superradiance [43]. The reflection coefficient R(ω) is completely specified by the BC at x = −L. We can point out three familiar cases. For a purely outgoing wave, to the left, we simply have Both of the latter two are conservative, they satisfy |R(ω)| = 1.
Conversely, R(ω) can be arbitrarily chosen and the BC at x = −L becomes automatically imposed. For instance, dissipation can be introduced by generalizing the latter reflectivities to with r ∈ [−1, 1] and |R(ω)| = |r| ≤ 1. The BC at
The time-dependent solution is the inverse of this transform, where β assumes any value β > α to ensure the integrand is always convergent along the path of integration. With these definitions, Eq. (1) is reduced to the ODE with source term where ψ 0 (x) = Ψ(0, x) andψ 0 (x) = ∂ t Ψ(0, x) are the initial conditions. Now, instead of pursuing the usual Green's function approach, we use a perturbative framework. The ODE (8), and BCs, can be rewritten in the Lippmann-Schwinger integral form is the Green's function of the free wave operator d 2 /dx 2 + ω 2 with BCs (2) and (3), and is the free-wave amplitude. The formal solution of Eq. (10) is the Dyson series which effectively works as an expansion in powers of V /ω 2 , so we expect it to converge rapidly for high frequencies.
Note that if we were to expand each term of the series with explicit use of (11) we'd get a panoply of powers of R(ω). Now, we may ask, is it possible to reorganize (13) and express it as a series in powers of R(ω)?
C. Resummation of the Dyson series and echoing structure We start by separating the Green's function (11) the open system Green's function, and the "reflection" Green's function. We can then write (10) as Now, in the same way as a Dyson series is obtained, we replace theΨ(ω, x ) in the third integral with the entirety of the rhs of Eq. (16) evaluated at x . Collecting powers of R(ω) yields where, for better clarity, we chose not to write the functions' arguments. If we repeat the process one more time, i.e. by replacing Eq. (16) withΨ in the last integration in (17), we get and a pattern starts to emerge. The first line does not contain any g r , the factor of R contains one g r arranged in all possible distinct ways with the g o 's, the factor of R 2 contains two g r 's also arranged in all possible ways, and so on and so forth. If we continue this process we end up with a geometric-like series in powers of R, with each term a Dyson series itself: the series stemming from the first line of (18), and the reflectivity terms, which can be re-arranged as, where x 0 := x, S k is the permutation group of degree k and 1 n!(k−n)! σ∈S k represents the sum on all possible distinct ways of ordering n g r 's and k−n g o 's, resulting in a total of |S k | n!(k−n)! = k n terms. For instance, for k = 3, n = 2, we have (22) which can be more eficiently obtained by mantaining the functions' arguments in the same position, but instead interchanging the relative positions of the g r 's and the g o 's.
Without a doubt we've increased the mathematical complexity of the problem. Nonetheless, Eq. (21) has special significance: it's the frequency amplitude of the n-th echo of the initial burst. There's no proper way to show this since there's no mathematical definition of an echo. However, with the following discussion and further application of this formalism to the Dirac delta potential, we hope to provide enough justification.
If R = 0 thenΨ =Ψ o , the open system waveform, where only g o participates. Conversely, when we don't have a perfectly transmissible boundary (R = 0), we get an additional infinite number of Dyson series, as stated in Eq. (19). TheseΨ n terms are expected to give a smaller contribution toΨ as n increases. This is mainly due to two features in (21).
• First, when |R(ω)| < 1, R n (ω) is obviously an attenuation factor with a greater impact at large n. It indicates n partial reflections at the boundary, as physically done by the n-th echo. Moreover, echoes have the distinctive feature of being spaced by the same distance for any pair of successive echoes. The fact thatΨ (n+1) has an additional factor of R(ω) thanΨ (n) , hence an, independent of n, phase difference of arg[R(ω)], indicates this.
• Moreover, the fact that the Dyson series starts at k = n. Since g o and g r are of the same order of magnitude, it is natural to expect that the series starting ahead (with less terms) has a smaller magnitude and contributes less toΨ than the ones preceding them. The additional term thatΨ n possesses when compared toΨ n+1 , and can thus be used to evaluate their amplitude difference, is given by Furthermore, latter echoes are seen to vibrate less than the first echoes. As stated before, since the Dyson series is basically an expansion on powers of V /ω 2 , by starting at k = n,Ψ n skips the high frequency contribution to the series until that point. This is intuitively due to high frequency signals trespassing the potential barrier more easily than lower frequency signals, which is the reason why high frequency behaviour predominates in the earlier echoes.

D. Inversion into the time domain
Finally, let us make use of the inverse Laplace transform (7) to obtain the time-dependent solution of wave equation (1). We start with the open system perturbation (20). The frequency dependent terms are the Green's functions g o , which have a pole at ω = 0 (Eq. (14)), and the source term I which does not have any pole (Eq. (9)). Thus, to keep the integrand convergent we should integrate above ω = 0. The frequency integral of the k-th term of (20) is where we've chosen β = 1 (> 0), in Eq. (7). The integrand is singular except when k = 1, due to the term iωψ 0 (x) in I(ω, x), that cancels the ω in the denominator. For this term, we have whereas for the term −ψ 0 (x) of I(ω, x), we have to integrate a simple pole at ω = 0, which vanishes for t < |x − x 1 |: The initial signalψ 0 (x 1 ) did not have enough time to travel to the point of observation x, i.e. these points are causally disconnected. Finally, integration of ψ 0 (x 1 ) andψ 0 (x 1 ) with (25) and (26), respectively, yields the first term of Ψ o (t, x), If both R(ω) and V (x) vanish, this is nothing but the final solution Ψ(t, x). The equation above reveals that the initial waveform separates in two halves, propagating in opposite directions, just like an infinite plucked string.
For k = 1, we start by defining the causal distance, involving k interaction points besides the point of observation x and the source point x k+1 , for an elapsed time t.
With this definition the argument of the exponential in (24) is simply iωs k−1 . This integration, for k = 1, is If I(ω, x k ) was independent of ω, the term in brackets would be (is k−1 ) k−1 I(x k ). But since I has the linear form (9), we can write the term inside brackets as (is k−1 ) k−1 I(−i k−1 s k−1 , x k ) . Putting everything together yields a Taylor-like expansion, Inversion of Eq. (21) follows the same lines. Instead of s k , it is useful to define s n,k := (x − x 1 ) + · · · + (x n−1 + x n ) + |x n − x n+1 | + · · · + |x k−1 − x k | − t , (31) in order to write the frequency integral, corresponding to inversion of the k-th term of (21) through (7), as where we replaced the Green's functions by their explicit forms (14) and (15). Now, we can't go further unless we know R(ω) in detail. More specifically, its poles and divergent behaviour at ±i∞, which specify the choice of contour.
For completeness, we present below the calculation for R given by Eq. (5): with σ(s n,k ) := (x σ(1)−1 − x σ(1) ) + · · · + (x σ(n)−1 + x σ(n) ) and the reflected initial waveform, present only in the first echo (due to the Kronecker delta δ n,1 ). The more critical reader may realize that this method is only useful if the explicit form of R(ω) is known. For instance, this is not the case for a wormhole system, where R stands for the reflectivity of the Schwarzschild potential which can only be extracted numerically. Thus, one may ask if it is also possible to express the reflectivity of a generic potential as a perturbative series. The answer is yes.

E. Reflectivity series
If we send a wave from +∞ (e −iωx ), the reflectivity will be the factor of the reflected wave, R(ω)e iωx . The source term that corresponds toΨ 0 = e −iωx can be inspected from (12), with g = g o (we're trying to extract the reflectivity, so it's only natural to consider purely outgoing BCs at both sides), and is formally given by (36) which non surprisingly corresponds to a source pulse located at x → ∞. The factor of e −iωl takes care of the phase difference. Now, with this source term, the solution, given by the Dyson series (20), at x → ∞ has the formΨ(ω, x → ∞) = e −iωx + R(ω) e iωx , with the reflection coefficient expressed in terms of the potential, as promised.
The above expression can also be used to compute the system's QNMs, which are the poles of R(ω). In fact, there's an ongoing discussion on whether the QNMs of the system with purely outgoing BCs at both sides (as in the above case) coincide with the ones where a mirror is introduced, thus replacing the outgoing BC at one side with other, more complex, BC. The mirror + potential system's QNMs should also be the poles of its "reflectivity". This concept, however, is not defined in the case both BCs are not purely outgoing, that is, if the system is only partially open. We can't simply take the initial wave asΨ 0 = e −iωx , but instead, where R(ω) is NOT the reflectivity of the system but the reflectivity associated with the non trivial BC at some x = −L. This is the correct form forΨ 0 since, by Eq.
(10), it should be the complete solution of the system when there's no potential barrier and also reduce to e −iωx when the mirror vanishes, that is, a free wave travelling to the left. The reader may find comfort in this definition by noting thatΨ 0 computed through Eq. (12) with I given by Eq. (36) does indeed recover expression (38). Naturally, the system's reflectivity is still defined by the factor multiplying the outgoing wave at +∞, R(ω) e iωx . With the source term (36) and g given by Eq. (11), we just need to evaluate expression (13) at x → ∞ to extract which reduces to R = R if the potential vanishes, and to Eq. (37) if R → 0, as expected. We reinforce that the R(ω) in the above expression corresponds to the reflectivity associated with the non-trivial BC at x = −L whereas, in (37), it is the reflectivity of the potential barrier. We use the same letter for both since the mirror at x = −L can be either due to a BC at this point, or a potential barrier, in this case computable through Eq. (37).
One can see that Eq. (39) does not diverge where Eq. (37) diverges, for arbitrary potential. In other words, the mirror+potential system and the completely open potential system do not share the same spectrum of QNMs.
In the next section, we'll apply this apparatus to a specific potential to explicitly see this difference.

III. THE DIRAC DELTA POTENTIAL
We now apply the previous formalism to the Dirac delta potential, with V 0 > 0.

A. Open system solution: Ψo
Instead of employing Eq. (30) straight ahead, it's interesting to first compute the frequency amplitude from Eq. (20). The k = 1 term corresponds to the freely propagating initial waveform,Ψ i (ω, x). For k > 1, the k − 1 delta functions collapse all the integrals except the integration in x k , which is in fact a k-independent integral: Relabeling x k → x and treating the sum as a geometric series, simplifies the above tõ with the reflectivity of the Dirac delta potential (40), which could be directly computed from Eq (37). It diverges at the QNM WithΨ o in hand we just have to apply Eq. (7) to get the time-dependent solution: and Ψ i (t, x) given by (27). Direct application of Eq. (30) would even be more straightforward: Instead of a geometric series, the infinite sum that factors out is the Taylor series of e V0(|x|+|x |−t) .
Before we move on, we should point out the following. When there are no interactions, V 0 = 0, the two latter terms of (45) cancel each other and, as expected, Ψ o (t, x) = Ψ i (t, x). Also note that, unlike conservative systems, the QNM excitation coefficient C V0 is not a constant. One may then ask in what conditions does Ψ o (t, x) decay with the QNM behaviour for t → ∞. We expect this to happen when I(ω, x) is sufficiently localized in space, corresponding to more "physical" sources. Even a decay I(ω, x) ∼ e −a|x| , for some a > 0 gives Ψ o (t → ∞, x) ∼ e −at . For a gaussian source 2a . Even if the gaussian is disperse (small a), which makes b assume large values, for t b the integrand will contribute little and C V0 (t) is essentially independent of t. In the very limit t → ∞, C V0 will just be the integral of a gaussian in the real line, with convergent (and known) value and hence Ψ o (t → ∞, x) = C V0 e −V0t .

B. Echoes: Ψn
To obtain Ψ(t, x) we still need to get Ψ n (t, x), according to Eq. (19). Since we haven't yet specified R(ω), we must start at the frequency amplitude and employ Eq. (21) with potential (40).
As in the previous case, the delta functions will collapse all integrals in the k-th term of expansion (21), except the one in x k , implying that the Green's functions product in the integrand will have the form g r (x, 0)g r (0, 0) · · · g r (0, 0)g o (0, 0) · · · g o (0, 0)g o (0, x k ), (47) for the identity permutation σ(j) = j.
(48) In fact, there are only 4 possible algebraically different outcomes for the pair of functions at both ends of the product, which make the sum 1 n!(k−n)! σ∈S k over the g's simplify to To get the first term, for example, we have two g r 's at the ends, leaving k −2 spots for the remaining n−2 g r 's to be organized. The remaining terms follow the same reasoning. Albeit not directly apparent, we are dealing with a total of k n terms, as pointed out after Eq. (21), since Now, similarly to what happened in Eq. (41), renaming x k → x will make the integral (21) independent of k and a geometric-like series factors out for every term in Eq. (49). For g r (x, 0)g r (0, x ), for instance, what factors out is with R δ given by Eq. (43), where we've used the identity for the power of a geometric series, Using the same identity for the remaining terms yields This variety of terms can be physically interpreted. The first one, with the product R n−1 δ R n , corresponds to a wave sent left, towards the mirror, and also received from the mirror, travelling to the right. Take the second echo, n = 2, for instance. It reflects first at the mirror, then at the delta, and again at the mirror, picking up a factor R δ R 2 .
The second one, with R n+1 δ R n , is the reverse situation. The wave is sent right, towards the delta, and then also received from the delta, but travelling to the left. This situation can only happen for x ∈ [−L, 0], when the observer is inside the cavity.
The same is true for the remaning terms, with R n δ R n . Here, one of two situations happen. The wave is sent into the mirror and then received from the delta, or first sent into the delta and then received from the mirror, reflecting, in either case, an equal number of times at the delta and at the mirror.
The echoes' amplitudeΨ n is similar, in form, toΨ 0 .
Besides the presence of R n , the difference lies in the order of the pole of the QNM (44), due to the powers of R δ . Inversion will result in derivatives of the integrand, evaluated at the QNM. Thus, the echoes besides vibrating and decaying with the delta QNM, have a slightly different behaviour which we'll see below to be of the polynomial sort.
To proceed with inversion, with the use of Eq. (7), we consider R(ω) given by Eq. (5), to get A few comments must be made. Interaction of the source with the delta potential is being accounted for in the first integral of Eq. (55) whereas reflection at the mirror is being accounted in the second integral, hence the Θ functions ensuring that there is enough time for the source to reach the delta and the mirror, respectively.
The factor (1 − δ n,0 ) vanishes for n = 0 and is 1 otherwise. It's easy to see that it only vanishes for Ψ 1 , the first echo, which instead possesses the term on the first line of Eq. (54), corresponding to the reflection at x = −L of the left-travelling intial waveform. This is the only surviving term in case V 0 → 0, when there is no cavity.
It is interesting to note that the integrals themselves do not depend on n, apart from the integration limits. The difference between echoes mostly lies in the order of the derivative on V 0 . For instance, whenψ 0 = 0, the derivative will only act on the QNM exponential factor assigning the said polynomial behaviour to the echoes' waveform: (56) Figure 3 shows a "time-lapse" of the complete waveform given by the sum of the open-system solution, Eq.
(45), with the first 3 echoes, described by Eq. (54), with initial condition and system parameters The complete sequence of events can be seen in video format at: https://youtu.be/XfJNwuwbvnA .

C. QNMs
To substantiate the discussion at the end of the first section, let us compute the full system's reflectivity. Usage of Eq. (39) with potential (40) yields which simplifies to by using the geometric series identity and definition of R δ , Eq. (43). It is easy to check that R δ → R δ if R → 0 and viceversa, if R δ → 0 then R δ → R. Moreover, if R = −1, corresponding to a perfect mirror, then we should expect everything to be reflected back, independently of the potential. In this limit we can also see that R δ = −1.
More interestingly, at the delta QNM (44), R δ → ∞, the dependence on R δ cancels to give R δ = − 1+2R R , which is finite for non-vanishing R. Thus, ω = −iV 0 is not a QNM of the mirror+delta system. The QNMs are instead implicitly given by R(ω n )R δ (ω n ) = 1 .
(61) Fig. 1 plots the frequencies that respect the above, for R given by Eq. (5), with r = 1 (Dirichlet BC at x = −L), which are well-approximated by the expression Unsurprisingly, the imaginary part grows in magnitude with |n|. This implies that, at sufficiently long times, the perturbation will decay with the fundamental mode ω ±1 (Fig. 2), even if the initial perturbation decays according to the pure-delta QNM (44) (as we show in Figure 3). We believe that this is the most convincing demonstration to date that the late-time decay is indeed governed by the QNMs of the composite system.

IV. FINAL REMARKS
We have shown that a proper re-summation of the Dyson series solution of the Lippman-Schwinger equation accounts for the presence of echoes in the waveforms of extremely compact objects (termed "ClePhOs" in the nomenclature of Refs. [2,28]). We recover previous results, obtained with a completely different approach [38], but our approach also provides a few more insights.
The key result Eq.(21) besides confirming the lower frequency content, decaying amplitude and constant distance of successive echoes, also explicitly relates the echoes waveform Ψ n with the initial conditions and sources incorporated into I(ω, x), the potential of the system V (x), and the reflectivity of the wall R(ω) (the Time evolution of the waveform using initial conditions (57) and parameters (58). The plot shows the decay with the fundamental mode of the system (with mirror on the left), ω±1 ≈ −0.00205 (in red), for large t. The early echoes decay in a way that is governed by the QNMs of the pure delta. The high-frequency component is filtered out and progressively the signal is described by the modes of the composite system at late times, as it should.
latter two function as the right and left sides of the lossy cavity, respectively).
With the hypothetical future discovery of echoes in gravitational wave signals, the echo amplitudeΨ n can be extracted up to experimental and numerical error. Together with the knowledge of V (x) and I(ω, x), this turns Eq. (21) into an equation for R(ω), which encodes the information we currently lack on the quantum structure at the event horizon.
Further developments of our formalism should include: application to other potentials besides the simple Dirac delta; a careful analysis of the convergence properties of Eq. (21); extension of our methods to more than one spatial dimension; check if superradiant amplification is observed in Eq. (21) if |R(ω)| > 1; implementation of the reflectivity series (37) and (39) to QNM computation; confirm the polynomial behaviour of echoes in other systems besides the Dirac delta potential and use this information to echo modelling. The initial profile (57) quickly gives way to two pulses traveling in opposite directions at t = 6, as described by Eq. (27); the left propagating pulse interacts with the (delta) potential at t = 10.5 and gives rise to a transmitted pulse and a reflected one (t = 16). The reflected pulse eventually reaches the boundary, at t = 20, and will cross the potential at around t = 36 giving rise to the first echo. The wave confined to the cavity (mirror+potential) will produce all subsequent echoes. At t = 58, after 2L = 20 time units the second echo emerges out of the cavity and at t = 74 a third echo is about to be produced. These snapshots were obtained by adding three "echoes," and coincides up to numerical error, with the waveform obtained via numerical evolution of the initial data. In the central panel, the red line shows −e V 0 (x−10) , confirming that the initial decay is described by the QNMs (44) of the pure delta function (no mirror). try Canada and by the Province of Ontario through the Ministry of Economic Development & Innovation. This article is based upon work from COST Action CA16104 "GWverse", supported by COST (European Cooperation in Science and Technology). This work was partially supported by FCT-Portugal through the project IF/00293/2013, by the H2020-MSCA-RISE-2015 Grant No. StronGrHEP-690904.