Limits on High-Frequency Gravitational Waves in Planetary Magnetospheres

High-frequency gravitational waves (HFGWs) carry a wealth of information on the early Universe with a tiny comoving horizon and astronomical objects of small scale but with dense energy. We demonstrate that the nearby planets, such as Earth and Jupiter, can be utilized as a laboratory for detecting the HFGWs. These GWs are then expected to convert to signal photons in the planetary magnetosphere, across the frequency band of astronomical observation. As a proof of concept, we present the first limits from the existing low-Earth-orbit satellite for specific frequency bands and project the sensitivities for the future more-dedicated detections. The first limits from Juno, the latest mission orbiting Jupiter, are also presented. Attributed to the long path of effective GW-photon conversion and the wide angular distribution of signal flux, we find that these limits are highly encouraging, for a broad frequency range including a large portion unexplored before.


INTRODUCTION
The successful detection of gravitational waves (GWs) by the Laser Interferometer Gravitational-Wave Observatory (LIGO) opens up a window to observe the Universe otherwise inaccessible [1].This has motivated a series of ongoing and to-be-launched projects to detect the GWs with a frequency ranging from ∼ 10 3 Hz to orders of magnitude below that.Yet, the GWs with a frequency above that could have also been produced in the early cosmological events such as preheating and high-temperature phase transition and the violent astronomical activities of small-scale objects, e.g., merging of primordial black holes and intercommutation of cosmic strings.Thus, detecting high-frequency GWs (HFGWs) is of high scientific value (for a review, see, e.g., [2]).
However, the detection of HFGWs has been significantly less explored than that of low-frequency GWs [3][4][5].Due to the shorter wavelength of HFGWs, this task is more challenging.One traditional wisdom is to employ the inverse Gertsenshtein effect [6][7][8][9][10][11], where the HFGWs are expected to convert to signal photons in an astronomical [12][13][14][15][16] or artificial  magnetic field.To compensate for the weakness of gravitational coupling, the magnetic field needs to be either strong or distribute broadly in space.Nonetheless, the existing proposals are subject to a variety of weakness, such as relatively short path for high-efficient conversion (e.g., neutron star (NS) [11]), large uncertainty of cosmic magnetic field strength [15], highly specific detection frequency band (for a recent effort to address this, see [38]), and narrow angular distribution of signal flux (especially for some laboratory experiments).
Alternatively, in this Letter, we propose to detect the HFGWs using the nearby planets such as Earth and Jupiter as a laboratory, where the GW-photon conversion is expected to occur in their planetary magnetosphere.Due to its relatively big size, the path for the effective conversion in such a laboratory is typically long.Par-ticularly, as to be shown, such an effective conversion can be achieved across the full electromagnetic (EM) frequency band of astronomical observation, ranging from radio waves to PeV photons.Moreover, as the detectors are positioned in the planetary magnetosphere, the stochastic signals can be detected in a wide range of directions.Combining these features creates a new operation space for detecting the HFGWs (for applying the planets to detect dark matter, see, e.g., [39][40][41][42][43][44][45]).
As a proof of concept, we consider the satellite-based detectors at low Earth orbit (LEO), with a bird view to the dark side of Earth.Both diffuse sky background and sunshine are expected to be occulted by the Earth then.We will present the first limits for some specific frequency bands and project the sensitivities for the future more-dedicated detection.It is important to note that the variety of detector designs, such as terrestrial versus satellite-based, bird view versus bottom view, etc., can have significant impacts on the sensitivities.Therefore, this study should be considered as a starting point for more systematic exploration of the opportunities presented by such a laboratory, rather than a full demonstration of the sensitivity potential of this strategy.

GW-PHOTON CONVERSION PROBABILITY
With a WKB approximation, the inverse Gertsenshtein effect is characterized by a mixing matrix [11,46] (see Supplementary Materials (SM) at Sec.A for details [47]) Here ∆ M = 1 2 κB t encodes the GW-photon mixing, with κ = (16πG) 1/2 and B t being the component of external magnetic field B transverse to the GW traveling direction.∆ γ ≈ ∆ vac + ∆ pla is the effective photon mass.∆ vac = 7αω/(90π)(B t /B c ) 2 denotes the QED vacuum effect, with α the fine-structure constant, ω the angular fre- e /e.∆ pla = −m 2 pla /(2ω) represents the plasma-mass contribution with m 2 pla = 4παn c /m c , where n c and m c are the number density and invariant mass of charged plasma particles.By diagonalizing this mixing matrix, one can obtain the GW-photon conversion probability in a homogeneous magnetic field [15,48] Here Θ = 1 2 arcsin(∆ M l osc ) and l osc = 2/(4∆ 2 M + ∆ 2 γ ) 1/2 are the GW-photon mixing angle and oscillation length, respectively.L is the travel distance of GWs in the magnetic field.For a general path from ℓ 0 to ℓ 1 , the conversion probability can be evaluated as [11,49]: To have a taste, we first evaluate P 0 , the conversion probability for a radial path from zero altitude to infinity over the planet equator.We consider the NSs also for reference.The magnetic fields of both can be modeled as a magnetic dipole, with the magnetic axis aligned with the rotation axis.Then they are transverse to the radial path, with B t = −B 0 (r 0 /r) 3 .Here r 0 is radius and B 0 is surface magnetic field strength.We show P 0 as a function of frequency for the Earth, Jupiter and two benchmark NSs in Fig. 1.Because of their difference in plasma density profile and external magnetic field strength, the curves demonstrate quite different features.FIG.1: P0 as a function of f = ω/2π.We show (r0, B0) for the planets and (r0, B0, 2π/ΩNS) for the NSs, where r0, B0 and ΩNS (spinning velocity) are in unit of km, Gauss and second −1 , respectively.The green region denotes the EM observation frequency band in astronomy, ranging from radio waves to PeV photons.
For the planets, the plasma density is described by a barometric formula, which yields an exponentially suppresses plasma mass as r increases.If ω is not too small, we can set ∆ γ ≈ ∆ vac and then obtain where ω tra is determined by r 0 ≈ l osc (r 0 ), corresponding to the transition point of planet curves in Fig. 1.Consisting with Fig. 1 which is full-calculation-based, P 0 is approximately a constant for ω ≲ ω tra and drops as ∝ ω −4/5 for ω ≳ ω tra .These relations can be qualitatively explained with Eq. ( 2).In the homogeneous magnetic field, an optimal probability P ≈ ∆ 2 M L 2 can be achieved for sinc(L/l osc ) → 1 or l osc ≳ L, a case dubbed as "coherent conversion" [50,51], for given L. When l osc becomes smaller than L, P is suppressed by l 2 osc /L 2 .Extending this criterion to the planets, we have coherent conversion near their surface for L ∼ r 0 ≲ l osc (r 0 ) or ω ≲ ω tra , where P 0 ∝ ∆ 2 M (r 0 )r 2 0 ∝ B 2 0 r 2 0 .P 0 for the Jupiter in this region is then ∼ 10 4 times larger than that of the Earth, as its B 0 and r 0 are both ten times larger.For r 0 > l osc (r 0 ), the coherent conversion is suppressed near the planet surface.It can only occur for r * ≈ l osc (r * ) > r 0 , with a reduced rate 0 ω −4/5 .As shown in Fig. 1, the frequency band for EM astronomical observations falls entirely into the range of near-surface coherent conversion for both Earth and Jupiter.
In comparison, due to the large strength of their external magnetic field, the NSs tend to have a suppressed l osc .The near-surface coherent conversion thus becomes difficult to achieve.For the two benchmark NSs in Fig. 1, it takes place between 10 2 − 10 11 GHz for the NS2, and hardly occurs for the NS1 [52].In the high-frequency limit, the vacuum effect dominates.P 0 is given by the 2nd formula in Eq. (4).At the low-frequency end, the plasma effect becomes dominant, yielding a sine-wiggling P 0 .With the Goldreich-Julian model for the NS plasma density [53], where decays more slowly than ∆ vac does as r increases.The coherent conversion then takes place at a larger r * , compared to the high-frequency case, yielding a more suppressed P 0 as ω decreases [54].Now let us consider a satellite-based detector positioned at altitude H and latitude λ O ′ , and define its instant spherical coordinate system with ẑ′ pointing to the planet's center (see SM at Sec.A for details).Then we are able to evaluate P (Ω ′ ) of HFGWs traveling to this detector in all directions, where Ω ′ = {θ ′ , ϕ ′ }.Fig. 2 displays the θ ′ dependence of P (Ω ′ ), with P (θ ′ ) = (2πP 0 ) −1 2π 0 P (Ω ′ )dϕ ′ .While this figure is drawn with f = 10 8 GHz, the features that it demonstrates are almost unchanged for the GW frequency range of interest.Here we neglect the plasma effect which is expected to be small.We split the photon incoming directions into planet-cone (PC) (θ ′ < θ c = arcsin[r 0 /(r 0 + H)]) and outer-space (OS) (θ ′ > θ c ) regions, based on whether the line of sight intersects with planet surface.For all sets of (H/r 0 , λ O ′ ), P (θ ′ ) peaks sharply at the PC edge and drops quickly away from it.For λ O ′ = π/2, it drops all the way to zero at θ ′ = 0, as B t vanishes in this direction for a detector right above the planet poles.As H/r 0 increases, with a cost of smaller PC, P (θ ′ ) tends to have a bigger value inside the PC due to a longer GW-photon conversion path.

SENSITIVITY ANALYSIS
The stochastic GWs are typically isotropic and stationary [55].For a detector in the planetary magnetosphere, we have the GW-converted photon flux (see, e.g., [56]) where ∆Ω denotes the detector field of view (FOV), and h c is the GW characteristic strain.Note, the angular distribution of GW-converted photons is determined by P (Ω ′ ) and only those falling into the detector FOV can contribute to Φ γ .The signal and background counts for a narrow frequency band ∆ω, a short observation time ∆t and an effective detector area A are then given by Here ⟨P ⟩ det = ∆Ω P (Ω ′ )dΩ ′ /∆Ω, essentially determined by the detector position and pointing direction {θ ′ det , ϕ ′ det }, denotes the average GW-photon conversion probability over its FOV.If the detector points to the planet center (θ ′ det = 0), ⟨P ⟩ det is ∝ ∆Ω P (θ ′ ) sin θ ′ dθ ′ / ∆Ω sin θ ′ dθ ′ and increases with ∆Ω inside the PC.The 95% upper limit on h c can be derived from s/ √ b ≈ 1.64 in the large-background limit: In this Letter, we mainly consider the LEO satellite which has a bird view to the dark side of Earth.With the detector FOV restricted to inside the PC, the dominant backgrounds vary from atmospheric thermal emissions in the infrared (IR) band [57] to cosmic photon albedo (see, e.g., [58]) throughout the optical -γ-ray region.We simplify the analysis by assuming a uniform background flux ϕ b .We will briefly discuss the Jupiter case also.Given the strong angular dependence of P (Ω ′ ) (see Fig. 2), an efficient detection requires proper design for satellite orbit and full optimization of detector performance.Below let us consider some specific examples.
We first consider the detection of HFGWs in the Earth's magnetosphere.Take the Suzaku mission [64] as an example.The Suzaku has a low inclination orbit (see Tab. I).It revolves at low latitude and has a performance relatively insensitive to season alternation.Eq. ( 7) thus applies approximately for the entire observation period T dark .Also, as the Suzaku FOV is small, ⟨P ⟩ det does not vary much for θ ′ det ≪ θ c .For θ ′ det = 0 we have ⟨P ⟩ det ≈ 5 × 10 −35 .As for the background flux in the dark side, it has been measured by the Suzaku to be ϕ b ≈ 6.3 ×10 −8 cm −2 s −1 keV −1 arcmin −2 for ω ∼ 0.5 − 10 keV [64].For one-year operation, we have This limit can be scaled to other missions with similar properties, following Eq.( 7).Notably, the satellite at a high inclination orbit has quite different properties.It scans over the high latitude region also, where the conversion probability varies more inside the PC and becomes bigger for a region extending from the PC edge to its inside, compared to the lowlatitude case (see Fig. 2).The sensitivity thus could be optimized by taking a detector with either a large FOV covering the whole PC or a small FOV but pointing to the region near the PC edge.Moreover, the observation of such a satellite in the dark side of Earth is sensitive to the season alternation.So we need to generalize Eq. ( 7) by binning the detector FOV and observation time in this case.The limit shall be obtained from the combined statistics (see SM at Sec. B for details).
Next, we consider the HFGWs conversion in the Jovian magnetosphere.Take Juno, the latest mission orbiting the giant, as an example.The Juno is settled at a highly elliptical polar orbit with a Perijove ∼ 5000 km [72] and a high inclination angle (∼ 90 TABLE I: Benchmark scenarios for sensitivity study.For the conservative case, we take the Suzaku-like low inclination orbit [59] with one year operation.The detector FOV and effective area are assumed to be the same as those of the existing missions including Nimbus [60] (IR), Hubble [61, 62] (ultraviolet (UV)-Optical), Voyager [63] (extreme UV (EUV)), Suzaku [64] (X-ray) and Fermi-LAT [65] (γ−ray).For the optimistic case, we take the Safir-2-like high inclination orbit [66] with ten year observation.A FOV covering the whole PC is considered.The effective area is set in a way such that the corresponding etendue for the IR, UV-optical, EUV and X-ray bands [67] is consistent with that of future missions [68-71].
Jupiter's cloud.So ∆t is ∼ 10 5 s for the ∼ 35 rounds of its prime mission.As a conservative estimate, we consider the Juno observation of aurora emissions by the Jovian Infrared Auroral Mapper (JIRAM) [73] and UV emissions by the Ultraviolet Spectrograph (UVS) [74], which set respectively ϕ b ∼ 4.2 × 10 4 cm −2 s −1 eV −1 arcmin −2 for ω ∼ 0.35 eV [75] and ∼ 34 cm −2 s −1 eV −1 arcmin −2 for ω ∼ 8 eV [76].We take θ ′ det ∼ 0.8θ c for ⟨P ⟩ det and assume ⟨P ⟩ det to vary little with time.Then, with ∆Ω ∼ 10 −3 sr and 3 × 10 −4 sr, A ∼ 2 cm 2 and 6 cm 2 and ∆ω ∼ 1.4 × 10 13 Hz and 0.9 × 10 15 Hz, we have for the JIRAM and the UVS, respectively.As a reference, let us consider the conversion of HFGWs in the NS magnetosphere.The NSs are remote.So we have ⟨P ⟩ det ∼ P 0 and ∆Ω ∼ πr 2 0 /d 2 , where d is their distance to the Earth.Take the Magnificent Seven (M7) X-ray dim isolated NSs as an example, which have d ∼ O(100) pc and B 0 ∼ 10 13 Gauss, and hence P 0 ∼ 10 −18 and ∆Ω ∼ 10 −30 sr.The PN and MOS of the XMM-Newton telescope have measured the flux of, e.g., J0420 (d ≈ 345 pc), to be ϕ b ∼ 10 17 − 10 19 cm −2 s −1 keV −1 arcmin −2 for ω ∼ 0.5 − 1 keV, where A ∼ 1500 cm 2 and ∆t ∼ 10 5 s.Such an intensity consists well with the background model of thermal surface emissions [77].By applying Eq. ( 7), we find Due to the extreme smallness of ∆Ω, this limit is about five orders of magnitude worse than that in Eq. ( 8).This highlights in part the merit of nearby planets in performing the task of detecting HFGWs.It is also informative to compare these limits with the ones recast from the existing laboratory experiments for axion detection [24].In these experiments, the magnetic field is typically confined inside a long straight pipe, with a small cross-sectional area.The GW-photon conversion is suppressed outside the small opening angle of long pipe.The signal flux in Eq. ( 5) is thus strongly limited by the detector geometry (see SM at Sec. C for details).Taking the CERN Axion Solar Telescope (CAST) experiment [78] as an example, we have the recast limit PROJECTED SENSITIVITIES To demonstrate the potential of detecting stochastic HFGWs with the LEO satellites, we define two benchmark scenarios, dubbed "conservative" and "optimistic", in Tab.I.The choice of orbit in these two scenarios echoes the discussions above on its impacts on the detection efficiency.Fig. 3 shows the projected 95% C.L. limits on the characteristic strain of HFGWs for a wide range of frequencies, including a large portion unexplored before.
In the γ-ray band, the limits are strengthened with a rate faster than ∼ f −1/2 , due to the reduction of albedo backgrounds [85].For the X-ray region, the background flux is taken from the Suzaku measurement.Its conservative limits thus represent the first limits from this mission.For the UV-optical band, we consider the cosmic background [79,80] with a frequency-dependent reflectance of atmosphere [81,82].The limits get nearly one order of magnitude stronger at the right edge of the UV-optical band, due to the ozone absorption of cosmic photons.Such a trend extends to the EUV band, where photodissociation and photoionization become important [92] and a reflectance of ∼ 10 −3 [84] is thus approximately taken for cosmic photons.As for the IR band, the backgrounds are modelled with a blackbody spectrum at 294K [57].These thermal radiations peak at f ∼ 10 13 Hz, and are quickly suppressed for f > 10 14 Hz.We also present the limits from the observations of Jupiter and NSs in some specific frequency bands.For the optimistic estimate of Jupiter limits, we take ∆Ω = 4 sr, A = 100 cm 2 and ∆t = 10 7 s.Notably, the energy density contribution of GWs with such a large stochastic magnitude at these frequencies is too signifi- FIG.3: Projected 95% C.L. upper limits on the characteristic strain hc of stochastic HFGWs, set by the LEO satellites (red).We estimate these limits with the Suzaku X-ray measurements [64] for band IV, and the IR [57], UV-optical [79][80][81][82], EUV [83,84] and γ-ray [85] backgrounds for bands I, II, III and V, respectively.The dashed and solid bottom lines for these bands correspond to the conservative and optimistic estimates defined in Tab.I.Moreover, we present the conservative and optimistic limits as orange bands, based on the Juno observation of aurora emissions [75] and UV spectra at Jupiter [76].We also present the limits set by the X-ray observations of M7 NSs [77] as cyan band.For completeness, we recast the limits of ALPS, OSQAR and CAST, i.e., three existing laboratory experiments, as gray bands.We denote the highly-uncertain radio limits of EDGES and ARCADE [2] as blue and yellow bands [86].The black dotted line reflects the Big-Bang-Nucleosynthesis (BBN) constraints on the radiation energy density in the early Universe [87], while the khaki dotted line represents the maximal strain values set by the critical energy density of the Universe today.As a reference, we present the HFGWs from the Hawking evaporation of highly spinning primordial black holes (PBHs) in the early Universe [88][89][90] as green solid lines, and from the binary mergers of PBHs after the BBN [91] as thin-solid and solid blue lines for two benchmark scenarios.While the former is subject to the BBN constraints, the latter is less restricted.
cant, compared to the BBN constraints and the cosmic critical value today.Nevertheless, the presented limits lay a foundation for further searching for astrophysical signals of, e.g., PBH mergers and advancing detector technologies capable of detecting the signals from the early Universe.

SUMMARY AND OUTLOOK
In this Letter, we have proposed to detect the stochastic HFGWs in the planetary magnetosphere.Due to the relatively long path for effective GW-photon conversion, the wide angular distribution of signal flux and a full coverage of the EM observation frequency band in astronomy, this strategy creates a new operation space.With the proof of concept presented, we can immediately see several important directions for next-step explorations.
Firstly, extend the detection of HFGWs from the PC to the OS.As indicated by Fig. 2, the GW-photon conversion probability right outside the PC can be much higher than that inside.Moreover, a detector oriented toward the OS may receive considerably less atmospheric thermal radiation.Secondly, extend the satellite-based detection to terrestrial observations, which may allow the sensitivities to cover the radio bands.Thirdly, extend the dedicated study to the Jupiter and even the Sun, given their stronger magnetic field, larger space for an effective GW-photon conversion, and the active and upcoming missions.We leave these explorations with refined analysis to a paper in preparation [54].

Supplementary Material
The Supplementary Materials contain additional calculations and explanations in support of the results presented in this Letter.Concretely, we discuss the GW-photon conversion in planetary magnetosphere in Sec.A, analyze the impacts of the orbit design for the LEO satellites on the HFGWs detection in Sec.B, and calculate the recast limits from the laboratory experiments in Sec. C.

A. GW-Photon Conversion in Planetary Magnetosphere
As gravity couples with EM dynamics via ∼ g µν T EM µν , the GWs can be converted to photons while propagating in an external magnetic field B. With a WKB approximation, which is valid here since the wavelength for the GW-converted photons expected to observe through astronomical telescopes is much shorter than the characteristic spatial scale (comparable to the planet size) at which the external magnetic field varies, the leadingorder equations of motion are given by [11,46] Here G × and G + denote the cross and plus polarization modes of GWs; A || and A ⊥ denote the polarization modes of photons parallel with and perpendicular to B T , the B component transverse to the GW propagation direction.In the 4 × 4 propagation matrix, denote an effective mass of A || and A ⊥ , respectively.Here encode the plasma and QED vacuum effects, and ∆ CM ∝ B 2 t [46] can induce Cotton-Mouton birefringence.m 2 pla = 4παn c /m c is plasma mass square, with n c and m c being the number density and invariant mass of charged plasma particles.Moreover, the ∆ R term, which is proportional to B − B T , couples A || and A ⊥ and can yield Faraday rotation.The GW-photon mixing, which is a core of this study, is governed by ∆ M = 1 2 κB t , where For the GW-photon conversion in planetary magnetosphere where ω of interest (see Fig. 3) is mostly larger than the cyclotron angular frequency ω c ≈ 10 7 (|B|/Gauss) Hz, ∆ CM and ∆ R are negligibly small compared to ∆ pla [46].Thus, the equations in Eq. ( 12) can be separated to two decoupled ones, for the polarization modes {A || , G × } and {A ⊥ , G + }, respectively.For simplicity, one can further neglect the small difference between ∆ vac,|| and ∆ vac,⊥ , and assume ∆ vac,|| = ∆ vac,⊥ = ∆ vac .The two propagating equations then share the same form and become where ∆ γ ≈ ∆ vac + ∆ pla .This gives exactly the mixing matrix in Eq. ( 1).For the formalism developed so far, we have neglected the plasma effect in the dielectric tensor of media, which has been carefully analyzed in the context of axionphoton conversion in the highly magnetized anisotropic plasma of NSs [100].The dielectric tensor is adjusted by terms proportional to m 2 pla /(γ 3 ω 2 ), where γ is Lorentz factor.These terms can mix the transverse and longitudinal modes of GW-converted photons, and significantly alter the flux transfer in the direction of incident GWs.To look into this effect in our case, we demonstrate in Fig. 4 the m 2 pla /ω 2 factor as a function of r/r 0 for different GW frequencies in the Earth's magnetosphere.γ ≈ 1 has been assumed for non-relativistic plasma in the planetary magnetosphere.For f ≳ 10 11 Hz, namely the GW frequencies considered for sensitivity analysis in Fig. 3, we find m 2 pla /ω 2 ≪ 1 for r/r 0 ≳ 1.02, which corresponds to an altitude H ≳ 130 km above the Earth's surface.m 2 pla /ω 2 > 1 only occurs for r → r 0 and relatively low GW frequencies, where however the relative contributions to the integrated GW-photon conversion probability for the LEO detectors tend to be suppressed (see, e.g., Fig. 5 below).Therefore, for the LEO detectors with H ≳ 600 km considered in our study, the plasma effect in the dielectric tensor of media can be safely neglected.Due to the smallness of gravitational coupling, we have ∆ M ≪ ∆ γ and hence l osc = 2/(4∆ 2 M + ∆ 2 γ ) 1/2 ≈ 2/∆ γ in general.This implies a weak mixing between the GWs and photons.Nonetheless, ∆ vac and ∆ pla have opposite signs, and a maximum mixing can be realized if the two cancel well such that ∆ γ ≪ ∆ M .In the planetary magnetosphere, such a cancellation might be possible only for a tiny range of the GW traveling path.Its contribution to the total conversion probability is typically negligible.Instead, the "coherent conversion", where the conversion probability is maximized for given magnetic field and GW path, tends to be realized for a longer part of the path and hence contribute more to the GW-photon conversion.
To have a taste, in the main text we first analyze the GW-photon conversion for a radial path of the planets from zero altitude and latitude to spatial infinity.Due to the difference in spinning velocity, the plasma density of planets is described by Barometric formula where H cor ≪ r 0 denotes the scale height and n c,0 denotes the surface plasma density, with H cor ≈ 8.4 ( 27) km and n c,0 ≈ 10 26 (10 29 ) fm −3 for the Earth (Jupiter) [93][94][95], while the plasma density of NSs is described by Goldreich-Julian model [53].As n Bar c is exponentially suppressed at r − r 0 ≫ H cor , the plasma mass contribution is mostly negligible along the path for the GW frequency of interest.Thus, we focus on the vacuum effect and assume ∆ γ = ∆ vac .P 0 then can be solved exactly as [96][97][98] where Γ(z) is gamma function and Γ(a, z) is incomplete gamma function.∆ M (r 0 ) and ∆ vac (r 0 ) are defined at the planet surface where r = r 0 .
In the low frequency limit, where ∆ vac (r 0 ) r 0 ≲ 1 and hence r 0 ≲ l osc (r 0 ), the coherent conversion can be achieved near the planet surface.Eq. 17 then gives the first formula in Eq. ( 4): Notably, the GW-photon conversion can be coherent also in the region well above the planet surface.log 10 [ r losc(r) ] is uniformly smaller than zero for r ≳ r 0 , as shown in Fig. 5. But, the contributions to P 0 from the high-altitude magnetosphere are relatively few since the magnetic field and (left axis, solid curves) and log 10 [ r losc(r) ] (right axis, dashed curves) as a function of r/r0 in the Earth's magnetosphere (ωtra ≈ 10 25 GHz), along the radial direction considered in Fig. 1.The black and red curves denote the low-frequency (ω = 10 GHz) and high-frequency (ω = 10 28 GHz) scenarios, respectively.
hence ∆ M (r) become weakened as the altitude increases.So, P 0 receives most of the contributions from the region near the Earth surface, yielding a sharp peak at r ∼ r 0 for the solid black curve in Fig. 5.
In the high frequency limit, where ∆ vac (r 0 ) r 0 ≳ 1 and hence r 0 ≳ l osc (r 0 ), the coherent conversion can not be achieved near the Earth surface.Eq. ( 17) then gives the second formula in Eq. ( 4): where the transition frequency is determined by ∆ vac (r 0 ) r 0 ≈ 1. Eq. ( 19) matches to Eq. ( 18) smoothly at ω = ω tra .In this limit however the coherent conversion can be realized in the magnetosphere at high altitude, where r = r * > r 0 .r * can be estimated by solving ∆ vac (r * )r * = 1.With ∆ M (r) = ∆ M (r 0 )(r 0 /r) 3 and ∆ vac (r) = ∆ vac (r 0 )(r 0 /r) 6 , we obtain r * = r 0 (∆ vac (r 0 )r 0 ) 1/5 .This implies that P 0 here receives more contributions at r ≳ r * .These features are shown in Fig. 5 for the Earth.Indeed, the contributions from the region near the Earth surface get suppressed in the high-frequency case, which leads to a relatively gentle peak for the solid red curve near r * /r 0 ≈ 10 −5 (ω/GHz) 1/5 .The features of the Earth and Jupitor curves in Fig. 1 can be well-explained with the formulae obtained by neglecting the plasma mass of photons in the planet magnetosphere, while these curves (together with Fig. 5) are drawn based on a full calculation.Actually, we can further demonstrate that such a negligence has limited impacts only on the accuracy of results, for the GW frequency range of interest.To see this point, we demonstrate in Fig. 6 the ratio of conversion probabilities P [H i , 800]/P vac [H i , 800] for a typical LEO detector, as a function of H i .From this figure, one can see that the P vac [H i , 800], calculated with the plasma term ∆ pla in ∆ γ being turned off, deviates from the full calculation P [H i , 800] by a factor ≲ 2 for f > 10 11 Hz.This deviation is maximized while H i → 0 and f → 10 11 Hz, as the plasma density gets exponentially enhanced with r → r 0 and meanwhile ∆ pla ∝ m 2 pla /ω becomes more important for small ω.Such a deviation however is too small to cause a significant impact on the sensitivity discussions in this proof-of-concept study.
For sensitivity demonstration, we consider a satellitebased detector.The instant coordinate system for this detector is displayed in Fig. 7.As the GW-photon conversion probability P (Ω ′ ) exhibits a stronger dependence on θ ′ than on ϕ ′ , it is convenient to show its θ ′ dependence in Fig. 2, by integrating out ϕ ′ .
As shown in Eq. ( 7), the limits on h c are sensitive to the satellite position, the detector pointing direction {θ ′ det , ϕ ′ det } and its FOV ∆Ω.Considering h c,95% ∝ ∆Ω 1/2 ⟨P ⟩ det , we show in Fig. 8 ∆Ω 1/2 ⟨P ⟩ det as a function of ∆Ω for different satellite latitudes and detector pointing directions θ ′ det .For each curve, the quantity increases monotonically with ∆Ω in general, and approximately scales as ∆Ω 1/2 for small FOV.As for the θ ′ dec dependence, the curves change more drastically for the high latitude satellite, in accordance with the angular dependence shown in Fig. 2. ∆Ω 1/2 ⟨P ⟩ det reaches its maximum when the satellite operates at λ O ′ = π/2, with the detector FOV covers the whole PC region.As the PC edge contributes more to the signal rate, a sensitivity of the same level could be also obtained for a detector with relatively small FOV but pointing to the PC edge.

B. Low Earth Orbits for Satellites
For a satellite-based detector, the signal and background photon fluxes are direction and time dependent.To evaluate the orbit quality accurately and hence optimize its design, one can bin the events w.r.t.direction and time while analyzing the sensitivity reach.Eq. ( 6) is then generalized to for the i-th angular and j-th time bin.Given the count of observed events n ij , the exclusion limit on h c is set by the test statistic [99] q which is Poisson-statistics-based and defined against the background + signal (h c ) hypothesis.The projected 95% C.L. upper limit on h c is then obtained by requesting q(h c,95% )| n=b = 2.7, where n ij is assumed to be the same as the predicted background count.In the Gaussian statistics, Eq. ( 22) is reduced to q(h c ) = i,j s 2 ij /(b ij + s ij ).It can be further simplified to q(h c ) = i,j s 2 ij /b ij in the limit of large background, where b ij ≫ s ij .This yields the criteria of s/ √ b = 1.64 for the one-bin analysis, which has been applied for the derivation of Eq. (7).To have a taste on the impacts of orbit design, let us consider two benchmark scenarios defined in Tab.II.In Scenario I, the orbit is Suzaku-like, characterized by a low inclination angle.As the detector is restricted to revolve in the low latitude region (|λ O ′ | ≲ 30 • ), the total time for its staying in the orbit of dark side, dubbed "dark orbit", varies little with the season alternation.This yields an overall fraction f dark ≈ 0.35 or equivalently effective observation time T dark ≈ 1.1 × 10 7 s for one-year revolution.In scenario II, the orbit is SAFIR 2-like, characterized by a high inclination angle.As the detector can reach the high-latitude region, the signal rate fluctuates more during each orbit period.The total time for its revolving in the dark orbit varies strongly with the season alternation.For one-year revolution, we have an overall fraction f dark ≈ 0.2 and effective observation time T dark ≈ 6.3 × 10 6 s.
To bin the event count in each benchmark scenario, we divide the dark orbit by a time interval ∆t = 120 s.Also, we consider the full PC region as the detector FOV and separate it into four concentric circular belts.We then evaluate q(h c )| n=b by taking a summation of s ij and b ij over the time-angular bins for 10-year revolution.Finally we display the test statistic q(h c )| n=b in Fig. 9 as a function of h c for these two scenarios (with A = 250 cm 2 ).The SAFIR 2-like orbit has a better sensitivity performance, although it spends less time in the dark orbit.This could be attributed to a larger GW-photon conversion probability off the PC center while the detector revolves in the high-latitude region (see Fig. 2).Also, as one merit of high-inclination orbit, the dark-orbit time is mainly contributed by the polar-night season, enabling a more efficient usage of the telescope resources.

C. Recast Limits from Laboratory Experiments
The limits obtained in the laboratory experiments for axion detection can be recast to constrain the HFGWs (for a review, see [2]).Since the magnetic field in these experiments is often confined inside a straight-long cylinder pipe, the conversion probability of stochastic HFGWs into photons exhibits a strong angular dependence.This will necessarily influence the estimate of the GW-converted photon flux in Eq. (5).Explicitly, let us consider the coherent conversion in an experimental facility shown in Fig. 10.For the central pixel of its detector, the angular integral in Eq. ( 5) reduces to an integration over θ, which can be put into two parts.The first one accounts for the signals re- arXiv:2305.01832v3 [hep-ph] 30 Mar 2024 quency and B c = m 2

2 FIG. 2 :
FIG.2: P (θ ′) as a function of θ ′ for a satellite-based detector, with f = 10 8 GHz.The color and line style denote the normalized altitude H/r0 and the latitudes λ O ′ respectively for the satellite.The vertical lines denote θc for the PC with different H/r0 values.

2 Ω G W = 1 B
fro m PB H ev ap or at ion S G W B s fr om P B H m er ge r h 0

FIG. 4 :
FIG.4: m 2 pla /ω 2 as a function of r/r0 in the Earth's magnetosphere, along the radial direction considered in Fig.1.

FIG. 6 :
FIG.6:P [H i ,800]  P vac [H i ,800] as a function of Hi, for a LEO detector over the equator.Here P [Hi, 800] and P vac [Hi, 800] are GW-photon conversion probabilities.They are calculated by integrating the radial-direction contributions from an altitude Hi < H to the detector location H = 800 km, with the plasma term ∆ pla in ∆γ being turned on and off respectively.

FIG. 7 : 2 FIG. 8 :
FIG. 7: Instant coordinate system for a satellite-based tector.{x, y, z} and {x ′ , y ′ , z ′ systems centered at the planet center O and the satellite location O ′ , respectively, with z axis and z ′ axis respectively oriented towards the planet north pole and center.λ O ′ denotes the latitude of the satellite.θc is the open angle of the PC that is tangential to the planet's surface.θ ′ and ϕ ′ are the polar and azimuthal angles the line of sight (solid green line) defined in the {x ′ , y ′ , z ′ } frame, respectively.

FIG. 10 :
FIG.10: Schematic for the GW-photon conversion in a straight-long cylinder pipe which has a length L0, a crosssectional A = πR 2 and an internal magnetic field B0 (with B0 = |B0|).The black dot denotes the central pixel of the detector at the right end of the pipe.The red lines denote the incoming GWs.These GWs travel in the pipe by a distance L, in a direction with 0 ≤ θ < π/2.θp ≈ R/L0 ≪ 1 denotes the opening angle of pipe.
•).For each orbit, the Juno spends a few hours around the Perijove in observing the