Theory for Anomalous Terahertz Emission in Striped Cuprate Superconductors

Recent experiments in the doped cuprates La$_{2-x}$Ba$_x$CuO$_4$ have revealed the emission of anomalous terahertz radiation after impulsive optical excitation. Here, we theoretically investigate the nonlinear electrodynamics of such striped superconductors and explore the origin of the observed radiation. We argue that photoexcitation is converted into a photocurrent by a second-order optical nonlinearity, which is activated by the breaking of inversion symmetry in certain stripe configurations. We point out the importance of including Umklapp photocurrents modulated at the stripe periodicity itself, which impulsively drive surface Josephson plasmons and lead to a resonant structure of outgoing radiation, consistent with the experiments. We speculate on the utility of the proposed mechanism in the context of generating tunable terahertz radiation.

Electronic phases of quantum matter are typically distinguished by their signature electromagnetic responses [1][2][3][4][5][6][7][8].For instance, in strongly anisotropic materials, the onset of superconductivity manifests as the appearance of Josephson plasmon (JP) edges in the terahertz reflectivity [9][10][11][12].Likewise, states with an incommensurate charge density wave (CDW) exhibit a characteristic quasiparticle gap in the low-frequency conductivity, accompanied by resonance at the pinning frequency [13][14][15].Recent experiment [16] in the doped single-layer cuprates La 2−x Ba x CuO 4 demonstrated a unique nonlinear optical signature of the "superstripe" phase, i.e., when both superconducting and stripe orders are present and intertwined [17][18][19][20][21][22].The discovered phenomenon can be summarized as follows: upon a strong optical pump pulse, outgoing radiation is observed at terahertz frequencies, with a spectrum peaked at the Josephson plasmon resonance (JPR).These findings are surprising because these materials are nominally centrosymmetric, and such radiation is not expected in the absence of a current or magnetic field bias.In this paper, we provide a theoretical interpretation of these experimental observations.Crucially, this nonlinear effect cannot be understood from the perspective of either the superconducting or stripe orders individually and requires combining both types of symmetry breaking simultaneously.
The phenomenon of "radiating stripes" in LBCO provides several important insights into the nature of this striped superconductor.The presence of a subharmonic optical response constrains the symmetries of the stripe order [23] and, in particular, implies inversion symmetry is broken [24][25][26][27][28]. Consequently, the impulsive optical excitation is now coupled to both bulk and surface Josephson plasmons, which are otherwise symmetry-odd infrared-active modes [29].We argue that the outgoing terahertz radiation originates from the surface Josephson plasmons and is sensitive to the stripe order, which outcouples these otherwise silent modes.From the observation that the terahertz radiation is sharply peaked at the Josephson plasmon resonance, we can infer that a photocurrent is being generated both at zero momentum and at momenta corresponding to the reciprocal vectors of the CDW lattice, yielding insight into the microscopic physics of the striped state.
The key finding of this paper can be summarized as follows.The role of the CDW order is that it i) makes the surface modes optically active [Fig.1] and ii) gives rise to a nonzero Umklapp photocurrent [24].We show that Umklapp photocurrents resonantly drive the surface plasmons, which exhibit a large density of states near the Josephson plasma resonance ω JPR ≃ 0.5 THz, thereby resulting in sharp in frequency radiation [Fig.2] consistent with the experiment of Ref. [16].We point out that the proposed mechanism is generic and potentially useful for designing platforms for the generation of tunable terahertz radiation [30][31][32][33][34].We remark that Umklapp currents carry a large momentum of the CDW reciprocal lattice so that nominally, the resulting emitted light is expected to be far away from the light cone and, thus, decay on a length scale of the order of a few wavelengths at most.Therefore, the proposed mechanism of resonant coupling to the surface collective modes provides a pathway towards detecting Umklapp currents with far-field optics.
Photocurrent generation-A crucial aspect of the experimental findings reported in Ref. [16] is that the frequency of the pump pulse Ω pm ≃ 375 THz is at least two orders of magnitude larger than ω JPR .This separation of energy scales suggests that the pump pulse drives the mobile electronic degrees of freedom that downconvert the large incoming frequency into a low-frequency photocurrent [8,[35][36][37][38][39][40].This is further bolstered by the expectation that inversion symmetry is broken for a realistic pattern of charge order [24].The photocurrent is characterized by a rank-three conductivity tensor: where indices a, b represent the Cartesian directions.Following Ref. [16], we assume that the pump-pulse width ∆Ω pm ≃ 5 THz is also much larger than ω JPR , in turn implying that the photocurrent appears featureless for ω ≲ ∆Ω pm , Ω pm : This should be contrasted to the experiment [16], where the observed radiation is sharply peaked in frequency at ω JPR .Intuitively, the broadband photocurrent acts like an impulsive drive to coherent Josephson plasma modes, both bulk and surface ones.Here we focus on the surface excitations because the scenario based on the bulk ones is inconsistent with the experimental data, as we discuss in the Supplemental Materials [41] and further elaborate on below.In the following discussion, we investigate a semiphenomenological model, where we assume the presence of such a broadband drive and study how it interacts with Josephson plasmons.
To further appreciate the importance of the surface excitations, we argue that the photocurrent generation occurs at the surface of the sample.Indeed, the frequency of the incoming photoexcitation is large so that the skin depth of how light penetrates the sample is small.From the equilibrium optical properties [16], we estimate it to be d ≈ 200 nm; this short electronic length scale is much shorter than any length scale associated with the collective excitations of the superconductor [42,43].To illustrate this explicitly, we consider the bulk polariton dispersion in the vicinity of q = 0: where ϵ ∞,ab is the high-frequency in-plane dielectric constant of the medium.Equation (2) allows defining the plasmon coherence length as l JP = c/(ω JPR √ ϵ ∞,ab ) ∼ 100 µm ≫ d.Hence, from the perspective of collective modes, one indeed can consider the photocurrent as a surface phenomenon.We remark that in the experiment, the beam-spot size is large d beam ∼ 500 µm ≳ l JP so that the pump can be approximately described as uniform along the interface.
For future reference, the electromagnetic response of layered materials, such as cuprates, is encoded into the anisotropic dielectric tensor ε(ω) = diag(ε ab , ε ab , ε c ). Motivated by the two-fluid model of anisotropic superconductors [41,[44][45][46][47], which accurately captures optical linear response of cuprate superconductors [45,46], we choose the following form: This form so far does not include the charge order -we will return to this below.The second term in Eq. (3) describes the reactive response of the superconducting fluid so that ε α (ω) ∼ ω −2 for ω → 0; ω ab and ω c = ω JPR are the in-plane and c-axis plasma frequencies.The strong anisotropy of cuprates implies ω ab ≫ ω c .We also expect the in-plane plasmons to be strongly damped γ ab ≫ γ c .The third term represents the normal fluid; γ α is a phenomenological damping parameter, which is proportional to the corresponding normal-fluid conductivity [41,45].
Unless specified otherwise, we use the following parameters [16,46,48]: ω c = 1 THz, γ c = 0.1 THz, ϵ ∞,c = 25, and ϵ ∞,ab = 4.Our conclusions do not rely on the choice of ω ab and γ ab as long as ω ab ≫ ω c .Surface Josephson plasmons-We begin by reviewing the properties of surface Josephson plasmons in the absence of CDW order [46,49,50].These are evanescent collective modes that are confined to the interface x = 0 [inset of Fig. 2], i.e., they decay both into the air and into the sample [46,51,52].Of the most physical interest for us below are the modes with the magnetic field pointing along the y-axis: Note that because of the translational invariance along the z-axis, we choose the same dependence on z inside and outside the sample.Here k a and k m are yet unknown wave vectors that depend on both ω and q z .We find them by solving the Maxwell equations in each media: 2. Emission spectrum.We find that the outgoing radiation can be sharply peaked around ωJPR only in the presence of a nonzero Umklapp photocurrent JU ̸ = 0, Eq. ( 13).This Umklapp component impulsively drives high momenta surface Josephson plasmons that can emit light, Fig. 1(b).Here we fixed A = 0.1ε∞,c.
To describe evanescent electromagnetic waves, we choose the roots with Re k a , Re k m > 0. By matching the Fresnel boundary conditions (E air,z = E mat,z and B air,y = B mat,y at x = 0), we obtain an implicit equation on the dispersion ω s (q z ) of the surface Josephson plasmons: The spectrum of these excitations is shown in Fig. 1(a).We find that ω s (q z ) quickly saturates at around ω JPR .
In other words, we expect a large density of states of these excitations near ω JPR .It is worth pointing out that the fact that the saturation occurs near the bulk plasmon resonance, ω s (q z ) → ω JPR for cq z ≳ ω JPR , is a consequence of the strong anisotropy ω ab ≫ ω c .For instance, in isotropic superconductors with ω ab = ω c and ε ∞,ab = ε ∞,c = 1, a similar saturation occurs but at a notably lower frequency ω JPR / √ 2. We finally remark that the two other surface Josephson plasmons exhibit a similar saturation but at much higher frequency.
An immediate problem we encounter is that the surface excitations lie outside the light cone and, therefore, cannot radiate.This issue is evaded by the charge order, which couples the mode with wave vector q z to the ones with wave vectors q z + nQ CDW , where n = 0, ±1, ±2, . . .and Q = Q CDW is the ordering wave vector along the z-axis.As such, the entire surface plasmon dispersion becomes backfolded to the reduced Brillouin zone defined by Q.Some modes back fold into the light cone [Fig.1] and, therefore, they can emit photons.This insight is crucial in explaining the experiment.In essence, the stripes act as a natural diffraction grating, which then out-couples the otherwise silent surface modes, much like a nano-fabricated corrugation would be engineered [30-34, 44, 53-56].
Here we describe the CDW order phenomenologically.Since the most dominant effect is expected to be due to the charge modulation along the z-axis, we neglect the fact that stripes have a non-trivial structure along the xand y-axes.Specifically, we assume that the CDW order parameter enters through the modulation of the c-axis plasmon frequency: This modulation, which we assume to be weak, modifies only the c-axis dielectric function, which we write as: where A captures the strength of the stripes order.For simplicity, we assume that the stripe period along the zaxis equals two lattice constants.This assumption is by no means crucial but allows us to make substantial analytical progress.Finally, to properly describe the momentum mixing, one should take into account the momentum dependence of the dielectric tensor (3).However, this dependence is expected to be nonessential at the scale of Q, as we show in [41], where we carefully consider the surface Josephson plasmons at large momenta.
To get the properties of surface Josephson plasmons in the medium with stripes, we proceed similarly as above but now take into account the Bragg mixing of the zaxis momenta.Specifically, we substitute the following evanescent-wave ansatz: where k a (q, ω) = q 2 − ω 2 /c 2 and ka = k a (q + Q, ω).
The wave vectors λ 1 and λ 2 , together with parameters γ 1 and γ2 that encode the mentioned mixing, are known functions of q z and ω, which are obtained by solving the Maxwell equations inside the sample with stripes [41].
The four remaining unknown coefficients α a , β a , α m , β m are related to each other via: B air,y = B mat,y and E air,z = E mat,z .These conditions, in turn, implicitly define the spectrum of surface plasmons through det M(q z , ω) = 0, where Schematic illustration of the modified by the stripes dispersion is shown in Fig. 1(b), where we reflect the back folding of high momenta surface modes into the reduced Brillouin zone.
An impulsive drive to the surface modes-Having established the structure of surface excitations, we turn to discuss the role of photocurrent.It lives on the surface of the material and results in the emission of radiation in both the air and sample [Fig.2].Motivated by the experiment [16], we turn to evaluate the frequency dependence of the electromagnetic field emitted into the air.We assume that J NL flows along the z-axis so that the magnetic field is oriented along the y-axis.Since J NL is confined to the surface of the superconductor, we can incorporate this drive through the generalization of the Fresnel formalism, where as above we solve Maxwell's equations inside the two media independently and then match the solutions using appropriate boundary conditions.One of them is the continuity of the tangential component of the electric field (follows from Faraday's law): E air,z = E mat,z at x = 0.The other equation is obtained from integrating the fourth Maxwell equation: The second term in Eq. ( 11) is parametrically small, so we neglect it and obtain: Given the possibility of having a non-zero Umklapp photocurrent [24], we write Eq. (1) as This insight that one can have J U ̸ = 0 is essential for understanding the experimental data, as we show below.
Using the form (13), we turn to compute the spectrum of outgoing radiation.To this end, we employ the same ansatz as in Eqs. ( 9)-( 10), except we now specialize on q z = 0.By invoking the derived boundary conditions, we arrive at: The coefficient α a , in turn, gives the amplitude of the emitted into the air radiation f (ω) -see Fig. 2.
In the absence of the Umklapp component J U = 0, we find that f (ω) displays a double peak structure, which comes from the splitting of the bulk JP resonance, encoded in the prefactors ε c /(ε 2 c − A 2 ).It is worth pointing out that since the surface excitations exhibit saturation around ω JPR , i.e., at the bulk resonance, it is not entirely clear that this splitting comes from the bulk rather than the surface.To resolve this question we consider in [41] isotropic superconductors, where the surface excitations are well separated from the bulk ones, and confirm that the double-peak structure originates from the bulk.We For ωJPR ≪ γc, we find that ESJP ∼ (ωJPR/γc) 3/2 for the surface scenario and EBJP ∼ (ωJPR/γc) 1/2 for the bulk one.
also further elaborate in [41] on the asymptotic behavior of f (ω) at both small and large frequencies.This doublepeak splitting was not observed in the experiment [16].
In addition, the resulting spectral function is too broad for realistic parameters to explain the experimental data.Most remarkably, provided J U ̸ = 0, we find the emission spectrum f (ω) becomes sharply peaked at around ω JPR .This peak originates from the fact that the Umklapp component now resonantly drives the surface plasmons at wave vectors around q z = Q; due to the back folding into the light cone, these plasmons can now radiate photons, as illustrated in Fig. 1(b).We further elaborate on the surface origin of this effect in [41].This result that the Umklapp photocurrent can give rise to a sharp emission provides a natural interpretation of the experimental data.
Bulk plasmon scenario-We now comment on the role of bulk polaritons.Since they are also being impulsively driven by J 0 ̸ = 0, Eq. ( 14), they can produce characteristic radiation peaked in frequency near ω JPR .While we cannot completely rule out the bulk scenario as a possible alternative interpretation of the experiment [16,41], there are two good reasons why the observed emission is dominated by the surface Josephson plasmons.The first one is that the resulting bulk emission spectrum turns out to be too broad for realistic parameters to explain the experiment [41].In contrast, in the surface scenario, the Umklapp photocurrent resonantly drives the surface modes resulting in sharp radiation without any fine-tuning of model parameters.The second argument invokes the experimental phenomenology that the observed radiation linewidth correlates with the out-ofplane CDW coherence length (once inversion symmetry is broken).This fact is inconsistent with the bulk scenario because the far-field reflectivity, which entirely determines the bulk emission properties, is known to be insensitive to the stripes [16,57].At the same time, the observed phenomenology can naturally be explained within the surface scenario since the corresponding emission is acutely sensitive to the CDW order (for an additional discussion, see [41]) -see also Fig. 3(a).
Finally, to further distinguish the two scenarios from each other, we theoretically predict the scaling of each response with temperature T .As T approaches T c from below the JPR softens to zero.As shown in Fig. 3, we find that the outgoing emission exhibits distinct behaviors for the two scenarios when ω JPR ≪ γ c : .
One important implication of this result is that the emission amplitude in the surface scenario E SJP is expected to become suppressed compared to the bulk one E BJP for T ≈ T c : E SJP /E BJP ∝ ω JPR /γ c .Since according to Eq. ( 14) both scenarios contribute to the outgoing radiation, we conclude that the surface scenario dominates at lower temperatures, while the bulk one might play a role in the immediate vicinity of T c .Conclusion-To summarize, our interpretation of the observed terahertz emission in striped superconductors consists of two related arguments.The first one is the downconversion of optical fields by fermionic quasiparticles from high-frequency high-intensity pump pulse to low-frequency regular and Umklapp photocurrents.Given that the emission is observed only in the striped phase, we expect that this downconversion arises from the CDW order [24].The second argument is that photocurrents impulsively drive collective modes of the sample.In particular, the Umklapp component resonantly couples to the surface Josephson plasmons, in turn resulting in spectrally sharp radiation at ω JPR and, thus, explaining the experiment [16].In the electrodynamics of striped superconductors, both stripe and superconducting orders combine and play significant roles.
For outlook, an intriguing open question is the microscopic origin of the photocurrents and whether it can be related to a pair density wave [17].The Umklapp photocurrent is interesting on its own, even without superconductivity, and might be relevant, for instance, in moiré systems [24,58].

Supplemental Material to "Theory for Anomalous Terahertz Emission in Striped
Cuprate Superconductors"

I. PLANE WAVES IN THE MEDIUM WITH STRIPES
Here we analyse plane waves in the sample with stripes.Because of the translational symmetry breaking along the z-axis, the structure of free harmonics contains the mixing between q z and q z + Q wave vectors; at the same time, along the x-axis, we can substitute an evanescent wave ansatz: where λ(ω, q z ) is yet unknown wave vector that depends on both ω and q z ; α and β are two coefficients that encode the strength of the mentioned mixing.Substituting this ansatz into the fourth Maxwell equation with the modified dielectric function, Eq. ( 8) of the main text, we obtain the electric field in the sample: Plugging this result into the third Maxwell equation, we get the following secular equation: where k m (q, ω) is determined implicitly through For consistency, when solving this equation, we choose the root with Re k m (q, ω) > 0. There are four eigenvalues of Eq. (S2): , where we defined km = k m (q z + Q, ω).We select the two of them with Re λ 1 , Re λ 2 > 0 that describe waves decaying into the sample.We choose these roots such that in the limit A → 0, we get λ 1 ≈ k m and λ 2 ≈ km .Equation (S2) also fixes the ratio between the amplitudes α and β:  1) does not depend on frequency ω, for ω ≲ ωJPR, we find that the amplitude of the radiation into the air displays a peak at ωJPR.This peak, though, is smeared out for realistic values of the quasiparticle damping γc.The inset illustrates that the photocurrent serves as a drive that emits radiation into the air and the superconductor.

II. BULK POLARITON SCENARIO
In this section, we consider a simplified phenomenological model, where we assume the presence of a broadband photocurrent, Eq. ( 1) of the main text, but we neglect stripes and do not take into account how they modify the dielectric properties of the sample.This is a bit artificial because the photocurrent generation itself is expected to occur due to the stripes [16].Nevertheless, this approach allows us to separate the response due to the bulk Josephson plasmons from the surface modes since the latter remain silent because they lie outside the light cone, as shown in Fig. 1(a) of the main text.In other words, the emission we compute below comes entirely due to the bulk excitations.
To get the emission properties in this model, we first solve the Maxwell equations in each media.In the air, we have an outgoing plane wave, which we write as: where B a (ω) encodes the amplitude of this plane wave.In the superconductor, the solution reads: where q m (ω) is yet an unknown wave vector that depends on ω, and B m (ω) is the amplitude of this medium harmonic.From the fourth Maxwell equation, we obtain the electric field in the sample: By plugging this into the third Maxwell equation, we get: q 2 m = ω 2 ε c (ω)/c 2 .Now, using the generalized Fresnel boundary conditions, derived in the main text, we finally obtain: This is the central result of this section, which we turn to discuss in detail.
Figure S1 shows the amplitude f (ω) of the emitted into the air radiation as a function of frequency ω.We find that even though the photocurrent is flat for ω ≲ ω JPR , f (ω) is peaked at ω JPR , i.e., the superconductor acts a low-frequency filter to the featureless drive J NL (ω).The reason for this peaked behavior is the enhanced due to the van Hove singularity density of states of bulk c-axis Josephson plasmons at ω = ω JPR .We now turn to discuss the dependence of f (ω) in more detail.First, at large frequencies ω ≳ ω JPR , where the photon dispersion inside the sample (bulk polariton) becomes approximately linear, the drive can emit light into both the air and the material.The fraction of how much it emits into the air is entirely determined by the relative speed of light c/c m = √ ε ∞,ab in the air to that in the superconductor.For instance, in case c = c m , half of the radiation would go into the air and the other half into the sample.Second, for small frequencies ω → 0, we find that f (ω) becomes suppressed and, in particular, f (0) = 0.The reason for this is that fundamentally, the superconductor can easily screen dc current, not allowing a low-frequency electric field to develop.We comment that the peak at ω JPR can be smeared out by finite γ c .For realistic parameters, this smearing effect is profound, and the resulting emission spectrum is too broad to explain the experiment [16].In contrast, the surface scenario provides sharp emission even if γ c is notable essentially because the Umklapp photocurrent resonantly drives the surface modes.
We turn to discuss an additional argument why the bulk scenario is inconsistent with the experimental phenomenology.In the experiment [16], the emission was detected in two LBCO samples corresponding to two different dopings.These samples have stripes with different relative strengths and out-of-plane coherence lengths, and the resulting emission was more intense and much narrower in frequency for the sample with stronger, more coherent stripes.The issue with the bulk scenario, however, is that the far-field reflectivity, which entirely determines the emission properties associated with the bulk modes, is known to be insensitive to the stripes [16,57].Therefore, the bulk plasmons are expected to produce similar radiation in the two samples, with the same lineshapes independent of the stripe correlation length and, thus, disagreeing with the experiment.In contrast, the emission from the surface Josephson plasmons is acutely sensitive to the CDW order, further justifying our expectation that the observed phenomenon comes from the surface modes.

III. SURFACE JOSEPHSON PLASMONS AT LARGE MOMENTA
In the main text, we considered bulk c-axis Josephson plasmons with flat dispersion, as encoded in momentumindependent electric permittivity (3).If one is interested in small momenta close to the light cone q ≃ ω/c, then this approach is well-justified.However, we also encountered in the main text the situation of Bragg mixing of surface modes with small momenta to those that have large momenta q ≃ Q CDW .For these latter modes, it might be essential to consider the momentum dependence of ε(ω, q), as we do in this section.On the phenomenological level, the leading correction to the flat caxis dispersion comes from the effects of quasiparticle compressibility, which, in turn, sets a very large Thomas-Fermi momentum scale q TF .We expect that in cuprates Q CDW ≲ q TF ≲ 2π/a 0 , implying that: i) Q CDW is still small enough so that one can essentially disregard momentum dependence of ε(ω, q) for q ≲ Q CDW ; ii) the lattice reciprocal momentum 2π/a 0 is large, and effects of compressibility cannot be ignored at such a momentum scale [59].In the literature, the dispersions of the usual surface plasmons and the surface Josephson plasmons in stacked 2D metals and anisotropic superconductors have been discussed in various geometries, including infinite half-space, thin film [60], and spherical particles [61].In our theory, we go beyond previous treatments of surface Josephson plasmons by including the effects of quasiparticle compressibility, which allows us to capture the dispersion of surface excitations at large momenta of the order of Q CDW .
We describe the electric response of superconductors within the two-fluid model [62], where the total electric current J = J n + J s is represented as a sum of the quasiparticle contribution J n and the one due to the superflow J s .For the normal component, we write Ohm's law but take into account the effects of electrochemical potential δµ = δρ/χ, where δρ represents the electric charge density and χ is the compressibility: For simplicity, we assume that the quasiparticle conductivity tensor σ = diag(σ ab , σ ab , σ c ) is frequency and momentum independent.For the superconducting part, we Dispersion ωs(qz) of the surface excitations for large momenta qz ≃ qTF, where the effects of quasiparticle compressibility become notable.In contrast to the behavior at low momenta near the light cone (inset), see also Fig. 1(a), where ωs(qz) displays a quick saturation near ωJPR, now ωs(qz) is no longer flat.In cuprates, we expect QCDW ≲ qTF ≲ 2π/a0.In this case, the stripes momentum QCDW is still quite small, ωs(QCDW) ≃ ωJPR, and the effects of compressibility are unimportant at this scale.At the same time, the lattice momentum 2π/a0 is large resulting in ωs(2π/a0) being substantially different from ωJPR.
Figure S2 shows the resulting dispersion ω s (q z ).While at low momenta, q z ≪ q TF , we reproduce the same behavior as in Fig. 1(a), i.e., ω s (q z ) shows a quick saturation near ω JPR , at larger momenta, ω s (q z ) becomes dispersive.In particular, for momenta q z ≳ q TF , we expect that ω(q z ) will be substantially different from ω JPR .The result in Fig. S2 implies that for the interpretation we offer in the main text to be consistent, in the experiment, we should have Q CDW ≲ q TF ≲ 2π/a 0 so that ω(Q CDW ) ≈ ω JPR .This latter requirement is essential from the perspective of back folding of the modes with momenta q z ≃ Q CDW , which should have frequency close to ω JPR .

IV. BULK VS SURFACE CONTRIBUTIONS TO THE SPECTRUM OF OUTGOING RADIATION
So far, we primarily studied strongly anisotropic superconductors with ω ab ≫ ω c .A unique feature of such materials is that the surface modes saturate at the frequency immediately below the bulk Josephson plasmon resonance ω JPR .This makes distinguishing whether the features in the spectral function of outgoing radiation arise from the surface or bulk modes difficult.To resolve this question, here we solve the problem for isotropic superconductors, where the surface plasmon frequency ω JPR / √ 2 is well separated from the bulk.Figure S3 shows the emission spectral function f (ω) in isotropic superconductors.For J U = 0, f (ω) displays a double-peak structure near the bulk resonance: the dip that appears at ω JPR corresponds to the bulk plasmon splitting in the striped phase due to the hybridization between zero and finite momentum bulk modes.Importantly, once J U ̸ = 0, f (ω) exhibits an additional sharp peak at the surface resonance.This confirms both that the sharp resonance in Fig. 2 is due to the surface Josephson plasmons and that surface Josephson plasmons can only be resonantly excited by a nonzero Umklapp photocurrent.We remark that the bulk features are not severely renormalized by J U .We conclude that the zero momentum photocurrent J 0 drives bulk modes and leads to a broad emission, while the Umklapp photocurrent J U primarily drives surface excitations and leads to a sharp peak at the surface plasma resonance -see the inset in Fig. S3 for a summary.

V. REFLECTION FROM A MEDIUM WITH BRAGG MIXING
Here we compute the reflection coefficient r p of the medium with stripes.Following the discussion in the main text, cf.Eqs. ( 9)- (10), we write the magnetic field as (q z = 0): Now, the coefficients α a , β a , α m , β m are related to each through the appropriate Fresnel boundary conditions, which in turn give an implicit equation on r p : det N = 0, where

(S25)
For A = 0, we reproduce the familiar expression: We compare the usual Josephson plasma edge appearing in the reflectivity spectra of a layered superconductor without stripes to reflectivity spectra renormalized by the stripes in Fig. S4 (left).We find that for A ̸ = 0, the Bragg mixing between zero momentum, q z = 0, and finite momentum, q z = Q CDW , gives rise to a plasma edge splitting.As shown in Fig. S4 (right), where following the preceding section, we consider isotropic superconductors and confirm that the new features in Fig. S4 (left) come from the bulk.In other words, far-field measurements, even in the presence of the CDW order, are sensitive to the bulk and not the surface excitations.However, up to relatively large stripe amplitudes, A ∼ 0.1ϵ ∞,c , the new features appear negligible and can be easily missed within the experimental precision.This is consistent with reflectivity experiments in stripes which have reported only a single plasma edge.For this reason, in our calculations for terahertz emission, we use the value of A = 0.1ϵ ∞,c as a reasonable upper limit to the stripe amplitude.
Notably, we find that the coupling to finite momentum modes can drastically enhance terahertz emission even though reflectivity probes are essentially insensitive to these modes up to large values of the stripe strength A.

scattering a b FIG. 1 .
FIG. 1.(a) Dispersion ωs(qz) of the surface Josephson plasmons.In strongly anisotropic superconductors, ωs quickly saturates at ωJPR.Dashed line corresponds to the light cone ω = cqz.(b) Schematic of the surface plasmon spectrum in the presence of a weak stripes order, which couples momenta qz to qz + QCDW.As such, the plasmon dispersion exhibits back folding to the reduced Brillouin zone, defined by the CDW wave vector QCDW.Notably, the surface plasmons that back fold inside the light cone can now radiate out.
FIG.2.Emission spectrum.We find that the outgoing radiation can be sharply peaked around ωJPR only in the presence of a nonzero Umklapp photocurrent JU ̸ = 0, Eq. (13).This Umklapp component impulsively drives high momenta surface Josephson plasmons that can emit light, Fig.1(b).Here we fixed A = 0.1ε∞,c.
FIG. S1.Bulk polariton scenario.While the photocurrent in Eq. (1) does not depend on frequency ω, for ω ≲ ωJPR, we find that the amplitude of the radiation into the air displays a peak at ωJPR.This peak, though, is smeared out for realistic values of the quasiparticle damping γc.The inset illustrates that the photocurrent serves as a drive that emits radiation into the air and the superconductor.

!
FIG.S3.Emission spectrum f (ω) in a system with ω ab = ωc = 1 THz, γ ab = γc = 0.1 THz, ε ∞,ab = ε∞,c = 1, A = 0.1.In such isotropic superconductors, bulk and surface resonances have notably different energies.This allows us to conclude that the homogeneous component J0 of the photocurrent drives the bulk plasmons only, while the Umklapp one JU affects the surface plasmons, in turn, providing a sharp resonant structure of f (ω).
FIG.S4.Reflectivity R(ω) in anisotropic (left) and isotropic (right) superconductors with stripes.Provided the amplitude A of the CDW order is nonzero, there occur new features in R(ω) near the bulk plasmon resonance ωJPR.At the same time, as indicated in the right panel, no features occur at the surface resonance, implying that the far-field photons are coupled only to the bulk plasmons.Parameters in the left (right) panel are the same as in Fig.2(Fig.S3).