The Echo Method for Axion Dark Matter Detection

: The axion is a dark matter candidate arising from the spontaneous breaking of the Peccei– Quinn symmetry, introduced to solve the strong CP problem. It has been shown that radio/microwave radiation sent out to space is backscattered in the presence of axion dark matter due to stimulated axion decay. This backscattering is a feeble and narrow echo signal centered at an angular frequency very close to one-half of the axion mass. In this article, we summarize all the relevant results found so far, including analytical formulas for the echo signal, as well as sensitivity prospects for possible near-future experiments.


Introduction
The nature of dark matter is one of the most intriguing puzzles of modern science. The lack of signals from thermally produced WIMP dark matter suggests the exploration of new routes. One possibility is hypothetical particles produced in the early universe by nonthermal means. Such is the case for sterile neutrinos [1], dark photons [2], and QCD axions [3][4][5] (or axion-like particles (ALPs) [6]), which have became very popular in the last few years.
The QCD axion is a hypothetical particle that was originally postulated as a solution to the strong CP problem of the Standard Model of particle physics, i.e., the puzzle of why the strong interactions conserve the P and CP symmetries [3]. The QCD axion and ALPs are leading dark matter candidates in the sub-eV mass range [7][8][9][10]. Their nonthermal production mechanisms include the realignment and decay of topological defects.
Several running experiments are carrying out the search for axions, and there are also many new interesting proposed ideas. See [11,12] for recent reviews. Most of these detection schemes are based on the axion to two-photon coupling [13].
Thanks to this coupling, a massive axion is allowed to decay spontaneously into two photons, and more interestingly, the decay rate is enhanced considerably in a photon background with photon energy close to half the axion mass. This stimulated axion decay has been discussed extensively in the recent literature [14][15][16][17][18][19][20]. In the axion rest frame, the photons ejected from the axion decay are forced to propagate along the trajectory of the incident photon. While one of the decaying photons is released together with the incident one, the other must be released exactly in the opposite direction due to the conservation of momentum. Reference [21] discussed the possibility of sending out to space a powerful beam of microwave (or radio) electromagnetic radiation as an axion decay stimulator and then "listening" for the incoming radiation, the echo. The power signal of the echo was estimated, giving promising results for radio astronomy searches. The atmosphere's transparency mainly constrains the frequency range for which the echo method is realistic.
At sea level, this corresponds approximately to the range 30 MHz-30 GHz, equivalent to axion masses between 2.5 × 10 −7 eV and 2.5 × 10 −4 eV. Of course, this is just to have an idea. At a certain altitude, the method could be extended even to infrared or visible light, depending on the telescope's frequency operation, location, weather, noise temperature, etc. A more detailed and exhaustive analysis of the echo signal can be found in [22].
The echo intensity has three relevant properties: (1) it depends strongly on the local axion dark matter momentum distribution; (2) it does not depend on the shape of the outgoing beam; (3) its frequency is shifted with respect to the outgoing beam's frequency by the axion's average momentum in the beam direction.
In this article, we provide a summary of the results found in [21,22]. In Section 2, we calculate, in a simple manner, the total power of the echo wave, while in Section 3, we specify how this power is spatially distributed. In Section 4, we discuss the sensitivity prospects of the method based on realistic radio astronomy equipment, and finally, in Section 5, we conclude.

Power in the Echo Wave
The first step of our discussion is to compute the power stored in the echo wave. Let A (0) ( x, t) be the outgoing beam vector potential and a( x, t) the ambient dark matter axion field. The interaction of the outgoing beam with the axion field results in a perturbative correction A (1) ( x, t) to the zero-order A (0) ( x, t). In the Coulomb gauge, this correction satisfies the inhomogeneous wave equation: We also neglected terms containing gradients of a. This is justified by the assumption of a nonrelativistic axion background. From now on, we will always ignore this kind of contribution.
For the time being, let us just consider one axion's momentum mode with energy density ρ( p ). The axion field can be written as: where E p = m 2 + p 2 . As E p = m + O(p 2 ), we use E p ≈ m throughout this article. Let us expand the source field A (0) , as well as the correction A (1) in Fourier modes as: Here, we assumed an incident field polarized in the directionê. A (1) is polarized necessarily in the directionk ×ê, except for small corrections of the order of p. We also defined ω(k) ≡ | k|. Keeping only terms relevant to the stimulated axion decay, i.e., the photon momentum modes k and q = p − k, Equation (1) becomes: By looking at the source term of the above equation, we see that when m > ω(k), the solution describes a wave that travels backwards with respect to the incident one. We call it the echo wave. The echo wave is excited when ω(q) is equal or very close to m − ω(k), i.e., when ω(k) ≈ m/2. Now, we look for a resonant solution using the ansatz A (1) (t, q ) = A(t, q )e −iω(q) t where A(t, q ) varies slowly in time with respect to A (1) (t, q ). Neglecting the second derivatives of A(t, q ) in Equation (5), we obtain the solution: where ( k, p ) = ω(k) + ω(| p − k|) − m. In the limit t 1, we can substitute sin( t/2)/ → πδ( ) in Equation (6). Integrating over all axion's momentum modes and assuming that the incident beam propagates mostly in a preferred direction (the direction of k is fixed), the echo power is found to be: where p is the axion's momentum projected in the direction of k. To obtain this result, we assumed that the energy is provided by a source emitting a constant power P 0 from t = 0 until t = t off . Let us perform the integral (7) in two scenarios: when the bandwidth of the beam δω is much bigger than the momentum bandwidth of the axion δp and in the opposite case. For δω δp, the axion distribution can be considered as if it were condensed at a single valuep . We then approximate ∂ρ/∂p = ρ δ(p −p ). We straightforwardly find: where: For δω δp, we have ∂P 0 /∂ω = P 0 δ(ω −ω). We find: where: For the DFSZ model of the QCD axion, an input power of 1 kW over a bandwidth of 1 kHz working for one hour provides a total echo power of P ∼ 10 −19 W m 10 −4 eV 2 , which is in principle not too difficult to detect with current radio astronomy technology. Unfortunately, as discussed in [21], due to the nontrivial axion's velocity distribution, this power spreads over a surface that eventually exceeds the detector's size. This effect is easy to visualize if we consider what happens for a single axion decay. If the decaying axion moves with a velocity not perfectly aligned with the incident photon (outgoing beam), the echo photon is released in a direction also different from the incident photon's trajectory. In fact, if the echo photon has energy Ω, the conservation of momentum implies Ω sin(χ) = mv ⊥ , where χ is the angle formed by the trajectories of both photons and v ⊥ is the component of the axion's velocity perpendicular to the incident photon momentum. From energy conservation, Ω is equal to m/2 plus small corrections; therefore, sin(χ) 2v ⊥ . As the axion's velocities are of the order of 10 −3 or smaller, we can write χ 2v ⊥ . It follows that, for instance assuming the isothermal sphere model [23] for the galactic halo, which has a velocity dispersion of 270 km/s, the echo spreads over approximately 10 6 km after 1 h. Even in the caustic ring model of the galactic halo [24], where the minimal axion transverse velocity we can achieve is about 5 km/s, the spread is at least 10 4 km after the same amount of time. The echo power signal is thus strongly dependent on the model for the axion phase-space distribution, as well as the size of the detection apparatus. The following section is devoted to analyzing these effects by computing the echo intensity as a function of the spacetime coordinates.

Echo Intensity
The total power in the echo wave, computed in the previous section, does not carry any information about how it is spatially distributed. As explained above, this spatial dependence is closely related to the axion's velocity distribution, and it especially depends strongly on the transverse components of the axion's velocity. This section is dedicated to characterizing these effects by computing the echo intensity as a function of the spacetime coordinates.

Echo of a Dish Antenna Beam
To develop a sense of how the transverse axion's velocities affect the signal, we perform a simple and useful computation of the echo intensity for some particular models of the axion's velocity distribution. Moreover, as an incident wave, we take the beam emitted by a parabolic dish antenna. Although the parabolic antenna case is the most realistic setup we can address in this work, our findings are limited by scenarios where analytical approximations are possible. These approximations will be explained throughout the text.
Given the axion field a(t, x) and outgoing beam magnetic field B (0) (t, x), the echo vector potential field A (1) (t, x) is determined by: We write the axion dark matter field, in a large volume V, as an expansion in momentum modes as: where p is the axion's momentum in the experiment's rest frame and φ p random phases. These phases are not known and are usually modeled as uniformly distributed random numbers [25][26][27]. In this article, we limit ourselves to considering the ensemble average: We leave a discussion of the statistics of the signal to future work. We normalize f a ( p ) as: We considered a parabolic antenna of radius R fed by a plane wave with electric field amplitude E 0 (ω), ω being the wave's frequency. The electric field emitted by the antenna is known analytically in the far-field zone, i.e., for distances of order ωR 2 or larger from the emission spot. Assuming that the center of the antenna is located at x = 0 and that the central component of the beam propagates alongẑ, the antenna's electric field can be written in spherical coordinates as: where J 1 (x) is the Bessel function of the first kind of order one. For Equation (16) to be a good approximation, we need most of the echo to be produced in the far-field zone.
Therefore, for an outgoing beam turned on at t = 0, our approximation is valid for t ωR 2 . We perform the integral (12) in the limit ωR 1, which is valid for any experimental setup we consider in this article. If the outgoing beam is emitted continuously from t = 0, the maximum echo spectral intensity at the emission spot is: where v ⊥ = p ⊥ /m is the axion's transverse velocity and Ω * = (m − p z )/2. The symbol means averaging over the axion's velocity distribution | f a ( v)| 2 .
We can see from the above result that the echo signal depends strongly on the local axion's velocity distribution. Indeed, it scales as the inverse of the momentum dispersion in the forward direction and the inverse of the average transverse velocity component. cannot be reduced to values smaller than 5 km/s. In this analysis, we saw how the axion's velocity distribution influences the echo signal.
Equation (17) shows the effects explicitly. For t < R v −1 ⊥ , the echo has not yet spread beyond the emission region, leading to the linear growth of the intensity, in agreement with the results of Section 2. After t = R v −1 ⊥ , the intensity saturates suddenly to a constant value, due to the transverse velocity effects. This abrupt change in the intensity behavior is a consequence of the approximations employed, the assumption that most of the echo is coming from the far-field zone. If we had also taken into account the near behavior of the outgoing beam, the transition to the saturated regime would have been smooth. Another limitation of our result is that Equation (17) gives us an approximated value for locations nearby the emission spot. A complete characterization of the outgoing beam would allow us to know the local values of the intensity in more detail. Unfortunately, there is no simple analytical expression for the antenna's emitted field in the near zone. Nevertheless, this issue is addressed in Section 3.2 by constructing a simple and well-behaved model of the outgoing beam in the near zone.

Paraxial Gaussian Beam
To gain further insights into the intensity as a function of the position and into the transition to the saturated regime, we take a simplified model for the outgoing beam. This beam features a simple shape in the near-field zone and matches the antenna beam's far-field zone behavior correctly. Given this beam, we use the paraxial approximation to obtain the local values of the echo intensity as well as its time evolution. We write the beam magnetic field as a transverse wave propagating in the z direction with small transverse corrections: For the spatial shape of the beam η( x, ω), we take the paraxial Gaussian beam used in laser physics: where (r, φ, z) are cylindrical coordinates. Here, w(z) is the effective radius of the beam and R(z) is the radius of curvature of the wavefronts. Their mathematical expressions are: where R = w(0) and z R = ωR 2 . The beam is turned on at t = 0. At a time t > 0, it extends from z = 0 to z = t. After the beam is turned off at t = t off , it extends from z = t − t off to z = t.
We describe the local dark matter axion field as a superposition of plane waves, in the same way as in Section 3.1 (see Equations (13)-(15)). In addition, we also assumed that f a ( p ) can be factorized in its forward and transverse parts. In the limit δp z t 1, after averaging over random phases, the spectral echo intensity at resonance can be written as: where the quantity T (t, x ⊥ ) has units of time and contains all the three-dimensional effects. In a one-dimensional treatment, we have simply T (t, x ⊥ ) = min(t, t off ) (see [22]), while in three dimensions, the behavior is more complicated. The lateral spread leads to the saturation of the intensity on a characteristic time scale, which is given below. On the other hand, we find null effects coming from the outgoing beam divergence. We first considered the case of a negligible axion velocity dispersion. In this case, all axions move with approximately the same velocity. We denote the axion's velocity in the plane perpendicular to the beam direction of propagation as v p . We analyzed the small velocity dispersion case in two scenarios: (a) t off < R/v p and (b) t off > R/v p .
In Scenario (a), for t < t off , T (t, x ⊥ ) grows linearly in time, featuring a transverse Gaussian profile with effective radius R. For t > t off , i.e., after the outgoing beam is turned off, the Gaussian keeps the maximum value reached at t = t off and moves rigidly with velocity v p . In Scenario (b), for t < R/v p , the intensity grows linearly in time with a Gaussian profile, just as in Scenario (a) before t off . Later, for R/v p < t < t off , the echo spreads along the direction of v p at a speed v p , forming a "sausage"-shaped profile (see Figure 1). The energy that in Scenario (a) accumulates at the emission spot now spreads laterally, making the echo intensity saturate at any coordinate x ⊥ when t > r cos(φ) v p , φ being the angle between x ⊥ and v p . The maximum saturated value is reached at locations where φ = 0 and R < r cos(φ) < v p t. The corresponding value for T within this region is: Finally, for t > t off , the "sausage" moves rigidly at speed v p , keeping a fixed length v p t off and a width 2R. Figure 1. Echo intensity form factor F(t, x ⊥ ) at resonance ω = ω * in the case of saturation due to the transverse average axion's velocity. The beam is turned on at t = 0 and turned off at mt off = 1.4 × 10 9 . The radius of the emitter is set to mR = 4000. The axion's velocity components are those of the caustic ring model, with minimum transverse velocity v x = 5 km/s and v z = 290 km/s. The velocity dispersion is 70 km/s in the z direction. The saturation time is mt sat ≈ mt a = 4.3 × 10 8 .
Next, we considered the case in which the axion transverse velocity dispersion δv ⊥ dominates over the transverse average velocity. Again, we divided our discussion into scenarios: (a) t off < R/δv ⊥ and (b) t off > R/δv ⊥ . In Scenario (a), the intensity grows at first linearly in time, keeping a Gaussian shape, as in the small velocity dispersion case. After t off , the maximum intensity stays constant at the value obtained at t = t off , and the Gaussian shape stays unchanged. Finally, for times larger than R/δv ⊥ , the value of the intensity decreases as (δv ⊥ t) −2 and the radial extension becomes δv ⊥ t, i.e., the energy spreads radially at velocity δv ⊥ . In Scenario (b), we have again the Gaussian profile with a linearly growing maximum value for t < R/δv ⊥ . For R/δv ⊥ < t < t off , the energy spreads radially, making the intensity saturate for all the positions satisfying r < δv ⊥ t (see Figure 2). The maximum value of T in this saturated regime is found at r = 0 and is: After t off , the intensity falls off as (δv ⊥ t) −2 .
The characteristic times t a and t b are not only the saturated values of T in the respective scenarios, but also the timescales over which the saturation is achieved. We can then define the saturation time: Figure 2. Echo intensity form factor F(t, x ⊥ ) at resonance ω = ω * in the case of saturation due to velocity dispersion effects. The beam is turned on at t = 0 and turned off at mt off = 7 × 10 5 . The radius of the emitter is set to mR = 100. The axion's velocity components are those of the isothermal model, with v x = 10 km/s and v z = 230 km/s, while the velocity dispersion is 270/ √ 3 km/s in each direction. The saturation time is mt sat ≈ mt b = 2.1 × 10 5 .
To confirm the results of the discussion above, we solve Equation (1) numerically in the paraxial limit, with the model beam Equations (18) and (19). We plot the intensity in the z = 0 plane for the caustic ring model and the isothermal sphere, corresponding to the small dispersion and large dispersion cases discussed above, respectively.
We write the echo intensity as: By comparing with Equation (20), we see that, in the limit δp z t 1, the echo intensity factor is F(t, x ⊥ ) = T (t, x ⊥ )/t sat . This means that, when the echo saturation is predominantly due to one of the transverse average axion's velocity or velocity dispersion, the maximum value of F tends to one for times larger than t sat , provided t off > t sat . Figure 1 shows the echo intensity factor at resonance, i.e., forω = (m + p z )/2, for parameters compatible with the caustic ring model. We can clearly see the echo's "hot spot" spreading in the direction of the axion average velocity until the beam is turned off at mt off = 1.4 × 10 9 . Thereafter, the hot spot stops spreading and starts traveling rigidly. Figure 2 shows the echo intensity factor at resonance in a case of saturation due to velocity dispersion effects for parameters compatible with the isothermal sphere model. The hot spot spreads out in all directions, while at the same time, the maximum intensity grows until the beam is turned off. After that, the maximum intensity quickly drops to zero. Finally, we would like to comment on the echo profile in the frequency space. Figure 3 shows F(t, 0 ) as a function of the beam frequency ω at different times. The intensity factor was multiplied by the spectral intensity inherited from the beam dĨ 0 /dω| ω=Ω , whose maximum value was normalized to one. The beam spectral intensity, normalized to its maximum value, is shown in black. The echo peaks at a frequencyΩ = (m − p z )/2, while the beam is centered atω = (m + p z )/2. The echo bandwidth is given by δωδp z /(2∆). It is then possible for the echo frequency band to be completely separated from the beam frequency band if: Such a configuration may be advantageous to separate the signal from the noise if a part of the latter comes from the outgoing beam itself or from backscattering in the atmosphere.

Echo Collected Power and Sensitivity
The power collected over a surface S laying in the z = 0 plane can be computed as: where S 0 is the effective cross-sectional area of the beam at z = 0, P 0 = I 0 S 0 is the power of the beam, and: We estimated the sensitivity of the echo method assuming that the echo power is collected by a dish antenna with radius R c , located in the plane z = 0 and concentric with respect to the emitter. In this configuration, and for R c > R, we have 1 The signal-to-noise ratio is given by Dick's radiometer equation: where P c is the collected power, T n the noise temperature, B the bandwidth of the signal given roughly by B ≈ min(δω, δp z /2)/(2π), and t m the measurement time.
We determined t m as the time at which the power collected by the receiver drops off. The echo spreads all over the receiver in a time t ⊥ . If at time t ⊥ , the beam is still on, the collected power will drop at t = t off . On the other hand, if the beam is turned off before t ⊥ , the power collected will drop off then at t = t ⊥ . In summary, the measurement time is given roughly by the maximum between t off and t ⊥ .
The solid blue lines of Figure 4 mark the parameter space where the echo method is sensitive, assuming a fixed amount of energy E = 10 MW y is spent to cover an octave in axion mass. We assumed this energy is delivered with the largest power compatible with the constraint t off ≥ (2δω) −1/2 . This choice maximizes the sensitivity for a given bandwidth. However, in some parts of the parameter space, too large a power would be required. We also assumed s/n = 5, T n = 20 K, R = 50 m, R c = 100 m and δω = δp z /2. For both the isothermal model (right panel) and the caustic ring model (left panel), we assumed minimal average axion transverse velocity: v p = 0 and v p = 5 km/s, respectively. For the isothermal sphere model, we used a velocity dispersion 220 km/s. Figure 4 shows a range of axion masses compatible with the working frequency of radio telescopes.
The dashed blue lines of Figure 4 correspond to the sensitivities using a power source of P 0 = 10 MW working for 1 y. The other parameters are the same as for the solid blue lines. Despite the fact that the amount of energy spent is the same as for the solid lines, the sensitivity is worse. It is likely possible, though, to approach the solid lines' sensitivity by setting up the experiment in a clever way. . Expected sensitivity of the echo method in the space of parameters (m, g) for the caustic ring (left) and isothermal model (right) of the galactic halo. These plots assume that the method consumes an energy of 10 MW y per a factor of two in the axion mass range. The difference between the solid and dashed lines is explained in the main text. The green regions are current bounds from axion dark matter haloscopes [31][32][33][34][35][36][37][38][39][40][41][42][43][44], and the grey region corresponds to bounds from the CAST experiment [45].

Conclusions
In this manuscript, we presented a description of the echo method for axion dark matter detection proposed in [21]. The echo intensity grows in time until it saturates after a characteristic time scale that depends on the transverse axion's velocity distribution. We classified this velocity distribution effect into two cases: small and large dispersion, depending on whether the transverse velocity dispersion is smaller or larger than the transverse average velocity. We also showed that it is possible to achieve a separation in frequency between the beam and echo bandwidths, thanks to the average velocity of the axion flow along the beam's direction relative to the lab frame. This may help reduce the noise from atmospheric backscattering and beam leakages.
We provided sensitivity estimates assuming a 100 m receiving dish and a power of 10 MW. These experimental parameters are attainable with currently available technology. For instance, the receiver could be a radio telescope such as the Green Bank and Effelsberg of FAST, while the power source could be a high-power klystron, as used for radar transmitters or particle accelerators.
Our sensitivity estimates agree with those of [21] when considering velocity distribution effects and also with the fact that the outgoing beam divergence does not play any role in the signal. The echo method appears an attractive approach to axion dark matter detection because it needs relatively old technology and because it is applicable over a wide range of axion masses.