Photon-Dark Photon Conversion with Multiple Level Crossings

Dark photons can oscillate into Standard Model (SM) photons via kinetic mixing. The conversion probability depends sensitively on properties of the ambient background, such as the density and electromagnetic field strength, which cause the SM photon to acquire an in-medium effective mass. Resonances can enhance the conversion probability when there is a level-crossing between the dark photon and background-dependent SM photon states. In this work, we show that the widely used Landau-Zener (LZ) approximation breaks down when there are multiple level-crossings due to a non-monotonic SM photon potential. Phase interference effects, especially when the dark photon mass is close to an extremum of the SM photon effective mass, can cause deviations from the LZ approximation at the level of a few orders of magnitude in the conversion probability. We present an analytic approximation that is valid in this regime and that can accurately predict the conversion probabilities in a wide range of astrophysical environments.


I. INTRODUCTION
Dark photons (DPs) are the gauge bosons of a hidden U (1) ′ symmetry that may be spontaneously broken (analogous to the Higgs mechanism) or broken explicitly (with the Stuckelberg action for a spin-1 field), yielding a non-zero dark photon mass.In general, DPs will mix with Standard Model (SM) photons at some level via the dimension-4 kinetic mixing operator.The additional terms in the Lagrangian capturing the effects of the DP field A ′ mixing with the SM photon field A are where F ′ and F are the dark and SM field strength tensors, ϵ ≪ 1 is the dimensionless coupling that determines the strength of kinetic mixing with the SM photon, and m A ′ is the DP mass.Kinetic mixing has attracted considerable interest as a portal to dark sectors that would be accessible at a wide range of energies since it is a marginal operator.In particular, ϵ could be generated by any number of mechanisms [1][2][3][4][5][6][7][8][9][10][11], for instance by loop diagrams with heavy matter fields charged under both U (1) ′ and the SM U (1), or alternatively from certain realizations of string theory.From a bottom-up effective field theory point of view, ϵ can be thought of as a free parameter accompanying a marginal operator.Astrophysical and cosmological observations provide some of the tightest constraints on the existence of DPs over a wide range of masses [12][13][14][15][16][17][18][19][20][21].These constraints often involve oscillations between DPs and SM photons leading to spectral imprints or new energy loss mechanisms.Furthermore, this conversion is resonantly enhanced if there is an energy level-crossing due to inmedium contributions to the SM photon effective mass.Large corrections to the SM photon mass arise from ambient free charges or strong background electromagnetic fields [22,23], and corrections to the DP effective mass can also arise due to a background of dark sector particles directly charged under the DP [24,25].
In-medium effective masses can vary considerably over spatial or temporal domains in astrophysical and cosmological environments.In the case of non-adiabatic transitions due to a level-crossing, the conversion probability is usually well-approximated by the Landau-Zener (LZ) formula [26,27], which is often employed in such calculations as the leading-order term in the stationary phase approximation.However, this formalism breaks down when applied to the special case of resonant conversions occurring near local minima or maxima of the SM photon's in-medium mass (which can depend on the photon frequency).This breakdown therefore occurs generically in environments where the density of SM particles is non-monotonic in space or time as traversed by SM photons.The conversion probability predicted by the LZ formula can deviate from the full solution to the Schrödinger equation (to leading order in ϵ) by several orders of magnitude.Therefore, there may be regions of DP parameter space at specific masses where existing constraints are subject to large corrections.In this work, we provide a simple analytic expression, Eq. ( 21), that is more generally applicable to this region of parameter space and discuss the potential impact on DP oscillation probabilities.We note that although we concretely work with DP oscillations in this paper, the formalism developed here can potentially be applied to axion-photon conversions as well as neutrino oscillations.
The rest of this paper is organized as follows.In Section II, we review various formalisms for computing transition probabilities for two-state systems.In Section III, we quantify the breakdown of the LZ approximation near the critical mass with an illustrative toy example.In Section IV, we examine the consequences of the breakdown of the LZ approximation in a few case examples arising in astrophysical and cosmological settings, such as neutron star magnetospheres and the solar chromosphere, as well as during the epoch of reionization.Concluding remarks follow in Section V.

II. PHOTON-DARK PHOTON OSCILLATION FORMALISM
A. Photon-dark photon oscillations in an arbitrary in-medium potential The kinetic terms in Eq. ( 1) can be made canonical by moving to the active-sterile basis described by the A a and A s fields, to leading order in ϵ.Using this in Eq. ( 1), we then have where F µν a and F µν s are the active and sterile field strengths and we have also included the coupling of A a to the SM current density J.This basis is often referred to as the interaction basis, since SM currents selectively couple to the active state A a .However, note that the mass matrix is non-diagonal in this case.As a result, although the sterile state is not sourced directly from SM currents, it arises indirectly from A a ↔ A s oscillations, which are the focus of this work.
In a medium, the free charges and electromagnetic fields in the background can source a potential that alters the dispersion relation of the active photons.We parameterize this effect by including a spatially-dependent inmedium mass m eff (x) for the visible field A a in Eq. ( 3).The in-medium mass-squared matrix is then given by In this Section, we remain agnostic about the specific form of m eff (x), leaving an exploration of specific examples to Section IV.We can track the propagation of the A a − A s system by solving the corresponding equation of motion from Eqs. (3) and (4).Following Ref. [28], we switch to Fourier space, where ω and k are the frequency and wavenumber of the field A µ , respectively, and assume that m eff (x) varies on scales much larger than k −1 .In this case, we can approximate this in-medium contribution as a constant on the scale of the de Broglie wavelength, such that the equation of motion for transverse modes takes the form of a standard wave equation, We note that an analogous dispersion relation of the form (ω 2 − M 2 L ) A L = 0 holds for longitudinal modes, but a dedicated study of their evolution is beyond the scope of this work.
To proceed, we expand Eq. (5) in the relativistic limit ω ≃ k ≫ m A ′ , m eff to obtain a linearized Schrödingerlike equation, as in Ref. [28].We choose to work with the spatial domain marked by the position z, but can equivalently use the temporal domain since we deal with the propagation of relativistic particles in this work.This gives i∂ z A = HA, where the total Hamiltonian H = H 0 + H 1 is split into diagonal and off-diagonal components, with Since ϵ ≪ 1, we can approximate H 1 ≪ H 0 and use the techniques of time-dependent perturbation theory to solve for the evolution of A. In particular, we switch to the interaction picture where , and U is defined to be such that z i marks the point at which we fix our initial condition A(z i ) = A int (z i ).Hence, the system evolves as A int (z) = e −i z z i dz ′ Hint(z ′ ) A(z i ), which in the Schrödinger picture is equivalent to The first factor of Eq. ( 9) is determined by the definition of H 0 in Eq. (6).To proceed, we evaluate H int using the Baker-Campbell-Hausdorff identity, where we have defined This form of H int in Eq. ( 10) can be used in the second factor of Eq. ( 9), after expanding to O(ϵ).Up to an irrelevant overall phase, this yields where we defined The probability of conversion between active and sterile states is then given by the square of the off-diagonal elements in the above expression, In vacuum, there is no in-medium contribution to the active photon, ∆ = 0, and the above integral can be performed analytically, yielding the standard result More generally, for ϵ ≪ 1 and ∆ ̸ = 0, one can numerically integrate Eq. ( 13) in order to calculate the in-medium conversion probability.In this case, the integrand in Eq. ( 13) is very oscillatory, making it difficult to achieve a high degree of numerical accuracy.Near resonances, however, it is possible to accurately approximate P Aa↔As analytically, as discussed in the next Subsection.

B. Stationary Phase Approximation and the Landau-Zener Formula
While the integrand in Eq. ( 13) typically exhibits highly oscillatory behavior, it varies slowly near stationary points z n , defined by Φ ′ (z n ) = 0, where the "prime" corresponds to a spatial derivative.Since the phase varies slowly in this region of coordinate space, its contribution to the integral is not cancelled out by other regions' contributions, where oscillations tend to interfere destructively and average to zero.From Eq. ( 11), it is evident that these stationary points occur when the resonance condition holds, ∆ = ∆ A ′ , analogous to level-crossings induced by matter effects within the context of neutrino oscillations [29,30].In this case, Eq. ( 13) can be evaluated analytically by use of the stationary phase approximation, in which Φ(z ′ ) and ∆ A ′ (z ′ ) in the integrand are Taylor expanded to second and zeroth order around z ′ ≃ z n , respectively.This gives where the sums are over all stationary points and in the first and second lines we have defined and respectively, where from Eq. ( 11) we have Φ ′′ (z n ) = ∆ ′ (z n ).In the second line of Eq. ( 15), the first term is simply the sum of individual probabilities from each resonance.The second term instead stems from the interference between two different resonances, which gives rise to what we will refer to below as phase effects.Such phase effects imprint oscillatory behavior into the conversion probability as a function of m A ′ and ω.
As shown schematically in Fig. 1, the conversion probability undergoes a large number of oscillations between any two resonances if |Φ nm | ≫ 2π, or equivalently where we have been explicit in regards to the spatial and frequency dependence.In this case, small variations in frequency within an observed resolution bandwidth could lead to variations in Φ nm that are much larger than 2π, such that the phase in the second term of Eq. ( 15) averages to zero within that frequency bin.A related possibility is if the spatial profile of m eff (z) varies between slightly different lines of sight within the angular resolution of the detector.Again, in the limit where the total acquired phase is large, small variations in the m eff (z) profile can lead to large phase variations within the resolution of the detector.This would potentially cause phase effects to average to zero within an angular bin, depending on the underlying spatial distribution of m eff .An additional factor in the loss of phase information is the coherence length, which determines the ability to maintain phase coherence between different mass eigenstates.This length is the distance over which wavepackets of different mass eigenstates seperate by a distance greater than the wavepacket width, and its dependence on the production process has been investigated extensively in the context of neutrinos [31,32].If the distance between resonance regions exceeds the coherence length, then different mass eigenstates decohere before reaching the subsequent resonance region, leading to the loss of phase information.For the purposes of this work, we assume that the states remain coherent, as we are working with relativistic DPs and do not assume that the source of SM photons or DPs is spatially or temporally localized; we will revisit these assumptions and their relevance to astrophysical environments in future work.
In the literature, phase effects are often assumed to average to zero.Ignoring the corresponding cross terms in Eq. ( 15), the conversion probability reduces to the standard LZ result [26,27] for a non-adiabatic two-level transition, Note that this approximation is valid only when the above expression remains less than unity, which would otherwise seemingly violate unitarity.In this work, we always deal with conversion probabilities P Aa↔As ≪ 1 such that multiple sequential conversions are highly suppressed.

C. Breakdown of Landau-Zener
If m A ′ happens to be near a local extremum of m eff (z), there exists a pair of resonant points z 1,2 such that z 1 ≲ z c ≲ z 2 and z 1 ≃ z 2 ≃ z c (see Fig. 1).For a particular value of m A ′ , such points coalesce near the position of the extremum, z c .We refer to this mass m A ′ = m eff (z c ) ≡ m c as the critical mass.At this critical point z c , the first two derivatives of Φ(z) both vanish due to the resonance appearing at the extremum of the m eff potential.Near z c , the contributions from the saddle points A n ∝ |Φ ′′ (z n )| −1 in Eq. ( 15) diverge.As a result, neither Eq. ( 15) or the LZ approximation in Eq. ( 19) are valid for m A ′ ≃ m c .
Instead, analytically determining the transition probability when m A ′ ≃ m c requires incorporating the cubic term in the Taylor expansion of Φ(z) around z c in Eq. ( 13), since Φ ′ (z c ) = Φ ′′ (z c ) = 0. We therefore define a dimensionless variable that quantifies the relative contribution of the quadratic term over the cubic term in the Taylor expansion.Hence, for ξ ≪ 1 the quadratic term becomes negligible and one should expand to include Φ ′′′ (z c ) in order to accurately evaluate the conversion probability near the critical mass.We incorporate this higher-order term by making use of well-known results for two "coalescing saddle points" [33,34], yielding which is a key result of this work and is only valid for values of m A ′ ≃ m c such that ξ ≲ 1.In Eq. ( 21), the entire expression is evaluated at the critical point z c , Ai(x) is the Airy function, and ζ ≡ σ(z 1 ) (2/|Φ ′′′ |) 1/3 Φ ′ .We note that this formula is valid for two resonant crossings near an extremum; incorporating additional resonance points is possible with the approximations of Ref. [34].

III. NON-MONOTONIC TOY POTENTIAL
In this Section, we begin to quantify the degree to which various treatments of DP-SM photon transition probabilities can differ.Motivated by the form of the Taylor expansion of m eff about its extremum, we consider the following quadratic toy model for the potential, where we have chosen to take the width of the peak near z c to be O(z c ) so that there is only one relevant length scale in the problem.For m A ′ < m c , there are two resonant level crossings where m A ′ = m eff .Because any potential will take a similar quadratic form near its extremum, we can use this toy model as a proxy to gain intuition for potentials in astrophysical systems, such as the ones discussed in the next Section.Conversion probabilities using Eq. ( 22) computed with different methods are shown in Fig. 2 as a function of As can be seen in the left panel of Fig. 2, when m A ′ is far from the critical mass m c , the approximation in FIG.2: Left: The conversion probability PA a ↔As as a function of δm (a dimensionless measure of the proximity of m A ′ to the critical mass mc) for the toy potential of Eq. ( 22) with mc zc = 2 × 10 4 , ω/mc = 10 2 , and ω zc = 2 × 10 6 .The region to the left of the dashed vertical line corresponds to ξ ≲ 1 (see Eq. ( 20)), indicating the expected breakdown of standard approximations.The various approximations, consisting of the LZ result of Eq. ( 19), adding the term that captures phase effects in Eq. ( 15), and our approximation in Eq. ( 21), are compared to a numerical evaluation of Eq. ( 13) and the vacuum conversion probability of Eq. ( 14).Right: The ratio between the full conversion probability computed numerically and the vacuum oscillation probability as a function of ω/mc for fixed mc zc = 2 × 10 4 and for different values of δm.For larger frequencies, the conversion probability reduces to its vacuum value as the transition becomes increasingly non-adiabatic.
Eq. ( 15) (labelled "phase effects") matches well with the numerical evaluation of Eq. (13) (labelled "numerical").In this case, the standard LZ result of Eq. ( 19) (labelled "LZ") accurately captures the typical value of the transition probability, averaged over the oscillatory features with varying m A ′ .As m A ′ approaches m eff from below, the two resonance points converge spatially, eventually merging at the critical point when δm = 0, corresponding to ξ → 0 in Eq. (20).As discussed in the previous Section, the approximations detailed in Eqs. ( 15) and (19) are no longer valid for ξ ≲ 1 (left of the vertical dashed line).Nevertheless, our approximation in Eq. ( 21) (labelled "our approx.")remains accurate and proves to be a reliable method of tracking the conversion probability near this critical mass.We emphasize that these approximations are substantially faster to evaluate numerically as compared with the full numerical solution of the Schrödinger equation, primarily due to the oscillatory nature of Φ(z).
In-medium resonances significantly amplify the conversion probability with respect to the vacuum value of Eq. ( 14).To encapsulate the enhancement due to such in-medium resonances, we first note that Eq. ( 16) can be rewritten as A n = 2π µ n , where we have introduced a dimensionless resonance enhancement parameter that is analogous to the adiabaticity parameter of Keldysh in Ref. [35].Since P Aa↔As ∼ ϵ 2 max A n in the LZ approximation, µ ≡ max µ n effectively quantifies the degree that any such resonance enhances the conversion probability over the vacuum value.Hence, for µ ≪ 1, we expect the in-medium modifications to the transition probability to be suppressed, such that P Aa↔As approaches the vacuum value of Eq. ( 14).
This form of the resonance enhancement parameter µ is to be expected.To see this, note that for a nonmonotonic potential, such as the one of Eq. ( 22), µ n ∼ L res /L vac , where L res is the separation between z n and a nearby resonant point, and L vac ∼ 1/∆ A ′ is the vacuum oscillation length.Hence, for µ ≪ 1, the potential as seen by the DP undergoes a spatially abrupt change compared to the oscillation length, resulting in an extremely non-adiabatic transition where it is valid to use the "sudden approximation."Since by the uncertainty relation the DP can only resolve length scales greater than 1/∆ A ′ ∼ L vac for a momentum difference between active and sterile states of δk = ∆ A ′ , it cannot effectively discern the presence of the potential over a distance of L res ≪ L vac .Conversely, for µ ≫ 1, the potential varies over long distances compared to the vacuum oscillation length, indicating a need to account for m eff using one of the other approximations detailed in the previous Section.
The validity of the sudden approximation (i.e., using the vacuum transition probability) for µ ≪ 1 is verified numerically in the right panel of Fig. 2, which compares numerical evaluations of P Aa↔As with Eq. ( 13) to the vacuum value in Eq. ( 14), as a function of ω (expressed in terms of the dimensionless quantity ω/m c ) and for various representative choices of the DP mass in the form of δm.As ω increases, the vacuum oscillation length also increases compared to the effective width of the m eff profile, such that the conversion probability progressively approaches its vacuum value.

IV. ASTROPHYSICAL AND COSMOLOGICAL ENVIRONMENTS
In this Section, we provide a few relevant examples of environments with non-monotonic effective photon masses where it is necessary to use Eq. ( 21) in order to accurately evaluate the DP conversion probability near the critical mass.It may be necessary to update some astrophysical and cosmological bounds on DPs for particular DP masses in light of these considerations.

A. Neutron Star Magnetospheres
Among astrophysical environments, neutron star (NS) magnetospheres are the most extreme in their variations in free charge density and electromagnetic field strength; they therefore pose an environment where m 2 eff can vary substantially, affecting conversion of DPs to SM photons [36].In particular, the rotating magnetic fields source electric fields that exert forces much larger than the gravitational binding energy on the NS surface charges.As a result, charges are pulled off of the surface and fill the magnetosphere of the NS.The charges then redistribute themselves in a way such that the Lorentz force on them cancels, producing a corotating plasma surrounding the NS.The charge density of this plasma can be approximated by the Goldreich-Julian model [37], which has been generalized to account for relativistic effects [38] n Here, r and θ are polar coordinates such that r is the distance from the center of the NS, θ is the polar angle with respect to the NS's rotation axis Ω, B is the magnetic field outside the NS (assumed to be dominated by its dipole component), and Ω is the rotational frequency.The functions and incorporate relativistic corrections, where we have defined r = r/r g and β = 2/5(r NS /r g ) 2 for a NS of radius r NS and Schwarschild radius of r g .Note that at large distances, lim r→∞ F 1 (r) = −1/2r and lim r→∞ F 2 (r) = −2/3.Following Ref. [39], we take the rotation axis to be aligned with ẑ, such that and where B 0 is the magnetic field strength at the NS surface, α B is the orientation angle of the magnetic field with respect to the rotation axis, and r LC = 1/Ω is the lightcylinder radius.
In the presence of large magnetic fields (with e /e is the critical value of the magnetic field), both the magnetic field and the plasma contribute to m eff for the propagating photon mode (that has its electric field parallel to the plane containing the propagation and NS magnetic field vectors) but with opposite signs [36,40,41].We can decompose the effective SM photon mass as m 2 eff = V B + V pl , where b = B(r)/B c , ω 2 pl = 4παn GJ e /m e is the plasma frequency, and qB is a fitting function that reproduces the correct b ≪ 1 and b ≫ 1 behavior [41] qB The frequency-dependent B-field contribution to the photon potential V B originates from non-linear vacuum polarization effects (analogous to the Euler-Heisenberg term in pure QED) and contributes effectively only when B(r) ≳ B c .The opposing signs of V B and V pl give rise to a non-monotonic profile for m 2 eff (r) .Note that the potential is non-monotonic regardless of the inclusion of relativistic corrections in Eq. (25).
The exact form of the SM photon potential m 2 eff (r) is sensitive to various NS parameters, such as the magnetic field strength, rotation speed, and frequency ω, resulting in a wide range of possible profiles.Since the electron density and magnetic field both scale as ∝ 1/r 3 far from the NS's surface, the potentials of Eq. ( 30) scale as V B ∝ −r −6 and V pl ∝ r −3 .As a result, m 2 eff < 0 when V B dominates at intermediate distances (which precludes the possibility of a resonance in this region) and m 2 eff > 0 FIG.3: Left: Contours of the critical mass mc in the environment of a neutron star, in the plane spanned by the rotation period P of the NS and frequency ω for a fixed value of the magnetic field B0/Bc = 10.The critical mass varies considerably for fixed neutron star properties, meaning that any conversion processes occurring over a range of frequencies would only be accurately described by the approximations developed in this work.Right: As in Fig. 2, but for the neutron star potential of Eq. ( 30) with rLC = 300 km, B0/Bc = 10, and ω = 0.08 eV.The similarity of this behavior with that of the toy model underscores the necessity of using the approximation near the critical mass.
when V pl dominates at large distances.Hence, the potential reaches its extremum at the turnover point r c when V pl (r c ) ∼ −V B (r c ).For B 0 ≲ B c , the corresponding critical radius and mass take simple analytical forms, As discussed in Section II C, when the dark photon is close to but slightly less than the critical mass, there exist two nearby resonant points r 1 ≃ r 2 such that r 1 ≲ r c ≲ r 2 .These resonances are physically realized only if the critical radius r c lies in the range r NS ≤ r c ≤ r LC , which is equivalent to The left panel of Fig. 3 shows that the critical mass m c spans many orders of magnitude for a representative range of possible values for the light-cylinder radius r LC and frequency ω, fixing B 0 /B c = 10 and using the full form for the potential m 2 eff .In the right panel, we show the conversion probability as a function of δm, using the various approximations discussed in Section II, for a particular choice of NS parameters.By comparing to Fig. 2, we note that our findings qualitatively resemble the probabilities for the toy example of the previous Section, further illustrating that Eq. ( 21) is required to accurately compute P Aa↔As for δm ≪ 1.

B. Intergalactic Medium
The free electron fraction of the intergalactic medium (IGM) has fluctuated considerably since recombination, diminishing during the dark ages and cosmic dawn before increasing again during reionization, all while the density of the Universe was redshifting as ∼ (1 + z) 3 .As a result, cosmic microwave background (CMB) photons experience large variations in their effective mass as they traverse the IGM, where ω p,e = 4πα n e /m e and ω p,HI = 4πα n HI /m e are contributions to the plasma frequency from the average number density of free electrons n e and electrons bound in neutral hydrogen n HI , respectively.As shown in the top-left panel of Fig. 4, the resulting SM photon potential has two local extrema independent of frequency near z ∼ 10, corresponding to just before and after reionization, as well as additional local maxima for frequencies significantly larger than the CMB blackbody temperature, i.e., ω/T ≳ few [17,43] .The corresponding values of the critical mass are roughly m c ∼ few × 10 −13 eV, m c ∼ 10 −14 eV, and m c ∼ 10 −11 eV × (T /ω) 2 .We note that the first two of these are subject to uncertainties in the ionization history.We postpone a more careful consideration of alternative parametrizations of reionization [44] as well as the incorporation of fluctuations in the plasma density (along the lines of Refs.[18,19]) to future study.Previous studies have derived strong constraints on DPs from requiring that the CMB spectrum remains a FIG.4: Top left: The spatially averaged evolution of m eff in the intergalactic medium as a function of redshift for different frequencies.Reionization induces extrema in the potential experienced by CMB photons.Top right: As in Fig. 2, but for the IGM potential shown in the top-left panel for a frequency of ω = 7 × 10 −4 eV and mc ≃ 2.5 × 10 −13 eV.Bottom: As in Fig. 2, but with two values of fixed dark photon masses and instead varying the frequency.Also shown as a grey vertical band is the frequency resolution of FIRAS, 24.6 GHz [42].For values of the dark photon mass near the critical mass, the frequency oscillations do not necessarily average out within the frequency resolution of the instrument, meaning that phase effects are potentially significant in determining the observed conversion probability.Note that in the bottom panels, the values of δm correspond to ξ > 1 where we expect our approximation is not necessary to compute the conversion probability.
blackbody, as resonant conversion from CMB photons to DPs could induce obervable spectral distortions in the frequency range probed by the Far InfraRed Absolute Spectrophotometer (FIRAS) [45].Accounting for both the mean value and fluctuations of the plasma density has led to strong constraints at the level of ϵ ≲ 10 −7 for 10 −15 eV ≲ m A ′ ≲ 10 −6 eV [17][18][19]21].Such studies have employed the LZ approximation discussed above, but this is expected to break down for m A ′ ≃ m c , necessitating the use of Eq. ( 21).This is demonstrated in the top-right panel of Fig. 4, which shows the conversion probability for DP masses near the critical point as com-puted using different approximation schemes, analgous to the toy example of Fig. 2. Given the form of the mean photon potential (i.e., ignoring plasma density fluctuations), DPs with 10 −13 eV ≲ m A ′ ≲ 10 −14 eV have three or more resonance points, with two critical points occurring around the time of reionization.It is noteworthy that both ω(t) and Φ ′′ (t) decrease steeply with the age of the Universe, meaning that resonances that occur later are more adiabtic.As a result of the ∼ 1/(ω 2 Φ ′′ ) scaling of the A n in Eq. ( 16), resonances occurring at later times therefore have a much greater contribution to the total conversion probability as computed with the LZ formal-FIG.5: Left: The non-monotonic plasma frequency ωp and temperature T profile above the solar photosphere.The solar transition region which separates the chromosphere and corona is also shown.Right: As in Fig. 2, but with the chromospheric potential depicted in the left panel.ism of Eq. ( 19).We therefore ignore the earliest level crossings that occur well before reionization in computing the conversion probability and include only the ones occurring around the time of reionization.
The conversion probability fluctuates sharply as a function ω and m eff due to the phase effects discussed in Section II B. Effects from variations in ω are potentially observable depending on the frequency resolution of the detector, which is shown for FIRAS as the vertical gray region in the bottom row of Fig. 4. Variations in m eff arise from the fact that different lines of sight trace slightly different ionization histories due to the patchy morphology of reionization, leading to shifts in the SM photon potential and conversion probability along slightly different lines of sight.We note that this opens the possibility for the constraints on DPs from CMB spectral distortions to form more of a "fog" (due to theoretical uncertainties in computing the conversion probabilities) rather than a sharp exclusion boundary in parameter space.We leave the question of whether or not these variations average out across different lines of sight to future work.We also note that these considerations may be relevant for the dark screening effect that was recently proposed in Ref. [46], since there the photons passing through halos containing a gas overdensity will have two resonance points near the maximum plasma frequency.

C. Solar chromosphere
In stellar chromospheres, the density of free electrons n c e exhibits a distinct peak located approximately 10 3 km above the photosphere, as depicted in Fig. 5 [47].This non-monotonic electron density translates to a non-monotonic m 2 eff = 4πα(n c e ) 2 /m e .Additionally, the form of the profile yields a small value of ξ for a wide range of frequencies if the DP mass is near the critical mass.For instance, for δm ∼ 0.01, ξ ≲ 1 for ω ≳ 0.1 eV.This implies the unavoidable need for the approximation of Eq. ( 21) to accurately assess the conversion probability shown in Fig. 5.This may be important to incorporate for DP searches involving the Sun, for instance strategies that could complement existing searches involving resonant level crossing in the solar atmosphere (e.g., Ref. [48]).We explore this possibility in future work.

V. CONCLUSION
In this paper, we have developed an accurate approximation for the solution to the Schrödinger equation for a two-state system (specifically, photons and dark photons) with multiple nearby resonances about the extremum of a potential.This approximation can be viewed as an extension of the LZ formula, which is widely used for computing the transition probability between photons and dark photons.Using a toy model, we find that there can be large corrections to the conversion probability for specific dark photon masses in a given potential (i.e., given some properties of the background environment that affect the propagation of photons).
We have highlighted the application of this formalism to various astrophysical systems where multiple resonances between dark photons and photons are possible.Some of these systems have been previously used to constrain the existence of dark photons through their observed spectral signatures, which makes it important to quantify the effects of multiple resonances on the con-straints.We have not attempted to update any of these constraints from specific astrophysical systems, and leave such an analysis to future work.We note that applying our formalism to neutron stars seems especially promising due to the large resonant enhancement factors over the vacuum conversion probability as well as the wide range of possible values of the critical mass.
There are several additional subtleties in the conversion between photons and dark photons that we have not considered here, such as the role of decoherence between different states, which may be particularly relevant when the dark photon is non-relativistic (e.g., in the case of dark photon dark matter).We note that much of the formalism here, along with the associated subtleties, may be transferred over to axion-photon conversions and neutrino oscillations.We leave consideration of all of these effects to future work.

FIG. 1 :
FIG.1:A schematic representation of a peaked m eff -profile (gray line) with the distance Lres between the resonances (black dots) marked for a specific dark photon mass.The critical point is marked by the black "×" at the top of the potential.Left: A depiction of how the conversion probability oscillates spatially, where Lvac is the relevant conversion length scale in vacuum.Right: Propagating photon and dark photon wave packets are depicted in red and blue, respectively, for the two mass eigenstates, with the coherence length L coh also shown.