# Screener3D: a Gaseous Time Projection Chamber for Ultra-low Radioactive
Material Screening

Haiyan Du  INPAC; Shanghai Laboratory for Particle Physics and Cosmology; Key
Laboratory for Particle Astrophysics and Cosmology (MOE),
School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai
200240, China  Chengbo Du  Yalong River Hydropower Development Company, Ltd.,
288 Shuanglin Road, Chengdu 610051, China  Karl Giboni  Ke Han
<EMAIL_ADDRESS>INPAC; Shanghai Laboratory for Particle Physics and
Cosmology; Key Laboratory for Particle Astrophysics and Cosmology (MOE),
School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai
200240, China  Shengming He  Liqiang Liu  Yalong River Hydropower Development
Company, Ltd., 288 Shuanglin Road, Chengdu 610051, China  Yue Meng  INPAC;
Shanghai Laboratory for Particle Physics and Cosmology; Key Laboratory for
Particle Astrophysics and Cosmology (MOE),
School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai
200240, China  Shaobo Wang  INPAC; Shanghai Laboratory for Particle Physics
and Cosmology; Key Laboratory for Particle Astrophysics and Cosmology (MOE),
School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai
200240, China  SPEIT (SJTU-ParisTech Elite Institute of Technology), Shanghai
Jiao Tong University, Shanghai, 200240, China  Tao Zhang  Li Zhao  INPAC;
Shanghai Laboratory for Particle Physics and Cosmology; Key Laboratory for
Particle Astrophysics and Cosmology (MOE),
School of Physics and Astronomy, Shanghai Jiao Tong University, Shanghai
200240, China  Jifang Zhou  Yalong River Hydropower Development Company, Ltd.,
288 Shuanglin Road, Chengdu 610051, China

###### Abstract

In experiments searching for rare signals, background events from the detector
itself are some of the major factors limiting search sensitivity. Screening
for ultra-low radioactive detector materials is becoming ever more essential.
We propose to develop a gaseous time projection chamber (TPC) with a
Micromegas readout for radio screening. The TPC records three-dimensional
trajectories of charged particles emitted from a flat sample placed in the
active volume of the detector. The detector can distinguish the origin of an
event and identify the particle types with information from trajectories,
which significantly increases the screening sensitivity. For $\alpha$
particles from the sample surface, we observe that our proposed detector can
reach a sensitivity higher than 100 $\mu$Bq$\cdot$m-2 within two days.

Charged Particle detector, surface $\alpha$ measurement, Ultra-low
Radioactive, Material Screening

##  I Introduction

Experiments searching for neutrinoless double beta decay DellOro:2016tmg  and
dark matter direct-detection signals Roszkowski:2017nbc ; Liu:2017drf  have
extremely low signal rates, if any at all. The background rate recorded by the
detector is essential to the search sensitivity. In underground laboratories
such as the China Jinping Underground Laboratory (CJPL) Cheng:2018lcf , the
influence of cosmic rays is reduced by orders of magnitude and is no longer
dominant. Background events result almost entirely from the detector itself
and the surrounding lab environment. Therefore, extensive material screening
campaigns to select materials with low radioactivities have been conducted by
many experiments (for example, Wang:2016eud ; Jiang:2018pic ; Abgrall:2016cct
; Akerib:2014rda ; Leonard:2007uv ; Alessandria:2011vj ). The majority of
efforts emphasize the bulk radioactivities of materials using screening
techniques such as $\gamma$-ray measurement with high purity germanium (HPGe)
detectors, inductively coupled plasma mass spectrometry (ICP-MS), and neutron
activation analysis (NAA). Note that these techniques have significantly wider
applications than only screening low-radioactive materials (for example,
HPGe_application ; ICP-MS_application ).

Surface radioactivity often differs from those in the bulk and may
significantly contribute to the background budget. The surface of detector
components may become contaminated with additional radioactivity through
machining, handling, or exposure to air. When directly facing a detector
active volume, $\alpha$ or $\beta$ particles emitted from the surface may
directly introduce background events. The surface background is particularly
critical for modular detectors where large surface areas face detector
modules. For example, $\alpha$ events from the supporting structure and
surface of a detector are the most prominent source in the Cryogenic
Underground Observatory for Rare Events (CUORE) bolometer array
Alduino:2016vtd . Even if surface contaminations are far from detectors,
radioactive impurities may emanate from the surface and travel to the detector
via the circulation of the liquid or gas detector medium. For example, in the
PandaX-4T experiment, 222Rn from the detector and inner surface of the
circulation pipe is a major background source in the search for dark matter
Zhang:2018xdp .

Surface radioactivity is often measured using semiconductor detectors,
scintillator arrays, gas ionization counters, etc. The detection areas of
commercial silicon semiconductor detectors are frequently smaller than 30 cm2
(e.g. OrtecWebsite ). An example of a scintillator array is BiPo-3, which
utilizes 40 scintillator modules for a total screening area of 3.6 m2 BiPo-3 .
The background rate achieved by the BiPo-3 detector is 0.9 $\mu$Bq$\cdot$m-2
for 208Ti and 1 $\mu$Bq$\cdot$m-2 for 214Bi BiPo-3 ; BiPo-3_background . The
extremely sensitive detector requires dedicated support; therefore, the cost
for measurement is high. Gas ionization counters can have large active areas
and particle identification capability using pulse shape analysis. Ultralo
1800, offered by XIA XIAWebsite , is one such detector with an $\alpha$
background of approximately 280 $\mu$Bq$\cdot$m-2. However, it fully utilize
the three-dimensional track recording capabilities of gas detectors. Proposals
such as BetaCage BetaCage  and $\mu$ time projection chamber ($\mu$-TPC) upic
are no longer active or are in the early stage of development. Our proposal
differs from $\mu$-TPC in terms of sample placement, readout modules, and
operating gas medium.

We propose to construct a large-area, high-efficiency, and low-cost gaseous
TPC to measure surface contaminations, which we call Screener3D. Screener3D
will be able to measure the energy and trajectories of $\alpha$ and $\beta$
particles with position-sensitive readout modules of the Micro-MEsh Gaseous
Structure (Micromegas MM ). In a typical argon TPC operating at atmospheric
pressure, $\alpha$ particles of 10 MeV travel approximately 10 cm in an almost
straight line. Trajectories of $\beta$ particles of the same energy are
significantly longer and meander. Along the particle trajectories, the gas
medium is ionized. The number and initial position of ionization electrons
carry the three-dimensional trajectory information of the original particle.
When combining energy and trajectory information, we can better distinguish
signals from the background and identify the origin of the signals. Samples
are placed in the TPC and the emitted charged particles are recorded with high
efficiency. For $\alpha$ particles, in particular, the energy can be measured
more accurately without energy loss in entrance windows. In this paper, we
provide an overview of the proposed design and sensitivity studies using
simulations.

##  II Design overview

The central part of Screener3D is a low-background, large-area, and high-
granularity gaseous TPC. This type of gaseous TPC is extensively used in
nuclear and particle physics (e.g. CAT-TPC ; PandaX-III ) All the materials
used to build the TPC will have low radioactivities to control the background
contributions from the detector itself. The designed readout area, which
determines the largest sample area that can be measured, is approximately 2000
cm2. The Micromegas readout will be based on 3 mm wide strips, a granularity
sufficient for $\alpha$ and $\beta$ tracks in the order of 10 cm. In addition
to the key detector performance specifications, Screener3D must be stable over
an extended period. The long-term stability is maintained with a gas
circulation and purification system and a real-time monitoring system. With
the auxiliary systems, we can measure samples for a prolonged period of up to
weeks to increase the measurement sensitivity. We describe the main components
of Screener3D and material screening plans in this section.

###  II.1 Gaseous TPC with Micromegas

Figure 1: Schematic design of Screener3D with the main components labeled.
(Top) The inside view with field cage drawer pulled out. Samples are placed at
the bottom of the drawer (on the cathode). (Bottom) The detector with front
cap closed and ready for measurement.

The Screener3D TPC primarily comprises a readout plane with Micromegas
modules, a cathode, a field cage, a gas medium, and an outer vessel (Fig. 1).
The rectangular readout plane on the top collects drift electrons after
amplification via avalanche. The cathode is at the bottom and provides a
negative high voltage. The field cage connects the readout plane and cathode
mechanically and aids in maintaining a uniform drift electric field in the
active volume, which is the cuboid space enclosed by the aforementioned
components. The active volume is 60 cm long, 40 cm wide, and 10 cm high, which
is filled with 1 bar of an argon and isobutane gas mixture during operation.
All the components are fixed inside the vessel with signal feedthroughs and
various ports for gas circulation and pumping.

During measurement, a thin slab sample, such as a sheet of metal, plastic
film, or silicon wafer, would be placed on the cathode plane. Radioactive
contaminations on the top surface of the sample emit $\alpha$ and $\beta$
particles into the active volume. We focus primarily on the identification of
$\alpha$ particles from samples since the energy deposition in the unit
distance is large and trajectory characteristics are easier to define. The
typical $\alpha$ emissions emanate from uranium and thorium decay chains.
Nuclides such as 232Th, 238U, 214Po, and 212Po emit $\alpha$ particles with
energies of 4.0, 4.1, 7.7, and 8.8 MeV, respectively. A few other isotopes
along the decay chains emit $\alpha$ particles in the energy range of 4 to 7
MeV.

The readout plane may consist of multiple Micromegas modules. In our current
design, six modules of approximately $20\times 20\,\rm{cm}^{2}$ are tiled
together to form a total readout area of 2400 cm2. In Fig. 1, the modules,
supported by an acrylic frame, are shown face-down. For each module, the
active area is split into $64\times 64$ diamond-shaped pads, each of which has
a diagonal length of 3 mm. The pads are inter-connected in the vertical or
horizontal direction and form the so-called X or Y readout strips. With 64 X
and 64 Y strips, the total number of electronics channels is 128 per module.
Flat signal cables fabricated from a Kapton-based flexible printed circuit
board (PCB) pass through the feedthroughs on the top of the vessel to send the
signals to readout electronics (not shown in the figure) outside the TPC.
Micromegas modules are frequently biased at a few hundreds of volts and the
voltage is also provided through the flat Kapton cable.

Fig. 1 (top) also shows the field cage and cathode as a drawer pulled out from
the acrylic frame. The sidewall of the field cage is created from 2 cm-thick
acrylic plates with a Kapton-based flexible PCB attached to the outside. There
are uniformly spaced copper strips on the PCB, and the adjacent strips are
connected with surface mount resistors to maintain a uniform electric field.
At the bottom of the drawer, a polished copper plate with a thickness of 0.2
cm is used as the cathode. The cathode is electrically connected to the bottom
strip of the field cage PCB. A uniform electric field is formed within the
field cage when a negative high voltage is provided to the cathode. The entire
setup is enclosed in an outer vessel created from 1 cm-thick stainless steel
(SS). For easy sample loading with the drawer-style field cage, the opening of
the vessel is at the front (Fig. 1 (bottom)). Ports for gas lines, a vacuum
pump, and feedthrough for the cathode high-voltage cable are added on the side
of the vessel.

Samples are placed on the cathode for measurement. In our current design shown
in Fig. 1, the front panel of the vessel is unmounted, the field cage is
pulled out, and the sample can be placed on the cathode. When closed again, we
pump down the vessel to $\mathcal{O}(1)$ mPa before flushing argon gas in. We
expect the entire sample loading process to be accomplished in a period of an
hour. Alternative designs, such as a stationary field cage with a slit opening
on the front, are also under consideration to further streamline the sample
loading process.

When samples are placed in the active volume, the counting efficiency is high
and the measurement of $\alpha$-particle energy is more accurate. An aspect of
concern may be the disturbance of electric potential lines when samples are on
the cathode. We simulated the effect of a thin Kapton sample using the finite
element analysis software COMSOL Multiphysics COMSOLWebsite  (Fig. 2). The
black wireframes in Fig. 2 represent the internal structure of the TPC, such
as the field cage, cathode, and readout modules. We supplied the cathode with
a 1 kV negative voltage and the readout plane with a 250 V negative voltage.
The Kapton sample 54.4 cm long, 34.4 cm wide, and 0.1 cm thick is shown as the
innermost rectangle in Fig. 2 (a). Fig. 2 (b) shows the side view of field
cage. The field lines inside the field cage exhibited minimal distortions. In
the X–Z plane , the electric field in the X-direction was 3.6% of that in the
vertical Z-direction at 2 mm from the edge of the sample, while the value was
less than 1% at 2 cm from the edge of the sample. The uniformity of the
electric field was almost the same for a 0.1 cm-thick copper or silicon sample
of the same size; thus, it was negligible for our measurement.

The distortion of the electric field at the edge of the sample became more
severe with increased sample thickness. For a 1 cm-thick copper sample, the
distortion became smaller than 1% when measured at approximately 7 cm from the
edge. Therefore, we can mitigate this problem using more aggressive fiducial
cuts. However, we should note that the main objective of Screener3D is to
measure thin millimeter-scale samples.

(a) Top view of the detector showing potential on the X–Y plane (b) Side view
of the detector showing potential on the X–Z plane

Figure 2: Electric field simulation using COMSOL, which demonstrated minimal
impact to the electric field of samples on the cathode. (a) Top view of the
detector showing the potential on the X–Y plane. The innermost black rectangle
shows the boundary of the sample, which had a slightly smaller footprint than
the active volume of the TPC. (b) Side view of the detector showing the
potential on the X–Z plane. Red lines with arrows indicate the electric field
lines. We can observe that the field lines were in almost perfect vertical
directions above the sample’s footprint.

###  II.2 Gas system and slow monitoring system

The energy response and resolution of gaseous TPC are highly sensitive to the
purity of gas medium, particularly electronegative impurities such as oxygen
and water. The impurities may attract and absorb drift electrons in the active
volume. They may also reduce the gain in the avalanche gap of the Micromegas.
For the stability of the prolonged measurement period of Screener3D, a gas
circulation and purification system (shown in Fig. 3) are designed. Before
flowing into the TPC, gas mixtures travels through a set of gas purifiers.
Given the different precision requirements of measurement, the purifier may be
a charcoal-based absorber (primarily to absorb radio-impurities such as
radon), a chemical-reaction-based getter (primarily to absorb electronegative
impurities), or a combination of both. A circulation loop is added for the
continuous purification of the working gas for an extended measurement of low-
radioactivity samples.

To quantify the long-term (in)stability of the detector, we monitor key
operation parameters, such as ambient temperature, detector temperature, gas
pressure, Micromegas bias voltage, leakage current, TPC drift voltage, and
event rates . We design and build a centralized slow monitoring system
slow_control  (SMS) with all the parameter values stored in a database and are
accessible in real-time on a web browser. With pre-defined safety ranges of
the key parameters, the SMS can also send alerts via email or phone text
messages. The SMS data would also enable us to reject or correct poor-quality
data based on the variation of operating conditions. The SMS system may also
store sample information, such as area, material, and screening time to serve
as a centralized catalog of measurements.

Figure 3: Schematic diagram of the Screener3D gas system. Red arrows denote
the filling loop. Green arrows denote the circulation loop. Blue arrows are
the loops shared by the filling and circulation. After the pressure-reducing
valve, the working pressure is limited to 1 bar.

###  II.3 Intrinsic background sources and mitigation

To improve the signal-to-background ratio in screening, we must minimize the
radioactive background from Screener3D itself. The inner part of the
Screener3D TPC will be built with materials of low radioactivity, screened
using HPGe and/or mass spectrometer techniques. The cathode will be fabricated
from high-purity oxygen-free copper particularly screened for low uranium and
thorium contaminations. The field cage is designed to have only low-
radioactivity acrylic facing the active volume to avoid background events from
the PCB and resistors. Surface contaminations of the Micromegas are more
challenging to measure and control. In this study, we used microbulk
Micromegas MM , which is created from Kapton and copper using the PCB
lithography technique and has been proven to have good performance and low
surface radioactivityMMbackground . Furthermore, the surfaces of detector
components, particularly those facing the active volume of the TPC, will be
thoroughly cleaned to reduce the emission of $\alpha$ particles. The inner
volume of the Screener3D will be always filled or flushed with argon gas to
avoid radon contamination from the air. In our sensitivity study later, we
assume all the detector components are created from low-radioactivity
materials and no secondary surface contaminations are introduced in
construction.

In addition to minimizing the radioactive emission from detector components,
we can further suppress the background with topological information from
particle tracks, which is the key feature of Screener3D. Tracks of $\alpha$
particles from a sample and various TPC components are shown in Fig. 4. The
high energy $\alpha$ tracks are mostly straight in 1 bar of argon. A large
blob at one end of the track highlights the $dE/dx$ increase when $\alpha$
particles stop. $\alpha$ events originating from the sample would have
starting points from the sample surface and an upward-going direction in our
detector geometry. For the contaminations of the field cage, $\alpha$ tracks
enter the active volume from the side and may travel upward or downward.
Surface $\alpha$ events from the readout plane point downward and can be
easily distinguished from those of the samples. During sample measurement, the
majority of the surface area of the cathode is covered by a sample and no
$\alpha$ particle from the cathode could enter the active volume excepted the
very edge (Fig. 4). With the distinguished track features, we can
significantly suppress $\alpha$ backgrounds from the detector and achieve a
higher sensitivity. In the following section, we demonstrate the suppression
power using simulation data.

Figure 4: Illustration of $\alpha$ particle trajectories in a gaseous TPC.
$\alpha$ particles travel mostly in a straight line in 1 bar of argon and 5%
isobutane gas with a prominent Bragg peak at the end. We utilize this feature
for background suppression in our detector.

##  III Expected screening sensitivity

The screening sensitivity of Screener3D depends on the background event rates
recorded by the detector. In this section, we describe the setting up of a
detector model in the simulation framework Geant4 Geant4  and major background
sources for $\alpha$ particle measurement. The energy deposition process and
gaseous detector response were implemented in the simulation. Signal selection
and background suppression efficiencies were first calculated using energy
cuts and particle identification. The unique tracking capability of the
gaseous detector was then exploited to distinguish $\alpha$ particles from
different origins and with different orientations to further reduce the
$\alpha$ background rate. We finally present the material screening
sensitivity for surface contamination of $\alpha$-emitting impurities.

###  III.1 Screener3D background

Figure 5: Geometry of main components of Screener3D as constructed in the Geant4 simulation. All the dimensions followed the schematic design as shown in Fig. 1. Table 1: Radioactivity level of bulk material in Screener3D  Material | Radioactivity
---|---
|  238U |  232Th |  40K |  60Co |  137Cs |  222Rn |  39Ar
| (mBq$\cdot$kg-1) | (mBq$\cdot$kg-1) | (mBq$\cdot$kg-1) | (mBq$\cdot$kg-1) | (mBq$\cdot$kg-1) | ($\mu$Bq$\cdot$m-3) | (mBq$\cdot$kg-1)
Gaseous argon | 0.0018 | 0.0004 |  |  |  | 10 | 103
Acrylic | 0.003 | 0.003 |  |  |  |  |
Oxygen-free copper | 0.38 | 0.51 | 4.00 | 0.20 | 0.16 |  |
Stainless | 1.70 | 2.74 | 13.95 | 1.03 | 2.36 |  |

Fig. 5 shows the detector geometry constructed in the Geant4 simulation Geant4
. The setup consisted of the sensitive volume, field cage, readout plane,
cathode, and SS vessel with dimensions identical to the detector’s conceptual
design. The sensitive volume in the center was a $60\times 40\times 10$ cm3
cuboid with 1 bar argon and 5% isobutane mixture. The readout plane is
expected to be instrumented using microbulk Micromegas, which has complicated
mechanical structures, such as a lithographed copper mesh, avalanche holes,
and underlining PCB patterns MM . The structure was substituted with a 0.1 mm-
thick copper. The readout plane on the top of the sensitive volume consisted
of 6 Micromegas modules and was 58.0 cm long and 38.5 cm wide. The cathode,
created from 2 mm-thick oxygen-free copper, had the same surface area and
faced the readout plane beneath the sensitive volume. A 4 cm-thick acrylic
frame surrounding the sensitive volume from four sides represented the
mechanical structure of the field cage and acrylic frame. We placed a
55$\times$35$\times$0.01 cm3 sample on the cathode in the simulation. All the
components described above were enclosed in a 1 cm-thick SS vessel with an
inner volume of $80\times 60\times 15$ cm3.

We simulated major background contributions from surface and bulk containments
of all components. The energy deposition process of $\alpha$, $\beta$, and
$\gamma$ particles in the TPC were simulated in Geant4. For the study of
$\alpha$ contamination on the sample surface, we expected the $\beta$ and
$\gamma$ background events from detector bulk material to have a marginal
impact. We confirmed this in our simulation first before focusing on the
$\alpha$ background from surfaces of detector components.

Major bulk background sources such as 238U, 232Th, and other radio-
contaminations in different materials are listed in Table 1. The radio-
contamination of oxygen-free copper and SS followed measurements of the PandaX
Collaboration Zhang:2018xdp . For acrylic, the radioactivities were assumed to
be the same as that used by Reference Juno_acrylic . The working gas medium
for $\alpha$ contamination measurement was argon. The radio-purity of
commercial argon gas varies significantly from different vendors or even
different batches of the same vendor. For 238U and 232Th in argon, we adopted
values from GERDA BUDJAS2009755 . Moreover, atmospheric argon contains
approximately 1 Bq$\cdot$kg-1 of cosmogenic 39Ar Agnes:2015ftt , which is a
concern for $\beta$ measurement but will not introduce background in the high-
energy region for $\alpha$ measurement.

222Rn in argon is listed separately from the parent 238U chain. Radon is a
noble gas and is frequently at the $\mu$Bq-mBq$\cdot$m-3 level in commercial
argon gas supply. The contamination can be reduced by at least three orders of
magnitude after passing through an activated carbon cold trap Simgen . We
assumed a concentration of 222Rn of 10 $\mu$Bq$\cdot$m-3 in our simulations.

Figure 6: Low-energy spectra of Screener3D. Main contributions in the energy
range were from $\gamma$ and $\beta$ events.

Fig. 6 shows the energy deposition of $\beta$ and $\gamma$ events in
Screener3D. Owing to the relatively small $dE/dx$ of $\beta$ and $\gamma$
particles, the typical energy deposition is less than 1 MeV, which is
distinctively different from that of $\alpha$ particles. Therefore, we defined
the region of interest (ROI) for $\alpha$ particle measurement as above 1 MeV.
The upper limit for $\alpha$ ROI was 10 MeV, which is above the highest
$\alpha$ energy of naturally occurring radioactivity of the decay chains of
uranium and thorium.

Only the $\alpha$-emitting surface directly facing the sensitive volume,
including argon gas, readout plane, field cage, and cathode contribute to the
$\alpha$ background. $\alpha$ particles with energy below 10 MeV travel less
than 20 $\mu$m in copper and less than 100 $\mu$m in acrylic. Therefore,
$\alpha$ contamination near the surface were our main focus. The concentration
of $\alpha$-emitting radio-isotopes, particularly decay chains of 238U and
232Th, on the surface is frequently higher than that in the bulk, owing to
secondary contamination. The surface radioactivities of acrylic, oxygen-free
copper, and Micromegas readout plane are shown in Table 2. We deduced a
conservative surface contamination level of acrylic by comparing two sets of
acrylic contamination levels with and without the surface cleaning from the
JUNO Collaboration Juno_acrylic . For copper, we adopted the maximum values
measured by the CUORE Collaboration in dedicated bolometer arrays CUORE_copper
. The radioactivity levels of the readout plane followed the results of
microbulk Micromegas measured using the BiPo-3 detector MMbackground ; PandaX-
III .

The depth profile of surface contaminations can be approximated using an
exponential curve,

$\rho={\rho_{0}e^{-\frac{x}{\lambda}}}+b,$ (1)

where $\rho$0 is the surface contamination per unit volume at the very
surface, $x$ is the distance into the bulk material, $\lambda$ is the
characteristic penetration depth of the contamination, and $b$ is the
contamination level in the bulk material. At the very surface, the
contamination level is $\rho_{0}+b$ CUORE_copper . $\lambda$ may vary
significantly with the contamination type and causes of contamination. The
depth profile can be verified and the exact value of $\lambda$ determined by
measuring the peak shape of $\alpha$ spectrum, but they are challenging
experimentally. In our simulation, we set $\lambda=0.1\,\mu$m, which is a
shallow distribution, but the effects can still be experimentally identified.
Additionally, we compared our results with $\lambda=1\,\mu$m, to study the
influence of the penetration depth of the contamination. For both $\lambda$
values, the total amount of the surface contaminations in the materials was
the same.

Table 2: Radioactivity level of surface in Screener3D Material | Radioactivity(mBq$\cdot$m-2)
---|---
|  238U |  232Th
Acrylic | 0.25 | 0.08
Oxygen-free copper | 1.30 | 1.30
Readout plane | 0.45 | 0.14

Fig. 7 shows the energy spectra of a copper sample on the cathode plane and
different detector components for both $\lambda$ values. The count rates of
the sample, field cage, readout plane, and cathode were scaled with the
surface area and the corresponding surface radioactivity level listed in Table
2. For the gas medium, the spectrum was determined by the contamination level
of radon and total argon mass. The bulk radioactivity level used in the
simulation is listed in Table 1.

In the ROI between 1 and 10 MeV, the total expected background rate was 82.4
(76.1) counts per day for $\lambda=0.1\,\mu$m ($1\,\mu$m). Approximately 46.0%
(45.4%) of the background was from the readout plane, 36.1% (35.8%) from the
cathode, 17.6% (18.5%) from the field cage, and 0.3% (0.3%) from the working
gas. The difference in $\alpha$ peak shapes of Fig. 7 (a) and Fig. 7 (b)
highlight the impact of depth profile. For relatively radio-active samples
with high contamination levels, we can potentially have a quantitative
estimation of the depth profile. In the two figures, the spectra from argon
gas are identical.

(a) Geant4-simulated energy spectra assuming $\lambda=0.1\,\mu$m (b)
Geant4-simulated energy spectra assuming $\lambda=1\,\mu$m

Figure 7: Geant4-simulated energy spectra of Screener3D. For a longer
penetration depth of $\lambda=1\,\mu$m, $\alpha$ particles are more likely to
lose partial energy before escaping the surface. Therefore, the tail of the
characteristic peak is more apparent in figure (b).

###  III.2 Detector response simulation

Detector response, including electron diffusion, energy resolution, readout
scheme, and electronics response, was added to the Geant4 simulation data to
produce mock detector data for downstream analysis. The detector response was
simulated using the REST framework REST , which was developed and used by
PandaX-III PandaX-III  and other gas TPC projects. If we use an $\alpha$
particle traveling in gas as an example, Geant4 recorded the differential
energy deposition $dE/dx$ along the trajectory. In REST, the differential
energy deposition was converted to the number of ionization electrons in gas.
The ionization electrons diffused while drifting to the readout plane. REST
interfaced with the Garfield program Garfield  to calculate diffusion
coefficients. For argon/isobutane(5%) mixture at 1 bar and under an electric
field of 100 V/cm, the transverse (longitudinal) diffusion coefficient was
0.048 (0.044) $\sqrt{\text{cm}}$. The diffusion of the electrons was then
calculated as $D_{ic}\sqrt{L_{d}}$, where $D_{ic}$ is the diffusion
coefficient and $L_{d}$ is the drifting distance. When the diffused electrons
reached the readout plane, the position and arrival time were recorded.
Following the strip-readout schemes of Micromegas in the actual detectors,
ionization electrons were grouped by strips, and pulses were generated strip
by strip (Fig. 8 (a)). The colored pulses represent signals from different
strips. The pulse amplitude in Fig. 8 is proportional to the number of
electrons collected per strip. The amplitudes were smeared according to a
Gaussian function to account for an energy resolution of 3% full width at half
maximum (FWHM) at 2.5 MeV in the detector response simulation. Pulse widths
were determined from electronics shaping and a 1 $\mu$s shaping time was used
(Fig. 8 (a)). The relative timing among pulses indicated the arrival time of
electrons. Under a relative constant drift velocity of approximately 30
mm/$\mu$s, the arrival time denoted the relative energy deposition position in
the drift direction. For a sampling rate of 50 MHz and a record length of 512
sample points per signal, we recorded approximately 10 $\mu$s of data per
signal window. For a drift length of 10 cm in our TPC, the signal window was
more than adequate, as shown in the figure. We also implemented triggers in
REST, but the trigger threshold effect was negligible for $\alpha$ events in
the ROI.

For comparison, Fig. 8 (b) shows an $\alpha$ event recorded by a gaseous
detector. The detector was equipped with a $20\times 20$ cm2 Micromegas and
filled with 1 bar of argon and 5% isobutane gas mixture. The height of the
drift volume was 10 cm, the same as in our simulation. Our mock pulses re-
produced the key features of detector pulses. Next, the mock pulses of each
simulated event were analyzed using the same procedure as detector data.

(a) Simulated waveforms recorded by Micromegas (b) Waveforms of an $\alpha$
event in a TPC measured with a microbulk Micromegas readout plane.

Figure 8: Simulated (a) and measured (b) waveforms of an $\alpha$ event. The
measurement was performed with 1 bar of argon and 5% isobutane gas mixture.

Energy was reconstructed using mock pulses and the resulting spectra are shown
in Fig. 9. Compared with Fig. 7, the $\alpha$ peaks were broadened owing to
the detector resolution. Since the effective area of the readout plane did not
cover the full area of the field cage, partial energy of background events
from the field cage were not be recorded. Therefore, the spectrum of
background events from the field cage shifted to a lower energy compared with
that in Fig. 7. Meanwhile, a small peak was observed in front of some
characteristic peaks from the readout plane.

(a) Simulated energy spectra with detector response assuming
$\lambda=0.1\,\mu$m (b) Simulated energy spectra with detector response
assuming $\lambda=1\,\mu$m

Figure 9: Simulated energy spectra of Screener3D with detector response added.

###  III.3 Track reconstruction

The main application of mock pulses is track reconstruction. As mentioned
earlier, the amplitude and timing of pulses of different strips contain all
the information we would collect from an actual detector. Therefore, this
information can be used to reconstruct the particle tracks and differential
energy loss along the tracks in a TPC. Track reconstruction was performed in
the X–Z and Y–Z planes, where Z denotes the drift direction and X and Y
represent transverse directions. In 1 bar argon, $\alpha$ particles would
rarely be scattered by a large angle and the tracks are mostly straight. Fig.
10 (a) shows an $\alpha$ track in the X–Z and Y–Z planes. The red trace
represents the trajectory and yellow dots represent the energy deposition
vertex. The size of yellow dots illustrates the relative amount of energy
deposition. The Bragg peak is prominent at the end of the trajectory on the
top left of the figures.

In Fig. 10 (b), we show the reconstructed track in the X–Z and Y–Z planes of
the same event as in Fig. 10 (a). Each red (or green) point in the figures
represents a triggered strip signal, and the size of the point represents the
amount of deposited energy. The X (Y) coordinate was determined from the
position of the strip. The Z-axis indicates the position of the ionized
electron relative to that of the first electron to reach the readout plane, as
determined using the drift velocity and pulse timing. For long, near-straight
$\alpha$ tracks, sorting by timing or X/Y position can reliably reconstruct
the true tracks. The black lines in Fig. 10 (b) represent a well-reconstructed
track. Track reconstruction can be challenging when the tracks are short. For
example, the length of the X–Z and/or Y–Z tracks may be very short because of
the projection angle. The tracks may also be short when $\alpha$ particles
only deposit partial energy in the active volume.

(a) Geant4-simulated tracks of an $\alpha$ particle (b) Hits recorded using
Micromegas strips and the reconstructed track

Figure 10: Geant4-simulated (a) and reconstructed (b) tracks of an $\alpha$
particle. In (a), the straight red line represents the $\alpha$ particle track
and yellow dots represent the relative size of energy deposition. In (b), each
red/green dot represents a triggered strip signal and the size of the dot
represents the amount of deposited energy in the strip. The black line is the
reconstructed track by connecting nearby hits.

###  III.4 Screener3D background suppression

In our ROI from 1 to 10 MeV for $\alpha$ measurement, effectively all the
$\beta/\gamma$ background were removed. The ROI cut also removed 11% of
$\alpha$ background and maintained 98.9% of the signal events. Alpha
background could be further suppressed using fiducial, angle, and hit-number
cuts. We explain the cuts in detail as follows.

We set the cut values to maximize the detector sensitivity $S_{d}$. For
measurements with a large number of backgrounds $B$, background fluctuation
follows a Poisson distribution and equals $\sqrt{B}$. Therefore, the detector
sensitivity is proportional to $\epsilon_{s}/\sqrt{\epsilon_{b}}$, where
$\epsilon_{s}$ is the signal detecting efficiency and $\epsilon_{b}$ is the
efficiency in which background events are maintained through selection cuts.
For each cut, we selected the cut values that maximized the
$\epsilon_{s}/\sqrt{\epsilon_{b}}$.

Figure 11: Distribution of tracks starting point. We created a cut at the
fiducial distance of 27 mm to effectively remove the background events from
the field cage and cathode. Figure 12: Distribution of track direction,
defined as the cosine of the angle between the track and vertical direction.
The distribution is for events that survived the fiducial cut. We created a
cut at 0 to effectively remove the background events from the readout plane
and gas. Figure 13: Distribution of number of hits of an event. We created a
cut at 10 (9 when $\lambda=1\,\mu$m) to remove the events with short tracks.

We devised two competing algorithms to identify the Bragg peak of a track and
thus the track’s origin and direction. A track in the X–Z or Y–Z plane was
split into two segments from the mid-point, and the energy in each segment was
calculated. The Bragg peak was more likely to be in the segment with larger
energy. The energy of each hit point and $dE/dx$ along the vertical direction
could also be used to locate Bragg peaks. When particle energy was higher than
1 MeV and the total number of hits was larger than 10, the identified
efficiency of the track’s direction was approximately 97%.

We identified the starting point of each track, and the distribution is shown
in Fig. 11. The X-axis in the figure is called the fiducial distance, which is
defined as the distance in the inward direction from the inner surface of the
field cage. The majority of the events from the field cage and cathode
originated from the edge of the sensitive volume. For those events, the
fiducial distance was less than 50 mm. We performed a fiducial cut at the
fiducial distance of 27 mm and discarded any events with a smaller fiducial
distance, i.e., closer to the edge of the active volume. After the fiducial
cut, the signal efficiency was at 96.2%. Less than 3% of the events from field
cage and cathode background remained. The cut was less effective for other
background sources, and the efficiencies are listed in Table 3.

Fig. 12 shows the distribution of $\cos(\theta)$, where $\theta$ is the angle
between the particle traveling and upward vertical directions. The
distribution is for events after the fiducial cut. The majority of $\alpha$
events from the sample had positive $\cos(\theta)$ values while events from
the readout plane had negative values, as shown clearly in the figure. Events
from the sample with negative $\cos(\theta)$ values were mostly short events
with incorrect track reconstruction. The same was true for events from the
readout plane but with positive $\cos(\theta)$ values. Events from the cathode
and field cage that survived the fiducial cut had incorrectly reconstructed
starting points and possibly orientations. Therefore, the distribution of
these events is nearly symmetrical in the figure. We created a cut at 0 and
rejected events with $\cos(\theta)<0$.

We also counted the number of hit strips to further reject background events.
The distribution of the hits number after the fiducial cut and angle cut is
shown in Fig. 13. We maintained only events that triggered more than 10 (9
when $\lambda=1\,\mu$m) strips. The cut rejected background events from the
readout plane most effectively since the remaining events that survived
previous cuts were primarily particles with short tracks. The peak structure
at the number of hits equal to 10 for events from the cathode was an artifact
from the angle cut, where we had two algorithms to calculate the angles for
tracks longer or shorter than 10 with different efficiencies.

Table 3: The $\alpha$ background rate after each cut ($\lambda=0.1\,\mu$m) Cut parameter | Energy | Fiducial | Angle | Hit number
---|---|---|---|---
| cut(%) | cut(%) | cut(%) | cut(%)
Field cage | 83.9 | 2.9 | 0.7 | 0.5
Cathode | 79.5 | 2.7 | 1.5 | 1.2
Readout plane | 95.1 | 80.3 | 7.1 | 1.2
Working gas | 92.0 | 77.1 | 40.5 | 28.9
Sample | 98.9 | 96.2 | 92.2 | 74.2
Table 4: The $\alpha$ background rate after each cut ($\lambda=1\,\mu$m) Cut parameter | Energy | Fiducial | Angle | Hit number
---|---|---|---|---
| cut(%) | cut(%) | cut(%) | cut(%)
Field cage | 83.5 | 2.8 | 0.7 | 0.6
Cathode | 78.7 | 2.7 | 1.3 | 0.9
Readout plane | 94.0 | 79.3 | 10.1 | 1.2
Working gas | 92.0 | 77.1 | 40.5 | 31.3
Sample | 97.0 | 94.3 | 87.5 | 67.9

Fig. 14 shows the signal and background spectra after all cuts. The
backgrounds were suppressed by nearly two orders of magnitude while
maintaining a high signal efficiency. For $\lambda=0.1\,\mu$m, the combination
of all the cuts reduced the $\alpha$ background rate from 82.3 to 0.9 counts
per day while maintaining 74.2% of the signals. Among the remaining background
events, 49.7% emanated from the readout plane, 36.5% from the cathode, 7.1%
from the field cage, and 6.7% from the gas. Assuming $\lambda=1\,\mu$m, the
background was suppressed from 76.1 to 0.8 counts per day while 67.9% of the
signal was preserved. Approximately 52.5$\%$ of the remaining background was
from the readout plane, 28.7$\%$ from the cathode, 10.2% from the field cage,
and 8.6% from the gas. Owing to more $\alpha$ events generating below the
surface and thus losing partial energy before they were emitted from the
surface, the total residual events of both the background and signal were
slightly fewer than those of $\lambda=0.1\,\mu$m.

(a) Energy spectra after all the cuts assuming $\lambda=0.1\,\mu$m (b) Energy
spectra after all the cuts assuming $\lambda=1\,\mu$m

Figure 14: Energy spectra after all the cuts. The total backgrounds were
suppressed by nearly two orders of magnitude while nearly 70% signals were
preserved.

###  III.5 Screener3D measurement sensitivity

Table 5: Sensitivities of Screener3D for total $\alpha$ measurement (the results in parentheses from $\lambda=1\,\mu$m) Measurement time | Background events | Sensitivity
---|---|---
(day) | (counts) | ($\mu$Bq$\cdot$m-2)
2 | 1.86 (1.59 ) | 97 (99)
7 | 6.50 (5.56) | 46 (47)

Assuming no events are observed from the samples statistically different from
the fluctuation of background, we can evaluate the sensitivity of Screener3D
for a particular screening period T in days. With a background rate of
$R_{bkg}$ counts per day, the total background for a measurement is
$B=R_{bkg}T$. We use the Poisson fluctuation of small background values to
calculate the upper limit on counts and subsequent surface background rate.
For small $B$ and a confidence level of $\Gamma$, the upper limit for counts
of signal $\hat{S}$ can be expressed as

$\sum_{k=\hat{S}+B}^{\infty}p_{B}(k)=1-\Gamma,$ (2)

in which $\hat{S}$ can be evaluated numerically or with a lower incomplete
Gamma function. Finally, the sensitivity can be calculated using the signal
efficiency, $\hat{S}$, and sample surface areas. Table 5 shows the sensitivity
of Screener3D with typical measurement periods. We can reach a sensitivity of
97 $\mu$Bq$\cdot$m-2 at the 90%confidence level (C.L.) for a two-day
measurement, assuming a typical sample such as copper and an exponential
surface contamination profile with $\lambda=0.1\,\mu$m. The sensitivity would
be higher than 50 $\mu$Bq$\cdot$m-2 when the measurement duration is extended
to one week. For a characteristic depth $\lambda=1\,\mu$m, the results are
very similar (Table 5).

Table 6: Sensitivities of Screener3D for $\alpha$ measurement from polonium (the data in parentheses from $\lambda$=1 $\mu$m) Nuclide | Measurement time | Background events | Sensitivity
---|---|---|---
| (day) | (counts) | ($\mu$Bq$\cdot$m-2)
212Po | 2 | 0.06 (0.04) | 75 (75)
214Po | 2 | 0.18 (0.1) | 75 (75)

Screener3D can also measure $\alpha$ rates of a specific nuclide, such as
212Po of the thorium decay chain and 214Po of the uranium decay chain. The
signature $\alpha$ peak of 212Po is at 8.8 MeV and we define a specific ROI
between 8 and 9 MeV. In the ROI, 7% of background events and 98% of the
signals are reserved after the energy and topology cuts (fiducial and angle
cuts). The final background rate is 0.03 (0.02 when $\lambda=0.1\,\mu$m)
counts per day. For 214Po in the ROI of 7 and 8 MeV, 97% of the background is
removed while 98% of the signals remains. The final background rate is 0.09
(0.05 when $\lambda=0.1\,\mu$m) counts per day. For the near-zero background
measurements at specific peaks, the Eq. 2 is no longer applicable. We use the
Feldman-Cousins approach FC  with zero background as an approximation for this
calculation. The sensitivity is 75 $\mu$Bq$\cdot$m-2 for two-day measurements
of both polonium isotopes.

##  IV Summary and outlook

Background events from detector surfaces may have a severe negative impact on
experiments searching for rare events. Screening methods for surface
contaminations have not been utilized as widely as those for bulk
contaminations (such as gamma and mass spectrometers). In this paper, we
propose a charged-particle detector using the TPC concept for surface
contamination screening. The TPC design concept and detector background
control are described in detail. In addition, Screener3D can significantly
reduce $\alpha$ background rate from itself through topological
characteristics of particle trajectories, a distinctive feature of the gaseous
detector. Eventually, we can have a detector with fewer than one $\alpha$
background count per day. With a large active area of approximately 2000 cm2,
the TPC we propose may attain a sensitivity higher than 100 $\mu$Bq$\cdot$m-2
for surface $\alpha$ contaminations with two days of measurements. If we only
count $\alpha$ particles from specific $\alpha$ peaks, the sensitivity can be
further increased.

Screener3D is currently designed primarily for surface $\alpha$ particle
screening, and the ROI is typically above 1 MeV. We can also adopt the
detector for surface $\beta$ counting in the energy range below 1 MeV.
However, in this energy range, $\gamma$-rays from the detector and nearby
laboratory environment are the dominating background. Therefore, a
$\gamma$-shielding facility is required and is currently under design. We are
also investigating track reconstruction algorithms to better reconstruct the
more meandering $\beta$ tracks, aiming to locate the starting position with
high precision. Screening capability for surface $\beta$ will be reported in a
future manuscript.

##  V Acknowledgements

This work was supported by the grant from the Ministry of Science and
Technology of China (No. 2016YFA0400302) and the grants from the National
Natural Sciences Foundation of China (No. 11775142 and No. U1965201). We
appreciate the support from the Chinese Academy of Sciences Center for
Excellence in Particle Physics (CCEPP).

## References

  * (1) S. DellOro, S. Marcocci, M. Viel et al, Neutrinoless double beta decay:2015 review. Adv.High Energy Phys 2016, 2162659 (2016). http://dx.doi.org/10.1038/10.1155/2016/2162659
  * (2) L. Roszkowski, E.M. Sessolo, and S. Trojanowski, WIMP dark matter candidates and searches-current status and future prospects. Rept. Prog. Phys 81(6), 066201 (2018). http://dx.doi.org/10.1088/1361-6633/aab913
  * (3) J. Liu, X. Chen, and X. Ji, Current status of direct dark matter detection experiments. Nature Phys 13(3), 212–216 (2017). http://dx.doi.org/10.1038/nphys4039
  * (4) J.-P. Cheng, K.-J. Kang, J.-M. Li et al, The China Jinping Underground Laboratory and its Early Science. Ann. Rev. Nucl. Part. Sci 67,231–251 (2017). http://dx.doi.org/10.1146/annurev-nucl-102115-044842
  * (5) N. Abgrall, I.J. Arnquist, F.T. Avignone III, et al, The Majorana Demonstrator radioassay program. Nucl. Instrum. Meth A, 828, 22–36 (2016). http://dx.doi.org/10.1016/j.nima.2016.04.070
  * (6) D.S. Akerib, H.M. Araújo, X. Bai, et al, Radiogenic and Muon-Induced Backgrounds in the LUX Dark Matter Detector. Astropart. Phys 62,33–46 (2015). http://dx.doi.org/10.1016/j.astropartphys.2014.07.009
  * (7) F. Alessandria, E.Andreotti, R.Ardito et al, CUORE crystal validation runs: results on radioactive contamination and extrapolation to CUORE background. Astropart. Phys 35, 839–849 (2012). http://dx.doi.org/10.1016/j.astropartphys.2012.02.008
  * (8) H. Jiang, L.P. Jia, Q. Yue, et al, Limits on Light Weakly Interacting Massive Particles from the First 102.8 kg x day Data of the CDEX-10 Experiment. Phys. Rev. Lett 120(24),241301 (2018). http://dx.doi.org/10.1103/PhysRevLett.120.241301
  * (9) D.S. Leonard, P. Grinberg, P. Weber, et al, Systematic study of trace radioactive impurities in candidate construction materials for EXO-200. Nucl.Instrum. Meth. A 591, 490–509 (2008). http://dx.doi.org/10.1016/j.nima.2008.03.001
  * (10) X. Wang, X. Chen, C. Fu et al, Material Screening with HPGe Counting Station for PandaX Experiment. JINST 11(12),T12002 (2016). http://dx.doi.org/10.1088/1748-0221/11/12/T12002
  * (11) X. Guan, L. Ge, G. Zeng, et al. Determination of gross $\alpha$ and $\beta$ activities in Zouma River based on online HPGe gamma measurement system. Nucl. Sci. Tech. 31(12),120 (2020). http://dx.doi.org/10.1007/s41365-020-00828-0
  * (12) L. Yin, Q. Tian, X. Shao, et al. ICP-MS measurement of uranium and thorium contents in minerals in China. Nucl. Sci. Tech. 27(1) 10 (2016). http://dx.doi.org/10.1007/s41365-016-0018-5
  * (13) C. Alduino, K. Alfonso, D.R. Artusa, et al, Measurement of the two-neutrino double-beta decay half-life of 130Te with the CUORE-0 experiment. Eur. Phys. J. C 77(1), 13 (2017). http://dx.doi.org/10.1140/epjc/s10052-016-4498-6
  * (14) H. Zhang, H.G. Zhang, A. Abdukerim, et al, Dark matter direct search sensitivity of the PandaX-4T experiment. Sci. China Phys. Mech. Astron 62(3), 31011 (2019). http://dx.doi.org/10.1007/s11433-018-9259-0
  * (15) Ortec detectors. https://www.ortec-online.com/products/radiation-detectors/
  * (16) P.Loaiza, A.S. Barabash, A. Basharina-Freshville, et al, The BiPo-3 detector. Appl. Radiat. Isot123, 28242294(2017). https://doi.org/10.1016/j.apradiso.2017.01.021
  * (17) A.S. Barabash, A. Basharina-Freshville, E. Birdsall, et al, The BiPo-3 detector for the measurement of ultra low natural radioactivities of thin materials. J. Instrum12, P06002(2017). https://doi.org/10.1088/1748-0221/12/06/p06002
  * (18) XIA Ultralo 1800 gas counter. https://xia.com/ultralo.html
  * (19) R. Bunker, et al. The BetaCage, an Ultra-Sensitive Screener for Surface Contamination. AIP Conference Proceedings 1549, 132–35 (2013). https://doi.org/10.1063/1.4818093
  * (20) H. Ito, et al, Development of an Alpha-Particle Imaging Detector Based on a Low Radioactive Micro-Time-Projection Chamber. Nucl. Instrum. Meth. A. 953, 163050 (2020). https://doi.org/10.1016/j.nima.2019.163050
  * (21) S. Andriamonje, D. Attie, E. Berthoumieux at el, Development and performance of Microbulk Micromegas detectors. Journal of Instrumentation 4, 02001(2010). https://doi.org/10.1088/1748-0221/5/02/p02001
  * (22) L. Yang, J. Xu, Q. Li, et al. Performance of the CAT-TPC based on two-dimensional readout strips. Nucl. Sci. Tech. 32, 85 (2021). https://doi.org/10.1007/s41365-021-00919-6
  * (23) X. Chen, F. Fu, J. Galan, et al. PandaX-III: Searching for neutrinoless double beta decay with high pressure 136Xe gas time projection chambers. Sci. China Phys. Mech. Astron. 60, 061011 (2017). https://doi.org/10.1007/s11433-017-9028-0
  * (24) COMSOL. https://www.comsol.com
  * (25) X. Yan, X. Chen, Y. Chen, et al, Slow Control System for PandaX-III experiment. JINST16(05),T05004 (2021). https://doi.org/10.1088/1748-0221/16/05/t05004
  * (26) J. Castel, I. Coarasa, T. Dafni, et al, Background assessment for the TREX dark matter experiment. Eur. Phys. J. C. 79, 782 (2019). https://doi.org/10.1140/epjc/s10052-019-7282-6
  * (27) S. Agostinelli, J. Allison, K. Amako, et al, (GEANT4 Collaboration), GEANT4: A Simulation toolkit. Nucl. Instrum. Meth. A. 506, 250 (2003). https://doi.org/10.1016/S0168-9002(03)01368-8
  * (28) D. Budjá, A.M. Gangapshev, J. Gasparro, et al, Gamma-ray spectrometry of ultra low levels of radioactivity within the material screening program for the GERDA experiment. Applied Radiation and Isotopes 67(5), 755 – 758(2009). http://dx.doi.org/10.1016/j.apradiso.2009.01.019
  * (29) P. Agnes, L. Agostino, I.F.M. Albuquerque, et al, Results From the First Use of Low Radioactivity Argon in a Dark Matter Search. Phys. Rev. D 93(8),081101, (2016). http://dx.doi.org/10.1103/PhysRevD.93.081101
  * (30) H. Simgen, and G. Zuzel, Analysis of the 222Rn concentration in argon and a purification technique for gaseous and liquid argon. Appl Radiat Isot. 67(5), 922-5(2009). http://dx.doi.org/10.1016/j.apradiso.2009.01.058. Epub 2009 Jan 25. PMID: 19251429
  * (31) C.Y. Cao, N. Li, X.Y Yang at el, A practical approach of high precision U and Th concentration measurement in acrylic. Nucl. Instrum. Meth. A 1004, 165377(2021). https://doi.org/10.1016/j.nima.2021.165377
  * (32) F. Alessandria, R. Ardito, D.R. Artusa at el, Validation of techniques to mitigate copper surface contamination in CUORE. Astroparticle Physics 45, 13-22(2013). https://doi.org/10.1016/j.astropartphys.2013.02.005
  * (33) J. Galan, X. Chen, H. Du, et al, Topological background discrimination in the PandaX-III neutrinoless double beta decay experiment. J. Phys. G: Nucl. Part. Phys. 47, 045108 (2020). https://doi.org/10.1088/1361-6471/ab4dbe http://dx.doi.org/10.1007/s11433-017-9028-0
  * (34) R. Veenhof, Garfield - simulation of gaseous detectors. http://garfield.web.cern.ch/garfield/
  * (35) G. J. Feldman and R. D. Cousins, Unified approach to the classical statistical analysis of small signals. Phys. Rev. D 57, 3873 (1998). https://doi.org/10.1103/PhysRevD.57.3873


# Effect of Earth-Moon’s gravity on TianQin’s range acceleration noise

Xuefeng Zhang<EMAIL_ADDRESS>TianQin Research Center for Gravitational
Physics & School of Physics and Astronomy, Sun Yat-sen University (Zhuhai
Campus), Zhuhai 519082, P. R. China  Chengjian Luo<EMAIL_ADDRESS>TianQin Research Center for Gravitational Physics & School of Physics and
Astronomy, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082, P. R. China
Lei Jiao  TianQin Research Center for Gravitational Physics & School of
Physics and Astronomy, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082,
P. R. China  Bobing Ye  TianQin Research Center for Gravitational Physics &
School of Physics and Astronomy, Sun Yat-sen University (Zhuhai Campus),
Zhuhai 519082, P. R. China  Huimin Yuan  TianQin Research Center for
Gravitational Physics & School of Physics and Astronomy, Sun Yat-sen
University (Zhuhai Campus), Zhuhai 519082, P. R. China  Lin Cai
<EMAIL_ADDRESS>Center for Gravitational Experiments, School of Physics,
MOE Key Laboratory of Fundamental Physical Quantities Measurement & Hubei Key
Laboratory of Gravitation and Quantum Physics, PGMF, Huazhong University of
Science and Technology, Wuhan 430074, P. R. China  Defeng Gu  TianQin Research
Center for Gravitational Physics & School of Physics and Astronomy, Sun Yat-
sen University (Zhuhai Campus), Zhuhai 519082, P. R. China  Jianwei Mei
TianQin Research Center for Gravitational Physics & School of Physics and
Astronomy, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082, P. R. China
Jun Luo  TianQin Research Center for Gravitational Physics & School of Physics
and Astronomy, Sun Yat-sen University (Zhuhai Campus), Zhuhai 519082, P. R.
China  Center for Gravitational Experiments, School of Physics, MOE Key
Laboratory of Fundamental Physical Quantities Measurement & Hubei Key
Laboratory of Gravitation and Quantum Physics, PGMF, Huazhong University of
Science and Technology, Wuhan 430074, P. R. China

###### Abstract

TianQin is a proposed space gravitational-wave detection mission using
circular high Earth orbits. The geocentric concept has raised questions about
the disturbing effect of the nearby gravity field of the Earth-Moon system on
the highly-sensitive intersatellite ranging measurements. Here we examine the
issue through high precision numerical orbit simulation with detailed gravity-
field models. By evaluating range accelerations between distant free-falling
test masses, the study shows that the majority of the Earth-Moon’s gravity
disturbances are not in TianQin’s detection frequency band above $10^{-4}$ Hz,
and hence present no showstoppers to the mission.

††preprint: APS/123-QED

##  I  Introduction

The current TianQin design assumes high Earth orbits with an orbital radius of
$10^{5}$ km [1]. The nearly equilateral-triangle constellation stands almost
vertically to the ecliptic. High precision laser ranging interferometry tracks
distance changes between well-protected test masses (TM) in separate drag-free
controlled satellites, within a preliminary frequency range of $10^{-4}-1$ Hz.
Proximity to the Earth has certain benefits such as lower launch cost, shorter
transfer duration, easier communication, availability of GNSS (Global
Navigation Satellite System), etc. Other geocentric mission concepts [2]
include OMEGA [3, 4], GEOGRAWI/gLISA [5], GADFLI [6], B-DECIGO [7], etc.
TianQin’s orbit is different from them in both the orbital radius and
orientation.

Space-based gravitational-wave (GW) detectors are subject to various
influences from the surrounding environment, including gravity field, thermal
radiation, plasma, magnetic field, solar-wind particles, galactic cosmic rays,
micrometeorites, etc., that are existential in outer space. Environmental
effects can strongly affect performance and lifetime of the sciencecraft.
Quite prominently, the space gravity-field environment, encompassing
gravitational perturbations from the central and other celestial bodies, plays
an important role. It is particularly the case for geocentric missions due to
their closeness to the Earth and the Moon. More specifically, the effects on
TianQin are two-fold. On large scales, the perturbations distort the nominal
equilateral triangle of the constellation. The resulting unequal and time-
varying arm-lengths have far-reaching implications on science payload design
and data processing strategies (e.g., [8, 9]). The distortion can be reduced
by orbit optimization [10, 11] and control to meet the stability requirements
of the science payloads. On small scales, the perturbations impinge on TMs’
geodesic motion under nearly pure gravity. Since space-based detectors
accurately measure arm-length variations between TMs, they respond not only to
GWs (radiation zone), but equally well to Newtonian gravity fields (near
zone), in targeted frequency bands. Appearing as environmental noise, the
latter should be avoided or mitigated.

Ideally, gravitational perturbations in space should only manifest as long-
term and slow changes in inter-spacecraft displacement measurements. If there
exists a proper separation in the frequencies of gravity-field fluctuations
and GWs, then the GW signals, superimposed on top of a smooth and slow-varying
background, can be extracted ([12], Sec 2.1.1). Therefore the GW detection
relies heavily on the “quietness” of the ambient gravity-field environment in
the measurement band.

The problem of environmental gravity disturbances were recognized early on in
designing ground-based detectors [13], and hence is not unique for space
missions, where the problem is thought to be much less severe. In ground-based
detectors, Newtonian or gravity-gradient noise caused by terrestrial gravity
fluctuations poses a limitation to sensitivity improvement below $\sim 10$ Hz
[14]. Multiple strategies have been developed to effectively mitigate such
noise, and the techniques have a major influence on designing next-generation
ground-based detectors.

If not handled properly in space GW detection, disturbing gravity fields may
induce excessive “orbital noise” that encroaches on the sensitivity curve,
causing a situation somewhat similar to galactic foreground noise in lower
frequencies [15]. The potential risk has drawn attentions for TianQin, and may
raise concern for other geocentric concepts as well. With regard to LISA [16],
the majority of the effect is expected to be out of the sensitive frequency
band because of its heliocentric yearly orbits and being placed far away from
the Earth-Moon system ($\sim 20^{\circ}$ trailing angle, $\sim 5\times 10^{7}$
km).

In general, gravity disturbances in space constitute an important potential
noise source for inter-spacecraft measurements. In this work, we aim to
determine the amplitudes and frequencies of the disturbances for TianQin’s
orbit, and quantitatively evaluate the impact on TianQin’s acceleration noise
requirement. The forward modeling takes into account a variety of main
gravitational perturbations including the gravity fields of the Earth (static
and tidal), the Moon, and the Sun, as well as other solar-system bodies. It
requires realistic and accurate orbit propagation that is also used in
performance assessment and data analysis of gravity mapping missions such as
GRACE [17], GRACE Follow-On [18], GOCE [19], and GRAIL [20]. However, for
TianQin, a problem with insufficiency of double precision arithmetic has
emerged owing to the high measurement accuracy requirement over the long
baseline. To tackle the issue, an earlier attempt was made in [21], where
analytical expansions of perturbed orbits were derived. Unfortunately, the
approach cannot handle complicated gravity field models, and only the Earth’s
static gravity field was considered without realistic Earth’s rotation
(precession, nutation, etc), the Earth’s tides, and third bodies. It motivated
us to take a fully numerical approach to be shown in this paper. For other
works regarding environmental magnetic and plasmic effects on TianQin, one can
refer to, e.g., [22, 23, 24].

This is our third paper of the concept study series on TianQin’s orbit and
constellation. It is based on the previous work of orbit optimization and
constellation stability [10, 11], and shifts the attention to small-scale
orbital motion through much refined simulation. The paper is organized as
follows. In Section II, three types of intersatellite observables are
analysed, and the range acceleration is chosen for evaluating the impact. In
Section III, we describe the high precision orbit propagator, detailed force
models, and orbital parameters used in the assessment. Section IV presents the
amplitude spectral density (ASD) results of the calculated range
accelerations. At the end, the conclusions are made in Section V.

##  II  Observables and Criteria

For the evaluation purpose, the numerical simulation should provide an
observable accuracy better than the instrumental measurement noise level. The
selectable intersatellite observables include the (instantaneous) range,
range-rate, and range acceleration. Mathematically, they are interchangeable
by differentiation and integration. But their numerical calculations require
different computational resources. Here we estimate the magnitudes of their
numerical ranges (numbers of significant digits required) for TianQin. First,
the range between two satellites is given by

$\rho=|\mathbf{r}_{2}-\mathbf{r}_{1}|,$ (1)

where $\mathbf{r}_{1,2}$ denotes the position vector of each satellite
relative to the Earth’s center. Taking the baseline $1.7\times 10^{8}$ m and
the displacement measurement noise $1\times 10^{-12}$ m/Hz1/2 [1], the
numerical representation of the range observable requires at least 20 digits,
exceeding the 16 digits of the double-precision format (64 bits). Second, the
range rate reads

$\dot{\rho}=\hat{\mathbf{e}}_{12}\cdot(\dot{\mathbf{r}}_{2}-\dot{\mathbf{r}}_{1}),$
(2)

with the unit vector
$\hat{\mathbf{e}}_{12}=(\mathbf{r}_{2}-\mathbf{r}_{1})/\rho$. The relative
velocities between the TianQin satellites is expected to be within $\pm 5$ m/s
[10]. Taking the range-rate measurement noise $5\times 10^{-14}$ m/s/Hz1/2
($\sim 2\pi f\times 10^{-12}$ m/Hz1/2 at the crossover frequency $f\sim
10^{-2}$ Hz of the displacement and residual acceleration noises [1]), the
dynamical range of $\dot{\rho}$ takes up about 15 digits. Third,
differentiating the range rate yields the equation for the range acceleration:

$\ddot{\rho}=\hat{\mathbf{e}}_{12}\cdot(\ddot{\mathbf{r}}_{2}-\ddot{\mathbf{r}}_{1})+\frac{1}{\rho}\left(|\dot{\mathbf{r}}_{2}-\dot{\mathbf{r}}_{1}|^{2}-\dot{\rho}^{2}\right),$
(3)

where, on the right-hand side, the first term represents projected
differential acceleration, and the second term centrifugal acceleration. The
gravitational acceleration of one TianQin satellite is in the order of
$10^{-2}$ m/s2. This is 13 order of magnitude greater than the residual
acceleration noise level of one TM, i.e., $1\times 10^{-15}$ m/s2 [1]. For
either $\dot{\rho}$ or $\ddot{\rho}$, if one takes into account that numerical
errors compounded over time may occupy 2-3 digits, and redundant numerical
accuracy another 1-2 digits, then the requirement would exceed 16 digits.
Therefore, the commonly used double-precision arithmetic is insufficient in
representing the intersatellite observables, and the associated roundoff error
becomes a bottleneck for precision improvement (cf. Fig. 1).

Among the three observables, the range acceleration appears more favorable for
taking up less digits in numerical computation. In the frequency domain,
acceleration and displacement can be easily converted. For evaluating gravity
disturbances in space, we henceforward adopt the range acceleration as the
main observable (cf. [25]), and directly compare its ASD with the
intersatellite residual acceleration noise requirement $\sqrt{2}\times
10^{-15}$ m/s2/Hz1/2 at $10^{-4}-10^{-2}$ Hz as the criteria, which is simply
$\sqrt{2}$ of the residual acceleration noise of a single TM (cf. [26]). Note
that this flat noise requirement is preliminary and expected to be relaxed
near $10^{-4}$ Hz in the future [1, 16].

##  III  Simulation and Force Models

The evaluation requires careful calculation and modeling of satellite orbits
and gravity fields. The accuracy of numerical integration must surpass the
noise requirement $\sqrt{2}\times 10^{-15}$ m/s2/Hz1/2 of the range
acceleration observable by at least one order of magnitude. The force modeling
should be sufficiently detailed and up-to-date to reflect as many significant
gravity disturbances as possible, particularly those that may enter the
detection band.

###  III.1 Quadruple Precision Orbit Propagation

There exist a few strategies to tackle the inadequacy of double precision. A
straightforward way is by extending to 34 significant digits with quadruple
precision arithmetic (128 bits). The potential downsides are low execution
speed and heavy programming workload. Following this “brute force” approach,
the TQPOP (TianQin Quadruple Precision Orbit Propagator) program based on
MATLAB has been developed so as to evaluate the range acceleration at
$<10^{-15}$ m/s2/Hz1/2 levels. The quadruple precision data type is applied to
all the necessary aspects of the program, including parameter inputs,
ephemeris data outputs, reference frame transformations, time conversion,
numerical integration, force models, etc. For the nearly circular high orbits,
the integrator uses the 8th-order embedded Prince-Dormand (DP87) method [27]
with a constant step size of 50 seconds (Nyquist frequency $10^{-2}$ Hz). The
algorithm provides a relative truncation error of $<10^{-20}$ (more than 20
significant digits) in both satellite positions and velocities. Thereby the
range acceleration error is estimated to be $<10^{-22}$ m/s2 and well below
$10^{-15}$ m/s2. The roundoff error due to finite digits is approximately
$10^{-33}$ m/s2/Hz1/2, and no longer poses a limiting factor (see Fig. 1),
which otherwise would overwhelm gravity field signals in the case of double
precision. To mitigate low efficiency of quadruple precision calculations,
great effort was made on optimizing code execution to have reduced the run
time significantly. Other quadruple precision orbit simulations can be found
in, e.g., [28, 29].

Figure 1:  The ASD of the range acceleration $\ddot{\rho}$ between two
satellites in circular orbits of the radius $10^{5}$ km and separated by
$120^{\circ}$ in phase. The orbits are integrated with a constant step size of
50 s and under the central force of the Earth’s point mass. The roundoff error
of the quadruple precision arithmetic is at the level of $10^{-33}$ m/s2/Hz1/2
(at $10^{-4}$ Hz). The curve tilts up toward low frequencies due to
accumulation of roundoff errors over time. For comparison, the roundoff error
of double precision is also shown, but at a much higher level of $10^{-15}$
m/s2/Hz1/2 (at $10^{-4}$ Hz), hence not sufficient for the accuracy
requirement.

###  III.2 Detailed Force Models

As the satellites are drag-free controlled, we only consider pure free-fall
orbits of the TMs in order to focus on the gravitational perturbative effects.
Excluding non-gravitational forces, the force models implemented are
summarized in Table 1. The types of the gravity field models are comparable
with those used in the Earth’s gravity field determination in satellite
missions such as GRACE [17, 30] and GOCE [19].

The solar system ephemeris uses DE430 [31] including all eight planets and the
Moon. The effect from the main belt asteroids is estimated not to enter LISA’s
detection band [32] , nor TianQin’s due to the shorter arm-length. Hence they
are not included in the simulation.

For the Earth’s orientation, the International Astronomical Union (IAU) 2006
precession and IAU 2000A nutation models [33] are used with the help of the
Standards Of Fundamental Astronomy (SOFA) software collection [34]. The
Earth’s polar motion adopts the IERS Earth Orientation Parameters (EOP) 14 C04
data series [35].

Earth’s non-spherical static gravity field is provided by the EGM2008 model
[36], following the recommendation of IERS (2010) [33]. The normalized
spherical harmonic coefficients ($\overline{C}_{nm}$, $\overline{S}_{nm}$) are
kept up to the 12th degree and order. High-degree terms decay rapidly with
increasing radius as $1/r^{n+1}$. Our numerical tests and perturbation
analysis shows that the effect from the 9th-degree gravity field has already
dropped below $10^{-15}$ m/s2/Hz1/2. The contribution from the 12th degree
sinks deeper to the level of $10^{-18}$ m/s2/Hz1/2. Hence we deem it safe to
truncate at the degree and order 12.

Temporal variations of the Earth’s gravity field are added as corrections to
the spherical harmonic coefficients. To model the Earth’s tidal effects, we
have followed IERS (2010) [33] and taken into account solid Earth tides (an-
elastic), ocean tides, solid Earth pole tide, and ocean pole tide, as
specified in Table 1. The widely-used ocean tide model FES2004 [37] includes
long-period ($\Omega_{1}$, $\Omega_{2}$, $S_{a}$, $S_{sa}$, $M_{m}$, $M_{f}$,
$M_{tm}$, $M_{sqm}$), diurnal ($Q_{1}$, $O_{1}$, $P_{1}$, $K_{1}$), semi-
diurnal ($2N_{2}$, $N_{2}$, $M_{2}$, $S_{2}$, $K_{2}$), and quarter-diurnal
($M_{4}$) waves. The coefficients up to the degree and order 10 are used.
Additionally, atmospheric tides are incorporated, though their effect is small
compared to the solid Earth and ocean tides. The associated model [38]
consists of the diurnal and semi-diurnal waves $S_{1}$ and $S_{2}$ in the
highest frequency constituents. The correction is made up to the degree 8 and
order 5. The non-tidal temporal gravity changes have been estimated to be
orders of magnitude smaller than the static gravity [39], and will be
discussed elsewhere.

Table 1:  The list of force models implemented in the simulation.  Models | Specifications
---|---
Solar system ephemeris | JPL DE430 [31]
Earth’s precession & nutation | IAU 2006/2000A [33]
Earth’s polar motion | EOP 14 C04 [35]
Earth’s static gravity field | EGM2008 ($n=12$) [36]
Solid Earth tides | IERS (2010) [33]
Ocean tides | FES2004 ($n=10$) [37]
Solid Earth pole tide | IERS (2010) [33]
Ocean pole tide | Desai (2003) [33]
Atmospheric tides | Biancale & Bode (2003) [38]
Moon’s libration | JPL DE430 [31]
Moon’s static gravity field | GL0660B ($n=7$) [20]
Sun’s orientation | IAU [40], Table 1
Sun’s $J_{2}$ | IAU [40], Table 1
relativistic effect | post-Newtonian [41]

Moon’s liberation varies about $\pm 8^{\circ}$, and is provided by DE430
([31], IIE). For the Moon’s static gravity field, we use GL0660B [20] up to
the degree and order 7, and the effect of the 7th degree and order alone is
below $10^{-20}$ m/s2/Hz1/2. The model was obtained from the GRAIL (Gravity
Recovery and Interior Laboratory) mission [20] with improved low-degree
harmonics. The lunar tide is not included, since the effect is quite small and
(semi-)monthly periodic, hence out of the detection band, owning to the Moon’s
tidal locking with the Earth. The values of the Sun’s oblateness $J_{2}$ and
orientation are taken from [40] (see also [31], Table 9). Moreover,
relativistic effect is added as post-Newtonian correction terms to the
equations of motion [41]. The effect is slow varying and expected to be
outside the detection band.

###  III.3 Orbital Parameters

The initial orbital parameters are given in Table 2. The integration lasts for
one observation window of three months [1], that is, from 06 Jun. to 04 Sep.
2004 for 90 days, when the orbital plane is facing the Sun within $\pm
45^{\circ}$. The year 2004 is chosen without preference but to take advantage
of the available EOP observation data, which is more accurate than prediction
in 2030s. Our tests have shown that the dominant spectral behavior does not
depend on a specific year chosen.

To make the simulation more realistic, we use the optimized initial orbital
elements in Table 2 that can meet TianQin’s constellation stability
requirement (e.g., the breathing angles within $60\pm 0.1^{\circ}$) for three
months [10, 11]. The optimization removes linear drift in the arm-lengths and
breathing angles, and prevents the nearly equilateral-triangle constellation
from having severe distortion. The initial eccentricities are set to zeros to
keep the orbits almost circular. Note that even if one starts with less
optimized initial orbital elements (e.g., the nominal values, $a=10^{5}$ km,
etc), the dominant spectral behavior of three months (cf. Fig. 2) will be
unaffected.

Table 2:  The initial orbital elements of the TianQin constellation in the J2000-based Earth-centered equatorial coordinate system at the epoch 06 Jun. 2004, 00:00:00 UTC for the evaluation purposes. Here $a$ denotes the semimajor axis, $e$ the eccentricity, $i$ the inclination, $\Omega$ the longitude of ascending node, $\omega$ the argument of periapsis, and $\nu^{\rm ini}$ the truly anomaly.  | $a$ | $e$ | $i$ | $\Omega$ | $\omega$ | $\nu^{\rm ini}$
---|---|---|---|---|---|---
SC1 |  $100000.0$ km | 0 | $74.5^{\circ}$ | $211.6^{\circ}$ | $0^{\circ}$ | $30^{\circ}$
SC2 |  $100009.5$ km | 0 | $74.5^{\circ}$ | $211.6^{\circ}$ | $0^{\circ}$ | $150^{\circ}$
SC3 |  $99995.0$ km | 0 | $74.5^{\circ}$ | $211.6^{\circ}$ | $0^{\circ}$ | $270^{\circ}$

##  IV  Spectral Results

It should be emphasized that the purpose of this work is to determine the
frequency-domain effects (especially $>10^{-4}$ Hz) of various gravity
disturbances on the range acceleration observable, and it concerns less about
the absolute accuracy of an integrated orbit, which may drift away from true
values over long time scales outside the frequency band of interest.

###  IV.1 Total Effect

The overall result of the range acceleration ASD is presented in Fig. 2 for
the arm SC1-SC2 using the models of Table 1 assembled together in the
simulation. The plots for the other two arms severely overlap with the first
one, hence not presented for clarity.

In the frequency domain, the gravity field signals are dominating below
$10^{-4}$ Hz, and roll off rapidly in amplitudes toward high frequencies, and
intersect with the lower end of the range acceleration noise requirement at
$1\times 10^{-4}$ Hz The steep fall-off roughly follows a power law of $\sim
1.7\times 10^{-15}$ m/s2/Hz1/2 $\times(0.1\ \textrm{mHz}/f)^{24}$ near
$f=1\times 10^{-4}$ Hz. The plot demonstrates that the effect of the gravity
field models in Table 1 does not enter the detection band $>10^{-4}$ Hz. Note
that the slanted part of the ASD curve ($<10^{-17}$ m/s2/Hz1/2 and $>10^{-4}$
Hz, marked by “Numerical error” in Fig. 2) is an artifact of numerical
interpolation of the EOP data.

The frequency-domain behavior somewhat resembles the one in the intersatellite
laser ranging measurement result of the GRACE Follow-On mission [18], which
also shows a steep fall-off, but at a higher frequency ($\sim 4\times 10^{-2}$
Hz) because of its low orbit altitude of approximately 500 km.

Figure 2:  The range acceleration ASD of two TianQin satellites SC1 and SC2,
calculated in quadruple precision and with the models of Table 1 and step size
50s using 90 days of data. The orbital period 3.6 days corresponds to
$3.2\times 10^{-6}$ Hz. The plots for SC1-SC3 and SC2-SC3 are nearly identical
to the one shown above.

###  IV.2 Effect Breakdown

Now we examine various contributions to the total range acceleration ASD of
the arm SC1-SC2. The result is presented in Fig. 3.

The Earth’s non-spherical static gravity field with degrees $n\geq 3$
dominates above $5\times 10^{-5}$ Hz in the total ASD, indicated by the
overlapping of the blue and red curves. The effect decreases rapidly toward
high frequencies and impinges on the noise requirement at $10^{-4}$ Hz. At the
high altitude of TianQin, high-degree harmonics are effectively attenuated.

The contribution from the Moon’s non-spherical static gravity field ($n\geq
2$) is minute and only sticks out in the low-frequency region. One may have
expected so since the Moon is a slowly rotating body. The same argument also
justifies the omission of the lunar tides in the simulation. Nevertheless, the
Moon’s point mass and its orbital motion play an important role, largely
accounting for the total ASD below $4\times 10^{-5}$ Hz.

The effects of relativity and Earth’s tidal gravity field (solid Earth,
oceanic, pole, and atmospheric, cf. Table 1) are considerably smaller than the
total effect, and both peak at low frequencies away from the detection band.
These low frequency disturbances show no significant coupling into high
frequencies, and do not induce pronounced range acceleration response above
$10^{-4}$ Hz.

Figure 3:  Components of the range acceleration ASD of two TianQin satellites
SC1 and SC2. The total ASD in Fig. 2 is duplicated in dotted red curve for
comparison.

###  IV.3 Model Errors

The models inevitably contain errors. To estimate their effect, a
straightforward way is to determine whether discrepancies between different
models can significantly alter the spectral result in Fig. 2.

For a cross-check, we test another set of gravity field models shown in Table
3 where several replacements are made to Table 1. The substitute models are
deemed less accurate than the corresponding more recent ones in Table 1, and
thus can mimic model errors (see, e.g., [42] ). In Fig. 4, both ASD results
show good agreement with each other, and the difference $\Delta\ddot{\rho}$ is
well below the noise requirement. In addition, the spectral behavior above
$10^{-15}$ m/s2/Hz1/2 is also confirmed by running the flight-qualified, open
source program GMAT [43] in double precision. Hence the overall frequency-
domain behavior appears to be robust, which instills more confidence in the
results.

Table 3:  The list of replacement force models to Table 1 used for spectrum comparison.  Models | Specifications
---|---
Solar system ephemeris | JPL DE405 [44]
Earth’s precession & nutation | IAU 1976/1980 [41]
Earth’s static gravity field | EGM96 ($n=12$) [45]
Moon’s libration | JPL DE405 [44]
Moon’s static gravity field | LP165P ($n=7$) [46]
Figure 4:  The range acceleration ASD of two TianQin satellites SC1 and SC2
using the models of Table 3 for replacement to Table 1 (Tab1+3, blue). The ASD
in Fig. 2 is duplicated in dashed red curve (Tab1) for comparison. The ASD of
their difference $\Delta\ddot{\rho}$ is marked by cyan curve.

##  V  Conclusion

The TianQin mission, to be deployed in a high Earth orbit, shares
technological similarities with low-Earth gravimetry missions using satellite-
to-satellite tracking. They diverge on a key point that the Earth’s gravity
field signals targeted in gravimetry missions become environmental noise in
TianQin’s GW detection. Hence TianQin must keep a safe distance from the Earth
by flying high enough, so as to push the Earth’s gravity field interference
out of the detection band. This work has been devoted to evaluating and
examining this type of effect, and two main conclusions can be drawn here.

1\. With the orbital radius of $10^{5}$ km for TianQin, the current models
show that the effect of the Earth-Moon’s gravity field dominates at low
frequencies, and that the amplitude rolls off rapidly toward high frequencies
and intersects with the range acceleration noise requirement ($\sqrt{2}\times
10^{-15}$ m/s2/Hz1/2) at $10^{-4}$ Hz, right on the lower end of the
preliminary detection band. To provide more context, the gravity field signals
from GRACE-FO laser ranging interferometry along a $\sim 200$ km baseline
falls off at about $4\times 10^{-2}$ Hz [18] with an orbit altitude of $\sim
500$ km.

2\. The high-precision numerical simulations help to rule out the majority of
the perturbing gravity sources for TianQin, including the Sun’s point mass and
its $J_{2}$, the solar system planets’ point masses (under their orbital
motion), the Earth’s static gravity (with its rotation), the Earth’s tidal
gravity changes (solid Earth, oceanic, pole, and atmospheric), Moon’s static
gravity, relativistic effect, etc (cf. Table 1). These effects are slowly
varying, not entering the detection band, and present no show-stoppers for
TianQin. The Newtonian gravity-field environment at a distance of $10^{5}$ km
from the Earth is expected to be fairly “quiet” for TianQin.

The results can provide useful inputs and guidelines to several aspects of the
mission concept studies, such as orbit selection, noise reduction, and data
processing. For future works, further refined gravity models will be explored
to identify other possible noise sources. On another note, the strong low-
frequency gravity field signals ($<10^{-4}$ Hz) illustrated in Fig 2, which
carry long-wavelength gravity information, may find potential applications in
geodesy and geophysics [47]. This may help to enrich TianQin’s secondary
science output.

###### Acknowledgements.

The authors thank Qiong Li, Xiaoli Su, Fan Yang, Hao Zhou, Vitali Müller,
Gerhard Heinzel, Bo Xu, Shuai Liu, Yi-Ming Hu, Hui-Zong Duan, Lin Zhu, Jinxiu
Zhang, Cheng-Gang Shao, Shan-Qing Yang, Hsien-Chi Yeh, and anonymous referees
for helpful discussions and comments. The work is supported by the National
Key R&D Program of China (No. 2020YFC2201202). X.Z. is supported by NSFC Grant
No. 11805287.

## References

  * J. Luo et al. [2016] J. Luo et al., Class. Quantum Grav. 33, 035010 (2016).
  * NAS [2012] _NASA Gravitational-wave mission concept study final report_ (2012), eprint https://pcos.gsfc.nasa.gov/physpag/GW_Study_Rev3_ Aug2012-Final.pdf.
  * Hiscock and Hellings [1997] B. Hiscock and R. W. Hellings, Bull. Am. Astron. Soc. 29, 1312 (1997).
  * Hellings et al. [2011] R. W. Hellings, S. L. Larson, S. Jensen, C. Fish, M. Benacquista, N. Cornish, and R. Lang, _A low-cost, high-performance space gravitational astronomy mission_ (2011), eprint https://pcos.gsfc.nasa.gov/studies/rfi/GWRFI-0007-Hellings.pdf.
  * Tinto et al. [2015] M. Tinto, D. DeBra, S. Buchman, and S. Tilley, Rev. Sci. Instrum. 86, 014501 (2015).
  * McWilliams [2011] S. T. McWilliams, _Geostationary Antenna for Disturbance-Free Laser Interferometry (GADFLI)_ , arXiv:1111.3708 [astro-ph.IM] (2011).
  * S. Kawamura, et al. [2018] S. Kawamura, et al., Int. J. Mod. Phys. D 27, 1845001 (2018).
  * Folkner et al. [1997] W. M. Folkner, F. Hechler, T. H. Sweetser, M. A. Vincent, and P. L. Bender, Class. Quantum Grav. 14, 1405 (1997).
  * Tinto and Dhurandhar [2014] M. Tinto and S. V. Dhurandhar, Living Reviews in Relativity 17, 6 (2014).
  * Ye et al. [2019] B. Ye, X. Zhang, M. Zhou, Y. Wang, H. Yuan, D. Gu, Y. Ding, J. Zhang, J. Mei, and J. Luo, Int. J. Mod. Phys. D 28, 1950121 (2019).
  * Tan et al. [2020] Z. Tan, B. Ye, and X. Zhang, Int. J. Mod. Phys. D 29, 2050056 (2020).
  * Barke [2015] S. Barke, PhD dissertation, Gottfried Wilhelm Leibniz Universität Hannover (2015).
  * Weiss and Muehlner [1972] R. Weiss and D. Muehlner, Quarterly Progress Report, Research Laboratory of Electronics (MIT) 105, 54 (1972).
  * Harms [2019] J. Harms, Living Reviews in Relativity 22, 6 (2019).
  * Nelemans [2009] G. Nelemans, Class. Quantum Grav. 26, 094030 (2009).
  * LIS [2017] _LISA Laser Interferometer Space Antenna, A proposal in response to the ESA call for L3 mission concepts_ , arXiv:1702.00786 (2017).
  * Tapley et al. [2004a] B. D. Tapley, S. Bettadpur, M. Watkins, and C. Reigber, Geophys. Res. Lett. 31, L09607 (2004a).
  * K. Abich et al. [2019] K. Abich et al., Phys. Rev. Lett. 123, 031101 (2019).
  * Floberghagen et al. [2011] R. Floberghagen, M. Fehringer, D. Lamarre, D. Muzi, B. Frommknecht, C. Steiger, J. Piñeiro, and A. da Costa, J. Geod. 85, 749 (2011).
  * Konopliv et al. [2013] A. S. Konopliv, R. S. Park, D.-N. Yuan, S. W. Asmar, M. M. Watkins, J. G. Williams, E. Fahnestock, G. Kruizinga, M. Paik, D. Strekalov, et al., J. Geophys. Res. Planets 118, 1 (2013).
  * Liu et al. [2019] S. Liu, L.-X. Liu, Z.-B. Tan, Y.-M. Hu, X. Zhang, J. dong Zhang, J. Zhang, J. Mei, and J. Luo (2019), internal technical report (unpublished).
  * Lu et al. [2020a] L.-F. Lu, Y. Liu, H.-Z. Duan, Y.-Z. Jiang, and H.-C. Yeh, Plasma Sci. Technol. 22, 115301 (2020a).
  * Su et al. [2020] W. Su, Y. Wang, Z.-B. Zhou, Y.-Z. Bai, Y. Guo, C. Zhou, T. Lee, M. Wang, M.-Y. Zhou, T. Shi, et al., _Analyses of residual accelerations for TianQin based on the global MHD simulation_ , arXiv:2004.00254 (2020).
  * Lu et al. [2020b] L.-F. Lu, W. Su, X. Zhang, Z.-G. He, H.-Z. Duan, Y.-Z. Jiang, and H.-C. Yeh, _Effects of the space plasma density oscillation on the inter-spacecraft laser ranging for TianQin gravitational wave observatory_ (2020b), submitted.
  * Müller [2017] V. Müller, PhD dissertation, Gottfried Wilhelm Leibniz Universität Hannover (2017).
  * M. Armano et al. [2016] M. Armano et al., Phys. Rev. Lett. 116, 231101 (2016).
  * Prince and Dormand [1981] P. J. Prince and J. R. Dormand, J. Comp. Appl. Math. 7, 67 (1981).
  * McCullough et al. [2015] C. McCullough, S. Bettadpur, and K. McDonald, Journal of Spacecraft and Rockets 52, 766 (2015).
  * Wöske et al. [2016] F. Wöske, T. Kato, M. List, and B. Rievers, Adv. Astr. Sci. 158, 2445 (2016).
  * Tapley et al. [2004b] B. D. Tapley, S. Bettadpur, J. C. Ries, P. F. Thompson, and M. M. Watkins, Science 305, 503 (2004b).
  * Folkner et al. [2014] W. M. Folkner, J. G. Williams, D. H. Boggs, R. S. Park, and P. Kuchynka, IPN Progress Report 42-196, Jet Propulsion Laboratory (2014).
  * Bronicki et al. [2018] D. Bronicki, B. Bolen, and S. Larson (2018), poster in the 12th International LISA Symposium.
  * Petit and Luzum [2010] G. Petit and B. Luzum, Technical Report 36, BUREAU INTERNATIONAL DES POIDS ET MESURES SEVRES (FRANCE) (2010).
  * [34] eprint http://www.iausofa.org/.
  * [35] eprint http://hpiers.obspm.fr/iers/eop/eopc04/.
  * Pavlis et al. [2012] N. K. Pavlis, S. A. Holmes, S. C. Kenyon, and J. K. Factor, J. Geophys. Res. 117, B04406 (2012).
  * Lyard et al. [2006] F. Lyard, F. Lefevre, T. Letellier, and O. Francis, Ocean Dyn. 56, 394 (2006).
  * Biancale and Bode [2006] R. Biancale and A. Bode, STR 06/01, GeoForschungsZentrum Potsdam (2006).
  * Gruber et al. [2011] T. Gruber, J. L. Bamber, M. F. P. Bierkens, H. Dobslaw, M. Murböck, M. Thomas, L. P. H. Van Beek, T. van Dam, L. L. A. Vermeersen, and P. N. A. M. Visser, Earth Syst. Sci. Data 3, 19 (2011).
  * Archinal et al. [2011] B. A. Archinal, M. F. A’Hearn, E. Bowell, A. Conrad, G. J. Consolmagno, R. Courtin, T. Fukushima, D. Hestroffer, J. L. Hilton, G. A. Krasinsky, et al., Celest. Mech. Dyn. Astr. 109, 101 (2011).
  * McCarthy [1996] D. D. McCarthy, Technical Report 21, Central Bureau of IERS - Observatoire de Paris (1996).
  * Loomis et al. [2012] B. D. Loomis, R. S. Nerem, and S. B. Luthcke, Journal of Geodesy 86, 319 (2012).
  * [43] eprint http://gmatcentral.org.
  * Standish [1998] E. M. Standish, IOM 312.F-98-048, Jet Propulsion Laboratory (1998).
  * F. Lemoine et al. [1998] F. Lemoine et al., Technical Paper NASA/TP-1998-206861, NASA, Greenbelt, MD (1998).
  * Konopliv et al. [2001] A. S. Konopliv, S. W. Asmar, E. Carranza, W. L. Sjogren, and D. N. Yuan, Icarus 150, 1 (2001).
  * Jiao and Zhang [2020] L. Jiao and X. Zhang (2020), poster in 2019 Annual Conference of HUST CGE and SYSU TRC.