Axion Dark Matter eXperiment: Run 1A Analysis Details

The ADMX collaboration gathered data for its Run 1A axion dark matter search from January to June 2017, scanning with an axion haloscope over the frequency range 645-680 MHz (2.66-2.81 ueV in axion mass) at DFSZ sensitivity. The resulting axion search found no axion-like signals comprising all the dark matter in the form of a virialized galactic halo over the entire frequency range, implying lower bound exclusion limits at or below DFSZ coupling at the 90% confidence level. This paper presents expanded details of the axion search analysis of Run 1A, including review of relevant experimental systems, data-taking operations, preparation and interpretation of raw data, axion search methodology, candidate handling, and final axion limits.


I. INTRODUCTION
An axion field is a consequence of the Peccei-Quinn solution to the strong-CP problem of particle physics [1][2][3].The Peccei-Quinn mechanism implies several dynamical processes that produce cold axions in the early Universe.These axions may constitute some fraction or all of the dark matter [4][5][6].An axion field serves as an excellent candidate for cold dark matter as it has very small primordial velocity dispersion, feeble couplings to itself and to the Standard Model fields.It describes a very weakly interacting massive particle, the axion, with lifetime much longer than the age of the Universe.Matter under these conditions is expected to collapse into large scale structures that seed the Universe's galaxies, including our own.For axions to saturate the ΛCDM model's dark-matter density, many numerical and analytical studies of QCD (quantum chromo-dynamics) prefer an axion mass in the 1 − 100µeV range [7][8][9][10][11].
Direct searches for axions use several techniques, based in either pure laboratory methods or by providing the requisite axions through external sources.The Axion Dark Matter eXperiment (ADMX) utilizes the axion haloscope technique [16], and consists of a cold microwave cavity threaded by a static magnetic field coupled to the ambient axion field via an inverse-Primakoff process.Axions are resonantly converted to photons in a cavity mode tuned to the axion field's frequencies, which are centered at f = ε/h where h is the Planck constant and ε is the expected total axion energy consisting of the rest mass energy and kinetic energy.The microwave haloscope has proven to be the most sensitive technique to search for axions as the dark matter over the favored range.Axion searches prior to the year 2017, including ADMX, have managed to approach KSVZ sensitivity, but not the stricter DFSZ model under the conservative assumption of a virialized galactic halo.ADMX entered as a DOE Generation 2 dark-matter search in 2014 with the explicit target of reaching DFSZ sensitivity over the frequency range 0.5−10 GHz (2−40µeV) [17], see Fig. 1, and more recently has been split into G-2 and EFR (extended frequency range) efforts.The first of these searches, referred to as "Run 1A," gathered data between January and June 2017 over the frequency range 645 − 680 MHz (2.66 − 2.81µeV) and found no axions in that range [18].This paper explains in more detail the operations and analyses that resulted in the Run 1A search reaching DFSZ limits.
The remainder of the this paper is structured as follows.Section II provides a brief overview of the ADMX apparatus, with focus on the axion converter-detector assembly, including the cavity and magnet, the radio frequency (RF) receiver chain, and the data acquisition system.Section III presents the data-taking operations, including search guidelines, data-taking cadence, scanning operations, logging of data, and candidate rescans.Section IV reviews the preparations performed on the raw experimental data prior to searching for axion sig-

Frequency [GHz]
FIG. 1. Mass and photon coupling parameter space of axionlike particles over the QCD axion dark-matter region.This plot is current as of early 2017, prior to publication of the Run 1A limits [18].Limits produced on the axion-photon coupling are opaque [19][20][21], and the ADMX generation 2 discovery potential as projected in early 2017 is displayed in translucent blue [22].Benchmark couplings from the KSVZ and DFSZ models are given by the diagonal red band's upper and lower edges respectively.Previous ADMX limits and sensitivity projections assume a local axion dark-matter density of ρa = 0.45 GeV/cc.
nals, paying particular attention to noise characterization, backgrounds in the integrated power spectra, and the classification of unfit data.Section V presents the statistical models used to search for the axion, how multiple observations are composed to form the grand statistical spectra, and the process of converting sensitivity spectra to limits.Section VI reviews the different axion dark-matter models searched for, including the model used during data-taking operations and the more sophisticated models of the final analysis.Section VII presents the procedure for handling axion candidates during Run 1A, including initial identification using the live analysis output, generation of synthetic axions, rescans and further procedures, residual candidates, and the reintroduction of rescan data for the construction of the grand spectra.Section VIII consolidates the cumulative findings into a grand spectrum, providing the net Run 1A axion limits.Section IX summarizes the findings and sets the stage for Run 1B [23,24].

II. THE HALOSCOPE APPARATUS
This section provides a summary of the ADMX haloscope apparatus, specifically the magnet, resonant cavity, cold receiver chain, warm receiver, and data acquisition (DAQ) systems.The focus is on components directly relevant for the axion search.These and other supporting components are discussed in more detail in the recent ADMX instrumentation paper of Run 1A and 1B [25].
The main ADMX apparatus is an Earth-bound tunable haloscope, designed to detect relic axions via microwave power emitted from a resonantly enhanced inverse-Primakoff process.The experimental apparatus, shown in Fig. 2, is divided into a graduated cold space, which at its coldest level (∼ 150 mK) is regulated by a dilution refrigerator (DR), and surrounded by increasingly warmer spaces leading up to room temperature.The cold space lies primarily within the main magnet housing.The outer steel housings of the magnet cryostat act as shielding against heat and radio frequency interference (RFI) leakage [26,27].Inside the magnet bore housing is the detector insert, containing the cavity, cold electronics, and sub-Kelvin cryogenics.Two more layers of heat/RFI shields separate the coldest elements, including the resonant cavity, from the 4.2 K main magnet.Power from the cavity is transmitted through a tunable antenna in the top plate and out of the cold space through an ultra-low noise transmission and amplification chain.This paper concentrates only on the receiver chain used in Run 1A, which tracked a TM 010 -like mode in the cylindrical cavity.In the warm space, the transmitted signal is then mixed down in frequency and digitized as a voltage time series and a frequency power spectrum over a 25 kHz bandwidth, which roughly matches the bandwidth of the TM 010 -like resonance.These digitized power spectra are where the axion signatures would be expected to appear.
The warm space also includes all equipment outside the magnet shield, and applies to the majority of the cryogenics and gas handling infrastructure.Also at room temperature is the DAQ infrastructure that monitors and controls the components of the warm and cold spaces.The DAQ monitors and controls many physical components of the experiment, such as running the DR, operation of the quantum-limited amplifiers and mechanical control of the cavity, monitoring of the cryogenics, cavity, magnets, and a complex of sensors for magnetic field strength, pressures, cryogen levels and temperatures at every stage of the insert.The lowest layer of the DAQ software is based on EPICS [25,28], which provides a uniform software interface for interaction with the instruments.The DAQ infrastructure also executes the numerous routines used in data taking and monitors its status using a live analysis.This sub-section reviews the primary components responsible for converting ambient axions into microwave photons: the main magnet and cavity.2. Cutaway model of the ADMX apparatus.The "Insert" (below and including the top conflat flange) is comprised of three main sections: "the cold space" (below 4 Kelvin), the "reservoir" (4-12 Kelvin), and a 12-50 Kelvin space.The cold space is primarily cooled by a DR below 1 Kelvin (K) and also contains the cavity and radiation shields to prevent components from being radiatively coupled.The "field-free region" containing cold electronics components sensitive to high magnetic fields is also in the cold space as it is thermally sunk to the DR.Sitting on top of the cold space is the helium reservoir.This portion is maintained at 4 K and houses a helium reservoir and other components that are thermally sunk to it.Above that is the warm space.This includes baffles, supporting structural components and a top conflat flange, above which are the only components visible during operations.The Insert sits inside a cryostat, which holds the superconducting magnet.The main magnet can generate fields exceeding 8 Tesla in the cavity.Labeled components will be elaborated on later in Section II.This figure is modified from a version that first appeared in [25].

Main Magnet
The static magnetic field is supplied by a superconducting magnet manufactured by Wang NMR of Livermore, CA [26,29].The magnet consists of niobiumtitanium windings immersed in liquid 4 He.During normal Run 1A data-taking operations, the field at the center of the solenoid was 6.8 Tesla, falling off to about 70% of this value at the center of the end plates of the cavity, as seen in Fig. 3.The superconducting coil is a 1.12 meter tall solenoid with a 60 cm inner diameter bore.The magnet winding is composed of four concentric superconducting solenoids, containing 99 km of copper-stabilized niobium-titanium wire wound around a stainless steel spool piece and potted in epoxy.The rated maximum field at the center of the magnet is 8.5 Tesla at a current of 248.96Amperes and a stored energy of 16.54 MJoules.The central field of the main magnet runs linearly with the current through the coils where I is the supplied current to the coils.The main magnet's current is continuously maintained by a power supply.Once cooled to operating temperature, the magnet requires a supply of approximately 2000 liters of liquid helium per month for continuous cooling during datataking operations.

Cavity
Axions passing through the field of the main magnet transition into microwave photons with some small probability.To capture and maximize this conversion, a highly conductive cavity is placed in the magnet bore at the center of the solenoid field.The ADMX Run 1A cavity is a 136 liter copper-plated stainless-steel cylinder with two copper tuning rods.The cylinder is made of a tube-section forming the cavity walls plus two removable end caps.The end caps form a low-resistance connection to the walls via a knife edge on the walls pressed firmly into the end plates.Both the stainless-steel cavity and the copper tuning rods are plated with OFHC (oxygen free high conductivity) copper with a minimum thickness of 0.075 mm and then annealed for 8 hours at 400 Celsius in vacuum.The annealing process further increases the grain-size of the copper crystals allowing for longer electron scattering path lengths as the copper is cooled and enters the anomalous skin depth regime.The quality factor Q of the tracked TM 010 -like mode in the unloaded cavity ranges from 20, 000 at room temperature to 160, 000 at cryogenic temperatures [29].The unloaded Q of the cavity primarily depends on the resistive skin losses of the copper plating.The oscillating electric fields of the cavity mode penetrate the resistive copper walls of the cavity.For more information on the construction of the Run 1A cavity, see Khatiwada et al. [25], Hotz [29].
The cumulative interaction strength between the ambient axions and electromagnetic fields over the interior of the cavity is given in natural units by the Lagrangian where V is the inner volume of the cavity, α is the fine structure constant, g γ is the dimension-less modeldependent coupling constant, f a is the Peccei-Quinn symmetry-breaking scale, a is the axion field, ⃗ E is the electric field strength, and ⃗ B is the magnetic field strength.For a highly resonant cavity threaded by a strong magnetic field, passing relic axions may be seen to act as an external driving force on cavity wave-forms with non-trivial Lagrangian L cav so long as the axionphoton coupling is feeble enough that the probability of being converted is small.The dynamics of the system become resonant around the cavity modes as the boundary fields and losses per cycle are small in comparison to the stored energy.
The power deposited into the cavity from a monotone source averaged over times much longer than the coherence time of the cavity mode, t coh = Q L /ν 0 , where ν 0 is the mode center frequency and Q L is the loaded quality factor, in SI units is found to be Here ϵ 0 is the permittivity of free space, c is the speed of light in vacuum, C mode is the form factor quantifying the mode and magnetic field alignment, and T ν0 (ν) is the mode envelope shape, which is expected to follow a Lorentzian form The cavity form factor parametrizes the overlap of the magnetic field and electric field in the cavity The magnetic field under is dominated by several orders of magnitude by the external main magnet under normal operating conditions, which has profile ⃗ B 0 .The electric field in the vicinity of a single cavity mode is dominated by the mode's form ⃗ E ξ where ξ is the multi-index of the mode wave-form, which will often be parameterized by modified cylindrical harmonics (n, l, m, X) where the last index X will be used to distinguish modes split by the axial-symmetry-breaking tuning rods.The form factor of the ξ mode is then well approximated by The main magnet's field is oriented vertically along the cavity's axis, though it does diverge some at either end, see Fig. 3, making the value of the form factor primarily dependent on the shape of the mode's axial electric field.For the ADMX cavity, there are, among others, transverse electric (TE) modes and transverse magnetic (TM) modes.The TM modes are the only category that contain an axial electric field desired for large form factors.The axial electric field for a TM nlm mode of an empty right-circular cylinder is where E amp (t) is the time dependent component of the field, J m is a cylindrical Bessel function, x ml is the l-th root of J m (x) = 0, R is the cavity radius, and d is the cavity height.For the rod-less cavity and magnet configuration used in ADMX, the chosen TM 010 -like mode maximizes the form factor C m .Two copper-plated stainless steel tuning rods are placed in the cavity in order to tune the axion-coupled modes.The rods run through the length of the cavity and are mounted on rotating alumina oxide rotary armatures.The rods are 0.05 m in diameter [25,29].The rods ideally create null boundary conditions in the electric field of the cavity mode, effectively shrinking the extent of the cavity in the horizontal plane, increasing the TM mode frequencies, and splitting modes by breaking the axial symmetry.The presence of the tuning rods also reduce the quality factor from that expected of an empty cavity by about a factor of two.Identifying the split TM 010 -like modes can be done heuristically or through simulation by increasing the thickness of the rods from zero (bare cavity) to their physical dimensions [30].The TM 010 -like mode branch tracked in Run 1A is identified with the name TM 010c in Lyapustin [27].Rotating the rods from near the wall of the cavity to near the center tunes the TM 010 -like mode from roughly 580 MHz to 890 MHz.The rods are moved by stepper motors located on top of the experimental insert.The stepper motors are connected to long G10 shafts that drive gearboxes on top of the cavity, which in turn are connected to the alumina shafts of the tuning rods.The gear boxes consist of two anti-backlash worm gear reductions, both geared down by 140:1 for a total reduction of 19600:1.Combined with the precision of the stepper motors, the angle of the rods can be stepped at the level of micro-radians.In practice, this stepping resolution allowed for tuning the cavity at 100-200 Hz per step during Run 1A operations.
The motion of the rods is most effective in tuning the TM modes, with TE frequencies remaining nearly constant.However, when the TM 010 -like mode crosses another mode, their waveforms mix and share energy.This mixing degrades both the form factor C nlmX and quality factor and severely complicates the flow of power.As such, the axion search becomes insensitive at these crossings.The two-rod system can be maneuvered to avoid many of these crossings, but not all.For Run 1A, one of the rods was put in the wall position while the other was left free to tune.Frequencies of the cold space circulators and the MSA (microstrip SQUID amplifier) coincided with a range of good form factor between two major mode crossings, as seen in Figs. 4, 8.The form factor is simulated over this range using FEM (finite element method) multi-physics software Comsol [25].
The cavity contains three microwave ports: one "weak" port on the bottom plate and two tunable "major" ports on the top plate.The ports allow for the extraction or insertion of RF power with the cavity.The weak port is a fixed, short antenna that is extremely under-coupled to the TM 010 -like mode.The weak port Errors (shaded regions) in the form factor and quality factor were computed to be 6% and 2.2% respectively, as discussed in Section IV E and presented in Table I.
can be used to inject signals into the cavity without degrading the loaded quality factor of the mode Q L or extracting signal power.The major ports consist of movable antennae that can be inserted or withdrawn via stepper motors attached to linear gear drives to vary the coupling strength.Like the tuning rods, the major ports' linear drives are actuated by stepper motors attached to long G10 shafts which turn 140:1 anti-backlash worm gear drives.The intention of having two major ports in Run 1A was to perform two simultaneous axion searches on two TM modes, a TM 010 -like and a TM 020 -like.However, for reasons covered in Khatiwada et al. [25], only the TM 010 -like mode was observed and the other port was uncoupled from the cavity.Only the TM 010 port and the connected receiver chain will be referenced from here on.A vector network analyzer (VNA) can be used to measure the TM 010 -like mode's frequency, quality factor, and the major port's coupling to the cavity.The VNA measures the swept transfer function (S 12 ) of the cavity by injecting a tone into the weak port of the cavity and measuring the same tone's amplitude transmitted through the major port.The frequency of the injected tone is swept across the band of the resonant mode.At frequencies far outside the central resonant frequency of the cavity, the injected tone is largely reflected by the cavity at the weak port or absorbed by the cavity walls so almost none is extracted at the major port.On resonance, the injected tone enters the cavity, excites the mode, and its power is extracted at the major port with far less loss.The output of the VNA is the swept response across the cavity mode.The expected response of an unmixed isolated mode is the Lorentzian distribution of Eqn. 4, from which the central frequency and Q-width can be extracted from a fit, seen in Fig. 5.
The VNA also determines the coupling properties between the major port and cavity via a swept reflection measurement (S 11 ).The VNA sends power towards the major port of the cavity through a circulator.A circulator allows for the directional injection of power along a transmission line.The power incident on the cavity reflects back with the form where β is the coupling strength parameter.The coupling strength can be expressed as where Q 0 is the cavity Q-factor with the major port uncoupled, and Q ext is the contribution to the quality factor from external losses such as the major port.For a given coupling, the response is total reflection off resonance and a dip on resonance where power is absorbed by the coupled cavity.Critical coupling occurs when β = 1 and all power passes into the cavity at the central frequency.A good impedance match is marked by a deep trough in the reflected baseline on resonance.The depth of an antenna is adjusted to maximize trough depth for the TM 010 -like mode as proxy for critical coupling.More precise analyses of coupling strength that can track overcoupling and other conditions have been implemented in subsequent runs [24].Conventionally, when the difference between the minima of the trough and the off-resonance baseline reaches −30 dB the antenna is considered critically coupled.Such a trough means that only 0.1 percent of the on-resonance incident power is reflected.External losses from the major port lower the mode quality factor as where Q L is the quality factor of the loaded cavity.During critical coupling, half of the power is lost to the major port, meaning Q loaded = Q free /2.Run 1A is configured to run at critical coupling, meaning at best only half of the axion power generated in the cavity is expected to leave through the major port.

B. Run 1A Cold Receiver Chain
The TM 010 major port connects to the RF receiver chain, as seen in Fig. 6.The chain begins with a low-pass filter and runs through the cold space where the quantum-limited electronics are housed.A bucking magnet surrounds the quantum amplifiers and other electronics sensitive to stray magnetic fields, actively canceling the field from the main magnet to tens of Gauss.Two hall probes are located in the field free region of bucking coil to confirm the field is cancelled to within a few Gauss [25].The electronics inside the field-free region are contained in an OFHC copper frame called the "squidadel" containing the cryogenic RF electronics and quantum amplifier package.This includes quantum noise limited amplifiers, circulators, switches, and temperature sensors.Physical and noise temperatures of the cryogenic electronics housed in the squidadel largely determine the noise temperature of the system, therefore the squidadel is kept thermalized to the dilution refrigerator mixing chamber, the coldest part of the system.
The small power emitted from the cavity passes through the first series of switches and circulators to the first stage quantum amplifiers and further to the HFET amplifier before passing into the warm space.The MSA boosts the signal with a characteristic gain of 20 − 25 dB, followed by an HFET amplifier in the 10 K space providing a boost of 30 dB.The squidadel is wired completely with copper coaxial cables whereas wiring from main port antenna to the first stage quantum amplifier is NbTi.Coaxial cables in the input chain are stainless steel.Characterizing the early-stage electronics is an extremely important procedure in determining the sensitivity of the experiment and is covered in Section IV as well as in [25].The total gain in the cold space is 50 − 55 dB.

C. Warm Receiver Chain
The remainder of the receiver chain runs through the successively warmer cryogenic spaces and vacuum feeds into the room temperature space.The total roomtemperature amplification is approximately 40 dB, bringing the overall expected power to the pico-watt level.This power is then directed to a variable-frequency superheterodyne receiver for digitization, as seen in Fig. 7.The power and gain of the receiver chain with the MSA switched out of the system were measured as a function of frequency every few weeks and found to change by a negligible amount, below the 1% level.
Treatment of the insert RF emissions into a digitized data set is performed in several steps, see Fig. 7.The first step is to mix the signal such that the mode center frequency is centered at 10.7 MHz via a local oscillator set to f = ν 0 + 10.7 MHz.The down-mixed voltage is bandpass filtered in a 30 kHz window centered at 10.7 MHz, which is expected to cover at least the full width at half maximum of the TM 010 -like mode.The analog signal now exists only in this vicinity of 10.7 MHz and is capable of being sampled quickly enough to resolve structures well below the expected total axion signal width.Time-wise sampling then occurs at a rate of 400 Mega-FIG.6.Schematic of cold RF receiver for Run 1A, first shown in Khatiwada et al. [25].Labeled components include circulators C1 and C2 to guide cavity bypass and reflection measurements, a low pass filter to limit incident thermal radiation at higher frequencies to be emitted by the receiver, switches (S * ), MSA, DC block, variable resistor (hot load), and the HFET (heterostructure field-effect transistor) amplifier.All components are connected with coaxial cabling.Switches allow the characterization of the system noise by toggling between a heated "hot load" and the main cavity (S1), and bypass the MSA (S2 and S3).Other pathways are for auxilliary measurements such as mode mapping (2) and MSA tuning (3), and are not used for axion observation.
samples/s with a 10-bit digitizer [25].To optimize the resolution/precision function of potential analyses, digitized samples are partitioned into bins 8 samples wide and averaged.This down sampling and averaging improves the signal to noise of samples by √ 8 ∼ 3 times and increasing the bit depth by log 4 (8) = 1.5.The resulting re-binned sample has an effective width of 25 MHz, well above the 2 × 10.7 MHz Nyquist rate of the central frequency.The re-binned time series are split into ≈ 10 ms blocks and temporarily stored to a circular RAM queue.No safety mechanisms exist to preserve the oldest unprocessed sequences, which are overwritten by the newest recordings.If overwritten, those data are lost and a flag is raised in the integration's metadata, indicating that scans processed are not necessarily consecutive.This becomes important for the high-resolution stored data.Once the 10 -ms series blocks are pulled off the RAM buffer, a Fast Fourier Transform (FFT) is performed on the concatenated series to 12 MHz resolution, near the Nyquist limit.With the well-resolved contributions about the resonance centered at 10.7 MHz, a digital mixer centers a 30 kHz band and shifts these contributions down to begin at 0 Hz, removing contributions above 30 kHz with a low-pass filter.We now have drastically down sampled the coherent frequency-space data set of a 10 ms sample at a rate resolving the full 30 kHz.The high resolution time series data is formed by performing and inverse FFT on each 10 ms sample and concatenating the time series in order.The medium resolution (MR) power spectral density, or power spectrum, of the scan is calculated from the modulus squared of each 10 ms sample, averaged over the 100 second scan (∼ 10 4 samples).The power spectrum array then has its edge bins removed, producing the 256 bin wide form of the raw spectrum saved for the MR analysis.The width of each bin the MR spectrum is ≈ 100 Hz, the Nyquist limit for each 10 ms sample.It is in this power spectrum that the axion dark-matter signal would appear as a localized excess.

D. Experiment Status Measurements
The warm receiver chain and digitizer are located within the larger structure of the ADMX DAQ system.A number of measurements are recorded by the DAQ to monitor the state of the insert, magnet, and ancillary cryogenics.Interpretation of these experimental conditions plays a critical role in performing the offline analysis.Here we briefly describe the measurement and interpretation of these experimental conditions.A more detailed account can be found in Khatiwada et al. [25].
The data recorded from the experiment can be divided into periodically sampled experimental state information and radio-frequency measurements taken in the course of the axion search.The experimental state information consists of status readings from temperature, pressure, field, and current sensors.
For temperatures above 1 K, an assortment of resistance sensors are used to read out temperatures at various thermal stages.For temperatures below 1 K, temperature sensors are sampled with a four-wire resistance measurement using a Lakeshore Alternating Current (AC) resistance bridge [25].The temperatures of the cavity and quantum electronics package are measured using Cernox temperature sensors and the temperature of the mixing chamber mounted to the cavity is measured using a Ruthenium Oxide sensor [25].During operations, it was observed that the Cernox sensors had a large magneto-resistance at temperatures below 1 K.With the magnet ramped to 6.8 T, The Cernox temperature sensors on the main cavity read 70% higher temperatures compared to the magnet ramped down, while the Ruthenium Oxide temperature sensor on the mixing chamber increased by 2%.Thus, in Run 1A, the temperature of the cavity was read by the Ruthenium Oxide temperature sensor mounted to the mixing chamber.Because the quantum electronics package was kept in a field-free region, Cernox temperature sensors located on the package did not suffer any appreciable effects from the magnetic field, and were used to measure the physical temperature of the quantum amplifier.
The main magnet state is captured by several sensors on the magnet's power supply as well as Hall probes that directly measure the magnetic field parallel to the probe wire, which are set in an azimuthal (vertical) orientation in areas of high magnetic field as well as in low-field/highsensitivity areas like the so-called field-free region.The power supply monitors the voltage and amperage being fed to the main magnet at a sampling rate of a few minutes.As stated by the manufacturer, the peak magnetic field for an empty bore can be calculated from Eqn. 1, though minimal impact from insert materials and bucking coil are expected near the cavity.The Hall probes record the magnetic field with a period of about an hour.The magnetic field throughout the cavity may be modeled using the combined data and then used to compute the form factor.

III. DATA-TAKING OPERATIONS
Data-taking operations for Run 1A occurred between January 18 and June 11, 2017.During that time, the TM 010 -like mode of the main cavity was tuned over the range of 645 − 680 MHz and observed for power excesses consistent with an axion DM signal.The identification and handling of candidates was also performed during this period.This section decomposes the data-taking process into its constituent parts and cadences down to the level of a single RF digitization, overviewed in the previous section.

A. Overall Structure of data-taking
The global imperative for data-taking in Run 1A was to cover the viable frequencies at the intersection of the cavity and receiver chain operational ranges.Once at sub-Kelvin temperatures, the MSA was established to be the limiting device in the receiver chain, providing the target gain over the range 645 − 680 MHz, quickly decaying for lower and higher frequencies.This range coincides with a continuous stretch of the TM 010 -like cavity mode with high form factor C, high quality factor Q, and without significant mode crossings.The frequency range is accessed through tuning rod configurations where one rod is held at the wall and the other is free to turn, as seen in Fig. 8 6 is used as the input of the warm space amplification and digitization chain.The raw data from a 100 s integration is recorded in two forms: as a power spectrum with ≈ 100 Hz wide bins and as a 25M element voltage time series.
dividually scanned to DFSZ sensitivity, investigated for candidates, and cleared of candidates before moving onto the next segment.The range is modularized for faster turn-around in the case of a prominent candidate, and for feedback on the state of the MSA.The scan process of each segment is structured in a sequence of up-tuning and down-tuning sweeps of the TM 010 -like mode over several weeks to provide uniform coverage.Multiple sweeps separated in time by days or weeks also provides an opportunity to sample possible periodic and transient behavior of the axion field.This cadence is only rarely broken in the case of a promising candidate.
Once the scans were complete to DFSZ sensitivity according to the live analysis, the procedure transitioned into candidate handling.DFSZ sensitivity is established based on criteria covered in sections Sec.V.The investigation of candidates and the handling decision tree was established in the Candidate Handling Protocol covered in Sec.VII.Once candidates for one segment are handled, the process begins again for the next segment.

B. Cadence of a data-taking Cycle
The process of scanning a segment is broken up into a smaller cadence called a data-taking cycle.The datataking cycle lasts for approximately 150 seconds and in-cludes operations of moving tuning rods, active cavity measurements, and the passive integration of RF chain emissions.The data-taking process at and below this level is automated and was often left in continuous and unhindered operation for days barring a necessary manual bias of the MSA.
At the head of the data-taking cycle are several active measurements of the cavity and receiver chain using the VNA.The injected sweep signals are far more powerful than a potential axion signal or impinging RFI, making them easier to see, while still conforming to the receiver operating parameters optimized for the passive integrations.The first measurement is a swept transmission S 12 from the weak port through the receiver including the cold and warm amplifiers.This measurement shows the cavity response, yielding a measurement of mode frequencies and quality factors.S 12 measurements take approximately 10 seconds each.A wide-band transmission scan is made every 10-th scan cycle for mode mapping purposes.The transmission measurement is followed by a reflection measurement taking approximately 20 seconds, where the signal is sent through the bypass line and is directed towards the cavity by a circulator (C 1 in Fig. 6).The reflection sweep near the resonance is mostly absorbed by the cavity, while off resonance the signal is reflected and passes back through the cold and warm amplification of the RF system.This measurement yields a wide-band measurement of system gain for noise calibrations, the coupling between the antenna, and the cavity mode of interest, as well as the Q L of the cavity mode of interest.Coupling of β = 1 is consistent with critical, or impedance-matched, coupling.The duration of transmission and reflection measurements have been reduced to less than a second each in subsequent runs by optimizing input attenuation and power [24].
The receiver integrates emissions from the cold RF for the remainder of the data-taking cycle, which comes to an observation duty factor of approximately 0.7.At the end of the data-taking cycle, there are two mechanical tuning processes to prepare for the next integration.The tuning rate of the rod can be estimated given a target sensitivity and regional values for the system temperature, quality factor, etc., or can be set manually to a given number of steps per cycle.The tuning rate of the warm stepper motor was set between 0.1 − 0.2 radians per cycle during the first passes of a section's bandwidth, which translates to 1 − 2 kHz per cycle.The main port antenna is also given the opportunity to tune to alter the coupling to the main cavity.Note that rod and antenna tuning must occur slowly to avoid overtaxing the dilution refrigerator.

C. Logging of Search Data
Each data-taking cycle is stored as an independent entry in an SQL database on the DAQ's main control computer.Each cycle contains its own unique serial number and timestamp when the cycle is initiated.The integrated power spectrum and a collection of markers are added to the cycle entry to more easily identify the conditions surrounding the experiment.Each entry contains the sum total information necessary to perform an axion search.

IV. DATA PREPARATION FOR AXION SEARCH
With data collected, preparations are then made to characterize the state of the experimental apparatus and assess the quality of the measurements to search for the axion.This includes the verification of the cavity magnetic field, characterization of the receiver chain as measured by its effective noise and gain properties and its persistent background structure.This section details how the measurements necessary to axion search are interpreted from their recorded state into actionable data.

A. Direct Measurements
Temperatures within the experiment are measured by arrays of sensors placed at every level of the insert and on the main magnet casing, both in areas of high magnetic field and in the field-free region.Voltage time series from the sensor leads are interpreted through EPICS and converted into temperature readouts.Sensors of differing makes and models were tested against one another in and out of the magnetic field in order to study the uncertainties and biases present during data-taking operations, as detailed in [25].The errors in temperature sensors are reflected in Table I, and are integrated in the final sensitivities and limits.
The magnetic field in the main cavity is computed from the current supplied to the main magnet and confirmed by the Hall probes placed throughout the insert.The maximum magnetic field in the solenoid center is computed from the current via Eqn. 1 and modeled in form by numerical simulation as indicated in [25], including the counter field induced by the bucking coil surrounding the field-free region.
Active transmission and reflection measurements taken during each data-taking cycle are analysed to assess the impedance match between the cavity and main port, to match the transmission function through the port to a Lorentzian distribution, and to extract the quality factor and mode central frequency.During Run 1A, the weak port was partially dislodged and became decoupled from the cavity during the insertion and commissioning process, marginalizing the effectiveness of the cavity transmission measurements.As a result, S 11 reflection measurements were used to assess coupling β and loaded quality factor Q L of the TM 010 -like cavity mode.Recall that the normalized reflection power spectrum is modeled by the form given in Eqn. 8.The logged S 11 measurement is fit to this form by a least-squares analysis of the swept spectrum and the parameter values are used to assess the loaded Q L and coupling β of the cavity-main-port state, see Fig. 5.

B. Axion Power in Cavity and Main Port Transmission
The figures central to the deposition of axion power into the cavity and transmission into the receiver can now be computed.The power-per-axion density deposited in the cavity at a single frequency is given by The form factor of the TM 010 -like mode C 010 is computed from Eqn. 6 using numerical models of the magnetic field and mode electric field distributions [25].The volume of the cavity is computed from [26] and is invariant to high precision during low temperature operations.The maximum magnetic field is computed from Eqn. 1, and has error calculated at the level of the power supply current stability.The numerical error of B 2 max × C 010 × V is available in Table I.
The loaded quality factor Q L is computed from the reflection measurement detailed in the previous sub-section and uses a rolling average and uncertainty from the previous 10 measurements.The portion of the power transmitted from the cavity into the main port, and past the first low-pass filter is modulated by the coupling strength factor β, and the low-pass filter transmission function T f ilter where the filter is near-transparent at the Run 1A frequencies (T f ilter ≈ 1).

C. Receiver Gain and Noise
The remainder of the receiver chain transmits the cavity output, but also imprints its own structure into the digitized power spectrum down to the scale ν o /2Q L ∼ 10 kHz (see [25,31]), an order of magnitude wider than the expected axion signal width of ≲ 1 kHz.Directly modeling the total transmission function prior to Run 1A proved to be too unreliable due to high variability in the response of devices and strong inter-device couplings.This subsection analyzes the receiver transmission heuristically to characterize its gain structure and noise.
The noise background is expected to be overwhelmingly thermal.Noise in the power spectrum is first contributed by the fluctuations about the mean photon occupation function of the cavity ⟨n γ (ν)⟩ = 1 e hν/kBT − 1 (13) where k B is the Boltzmann constant, T is the physical temperature of the thermal source, h is Plank's constant, ν is the excitation frequency.The expected power spectrum per frequency of the cavity is derived from the mean energy spectrum where Z is the canonical partition function.Note that the first term of the mean energy spectrum is the vacuum contribution to the cavity energy and the second term is the contribution from non-trivial occupation.The emission of photons out of the cavity, ignoring reflections and attenuation for now, is then given by the free flow of power out of the strong port.The expected power spectrum per frequency is In the limit of the photon energy much less than the bath temperature hν ≪ k B T , the expected spectrum density flattens and the power spectrum becomes where b is the integrated bandwidth and ν is taken as the center frequency.
The distribution of power fluctuations can be found by looking at the occupation probability of individual states, given by the density function where r is the occupation number of the state(s) with frequency ν.The emission rate of photons over the bandwidth b is expected to be ∼ k B T b/hν.This comes to ∼ 16 bandwidth-emissions for a typical Run 1A temperature of T = 500 mK at a center frequency ν = 660 MHz.Each short continuous δt ≈ 10 ms integration (bandwidth b = 1/t ≈ 100 Hz) would then produce an MR emission power spectrum in an exponential distribution of similar shape to Eqn. 17.The rate parameter of an exponential distribution is identified with the variance square root λ = µ = σ where the mean of the power spectrum is already known to be µ = ⟨P n ⟩ = k B T b from Eqn. 16.
The ∆t ≈ 100 second integration period of each data taking cycle is made up of n∆t/δt ≈ 10 4 MR power spectra, which are considered identically distributed over that time scale and independent as the integration time exceeds the timescale of power fluctuations δt > 1/2b.Averaging those spectra will shift the distribution shape from an exponential to a normal shape in the large-n limit according to the central limit theorem, leaving the mean unchanged but re-scaling the standard deviation by a factor of 1/ √ n.Therefore, we find the distribution of emitted power fluctuations to have standard deviation which matches the Johnson spectra [32,33].
The width of the noise distribution is proportional to the thermal temperature, as is the mean thermal power.One can construct a bin-wise signal-to-noise ratio as a proxy to the sensitivity of the instrument to detect a signal of known power The distribution of power fluctuations through the remainder of the receiver chain is expected to remain thermal, therefore a system noise temperature is used to characterize the net distribution of noise in the receiver.The details of computing the system temperature can be found in [25], but we explain them briefly here as they will enter into the gain structure characterization in the next sub-section.
Amplifiers and attenuators in the receiver chain modify their inputs, both signals and noise powers, as well as inject noise power according to their own blackbody spectrum.A simple composition model for the power emitted at the end of the chain is where P ni and G i respectively are the noise power and the gain of the i-th component of the receiver.Assuming thermal power spectra at each stage (⟨P ni ⟩ = k B T ni b), the total system noise can be computed as a temperature by dividing out the total gain of the receiver If the gain of the first components are large, then one can effectively truncate the series after the first several terms.Note that the contributions of the power spectrum fluctuations also follow this relation The gain of the MSA varied with both its center frequency and physical temperature throughout the main data run, as did the HFET amplifier to a lesser degree, requiring periodic calibration to optimize the system temperature.Contributions from amplifiers and attenuators downstream of the HFET were measured to be effectively constant.Also during the run, the switch for the heated load malfunctioned, so the primary noise calibration came from an "on-off resonance" method, described below and in Du et al. [18], Khatiwada et al. [25].
The physical temperature of the cavity in Run 1A (∼ 150 mK) was significantly different than the physical temperature of the milli-Kelvin electronics (∼ 300 mK).This difference implied that the relative thermal power on and off resonance encoded sufficient information to determine the system noise in the same way a heated load measurement does.The method for calibrating T sys during Run 1A then became: 1. Use a "hot-cold" measurement as the primary calibration to determine T sys absolutely at a few fixed frequencies, as well as the noise contributions from the RF components downstream of the MSA 2. Use a faster "SNR-Improvement" measurement to determine T sys on every digitization The ideal hot-cold method measures the power at the receiver in a bandwidth of multiple Q-widths about the resonant frequency of the cavity while the antenna is critically coupled.The powers measured on and off resonance are then compared to determine the noise contributions of the MSA amplifier from: R = T attenuator + T amps T cavity + T amps where T attenuator is the physical temperature of attenuator A as seen in Fig. 6, T cavity is the physical temperature of the cavity, T amps is the noise contribution from all the electronics of the receiver chain, and R is the ratio of power off-resonance to the power on-resonance.Temperatures for the cavity and attenuator were taken from recordings of the Ruthenium Oxide temperature sensors closest to it [25].The R factor can be determined by a model of the cold RF system and fitting a model to the power spectrum as a function of frequency, discussed in [18,25].Parameters had to be changed slightly to accommodate the typical conditions for this in-situ measurement for Run 1A.Off-resonance noise power is dominated by the attentuator (physical temperature typically ∼ 300 mK) and also contains a contribution from the receiver noise temperature.On-resonance noise power is the sum of the cavity (physical temperature typically ∼ 150 mK) and the receiver.The transition between the two creates a dip of order 20% in power seen at the cavity resonance in Fig. 9.This measurement yielded a system noise temperature typically of ∼ 500 mK, but was found to vary substantially throughout the run due to different gains of the MSA, which required frequent re-optimization as discussed in Khatiwada et al. [25].
Example power measurements used to calibrate the system noise temperature Noise power on versus off resonance acts as an effective hot-cold load.The asymmetry of the shape is a result of interactions between RF components.This figure first appeared in [18].

D. Receiver Shape Removal
As was seen in the previous sub-section, the spectra from a well-equilibrated receiver chain sans persistent signal is expected to be flat, characterized by a single system temperature.This was rarely the case in practice for Run 1A.There are multiple causes of a non-flat power spectrum by experimental hardware, including frequency dependant gain variations before and after mixing down the target frequency, and frequency-dependant noise variations.The last of these is expected to be suppressed as the early receiver chain is nearly homogeneous and stable in temperature over a data taking cycle.
It is crucial to know the structure of the gain so that it is not confused with potential axion signals.Removing the structure and flattening the receiver power spectrum allows for a more straightforward interpretation of the data as Johnson spectra dominant.The flattened spectra also has computational advantages for the analysis as the gain response of an incoming axion signal is similarly flattened, producing a convolutional filter form to the optimal signal search, a topic that will be covered in Section V.The operational goal of receiver shape removal is to completely remove the frequency dependant receiver chain response while retaining potential signals and thermal components.
Finding the true gain response of this complex system is a difficult task to perform from status measurements of components alone, and has met with limited success [31].Techniques rooted in purely heuristic fitting of the power spectra have had far greater success [31,34,35].Run 1A used filters of Savitzky-Golay (SG) type and a low order polynomial fits for the off-line and live analyses respectively.

Polynomial Filter
A sixth-order polynomial is used by a live analysis process run in parallel with data-taking operations to fit the receiver structure for computational speed.The live analysis is responsible for providing a real-time feed of the experiment's sensitivity and for identifying potential candidates.
The polynomial fit works well for spectra that are already nearly flat, but its quality degrades significantly over spectra containing larger and more complicated variations, as seen in Fig. 10.There were multiple periods during Run 1A, either due to the operating state of the MSA or other elements, where the receiver structure was highly perturbed.A low-order polynomial fit does a poor job of conforming to large scale fluctuations at multiple harmonics [30].This shortcoming resulted in some background-subtracted scans with residual structure above the noise floor, resulting in an excess of axion candidates that impacted rescan operations.The candidate procedure will be discussed more in Sec.VII.

Savitzky-Golay Filter
The SG filter was used to fit the receiver structure after data-taking operations were complete to conduct the analysis that provided the limits in Du et al. [18].The SG filter proved to be much more versatile over the range of structures experienced in Run 1A, as seen in Fig. 10.The parameters used for the SG filter were d = 4 for the spline order and L = 121 for the size of the box window function.These values were chosen as they produced adequate fits to a wide range of observed backgrounds with minimal impact to potential axion signal shapes.
The general thinking is that you want to choose filter parameters that fit typical background shapes well, but are unable to fit signal shapes well.I'm not sure I can give why precisely those numbers were used, but my recollection is that they were tuned to minimize our 'fudge factor', the inefficiency (as measured on synthetic signals) introduced by background subtraction.
The SG filter has been found to attenuate axion signals and imprint small negative correlations between processed spectrum bins [34].These anti-correlations become problematic for scans repeated over the same frequency range under the same conditions, such as may occur for a re-scan, as seen in Fig. 11.Accumulation of anti-correlations will occur very rarely under normal scanning operations.No new candidates were produced under the SG filter.

Preparation of a Spectrum
We may now prepare the fitted scan for axion search.Using the findings of the previous sub-section, we can see that dividing the power spectrum by the background Anomolous Background FIG. 10.Sample fits to two observed receiver spectra by a six-order polynomial and Savitzky-Golay windowed spline.The two raw spectra represent the typical (left) and anomalous (right) imprints of receiver structure.The polynomial is adequate to fit small and low harmonic structure to the noise floor, but fails to respond to large structure.Savitzky-Golay is seen to provide a more consistent fit over the range of experienced structures in Run 1A, though artifacts of the spline fit can be seen in the form of fluctuations for scales at and below the window size.Such fluctuations lead to small negative correlations in the background-subtracted spectra.shape produces a dimensionless mean-normalized spectrum with unit mean and, in the absence of residual background structure, Johnson-distributed fluctuations with uniform variance over the entire frequency band.We are only concerned with power excesses above the mean as our axion signals are of narrow bandwidth, making the relevant spectra The dimensionless fluctuations are of random normal distribution with fractional width relative to the mean calculated to be which is expected to apply across the spectra as long as the receiver components up to the first stage amplifier are in thermal equilibrium and the gain structure of the receiver thereafter is much wider than the scan bandwidth.This means one can compare the fluctuation statistics across an MR scan, and between MR scans so long as the scan times match.The dimensionless fluctuations show overwhelming normal distribution shape, as seen in Fig. 12.
To find the unit-full power excesses, recall that the mean of O P is identified with the average thermal power out Multiplying the fluctuations by the mean power produces the normalized fluctuation power [18] δ w = δ × P .
These are the power fluctuations against which the search for axion emissions from the cavity is performed.
The SG filter was able to fit the background well enough so as to not change the noise width more than a factor of two from the theoretical expectation of Eqn. 25 for 96.6% of usable Run 1A scans, compared to 90.6% for the polynomial filter, providing a more complete sampling of the noise distribution, which was found to be overwhelmingly Gaussian normal in shape, as seen in Fig. 12.

E. Errors and Uncertainties
This sub-section provides the systematic uncertainties for parameters used in the analysis, summarized in Table I.First, the uncertainty on the quality factor was quantified by repeatedly measuring the quality factor in a narrow range of frequencies, in this case, from 645 − 647 MHz, where the quality factor was not expected to change much according to models.The fractional uncertainty on the quality factor in this range was determined to be 2.2%.The fractional uncertainty on the main port coupling was also computed over the same frequency range, and determined to translate to a transmitted power uncertainty of 0.55%.The uncertainty of the individual temperature sensors above 1 K were taken from their data sheets and the sensors below the 1 K stage were computed as indicated in Sec.II D and [25].The FIG. 12. Histograms of normalized power fluctuations for all scans used in Run 1A MR analysis after background subtraction using six-order polynomial (orange) and Savitzky-Golay (blue) filters, with green indicating overlap.The fluctuations δ are normalized to the thermal spectrum width.The Gaussian curve (red line) fitted to expected distribution of fluctuations shows the distribution is well fit at the center but is slightly deformed from a purely thermal spectrum, containing notable excess counts by ±3σ for the six-order polynomial background subtraction.The Savitzky-Golay background subtraction technique fares somewhat better, with excesses appearing closer to ±4σ from center.Injected power in the form of RFI is responsible for some portion of the positive power excesses, while much of the remaining large power excesses can be attributed to sub-optimal background fitting.
combined error on the system temperature T sys is calculated at 7.1% by adding in quadrature the individual error contributions from the components in Eqn.21, being dominated by the cavity (Scientific Instruments RO600 Ruthenium-Oxide sensor on the mixing chamber) and first stage amplifier (Lakeshore CX-1010 Cernox sensor on the squidadel).The uncertainty of the receiver noise temperature, and therefore the system noise, was primarily given by two contributions: the uncertainty in fit for the HFET noise, and the uncertainty in fit of the MSA noise.The total fractional uncertainty for the system noise amounted to 7.5%.Further details can be found in [25] and the supplementary material of [18].The final systematic uncertainty of ±13%, as shown in Table I, was computed simply by adding all listed uncertainties in quadrature as the errors are assumed to be independent from one another.

F. Data Cuts
Not all data taken during Run 1A operations is qualified to enter into the axion search.Scans taken as part of  I. Dominant sources of systematic uncertainty.The uncertainties were added in quadrature to attain the uncertainty on the total axion power from the cavity, shown in the bottom row.For the first entry, B is the magnetic field, V is the cavity volume, and C is the form factor.The RF model fit and temperature sensor calibration discussed in [18].The total uncertainty from these independent sources was computed from the quadrature sum and is incorporated into the final limit computations in Section VIII.
• SAG tests, • un-subtract-able RFI, • mode navigation errors, • abnormal cryogenics conditions, • and incomplete logs are flagged for omission.These flags are either explicitly written to the scan or logged by hand in an electronic log that is updated throughout experimental operations, upgrades, and day-to-day upkeep.Further cuts to the data were made during both the live and off-line analyses due to derived knowledge of the experiment being in a poor, or poorly understood, state.This subsection presents the conditions used to cut data from the axion search analysis and its impact on the number of viable scans.The Run 1A data set consisted of a total of 173,048 scans and raw integrated power spectra: 138,680 for preliminary scans and 34,396 for leveling and rescans.After implementing the analysis cuts described above and itemized in Table II, 78,958 scans remained for the axion search analysis.Motivation for these cuts were as follows.First, quality factors lower than 10, 000 or greater than 70, 000 were omitted from the analysis as they were determined to be either compromised by a mode crossing or non-physical and the result of a poor fit to the reflection scan.System noise temperatures below 0.2 K and above 5.0 K were excluded as these were likely a poor fit by the SNRI or other temperature fitting mechanism described above.Temperatures below the lowest physical temperature in the experiment of 0.15 K are flagged as non-physical.Additionally, the SG fit to the power spectra background in the offline analysis was required to have a fractional standard deviation relative to the radiometer equation between 0.5 − 1.2 among 95% of the least deviant points in the spectra.This proved sufficient to reject poor fits while retaining potential axion signals.

V. POWER EXCESS SEARCH
This section presents the statistical method used to search for persistent signals such as the axion in the MR integrated power spectra.Likelihood functions for power excesses in a scan, the axion signal, and conditions of the apparatus are first established, followed by searches on single data-taking cycles, which are then combined into a "grand search spectrum".

A. Statistical Modeling
It was established in the previous section that an integrated spectra's background noise is dominated by a thermal spectrum from the receiver chain where the receiver is expected to be in local equilibrium over the course of a single data-taking cycle.Under the integration and sampling rates set in the MR data as discussed in Section IV, the ∼ 10 4 -sampled Poisson-distributed thermal background power spectra when co-added ideally form a raw power spectrum with random normal distribution of mean power set by the product of the system temperature and the receiver's total gain T sys G tot,ν , and distribution width of σ = µ/ √ b∆t.Once gain structures and mean power have been removed, the remainder of the spectrum optimally consists only of the random-normal thermal fluctuations and externally induced power excesses.
The signal from the local axion field, broadly speaking, would present itself as a power excess of magnitude set in Eqn. 3 and an expected Q-width Q a ∼ 10 6 set by the virial velocity of the Milky Way.Note that the axion signal is present as a power excess as opposed to a deficiency as the conversion of a thermal photon to an axion via an inverse-Primakoff or other process is suppressed by an additional factor of the ratio between local photon and axion occupations n γ /n a .This implies that the statistical search for an axion signal can be straightforwardly conducted as a one-sided p-value test with underlying random normal uncertainties.
The likelihood function for an axion model hypothesis test is then taken to be where H a is the axion model hypothesis, D s is the data set, S is the action of the data and model hypothesis, λ is the data uncertainty measure, i is the index over observations, and N is the number of observations.It has already been shown that the integrated power spectra observations are nearly bin-wise independent in their noise background.The action of the data and model hypothesis will be taken over the power outputs from the cavity where P a is the expected transmitted axion power given by Eqn. 3 and P s is the recorded power spectrum.The expected axion power density P a is dependent on the transmission function of the cavity mode T νo , the cavity coupling β and quality factor Q L , magnetic field B, the observed frequency ν, the local axion density ρ a , and the overall shape of the axion frequency distribution function above the rest mass f (ν − ν a ).Note that the variance about the axion power is essentially vanishing in comparison to the expected power spectrum due to the high occupation values in the relic axion condensate, a topic that is discussed in Sections II, VI.The action is modulated by an uncertainty measure where σ P is the uncertainty in the recorded power spectrum and is constant over the scan.The uncertainty in a scan's power, is given in terms of the expected noise power uncertainty as computed in Eqn.25.The action and uncertainty measures form a χ 2 -like kernel for the likelihood, which takes the form where i is the index over bins in the digitized spectra, T ν is the cavity mode transmission shape, ⟨P ⟩ tot is the total power of the model axion signal, and p a is the probability distribution function of relic axions.To better track the coupling parameter, let us expand the model axion power as a fraction of the DFSZ benchmark model ⟨P ⟩ tot = A × P DF SZ , where A is the fraction of the model's power relative to the benchmark and P DF SZ is computed using Eqn. 3 sans the mode transmission shape.Refactoring the χ 2 figure in powers of A where the zeroth, first, and second power terms may be referred to as respectively.There are several priors that need to be declared in order to parameterize the space of incident conditions and hypothesis to be tested.This analysis operates under the following assumptions: • Relic axions dominate the local axion distribution.
• Relic axions make up 100% of the DM.
• The likelihood of the axion rest mass m a is uniform over the Run 1A frequency range.
• The likelihood of the axion-photon coupling g aγγ is uniform over all values covered by this search.
• Incident conditions of the experiment state are random normally distributed about a data-taking cycle.
The parameters of interest for the posterior distribution function are the intrinsic axion mass and the axionphoton coupling strength.One can marginalize over all other parameters above to reduce the likelihood function to be over only the parameters of interest.The tests here are organized first over axion masses, then the local distribution of the axion field, and lastly coupling strength.
Using axion power as proxy for the coupling, the posterior probability function for a given axion mass and local distribution reduces to a Gaussian with mean estimated by and width where .
The expected value of the axion power given an instance of data is computed from the posterior mean µ A .The sensitivity of the data to power excursions of the specified shape is set by the posterior width σ A .
The expected power's deviation from the null hypothesis of no axion (P a = 0) in units of the distribution width gives a measure of the significance for the potential fit and is the metric used to determine candidacy of axion signals.One may also construct confidence intervals over the posterior distribution from a data instance, though the significance of such statistics will not be made in this section.In Section VII we cover the procedure of identifying candidates and delineating their causes from statistical fluctuation, to RFI, to the relic axion field.The creation of limits after candidates are deferred to Section VIII.The next sub-section computes the above statistics on a single data cycle.

B. Single Scan Analysis
A single data-taking cycle is the natural choice over which to perform the axion search analysis outlined above.This subsection provides the computational details in searching for the axion in each RF integration.The following subsection presents the means by which the analyses are combined, extending the search over the entire scanned range.
The axion mass parameter is used to serialize the tests on a single scan.Axion masses are only sensible to be tested if there is an overlap between the power spectrum and the assumed axion distribution line shape, so we restrict ourselves to masses less than the uppermost edge of the scan frequency range (ν a < ν max ) and greater than the lower-most edge of the scan frequency range minus the width of the 99 % axion line shape (ν a > ν min −∆ν 99% ).The mass tests may be formed anywhere along this continuous segment of parameter space, though the utility of sampling below the MR bandwidth is marginal.The sampling of masses in the search for candidates is performed on the resolution of the MR bins, ∼ 100 Hz, with a starting point of 645 MHz.The final limits in Section VIII will be averaged over many adjacent mass samples.Now consider how each of the terms in Eqn.33 are calculated.The first term involves only the scan data and need be calculated only once for each scan regardless of the hypothesized axion mass.Its computation scales as the number of bins in the scan.The other two terms of the test statistic contain factors of the transmitted axion lineshape into the scan and are refactored for more efficient calculation over multiple axion masses.For both terms, the Lorentzian cavity transmission function T νo is factored out of the emitted axion power P emitted = T νo P a and is instead used to modulate the noise error estimate σ w = σ/T νo and the power spectrum P w = P emitted /T νo so that they now vary over the spectrum, see Fig. 13.The bandwidth of each digitized power spectrum is on the order of the Q-width of the TM 010 -like mode, less than 1 : 10 4 of the central frequency.It is by this same fraction that the supposed axion line shape will change its width across mass tests spanned by a single scan.To reduce the computational costs of generating a unique axion line shape for each probed mass, the line shape is generated once for the scan's central frequency (p a,i = bi dνf (ν a − ν)) and is approximated as translationally invariant over the scan.The factor of frequency in the axion power equation is also held constant at the scan's center point.The second term of Eqn.33 that computes the inner product of the power spectrum and line shape can now be phrased as a discretized convolution.Given proper zero-padding to the line shape and power spectrum, the intersection can be computed over the whole scan using the (discretized) convolution theorem where the Fourier transform is discrete and implemented using the fast Fourier algorithm.The third term of Eqn.33 is the overlap of the squared axion line shape with the power spectrum's modulated error, and can also be phrased as a discrete convolution, using similar zero-padding as the second term, Now one can compute each of the statistics in the previous sub-section to compute quantities such as the fit significance of axion hypotheses, Fig. 14.This use of convolutions to evaluate the data's senstivity to a given axion signal shape allows one to think of the search technique as an optimal filter.

C. Grand Spectrum
This sub-section integrates the statistics of single scan analyses into a grand spectrum of tests.This will be accomplished by deriving the arithmetic rules for the combination of single scan statistics for each figure of interest, then applying them to the set of viable scans.
Recall that the independence of each measurement allows the total likelihood function of Eqn.28 statistic to be decomposed as a product of single scan likelihoods.The total χ 2 statistic therefore can also be decomposed into a sum of single scans where N s is the set of viable scans.Further, the decomposition of the statistic in powers of the axion signal strength also factorize to where the i s index runs over the bins of scan s.Recall that we have organized all the tests according to prescribed axion rest masses and not bin frequency, minimizing the complications from scans with mismatched ranges and bin boundaries.The expectation and uncertainty statistics of interest to the grand spectrum are then seen to have the following addition rules where h s are the expectation value weights given by One can also subtract scans from the set by reversing the sign of the arithmetic operation.These are used to identify candidates, and ultimately provide the exclusion limits in Sections VII and VIII.

VI. AXION SIGNAL MODELS
An axion signal emitted from the cavity is stimulated from the oscillations in the axion field passing through during the integration phase of a data-taking cycle.This section presents the models used to emulate the signal generated by the ambient axion dark-matter distribution.This analysis takes a similar approach to that of Foster et al. [36] for modeling the response from a classical axion distribution, where the net axion field passing through the cavity during the observation is taken as a superposition of plane waves, which is considered accurate for a classical field on distances much shorter than the curvatures of an axion's orbit.The axion field inside the cavity is nearly homogeneous as its extent is on the order of the axion's Compton length, much less than the de Broglie length that sets the spatial fluctuations.The total axion field in the cavity is then where ρ DM is the local dark-matter density, N a = O(V ρ DM /m a ) is the number of axions passing through the cavity at any one time, E i is the energy of the i-th axion as measured from the cavity's rest frame, and ϕ i is the phase of the wave at the start of observation (t = 0).Most axions passing through the cavity will be bound to the Milky Way halo, having speed less than the escape speed of ≈ 560 km/s relative to the galactic center [37].
The motion of the cavity, also in orbit, has been measured to be of speed 230 ± 5 km/s [38] co-rotating with spin of the Milky Way.The motion of the majority of axions is therefore of the order ∼ 10 −3 c and the energy of each axion is well approximated by the lowest order kinematic expansion of the energy E = m a +m a ∆v 2 /2+O(∆v 4 /c 4 ) where ∆v is the relative speed of the axion to the cavity.The modulus squared of the axion field produces the power spectrum of the axion field as used in Eqn. 3 is then found to be where δν = ν − ν 0 is the frequency relative to the rest mass, and f DM (δν) is the local frequency distribution function of the axion halo.Several distribution functions were used in the search for axions during Run 1A, which are detailed in the following subsections.
A. Top-Hat Model (Live Analysis) The first and simplest distribution model is a top-hat function the width of seven MR bins (∼ 700 Hz).This model roughly reproduces the width of the local axion distribution, which at 650 MHz is expected to have an overall width in proportion to the expected virialaized speed squared w ≈ 650 MHz v 2 /2c 2 ∼ 700 Hz, and is the most robust of the models to detect power excesses at the expected width.The top-hat model is used only during the live analysis that occurs in parallel to datataking operations, informing on the overall health of data and identifying axion candidates.

B. Isotropic Isothermal Sphere
The standard halo model (SHM) distribution is the most common shape used by axion searches and direct dark-matter searches in general.The SHM is based on the assumption that the MW halo is given by a thermalized pressure-less self-gravitating sphere of particles.More specifically, we use the truncated isothermal sphere model [39], which is constructed with finite mass and has a cutoff at the halo escape speed.
In the frame of the galactic center, the velocity distribution at the solar radius takes on the near-Maxwell-Boltzmann form where σ is the dark matter velocity dispersion.The approximation of a full Maxwell-Boltzmann distribution is adequate for the MR analysis and retains an analytic form when boosted from the galactic frame by ⃗ v lab into the cavity frame where ⃗ v lab is the velocity of the lab relative to the galactic center.The motion of the lab during Run 1A was set to the orbital velocity of the sun around the galactic center.One could also incorporate the shifting motions of the Earth orbiting about the Sun and the Earth's spin, however these changes would only impact the overall width of the signal by less than 10% [34], and therefore were ignored during the analysis.

C. N-Body Model
The third line shape used in the Run 1A analysis is rooted in the highly detailed Romulus25 cosmological simulation [18,40,41], where the halos of MW-like galaxies were found to be notably denser and have narrower line shape when sampled from the reference of a Sunlike orbit.The signal shape generated takes a Maxwell-Boltzmann-like form parameterized by three constants.The N-body inspired line shape takes the form where ν 0 is the rest frame frequency, the exponent parameters are set to α = 0.36, β = 1.39, and the distribution temperature is given by T = 4.7 × 10 −7 .This shape has a width that is 1.8 times narrower than the SHM line shape.The Romulus25 analysis also revealed FIG.16.The decision tree for live candidate handling during Run 1A.The tree is entered into at the end of initial sensitivity scans of each 5-10 MHz section and traversed based on the identification and classification and candidate axion signals.Candidates are initially identified using the live analysis parameters.
a local expected DM density of ρ DM ≈ 0.6 GeV/cc, a notable increase from 0.45 GeV/cc used in previous ADMX analyses.

VII. CANDIDATE HANDLING
The processing and classification of candidate signals provides the means to claim detection of axion dark matter or to place limits on its mass and coupling.This section details the protocol used to identify, test, and re-test candidates through several filtering mechanisms unique to axion-photon conversion to robustly classify power excesses observed in the MR data.The decision protocol used for candidate handling in Run 1A is summarized in the decision tree of Fig. 16.

A. Live Analysis
A data analysis was run in parallel to data-taking operations, referred to as the live analysis, to inform the operators on the in situ sensitivity of the data and to identify candidate axion signals.
The live analysis operates on MR data under the chosen options to model the raw scan background using the six-order polynomial, and filters the prepared spectra using the top-hat signal model.These options were chosen for their low computational cost and robustness.Also, axion mass tests were made more sparsely to further reduce the computational cost, with separation of seven FIG. 17. Candidate search results (Level 1) for each segment.Candidate status was determined from the significance measure ΣA using a +3σ threshold.The number of candidates observed above the threshold can be found in Table III, and overall were seen at a significantly higher rate than expected from purely thermal fluctuations.The top-hat shape was used to determine candidate status via the live analysis run in parallel with data taking operations.The Maxwell-Boltzmann and N-body line shapes were implemented in the offline analysis and found no additional candidates.
MR bin widths to reflect the size of the top-hat signal filter.
The significance statistic of Eqn.37 was used to search for candidate axion masses.The initial set of candidates were identified with a threshold of +3σ, over which a particular mass test is considered to be a candidate.The expected number of candidates per 10 MHz segment in this mass range is ≈ 15.During the candidate handling procedure, an excess of persistent and non-persistent candidates were found in segment 3.This excess of nonpersistent candidates was found to be the result of inadequate background modeling by the six-order polynomial fit in the presence of significant background structure.Residual background structure larger than the noise floor were seen as power excesses (and deficits) by the analysis with enough regularity to skew the rate at which candidates were identified using a prescribed threshold.Also, a deficit of candidates was found in segment 1 by the live analysis.Possible systematic errors such as errant background subtraction or poorly fitted cavity Lorentzian were investigated but not found be the cause.

B. Candidate Procedure
Once candidates masses in a segment are identified, a targeted re-scan of the immediate area around each candidate frequency is conducted to the same sensitivity as the initial observations.The re-scan data was again analysed using the live analysis procedure to test for per-sistence using the same +3σ threshold.Candidates that did not persist were classified as random thermal variations.Persistent candidates were again rescanned, but to an SNR 5/3 that of the initial scans, corresponding to a slowdown in scan speed of 9/25.The re-scans are analysed for persistence using the same live analysis procedure with an updated threshold of 5σ.Candidate masses that do not pass this threshold are classified as random thermal variations and removed from further consideration.No thermal candidates are expected to pass this second rescan stage.
Remaining candidates are put through further persistence tests that also test for axion and dark-matter specific qualities.The first process is a further rescan of the candidate, again at slowed rate, and a test not only of the signal's persistence over the rescan but also of the signal's response to the Lorentzian envelope of the coupled mode T ν0 × |a| 2 , which confirms that the cavity is the origin of the power excess.The test is performed using the live analysis procedure, checking for significance that scales like 1/T ν0 .Persistent candidates are also tested at this point for a shape that matches either a Maxwellian distribution such as described in the previous section or a non-thermalized shape with dominant fine structure such as those proposed in Banik and Sikivie [42].This is performed both with the MR data as well as the full time series, mixed to resolution at the Hz level.Given the large range of proposed signals, few candidates are ruled out by this step unless they are un-physically wide.
The penultimate test in the candidate protocol is an RFI search at the remaining candidate frequencies.Using an exterior antenna and signal analyser, ambient signals are searched for about the ADMX insert and DAQ.Candidates found to have a counterpart external RF signal at the same frequency and shape, or at intermediate frequencies of the receiver, are classified as RFI.
The final test in the Run 1A candidate handling protocol is a ramping of the main magnet while taking data about a single persistent candidate.The data is then analyzed for a ∝ B 2 response in the emitted signal power particular to the inverse-Primakoff process.A signal exhibiting this response is categorized as a robust candidate for axion dark matter.Table III shows the number of candidates for each segment at each stage of the decision tree.No candidates survived past the fourth round.

C. Synthetic Axions
Run 1A saw the beginnings of a candidate injection system where artificial signals were injected into the power spectrum, using both hardware and software methods, in order to test instrument sensitivity.
Hardware axions were injected into the cavity via the weak port through with the Synthetic Axion Generator (SAG).The structure of the SAG and its placement in the RF system has been presented in Khatiwada et al. [25].
Injection frequencies were set prior to the scanning of each segment.Data from data-taking cycles with an injection were flagged and assessed independent of the main axion search.Blind injection did not occur until Run 1B [24,25].A second complete data-taking cycle at the same tuning rod configuration was taken immediately following an injection so as to not mask the data in the vicinity of the synthetic's central frequency.Only the second scan was used in the axion search analysis presented here.
Software injections of axion signals were also used to test the performance of the analysis.Simulated axion signals were imposed directly onto the raw power spectra in the forms of the Maxwellian and N-body inspired lineshapes, modulated by the fitted Lorentzian response of the cavity.We injected 25,000 software-simulated signals into the data set with couplings varying between DFSZ and 10 times KSVZ, ran through the analysis process, and evaluated the resulting candidate power to determine the systematic uncertainty associated with the background subtraction.Figure 18 shows the effect of the injected signals in both the background-subtracted spectra and the final filtered and combined spectrum.

VIII. EXCLUSION LIMITS
The previous section showed how the candidate procedure identified and classified candidates in the grand axion search analysis.Each of these greater-than-+3σ candidates were classified as non-axion in origin, therefore making the observations consistent at that level with the null-signal hypothesis under the live analysis top-hat line shape.Analyses using the Maxwellian and the Nbody axion lineshapes were performed after the run's data-taking operations.They both found no new candidates at the same threshold.
Given the null outcome, upper limits can then be formed on the axion power P a , and therefore coupling g aγγ under the priors presented in Sec.V for each axion mass tested.The one-sided 90% limits on the powerderived coupling are shown in Fig. 19 for both the SHM and N-body line shapes.

IX. SUMMARY AND CONCLUSIONS
This paper elaborates on the background and methods used by the ADMX collaboration to gather its Run 1A data and perform the axion dark matter search first reported on in Du et al. [18], which produced the first axion-photon coupling limits at DFSZ sensitivity.The ADMX apparatus was reviewed, with special concentration on the microwave cavity, receiver, and the magnet system that enabled resonant conversion of an ambient coherent axion field.Data-taking operations were explicated from their global structure down to the individual measurements of a single data-taking cycle.Interpretations of the raw data and its preparation for axion search were then presented, including the classification of the cavity and receiver transmission properties, noise temperatures, and structure of the receiver chain emitted noise power in the limit of many photon emissions per observation.The statistical basis for the axion search was then presented, both at the level of a single observation and the level of the wider Run 1A data set, and co-added into a grand spectrum spanning the observation range.The space of lineshapes for the local axion distribution and their usage were presented for both the live analysis and the off-line analyses.The axion signal candidate handling procedure was outlined, with candidates observed in the Run 1A data being catalogued either as random thermal fluctuations or interference external to the insert.No axion candidates persisted through the candidate-search process, making the data consistent with the null axion hypothesis and allowing detection limits to be formed over the observed range, with the exception of a set of RFI signals that could not be reliably removed from the data due to their strong timedependent properties.
Several updates and improvements have been made to the instrumentation and analysis that have since been used for Run 1B [23,24] and Run 1C still in progress [43].Areas highlighted in the Run 1A analysis are as follows: an improved modeling of the receiver structure during live analysis to a Padè filter from the six-order polynomial, which drastically reduces the number of candidates found in the initial search, bringing the number to near that expected from purely thermal contributions; synthetic axion injection, which was tested but not thoroughly integrated in Run 1A, has been fully incorporated into the data acquisition process using both hardware injected blinded and un-blinded candidates and software candidates that have been used to good effect to improve confidence in a search's sensitivity to the axion; and last to be mentioned here, the MSA quantum-limited amplifier of Run 1A has been replaced by JPAs (Josephson Parametric Amplifiers) in Runs 1B and 1C and have proven to be more stable and tunable.
Lastly, the possibility of fine structure existing in the local axion energy distribution presents a distinct opportunity to improve the observation's sensitivity and has motivated the collaboration to store the full time series of each observation.The resolution of this data extends to the tens of milli-hertz, the optimal level between the finest motivated fine structure features [42] and smearing due to orbital modulation from the Earth's spin and motion around the solar system center of mass.This high resolution data is significantly more sensitive to such fine structures.Therefore, as ADMX brings in more data at DFSZ sensitivity or better, searches for fine structure will play an increasingly important role in furthering the discovery potential for axion dark matter at increasingly smaller couplings or proportions of the total dark-matter density.

X. ACKNOWLEDGEMENTS
This work was supported by the U.S. Department of Energy through Grants No. DE-SC0009800, No.

FIG. 3 .
FIG. 3. Cutaway model of the ADMX main magnet's field, simulated in the absence of the rest of the apparatus using the COMSOL multi-physics simulation software.The cutaway model of the cavity is included for reference.The size of the arrow indicates the strength of the field at that point in addition to the color scale.

FIG. 4 .
FIG.4.ADMX main cavity and TM010-like simulated mode properties: (a) view of the interior of the main cavity; (b) TM010-like mode electric field in side view; (c) simulated mode form factor (blue dashed) and measured cyogenic loaded quality factor (red solid) over the Run 1A tuning path.Errors (shaded regions) in the form factor and quality factor were computed to be 6% and 2.2% respectively, as discussed in Section IV E and presented in TableI.

FIG. 5 .
FIG. 5. Sample transfer functions of the resonant cavity around well-coupled mode (TM010-like).The mode is centered at 660.05 MHz.(Left) (top) Illustration of the mode transmission function S12 from the weak port through the main port with a Lorentzian envelope about the mode center (Eqn.4).(bottom) Example S12 scan with fit Lortentzian with the cavity at sub-Kelvin temperature.(Right) (top) Illustration of the mode reflection function S11 in/out of the main port with envelope of the form in Eqn. 8. (bottom) Example S11 scan with the best fit envelope with cavity at sub-Kelvin temperature.

FIG. 8 .
FIG. 8. Tuning configurations of Run 1A: (a) The electric field strength component ⃗ E • ẑ for an example rod configuration for Run 1A included one rod being placed at the wall while the other turned from 107 − 141 degrees; (b) measured cyrogenic QL of the TM010-like over the rod tuning range; (c) simulated transmission (S12) spectra of the cavity, highlighting modes between 500 − 1000 MHz over the rod tuning range.The target TM010-like mode over the Run 1A scanning range is bracketed in black.

5 FIG. 11 .
FIG. 11.Scaling of fractional fluctuation distribution width with number of co-added signals for randomly generated simulated receiver spectra with typical shape (left panel) and anomalous receiver shapes (right panels).Solid lines indicate the average width, while the shaded band gives the one sigma uncertainty.The ideal variance scaling of co-added uncorrelated noise is σ 2 /µ 2 ∝ 1/Nscans (top panels dashed line).Fractional differences (lower panels) are relative to the ideal scaling.

FIG. 14 .
FIG.14.Fit significance ΣA using the sample scan from Fig.13for each line shape model used in the live and offline analyses.Different line shape models are discussed in Section VI.

FIG. 18 .
FIG. 18. (Upper) Series of background-subtracted single scans with synthetic axion signals with the N-body inspired signal shape [40], one at KSVZ coupling and one at DFSZ coupling.The KSVZ signal is easily visible in these individual spectra; the DFSZ signal, being a factor of 7 smaller, is not.(Lower) Same data after the individual scans have been optimally filtered and combined.Both KSVZ and DFSZ signals are visible with high SNR. Figure first appeared in Du et al. [18].
. The 35 MHz range is divided into four segments, observed chronologically as: 650−660 MHz; 660−670 MHz; 670 − 680 MHz; and 645 − 650 MHz.Each segment is in- FIG. 7. Schematic of the Run 1A digitizer chain.Note that Output (1) of the cold RF receiver shown in Fig.
receiver studies,

TABLE II .
Conditions and outcomes for data cuts in the offline analysis.Axion search data collected in these conditions were considered irrelevant or of poor quality.Conditions were applied in top-to-bottom order, making a condition's number of scans removed dependant on the conditions listed above it.
Each line shape used in Run 1A is represented both in its smooth form (dashed) and as a binned sample (solid) where a bin edge is aligned with the axion rest mass.The width of the top-hat line shape was equivalent in width to seven MR bins and was not varied throughout the scan range.The widths of the Maxwell-Boltzmann and Nbody inspired line shapes were varied on a scan-by-scan basis as described in Section V.

s) through region to DFSZ sensitivity Targeted rescan candidates at original SNR Rescan persistent candidates with 5/3- enhanced SNR Check signal shape and response to mode's Lorentzian Search for external/internal RFI Ramp main magnet to check for B-squared response
pe and shape check First sweep(

TABLE III .
Number of candidates present at each stage of the candidate handling decision tree.The identified RFI in the region 660.16 − 660.27MHz was unfortunately found to change its intensity and shape with time and was unable to be subtracted from the power spectrum.Axion tests overlapping with these signals were omitted when calculating the final exclusion limits.