Observation of the Crab Nebula with the HAWC Gamma-Ray Observatory

TL;DR

This paper reports the detection and spectral analysis of the Crab Nebula using the HAWC gamma-ray observatory, demonstrating HAWC's high sensitivity and wide field-of-view in the TeV energy range, and establishing it as a new tool for sky surveys.

Contribution

First detailed spectral measurement of the Crab Nebula with HAWC, showcasing its capabilities and establishing a baseline for future TeV sky surveys.

Findings

Crab spectrum well-fit by a curved power law between 1 and 37 TeV.

HAWC's sensitivity exceeds previous observatories at energies above 10 TeV.

Systematic flux uncertainty estimated at ±50% in the 1-37 TeV range.

Abstract

The Crab Nebula is the brightest TeV gamma-ray source in the sky and has been used for the past 25 years as a reference source in TeV astronomy, for calibration and verification of new TeV instruments. The High Altitude Water Cherenkov Observatory (HAWC), completed in early 2015, has been used to observe the Crab Nebula at high significance across nearly the full spectrum of energies to which HAWC is sensitive. HAWC is unique for its wide field-of-view, nearly 2 sr at any instant, and its high-energy reach, up to 100 TeV. HAWC's sensitivity improves with the gamma-ray energy. Above $\sim$1 TeV the sensitivity is driven by the best background rejection and angular resolution ever achieved for a wide-field ground array. We present a time-integrated analysis of the Crab using 507 live days of HAWC data from 2014 November to 2016 June. The spectrum of the Crab is fit to a function of the form $\phi(E) = \phi_0 (E/E_{0})^{-\alpha -\beta\cdot{\rm{ln}}(E/E_{0})}$. The data is well-fit with values of $\alpha=2.63\pm0.03$, $\beta=0.15\pm0.03$, and log$_{10}(\phi_0~{\rm{cm}^2}~{\rm{s}}~{\rm{TeV}})=-12.60\pm0.02$ when $E_{0}$ is fixed at 7 TeV and the fit applies between 1 and 37 TeV. Study of the systematic errors in this HAWC measurement is discussed and estimated to be $\pm$50\% in the photon flux between 1 and 37 TeV. Confirmation of the Crab flux serves to establish the HAWC instrument's sensitivity for surveys of the sky. The HAWC survey will exceed sensitivity of current-generation observatories and open a new view of 2/3 of the sky above 10 TeV.

Tables

Table 1. Table 1: Steps in the HAWC Event Reconstruction. The best core and direction reconstructions are applied with gradual narrowing of the hits used.

Step	Description	Hit Selection
1	Calibration
2	Hit Selection
3	Center-of-Mass Core Reconstruction	Selected Hits
4	SFCF Core - First Pass	Selected Hits
5	Direction - First Pass	Selected Hits
6	SFCF Core - Second Pass	Selected Hits within 50 ns of First-pass Plane
7	Direction - Second Pass	Selected Hits within 50 ns of First-pass Plane
8	Compactness	Selected Hits within 20 ns of Second-pass Plane
9	PINCness	Selected Hits within 20 ns of Second-pass Plane

Table 2. Table 2: Cuts used for the analysis. The definition of the size bin ℬ ℬ \mathcal{B} is given by the fraction of available PMTs, f hit subscript 𝑓 hit f_{\rm{hit}} , that record light during the event. Larger events are reconstructed better and ψ 68 subscript 𝜓 68 \psi_{68} , the angular bin that contains 68% of the events, reduces dramatically for larger events. The parameters 𝒫 𝒫 \mathcal{P} and 𝒞 𝒞 \mathcal{C} (Section 2.6 ) characterize the charge topology and are used to remove hadronic air shower events. Events with a 𝒫 𝒫 \mathcal{P} less than indicated and a 𝒞 𝒞 \mathcal{C} greater than indicated are considered photon candidates. The cuts are established by optimizing the statistical significance of the Crab and trend toward harder cuts at larger size events. The number of excess events from the Crab in each ℬ ℬ \mathcal{B} bin per transit is shown as well.

$ℬ$	$f_{\mathrm{hit}}$	${\psi }_{68}$	$𝒫$ Maximum	$𝒞$ Minimum	Crab Excess Per Transit
1	6.7 - 10.5%	1.03	2.2	$>$7.0	68.4 $\pm$ 5.0
2	10.5 - 16.2%	0.69	3.0	9.0	51.7 $\pm$ 1.9
3	16.2 - 24.7%	0.50	2.3	11.0	27.9 $\pm$ 0.8
4	24.7 - 35.6%	0.39	1.9	15.0	10.58 $\pm$ 0.26
5	35.6 - 48.5%	0.30	1.9	18.0	4.62 $\pm$ 0.13
6	48.5 - 61.8%	0.28	1.7	17.0	1.783 $\pm$ 0.072
7	61.8 - 74.0%	0.22	1.8	15.0	1.024 $\pm$ 0.053
8	74.0 - 84.0%	0.20	1.8	15.0	0.433 $\pm$ 0.033
9	84.0 - 100.0%	0.17	1.6	3.0	0.407 $\pm$ 0.032

Table 3. Table 3: Summary of primary contributions to HAWC systematic error in measuring photon fluxes. The different effects are described in the text. Systematics in the overall flux, the spectral index of sources, and the energy scale are shown. The systematics claims are conservative and are likely to improve with more understanding and better modeling.

Systematic	Overall Flux	Spectral Index	log10(E)
Charge Resolution/ Relative Quantum Efficiency	$\pm$ 20%	$\pm$ 0.05	0.1
PMT Absolute Quantum Efficiency	$\pm$ 15%	$\pm$ 0.05	0.1
Time Dependence, PMT Layout and Crab Optimization	$\pm$ 10%	$\pm$ 0.1	0.1
Angular Resolution	$\pm$ 20%	$\pm$ 0.1
Late Light Simulation	$\pm$ 40%	$\pm$ 0.15	0.15
Total Flux	$\pm$ 50%	$\pm$ 0.2	0.2

Equations

$$S_{i}=S(A,\vec{x},\vec{x}_{i})=A\Big{(}\frac{1}{2\pi\sigma^{2}}e^{-{|\vec{x}_{i}-\vec{x}|^{2}}/{2\sigma^{2}}}+\frac{N}{(0.5+{|\vec{x}_{i}-\vec{x}|}/{R_{m}})^{3}}\Big{)}$$

$$\mathcal{C}=\frac{{\mathrm{N}}_{\mathrm{hit}}}{\mathrm{CxP}{\mathrm{E}}_{40}}$$

$$\mathcal{P}=\frac{1}{N}\sum _{i=0}^N\frac{({\zeta }_i-⟨{\zeta }_i⟩{)}^2}{{{\sigma }_{{\zeta }_i}}^2}$$

$$\rho (\psi )=\alpha G_1(\psi )+(1-\alpha )G_2(\psi )$$

$$G_i(\psi )=\frac{1}{2\pi {\sigma }_i^2}{e}^{-\frac{{\psi }^2}{{2{\sigma }_i}^2}}$$

$$\ln (\mathcal{L}(\vec{a}))=\sum _{\mathcal{B}=1}^9\sum _{p=1}^N\ln \left(\frac{{(B_p+S_p(\vec{a}))}^{M_p}{e}^{B_p+S_p(\vec{a})}}{M_p!}\right)$$

$$\varphi (E)={\varphi }_0(E/E_0{)}^{-\alpha -\beta \cdot \ln (E/E_0)}$$

$$\varphi (x)=\{\begin{array}{ll}0& \text{if }E\ge E_{\mathrm{high}}\\ {\varphi }_0{E}^{-\alpha },& \text{otherwise}\end{array}$$

$$\varphi (x)=\{\begin{array}{ll}0& \text{if }E\le E_{\mathrm{low}}\\ {\varphi }_0{E}^{-\alpha },& \text{otherwise}\end{array}$$

Full text

We present a time-integrated analysis of the Crab using 507 live days of HAWC data from 2014 November to 2016 June. The spectrum of the Crab is fit to a function of the form $\phi(E)=\phi_{0}(E/E_{0})^{-\alpha-\beta\cdot{\rm{ln}}(E/E_{0})}$. The data is well-fit with values of $\alpha=2.63\pm 0.03$, $\beta=0.15\pm 0.03$, and log${}_{10}(\phi_{0}~{}{\rm{cm}^{2}}~{}{\rm{s}}~{}{\rm{TeV}})=-12.60\pm 0.02$ when $E_{0}$ is fixed at 7 TeV and the fit applies between 1 and 37 TeV. Study of the systematic errors in this HAWC measurement is discussed and estimated to be $\pm$50% in the photon flux between 1 and 37 TeV.

Confirmation of the Crab flux serves to establish the HAWC instrument’s sensitivity for surveys of the sky. The HAWC survey will exceed sensitivity of current-generation observatories and open a new view of 2/3 of the sky above 10 TeV.

gamma rays: observations — surveys — acceleration of particles — pulsars: individual (Sn-1054) — ISM: individual (Crab Nebula)

††journal: ApJ

\AuthorCollaborationLimit

=500

1 Introduction

The Crab Pulsar Wind Nebula (the Crab Nebula or the Crab) occupies a place of special distinction in the history of high-energy astrophysics. It was the first high-confidence TeV detection in 1989 using the Whipple telescope (Weekes et al., 1989) and is the brightest steady source in the Northern sky above 1 TeV. It has been observed with imaging atmospheric Cherenkov telescopes (IACTs) since (Tanimori et al., 1998; Aharonian et al., 2004, 2006; Celik, 2008; Aleksić et al., 2015). The first observation using a ground array was the 2003 Milagro detection (Atkins et al., 2003), and the signal was subsequently seen in other ground arrays (Amenomori et al., 2009; Bartoli et al., 2015).

The TeV emission arises from inverse-Compton (IC) up-scattering of low-energy ambient photons by energetic electrons accelerated in shocks surrounding the central pulsar (Atoyan & Aharonian, 1996). Photons from synchrotron emission of the electrons themselves are likely the dominant IC target with sub-dominant contributions from the cosmic microwave background and the extragalactic background light (Martín et al., 2012). Despite rare flaring emission below 1 TeV (Tavani et al., 2011; Abdo et al., 2011), and a potential TeV flare (Aielli et al., 2010), the Crab is generally believed to be steady at higher energies (Abramowski et al., 2014; Aliu et al., 2014; Bartoli et al., 2015). Consequently, the Crab Nebula has been adopted as the reference source in TeV astronomy and is a reliable beam of high-energy photons to use for calibrating and understanding new TeV gamma-ray instruments.

The High Altitude Water Cherenkov (HAWC) observatory is a new instrument sensitive to multi-TeV hadron and gamma-ray air showers, operating at latitude of +19∘N at an altitude of 4,100 meters in the Sierra Negra, Mexico. HAWC consists of a large 22,000 $\rm{m}^{2}$ area densely covered with 300 Water Cherenkov Detectors (WCDs), of which 294 have been instrumented. Each WCD consists of a 7.3-meter diameter, 5-meter tall steel tank lined with a plastic bladder and filled with purified water. Figure 1 shows a schematic of the WCD and an overhead view of the full instrument. At the bottom of each WCD, three 8-inch Hamamatsu R5912 photomultiplier tubes (PMTs) are anchored in an equilateral triangle of side length 3.2 meters, with one 10-inch high-quantum efficiency Hamamatsu R7081 PMT anchored at the center.

A high-energy photon impinging on the atmosphere above HAWC initiates an extensive electromagnetic air shower. The resulting mix of relativistic electrons, positrons and gamma rays propagates to the ground in a thin tortilla of particles at nearly the speed of light. Energetic particles that reach the instrument can interact in the water and produce optical light via Cherenkov radiation. The high altitude of HAWC sets the scale for the photon energy that can be detected. At HAWC’s altitude, the shower from a 1 TeV photon from directly overhead will have about 7% of the original photon energy left when the shower reaches the ground. The fraction of energy reaching the ground rises to $\sim$28% at 100 TeV. The detector is fully efficient to gamma rays with a primary energy above $\sim$1 TeV. Lower-energy photons can be detected when they fluctuate to interact deeper in the atmosphere than typical.

The voltages on the HAWC PMTs are chosen to match the PMT gains across the array. PMT pulses are amplified, shaped, and passed through two discriminators at approximately 1/4 and 4 PEs (Abeysekara et al., 2016) and digitized. The length of time that PMT pulses spend above these thresholds (time-over-threshold or ToT) is used to estimate the total amount of charge collected in the PMT. Noise arises from a number of sources including PMT afterpulsing, fragments of sub-threshold air showers, PMT dark noise, and other sources. The 8-inch PMTs have a hit rate (a hit being each time the PMT signal crosses the 1/4 PE threshold) due to the combined effect of these sources of 20–30 kHz and the 10-inch PMTs have a hit rate of 40–50 kHz.

The data from the front-end electronics is digitized with commercial time-to-digital converters (TDCs) and passed to a farm of computers for real-time triggering and processing. Events are preserved by the computer farm if they pass the trigger condition: a simple multiplicity trigger, requiring some number, $N_{\rm{thresh}}$, PMTs hit within 150 ns. Hits 500 ns prior to a trigger and up to 1000 ns after a trigger are also saved for reconstruction. During the operation of HAWC, $N_{\rm{thresh}}$ has varied between 20–50. The trigger rate at the time of writing, due primarily to hadronic cosmic-ray air showers, is $\sim$24 kHz with $N_{\rm{thresh}}=28$.

The reconstruction process involves determining the direction, the likelihood for the event to be a photon, and the event’s size. A first-look reconstruction is applied at the HAWC site. In this analysis, all the data has been reconstructed again (the fourth revision, or Pass 4, of the reconstruction process) off-site in order to have a uniform dataset and the best calibrations available. The chief background to gamma-ray observation is the abundant hadronic cosmic-ray population. Individual gamma-ray-induced air showers can be distinguished from cosmic-ray showers by their topology and the presence of deeply penetrating particles at the ground.

The strength of HAWC over the IACT technique is that photon showers may be detected across the entire $\sim$ 2 sr field-of-view of the instrument, day or night, regardless of weather conditions. As such, HAWC is uniquely suited to study the long-duration light curve of objects and to search for flaring sources in real time. Additionally, since sources are observed on every transit, HAWC obtains thousands of hours of exposure on each source, greatly improving the sensitivity to the highest-energy photons.

Section 2 outlines the algorithms by which the direction, size, and type (photon or hadron) of each shower is determined. Section 3 describes the identification of the gamma-ray signal from the Crab Nebula. The fit to the Crab energy spectrum, including a treatment of systematic errors, is described in Section 4. Finally, a discussion of the result is presented in Section 5, including a comparison to prior spectra measured by peer experiments and a computation of the sensitivity of the HAWC instrument, anchored in the agreement of the HAWC measurement to other experiments.

2 Air Shower Reconstruction

Events from the detector are reconstructed to determine the arrival direction of the primary particle and the size of the resulting air shower on the ground, a proxy for the primary particle’s energy. Table 1 summarizes the steps in reconstruction of HAWC events as explained below.

To illustrate the event reconstruction, Figure 2 shows a strong gamma-ray candidate from the Crab Nebula. In Section 2.1, the simulation is briefly described as it is key to evaluating the reconstruction process. Section 2.2 describes the calibration, by which the time and light level in individual PMTs are determined. Section 2.3 discusses the selection of PMT signals for reconstruction and the event size measurement. The direction reconstruction occurs in two steps, first the core reconstruction, described in Section 2.4, and then the direction determination, described in Section 2.5. The air shower core, the dense concentration of particles along the direction of the original primary, is needed to make the best reconstruction of the air shower’s direction since the air shower arrival front is delayed from a pure plane, depending on the distance from the core. The identification of photon candidates is presented in Section 2.6. The directional fit is iterated to suppress noise and this iteration is explained in Section 2.7.

2.1 Simulation

The HAWC instrument is modeled using a combination of community-standard simulation packages and custom software. The CORSIKA package (v7.4000) is used for simulation of air showers, propagating the primary particles through the atmosphere to the ground (Heck et al., 1998). At ground level, a GEANT 4 simulation (v4.10.00) of the shower particles is used to propagate the ground-level particles through the HAWC tanks and to track the Cherenkov photons to the faces of the PMTs (Agostinelli et al., 2003).

The response of the PMTs and the calibration are approximated with a custom simulation that assumes that recorded light is faithfully detected with some efficiency and an uncertainty in the logarithm of the total charge recorded. Decorrelated single PE noise is added. The absolute PMT efficiency for detecting Cherenkov photons is established by scaling the simulated PMT response to vertical muons to match data. Most muons passing through HAWC are minimum ionizing with nearly constant energy loss. Vertical muons, therefore, are a nearly constant light source and convenient for establishing to total PMT efficiency. Simulated events are subsequently reconstructed by the same procedure as experimental data to study the performance of the algorithms.

2.2 Calibration

The first step in the reconstruction process is calibration, the processes by which true time and light level in each PMT are estimated from the TDC-measured threshold-crossing times of each PMT (Lauer, 2013; Solares et al., 2016).

The calibration associates the measured ToT in each PMT with the true number of PEs. To give a sense of scale, the ToT for a single PE crossing the low-threshold discriminator (about 1/4 PE) is $\sim$100 ns. Above a few PEs, the higher-threshold (about 4 PEs) ToT is used for charge assignment and a high-threshold ToT of 400 ns roughly corresponds to a charge of $10^{4}$ PEs. The time scale for these ToTs is determined by the shaping of the front-end electronics and is chosen to be longer than the characteristic arrival time distribution of PEs during an air shower so as to integrate the whole air shower arrival into one PMT hit.

In addition to the PE measurement, the calibration system accounts for electronic slewing of the PMT waveform: Lower-PE waveforms cross threshold later than contemporaneous high-PE waveforms.

Subsequent reconstruction algorithms treat all PMTs, the 8-inch and 10-inch, as identical, despite the larger size and greater efficiency of the 10-inch PMTs. To accommodate this, an effective charge $Q_{\rm{eff}}$ is defined. For $Q_{\rm{eff}}$, PE values from the central 10-inch PMTs are scaled by a factor of 0.46 to place them on par with the 8-inch PMTs.

Finally, each PMT has a single calibrated timing offset that accounts for the different cable lengths and any other timing delays that may differ from PMT to PMT. These delays are measured to within a few ns by the calibration system and are refined to sub-ns precision by iterated fits to hadronic air showers. Since the hadronic background is isotropic, the point of maximum cosmic-ray density is overhead and the PMT timing pedestals are chosen to ensure this is true. A final small ($\sim$0.2∘) rotation of all events is performed to ensure that the Crab Nebula appears in its known location.

Throughout the analysis, the Crab is assumed to be at a location of 83.63∘ right ascension and 22.01∘ declination, in the J2000.0 epoch, taken from Aharonian et al. (2004). While the pulsar position is known more precisely (e.g. Comella et al. (1969)), this precision is sufficient for use in HAWC.

2.3 Hit Selection and Event Size Bins

As described in Section 1, the HAWC DAQ records 1.5 $\mu$s of data from all PMTs that have a hit during an air shower event. A subset of these hits are selected for the air shower fit. To be used for the air shower fit, hits must be found between -150 and +400 ns around the trigger time. Hits are removed if they occur shortly after a high-charge hit under the assumption that these hits are likely contaminated with afterpulses. Additionally, hits are removed if they have a pattern of TDC crossings that is not characteristic of real light; they cannot be calibrated accurately. Finally, each channel has an individual maximum calibrated charge, typically a few thousand PEs, but no more than 104 PEs, above which the PMTs are not used. Above $\sim$104 PEs, corresponding to a ToT of $\sim$400 ns, prompt afterpulsing in the PMTs can artificially lengthen the ToT measurement giving a false measurement. Channels are considered available for reconstruction if they have a live PMT taking data which has not been removed by one of these cuts.

The angular error and the ability to distinguish photon events from hadron events is strongly dependent on the energy and size of events on the ground. We adopt analysis cuts and an angular resolution description that depends on this measured size. The data is divided into 9 size bins, $\mathcal{B}$, as outlined in Table 2. The size of the event is defined as the ratio of the number of PMT hits used by the event reconstruction to the total number of PMTs available for reconstruction, $f_{\rm{hit}}$. This definition allows for relative stability of the binning when PMTs are occasionally taken out of service.

For this analysis, events are only used if they have more than 6.7% of the available PMTs seeing light. Since typically 1000 PMTs are available, typically a minimum of 70 PMTs is needed for an event. This is substantially higher than the trigger threshold. The data between the trigger threshold and the threshold for $\mathcal{B}=1$ in this analysis consists of real air showers, and techniques to recover these events and lower the energy threshold, beyond what is presented here, are under study.

Figure 3 shows the distribution of true energies as a function of the $\mathcal{B}$ of the events. The distribution of energies naturally depends heavily on the source itself, both its spectrum and the angle at which it culminates during its transit. A pure power-law spectrum with a shape of $E^{-2.63}$ and a declination of 20∘ was assumed for this figure. As $\mathcal{B}$ is a simple variable — containing no correction for zenith angle, impact position, or light level in the event — the energy distribution of $\mathcal{B}$ bins is wide. Section 5.3 discusses planned improvements to this event parameter that will measure the energy of astrophysical gamma rays better.

Bin $\mathcal{B}=9$ bears particular attention. It is an “overflow” bin containing events which have between 84% and 100% of the PMTs in the detector seeing light. Typically, a 10 TeV photon will hit nearly every sensor and the $\mathcal{B}$ variable has no dynamic range above this energy. This limit is not intrinsic to HAWC and variables that utilize the light level seen in PMTs on the ground, similar to what was used in the original sensitivity study (Abeysekara et al., 2013), have dynamic range above 100 TeV. These variables, not used in this analysis, will improve the identification of high-energy events. This is discussed farther in Section 5.3.

2.4 Core Reconstruction

In an air shower, the concentration of secondary particles is highest along the trajectory of the original primary particle, termed the air shower core. Determining the position of the core on the ground is key to reconstructing the direction of the primary particle. In the sample event, Figure 2, the air shower core is evident in Figure 2a. The image is an overhead view of the HAWC detector with circles indicating the WCD location and the PMTs within the WCDs. The colors indicate the amount of light (measured in units of PEs) seen in each PMT. The air shower core is evident as the point of maximum PE density.

The PE distribution on the ground is fit with a function that decreases monotonically with the distance from the shower core. The signal in the $i$th PMT, $S_{i}$, is presumed to be

[TABLE]

where $\vec{x}$ is the core location, $\vec{x}_{i}$ is the location of the measurement, $R_{m}$ is the Molière radius of the atmosphere, approximately 120 m at HAWC altitude, $\sigma$ is the width of the Gaussian, and $N$ is the normalization of the tail. Fixed values of $\sigma=10$ m and $N=5\cdot 10^{-5}$ are used. This leaves three free parameters, the core location and overall amplitude, $A$.

The functional form used in this algorithm, termed the Super Fast Core Fit (SFCF), is a simplification of a modified Nishimura-Kamata-Greisen (NKG) function (Greisen, 1960) and is chosen for rapid fitting of air shower cores. The NKG function has an additional free parameter, the shower age, and involves computationally intensive power law and gamma function evaluation. The SFCF hypothesis in Equation 1 is similar but numerical minimization can converge faster because: the function is simpler, the derivatives are computed analytically, and the lack of a pole at the core location.

Figure 2b shows the recorded charge in each PMT as a function of the PMT’s distance along the ground to the reconstructed shower core. The fit for this event is shown along with the PINCness moving average from Section 2.6. While the full NKG function would describe the lateral distribution better, the SFCF form allows rapid identification the center of showers and this is sufficient for the present analysis. Cores can be localized to a median error of $\sim$2 meters for large events ($\mathcal{B}=8$) and $\sim$4 meters for small events ($\mathcal{B}=3$) for gamma-ray events with a core that lands on the HAWC detector. The error in reconstructing the shower core increases as the core moves further from the array. For example, a shower with a core that is 50 meters from the edge of the array will have an error in the location of the core of $\sim$35 meters.

2.5 Direction Reconstruction

To first order, the air shower particles arrive on a plane defined by the speed of light and direction of the primary particle. In fact, the shower front has a slight conical shape centered at the air shower core. Several effects lead to this shape. First, particles far from the core arrive due to multiple scattering and longer travel distances. Second, the multiple scattering of particles at the edges of the shower leads to a broader arrival time distribution than at the core. Since the number of particles decreases with increasing distance from the shower core there are fewer opportunities to sample from the particle arrival time distribution. This decrease in sampling leads to a delay in the measured arrival time of the shower.

The conical shape of the air shower front can be easily seen with the sample event in Figure 2. Figure 2c shows, for each PMT in the sample event, the calibrated time each PMT saw light. The color trend is due to the inclined direction of the air shower. Taking out this inclination, we can see the curved shower front. The event in Figure 2 was chosen for display because — taken from a high signal-to-background sample of data — it is very likely a photon from the Crab. If we assert the origin of this particle is the Crab Nebula, we know the air shower plane precisely and can make Figure 2d: We adjust the times of each PMT hit assuming a pure planar air shower originating at the Crab Nebula. Figure 2d shows the plane-corrected time as a function of the PMT distance from the core of the air shower along the ground. A pure planar shower would be a horizontal line on this figure. The delay of particles far from the core, is evident.

In the present analysis, we use the reconstructed core location to correct for this effect. A combined curvature/sampling correction — a function of the distance of hits from the shower core and the total charge recorded in the PMT — is utilized for this correction. The curvature/sampling correction is based on a combination of simulation and Crab observations. The rough functional form is tabulated using gamma-ray simulation. The simulation-optimized curvature/sampling correction yields a measured angular resolution approximately a factor of 2 worse than predicted from simulation. The origin of this discrepancy is likely due to some oversimplification of the electronics simulation. Repeated fits to the Crab have yielded a modification to the curvature/sampling correction that is a simple quadratic function of the distance between a hit and the shower core. While the origin of the discrepancy is under investigation, it amounts to a relatively small correction, approximately 2 ns/100 meters. Nonetheless, the improvement in the angular resolution is nearly a factor of two for all $\mathcal{B}$. Remaining disagreement between the simulated and measured angular resolution is adopted as a systematic error.

After correcting for the sampling and curvature, the angular fit is a simple $\chi^{2}$ planar fit and has been described before (Atkins et al., 2003).

2.6 Photon/Hadron Separation

Hadronic cosmic rays are the most abundant particles producing air showers in HAWC and constitute the chief background to high-energy photon observation. The air showers produced by high-energy cosmic rays and gamma rays differ: gamma-ray showers are pure electromagnetic showers with few muons or pions. Conversely, hadronic cosmic rays produce hadronic showers rich with pions, muons, other hadronic secondaries, and structure due to the showering of daughter particles created with high transverse momentum. In HAWC, these two types of showers appear quite different, particularly for showers above several TeV.

Figure 4 shows the lateral distributions for two showers, an obvious cosmic ray (left) and a strong photon candidate (right) from the Crab Nebula. The effective light level $Q_{\rm{eff}}$ falls off for hits further from the shower core in both showers, but in the hadronic shower there are sporadic high-charge hits far from the air shower’s center. This clumpiness is characteristic of hadronic showers and arises from a combination of penetrating particles (primarily muons) and hadronic sub-showers which are largely absent in photon-induced showers.

Two parameters are used to identify cosmic-ray events. The first parameter, compactness, was used in the sensitivity study (Abeysekara et al., 2013). The variable ${\rm{CxPE}}_{40}$ is the effective charge measured in the PMT with the largest effective charge outside a radius of 40 meters from the shower core. We then define the compactness, $\mathcal{C}$, as

[TABLE]

where ${\rm{N_{hit}}}$ is the number of hit PMTs during the air shower. ${\rm{CxPE}}_{40}$ is typically large for a hadronic event, so $\mathcal{C}$ is small.

In addition to the largest hit outside the core, the “clumpiness” of the air shower is quantified with a parameter $\mathcal{P}$, termed the PINCness of an event (short for Parameter for Identifying Nuclear Cosmic-rays). $\mathcal{P}$ is defined using the lateral distribution function of the air shower, seen in Figure 4. Each of the PMT hits, $i$, has a measured effective charge $Q_{\rm{eff},i}$. $\mathcal{P}$ is computed using the logarithm of this charge $\zeta_{i}={\rm{log}}_{10}(Q_{\rm{eff},i})$. For each hit, an expectation is assigned $\langle\zeta_{i}\rangle$ by averaging the $\zeta_{i}$ in all PMTs contained in an annulus containing the hit, with a width of 5 meters, centered at the core of the air shower.

$\mathcal{P}$ is then calculated using the $\chi^{2}$ formula:

[TABLE]

The errors $\sigma_{\zeta_{i}}$ are assigned from a study of a sample strong gamma-ray candidates in the vicinity of the Crab.

The $\mathcal{P}$ variable essentially requires axial smoothness. Figure 4 shows the moving average $<\zeta_{i}>$ for two sample events. The hadronic event in Figure 4 is “clumpy” and has several hits that differ sharply from the moving average yielding a large $\mathcal{P}$.

The $\mathcal{C}$ and $\mathcal{P}$ variables are well-modeled in simulation. Figure 5 shows the measured distribution of $\mathcal{C}$ in the vicinity of the Crab Nebula (the Crab region) and in an annular reference region around the Crab (the background region). The background region is scaled to have the same solid angle as the Crab region. The distributions in the vicinity of the Crab are made of a combination of hadronic cosmic-rays and true photons from the Crab. Figures 5c and 5d show these distributions with the background distribution subtracted. The subtraction yields the data-measured distribution of $\mathcal{C}$ for gamma rays from the Crab. Figure 6 is a comparable figure for $\mathcal{P}$. Figures 5 and 6 are compared to a simulation prediction from the final fitted flux from Section 4; the simulation agreement is evident.

2.7 Noise and Fit Refinement

HAWC’s outer 8-inch PMTs individually trigger at some 20–30 kHz and the 10-inch central PMTs at 40-50 kHz. Of this random noise, roughly half (with large uncertainties) are believed to be due to real shower fragments and roughly half due to non-shower sources like radioactive contaminants in the PMT glass or PMT afterpulsing. Approximately 1 kHz of the noise is “high PE” noise from individual accidental air shower muons (10–200 PEs) and can be correlated between the PMTs within a WCD that the muon hits.

This noise can bias the air shower reconstruction and the muon noise has the potential, if not removed, to confuse the identification of gamma-ray showers because a single muon is enough to indicate that a shower is of hadronic origin.

In order to achieve the best direction reconstruction and to avoid falsely rejecting true photons, the SFCF core reconstruction and the plane fit are each performed twice as outlined in Table 1. During the first pass, all selected hits are used to locate the core and initial direction. After this first “rough” fit, hits that are more than $\pm$50 ns from the curvature/sampling corrected air shower plane are removed and the shower is fit a second time.

The computation of photon/hadron separation variables is done with hits that are within $\pm$ 20ns from the curvature/sampling-corrected air shower plane. With these cuts, only some $\sim$4% of gamma ray events will have an accidental muon contributing to the photon/hadron variable computation and risk being falsely rejected.

3 Crab Nebula Signal

Once individual events are reconstructed, the identification and characterization of sources proceeds. Section 3.1 discusses the first 553 days of data-taking. Section 3.2 describes how the residual cosmic-ray background, after gamma/hadron separation cuts, is estimated in the vicinity of the Crab. Section 3.3 describes the validation of the angular resolution. Section 3.4 describes the optimization of photon/hadron discrimination cuts.

3.1 Dataset

We consider here data taken by HAWC between 2014 November 26 and 2016 June 2, a total elapsed time of 553 days. The detector was not taking data for a cumulative time of 40 days for various operational reasons. Additionally, a further 7 days of data was rejected due to trigger rate instability. This yields a total livetime of 506.7 days, an average fractional livetime of 92%. Figure 7 shows the fractional livetime achieved in blocks of 10 days. With the exception of one period of extended downtime (due to a failure of the power transformer at the site), during no single 10-day block was the detector live less than 75% of the time.

The occasional downtime does not heavily bias the exposure. Figure 7 shows the exposure of the instrument, measured as the fractional deviation of the number of shower events observed as a function of reconstructed right ascension. The anisotropy in cosmic-ray arrival direction is subdominant at $\sim$10*-3* (Abdo et al., 2009) and we can conclude that the exposure is flat to within $\pm$2%. During the 507-day livetime, the Crab is visible above a zenith angle of 20∘ for $\sim$1400 hours and above a zenith angle of 45∘ for $\sim$3200 hours.

3.2 Background Estimation

Even with strict photon/hadron discrimination cuts, the data is still dominated by hadronic cosmic-ray events. Fortunately, the directions of hadronic cosmic rays are randomized by magnetic deflection in transit from their sources and the population of cosmic rays is isotropic to a few parts in 103 (Abdo et al., 2009). Gamma-ray sources, by comparison, appear as localized “bumps” on this smooth background. In order to identify gamma-ray sources, this background contamination must be estimated.

We utilize an algorithm developed for analysis of Milagro data (Abdo et al., 2012) called direct integration. Due to the strong cosmic-ray rejection capability of HAWC, the number of background events observed in the highest-energy $\mathcal{B}$ bins is very sparse. To compensate for the sparseness of the background, the direct integration background computation has been averaged by $0.5^{\circ}$. This averaging has been studied using simulated data sets and does not adversely affect the background estimate.

Knowing the background at each pixel in the sky, we can evaluate the likelihood that there is a gamma-ray source at a specific location and the photon flux from each source.

3.3 Point Spread Function

The point spread function, $\rho(\psi)$, describes how accurately the directions of gamma-ray events are reconstructed. Here, $\psi$ is the space-angle difference between the true photon arrival direction and the reconstructed direction. To a good approximation, the point spread function of HAWC is the sum of two 2-dimensional Gaussians with different widths.

[TABLE]

where $G_{i}$ is a Gaussian distribution with width $\sigma_{i}$

[TABLE]

which is normalized to unity across the unit sphere.

Figure 8 exhibits the measured angular resolution in HAWC data in two size bins $\mathcal{B}=3$ and $\mathcal{B}=8$. The solid-angle density of recorded events $dN/d\Omega$ in the vicinity of the Crab is shown as a function of $\psi^{2}$. Bins of $\psi^{2}$ have constant solid angle (in the small-angle approximation), so any remaining cosmic-ray background shows up as a flat component and the gamma rays are evident as a peak near $\psi^{2}=0$. The improvement in angular resolution for larger events is clear.

Fits to this functional form of Equation 4 can have highly coupled parameters. It is more useful and traditional to quantify the resulting fits with the 68% containment radius, $\psi_{68}$, the angular radius around the true photon direction in which 68% of events are reconstructed. Figure 9 shows $\psi_{68}$, for each $\mathcal{B}$ of the analysis, measured on the Crab and predicted from simulation. At best, events are localized to within 0.17∘, the best angular resolution achieved for a wide-field ground array.

Knowing the angular resolution is critical to subsequent steps of the analysis. Figure 8 indicates that the simulated angular resolution is in good agreement with measurements of the Crab Nebula. This is important because the angular resolution of HAWC for objects at declinations above and below the Crab will differ. While the measured PSF at the position of the Crab cannot be easily extrapolated to other declinations, the simulation can be used to predict the shape of the PSF at any declination. Therefore, the data-simulation agreement shown in Figure 8 is an important verification step.

3.4 Cut Selection and Gamma-Ray Efficiency

The two parameters described in Section 2.6, the compactness, $\mathcal{C}$, and PINCness, $\mathcal{P}$, are used to remove hadrons and keep gamma rays. Events are removed using simple cuts on these variables and the cuts depend on the size bin, $\mathcal{B}$, of the event. The cuts are chosen to maximize the statistical significance with which the Crab is detected in the first 337 days of the 507-day dataset. Concerns of using the data itself for optimizing the cuts are minimal with a source as significant as the Crab.

Table 2 shows the cuts chosen for each $\mathcal{B}$ bin. The rates of events across the entire sky going into the 9 bins, after hadron rejection cuts, vary dramatically, from $\sim$500 Hz for $\mathcal{B}$=1 to $\sim$0.05 Hz for $\mathcal{B}$=9. Figure 10 shows the predicted efficiency for gamma rays (from simulation) along with the measured efficiency for hadronic background under these cuts. The efficiency of photons is universally greater than 30% while keeping, at best, only 2 in $10^{3}$ hadrons. The efficacy of the cuts is a strong function of the event size, primarily because larger cosmic-ray events produce many more muons than gamma-ray events of a similar size.

The limiting rejection at high energies is better than predicted in the sensitivity design study (Abeysekara et al., 2013). The original study was conservative in estimating the rejection power that HAWC would ultimately achieve. With more than a year of data, we now know the hadron rejection of the cuts and can accurately compute the background efficiency.

4 Spectral Fit

Knowing the angular resolution and the background in each $\mathcal{B}$, the energy spectrum of the Crab Nebula may be inferred from the measured data. Section 4.1 describes the likelihood fit to the data. Section 4.2 describes the resulting measurement, and Section 4.3 describes the systematic errors to which this measurement is subject.

4.1 Likelihood Analysis

The HAWC data is fit using the maximum likelihood approach to find the physical flux of photons from the Crab (Wilks, 1938; Younk et al., 2016). In this approach, the likelihood of observations is found under two “nested” hypotheses where some number of free parameters are fixed in one model. This approach can be used to conduct a likelihood ratio test by forming a test statistic, TS, that indicates how likely the data is under a pure background hypothesis or to test the improvement of having additional free parameters in the functional form of the hypothesis spectrum.

The likelihood function is formed over the small (on the scale of the angular resolution) spatial pixels within 2 degrees of the Crab. Each pixel, $p$ has an expected number of background events of $B_{p}$ and, for a specific flux model, an expected number of true photons $S_{p}(\vec{a})$, where $\vec{a}$ denotes the parameters of our spectral model of the Crab. The predicted photon counts fall off from the source according the assumed point spread function. The likelihood $\mathcal{L}(\vec{a})$ is then the simple Poisson probability of obtaining the measured events in each pixel, $M_{p}$ under the assumption of the flux given by $\vec{a}$. The $\mathcal{B}$ dependence of each term in Equation 6 is suppressed.

[TABLE]

Specifically, we fit a differential photon flux $\phi(E)$ of the log parabola (LP) form:

[TABLE]

Here, $\phi_{0}$ is the flux at $E_{0}$, $\alpha$ is the primary spectral index and $\beta$ is a second spectral index that governs the changing spectral power across the energy range of the fit. In this formulation, $E_{0}$ is not fitted but is chosen to minimize correlations between the free parameters in the fit. When fit to an LP function, $E_{0}=7~{}{\rm{TeV}}$ produces good results.

4.2 Fit Results

We find the parameters for $\vec{a}$ that maximize the likelihood function under signal and background hypotheses and quantify the error region of $\vec{a}$ using Wilks’ Theorem (Wilks, 1938). Figure 11 shows the corresponding spaces of $\alpha$, $\beta$ and $\phi_{0}$ for the LP fit that are consistent with HAWC data at 1 and 2$\sigma$. The maximum likelihood occurs at $\alpha=2.63\pm 0.03$, $\beta=0.15\pm 0.03$, and log${}_{10}(\phi_{0}~{}{\rm{cm}^{2}}~{}{\rm{s}}~{}{\rm{TeV}})=-12.60\pm 0.02$. At this best flux the TS, compared to the background-only hypothesis, is 11225, a more than 100$\sigma$ detection.

The TS between an unbroken power-law hypothesis (with $\beta=0$) and the full LP fit is 142, so the spectrum is inconsistent with an unbroken power law at 12$\sigma$.

We quantify the energy range of this fit two ways. First, we take the spectral fit solution and compute the lower 10% quantile of true energy for $\mathcal{B}$=1 and the upper 90% quantile of true energy for $\mathcal{B}$=9. These are 375 GeV and 85 TeV respectively. This is the energy range over which, under the fitted hypothesis, most of the measured photons are expected to lie.

A more conservative approach can be made focusing on the lowest and highest energies where HAWC data could definitively reveal a sharp cutoff in the spectrum. To do this, we separately fit functions of the forms:

[TABLE]

and

[TABLE]

to find the highest $E_{\rm{low}}$ and lowest $E_{\rm{high}}$ that are, at 1$\sigma$, inconsistent with the HAWC observation. With this approach, we believe that we have positive detection of photons from the Crab between 1 and 37 TeV. This is not to say that higher or lower energy photons cannot be a part of the HAWC observation, but using the event size $\mathcal{B}$ to measure the energy of photons limits the dynamic range of the observation. Other sources at other declinations may yield different answers.

4.3 Systematic Errors

Table 3 summarizes the major systematic errors contributing to the measurement of the Crab spectrum with HAWC. These systematic errors have been investigated by computing the spectrum from the Crab under varying assumptions to study the stability of the results under perturbation of the assumptions.

For spectral measurements, a systematic error in three quantities is shown: the overall flux, the spectral index measured and the energy scale. The errors are summed in quadrature to arrive at a total systematic error. In addition to these systematic errors, a systematic error in the absolute pointing of the instrument has been studied.

4.3.1 Charge Resolution and Relative Quantum Efficiency

The charge resolution is a quantity that captures by how much individual PMT charge measurements can vary, for fixed input light, and is estimated to be 10-15% from studies using the HAWC calibration system. Additionally, PMTs vary in their photon detection efficiency by 15–20%. These factors are not simple numbers, but vary for different light levels in the detector and can change with the arrival time distribution during air showers. Varying these assumptions, the $\mathcal{P}$ and $\mathcal{C}$ change and impact the event passing rates, impacting the spectrum.

4.3.2 PMT Absolute Quantum Efficiency

PMTs have an efficiency for converting photons impinging on their surface into PEs detected by the PMT, typically between 20-30%. Of course a single “efficiency” number vastly simplifies the situation: the efficiency is divided between the efficiency for producing a PE and the efficiency for collecting a PE, varies across the face of the PMT, and is wavelength dependent. Additionally, the absorption of the water itself is wavelength dependent. Much of this is modeled, but the simulation carries uncertainties in the treatment and is difficult to validate. The calibration system, in particular, cannot yield the absolute PMT efficiency because it requires establishing the efficiency of the calibration system’s optical path to the PMTs much more precisely than is known. Furthermore, the laser for the calibration system is green light and must be extrapolated for application to blue Cherenkov light.

Instead, the absolute efficiency is established by selecting vertical muons in HAWC tanks by their timing properties. Vertical muons are typically minimum ionizing with a relatively constant energy loss. The simulated response to vertical muons is scaled to match data. Nevertheless, we estimate a $\pm 10\%$ uncertainty in the absolute PMT efficiency which propagates to the errors in the spectra of sources.

4.3.3 Time Dependence, PMT Layout and Crab Optimization

The HAWC instrument has changed over time. The main change is that PMTs or channels are occasionally removed during maintenance. Repaired channels have been re-calibrated and the calibration constants have been changed occasionally for other reasons. The event size bins, $\mathcal{B}$, are based on fractions of available PMTs to mitigate the impact of the varying numbers of PMTs. Furthermore, a single simulation with a single representative PMT layout is used to model the detector and this simplification results in a corresponding systematic error.

Different detector layouts were simulated to bound the impact of sporadic tubes being added or removed. As confirmation, the passing rate of background cosmic rays through the photon/hadron discrimination cuts was studied and shows comparable drifts to the simulation studies.

Additionally, the cuts used in the analysis were established by maximizing the statistical significance of the observations of the Crab Nebula during the first 337 days of data. While strictly not an a priori measurement, the Crab is strong enough in HAWC data that there is very little bias to the final measurement from this optimization process (and no bias for other weaker sources).

To investigate any potential bias, the data was divided in two pieces, a 337-day and a 170-day dataset. The 337-day dataset corresponds to the period over which the cut optimization was done and is the only dataset that could have an over-optimization bias. The Crab spectrum was then measured separately in each dataset. The fitted spectra differ by $\pm$10% in the flux and $\pm$0.1 in the spectral index, similar to what is expected from the varying number of PMTs.

It is unclear whether the different Crab spectra in the 337-day and 170-day datasets are due to overtuning the Crab or the changing detector later in the data-taking. They are similar size effects. Whatever the origin of the effect, it is a sub-dominant, but non-negligible, systematic error.

4.3.4 Angular Resolution

The chief uncertainties in the angular resolution arise from a mismatch between the data and the simulation and spectral dependence of the angular resolution. The impact of angular resolution has been studied by reconstructing the Crab spectrum under different angular resolution hypotheses.

4.3.5 Late Light Simulation

The single largest source of systematic error is how late light in the air shower is treated. Simulation suggests that the arrival time distribution of PEs at the PMTs should be well within $\sim$10 ns. Nevertheless, the distributions of $\mathcal{C}$ and $\mathcal{P}$ in background cosmic rays, as well as the raw PE distributions themselves, suggest some mis-modeled effect above about 50 PEs.

Dedicated studies of late light, using the the calibration system, have been seen to extend ToT measurements, thereby distorting the measured charge in PMTs, but an arrival time distribution much wider than expected from simulation is needed to explain the data.

Efforts to understand this systematic are aimed at measuring the entire PMT waveform for a sample of PMTs to better understand the arrival time distribution of PEs without requiring simulation. It is likely that this systematic error will be better understood in the future, but currently it dominates.

4.3.6 Absolute Pointing

The absolute pointing error is estimated to be no more than $0.1^{\circ}$ for sources that transit above 45∘ in HAWC. It is estimated to be no more than $0.3^{\circ}$ for higher-inclination sources.

The absolute pointing of the HAWC instrument is impacted by the timing calibration as discussed in Section 2.2. Each PMT has a calibrated offset to account for different cable lengths and other timing delays. These offsets are established coarsely by repeated reconstructions to force the peak of maximum cosmic-ray density to be overhead. With this correction, the location of the Crab is within 0.2∘ of its true location. A final detailed alignment is performed to put the Crab Nebula in its known location. The Crab itself must be in the correct location, by construction.

The absolute pointing error on other sources has been studied two ways. First, the Crab location has been fit using only events in bands of reconstructed zenith angle. The Crab location drifts by no more than $0.1^{\circ}$ up to a zenith angle of 45∘. Above that inclination, the Crab is weakly detected and we cannot independently demonstrate better than $0.3^{\circ}$ absolute pointing error. Furthermore, the Crab location has been reconstructed separately using data from each of the 9 $\mathcal{B}$ bins and they agree to within $0.1^{\circ}$.

Finally, other bright known sources, the blazars Markarian 421 and Markarian 501, agree with their known locations to within 0.1∘.

5 Discussion

5.1 Comparison to Other Experiments

Figure 13 shows the Crab spectrum measured with HAWC between 1 and 37 TeV compared to the spectrum reported by other experiments. It is consistent with prior measurements within the systematic errors of the HAWC measurement.

A number of improvements to the HAWC measurement — most notably the inclusion of a proper energy reconstruction — will reduce the systematic errors and increase the dynamic range of our measurements. The observation of the Crab validates this analysis for subsequent application to the sources across the rest of the sky.

5.2 Performance Figures

Figure 14 shows the effective area of HAWC, in this analysis, to photons arriving within 13∘ from overhead. The effective area is defined as the geometrical area over which events are detected, convolved with the efficiency for detecting events. The exact conditions for a photon to be considered detected are complicated for this analysis because, since we are performing a likelihood fit of all pixels in the vicinity of the Crab, photons that are poorly reconstructed play some role in the analysis. In order to have a well-formed effective area, we consider only photons defined within the 68% containment radius, from Table 1. For comparison, Figure 14 includes the progression of cuts, from the effective area without any photon/hadron discrimination without a strong angular accuracy cut to the full analysis cuts. The effective area can exceed the geometrical area of the instrument (about 2$\times 10^{4}$ m2) because events with a core location off the detector occasionally pass the imposed cuts.

This computation of the effective area is lower than in (Abeysekara et al., 2012). For this analysis, developed for steady, multi-TeV sources, we employ tight angular cuts and have a higher energy threshold than used in the initial design study. Furthermore, the re-triggering, defined by the lower edge of $\mathcal{B}=1$, limits the effective area below 1 TeV.

Figure 2 shows the computed differential sensitivity of HAWC to sources at the declination of the Crab utilizing the procedure of (Abeysekara et al., 2013) with the analysis presented here. A point source of differential photon spectrum $E^{-2.63}$ is simulated and fitted using the full likelihood fit. The flux required to be detected at 5$\sigma$ 50% of the time is shown for each bin, $\mathcal{B}$. The lines for each $\mathcal{B}$ are shown with a width corresponding to the width required to contain 68% of the events under the $E^{-2.63}$ hypothesis. A correction is made to adjust the $\mathcal{B}$ separation to a quarter decade in true energy, and the result is fitted. The sensitivity prediction from (Abeysekara et al., 2013) suffered from uncertain background at the highest energies. Now, with more than a year of data, we know the background precisely and can set the cuts appropriately. This has resulted in a more accurate (and more sensitive) analysis above 10 TeV. Below about 1 TeV, for a number of reasons, the sensitivity is somewhat worse than predicted in the original study. The background is larger than the original simulation-only prediction. Furthermore, in the current analysis we employ a relatively high cut (defined by $\mathcal{B}=1$) so that improperly modeled noise can be neglected.

5.3 Anticipated Improvements

The main limitation of this analysis is the reliance on the number of PMTs, used for the definition of $\mathcal{B}$, to simultaneously constrain the energy of photons, the angular resolution, and photon/hadron efficiency. Figure 3 shows that this is a poor energy estimation with each bin $\mathcal{B}$ spanning roughly an order of magnitude of energy. More critically, an overhead $\sim$10 TeV photon can trigger nearly every PMT in HAWC if the core lands near the center of the detector. Consequently $\mathcal{B}=9$ is an overflow bin of everything above 10 TeV. These limitations can be removed with an event parameter that accounts for the light level in the event and the specific geometry and inclination angle of events. Approaches like this are under development. The planned deployment of a sparse “outrigger” array should further increase the sensitivity to photons above 10 TeV (Sandoval, 2016).

Additionally, the principal systematic error (the modeling of late light) is conservatively estimated here and is being studied using the calibration system. It is likely that the effects of late light will be better modeled in the future.

Finally, the threshold for this analysis is established by including only events where more than 6.7% of the PMTs detect light. The typical number of live, calibrated PMTs is $\sim$1000, corresponding to a threshold of $\sim$70 PMTs. Events with 20–30 PMTs could be reconstructed if the noise could be confidently identified. A relatively high event size threshold is used in this analysis to reduce its dependence on the modeling of noise hits. Planned improvements in the modeling should lower the energy threshold of the spectrum analysis in future studies.

The HAWC instrument is performing well with survey sensitivity exceeding current-generation instruments above 10 TeV, sensitivity which HAWC maintains across much of its field-of-view. The all-sky survey conducted by HAWC probes unique flux space and reveals the highest-energy photon sources in the northern sky. Understanding the Crab gives confidence in the survey results.

We acknowledge the support from: the US National Science Foundation (NSF); the US Department of Energy Office of High-Energy Physics; the Laboratory Directed Research and Development (LDRD) program of Los Alamos National Laboratory; Consejo Nacional de Ciencia y Tecnología (CONACyT), México (grants 271051, 232656, 260378, 179588, 239762, 254964, 271737, 258865, 243290, 132197); L’OREAL Fellowship for Women in Science 2014; Red HAWC, México; DGAPA-UNAM (grants RG100414, IN111315, IN111716-3, IA102715, 109916); VIEP-BUAP; PIFI 2012, 2013, PROFOCIE 2014, 2015; the University of Wisconsin Alumni Research Foundation; the Institute of Geophysics, Planetary Physics, and Signatures at Los Alamos National Laboratory; Polish Science Centre grant DEC- 2014/13/B/ST9/945; Coordinación de la Investigación Científica de la Universidad Michoacana.