In the milliHertz frequency band stochastic gravitationalwave background can be composed of both astronomical and cosmological sources both can be anisotropic. Numerically depicting these anisotropies can be critical in revealing the underlying properties of their origins. For the first time we perform a theoretical analysis of the constraining ability of TianQin on multiple moments of the stochastic background. First we find that with a oneyear operation for a background with a signaltonoise ratio of 16 TianQin can recover the multiple moments up to l4. We also identified a unique feature of the stochastic background sky map which is the mirror symmetry along the fixed orbital plane of TianQin. Thirdly we explain the difference in anisotropy recovering ability between TianQin and LISA by employing the criteria of the singularity of the covariance matrix which is the condition number. Finally we find that since the different data channel combinations correspond to different singularities certain combinations might have an advantage in stochastic background mapmaking. We believe that the findings of this work can provide an important reference to future stochastic background analysis pipelines. It can also serve as a guideline for designing better gravitationalwave detectors aiming to decipher anisotropies in the stochastic background.
['General Relativity and Quantum Cosmology', 'Astrophysics - Cosmology and Nongalactic Astrophysics']
Mapping Anisotropies in the Stochastic GravitationalWave Background with TianQin
{'Mapping Anisotropies in the Stochastic Gravitational-Wave Background with TianQin': 'Zhi-Yuan Li, 1 Zheng-Cheng Liang, 1, ∗ En-Kun Li, 1 Jian-dong Zhang, 1 and Yi-Ming Hu 1, † 1 MOE Key Laboratory of TianQin Mission, TianQin Research Center for Gravitational Physics & School of Physics and Astronomy, Frontiers Science Center for TianQin, CNSA Research Center for Gravitational Waves, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082, China (Dated: September 18, 2024)  \nIn the milli-Hertz frequency band, stochastic gravitational-wave background can be composed of both astronomical and cosmological sources, both can be anisotropic. Numerically depicting these anisotropies can be critical in revealing the underlying properties of their origins. For the first time, we perform a theoretical analysis of the constraining ability of TianQin on multiple moments of the stochastic background. First, we find that with a one-year operation, for a background with a signal-to-noise ratio of 16, TianQin can recover the multiple moments up to l = 4. We also identified a unique feature of the stochastic background sky map, which is the mirror symmetry along the fixed orbital plane of TianQin. Thirdly, we explain the difference in anisotropy recovering ability between TianQin and LISA, by employing the criteria of the singularity of the covariance matrix (which is the condition number). Finally, we find that since the different data channel combinations correspond to different singularities, certain combinations might have an advantage in stochastic background map-making. We believe that the findings of this work can provide an important reference to future stochastic background analysis pipelines. It can also serve as a guideline for designing better gravitational-wave detectors aiming to decipher anisotropies in the stochastic background.', 'I. INTRODUCTION': "A stochastic gravitational-wave (GW) background (SGWB) consists of GWs that, when individually undetectable, collectively form a background [1]. Various mechanisms contribute to the SGWB, which can be categorized into astrophysical and cosmological origins [2-4]. Astrophysical origins primarily include double white dwarfs (DWDs) [5, 6], massive binary black holes (BBHs) [7], stellar-mass BBHs [8, 9], binary neutron stars (BNSs) [10, 11], extreme-mass-ratio inspirals (EMRIs) [12-14], and GW bursts [15]. Cosmological origins mainly contain inflation [16], first-order phase transitionss (PTs) [17], and cosmic defects [18, 19]. SGWB from either astrophysical origin [5, 7, 12, 20] or cosmological origin [21-27] can exhibit the potential for anisotropy [28]. The anisotropy can offer insights into the distribution of compact binaries within our Galaxy, the history of galaxy mergers [29-31], as well as the priceless primordial information of the early Universe [32-37]. For example, the structure of our Galaxy can be mapped using the lower multipole moments with l ≤ 4 [38, 39]. Conversely, the higher multipole moments, with l ≳ 100, are more suitable for detecting the subtle imprints of primordial fluctuations that are preserved in cosmological backgrounds [4, 34].  \nThe diverse origins of the SGWB allow it to span multiple frequency bands. At frequencies in the hundreds of hertz, ground-based detectors, despite not yet detecting SGWBs, have set an upper limit on the dimensionless energy density [40]. In the nano-hertz (nHz) frequency band, recent findings from Pulsar Timing Ar-  \nrayss (PTAs) have provided compelling evidence for the existence of SGWB [41-44]. Moving to the milli-hertz (mHz) frequency band, space-borne GW detectors, such as TianQin [45] and Laser Interferometer Space Antenna (LISA) [46], are expected to detect a Galactic foreground [47, 48].  \nWhen conducting GW detection, space-borne detectors face the challenge of canceling laser noise due to their motion. This issue can be addressed through time delay interferometry (TDI) techniques, which enable the formation of various types of TDI channels [49-53]. Among all the TDI channels, the unequal-arm Michelson, typically represented by the XYZ channel set, is the most widely utilized. This set can be further recombined to form the orthogonal AET channel set, with the T channel specifically designed to be insensitive to GWs. By employing the T channel as a noise monitor, the SGWB signal can be extracted through the correlation of data from other channels, which is often referred to as the null-channel method [54-62].  \nFrom the SGWB signal, one can derive a distribution of GW power across the sky, know as a sky map, to uncover the anisotropy of the SGWB. This process is termed 'map-making'. Initially, the SGWB signal produces a dirty map, which reflects the GW sky map as observed by detectors. The dirty map is susceptible to alterations due to variations in the detector's beam pattern. To refine the dirty map into a pure GW sky map, or clean map, it is essential to solve the deconvolution problem [63-65]. Thrane et al. are among the first to tackle this challenge [66], employing maximum likelihood estimates for the dirty map and then deconvolving it to recover a clean map [67-71]. Building on this foundation, Renzini et al. have introduced a frequentist maximum likelihood method to map GW power [72], conducting an all-sky analysis based on the complete likelihood of the  \ndata. In addition to the maximum likelihood method, a Bayesian spherical harmonic method has recently been proposed [73-75].  \nIn this paper, we harness maximum likelihood estimation to perform map-making with TianQin. In order to quantify the ability of TianQin, we simulate a Gaussian, stationary, and unpolarized point source into the data, both in the absence and presence of noise. To furnish a thorough analysis of the map-making process, we also examine the influence of several key factors, including the detector design, the spectral characteristics of the SGWB, the selected frequency range, and the specific TDI channel sets employed for the data analysis.  \nThe structure of the paper is organized as follows: we review the formalism of the SGWB detection in Sec. II. The maximum likelihood method and the deconvolution process for map-making are detailed in Sec. III. Based on the knowledge from the previous sections, we present the set up and results of map-making in IV. Finally, we give a brief summary in Sec. V.", 'A. Energy spectral density of stochastic background': "An SGWB is composed of a large number of GWs. Consequently, in transverse-traceless coordinates, the metric perturbations h ( t, ⃗x ) associated with the SGWB can be represented as a superposition of sinusoidal plane waves, characterized by frequency f and the unit vector ˆ n indicating the line-of-sight direction on the sky:  \nh ij ( t, ⃗x ) = ∑ P =+ , × ∫ ∞ -∞ d f ∫ S 2 d ˆ Ω ˆ n e P ij (ˆ n ) ˜ h P ( f, ˆ n ) × e i2 πf [ t +ˆ n · ⃗x ( t ) /c ] , (1)  \nwhere ⃗x denotes the location at which the GW measurement is conducted, c is the speed of light. The symbol e P ij (ˆ n ) signifies the GW polarization tensor corresponding to the polarization P , while ˜ h P ( f, ˆ n ) denotes the Fourier amplitude.  \nIn this work, we focus on the anisotropic, Gaussian, stationary, and unpolarized SGWB. Thus, the Fourier amplitude is treated as a zero-mean random variable, characterized by the one-sided power spectral density (PSD) P h :  \n⟨ ˜ h P ( f, ˆ n ) ˜ h ∗ P ' ( f, ˆ n ' ) ⟩ = 1 4 δ ( f -f ' ) δ PP ' δ 2 (ˆ n -ˆ n ' ) P h ( f, ˆ n ) , (2)  \nwhere the symbol * denotes the conjugate. The factor of 1/4 pertains to the definition of one-sided PSD, encompassing the summation over both polarizations. Given the assumption that the direction and intensity of the SGWB are independent, the P h ( f, ˆ n ) can be factorized  \ninto the spectral shape ¯ H ( f ) and angular distribution P h (ˆ n ):  \nP h ( f, ˆ n ) = ¯ H ( f ) P h (ˆ n ) . (3)  \nThe angular distribution P h (ˆ n ) can be expanded in terms of a set of basis:  \nP h (ˆ n ) = ∑ α P α e α (ˆ n ) , (4)  \nwhere e α (ˆ n ) can be spherical harmonic or pixel basis. With the definitions provided above, the dimensionless energy density spectrum Ω gw can be derived as follows [66]:  \nΩ gw ( f ) = 1 ρ c d ρ gw d(ln f ) = 2 π 2 f 3 3 H 2 0 ∫ S 2 d ˆ Ω ˆ n P h ( f, ˆ n ) , (5)  \nwhere the critical energy density ρ c = 3 H 2 0 c 2 / 8 πG , incorporating the gravitational constant G and the Hubble constant H 0 . The term d ρ gw represents the GW energy density contained within the frequency interval [ f , f +d f ]. Generally speaking, the energy density spectrum can be characterized by a power-law form with an index α and a reference frequency f r :  \nΩ gw ( f ) = Ω r ( f f r ) α , (6)  \nwhere Ω r denotes the reference energy density.  \nThus far, we has focused solely on the intrinsic anisotropy of the SGWB. However, the SGWB can also exhibit a degree of kinematic Doppler anisotropy, analogous to observations made in the cosmic microwave background (CMB) [76-79]. The Doppler effect can cause shifts in both the frequency and direction of the SGWB. By defining the unit vector ˆ v to indicate the direction of relative motion between the SGWB rest frame and the moving frame, and ˆ n D to signify the direction of the SGWB within the moving frame, the coefficient of frequency change due to the Doppler effect can be determined by  \nD = √ 1 -β 2 1 -β ξ , (7)  \nand the direction shift can be expressed as  \nˆ n = ˆ n D + ˆ v [( γ -1) ξ -γβ ] γ (1 -βξ ) . (8)  \nHere, β = | ⃗ v | /c , ξ = ˆ n D · ˆ v . Based on these definitions, the transformation of the energy density between the two frames can be written as:  \nΩ D GW ( f ) = ∫ S 2 d ˆ Ω ˆ n D 4 Ω GW ( D -1 f, ˆ n + ˆ v [( γ -1) ξ -γβ ] γ (1 -βξ ) ) . (9)  \nFurther details regarding this topic can be found in Ref. [77].", 'B. Detector noise': "Among the various space missions proposed for GW detection, this paper primarily focuses on TianQin [45] and LISA [46]. TianQin will consist of three identical drag-free satellites arranged in an equilateral triangular constellation orbiting the Earth. The orbital plane of TianQin is designed to always direct toward to J0806 ( θ s = -4 . 7 · , ϕ s = 120 . 5 · ), with each satellite separated approximately 1 . 7 × 10 5 km. Operating on a 'three months on + three months off' mode, TianQin will have an operational duration of six months each year. LISA, on the other hand, is designed to orbit the Sun, trailing Earth by about 20 · , with three satellites separated about 2 . 5 × 10 6 km apart. For the sake of convenience, shorthand notations will be employed in figures and equations: TQ for TianQin, and LS for LISA.  \nGiven the motion of space-borne detectors, phase noise becomes the dominant noise in GW detection. Fortunately, the phase noise can be depressed by several orders of magnitude through the implementation of TDI channels [59, 80-83]. Despite such alleviation, other types of noise, including position noise and acceleration noise, persist as significant considerations. For the conventional TDI channel set XYZ, which utilizes each satellite in harmony with its neighboring arm, the auto- and cross-PSD manifest as follows [84, 85]:  \nN II = 4 L 2 sin 2 [ f f ∗ ][ S p ( f ) + 2 ( cos 2 [ f f ∗ ] +1 ) S a ( f ) (2 πf ) 4 ] N IJ = -2 L 2 sin 2 [ f f ∗ ] cos [ f f ∗ ]( S p ( f ) + 4 S a ( f ) (2 πf ) 4 ) , (10)  \n̸  \nwhere I , J label the X/Y/Z channels, with I = J . S p and S a represent the position noise and the acceleration noise, respectively. Given the arm length L , the characteristic frequency f ∗ = c/ 2 πL . For TianQin and LISA, the characteristic frequencies f TQ ∗ and f LS ∗ are estimated to be approximately 0.28 Hz and 0.019 Hz.  \nBased on the channel set XYZ, one can further construct channel set AET by  \nA = 1 √ 2 (Z -X) , E = 1 √ 6 (X -2Y + Z) , T = 1 √ 3 (X + Y + Z) . (11)  \nWith the above construction, the auto-PSD of channel set AET can be derived from the the auto- and cross-PSD of channel set XYZ, as shown in Eq. (10) [85]. Furthermore, the cross-PSD of channel set AET is identically zero, thereby ensuring that the three channels are noiseorthogonal.", 'C. Detector response': 'In addition to detector noise, detector response is another important aspect to GW detection. The detector motion can lead to variations in the response to GWs over time. Hence, a long data stream needs to be divided into numerous segments. These segments should be long enough to cover the sensitive bands of the detector after transformation from the time domain to the frequency domain. Meanwhile, they should also be short enough to ensure that the response within each segment can be regarded as stationary, allowing the application of a short-term Fourier transform [1].  \nWithin a time interval τ wherein the response of channel I can be considered stationary, the frequency domain SGWB signal  \n˜ h I ( f, τ ) = ∑ P ∫ S 2 d ˆ Ω ˆ k F P I ( f, τ, ˆ n ) ˜ h P ( f, ˆ n ) e i2 πf ˆ n · ⃗x ( τ ) /c , (12)  \nwhere F P I denotes the response function. For the channel set XYZ,  \nF P I ( f, τ, ˆ n ) =2 sin 2 ( f f ∗ )[ ˆ u a ( τ )ˆ u b ( τ ) T ( f, ˆ n , ˆ u ( τ )) -ˆ v a ( τ )ˆ v b ( τ ) T ( f, ˆ n , ˆ v ( τ )) ] e P ab (ˆ n ) , (13)  \nwhere ˆ u and ˆ v are the orientations of arms associated with the channel I . The transfer function  \nT [ f, ˆ u ( τ ) , ˆ n ] = 1 2 [ sinc [ f 2 f ∗ (1 + ˆ n · ˆ u ( τ )) ] e -i f 2 f ∗ [3 -ˆ n · ˆ u ( τ )] +sinc [ f 2 f ∗ (1 -ˆ n · ˆ u ( τ )) ] e -i f 2 f ∗ [1 -ˆ n · ˆ u ( τ )] ] , (14)  \nwith sinc( x ) = sin x/x .  \nGiven the rotational symmetry of the X/Y/Z channels, an inherent equivalence exists among them in the SGWB detection. However, the A/E/T channels exhibit different responses: the A/E channels are responsive to GWs, whereas the T channel is typically less responsive unless it captures high-frequency GWs. Consequently, the T channel serves as a tool for monitoring detector noise [55, 56]. Furthermore, the response for A/E/T channels is not merely a linear superposition of the response for X/Y/Z channels; it also contains some exponential terms that significantly impact the response for the T channel. The specific form can be found in Eq. (18) of Ref. [85].  \nBy incorporating the response function and the varying arrival times of the signal at each channel, one can further define the antenna pattern [86]:  \nγ IJ ( f, τ, ˆ n ) = 1 2 ∑ P F P I ( f, τ, ˆ n ) F P ∗ J ( f, τ, ˆ n ) e i2 πf ˆ n ∆ ⃗x/c ,  \n(15)  \nFIG. 1. ORF of single channels for TianQin. The blue, orange, and green lines are denote the X/Y/Z channels, the A/E channels, and the T channel, respectively.  \n<!-- image -->  \nwhere the factor of 1 / 2 arises from the average of polarization. The separation vector between channels ∆ ⃗x = ⃗x I -⃗x J . Notably, when I = J , the antenna pattern refers to the auto-correlation of individual channels. Integrating the antenna pattern over all sky yields the overlap reduction function (ORF), which quantifies the reduction in correlation due to the non-parallel alignment and time delay between channels:  \nΓ IJ ( f, τ ) = 1 4 π ∫ S 2 d 2 ˆ Ω ˆ n γ IJ ( f, τ, ˆ n ) , (16)  \nwhere the factor of 1 / 4 π ensures that the ORF is properly scaled.  \nIn Fig. 1, we present the ORF for individual channels of TianQin. Given the symmetry inherent in the channel set, the ORFs for the X/Y/Z channels and the A/E channels are identical, respectively. For the detection channels X/Y/Z and A/E, the ORF scales with f 2 , while the ORF of the null channel T scales with f 8 . At lower frequencies, the ORF of the null channel is markedly inferior compared to that of the detection channels. Nonetheless, as the frequency increases, the ORF of the null channel can dramatically rise, thereby transitioning into a detection role. Fig. 2 further illustrates the corresponding antenna pattern when TianQin is positioned at perihelion, within the ecliptic coordinates [87]. For clarity, the antenna pattern is standardized relative to its apex, termed the hot spot. At the low-frequency threshold ( f = f TQ ∗ / 3), the hot spots of the antenna pattern for detection channels converge with the fixed orientation of TianQin, designated by the coordinates (lon, lat)= (120 . 5 · , -4 . 7 · ). As the frequency ascends to 3 f TQ ∗ , the hot spots migrate towards the lateral directions where TianQin is oriented. For the null channel, the antenna pattern demonstrates its minimal intensity in the direction TianQin is facing. Moreover, the overall intensity  \ndistribution across the celestial sphere remains relatively consistent despite the frequency surge.', 'III. MAP-MAKING ANALYSIS': "This section is bifurcated into two pivotal aspects: the initial part explores the utilization of the maximum likelihood method to ascertain the optimal estimator for the SGWB's dirty map, while the subsequent part delves into the application of singular-value decomposition (SVD) to transform a dirty map into a clean map.", 'A. Maximum likelihood method': "The data d I is the linear addition of the SGWB signal h I and the noise n I . In the frequency domain, we have  \n˜ d I ( f, τ ) = ˜ h I ( f, τ ) + ˜ n I ( f, τ ) . (17)  \nThe correlation measurement can be constructed by correlating two data:  \nD IJ ( f, τ ) = 2 τ d I ( f, τ ) d ∗ J ( f, τ ) , (18)  \nThe factor of 2 is consistent with the definition of onesided PSD. Combined with Eqs. (2), (3), and (12), the expectation of the correlation measurement is given by [85, 88]:  \n⟨ D IJ ( f, τ ) ⟩ = ¯ H ( f ) ∫ d 2 ˆ Ω ˆ n γ IJ ( f, τ, ˆ n ) P h (ˆ n )+ N IJ ( f, τ ) . (19)  \nTo ensure the expectation value attributes pure signals, one can further define the measurement as follows [89]:  \nC IJ ( f, τ ) = D IJ ( f, τ ) -N IJ ( f, τ ) . (20)  \nIn this work, we utilize the pixel basis grounded in HEALPix [90]. HEALPix, a grid-based methodology, is employed for the analysis and visualization of celestial data. It is uniquely adept at managing spherical data, characterized by its capacity to uniformly segment the sphere into equal-area and iso-latitude regions. Within this paradigm, the angular distribution  \nP h (ˆ n ) = P ˆ n ' δ (ˆ n , ˆ n ' ) . (21)  \nThen, the expectation of and variance of the above measurement C IJ [48]  \nµ IJ ( f, τ ) = ⟨ C IJ ( f, τ ) ⟩ = ∑ ˆ n ¯ H ( f ) γ IJ ( f, τ, ˆ n ) P ˆ n σ 2 IJ ( f, τ ) = ⟨ C 2 IJ ( f, τ ) ⟩ - ⟨ C IJ ( f, τ ) ⟩ 2 = 1 + δ IJ 2 T ∆ f N II ( f, τ ) N JJ ( f, τ ) W IJ ( f, τ ) , (22)  \nFIG. 2. Antenna pattern of different channels for TianQin, which is plotted on a Mollweide projection of the sky in ecliptic coordinates. In this representation, the blue star signifies the pointing direction of TianQin ([lon, lat]= [120 . 5 · , -4 . 7 · ]). Top and bottom panels correspond to frequencies of f TQ ∗ / 3 and 3 f TQ ∗ , respectively.  \n<!-- image -->  \nwhere ∆ f denotes the frequency resolution, the correction function W IJ is specifically designed for scenarios involving a strong SGWB:  \nW IJ ( f, τ ) = 1 + µ II ( f, τ ) N JJ ( f, τ ) + µ JJ ( f, τ ) N II ( f, τ ) N II ( f, τ ) N JJ ( f, τ ) + µ II ( f, τ ) µ JJ ( f, τ ) + (1 -δ IJ ) µ 2 IJ ( f, τ ) N II ( f, τ ) N JJ ( f, τ ) .  \n(23)  \nIn terms of Eq. (22), the signal-to-noise ratio (SNR) is expressed as follows [86]:  \nρ IJ = √ √ √ √ ∑ f,τ µ 2 IJ ( f, τ ) σ 2 IJ ( f, τ ) , (24)  \nand the log-likelihood for the P ˆ n is given by  \nln [ L ( n )] = const . -∑ fτ ∣ ∣ C IJ ( f, τ ) -¯ H ( f ) γ IJ ( f, τ, ˆ n ) P ˆ n ∣ ∣ 2 2 σ 2 IJ ( f, τ ) . (25)  \nThe process of maximizing the likelihood function is fundamentally aligned with optimizing the SNR. Consequently, the estimator that maximizes the likelihood for the clean map P ˆ n can be derived as [66]:  \nˆ P ˆ n = ( F -1 ) IJ ˆ n ˆ n ' X IJ ˆ n ' , (26)  \nwhere the dirty map represents the GW sky as observed by the detectors:  \nX IJ ˆ n ' = ∑ fτ γ ∗ IJ ( f, τ, ˆ n ' ) ¯ H ( f ) σ 2 IJ ( f, τ ) C IJ ( f, τ ) , (27)  \nand the covariance matrix of the dirty map  \nF IJ ˆ n ˆ n ' = ∑ fτ γ IJ ( f, τ, ˆ n ) ¯ H 2 ( f ) σ 2 IJ ( f, τ ) γ ∗ IJ ( f, τ, ˆ n ' ) . (28)  \nIt is crucial to note that the inverse of F IJ ˆ n ˆ n ' serves as the covariance matrix for the clean map:  \n⟨ ˆ P ˆ n ˆ P ∗ ˆ n ' ⟩ - ⟨ ˆ P ˆ n ⟩⟨ ˆ P ∗ ˆ n ' ⟩ ≈ ( F -1 ) IJ ˆ n ˆ n ' . (29)  \nHence, F IJ ˆ n ˆ n ' is commonly referred to as the Fisher information matrix (FIM). Furthermore, when considering multiple channel pairs { IJ } , the aggregation of these individual contributions leads to the formation of the total dirty map and total FIM:  \nX tot ˆ n = ∑ IJ X IJ ˆ n F tot ˆ n ˆ n ' = ∑ IJ F IJ ˆ n ˆ n ' . (30)", 'B. Deconvolution': "Eq. (26) indicates the possibility of recovering the clean map by inverting the FIM F IJ ˆ n ˆ n ' . However, the FIM is often a singular matrix, necessitating regularization for inversion. A common method to invert the FIM involves a SVD regularization scheme. Given the Hermitian nature of the FIM, its SVD representation takes the form  \nF = U Σ U ∗ , (31)  \nwhere U is a unitary matrix, and Σ is a diagonal matrix containing positive real eigenvalues s i of the FIM. When arranging the diagonal elements of Σ in descending order, a threshold s min can be selected to condition the matrix. Values below this threshold can be typically replaced with infinity or the smallest eigenvalue above the cutoff, resulting in the matrix Σ ' . Following this procedure, one can directly derive the regularized F ' using the resulting matrix Σ ' :  \nF ' = U Σ ' U ∗ , (32)  \nwith its inverse  \nF '-1 = U Σ '-1 U ∗ . (33)  \nBy multiplying the inverted-regularized FIM with the dirty map, the clean map can be derived:  \nˆ P ˆ n = ( F '-1 ) tot ˆ n ˆ n ' X tot ˆ n ' . (34)  \nThe process of discarding eigenvalues during SVD can lead to the loss of information associated with pixels of weak intensity, which can result in the generated clean map differing from the true sky map. In addition, when solving the equation represented by Eq. (34), round-off errors are introduced. These errors are inherently related to the condition number of the FIM.  \nWithin the framework of a system of equations denoted as Ax = b , the condition number κ ( A ) associated with the coefficient matrix A functions as an indicator of singularity, illustrating its susceptibility to perturbations in the input data. It is defined as follows [69, 91]:  \nκ ( A ) = ∥ A ∥∥ A -1 ∥ , (35)  \nwhere ∥ ... ∥ denotes the norm of matrix. In numerous applications, the 2-norm is often preferred due to its direct relationship with the eigenvalues of the matrix:  \n∥ A ∥ 2 = √ λ max ( A ∗ A ) , (36)  \nwhere λ max labels the maximum eigenvalue of the matrix. Given that the FIM usually exhibits off-diagonal elements with negligible imaginary parts, it can be accurately approximated as a symmetric, positive definite matrix. Under this assumption, the condition number κ ( A ) of the matrix simplifies to the ratio of its maximum to minimum eigenvalues:  \nκ ( A ) = λ max ( A ) /λ min ( A ) . (37)  \nThe condition number measures the extent to which the resulting clean map ˆ P ˆ n fluctuates in response to changes in the dirty map X . A high condition number suggests that the clean map is highly sensitive to even slight alterations in the dirty map, making it difficult to obtain accurate results. Conversely, a low condition number implies that changes in the dirty map have a minimal impact on the clean map, indicating that the FIM is 'well-behaved' and thus more effective at recovering the clean map.", 'IV. SET UP AND RESULTS': 'In this section, we will demonstrate the capability of TianQin in recovering the sky map for the SGWB. We narrow our attention to the pixel-based decomposition, with a specific emphasis on point sources.  \nFIG. 3. Auto-correlation data for the X channel of TianQin, with the top and bottom panels involving spectral indices α for SGWB at 2/3 and 3. Within each panel, the unfolded data is symbolized by the gray line, whereas the folded data is illustrated by the red line. Additionally, the black and green lines denote the noise and signal PSDs, respectively, corresponding to an SNR of approximately 16 over the course of 1-year operation.  \n<!-- image -->', 'A. Data simulation': "We start by generating random Gaussian SGWB signals, of which the covariance is obtained by summing contributions from each pixel across the all sky:  \nS in IJ ( f, τ ) = 4 π N pix ∑ ˆ n γ IJ ( f, τ, ˆ n ) ¯ H in ( f ) P in ˆ n , (38)  \nwhere the intensity of SGWB is determined by the spectral shape ¯ H in ( f ) and the angular distribution P in ˆ n . N pix denotes the number of pixels. Increasing the N pix enhances resolution but also escalates computational demands. To strike a balance, we opt for N in pix = 768 for both injecting and subsequently outputting the signal. In addition, to handle the computational load associated with analyzing the real-time response of the detector, which is subject to temporal variations, we divide the data into discrete time segments. We assume that within each brief interval, the detector's response  \nremains consistent. To encapsulate the full segment's response, we employ the response at the midpoint of the interval as a transient representation. For TianQin, we select segments of 3600 seconds, where the directional change in response is minimal [92]. This allows a compact period of 3.64 days to be divided into 87 segments. On the other hand, we simulate the random Gaussian frequency-domain noise based on the noise PSD provided in Eq. (10). By combining the generated signal with the noise, we derive the correlation measurement D IJ ( f ) as described in Eq. (18).  \nNext, we examine data gathered during the operational cycle of TianQin, which entails a three-month on phase followed by a subsequent three-month off phase. Incorporating a year's worth of data thus yields 50 periods, and these 50 periods can be then aggregated into a single period for analysis. By aggregating multiple periods into a single unit, one can significantly reduce computational effort. While this method may result in the loss of some low-frequency data, it is deemed acceptable as the data lies outside the sensitive frequency range of TianQin. The aforementioned data folding [93, 94] process is equivalent to rearranging the time summation in Eq. (27), maintaining the formulaic expression without alteration. Fig. 3 presents a sample of the data from the auto-correlated X channel for SGWB signal with different spectral indices α . The gray line represents the unfolded data from a single period, while the red line showcases the folded data after averaging across the 50 periods. For comparison purposes, we also plot the PSDs of the noise and SGWB signal using black and green lines, respectively. The SNRs for the signals with α values of 2/3 and 3 are approximately 16, achieved over the course of TianQin's 1-year operation.", 'B. Map-making': 'The process now pivots towards map-making. Unless otherwise specified, we will employ a Gaussian, stationary, unpolarized point source expanded up to multipole moments of l up to 4, with the signal PSD illustrated in Fig. 3. As depicted in Fig. 4, we initiate our analysis by injecting a pure signal with α = 2 / 3 into the XYZ channel set for data simulations. The left panel displays the resulting map after injection, with SGWB signal located at [lon, lat]= [5 π/ 3 , π/ 6]. It is noteworthy that, to ensure a non-negative distribution of spherical harmonics, Clebsch-Gordan coefficients are employed [73, 95]. The middle and right panels exhibit the clean maps recovered by truncating the data to f TQ ∗ / 3 and 3 f TQ ∗ , respectively. For straightforward comparison, these maps have been normalized to the brightest point of the injected map.  \nThe fixed orbital orientation of TianQin renders it incapable of distinguishing between two signals that are symmetrically aligned with the orbital plane. As a result, upon injecting a point signal, two point signals are symmetrically recovered relative to the orbital plane (marked  \nby the orange dashed lines), with the intensity of these recovered signals approximating half of the original injected signal. This mirror symmetry poses a substantial hurdle in precisely ascertaining the true orientation of the signal. We will further explore this issue in Sec. IV C by analyzing it through the lens of the Doppler effect.  \nFurthermore, a noteworthy trend in the antenna pattern γ IJ manifests as the increasing prominence of the imaginary component with rising frequency. This evolution is notably significant as it amplifies the process of inverting the FIM, a pivotal procedure for capturing finescale details 1 . As a result, determining the optimal frequency truncation becomes imperative. Under the noisefree premise, as depicted in the middle and right panels of Fig. 4, the distinction between a cutoff frequency of f TQ ∗ / 3 and 3 f TQ ∗ is negligible. To rigorously evaluate the effects of frequency truncation, we utilized the JensenShannon divergence (JSD) [96], a metric that quantifies the divergence between two distributions. A minimum value for JSD of 0 signifies that the two distributions are precisely identical, while a maximum value of ln 2 indicates that the two distributions are entirely distinct. With a computed JSD value of 0.017 nat, we infer that truncating the frequency to f TQ ∗ / 3 suffices for noise-free conditions.  \nThe effect of noise on map-making is a substantial concern. In Fig. 5, we further compare the results of mapmaking with and without noise in the top and bottom panels, where the data is truncated at 3 f TQ ∗ . For comparative analysis, the left and right panels refer to two different spectral shapes, α = 2 / 3 and α = 3, respectively. In the presence of noise, the capability to map the SGWB with α = 2 / 3 is significantly hindered, resulting in an increased JSD of 0.127 nat. Conversely, when α = 3, the recovered map becomes more accurate, with a JSD of 0.062 nat. This trend aligns with the observations in Fig. 4, where an increase in the cutoff frequency led to the gradual absorption of the high-frequency components of the signal. The transition from α = 2 / 3 to α = 3 noticeably intensifies the high-frequency components, which in turn facilitates the results for map-making.', 'C. Doppler-induced effects': "As previously discussed, TianQin faces a challenge due to mirror symmetry in map-making, as shown in Fig. 4 and Fig. 5. Eq. (9) implies that Doppler shifts can serve as a useful complementary tool for investigating the frequency profile of the SGWB and its inherent rest-frame anisotropies. We will explore strategies to tackle this symmetry issue in the context of Doppler effects.  \nCoupled with the Earth's velocity of 3 × 10 4 m / s, the space-borne detector's motion towards the signal direction leads to a frequency increase in the real signal and a  \nFIG. 4. Normalized sky map of a pure SGWB signal. The left panel illustrates the injected map for the SGWB signal which corresponds to the top panel of Fig. 3, positioned at [lon, lat]= [5 π/ 3 , π/ 6]. The middle and right panels refer to the recovered map, with the data truncated at f TQ ∗ / 3 and 3 f TQ ∗ , respectively. The orange dashed line indicates the orbital plane of TianQin, serving as a visual reference to highlight the symmetry observed in the recovered maps.  \n<!-- image -->  \nFIG. 5. Recovered sky map for the injected SGWB signal, as referenced in Fig. 3. The top and bottom panels correspond to the pure-signal and signal-plus-noise scenarios, respectively. In both cases, the cutoff frequency of the data is set to 3 f TQ ∗ . The orange dashed line signifies the orbit plane of TianQin.  \n<!-- image -->  \ndecrease in the frequency of the symmetric signal. This variation in frequency is instrumental in determining the precise location of a point source. Considering the point source shown in Fig. 4, we calculate the Doppler shift using Eq. (7). Our results indicate that the Doppler effect can cause frequency alterations and directional shifts of approximately 6 × 10 -5 and 2 × 10 -5 , respectively. However, when comparing the top and bottom panels of Fig. 5, it becomes evident that noise-induced errors reach a level of 10 -1 . Consequently, we conclude that the effects caused by the Doppler shift are not strong enough to overcome the mirror symmetry for the mapmaking process.", 'D. Singular matrix analysis': "The condition number of the FIM plays a critical role in the map-making process. Given that variations in the FIM can occur across different channels, it is essential to conduct a comparative analysis. To this end, we plan to utilize both the XYZ and AET channel sets for this study. Our goal is to evaluate and understand the differences in their performance and effectiveness within the mapmaking process.  \nIn Fig. 6, we calculate the condition number of the FIM using Eq. (37). The eigenvalues of the matrix are plotted in descending order, and we vary the number of eigenvalues retained for comparative analysis. It is important to recognize that while retaining a higher count  \n=2/3  \nFIG. 6. Condition number of FIM for different channel sets of TianQin and LISA. In the top panel, we display the condition number while holding the spectral index α constant, investigating the influence of altering the cutoff frequency and selecting different channel sets. The bottom panel maintains a constant cutoff frequency at 3 f ∗ , which allows us to examine the impact of different spectral indices on the condition number.  \n<!-- image -->  \nof eigenvalues incorporates more information, it generally inflate the condition number of the FIM. Therefore, our analysis focuses on the condition numbers associated with the first 192 eigenvalues.  \nThe top panel compares the condition numbers for TianQin and LISA at cutoff frequencies of f ∗ / 3 and 3 f ∗ , with a fixed spectral index α of 2/3. The condition number for LISA consistently remains lower than that for TianQin, especially at higher eigenvalue counts. The discrepancy can be attributed to the differing operational designs of the detectors: LISA is capable of surveying various directions in the sky throughout the year, offering dynamic coverage that enhances its condition number performance. In contrast, TianQin's operational focus remains fixed on J0806, leading to a stable and less varying antenna pattern. This stability significantly inflates the condition number of the FIM. Despite these differences, both detectors share certain characteristics in the  \nmap-making process. As expected, there is a notable correlation between an increase in the cutoff frequency and a decrease in the condition number, regardless of whether the XYZ or AET channel set is used. Interestingly, when the cutoff frequency is set at f ∗ / 3, the XYZ channel set typically exhibits a higher condition number compared to the AET channel set. However, the result dynamic changes when the cutoff frequency is elevated to 3 f ∗ , with the AET channel set then showing a higher condition number. The selection of both cutoff frequency and channel set plays a crucial role in influencing the condition number of the FIM for map-making.  \nThe results presented above are predicated on a specific spectral index, α = 2 / 3. Should this index vary, the outcomes would necessarily adjust to accommodate such changes. To illustrate this point, the bottom panel keeps the cutoff frequency constant at 3 f ∗ , and compares the results for α = 2 / 3 with those for α = 3. A distinct pat-  \nis observed: with α = 2 / 3, both TianQin and LISA tend to show a higher condition number for the XYZ channel set compared to the AET channel set; however, the situation is reversed when α = 3. This observation underscores that the spectral shape of the SGWB has a significant impact on the accuracy and reliability of the map-making process.", 'V. SUMMARY': 'In this paper, we have explored the map-making for the anisotropic SGWB using TianQin. We utilized the maximum likelihood method to construct maps from simulated data. By using a point source as a case study, we reconstructed its clean map up to multipole moments of order l ≤ 4 both in the absence and presence of noise, achieving an SNR of 16. We observed that the inclusion of noise necessitates the use of higher-frequency data to achieve satisfactory map-making results. A notable limitation of TianQin lies in distinguishing between a point source and its mirror position relative to the orbital plane. To tackle this challenge, we investigated the potential of harnessing Doppler-induced effects to break this symmetric degeneracy. However, the theoretical deviation in the signal between two symmetrical directions, induced by the Doppler effect, is considerably smaller than the statistical error introduced by the noise. Distinguishing signals originating from symmetrical directions remains a tough challenge.  \nRegarding both TianQin and LISA, we further investigated various factors that could influence the mapmaking process by analyzing the condition number of the FIM with respect to three key aspects: the cutoff frequency of the data, the choice of channel set, and the spectral profile of the SGWB. Our results indicated that adopting a sufficiently high cutoff frequency, specifically  \n- [1] J. D. Romano and N. J. Cornish, Living Rev. Rel. 20 , 2 (2017), arXiv:1608.06889 [gr-qc].\n- [2] N. Christensen, Rept. Prog. Phys. 82 , 016903 (2019), arXiv:1811.08797 [gr-qc].\n- [3] J. D. Romano (2019) arXiv:1909.00269 [gr-qc].\n- [4] N. Bartolo et al. (LISA Cosmology Working Group), JCAP 11 , 009 (2022), arXiv:2201.08782 [astro-ph.CO].\n- [5] S.-J. Huang, Y.-M. Hu, V. Korol, P.-C. Li, Z.-C. Liang, Y. Lu, H.-T. Wang, S. Yu, and J. Mei, Phys. Rev. D 102 , 063021 (2020), arXiv:2005.07889 [astro-ph.HE].\n- [6] V. Korol, E. M. Rossi, P. J. Groot, G. Nelemans, S. Toonen, and A. G. A. Brown, Mon. Not. Roy. Astron. Soc. 470 , 1894 (2017), arXiv:1703.02555 [astro-ph.HE].\n- [7] H.-T. Wang et al. , Phys. Rev. D 100 , 043003 (2019), arXiv:1902.04423 [astro-ph.HE].\n- [8] S. Liu, Y.-M. Hu, J.-d. Zhang, and J. Mei, Physical Review D 101 (2020), 10.1103/physrevd.101.103027.\n- [9] H. Wang, I. Harry, A. Nitz, and Y.-M. Hu, Phys. Rev. D 109 , 063029 (2024), arXiv:2304.10340 [astro-ph.HE].  \nat 3 f ∗ , can improve the mapping capabilities. Notably, at this cutoff frequency, the XYZ channel set can outperform the AET channel set in terms of map-making efficacy when the spectral index α is set at 2/3. Additionally, the map-making process for an SGWB with a spectral index α of 3 consistently yielded better results compared to an index α of 2/3.  \nIt is noteworthy that the injection of a point source represents an idealized scenario, selected because it provides a straightforward example to analyze the detection capabilities for anisotropic SGWB. Furthermore, we focused on single space-borne detectors, underscoring the critical necessity for accurate noise level estimation to facilitate effective SGWB detection. Failure to do so can significantly undermine the accuracy of detection [97]. A strategy to address this issue entails crosscorrelating data from multiple detectors [92, 98]. Unlike ground-based detectors, the complex relative motion of space-borne detectors substantially amplifies the computational demands for determining antenna patterns during the map-making process. Research into this topic is ongoing and will be expounded upon in forthcoming publications.', 'ACKNOWLEDGMENTS': "This work has been supported by the National Key Research and Development Program of China (No. 2020YFC2201400), the Natural Science Foundation of China (Grants No. 12173104), the Natural Science Foundation of Guangdong Province of China (Grant No. 2022A1515011862), and the Guangdong Basic and Applied Basic Research Foundation(Grant No. 2023A1515030116). Z.C.L. is supported by the Guangdong Basic and Applied Basic Research Foundation (Grant No. 2023A1515111184).  \n- [10] M. Surace, K. D. Kokkotas, and P. Pnigouras, Astron. Astrophys. 586 , A86 (2016), arXiv:1512.02502 [astroph.CO].\n- [11] D. Talukder, E. Thrane, S. Bose, and T. Regimbau, Phys. Rev. D 89 , 123008 (2014), arXiv:1404.4025 [gr-qc].\n- [12] H.-M. Fan, Y.-M. Hu, E. Barausse, A. Sesana, J.-d. Zhang, X. Zhang, T.-G. Zi, and J. Mei, Phys. Rev. D 102 , 063016 (2020), arXiv:2005.08212 [astro-ph.HE].\n- [13] T.-G. Zi, J.-D. Zhang, H.-M. Fan, X.-T. Zhang, Y.-M. Hu, C. Shi, and J. Mei, Phys. Rev. D 104 , 064008 (2021), arXiv:2104.06047 [gr-qc].\n- [14] C.-Q. Ye, H.-M. Fan, A. Torres-Orjuela, J.-d. Zhang, and Y.-M. Hu, (2023), arXiv:2310.03520 [gr-qc].\n- [15] Z. Wu, Y.-M. Hu, H.-M. Fan, and I. S. Heng, Phys. Rev. D 109 , 103004 (2024), arXiv:2308.15354 [gr-qc].\n- [16] A. H. Guth and S. Y. Pi, Phys. Rev. Lett. 49 , 1110 (1982).\n- [17] C. J. Hogan, Phys. Lett. B 133 , 172 (1983).\n- [18] T. W. B. Kibble, J. Phys. A 9 , 1387 (1976).\n- [19] P. Auclair et al. (LISA Cosmology Working Group), Living Rev. Rel. 26 , 5 (2023), arXiv:2204.05434 [astroph.CO].\n- [20] S. Liu, Y.-M. Hu, J.-d. Zhang, and J. Mei, Phys. Rev. D 101 , 103027 (2020), arXiv:2004.14242 [astro-ph.HE].\n- [21] X. Siemens, V. Mandic, and J. Creighton, Phys. Rev. Lett. 98 , 111101 (2007), arXiv:astro-ph/0610920.\n- [22] D. G. Figueroa, M. Hindmarsh, and J. Urrestilla, Phys. Rev. Lett. 110 , 101302 (2013), arXiv:1212.5458 [astroph.CO].\n- [23] V. Mandic, S. Bird, and I. Cholis, Phys. Rev. Lett. 117 , 201102 (2016), arXiv:1608.06699 [astro-ph.CO].\n- [24] K. Saikawa and S. Shirai, JCAP 05 , 035 (2018), arXiv:1803.01038 [hep-ph].\n- [25] C. Caprini and D. G. Figueroa, Class. Quant. Grav. 35 , 163001 (2018), arXiv:1801.04268 [astro-ph.CO].\n- [26] P. Auclair et al. , JCAP 04 , 034 (2020), arXiv:1909.00819 [astro-ph.CO].\n- [27] X. Wang, F. P. Huang, and X. Zhang, JCAP 05 , 045 (2020), arXiv:2003.08892 [hep-ph].\n- [28] B. Allen and A. C. Ottewill, Phys. Rev. D 56 , 545 (1997), arXiv:gr-qc/9607068.\n- [29] G. Cusin, I. Dvorkin, C. Pitrou, and J.-P. Uzan, Phys. Rev. Lett. 120 , 231101 (2018), arXiv:1803.03236 [astroph.CO].\n- [30] G. Capurri, A. Lapi, C. Baccigalupi, L. Boco, G. Scelfo, and T. Ronconi, JCAP 11 , 032 (2021), arXiv:2103.12037 [gr-qc].\n- [31] N. Bellomo, D. Bertacca, A. C. Jenkins, S. Matarrese, A. Raccanelli, T. Regimbau, A. Ricciardone, and M. Sakellariadou, JCAP 06 , 030 (2022), arXiv:2110.15059 [gr-qc].\n- [32] M. Geller, A. Hook, R. Sundrum, and Y. Tsai, Phys. Rev. Lett. 121 , 201303 (2018), arXiv:1803.10780 [hepph].\n- [33] A. C. Jenkins and M. Sakellariadou, Phys. Rev. D 98 , 063509 (2018), arXiv:1802.06046 [astro-ph.CO].\n- [34] Y. Li, F. P. Huang, X. Wang, and X. Zhang, Phys. Rev. D 105 , 083527 (2022), arXiv:2112.01409 [astro-ph.CO].\n- [35] S. Wang, V. Vardanyan, and K. Kohri, Phys. Rev. D 106 , 123511 (2022), arXiv:2107.01935 [gr-qc].\n- [36] S. Profumo and F. Yang, Phys. Rev. D 108 , 103528 (2023), arXiv:2306.07454 [astro-ph.CO].\n- [37] F. Schulze, L. Valbusa Dall'Armi, J. Lesgourgues, A. Ricciardone, N. Bartolo, D. Bertacca, C. Fidler, and S. Matarrese, JCAP 10 , 025 (2023), arXiv:2305.01602 [gr-qc].\n- [38] C. Ungarelli and A. Vecchio, Phys. Rev. D 64 , 121501 (2001), arXiv:astro-ph/0106538.\n- [39] K. Breivik, C. M. F. Mingarelli, and S. L. Larson, Astrophys. J. 901 , 4 (2020), arXiv:1912.02200 [astro-ph.GA].\n- [40] R. Abbott et al. (KAGRA, Virgo, LIGO Scientific), Phys. Rev. D 104 , 022004 (2021), arXiv:2101.12130 [gr-qc].\n- [41] G. Agazie et al. (NANOGrav), Astrophys. J. Lett. 951 , L8 (2023), arXiv:2306.16213 [astro-ph.HE].\n- [42] H. Xu et al. , Res. Astron. Astrophys. 23 , 075024 (2023), arXiv:2306.16216 [astro-ph.HE].\n- [43] J. Antoniadis et al. (EPTA, InPTA:), Astron. Astrophys. 678 , A50 (2023), arXiv:2306.16214 [astro-ph.HE].\n- [44] D. J. Reardon et al. , Astrophys. J. Lett. 951 , L6 (2023), arXiv:2306.16215 [astro-ph.HE].\n- [45] J. Luo et al. (TianQin), Class. Quant. Grav. 33 , 035010 (2016), arXiv:1512.02076 [astro-ph.IM].\n- [46] P. Amaro-Seoane et al. (LISA), (2017), arXiv:1702.00786  \n[astro-ph.IM].  \n- [47] M. Benacquista and K. Holley-Bockelmann, Astrophys. J. 645 , 589 (2006), arXiv:astro-ph/0504135.\n- [48] Z.-C. Liang, Y.-M. Hu, Y. Jiang, J. Cheng, J.-d. Zhang, and J. Mei, Phys. Rev. D 105 , 022001 (2022), arXiv:2107.08643 [astro-ph.CO].\n- [49] M. Tinto and J. W. Armstrong, Phys. Rev. D 59 , 102003 (1999).\n- [50] D. A. Shaddock, Phys. Rev. D 69 , 022001 (2004), arXiv:gr-qc/0306125.\n- [51] M. Vallisneri, Phys. Rev. D 72 , 042003 (2005), [Erratum: Phys.Rev.D 76, 109903 (2007)], arXiv:gr-qc/0504145.\n- [52] M. Muratore, D. Vetrugno, and S. Vitale, Class. Quant. Grav. 37 , 185019 (2020), arXiv:2001.11221 [astro-ph.IM].\n- [53] M. Tinto and S. V. Dhurandhar, Living Rev. Rel. 24 , 1 (2021).\n- [54] M. Tinto, J. W. Armstrong, and F. B. Estabrook, Phys. Rev. D 63 , 021101 (2001).\n- [55] C. J. Hogan and P. L. Bender, Phys. Rev. D 64 , 062002 (2001), arXiv:astro-ph/0104266.\n- [56] M. R. Adams and N. J. Cornish, Phys. Rev. D 82 , 022002 (2010), arXiv:1002.1291 [gr-qc].\n- [57] T. L. Smith, T. L. Smith, R. R. Caldwell, and R. Caldwell, Phys. Rev. D 100 , 104055 (2019), [Erratum: Phys.Rev.D 105, 029902 (2022)], arXiv:1908.00546 [astro-ph.CO].\n- [58] G. Boileau, N. Christensen, R. Meyer, and N. J. Cornish, Phys. Rev. D 103 , 103529 (2021), arXiv:2011.05055 [grqc].\n- [59] M. Muratore, D. Vetrugno, S. Vitale, and O. Hartwig, Phys. Rev. D 105 , 023009 (2022), arXiv:2108.02738 [grqc].\n- [60] J. Cheng, E.-K. Li, Y.-M. Hu, Z.-C. Liang, J.-d. Zhang, and J. Mei, Phys. Rev. D 106 , 124027 (2022), arXiv:2208.11615 [gr-qc].\n- [61] C. Gowling, M. Hindmarsh, D. C. Hooper, and J. Torrado, JCAP 04 , 061 (2023), arXiv:2209.13551 [astroph.CO].\n- [62] J. Alvey, U. Bhardwaj, V. Domcke, M. Pieroni, and C. Weniger, Phys. Rev. D 109 , 083008 (2024), arXiv:2309.07954 [gr-qc].\n- [63] N. J. Cornish, Class. Quant. Grav. 18 , 4277 (2001), arXiv:astro-ph/0105374.\n- [64] A. Taruya and H. Kudoh, Phys. Rev. D 72 , 104015 (2005), arXiv:gr-qc/0507114.\n- [65] A. Taruya, Phys. Rev. D 74 , 104022 (2006), arXiv:grqc/0607080.\n- [66] E. Thrane, S. Ballmer, J. D. Romano, S. Mitra, D. Talukder, S. Bose, and V. Mandic, Phys. Rev. D 80 , 122002 (2009), arXiv:0910.0858 [astro-ph.IM].\n- [67] A. I. Renzini and C. R. Contaldi, Mon. Not. Roy. Astron. Soc. 481 , 4650 (2018), arXiv:1806.11360 [astro-ph.IM].\n- [68] J. Suresh, A. Ain, and S. Mitra, Phys. Rev. D 103 , 083024 (2021), arXiv:2011.05969 [gr-qc].\n- [69] D. Agarwal, J. Suresh, S. Mitra, and A. Ain, Phys. Rev. D 104 , 123018 (2021), arXiv:2105.08930 [gr-qc].\n- [70] A. I. Renzini, J. D. Romano, C. R. Contaldi, and N. J. Cornish, Phys. Rev. D 105 , 023519 (2022), arXiv:2107.02292 [gr-qc].\n- [71] L. Xiao, A. I. Renzini, and A. J. Weinstein, Phys. Rev. D 107 , 122002 (2023), arXiv:2211.10010 [gr-qc].\n- [72] C. R. Contaldi, M. Pieroni, A. I. Renzini, G. Cusin, N. Karnesis, M. Peloso, A. Ricciardone, and G. Tasinato, Phys. Rev. D 102 , 043502 (2020), arXiv:2006.03313", '[astro-ph.CO].': "- [73] S. Banagiri, A. Criswell, T. Kuan, V. Mandic, J. D. Romano, and S. R. Taylor, Mon. Not. Roy. Astron. Soc. 507 , 5451 (2021), arXiv:2103.00826 [astro-ph.IM].\n- [74] A. I. Renzini, B. Goncharov, A. C. Jenkins, and P. M. Meyers, Galaxies 10 , 34 (2022), arXiv:2202.00178 [gr-qc].\n- [75] S. Rieck, A. W. Criswell, V. Korol, M. A. Keim, M. Bloom, and V. Mandic, (2023), arXiv:2308.12437 [astro-ph.IM].\n- [76] A. Kosowsky and T. Kahniashvili, Phys. Rev. Lett. 106 , 191301 (2011), arXiv:1007.4539 [astro-ph.CO].\n- [77] G. Cusin and G. Tasinato, JCAP 08 , 036 (2022), arXiv:2201.10464 [astro-ph.CO].\n- [78] N. M. J. Cruz, A. Malhotra, G. Tasinato, and I. Zavala, (2024), arXiv:2402.17312 [gr-qc].\n- [79] L. Heisenberg, H. Inchausp'e, and D. Maibach, (2024), arXiv:2401.14849 [gr-qc].\n- [80] M. Tinto, S. Dhurandhar, and D. Malakar, Phys. Rev. D 107 , 082001 (2023), arXiv:2212.05967 [gr-qc].\n- [81] M. Tinto, F. B. Estabrook, and J. W. Armstrong, Phys. Rev. D 65 , 082003 (2002).\n- [82] P.-P. Wang, W.-L. Qian, Y.-J. Tan, H.-Z. Wu, and C.-G. Shao, Phys. Rev. D 106 , 024003 (2022), arXiv:2205.08709 [gr-qc].\n- [83] M. Otto, Time-Delay Interferometry Simulations for the Laser Interferometer Space Antenna , Ph.D. thesis, Leibniz U., Hannover (2015).\n- [84] T. A. Prince, M. Tinto, S. L. Larson, and J. W. Armstrong, Phys. Rev. D 66 , 122002 (2002), arXiv:grqc/0209039.\n- [85] Z.-C. Liang, Z.-Y. Li, J. Cheng, E.-K. Li, J.-d. Zhang, and Y.-M. Hu, Phys. Rev. D 107 , 083033 (2023), arXiv:2212.02852 [astro-ph.GA].\n- [86] Z.-C. Liang, Z.-Y. Li, E.-K. Li, J.-d. Zhang, and Y.-M.  \n- Hu, Phys. Rev. D 110 , 043031 (2024), arXiv:2307.01541 [gr-qc].\n- [87] X.-C. Hu, X.-H. Li, Y. Wang, W.-F. Feng, M.-Y. Zhou, Y.-M. Hu, S.-C. Hu, J.-W. Mei, and C.-G. Shao, Class. Quant. Grav. 35 , 095008 (2018), arXiv:1803.03368 [grqc].\n- [88] L. Tsukada, S. Jaraba, D. Agarwal, and E. Floden, Phys. Rev. D 107 , 023024 (2023), arXiv:2208.14421 [astroph.IM].\n- [89] Z.-C. Liang, Z.-Y. Li, E.-K. Li, J.-d. Zhang, and Y.-M. Hu, Results Phys. 63 , 107876 (2024), arXiv:2403.18709 [gr-qc].\n- [90] K. M. G'orski, E. Hivon, A. J. Banday, B. D. Wandelt, F. K. Hansen, M. Reinecke, and M. Bartelman, Astrophys. J. 622 , 759 (2005), arXiv:astro-ph/0409513.\n- [91] D. R. Kincaid and E. W. Cheney, Numerical analysis: mathematics of scientific computing , Vol. 2 (American Mathematical Soc., 2009).\n- [92] Z.-C. Liang, Z.-Y. Li, E.-K. Li, J.-d. Zhang, and Y.-M. Hu, (2024), arXiv:2409.00778 [gr-qc].\n- [93] A. Ain, P. Dalvi, and S. Mitra, Phys. Rev. D 92 , 022003 (2015), arXiv:1504.01714 [gr-qc].\n- [94] A. Ain, J. Suresh, and S. Mitra, Phys. Rev. D 98 , 024001 (2018), arXiv:1803.08285 [gr-qc].\n- [95] G. B. Arfken, H. J. Weber, and F. E. Harris, in Mathematical Methods for Physicists (Seventh Edition) , edited by G. B. Arfken, H. J. Weber, and F. E. Harris (Academic Press, Boston, 2013) seventh edition ed., pp. 773814.\n- [96] M. L. Men'endez, J. Pardo, L. Pardo, and M. Pardo, Journal of the Franklin Institute 334 , 307 (1997).\n- [97] M. Muratore, J. Gair, and L. Speri, Phys. Rev. D 109 , 042001 (2024), arXiv:2308.01056 [gr-qc].\n- [98] Y. Hu, P.-P. Wang, Y.-J. Tan, and C.-G. Shao, Phys. Rev. D 107 , 024026 (2023)."}