New Spectral Evidence of an Unaccounted Component of the Near-infrared Extragalactic Background Light from the CIBER

TL;DR

This study presents the first measurement of the near-infrared extragalactic background light spectrum from the CIBER experiment, revealing a potential new diffuse component beyond known galaxy contributions.

Contribution

It provides new spectral measurements of the EBL in the 0.8-1.7 um range and highlights the need for an additional diffuse component not explained by existing models.

Findings

Measured EBL brightness at 1.4 um as 42.7+11.9/-10.6 nW/m2/sr using the Kelsall model

Using the Wright model results in unphysical EBL values below galaxy counts

Data suggests the existence of a new diffuse component beyond known galaxy light

Abstract

The Extragalactic Background Light (EBL) captures the total integrated emission from stars and galaxies throughout the cosmic history. The amplitude of the near-infrared EBL from space absolute photometry observations has been controversial and depends strongly on the modeling and subtraction of the Zodiacal light foreground. We report the first measurement of the diffuse background spectrum at 0.8-1.7 um from the CIBER experiment. The observations were obtained with an absolute spectrometer over two flights in multiple sky fields to enable the subtraction of Zodiacal light, stars, terrestrial emission, and diffuse Galactic light. After subtracting foregrounds and accounting for systematic errors, we find the nominal EBL brightness, assuming the Kelsall Zodiacal light model, is 42.7+11.9/-10.6 nW/m2/sr at 1.4 um. We also analyzed the data using the Wright Zodiacal light model, which results in a worse statistical fit to the data and an unphysical EBL, falling below the known background light from galaxies at <1.3 um. Using a model-independent analysis based on the minimum EBL brightness, we find an EBL brightness of 28.7+5.1/-3.3 nW/m2/sr at 1.4 um. While the derived EBL amplitude strongly depends on the Zodiacal light model, we find that we cannot fit the spectral data to Zodiacal light, Galactic emission, and EBL from solely integrated galactic light from galaxy counts. The results require a new diffuse component, such as an additional foreground or an excess EBL with a redder spectrum than that of Zodiacal light.

Tables

Table 1. Table 1: CIBER/LRS observed sky coordinates

			J2000		Ecliptic		Galactic		Solar Elong.
Field	Time	Altitude	$RA$	$Dec$	$\lambda$	$\beta$	$l$	$b$	$ϵ$
	[s]	[km]	[deg]	[deg]	[deg]	[deg]	[deg]	[deg]	[deg]
(2nd flight #36.265)
Vega	88–98	152–168	279.24	41.28	286.10	64.20	69.95	20.06	115.77
SWIRE/ELAIS-N1	109–190	185–283	242.84	54.65	209.24	72.38	84.39	44.67	93.19
NEP	210–282	299–329	270.87	66.45	342.70	89.64	96.26	29.46	90.22
Elat10	305–316	328–326	226.97	-2.55	225.25	14.46	356.29	45.71	115.69
Elat30	328–348	323–315	220.86	19.80	211.14	33.75	23.00	63.36	100.35
Bootes-A	358–432	310–243	218.51	34.75	200.68	46.59	58.42	66.79	91.39
Bootes-B	437–517	237–105	217.30	33.25	200.36	44.80	55.08	68.06	91.21
(3rd flight #36.277)
Lockman Hole	125–172	202–265	161.43	58.21	135.42	45.49	149.41	51.97	118.83
SWIRE/ELAIS-N1	194–239	284–315	242.81	54.59	209.32	72.32	84.31	44.71	105.65
NEP	256–309	320–324	270.63	66.28	311.48	89.62	96.06	29.56	89.76
Elat30	330–356	319–306	236.98	9.57	234.38	29.25	18.64	44.89	124.12
Bootes-B	375–425	296–244	217.30	33.27	200.35	44.82	55.13	68.06	132.31
Bootes-A	451–467	211–187	218.38	34.69	200.58	46.49	58.33	66.91	130.73
(4th flight #40.030)
DGL	160–225	272–401	251.99	68.76	154.03	82.84	100.27	36.19	88.59
NEP	240–300	425–505	270.80	66.15	307.86	89.48	95.91	29.49	90.32
Lockman Hole	315–365	520–558	161.23	58.58	135.01	45.74	149.07	51.65	69.33
Elat10	380–430	566–577	190.46	8.12	186.37	11.6	295.67	70.86	110.52
Elat30	445–495	577–562	193.08	28.03	179.76	30.70	110.75	89.08	102.31
Bootes-B	508–563	555–509	217.24	33.28	200.28	44.80	55.18	68.11	113.94
SWIRE/ELAIS-N1	641–696	395–275	242.76	54.68	209.05	72.38	84.45	44.71	102.07

Table 2. Table 2: Verification of the ISL estimate using the source counts data

Field	Totala,ba,bfootnotemark:	IGLaaIntegrated flux for $J>$17 mag in nWm-2sr-1 at 1.25 $\mu$m.	ISLaaIntegrated flux for $J>$17 mag in nWm-2sr-1 at 1.25 $\mu$m.	ISLaaIntegrated flux for $J>$17 mag in nWm-2sr-1 at 1.25 $\mu$m.	Reference
	(1)	(2)	(1)–(2)	Model
Bootes	5.45	4.31$\pm$1.02	1.14$\pm$1.02	1.93	Gonzalez et al. 2010
NEP	9.64	4.31$\pm$1.02	5.33$\pm$1.02	5.52	Oi et al. 2014

Table 3. Table 3: Summary of the EBL results

Wavelength	Nominal EBLaaMean value $\pm$ Statistical error + Systematic error(upper/lower)		Minimum EBLaaMean value $\pm$ Statistical error + Systematic error(upper/lower)
[$\mu$m]	[nWm-2sr-1]		[nWm-2sr-1]
0.80	41.1$\pm$4.8 $+15.7/-15.3$		-
0.83	30.6$\pm$4.4 $+15.8/-14.7$		-
0.86	23.1$\pm$4.2 $+15.1/-14.4$		-
0.90	29.6$\pm$3.2 $+14.1/-14.0$		-
0.95	34.2$\pm$3.3 $+13.8/-13.7$		-
1.00	31.2$\pm$3.5 $+13.4/-13.5$		-
1.05	31.3$\pm$3.2 $+13.1/-12.9$		7.9$\pm$3.8 $+1.5/-0.6$
1.11	36.7$\pm$2.8 $+12.6/-12.6$		15.4$\pm$3.4 $+1.1/-0.8$
1.18	32.5$\pm$2.2 $+13.3/-12.2$		12.9$\pm$3.2 $+4.0/-0.9$
1.25	38.7$\pm$2.1 $+12.6/-11.6$		20.1$\pm$3.1 $+3.9/-1.1$
1.33	41.7$\pm$2.0 $+11.3/-10.9$		24.3$\pm$2.8 $+2.7/-1.3$
1.42	42.7$\pm$2.3 $+11.7/-10.3$		28.7$\pm$3.0 $+4.1/-1.5$
1.51	37.8$\pm$2.9 $+16.4/-10.1$		23.8$\pm$3.5 $+9.9/-1.4$
1.60	41.1$\pm$3.3 $+15.6/-9.3$		27.1$\pm$3.8 $+10.0/-1.5$
1.70	35.7$\pm$2.8 $+13.3/-8.4$		24.8$\pm$3.0 $+7.2/-1.4$

Table 4. Table 4: Total Error Budget for the EBL measurement

Error term	Erroraain nWm-2sr-1 at 1.25 $\mu$m, or percent	Comment
Raw sky	$\pm$1	Detector noise, dark current subtraction
	3% ($\pm$12 at NEP)	Absolute calibration error
Airglow	1	Model fitting error, maximum at SWIRE
	$+$4	In case of no subtraction of the slow airglow component
ZL template	1%	Variation of sky difference spectra
	0.4%	ISL and DGL modeling error
	(1%)bbTotal integrated source counts including IGL and ISL.	Spectral response calibration error
ISL	10% ($\pm$5 at NEP)	Estimated from simulation (Poisson noise, aperture correction,
		and magnitude-cut effect)
	$\pm$2	Modeling error (scaling error to 2MASS counts at brighter magnitudes,
		and uncertainty of galaxy counts at fainter magnitudes)
DGL	20% ($\pm$4 at NEP)	Uncertainty of DGL/100 $\mu$m ratio (Arai et al., 2015)
EBL field mean	$\pm$2	Statistical error and field-to-field variation
(ZL model)	$+$12/$-$11	Systematic uncertainty including absolute calibration,
		airglow modeling error, and ISL and DGL modeling errors
EBL field mean	$\pm$3	Statistical error, field-to-field variation, and scaling error to IGL at 0.8–1.0 $\mu$m
(Minimum)	$+$4/$-$1	Systematic uncertainty including ZL template error,
		spectral calibration error, and airglow modeling error

Equations

$$I_{obs}=I_{sky}+A_t{e}^{-t/\tau }+A_h{e}^{-h/H}$$

Full text

cosmology: observations — dark ages, reionization, first stars — diffuse radiation — infrared: diffuse background — infrared: general — zodiacal dust

1 Introduction

The extragalactic background light (EBL) is the integrated intensity of all photons emitted and absorbed along a line of sight throughout the cosmic history. While the early objects are too faint to detect individually, the EBL as the integrated light of such objects is potentially a powerful probe to search the spectral signatures. As a result, the EBL is an important observable for understanding galaxy formation and evolution (Hauser & Dwek, 2001; Madau & Pozzeti, 2000). Of particular interest, redshifted ultraviolet radiation from primordial galaxies and black holes during cosmic reionization has been proposed as a measurable component of the near-infrared EBL. Such sources will contain a distinctive spectral Lyman break due to absorption by neutral hydrogen in the intergalactic medium and possibly a redshifted Lyman-alpha emission feature, leading to a spectral discontinuity at $\sim$1 $\mu$m for redshifts $z>$ 7 (Santos et al., 2002; Salvaterra et al., 2003; Fernandez & Komatsu, 2006; Yue et al., 2013a).

Measurements of the EBL in the range 1.25 to 5 $\mu$m with Cosmic Background Explorer (COBE) (Hauser et al., 1998; Wright & Reese, 2000; Cambrèsy et al., 2001; Levenson et al., 2007; Sano et al., 2015, 2016), Infrared Telescope in Space (IRTS) (Matsumoto et al., 2005, 2015) and AKARI (Tsumura et al., 2013c) have been difficult to reconcile with the integrated galaxy light (IGL) predicted from deep galaxy counts and contemporary models of galaxy evolution, exhibiting significantly larger surface brightnesses than predicted. Interpretation of this EBL excess has remained controversial, with some identifying it with reionizing sources and others with local foregrounds.

A significant challenge to the measurement of the EBL near the earth is to accurately subtract the foreground zodiacal light (ZL), which is sunlight scattered by interplanetary dust in the plane of the planets. The ZL has been modeled on the basis of the COBE all-sky observations (Kelsall et al., 1998), but since the ZL foreground dominates the EBL by factors $\gtrsim 20$ even a few percent uncertainty in the modeled ZL brightness changes the residual EBL level by large factors. Over the 1.4–4 $\mu$m range reported by the IRTS (Matsumoto et al., 2015), the EBL spectrum is similar to the ZL spectrum (Dwek et al., 2005). A promising approach to search for a reddening in the EBL spectrum, which would break the spectral degeneracy with ZL, is to extend the observation wavelengths to shorter than 1.25 $\mu$m.

In this paper, we report the first spectral measurement of the EBL in the unexplored wavelength range from 0.8 to 1.7 $\mu$m by a sounding rocket experiment, Cosmic Infrared Background Experiment (CIBER) (Zemcov et al., 2013; Bock et al., 2013; Tsumura et al., 2013a; Korngut et al., 2013).

2 CIBER experiment

2.1 The LRS instrument

One of the science instruments of CIBER, the Low Resolution Spectrometer (LRS) is designed to measure the absolute spectrum of the near-infrared sky at wavelengths 0.8–2 $\mu$m. A full description of the LRS is given in Tsumura et al. 2013a. The LRS consists of a 5-cm refractive telescope with five spectral slits providing 2.8-arcminute$\times$5.5-degree fields on the sky. A prism disperses the incident light perpendicular to the slits and the spectra are imaged with a 256$\times$256 HgCdTe detector array. Depending on the wavelength, a spectral resolving power of $\lambda$/$\delta\lambda$=15–30 is achieved.

2.2 CIBER flights

The CIBER experiment has flown four times. The initial experiment flew on 10:45 UTC February 25, 2009 at White Sands Missile Range (WSMR), New Mexico, USA, on a two-stage Terrier-Black Brant IX rocket 111 National Aerospace Administration (NASA) Sounding Rockets Program Office, The NASA Sounding Rocket Program Handbook (Wallops Island, VA: NASA). Scattered thermal emission from the rocket skin (Zemcov et al., 2013) affected the LRS data, so these are not used in this EBL analysis, though we did identify a new spectral feature in the ZL from the data taken during this flight (Tsumura et al., 2010).

Following the first experiment, we refurbished the instrument to improve the baffling performance to scattered thermal emission. The second experiment launched at 4:50 UTC July 11, 2010 from WSMR on a Terrier-Black Brant IX rocket. As a result of the instrument improvements and performance of the vehicle, this flight generated the highest quality data of all the CIBER flights.

For the next flight, the experiment was modified to measure polarization of the ZL. The third experiment flew on 9:00 UTC March 22, 2012 from WSMR on a Terrier-Black Brant IX rocket. Though we successfully measured the near-IR ZL polarization with these data, the photometric accuracy was not as good as that of the second flight due to strong airglow contamination. We exclude this data set for this EBL analysis, though we plan to present ZL polarization results in a future paper.

For the fourth and final experiment, we removed the polarizer and modified the LRS to match the second flight configuration. The fourth flight launched 3:05 UTC June 6, 2013 from Wallops Flight Facility (WFF), Virginia, USA, on a four-stage Black Brant XII rocket, which could carry the payload to much higher altitudes than the two-stage rocket used in the previous three flights. In this experiment we obtained longer observation times on all science fields, and could evaluate residual airglow emission from the upper atmosphere over long time scales.

2.3 Observations

Table 1 gives a summary of the science fields at varying equatorial, ecliptic and Galactic coordinates, and solar elongations in the three CIBER flights. The coordinates were determined from the positions of 2MASS (Two Micron All Sky Survey) stars matched to the positions of stars observed with the LRS (Kim et al., 2017).

For all the CIBER experiments, most of the observation time was devoted to observing the core science fields, SWIRE/ELAIS-N1, NEP, and Bootes-A/B, where the Spitzer and AKARI satellites carried out deep galaxy surveys and relevant ancillary data from ground-based telescopes are available (Bock et al., 2013). In order to investigate the ZL foreground, we observed the skies at low ecliptic latitudes, Elat-10 and Elat-30, where the ZL dominates the sky brightness. The Lockman hole, the lowest cirrus region, was also observed in the third and fourth flights, and a relatively high-cirrus region to investigate the diffuse Galactic light (DGL) was observed in the fourth flight. However, these additional fields were also useful for airglow investigation, because these observations were done at low altitudes in the early phase of the flight.

To minimize contamination by atmospheric emission and sun light, we set the launch date and time window so that the earth and sun avoidance angles of the science fields were greater than 40 and 115 degrees, respectively. Moreover, in order to compare our data directly with the COBE/DIRBE data with little uncertainty, we added another constraint on the solar elongations to be in the range covered by the COBE/DIRBE observations, 64–128 degrees (Hauser et al., 1998).

The actual coordinates of the observation were slightly displaced from the targeted coordinates of the rocket’s attitude control system, due to a modest misalignment between it and the LRS telescope. The differences are at most 0.5 degrees or one tenth of the field-of-view (FOV) of the LRS instrument. The pointing stability is better than 10 arcsec for an exposure time of 50 seconds, which is negligibly small compared to the spatial resolution of the LRS (84 arcsec). Accurate image registration using known 2MASS stars has been performed (Kim et al., 2017). The central LRS coordinates of each field are given in Table 1.

3 Data reduction & analysis

3.1 Image processing

The LRS instrument employs a 256$\times$256 HgCdTe focal plane array (FPA) to record the intensity of the dispersed light from a prism through five slits. Figure 1 is an example of the dispersed light image on the FPA with the photocurrent signal normalized to the peak intensity. Sky spectra are calculated from such spectral images using the following process.

The signal output of the FPA is sampled with multiple reads, recorded for each observation at a frame rate of 4 Hz. The intensity of the image in Figure 1 is directly proportional to the detector signal current, which is calculated from the slope of the accumulated charge in each pixel by line fits to time series data. In this image, stars are easily distinguished as bright spots dispersed horizontally. We masked all the detector pixels illuminated by stars brighter than about 12 mag in $J-$band and measured the average spectra of the diffuse sky emission using unmasked pixels. The star mask for each slit was generated by 2$\sigma$ clipping each row after summing the row pixels to improve the signal-to-noise. The source extraction is complete for stars brighter than 12 mag (Arai et al., 2015). In this analysis we also exclude the data from the centermost of the five slits as it exhibits anomalous behavior in the early part of the charge integration ramp, which arises from a transient behavior of the readout circuit at corresponding array pixels.

The spectrometer leaves a border of dark pixels that we use to continuously monitor the dark current. To obtain the net photocurrent, the dark current estimated from the masked area is subtracted from the illuminated area of the FOV using an image scaling from shutter-closed data. The dark current difference between the illuminated area and the masked area is estimated to be $\sim$0.03 e*-s-1*, which corresponds to a negligible surface brightness of $\sim$0.7 nWm*-2sr-1* (Arai et al., 2015).

3.2 Generation of sky spectra

The photocurrent spectrum for each slit is calculated by averaging all the photocurrent signals in spectral columns, the width of each bin corresponding to the spectral resolution, after removing the pixels of bright stars. The dispersion relation between pixel position and wavelength was measured in the laboratory (Tsumura et al., 2013a). The photocurrent spectra are converted to sky brightness units by a conversion factor for each wavelength determined by a calibration with a diffuse light source in the laboratory (Zemcov et al., 2013; Arai et al., 2015). Finally, the sky spectra for the four outer slits are averaged to obtain the final sky spectra. The final spectral binning from the native wavelength sampling of each pixel is performed so that each spectral bandwidth is equivalent to the spectral resolution, which is determined by the slit-apodized PSF and the prism dispersion, undersampling the spectrum, but reducing correlation between data points. Sky spectra from the second flight are shown in Figure 2.

3.3 Environmental emission

3.3.1 Thermal emission from the rocket skin

During flights we observe scattered thermal emission from the rocket skin corresponding to a temperature of $T\sim$400 K (Zemcov et al., 2013). This thermal emission was clearly observed in the first flight, and its origin was identified in laboratory tests (Tsumura et al., 2013a). For the second to fourth flights, the thermal emission was greatly reduced by improving the baffle performance, but still significant at wavelengths longer than 1.8 $\mu$m. For this work we used only unaffected data at wavelengths shorter than 1.7 $\mu$m, at which the brightness of the thermal emission is estimated from the blackbody fit to the data to be lower than 2 nWm*-2sr-1*, similar to the statistical error of the EBL measurement described later.

3.3.2 Airglow emission

We observe emission from atmospheric airglow and outgassing from the payload at altitudes lower than 250 km. As detailed in Figure 3, the emission spectrum of airglow was derived by splitting a field observation into two measurements, and taking the difference to isolate a time- or altitude-dependent signal. The initial field early in the observation sequence of each flight, e.g., the SWIRE/ELAIS-N1 field for the second flight, shows a strong broad peak around 1.6 $\mu$m, which is plausibly due to vibration transitions of OH molecules smoothed by the spectral resolution of the LRS. It should be noted that the airglow spectra are quite different from the spectra of either the isotropic background or the ZL to be presented in later sections.

We modeled the airglow emission based on the time and altitude dependence of the intensity in the fourth flight and a template of the emission spectrum, as described in Appendix A. For the second flight data, the difference between the spectra including and excluding a slow-decay airglow component gives a conservative estimate for the airglow systematic error. For the fourth flight data, the modeling error of the airglow is accounted in the statistical error of the sky brightness.

3.4 Measured sky spectra

Figure 4 shows the measured diffuse sky spectra of the five science fields after masking detected stars. Foreground emission must be subtracted from the measured absolute sky brightness, including the integrated star light (ISL), DGL, and ZL. The measured sky levels are similar to the photometric data of COBE/DIRBE instrument and IRTS in low-foreground regions(Kelsall et al., 1998; Matsumoto et al., 2005), but with deeper point-source removal due to CIBER’s higher angular resolution. The absolute photometric error of the sky brightness is described in the following sections (see Arai et al. 2015 for further details).

4 Photometric calibration

4.1 Calibration in the laboratory

We calibrated the LRS response to surface brightness in the laboratory before and after the flights using a range of light sources coupled to an integrating sphere to generate diffuse light. The radiance from the sphere was absolutely calibrated using an absolutely calibrated spectrograph provided by the National Institute of Standards and Technology (Zemcov et al., 2013; Tsumura et al., 2013a; Arai et al., 2015). The spectral calibration was performed using a tunable laser and a monochrometer with an accuracy of $\sim$1 nm for all four flights (Arai et al., 2015). The calibration factor, which is the conversion between the mean pixel photocurrent and surface brightness at each wavelength, is consistent within a 3% rms variation for the second and third flights, which we adopt as the overall uncertainty of the surface brightness calibration (Arai et al., 2015).

4.2 In-flight calibration from stars

The pre-flight calibration of the surface brightness for the fourth flight does not match the calibration determined for the second and third flights222Though the third flight employed a polarizing grid over four of the five slits, the central slit remained clear and measurements derived from it can be directly compared against the second and fourth flight calibrations.. It is reasonable to assume that the calibration does not change much through all the flights, because the optical system is common and there is no evidence that the intrinsic responsivity of PICNIC detectors changes appreciably with time. It appears that there was some failure in the calibration system used for the pre-flight calibration for the fourth flight, a situation which we did not identify until after the fourth flight, which did not allow payload recovery due to the higher altitude. To utilize the fourth flight data, we apply the previous calibration and verify that the LRS-observed fluxes of 2MASS stars are consistent with the catalogued fluxes as described below.

Though the rocket’s attitude control system (ACS) aligns the LRS pointing to within $15^{\prime}$ of each target, we require a finer understanding of the instrument pointing to proceed. We match bright detected stars to the 2MASS catalogue (Skrutskie et al., 2006) to generate the astrometric registration. To avoid mismatches arising from source confusion, we exclude stars fainter than 10 mag in $J$-band, stars that are next to brighter 2MASS stars within 500 arcsec, stars that have negative $J$–$H$ color, variable stars, and stars that have non-continuous spectra that could be due to instrumental and/or source confusion noise. Following these cuts, 105 LRS-detected sources across all fields in all flights are cross-matched with stars appearing in the 2MASS catalog. Details of the astrometric registration and star spectra are described in a separate paper (Kim et al., 2017).

A comparison between the star fluxes from the 2MASS catalog and the CIBER observations in all flights is shown in Figure 5. A tight correlation, with a slope consistent with unity and with an offset consistent with zero, are apparent in both the $J$ and $H$ bands after correction of the slit apodization. This indicates that most of the flux detected by the LRS is identical to the flux predicted from the 2MASS catalog. Mean ratios of the LRS measured fluxes to the 2MASS catalog fluxes in $J$ and $H$ bands are 0.45 and 0.46, respectively, and the difference from unity can be explained by the apodizing effect of the slit mask (Kim et al., 2017). Our Monte-Carlo simulation of the stellar photometry, using actual slit function of the LRS and taking the source confusion and instrumental noise into account, predicts the flux range of LRS measured fluxes well (Kim et al., 2017). The measured slope between CIBER and 2MASS is 1.01(1.04) for the second and third flights, and 1.03(1.01) for the fourth flight at $J$($H$) band, validating our assumptions for this flight’s calibration.

We could also check the calibration over the entire spectral band of the LRS by star observations. Several stars were accidentally observed in common over more than one flight, and the spectra showed identical spectral shape within the measurement error, indicating consistency of the calibrations among all flights (Kim et al., 2017).

5 Astrophysical foreground analysis

5.1 Integrated star light

To estimate the surface brightness of the ISL, we begin with the $J$, $H$ and $K$ band flux of the 2MASS sources (Skrutskie et al., 2006) within the LRS’s field of view (FOV), taking into account the survey completeness and the effective slit width (Arai et al., 2015). The FOV of the LRS is determined from the slit width and the point spread function measured in the laboratory using monochromatic light sources (Zemcov et al., 2013; Tsumura et al., 2013a; Arai et al., 2015). To simulate the completeness magnitude and calculate the continuous spectrum of ISL, we carried out additional simulations by the same method, but using a model of star counts based on stellar population synthesis, TRILEGAL (Girardi et al., 2005). The simulated counts of bright stars in the 2MASS bands agree with the 2MASS counts up to the detection limit of the LRS ($J\sim$13 mag) as shown for the NEP field in Figure 6, and with our observed star counts within Poisson errors as shown in Figure 7. The results in the other fields are similar. We also confirm that the ISL spectrum from the star count model is consistent with the $JHK$ color of the ISL from 2MASS for 9 $<J<$ 17 to within 2% accuracy.

We extend the ISL calculation to magnitudes fainter than the limiting magnitude of 2MASS $AB$=17 mag using the same star count model, but scaled by a factor of 1.2 to fit the 2MASS counts in the magnitude range from $J$=13 to 17 mag. The validity of the $\times$1.2 scaled model at fainter magnitudes than 17 mag was checked using deep star counts complete to $AB$=20 mag in $J$-band towards Bootes and NEP with NEWFIRM (Gonzalez et al., 2010) and WIRCam/CFHT (Oi et al., 2014), respectively. The raw source counts (including galaxies) towards Bootes are compared with the model star counts and also with well-established galaxy count data (Keenan et al., 2010) as shown in Figure 8. We confirm that the sum of the model star counts and the galaxy counts can account for the observed total counts. In Table 2, the model prediction for the stars fainter than the 2MASS detection limit is compared with the deep counts for Bootes and NEP. From this result, the uncertainty of ISL at faint magnitudes is estimated to be 1.5 nWm*-2sr-1* in $J$-band, dominated by uncertainties of the galaxy counts (Domínguez et al., 2011) in separating stars from galaxies.

The resulting ISL spectrum at the SWIRE/ELAIS-N1 field is shown in Figure 4. The ISL spectra in the other fields are identical to this even though the brightness varies by factor of 3. The total uncertainties of ISL at our observation fields are estimated to be 5–10%, which consists of the Poisson error, slight color difference between the star count model and the 2MASS counts, and uncertainty in the contribution of galaxies at fainter magnitudes (Keenan et al., 2010). The total systematic uncertainty of the ISL estimate propagated to the final EBL estimate is described in Section 6.5.

5.2 Diffuse Galactic light

The DGL, starlight scattered by interstellar dust, correlates with the far-infrared brightness of interstellar dust emission in optically thin regions (Brandt & Draine, 2012). The DGL spectrum was previously measured with CIBER by correlating the brightness over each field with the 100 $\mu$m intensity (Arai et al., 2015). The field-averaged DGL spectrum exhibits a continuum approximated by a model spectrum as shown in Figure 4. The LRS-measured DGL spectrum is consistent with the result in the visible with Pioneer10/11 (Matsuoka et al., 2011), with COBE/DIRBE (Sano et al., 2015), and with AKARI at wavelengths longer than 1.8 $\mu$m (Tsumura et al., 2013b). All the DGL measurements fit a DGL model spectrum from Zubko et al. (2004) with a scaling by a factor of 2 (Arai et al., 2015; Sano et al., 2015). In this work, we use the best-fit scaled Zubko et al. DGL model to estimate the DGL brightness according to the 100 $\mu$m brightness at each science field. We evaluate the uncertainty of our DGL estimate as 20% from the model scaling error to the measured DGL brightness, which does not appreciably affect the accuracy of the EBL measurement (Arai et al., 2015). The propagation of the DGL uncertainty to the EBL estimate is described in Section 6.5.

5.3 Zodiacal light

To obtain the final EBL spectrum, we subtract a smooth ZL component using two models based on COBE/DIRBE all-sky measurements (Kelsall et al., 1998; Wright, 2001). To confirm the applicability of the ZL model, we correlate the ISL- and DGL-subtracted surface brightness at DIRBE $J$-band with the Kelsall ZL model (Kelsall et al., 1998), shown in Figure 9. We observe a tight correlation between the data and the ZL model, consistent with the estimated uncertainties from statistical, systematic, and modeling errors. We also carry out the same test using the alternative Wright ZL model (Wright, 2001). The Wright ZL model gives a poorer fit to the data, but due to the measurement’s systematic errors, this model cannot be excluded.

To apply the ZL model calculated at $J$-band to the other LRS wavelengths, we scale a spectral template of the ZL. This template is computed from the ISL- and DGL-subtracted spectra by forming the pair-wise difference of the spectra, normalizing the spectra, and then computing the mean of all such spectral differences in a given flight with weights calculated from the statistical errors. The field differencing scheme effectively cancels out isotropic components of the emission. Most of the data points in the resulting ZL template spectra are consistent with each other within the uncertainties, as shown in Figure 10. Some differences between the flights beyond the expected errors may be due to systematic errors in the gain calibration of the LRS. In this work, therefore, we apply the ZL template measured from a given flight to its associated flight data to avoid increasing the effect of calibration error in the ZL subtraction. No ZL model uncertainty is explicitly propagated to the derived EBL brightness, though the absolute instrument calibration error is taken into account. Rather, the ZL modeling error is captured through the comparison of the Kelsall and Wright models.

6 Extragalactic background result

6.1 Nominal EBL spectrum

The nominal EBL spectrum computed from the combined second and fourth flight data after subtracting the Kelsall ZL model template spectrum is shown in Figure 11 and summarized in Table 3. The 1$\sigma$ statistical uncertainty is indicated by the error bars on the data points, while the total systematic uncertainty of the nominal EBL spectrum is indicated by an error band. The nominal near-IR EBL surface brightness is similar to or slightly lower than the DIRBE (Cambrèsy et al., 2001; Sano et al., 2015, references of the other previous works are therein.) and IRTS (Matsumoto et al., 2005, 2015) measurements derived using the Kelsall ZL model, but higher than the DIRBE measurement derived using the Wright ZL model (Wright, 1998; Levenson et al., 2007). We find a nominal EBL surface brightness of 42.7${}^{+11.9}_{-10.6}$ nWm*-2sr-1* at 1.4 $\mu$m. This is out of allowed range of the EBL of 15$\pm$5 nWm*-2sr-1* at 1.4 $\mu$m derived from HESS (High Energy Spectroscopic System) gamma-ray absorption spectra (H.E.S.S. Collaboration, 2013). Our derived EBL spectrum, which decreases toward visible wavelengths from a peak around 1.5 $\mu$m, cannot be attributed solely to residual ZL foreground, as the spectrum of the EBL is redder than the ZL component shown in Figure 10.

6.2 EBL with Wright ZL model

We must evaluate our nominal result given the large uncertainties of the ZL models. We also calculate the EBL assuming the Wright ZL model with a ”very strong no zodi principle” (Wright, 2001), and confirm the EBL is 7.0${}^{+13.1}_{-12.6}$ nWm*-2sr-1* in an equivalent DIRBE band at 1.25 $\mu$m, which is consistent with Levenson’s et al. EBL of 21.3$\pm$15.1 nWm*-2sr-1* within the errors. However, the derived EBL is negative for $\lambda$ $<$ 1.0 $\mu$m and less than the IGL for $\lambda$ $<$ 1.3 $\mu$m. The EBL derived from the Wright model is unphysical, it is also redder than the ZL spectrum. As previously noted, the Wright model does not fit the observed field difference data as well as the Kelsall model.

6.3 Minimum EBL spectrum

It is instructive to consider the EBL spectrum that results from making the most aggressive assumptions about the surface brightnesses of the ZL component. To obtain a conservative minimum estimate of the 1.0–1.8 $\mu$m EBL, we set the 0.8–1.0 $\mu$m signal levels to the IGL brightness as follows. We subtracted the ZL spectrum with an arbitrary scaling to eliminate any excess EBL at 0.8–1.0 $\mu$m and to match the IGL from galaxy counts (Domínguez et al., 2011). The scaling would be interpreted as an isotropic ZL component, but its physical origin is unclear.

The mean ”minimum EBL” spectrum is shown in Figure 11 (also in Table 3), with errors computed from the field-to-field variance and the statistical error of the derived EBL. The ZL brightness subtracted to obtain the minimum EBL is $\sim$5% higher than the Kelsall ZL model brightness at high ecliptic latitudes. Due to the blue color of the ZL spectrum, the derived EBL remains large at longer wavelengths: at 1.4 $\mu$m, the minimum EBL brightness is 28.7${}^{+5.1}_{-3.3}$ nWm*-2sr-1*. This exceeds the EBL upper limits from high-energy gamma-ray observations (H.E.S.S. Collaboration, 2013). Our result indicates that underestimating the ZL foreground (Dwek et al., 2005) cannot eliminate an excess above IGL near 1.5 $\mu$m. If the minimum EBL calculation is the true EBL level, the EBL must have a very red color increasing from 1.0 $\mu$m to 1.7 $\mu$m that is quite different from the IGL spectrum.

6.4 Field-to-field Variation of the EBL spectra

The variation of the EBL spectra over the five fields observed in the second flight using the Kelsall ZL model are shown in Figure 12. The brightness difference between the fields is only a few nWm*-2sr-1* at most wavelengths, except at short wavelengths which have relatively large statistical errors. The uncertainty-weighted mean of the spectra is adopted as the nominal EBL result, and the field-to-field variation is accounted as the statistical error of the final EBL result shown in Figure 11.

6.5 Systematic uncertainties of the EBL result

6.5.1 Airglow subtraction uncertainties

As described in Appendix A, we estimate airglow contamination using a double-exponential function. The uncertainty in the fit parameters is propagated through the analysis and is accounted as part of the total statistical uncertainty of the sky brightness measurement. Although we subtract both fast and slow decay airglow components from the second flight data using a model based on the fourth flight data, the evidence that the second flight data require a slowly decaying component is weak. We account for this possibility as the major systematic error in the airglow subtraction; difference between the EBL estimates with and without subtraction of the slow airglow component.

6.5.2 ISL and DGL Uncertainties

We determine the susceptibility of the EBL measurement to differences in the assumed ISL by changing the ISL amplitude by $\pm$10% from its nominal value. The ZL template spectrum is affected as shown in Figure 13. The total effect on the EBL estimate, including the change from the modified ZL template, is shown in the left panel of Figure 14. We estimate the ISL intensity error from uncertainties in the galaxy counts described in Section 5.1 to be 10%. The spectral uncertainty of the ISL is negligible compared to the absolute intensity error. We assign a systematic error of 5 nWm*-2sr-1* at 1.25 $\mu$m from the uncertainties associated with ISL.

In a similar fashion, we estimate the error in the EBL measurement from the DGL model by varying the DGL brightness $\pm$20% from its nominal value, as shown in the right panel of Figure 14. Again, the effect of the DGL on the ZL template spectrum is propagated through to the EBL estimate, as shown in Figure 13. In this case, the systematic error from DGL propagated to the mean EBL brightness is approximately 4 nWm*-2sr-1* at 1.25 $\mu$m.

6.5.3 Comparison of the EBL Results between Different Flights

EBL measurements from different flight data should give the same results, though data quality and systematic errors do vary from flight to flight. As discussed in Section 4, since the fourth flight laboratory calibration is not reliable, we use the third flight calibration factor for the fourth flight data. Note that the difference between the calibration factors of the second and third flights is only 1% in the entire LRS bands.

We compare the EBL results for each of the three flights with one another In Figure 15. As we expected, the fourth flight result agrees with the second flight result within the error bars, except for a dip in the fourth flight data at 1 $\mu$m, where the measurement error is large due to presence of the airglow systematic (Appendix A.2). The signal-to-noise ratio of the third flight data is comparably smaller because only the non-polarized channel, which corresponds to one of five slits on focal plane mask, is available for measuring the total brightness. The resulting spectrum has errors a factor of $2.2$ larger than the others, and we observe a large discrepancy from the other flights data near 1 $\mu$m.

6.5.4 Total Error Budget

The total error budget of our EBL measurement referenced to 1.25 $\mu$m is enumerated in Table 4. The statistical and systematic errors of the EBL estimate at other wavelengths are shown graphically in Figure 11. The ZL model error is not included in the final EBL error, but is bracketed by the difference between the two published ZL models, as shown in Figure 15. The different sources of statistical and systematic errors that arose in the data reduction and foreground subtraction are summed quadratically to produce the final EBL errors.

7 Discussions

If the EBL excess arises mainly from high-redshift objects, it obviously conflicts with the gamma-ray transparency of the intergalactic medium (IGM) corresponding to the number of EBL photons for the electron-positron pair creation, as discussed by many authors on previous EBL results (Aharonian et al., 2006; Raue et al., 2009; H.E.S.S. Collaboration, 2013). Therefore, the EBL excess would be attributed to an unaccounted component in the local universe where the gamma-ray absorption length in the line of sight is relatively short, unless new physical mechanisms which reduce the intergalactic gamma-ray absorption (de Angelis et al., 2008; Essey et al., 2010) are particularly required. In fact, an overall EBL spectrum, inferred from our nominal result and previous near-infrared data, is simply analogous to that of the redder galaxies at low redshifts; it exhibits a broad peak around $\sim$1.5 $\mu$m and the Rayleigh-Jeans tail at longer wavelengths. While there is no room to account a large number of missing galaxies in the local universe, unresolvable diffuse light of old stellar populations in galaxy halo and/or intergalactic medium may be source of the EBL excess having redder color. However, it should be noted that faint galaxy wings of resolved galaxies have little contribution to the EBL (Arendt et al., 2010; Donnerstein, 2015).

Our minimum EBL result exhibits a significant excess over the IGL level, but it is marginally consistent with the EBL limit obtained from the gamma-ray observations, and such a relatively small discrepancy could be due to an unaccounted source at high redshift. The gamma-ray limits on the EBL assume a specific spectral shape for the EBL. For a different spectral shape and a different EBL photon energy density as a function of redshift, the gamma-ray limits may be higher in a limited wavelength range as to accord with our minimum EBL result.

Intrahalo light (IHL), produced by tidally stripped stars in low-redshift galaxies, contributes significantly to the near-infrared background fluctuations (Cooray et al., 2012; Mitchell-Wynne et al., 2015), indicated by the excess clustering power on arcminutes scales (Kashlinsky et al., 2005; Thompson et al., 2007; Matsumoto et al., 2011; Zemcov et al., 2014). The IHL contribution to the EBL excess in the LRS bands, from the model fitting to the fluctuations measured with the CIBER Imager instrument, is 7.0 and 11.4 nWm*-2sr-1* at 1.1 and 1.6 $\mu$m, respectively, based on a relative EBL fluctuation amplitude of \delta$$I/I$$\sim20%, which is model dependent and could vary from 10-30% (Zemcov et al., 2014). Our measured EBL excess over the IGL is 25.6$\pm{12.9}$ and 29.7${}^{+15.9}_{-9.9}$ nWm*-2sr-1* at 1.1 and 1.6 $\mu$m for the nominal EBL, respectively, and is 4.4${}^{+3.6}_{-3.5}$ and 15.6${}^{+10.7}_{-4.1}$ nWm*-2sr-1* at 1.1 and 1.6 $\mu$m for the minimum EBL, respectively. Therefore, IHL may account for some of the EBL excess. More precise measurements are needed to compare the color between the intensity and fluctuations of the EBL.

First-light galaxies during the epoch of reionization (EoR), including PopII/PopIII stars (Cooray & Yoshida, 2004; Fernandez & Komatsu, 2006; Inoue et al., 2013) and direct collapse blackholes (Yue et al., 2013b), predict a distinctive spectral shape of the EBL excess, which peaks at $\sim$1.5 $\mu$m analogous to the Lyman break. Such objects with \delta$$I/I$$\sim5% can contribute effectively more to the mean EBL intensity than the fluctuations. However, an EBL intensity from the EoR as high as $\sim$10 nWm*-2sr-1* has difficulties due to the over-enrichment of the intergalactic medium by metals and an overproduction of the X-ray background (Madau & Silk, 2005). For a currently favored star-formation rate density at $z>$ 7 (Bouwens et al., 2012), the high-$z$ contribution to the EBL is at the level of 0.5–1 nWm*-2sr-1*.

Dark stars, hypothetical objects powered by annihilation of accreted weakly interacting massive particles (WIMPs) in the early universe, have been proposed to produce an EBL peaking at $\sim$2 $\mu$m (Maurer et al., 2012). Certain dark star models can be made consistent with the observed EBL excess (because their overall normalization is not constrained), but would be in conflict with the gamma-ray limits.

Even though the minimum EBL result is obtained by conservative subtraction of the ZL foreground, a spatially isotropic ZL component could still be present (Hauser et al., 1998). However, if this is the case, the dust grains must have an extraordinary spectral reflectance (Matsuura et al., 1995). Because the EBL residuals from our data do not have a ZL spectrum, an isotropic ZL component cannot explain our results. A definitive conclusion on the origin of the near-infrared EBL would be most easily achieved by an instrument in deep space (Matsuura, 2002; Matsuura et al., 2013) where the ZL foreground is absent (Weinberg et al., 1974; Hanner et al., 1974; Matsuoka et al., 2011)333The background observations with Pioneer10/11 were not accurate enough to significantly detect the optical EBL., or by precise multi-band fluctuation measurements (Lanz et al., 2014; Shirahata et al., 2016).

8 Summary

We achieved the first direct measurement of the diffuse background spectrum at 0.8–1.7 $\mu$m with the CIBER experiment. The EBL spectrum is derived using model-based foreground assessments of the ISL, DGL and ZL components. The residual absolute surface brightness, 42.7${}^{+11.9}_{-10.6}$ nWm*-2sr-1* at 1.4 $\mu$m, exceeds both the IGL amplitude and upper limits from the gamma-ray transparency of IGM. These measurements are consistent with previous results from the satellite observations at $>$1.25 $\mu$m within the measurement errors, if we compare to those analyses assuming the Kelsall ZL model.

If we assume the Wright ZL model, the derived EBL spectrum is unphysical, becoming less than the IGL at $\lambda<$1.3 $\mu$m and negative at $\lambda<$1 $\mu$m. For both ZL models, the derived EBL is redder than the ZL spectrum at $\lambda<$1.4 $\mu$m. The variation we observe between fields more closely matches the Kelsall ZL model than the Wright ZL model.

The minimum EBL levels at $>$1.25 $\mu$m are estimated using a model-independent method of the ZL subtraction by scaling the ZL template so that the average surface brightness of the EBL for 0.8 $<\lambda<$ 1 $\mu$m equals the IGL. The minimum EBL brightness, e.g., 28.7${}^{+5.1}_{-3.3}$ nWm*-2sr-1* at 1.4 $\mu$m, significantly exceeds the IGL level, and its spectrum exhibits a red color quite different from the blue continuum of the ZL.

These results could be explained by a new foreground component, such as color variations in the ZL, or by a new EBL component due to unaccounted emission from galaxies. Further observational studies in wider wavelength range by future space missions are needed to improve the accuracy of the EBL measurements, due to the uncertainties with the ZL subtraction. For the near term, we plan to extend the range of these absolute spectral measurements to include optical wavelengths.

We acknowledge the support of Grant-in-Aid from Japan Society for the Promotion of Science (KAKENHI 20.34, 18204018, 19540250, 21340047, 21111004, 24111717, 26800112 and 15H05744) and NASA APRA research grant which funded US contributions to the CIBER instrument (NNX07AI54G, NNG05WC18G, NNX07AG43G, NNX07AJ24G, and NNX10AE12G). Initial support was provided by an award to J.B. from the Jet Propulsion Laboratory’s Director’s Research and Development Fund. We thank the dedicated efforts of the sounding rocket staff at NASA Wallops Flight Facility and White Sands Missile Range, and also thank Dr. Allan Smith, Dr. Keith Lykke, and Dr. Steven Brown (NIST) for the laboratory calibration of LRS. A.C. acknowledges support from an NSF CAREER award AST-0645427 and NSF AST-1313319. H.M.L. was supported by NRF grant No.2012R1A4A1028713. S.M. thanks to Dr. Fumihiko Usui (Kobe University, CPS) and Dr. Takafumi Ootsubo (The University of Tokyo) for their help in coordinate calculations and for our discussions on the ZL models. This publication makes use of data products from the 2MASS, which is a joint project of the University of Massachusetts and the Infrared Processing and Analysis Center/California Institute of Technology, funded by the National Aeronautics and Space Administration and the National Science Foundation.

Appendix A Airglow subtraction method

A.1 A Simple Airglow Model with Exponential Functions

As shown in Figure 16, the airglow brightness at a given wavelength can be approximately fit by exponential functions modeled as the sum of time- and altitude-dependent components:

[TABLE]

where the time constant is $\tau$, the scale height is $H$, $t$ is the elapsed time after launch, $h$ is the payload altitude, and$I_{obs}$ and $I_{sky}$ are the observed raw and astronomical brightnesses, respectively. The amplitude parameters, $A_{t}$ and $A_{h}$, give the brightness of the time-dependent and altitude-dependent airglow components at a given wavelength at $t=0$ and $h=0$, respectively. The time- and altitude-dependent components are indistinguishable during the ascent phase. As a result, we fit the altitude-dependent component to the last science field (Bootes-A of the second flight) of the descent phase when the time-dependent component is negligible. We then fit the time-dependent component to the first science field (SWIRE of the second flight) after subtracting an altitude-dependent component. The model curves shown in Figure 16 are obtained by fitting to the data with $\tau$=20 s and $H$=22 km, which are fixed for all wavelengths to the values determined at the peak wavelength, i.e., 1.6 $\mu$m.

Using this result, we can successfully subtract the airglow component to negligible levels compared to the sky brightness. Airglow suppression is not necessary for fields near the apogee of the flight. The SWIRE data showed residual time dependence early in the observations even after the airglow model was subtracted, but the final quarter of the SWIRE data appears stable. The altitude dependence at Bootes-B showed anomalous behavior, which may be due to complex layering in atmospheric OH at lower altitudes, and we did not use these data for astrophysical analysis. For the second flight the final high-quality data set consists of five sky spectra at SWIRE, NEP, Elat10, Elat30, and Bootes-A in 0.8 $<\lambda<$ 1.7 $\mu$m.

A.2 Existence of slow decay airglow component

The above airglow model works only for a fast decay component with a time constant shorter than the integration time of a field. Although we can determine the decay time constant of the time-dependent airglow component from the second flight data alone, there may exist an unaccounted component with a time constant exceeding the integration time, i.e., $\tau>$100 s. We observe such a slow airglow component in the fourth flight data, which is sufficiently long for monitoring very slow signal variations. The time variation of the fourth flight at 1.6 $\mu$m (corresponding to the maximum brightness of the airglow) is shown in Figure 17. The fast decay airglow component is prominent early in the observations. Owing to the high altitude of the fourth flight, an altitude-dependent airglow component is not observed during the flight, so these data are used to constrain slow signal variations.

Since the sky brightness varies from field to field, the slowly decaying airglow component is not easily observed in the raw brightness data shown in Figure 17. In order to extract the time variation of the airglow from the sky brightness, we use the ratio of the sky brightness at 1.6 $\mu$m, the maximum airglow wavelength, to the sky brightness at 0.9 $\mu$m, the minimum airglow wavelength (see Figure 3). The 1.6-$\mu$m to 0.9-$\mu$m brightness ratios are nearly constant when the ZL dominates the total sky brightness. Even assuming an EBL twice as bright as the known IGL brightness, the color variations of the total sky brightness due to the differing contributions of ISL, DGL and EBL to the total sky brightness are estimated to be less than 2% rms, corresponding to a sky brightness of about 6 nWm*-2sr-1* at 1.6 $\mu$m. Therefore, any time variation of the 1.6 $\mu$m to 0.9 $\mu$m ratio exceeding the 2% level must be attributed to non-astrophysical emission.

As shown in Figure 18, the measured 1.6 $\mu$m to 0.9 $\mu$m brightness ratio for the fourth flight exhibits a slow component (in addition to the fast component seen early in the observation), and the ratio of the SWIRE field reaches an asymptotic brightness similar to the predicted value for astronomical sources of 0.62. The slow-decay component can be fit by an exponential function with a time constant of $\tau_{slow}$=144$\pm$19 s, while the fast-decay component can be fit with $\tau_{fast}$=17.5$\pm$0.2 s. Figure 19 shows the relative sky brightness at various wavelengths normalized by the sky brightness at 0.9 $\mu$m. These data are fit by exponential functions with a time constant fixed to $\tau_{slow}$=144 s, setting the amplitudes free. The spectral shape of the determined amplitude parameters with a prominent broad peak around 1.6 $\mu$m is similar to that of the fast decay airglow component, except for a feature around 1 $\mu$m, as shown in Figure 20. The 1 $\mu$m feature is not instrumental effect but real signal, because it is seen in the FOV area of the multiple slits but not in the masked area of the detector array.

A.3 Application of the Airglow Model

This slow decay airglow component may also be present in the second flight data, but unobservable due to the short time scale of the flight. The airglow signal has a large derivative early in the observation when the time-dependent airglow dominates, so are best compared using the 1.6-$\mu$m to 0.9-$\mu$m brightness ratio as shown in Figure 21. After adjusting for a time offset related to the time of opening the shutter door since ignition of the rocket, the shapes of the decays in the two flights are very similar. This result suggests that the second flight may be contaminated by the same airglow seen in the fourth flight.

The fast decay airglow spectra in the two flights are similar, as shown in Figure 20. Because of this, we fit the slow decay airglow spectrum in the second flight using the fixed time constant measured in the fourth flight, but with free amplitude parameters at each wavelength. The derived spectrum of the amplitude parameters is consistent with the fast decay airglow spectrum, as shown in Figures 20 and 22. For simplicity we assume that the slow decay airglow component has same spectral shape as the fast decay airglow component in both flights. From this modeling we can calculate the airglow intensity at a given wavelength and time using our knowledge of the 1.6 $\mu$m to 0.9 $\mu$m brightness ratio and the spectral template.

Our nominal result of the sky brightness measurement is obtained by subtracting the modeled airglow brightness. The fitting error of the exponential function is accounted as statistical error of the sky brightness measurement. For the second flight the full brightness difference with and without subtraction of the slow-decay component is conservatively taken as the systematic error associated with airglow subtraction.