Space Telescope and Optical Reverberation Mapping Project. V. Optical Spectroscopic Campaign and Emission-Line Analysis for NGC 5548

TL;DR

This study presents a detailed optical spectroscopic monitoring of NGC 5548, revealing complex emission-line responses, a larger broad-line region radius than previously inferred, and significant variability in emission-line lags over six months.

Contribution

It provides new insights into the structure and dynamics of the broad-line region through high-cadence spectroscopic data and highlights variability in emission-line responses over time.

Findings

Hβ and He II lags measured at 4.17 and 0.79 days

Broad-line region radius is ~50% larger than optical-only estimates

Emission-line responses changed significantly during the campaign

Abstract

We present the results of an optical spectroscopic monitoring program targeting NGC 5548 as part of a larger multi-wavelength reverberation mapping campaign. The campaign spanned six months and achieved an almost daily cadence with observations from five ground-based telescopes. The H$\beta$ and He II $\lambda$4686 broad emission-line light curves lag that of the 5100 $\AA$ optical continuum by $4.17^{+0.36}_{-0.36}$ days and $0.79^{+0.35}_{-0.34}$ days, respectively. The H$\beta$ lag relative to the 1158 $\AA$ ultraviolet continuum light curve measured by the Hubble Space Telescope is roughly $\sim$50% longer than that measured against the optical continuum, and the lag difference is consistent with the observed lag between the optical and ultraviolet continua. This suggests that the characteristic radius of the broad-line region is $\sim$50% larger than the value inferred from optical data alone. We also measured velocity-resolved emission-line lags for H$\beta$ and found a complex velocity-lag structure with shorter lags in the line wings, indicative of a broad-line region dominated by Keplerian motion. The responses of both the H$\beta$ and He II $\lambda$4686 emission lines to the driving continuum changed significantly halfway through the campaign, a phenomenon also observed for C IV, Ly $\alpha$, He II(+O III]), and Si IV(+O IV]) during the same monitoring period. Finally, given the optical luminosity of NGC 5548 during our campaign, the measured H$\beta$ lag is a factor of five shorter than the expected value implied by the $R_\mathrm{BLR} - L_\mathrm{AGN}$ relation based on the past behavior of NGC 5548.

Tables

Table 1. Table 1 Instrument characteristics and data-reduction parameters for all telescopes

Telescope	Instrument	Number of	Median	Wavelength	Wavelength	Pixel	Median	[O III]
		Epochs	Seeing	Dispersion	Coverage	Scale	SNR	$F_{\mathrm{var}}$
			(″)	(Å pixel-1)	(Å)	(″ pixel-1)		(%)
MDM	Boller & Chivens CCD Spectrograph	143	1.7	1.25	4225$-$5775	0.75	118	0.62
Lick	Kast Double Spectrograph	35	1.5	1.02	3460$-$5500	0.43	194	0.32
Asiago	Boller & Chivens CCD Spectrograph	21	4.0	1.00	3250$-$7920	1.00	160	0.27
APO	Dual Imaging Spectrograph	13	1.4	1.00	4180$-$5400	0.41	160	0.28
WIRO	WIRO Long Slit Spectrograph	6	2.1	0.74	5599$-$4399	0.52	217	0.47

Table 2. Table 2 Flux measurements for continuum and emission lines.

HJD$-$2,450,000	Telescope	$F_{5100}$	$F_{H\beta }$	$F_{H\beta ,\mathrm{SD}}$	$F_{\mathrm{He}\mathrm{II},\mathrm{SD}}$
6663.00	MDM	10.766 $\pm$ 0.075	726.012 $\pm$ 4.985	710.187 $\pm$ 4.897	21.720 $\pm$ 2.606
6663.65	Asiago	10.921 $\pm$ 0.040	741.771 $\pm$ 3.586	—	—
6664.03	MDM	11.154 $\pm$ 0.075	732.511 $\pm$ 5.156	715.057 $\pm$ 5.061	28.154 $\pm$ 3.222
6665.02	MDM	10.788 $\pm$ 0.075	724.537 $\pm$ 4.946	709.473 $\pm$ 4.860	25.135 $\pm$ 2.451
6667.02	MDM	10.872 $\pm$ 0.076	735.001 $\pm$ 5.393	711.347 $\pm$ 5.288	37.008 $\pm$ 3.964

Table 3. Table 3 Statistics for HST and MDM Light Curves

Emission Component	Epochs	Mean Flux	rms Flux	$F_{\mathrm{var}}$	$R_{\max }$
$F_{\lambda }$ (1158 Å)	171	043.48 $\pm$ 0.86	11.14 $\pm$ 1.21	0.255	4.07 $\pm$ 0.18
$F_{\lambda }$ (5100 Å)	133	011.96 $\pm$ 0.07	00.80 $\pm$ 0.09	0.066	1.33 $\pm$ 0.01
H$\beta$	133	738.49 $\pm$ 2.40	28.29 $\pm$ 3.38	0.038	1.22 $\pm$ 0.01
He II $\lambda$4686	133	078.71 $\pm$ 2.95	35.14 $\pm$ 4.17	0.444	7.48 $\pm$ 0.91
$F_{\lambda }$ (1158 Å, T1)	051	035.85 $\pm$ 1.79	12.61 $\pm$ 2.54	0.351	3.31 $\pm$ 0.15
$F_{\lambda }$ (1158 Å, T2)	120	046.72 $\pm$ 0.80	08.66 $\pm$ 1.13	0.184	2.36 $\pm$ 0.08
$F_{\lambda }$ (5100 Å, T1)	067	011.31 $\pm$ 0.08	00.67 $\pm$ 0.12	0.059	1.27 $\pm$ 0.01
$F_{\lambda }$ (5100 Å, T2)	067	012.51 $\pm$ 0.04	00.37 $\pm$ 0.06	0.029	1.13 $\pm$ 0.01
H$\beta$ (T1)	067	725.10 $\pm$ 3.80	30.47 $\pm$ 5.36	0.041	1.20 $\pm$ 0.01
H$\beta$ (T2)	067	750.14 $\pm$ 2.38	20.11 $\pm$ 3.33	0.026	1.13 $\pm$ 0.01
He II $\lambda$4686 (T1)	067	065.84 $\pm$ 4.16	33.61 $\pm$ 5.89	0.507	5.95 $\pm$ 0.73
He II $\lambda$4686 (T2)	067	089.90 $\pm$ 3.75	32.71 $\pm$ 5.34	0.362	4.07 $\pm$ 0.37

Table 4. Table 4 Rest-Frame Emission-Line Lags

Light Curves	${\tau }_{\mathrm{peak}}$	${\tau }_{\mathrm{cen}}$	${\tau }_{\mathrm{cen},\mathrm{T1}}$	${\tau }_{\mathrm{cen},\mathrm{T2}}$	${\tau }_{\text{JAVELIN}}$	${\tau }_{\text{JAVELIN},\mathrm{T1}}$	${\tau }_{\text{JAVELIN},\mathrm{T2}}$
H$\beta$ vs. $F_{\lambda }$(1158 Å)	${6.14}_{-0.98}^{+0.74}$	${6.23}_{-0.44}^{+0.39}$	${7.62}_{-0.49}^{+0.49}$	${5.99}_{-0.75}^{+0.71}$	${6.56}_{-0.49}^{+0.48}$	${6.91}_{-0.63}^{+0.64}$	${7.42}_{-1.07}^{+0.97}$
H$\beta$ vs. $F_{\lambda }$(1367 Å)	${5.90}_{-0.74}^{+0.25}$	${5.89}_{-0.37}^{+0.37}$	${7.24}_{-0.48}^{+0.49}$	${5.99}_{-0.82}^{+0.76}$	${6.12}_{-0.47}^{+0.46}$	${6.52}_{-0.57}^{+0.60}$	${7.11}_{-1.06}^{+1.03}$
H$\beta$ vs. $F_{\lambda }$(5100 Å)	${4.42}_{-0.25}^{+0.98}$	${4.17}_{-0.36}^{+0.36}$	${4.99}_{-0.47}^{+0.40}$	${3.10}_{-0.80}^{+0.77}$	${3.84}_{-0.59}^{+0.57}$	${5.15}_{-0.69}^{+0.68}$	${4.78}_{-1.17}^{+1.13}$
H$\beta$ vs. V band	${3.93}_{-0.98}^{+0.98}$	${3.79}_{-0.34}^{+0.37}$	${3.82}_{-0.47}^{+0.57}$	${4.13}_{-0.58}^{+0.55}$	${3.54}_{-0.46}^{+0.45}$	${4.89}_{-0.71}^{+0.66}$	${4.05}_{-0.78}^{+0.93}$
He II vs. $F_{\lambda }$(1158 Å)	${2.46}_{-0.25}^{+0.49}$	${2.69}_{-0.25}^{+0.24}$	${3.71}_{-0.38}^{+0.39}$	${3.19}_{-0.35}^{+0.36}$	${2.65}_{-0.27}^{+0.27}$	${3.27}_{-0.35}^{+0.35}$	${2.99}_{-0.26}^{+0.25}$
He II vs. $F_{\lambda }$(1367 Å)	${2.21}_{-0.25}^{+0.25}$	${2.45}_{-0.24}^{+0.25}$	${3.43}_{-0.43}^{+0.36}$	${3.16}_{-0.33}^{+0.29}$	${2.41}_{-0.26}^{+0.25}$	${3.04}_{-0.36}^{+0.35}$	${2.79}_{-0.25}^{+0.25}$
He II vs. $F_{\lambda }$(5100 Å)	${0.49}_{-0.25}^{+0.25}$	${0.79}_{-0.34}^{+0.35}$	${1.21}_{-0.36}^{+0.28}$	${0.85}_{-0.36}^{+0.36}$	${0.16}_{-0.37}^{+0.37}$	${1.13}_{-0.48}^{+0.51}$	${0.85}_{-0.35}^{+0.38}$
He II vs. V band	${0.49}_{-0.49}^{+0.25}$	${0.50}_{-0.26}^{+0.34}$	${0.40}_{-0.39}^{+0.43}$	${1.46}_{-0.27}^{+0.34}$	${0.44}_{-0.23}^{+0.23}$	${0.63}_{-0.36}^{+0.37}$	${0.91}_{-0.21}^{+0.23}$
$F_{\lambda }$(5100 Å) vs. $F_{\lambda }$(1158 Å)	${1.97}_{-0.49}^{+0.25}$	${2.23}_{-0.26}^{+0.31}$	${2.55}_{-0.33}^{+0.28}$	${2.77}_{-0.45}^{+0.41}$
H$\beta$ ${}_{\mathrm{full},\mathrm{MDM}}$ vs. $F_{\lambda }$(1158 Å)	${5.90}_{-0.49}^{+0.74}$	${6.26}_{-0.37}^{+0.38}$
H$\beta$ ${}_{\mathrm{full},\mathrm{all}\mathrm{sites}}$ vs. $F_{\lambda }$(1158 Å)	${5.65}_{-0.25}^{+0.74}$	${6.49}_{-0.37}^{+0.37}$

Table 5. Table 5 Rest-Frame H β 𝛽 \beta Velocity-Resolved Lags

Wavelength	$v_r$	${\tau }_{\mathrm{Full}}$	$r_{\max ,\mathrm{Full}}$	${\tau }_{\mathrm{T1}}$	$r_{\max ,\mathrm{T1}}$	${\tau }_{\mathrm{T2}}$	$r_{\max ,\mathrm{T2}}$
(Å)	(km s-1)	(days)		(days)		(days)
4743.82$-$4751.92	7000	0${1.74}_{-0.60}^{+0.59}$	0.61	0${2.05}_{-0.94}^{+1.04}$	0.79	${2.13}_{-0.66}^{+0.59}$	0.42
4751.92$-$4760.03	6500	0${2.34}_{-0.49}^{+0.49}$	0.72	0${2.57}_{-0.71}^{+0.60}$	0.88	${3.08}_{-0.77}^{+0.98}$	0.52
4760.03$-$4768.13	6000	0${2.94}_{-0.48}^{+0.45}$	0.73	0${3.41}_{-0.70}^{+0.63}$	0.88	${3.48}_{-0.66}^{+0.71}$	0.62
4768.13$-$4776.23	5500	0${3.18}_{-0.46}^{+0.38}$	0.75	0${3.25}_{-0.54}^{+0.57}$	0.83	${3.46}_{-0.62}^{+0.72}$	0.62
4776.23$-$4784.33	5000	0${3.92}_{-0.47}^{+0.50}$	0.73	0${3.12}_{-0.53}^{+0.57}$	0.85	${4.95}_{-0.94}^{+1.04}$	0.56
4784.33$-$4792.43	4500	0${4.69}_{-0.49}^{+0.56}$	0.73	0${3.94}_{-0.49}^{+0.63}$	0.85	${5.88}_{-0.70}^{+0.63}$	0.65
4792.43$-$4800.53	4000	0${6.19}_{-0.68}^{+0.63}$	0.63	0${7.00}_{-0.85}^{+0.74}$	0.83	${6.26}_{-1.29}^{+1.34}$	0.53
4800.53$-$4808.63	3500	0${7.26}_{-0.61}^{+0.62}$	0.68	0${8.45}_{-0.87}^{+0.94}$	0.79	${6.93}_{-1.57}^{+2.34}$	0.51
4808.63$-$4816.73	3000	0${8.69}_{-0.67}^{+0.64}$	0.63	${10.90}_{-0.89}^{+0.89}$	0.81	${6.40}_{-1.48}^{+2.71}$	0.42
4816.73$-$4824.84	2500	0${8.70}_{-1.13}^{+1.29}$	0.53	${10.54}_{-0.96}^{+1.55}$	0.77	${4.92}_{-0.86}^{+1.65}$	0.42
4824.84$-$4832.94	2000	0${7.63}_{-1.01}^{+0.87}$	0.58	${10.46}_{-1.07}^{+1.24}$	0.77	${4.31}_{-0.50}^{+0.62}$	0.66
4832.94$-$4841.05	1500	0${5.91}_{-0.53}^{+0.58}$	0.68	0${8.71}_{-0.87}^{+0.78}$	0.80	${5.36}_{-0.74}^{+0.62}$	0.70
4841.05$-$4849.15	1000	0${4.65}_{-0.45}^{+0.48}$	0.71	0${5.79}_{-0.53}^{+0.63}$	0.89	${5.36}_{-0.60}^{+0.55}$	0.74
4849.15$-$4857.25	0500	0${5.27}_{-0.49}^{+0.38}$	0.74	0${5.80}_{-0.65}^{+0.61}$	0.90	${6.11}_{-0.67}^{+0.51}$	0.76
4857.25$-$4865.35	0000	0${6.01}_{-0.49}^{+0.37}$	0.79	0${6.50}_{-0.60}^{+0.52}$	0.90	${6.59}_{-0.76}^{+0.57}$	0.71
4865.35$-$4873.45	0500	0${6.40}_{-0.39}^{+0.47}$	0.75	0${8.06}_{-0.70}^{+0.67}$	0.84	${6.43}_{-0.71}^{+0.64}$	0.67
4873.45$-$4881.56	1000	0${6.99}_{-0.46}^{+0.38}$	0.77	0${8.86}_{-0.39}^{+0.55}$	0.90	${6.87}_{-0.59}^{+0.53}$	0.74
4881.56$-$4889.66	1500	0${7.77}_{-0.37}^{+0.44}$	0.80	0${9.80}_{-0.46}^{+0.41}$	0.92	${7.05}_{-0.52}^{+0.57}$	0.73
4889.66$-$4897.76	2000	0${8.59}_{-0.39}^{+0.40}$	0.81	${10.65}_{-0.55}^{+0.53}$	0.90	${6.38}_{-0.68}^{+0.75}$	0.65
4897.76$-$4905.86	2500	${10.18}_{-0.67}^{+0.66}$	0.74	${11.13}_{-0.68}^{+0.67}$	0.87	${6.47}_{-2.38}^{+2.01}$	0.27
4905.86$-$4913.96	3000	${10.19}_{-0.92}^{+0.98}$	0.67	${11.18}_{-0.61}^{+0.71}$	0.86	${2.58}_{-0.98}^{+1.21}$	0.21
4913.96$-$4922.06	3500	0${8.23}_{-0.96}^{+1.03}$	0.66	0${9.82}_{-0.96}^{+0.86}$	0.81	${3.78}_{-1.46}^{+2.74}$	0.27
4922.06$-$4930.16	4000	0${7.25}_{-0.91}^{+0.87}$	0.65	0${8.80}_{-0.70}^{+0.71}$	0.87	${4.06}_{-0.97}^{+1.57}$	0.50
4930.16$-$4938.27	4500	0${6.01}_{-0.59}^{+0.55}$	0.62	0${7.01}_{-0.72}^{+0.81}$	0.88	${5.52}_{-0.70}^{+0.64}$	0.62
4938.27$-$4946.38	5000	0${5.04}_{-0.53}^{+0.59}$	0.47	0${6.63}_{-0.56}^{+0.61}$	0.88	${4.85}_{-0.68}^{+0.69}$	0.58
4946.38$-$4954.48	5500	0${4.00}_{-0.57}^{+0.56}$	0.40	0${5.22}_{-0.83}^{+0.84}$	0.80	${4.55}_{-0.64}^{+0.72}$	0.59
4954.48$-$4962.58	6000	0${3.41}_{-0.60}^{+0.63}$	0.43	0${3.16}_{-1.30}^{+1.05}$	0.83	${4.89}_{-0.70}^{+0.65}$	0.62
4962.58$-$4970.68	6500	0${3.20}_{-0.53}^{+0.61}$	0.53	0${3.46}_{-0.72}^{+0.74}$	0.85	${4.30}_{-0.96}^{+0.87}$	0.53

Table 6. Table 6 Responsivity, Slope of the Baldwin Relation, and Flux Deficit During the Anomaly for All Emission Lines

Line ID	${\eta }_{\mathrm{eff}}$	$\beta$	$F_{\mathrm{lost}}$
Ly$\alpha$	$0.30\pm 0.01$	$-0.73\pm 0.02$	9%
Si IV$+$O IV]	$0.45\pm 0.01$	$-0.58\pm 0.03$	23%
C IV	$0.25\pm 0.01$	$-0.75\pm 0.01$	18%
He II$+$O III]	$0.58\pm 0.04$	$-0.48\pm 0.04$	21%
H$\beta$	$0.15\pm 0.01$	$-0.85\pm 0.02$	6%

Table 7. Table 7 H β 𝛽 \beta Responsivity and AGN Optical Continuum Flux Density for NGC 5548 Over 25 Years

Year	${\eta }_{\mathrm{eff}}$	Mean $F_{\mathrm{cont}}$	rms $F_{\mathrm{cont}}$
1989	0.56 $\pm$ 0.04	06.54	1.27
1990	0.84 $\pm$ 0.03	03.79	0.91
1991	0.95 $\pm$ 0.09	06.06	0.92
1992	0.94 $\pm$ 0.05	03.34	1.17
1993	0.43 $\pm$ 0.04	05.69	0.87
1994	0.74 $\pm$ 0.04	06.40	1.11
1995	0.68 $\pm$ 0.04	08.71	1.01
1996	0.54 $\pm$ 0.03	07.07	1.52
1997	0.80 $\pm$ 0.07	04.73	0.89
1998	0.51 $\pm$ 0.02	10.05	1.44
1999	0.41 $\pm$ 0.04	08.48	1.82
2000	0.65 $\pm$ 0.06	03.59	1.20
2001	1.00 $\pm$ 0.12	03.65	0.86
2014	0.59 $\pm$ 0.01	07.44	0.50

Table 8. Table 8 H β 𝛽 \beta Rest-Frame Line Widths

Segment	Spectrum	${\sigma }_{\mathrm{line}}$	FWHM
		(km s-1)	(km s-1)
Full	rms	4278 $\pm$ 671	10161 $\pm$ 587
T1	rms	4155 $\pm$ 513	10861 $\pm$ 739
T2	rms	4856 $\pm$ 731	09103 $\pm$ 1279
Full	Mean	3691 $\pm$ 162	9496 $\pm$ 418
T1	Mean	3983 $\pm$ 150	9612 $\pm$ 427
T2	Mean	3939 $\pm$ 177	9380 $\pm$ 158

Table 9. Table 9 H β 𝛽 \beta Line Measurements, M BH subscript 𝑀 BH M_{\mathrm{BH}} Estimates, and Continuum Flux Densities

Segment	${\sigma }_{\mathrm{line}}$	${\tau }_{H\beta -\mathrm{opt}}$	Virial Product	$M_{\mathrm{BH}}$	$F_{5100,\mathrm{total}}$	$F_{5100,\mathrm{AGN}}$
	(km s-1)	(days)	(${10}^7$ $M_{\odot }$)	(${10}^7$ $M_{\odot }$)	(${10}^{-15}\mathrm{erg}{s}^{-1}{\mathrm{cm}}^{-2}{\AA }^{-1}$)
Full	4278 $\pm$ 671	${4.17}_{-0.36}^{+0.36}$	${1.49}_{-0.49}^{+0.49}$	${6.66}_{-2.17}^{+2.17}$	11.96 $\pm$ 0.07	7.44 $\pm$ 0.50
T1	4155 $\pm$ 513	${4.99}_{-0.47}^{+0.40}$	${1.68}_{-0.45}^{+0.44}$	${7.53}_{-1.99}^{+1.96}$	11.31 $\pm$ 0.08	6.79 $\pm$ 0.46
T2	4856 $\pm$ 731	${3.10}_{-0.80}^{+0.77}$	${1.43}_{-0.56}^{+0.55}$	${6.38}_{-2.53}^{+2.49}$	12.51 $\pm$ 0.04	7.99 $\pm$ 0.45

Equations

$$L(v_r,t)={\int }_0^{\infty }\Psi (v_r,\tau )C(t-\tau )d\tau ,$$

$$M_{\mathrm{BH}}=f\frac{c\tau \Delta V^2}{G}.$$

$$F_{\mathrm{var}}=\frac{\sqrt{{\sigma }^2-⟨{\delta }^2⟩}}{⟨f⟩},$$

$${\mathrm{e}\text{-}\mathrm{rms}}_{\lambda }=\sqrt{\frac{1}{N-1}\sum _{i=1}^N[(F_{\lambda ,i}-⟨F_{\lambda }⟩{)}^2-{\delta }_{\lambda ,i}^2]},$$

$$F_{\mathrm{H}\beta ,\mathrm{MDM}}=\varphi F_{\mathrm{H}\beta ,\mathrm{t}}$$

$$F_{5100,\mathrm{MDM}}=\varphi F_{5100,\mathrm{t}}+G.$$

$$\mathrm{EW}=\frac{F_{\mathrm{line}}}{F_{\mathrm{cont}}},$$

$$\log F_{\mathrm{line}}=A+{\eta }_{\mathrm{eff}}[\log F_{\mathrm{cont}}].$$

$$\log \mathrm{E}{\mathrm{W}}_{\mathrm{line}}=B+\beta [\log F_{\mathrm{cont}}],$$

$$\sigma^{2}_{\mathrm{line}}=\Bigg{(}\frac{c}{\lambda_{0}}\Bigg{)}^{2}\Bigg{(}\frac{\sum\lambda^{2}_{i}S_{i}}{\sum S_{i}}-\lambda^{2}_{0}\Bigg{)},$$

$$\frac{{\tau }_{\mathrm{H}\beta -\mathrm{opt}}}{{\tau }_{\mathrm{H}\beta -\mathrm{UV}}}=1-\frac{{\tau }_{\mathrm{opt}-\mathrm{UV}}}{{\tau }_{\mathrm{H}\beta -\mathrm{UV}}}.$$

$$R_{\lambda }\propto M_{\mathrm{BH}}^{2/3}{l}^{1/3}{\lambda }^{4/3},$$

$$R_{\mathrm{BLR}}\propto M_{\mathrm{BH}}^{1/2}{l}^{1/2}.$$

$$\frac{{\tau }_{\mathrm{opt}-\mathrm{UV}}}{{\tau }_{\mathrm{H}\beta -\mathrm{UV}}}\propto M_{\mathrm{BH}}^{1/6}{l}^{-1/6}{\lambda }^{4/3}.$$

Full text

77footnotetext: NSF Postdoctoral Research Fellow3131footnotetext: NASA Postdoctoral Program Fellow4343footnotetext: Einstein Fellow4444footnotetext: Sagan Fellow4848footnotetext: Carnegie-Princeton Fellow, Hubble Fellow5252footnotetext: Packard Fellow

1. Introduction

Broad emission lines are among the most striking features of quasars and active galactic nuclei (AGNs). These Doppler-broadened lines are emitted by gas occupying the broad-line region (BLR), which is located within several light-days to light-months of the central supermassive black hole (SMBH; e.g., Antonucci & Cohen, 1983; Clavel et al., 1991; Peterson et al., 1998, 2004; Bentz et al., 2009b; Grier et al., 2013). The geometry and kinematics of the BLR play a significant role in AGN research because these properties can be used to infer the mass of the central black hole (e.g., Gaskell & Sparke, 1986; Clavel et al., 1991; Kaspi et al., 2000; Denney et al., 2006, 2010; Pancoast et al., 2014). Additionally, it is possible that infalling BLR gas may fuel SMBH accretion (e.g., Peterson, 2006; Gaskell & Goosmann, 2016) and outflowing gas may be part of disk winds that carry away angular momentum from the disk and provide energy and momentum feedback to the host galaxy (e.g., Emmering et al., 1992; Murray & Chiang, 1997; Kollatschny, 2003; Leighly & Moore, 2004). Understanding the dynamical state and physical conditions of gas in the BLR is of key importance in completing our understanding of the AGN phenomenon.

Owing to its small angular size, the BLR is currently impossible to resolve spatially even for the closest AGNs. An alternative method to study this region is to resolve it in the time domain using reverberation mapping (RM), a technique that leverages the variable nature of quasars and Seyferts (Blandford & McKee, 1982; Peterson, 1993, 2014). AGNs exhibit stochastic flux variations, possibly because of inhomogeneous accretion and thermal fluctuations in the accretion disk (Czerny et al., 1999; Collier & Peterson, 2001; Czerny et al., 2003; Kelly et al., 2009; Kozłowski et al., 2010; MacLeod et al., 2010). Photons from the central engine ionize the BLR gas, which then echoes continuum-flux variations with a light-travel time lag, $\tau$. The emission-line flux $L(v_{r},t)$ at time $t$ and line-of-sight velocity $v_{\mathrm{r}}$ is related to the ionizing continuum by

[TABLE]

where $C(t-\tau)$ is the continuum emission at an earlier time $t-\tau$, and $\Psi(v,\tau)$ is the transfer function that maps the continuum light curve to the time-variable line profile (Blandford & McKee, 1982).

The transfer function—also known as the velocity-delay map—encodes important information about the BLR’s geometry and kinematics. There has been tremendous effort by many groups to recover velocity-delay maps (Rosenblatt & Malkan, 1990; Horne et al., 1991; Krolik et al., 1991; Ulrich & Horne, 1996; Bentz et al., 2010a; Pancoast et al., 2011; Li et al., 2013; Grier et al., 2013; Pancoast et al., 2014; Skielboe et al., 2015) and velocity-resolved line lags (e.g., Kollatschny, 2003; Bentz et al., 2009b; Denney et al., 2010; Barth et al., 2011; Du et al., 2016a). In order to obtain $\Psi(v,\tau)$, RM campaigns must have a combination of high cadence, long duration, high photometric precision, and high signal-to-noise ratios (SNR), which is often not achievable by ground-based programs. More typically, RM campaigns are able to only measure the mean emission-line lag $\tau$, which represents the response-weighted mean light-travel time from the ionizing continuum to the BLR.

Assuming that the broad-line width is a result of the virialized motion of gas within the black hole’s potential well, the emission-line lag and gas velocity dispersion inferred from the line width ($\Delta V$) can be used to infer the black hole (BH) mass using

[TABLE]

Here, $c\tau=R_{\mathrm{BLR}}$ is the characteristic radius of the BLR, and f is a dimensionless calibration factor of order unity that accounts for the unknown BLR geometry and kinematics. Ground-based RM campaigns have produced BH mass measurements for $\sim$60 local AGNs to date (see Bentz & Katz 2015 for references and a recent compilation). RM is also starting to be used for objects at cosmological distances (Kaspi et al., 2007; King et al., 2015; Shen et al., 2016) with the aims of studying the UV continuum and emission lines and calibrating BH masses at high redshifts.

The ionizing continuum is emitted at wavelengths $\textless$ 912 Å and is generally unobservable due to the Lyman limit of the host galaxy. Given this limitation, the far-ultraviolet (UV) continuum at $\lambda\approx 1100$–1500 Å should be used to derive emission-line lags because it is close in wavelength to the ionizing continuum and should therefore serve as an accurate proxy. However, wavelengths shorter than $\sim$3200 Å are inaccessible from the ground, so the rest-frame optical continuum is often used as a proxy for the ionizing source in low-redshift AGNs. Although the far-UV and optical continua have been shown to vary almost simultaneously in some cases (e.g., Clavel et al., 1991; Reichert et al., 1994; Korista et al., 1995; Wanders et al., 1997), more recent high-cadence studies have found that the optical continuum can lag the UV continuum by up to a few days (Collier et al., 1998; Sergeev et al., 2005; McHardy et al., 2014; Shappee et al., 2014; Edelson et al., 2015; Fausnaugh et al., 2016). This can significantly affect the measured broad-line lag if the BLR has a characteristic radius on the order of light days. The variable optical continuum has also been shown to have smoother features and smaller amplitudes than its UV counterpart (e.g., Peterson et al., 1991; Dietrich et al., 1993; Stirpe et al., 1994; Santos-Lleó et al., 1997; Dietrich et al., 1998; Shappee et al., 2014; Fausnaugh et al., 2016). These differences between the UV and optical continua suggest that the optical continuum is not fully interchangeable with the ionizing source for determining reverberation lags.

Furthermore, a long-standing assumption in RM is that the source of the ionizing photons in a typical Seyfert galaxy is physically much smaller than the BLR (about a factor of 100; e.g., Peterson, 1993; Peterson & Horne, 2004). This assumption implies that the disk size can be neglected when determining $R_{\mathrm{BLR}}$ from RM data. However, Fausnaugh et al. (2016) have shown that the optically emitting portion of the accretion disk has a lag similar to that of the inner portion of the BLR. If we assume a model in which the measured lags are purely dependent on the radial distance from the ionizing source, then the emission-line lags measured using the optical continuum may significantly underestimate the BLR characteristic radius. Since most RM campaigns use only optical data, it is imperative that we understand the systematic effects of using the optical rather than the UV continuum in RM studies and the relevant implications for BH mass estimates.

To this end, we present the results of a six-month ground-based RM program monitoring the galaxy NGC 5548 (redshift $z=0.0172$). This paper is the fifth in a series describing results from the AGN Space Telescope and Optical Reverberation Mapping (AGN STORM) campaign, the most intensive multi-wavelength AGN monitoring program to date. The campaign is centered around 171 epochs of daily cadence observations using the Cosmic Origins Spectrograph on the Hubble Space Telescope (HST). Concurrent with the HST program were four months of Swift observations and six months of ground-based photometric and spectroscopic observations. First results of the HST, Swift, and ground-based photometry programs were presented by De Rosa et al. (2015), Edelson et al. (2015), and Fausnaugh et al. (2016) (Papers I–III, respectively). Goad et al. (2016) (Paper IV) explore the anomalous behavior of the UV continuum and broad emission-line light curves observed during a portion of this campaign. This paper focuses on the ground-based spectroscopic data and emission-line analysis.

NGC 5548 is one of the best-studied Seyfert galaxies and has been the subject of many past RM programs. Most notably, it was the target of a 13-year campaign carried out by the AGN Watch consortium (Peterson et al., 2002, and references therein), which was initially designed to support UV monitoring of NGC 5548 carried out by the International Ultraviolet Explorer (IUE; Clavel et al., 1991). Individual years of this campaign achieved median sampling cadences of 1–3 days for spectroscopic observations. Subsequently, NGC 5548 was monitored in programs described by Bentz et al. (2007), Denney et al. (2009), Bentz et al. (2009b), and De Rosa et al. (AGN12, results in preparation) with campaign durations of 40 days, 135 days, 64 days, and 120 days (respectively), and each with a median sampling cadence of $\sim 1$ day. A more recent RM program described by Lu et al. (2016) monitored this AGN for 180 days with a median spectroscopic sampling of $~{}\sim 3$ days. The 2014 AGN STORM campaign’s combination of daily cadence, six-month duration, and multi-wavelength coverage makes it the most intensive RM campaign ever conducted.

There are two primary goals of the present work. The first is to compare the H$\beta$ emission-line lag measured against simultaneously observed far-UV and optical continua in order to understand the effects of substituting the optical continuum for the ionizing continuum in reverberation measurements. The second goal is to examine in detail the responses of the optical emission lines to continuum variations and compare them to those of the UV lines, which will provide a more complete picture of the structure and kinematics of the BLR than previous studies that used only optical data.

We describe the spectroscopic observations and reductions in Section 2. Section 3 details our procedures for flux and light-curve measurements. In Section 4, we present our analysis of emission-line lags, line responses, line profiles, and BH mass measurements. We discuss the implications of our results and compare our measurements to those from previous campaigns in Section 5. Section 6 summarizes our findings. We quote wavelengths in the rest frame of NGC 5548 unless otherwise stated.

2. Observations and Data Reduction

Spectroscopic data were obtained from five telescopes: the McGraw-Hill 1.3-m telescope at the MDM Observatory, the Shane 3-m telescope at the Lick Observatory, the 1.22-m Galileo telescope at the Asiago Astrophysical Observatory, the 3.5-m telescope at Apache Point Observatory (APO), and the 2.3-m telescope at the Wyoming Infrared Observatory (WIRO). Observations at MDM were carried out with a slit width of 5″ oriented in the north-south direction, and spectra at the other telescopes were taken with a 5″-wide slit oriented at the parallactic angle (Filippenko, 1982). The optical spectroscopic monitoring began on 2014 January 4 (UT dates are used throughout this paper) and continued through 2014 July 6 with approximately daily cadence.

Table 1 lists the properties of the telescopes and instruments used to obtain spectroscopic data, and Figure 1 shows the mean spectrum constructed using data from Asiago, which obtained the only spectra that cover the full optical wavelength range. MDM contributed the largest number of spectra with 143 epochs. The 35 epochs of Lick spectra were obtained by several groups of observers who used slightly different setups and calibrations. The Kast spectrograph (Miller & Stone, 1993) at Lick Observatory has red-side and blue-side cameras, but since the red-side setup was very different for each group, we present only the blue-side data here. Asiago, APO, and WIRO contributed 21, 13, and 6 epochs of spectra, respectively. Our analysis focuses primarily on the MDM dataset for homogeneity.

Data-reduction procedures included bias subtraction, flat fielding, and cosmic ray removal using the L.A. Cosmic routine (van Dokkum, 2001). The one-dimensional spectra were extracted from a 15″-wide region centered on the AGN and with consistent background sky apertures for all observations. We used optimally-weighted extractions for the stellar spectra (Horne, 1986) but unweighted extractions for the AGN spectra. This is because the optimal extraction method requires the spatial profile of the target to be a smooth function of wavelength, and tends to truncate the peaks of strong emission lines such as [O III] that have different spatial extents from the surrounding continuum.

The data were wavelength calibrated using night sky lines and flux calibrated using standard stars. Our most frequently used flux standard stars were Feige 34, BD 332642, and HZ 44. For nights when multiple exposures were taken, we aligned the flux-calibrated one-dimensional spectra by applying small wavelength shifts to each spectrum before combining them. We do not expect significant differential atmospheric refraction (Filippenko, 1982) because of the large slit width use for our observations.

For the MDM data, the first 133 epochs were flux calibrated using Feige 34, while the last 10 epochs, taken from 2014 June 20 to 2014 June 30, were flux calibrated with BD 332642. This caused spurious changes in the shape of some emission-line features, so we use only the first 133 MDM epochs for our present analysis.

2.1. Spectral Flux Calibrations

To place the instrumental fluxes on an absolute flux scale, we measured the narrow [O III] $\lambda$5007 line flux from spectra taken under photometric conditions and scaled all other nightly spectra to have the same [O III] flux. There were 21 epochs identified as having been observed under photometric conditions by the MDM observers. We determined the flux of the [O III] line ($\lambda_{\mathrm{observed}}=5093$ Å) by first subtracting a linear fit to continuum windows on either side of the line, then integrating over a fixed wavelength range. We used the rest-frame wavelength ranges 4976.5$-$4948.0 Å and 5027.7$-$5031.6 Å to fit the continuum and integrated over the range 4980.5$-$5026.7 Å for the line flux. The 2$\sigma$ outliers from this set of [O III] flux measurements were discarded, the mean was recomputed, and this process was repeated until there were no more 2$\sigma$ outliers, which resulted in a total of 16 final photometric spectra. The mean spectrum of these 16 epochs has an [O III] $\lambda$5007 line flux of $(5.01\pm 0.11)\times 10^{-13}$ erg s*-1* cm*-2*, which represents our best estimate of the true [O III] flux for NGC 5548 during this campaign and is not expected to vary over a six-month period. For comparison, Peterson et al. (2013) found the [O III] flux in NGC 5548 to be $(4.77\pm 0.14)\times 10^{-13}$ erg s*-1* cm*-2* in their 2012 monitoring campaign, and the difference is within the range of total [O III] variability observed for NGC 5548 over the course of 21 years (see Peterson et al., 2013).

In addition to the intrinsic variability of the AGN, many other factors contribute to nightly variations in the spectra. These include changes in transparency due to clouds, changes in seeing conditions, inconsistent instrument focus, and miscentering of the AGN in the slit during observations. We used the flux scaling method described by van Groningen & Wanders (1992) to align the nightly spectra and place them on a consistent flux scale. For each spectrum in the dataset, the algorithm looks for a combination of wavelength shift, multiplicative scale factor, and Gaussian kernel convolution that minimizes the residual between each individual spectrum and a reference spectrum over a region containing the narrow [O III] line.

We constructed a separate reference spectrum for each telescope by averaging the highest SNR spectra in each dataset, then broadened the reference spectrum so that the [O III] line width matches the broadest [O III] line width in the dataset. This extra broadening of the reference spectrum helps to reduce the [O III] residuals from spectral scaling (Fausnaugh, 2016). We then scaled each spectrum to have the same [O III] flux as the photometrically calibrated mean MDM spectrum. This brings all spectra to a common flux scale after spectral scaling.

To assess the accuracy of spectral scaling, we estimated the intrinsic fractional variability of the residual [O III] $\lambda$5007 light curve after correcting for random measurement errors,

[TABLE]

where $\sigma^{2}$ is the [O III] flux variance, $\langle\delta^{2}\rangle$ is the mean-square value of the measurement uncertainties determined from the nightly error spectra produced by the data reduction pipeline, and $\langle f\rangle$ is the unweighted mean flux. The $F_{\mathrm{var}}$ for the [O III] $\lambda$5007 light curve gives a good estimate of the residual flux-scaling errors (Barth & Bentz, 2016), and the value for each telescope is listed in the last column of Table 1. We found $F_{\mathrm{var}}$ to be between 0.27% and 0.62% for all telescopes, which means there is an additional scatter of less than 1% in the [O III] light curve above the measurement errors. These $F_{\mathrm{var}}$ values are consistent with or better than the best values typically obtained in ground-based campaigns. For example, Barth et al. (2015) found $F_{\mathrm{var}}$ values ranging from 0.5% to 3.3% for individual AGNs in the 2011 Lick AGN Monitoring Project.

Figure 2 shows (in black) the mean and root-mean-square (rms) residual spectra for the MDM dataset. The rms spectrum indicates the degree of variability at each wavelength over the course of the campaign. Both the broad H$\beta$ and He II $\lambda$4686 emission lines exhibit strong variations, and the H$\beta$ rms profile appears to have multiple peaks. Traditionally, the rms spectrum is constructed such that the value at each wavelength is taken to be the standard deviation of fluxes from all epochs, but this does not take into account Poisson or detector noise, which may bias the rms profile by a small amount (Barth et al., 2015). Park et al. (2012b) suggest using the SNR for each spectrum as the weight for that spectrum in calculating the rms, or using a maximum-likelihood method to obtain the rms. We adopt a simpler approach that uses the excess variance as a way to exclude variations that are not intrinsic to the AGN. This “excess rms” value at each wavelength is defined as

[TABLE]

where $N$ is the total number of spectra in the dataset, $\langle F_{\lambda}\rangle$ is the mean flux at each wavelength, and $F_{\lambda,i}$ and $\delta_{\lambda,i}$ are the wavelength-specific fluxes and associated measurement uncertainties from individual epochs, respectively. This method estimates the degree of variability above what is expected given the measurement uncertainties and pixel-to-pixel noise.

3. Spectroscopic Flux Measurements

The 5100 Å continuum flux density was determined by averaging the flux over the rest-frame wavelength range 5070$-$5130Å. The H$\beta$ line fluxes were measured from the scaled spectra using the same method as for [O III] $\lambda$5007, where we subtracted a linear fit to the surrounding continuum (wavelength windows 4483.0$-$4542.0 Å and 5033.5$-$5092.5 Å) and integrated across the line profile (4748.4$-$4945.1 Å). The uncertainty in each measurement is a combination of Poisson noise and residuals from spectral scaling. We computed the spectral scaling uncertainty by multiplying each flux measurement by the [O III] $F_{\mathrm{var}}$ value for that dataset, then adding this value in quadrature to the Poisson noise to obtain the final flux uncertainty for each measurement. There is an additional source of spectral scaling uncertainty from slight differences in the overall spectral shape from night to night. This effect is likely small for H$\beta$ because it is very close to the [O III] $\lambda$5007 line that anchors the spectral scaling.

Spectrophotometric calibrations of the reference spectra, as described in the previous section, converted all instrumental fluxes to absolute fluxes, which means that measurements from all telescopes should now be on the same flux scale. However, light curves from different observing sites may be offset from each other owing to aperture effects (Peterson et al., 1995, 1999). While our observations were standardized to have the same 5″ $\times$ 15″ aperture size, significant differences in image quality between observing sites could still cause flux offsets.

To intercalibrate the H$\beta$ light curves, we used data points from each non-MDM telescope ($F_{\mathrm{H\beta,t}}$) that are nearly contemporaneous with MDM observations ($F_{\mathrm{H\beta,MDM}}$) and performed a least-squares fit to the equation

[TABLE]

to find the scale factor $\phi$ that puts each line light curve on the same flux scale as the MDM data. For the continuum intercalibration, we also include an additive shift $G$ to account for the differences in the host-galaxy flux admitted by different apertures:

[TABLE]

The scale factors for the Lick, Asiago, APO, and WIRO light curves are $\phi=[0.961,~{}0.963,~{}1.037,~{}0.918]$, and the shift constants are $G=$[$-0.155$, $-0.640$, $-0.041$, $0.024$] in units of $10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$. The combined continuum and H$\beta$ light curves are shown in Figure 3 (THJD = HJD $-$ 2,450,000), and the 5100 Å continuum and H$\beta$ fluxes are listed in Table 2.

We attempted to measure the He II $\lambda$4686 flux from the nightly spectra. However, this line is very weak and also heavily blended with the broad H$\beta$, as shown in Figure 2. Thus, we were unable to obtain a He II light curve using the linear interpolation method to remove the continuum.

3.1. Spectral Decomposition

To more accurately remove the continuum underlying the emission lines and to deblend the broad emission features from each other, we employed the spectral decomposition algorithm described by Barth et al. (2015). The components fitted in this procedure include narrow [O III], broad and narrow H$\beta$, broad and narrow He II, Fe II emission blends, the stellar continuum, and the AGN continuum. The host-galaxy starlight was modeled with an 11 Gyr, solar metallicity, single-burst spectrum from Bruzual & Charlot (2003). For the Fe II model component, we tested three different templates from Boroson & Green (1992), Véron-Cetty et al. (2004), and Kovačević et al. (2010). The Fe II templates were broadened by convolution with a Gaussian kernel in velocity. The free fit parameters for Fe II include the velocity shift relative to broad H$\beta$, the broadening kernel width, and the flux normalization of the broadened template spectrum. The Boroson & Green (1992) and Véron-Cetty et al. (2004) templates are monolithic and require only one flux normalization parameter, whereas the Kovačević et al. (2010) template has five components that can vary independently in flux. The Kovačević et al. (2010) template achieves the best fit to the nightly spectra, presumably a result of the larger number of free fit parameters due to the multi-component Kovačević et al. (2010) template.

We made several modifications to the spectral fitting procedures used by Barth et al. (2015). First, because of the complex line profiles, we used sixth-order Gauss-Hermite functions (van der Marel & Franx, 1993) to fit the broad and narrow H$\beta$ and narrow [O III] lines instead of fourth-order functions. Second, there is significant degeneracy between the weak Fe II blend and the continuum flux in the nightly fits. Since the Fe II fit is poorly constrained and sometimes varied drastically from night to night, the continuum model flux also varied significantly as a result, which in turn introduced noise to the broad H$\beta$ fit component. To address this issue, we constrained the Fe II flux to lie within 10% of the value from the fit to the mean spectrum (Barth et al., 2013). We also fixed the Fe II redshift to that of the mean spectrum and constrained the Fe II broadening kernel to be within 5% of its value from the mean spectrum fit. The He I $\lambda$4922 and $\lambda$5016 lines are very weak and are heavily blended with broad H$\beta$, making it impossible to constrain their fit parameters. We therefore do not fit for these components in our model.

The broad He II $\lambda$4686 component has very low amplitude compared to the other fit components and it is blended with the blue wing of broad H$\beta$. It is also highly variable, as demonstrated by the broad bump in the rms spectrum. This made it difficult to fit the He II broad-line profile accurately, and the width varied significantly from night to night when fitted as a free parameter. Since the He II $\lambda 1640$ and $\lambda 4868$ lines are expected to form under the same physical conditions and should thus have similar widths, we used fits to the $\lambda 1640$ line in concurrent HST spectra to constrain the $\lambda 4686$ line width.

The He II $\lambda 1640$ line was modeled with five Gaussian components (De Rosa et al., in prep.), and we took the three broadest components to represent the broad He II $\lambda 1640$ line profile. For each MDM spectrum, the He II $\lambda 4686$ broad-line full width at half-maximum intensity (FWHM) was allowed to vary within 3 Å of the He II $\lambda 1640$ FWHM measured from the closest HST epoch. The first 23 epochs from the MDM campaign do not have corresponding HST spectra, so for each of these “pre-HST” epochs, we found the three epochs from later in the campaign with the closest matching 5100 Å continuum flux density. We then used the weighted mean of the broad He II $\lambda 1640$ widths from these three nights as the width constraint for the pre-HST epoch, where the weight were determined by how closely the 5100 Å fluxes of the later epochs matched that of the pre-HST epoch. The He II $\lambda 1640$ line width was highly variable during the HST campaign, and the model FWHM widths used to constrain the spectral decomposition have a mean of 48 Å, with a minimum of 28 Å and maximum of 59 Å.

We applied spectral decomposition to the data from all telescopes, but since the MDM dataset is the largest and has the highest data quality and consistency, we use this dataset for all subsequent analysis. Figure 4 shows the fit components for the mean MDM spectrum, where the black spectrum is the data and the red spectrum is the sum of all the model components. The model does not fit the detailed structure of the broad H$\beta$ line well, especially in the line core. To prevent this from impacting our measured H$\beta$ fluxes, we subtracted all the other well-modeled fit components except the broad and narrow H$\beta$ components from the full spectrum, then obtained the H$\beta$ line flux by integrating over the same wavelength range used to measure the flux without spectral decomposition. The He II $\lambda$4686 flux was taken to be the total flux in the broad- and narrow-line models for each night. The narrow H$\beta$ and He II $\lambda$4686 line fluxes from fits to the mean spectrum are $48.4\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}$ and $8.5\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}$ (respectively), with uncertainties of $\sim 2$% from the overall photometric scale of the data. The ratio of the narrow H$\beta$ flux to the [O III] $\lambda 5007$ flux is $F_{\mathrm{H\beta}}/F_{\mathrm{[O~{}III]}}=0.099\pm 0.002$, which is in good agreement with the value of $F_{\mathrm{H\beta}}/F_{\mathrm{[O~{}III]}}=0.110\pm 0.010$ found by Peterson et al. (2004).

The red spectrum in Figure 2 shows the rms of the MDM spectra after subtracting the AGN continuum and stellar continuum models from individual spectra so that only the emission-line components remain. This rms spectrum is expected to be a more accurate representation of the emission-line variability than the rms of the full spectra (Barth et al., 2015). We show the rms here and not the excess rms defined by Equation 4 because, for parts of the spectra dominated by continuum emission, the continuum-subtracted flux could be lower than the total flux uncertainties and the e-rms would be undefined. Thus, we use the excess rms only for the full spectrum and not for individual fit components.

Panels (a) and (c) of Figure 5 show the H$\beta$ mean and rms fluxes as a function of line-of-sight velocity ($v_{\mathrm{r}}$) after spectral decomposition, and panel (b) shows the difference between the T1 and T2 mean fluxes. The rms flux has a statistical uncertainty of $\sim 12\%$ and the [O III] residuals are much lower compared to the case with no spectral decomposition (Fig. 2). The rms profile still has jagged features, which likely reflect real variability across the broad emission line.

Figure 6 shows the 1158 Å UV continuum light curve from Paper I, the MDM optical 5100 Å continuum light curve, the V-band photometric light curve from Paper III, and the MDM H$\beta$ and He II $\lambda$4686 emission-line light curves. The He II light curve reaches a flat-bottomed minimum near THJD = 6720. This is because the He II flux includes contributions from the broad- and narrow-line components, so when the broad-line flux is near zero, the total He II line flux stays at a minimum value equal to the narrow-line flux. Light-curve statistics that quantify the variability of NGC 5548 during the monitoring period are given in Table 3. $F_{\mathrm{var}}$ is as defined in Equation 3 and $R_{\mathrm{max}}$ is the ratio between the maximum and minimum fluxes.

3.2. Host-Galaxy Flux Removal

We measured the host-galaxy contribution to the continuum using an “AGN-free” image of NGC 5548 generated by Bentz et al. (2013) after performing two-dimensional surface brightness decomposition on HST images of the galaxy. We found that the amount of starlight expected through a 5″$\times$ 15″ aperture with a slit position angle of $0\degree$ is $F_{\mathrm{5100,gal}}=(4.52~{}\pm~{}0.45)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$. Subtracting this from the mean continuum flux density of $F_{\mathrm{5100}}=(11.96~{}\pm~{}0.07)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$ gives a mean AGN flux of $F_{\mathrm{5100,AGN}}=(7.44~{}\pm~{}0.50)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$, which is consistent with the value of $F_{\mathrm{5100,AGN}}=(7.82~{}\pm~{}0.02)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$ measured from the power-law component of the spectral decomposition for the mean MDM spectrum.

3.3. Anomalous Emission-Line Light-Curve Behavior

RM analyses typically assume that the emission-line light curve responds linearly to continuum variations and is a lagged, scaled, and smoothed version of the continuum light curve. However, this does not appear to be the case for a portion of our campaign. As described in Paper I and Paper IV, there are significant differences between the UV continuum and emission-line light curves after THJD = 6780. The continuum flux increased while the emission-line fluxes either decreased or remained roughly constant in a suppressed state for the remainder of the campaign.

The H$\beta$ emission line shows similar behavior to that of Ly$\alpha$ and C IV, in that there is a marked difference in the line response between the first and second halves of the campaign. This “decorrelation” phenomenon is illustrated in Figure 7(a), where the H$\beta$ and UV continuum light curves trace each other well in the first half of the campaign, but the continuum flux continues to trend upward beyond THJD = 6740, while the H$\beta$ flux begins to fall. For the remainder of the campaign, the H$\beta$ flux remains in a suppressed state and the light curve does not follow the continuum light curve well in that prominent features in the continuum light curves (e.g. THJD $\approx$ 6770, 6785) are not present in the emission-line light curve. The He II $\lambda$4686 light curve behaves similarly and decorrelates from the continuum at around THJD = 6760 (Figure 7b). This phenomenon is also apparent when comparing the H$\beta$ light curve to the 5100 Å continuum, as shown in Figure 7(c). Though there are small differences between the light curves in T1 (such as around THJD 6685), the overall correlation in T1 is significantly better than in T2.

Owing to this change in the emission-line response, we followed the procedures presented in Paper I for determining the UV emission-line lags, and divided the 5100 Å continuum and optical emission-line light curves into two subsets, labeled T1 and T2, to examine the lag of each segment separately. The subsets are divided at THJD = 6747 and each has 67 epochs. The 1158 Å and V-band continuum light curves were also separated into two segments at THJD = 6747. Note that this is a different dividing epoch from the one used in Paper I. The characteristics of the half-campaign light curves are given in the bottom portion of Table 3. Figure 5 shows the MDM mean and rms spectra for T1 and T2 in gray and orange, respectively. While the three mean spectra look almost identical, the T1 and T2 rms spectra are significantly different, which indicates changes in the amount of variability between the two campaign halves.

4. Data Analysis

In the following sections, we examine properties of the BLR by measuring the emission-line responses to continuum variations. We also discuss the anomalous behavior of the emission-line light curves observed during this campaign and BH mass measurements using this dataset.

4.1. Emission-Line Lags

We measured the H$\beta$ and He II $\lambda$4686 lags relative to both the 5100 Å continuum and the 1158 Å continuum. All light curves were detrended by subtracting a linear least-squares fit to the data to remove long-term trends that may bias lag calculations (Welsh, 1999). In this case, we found very weak trends for all the light curves, and detrending has a very small effect ($\sim 0.01$ days) on the measured lags. We computed the cross-correlation coefficient $r$ for lags between $-20$ and 40 days in increments of 0.25 days using the interpolated cross-correlation function (ICCF; White & Peterson, 1994). Two lag estimates were made for each light-curve pair—the value corresponding to $r_{\mathrm{max}}$ ($\tau_{\mathrm{peak}}$) and the centroid of all values with $r>0.8r_{\mathrm{max}}$ ($\tau_{\mathrm{cen}}$). Estimates for the final $\tau_{\mathrm{peak}}$ and $\tau_{\mathrm{cen}}$ values and their uncertainties were obtained using Monte Carlo bootstrapping analysis (Peterson et al., 2004), where many realizations of the continuum and emission-line light curves were created by randomly choosing $n$ data points with replacement from the observed light curves, where $n$ is the total number of points in the dataset. If a data point is picked $m$ times, then the uncertainty on that point is decreased by a factor of $m^{1/2}$. Each value is then varied by a random Gaussian deviate scaled by the measured flux uncertainty. We constructed $10^{3}$ realizations of each light curve and computed the cross-correlation function (CCF) for each pair of line and continuum light-curve realizations to create a distribution of $\tau_{\mathrm{peak}}$ and $\tau_{\mathrm{cen}}$ values. The median value from each distribution and the central 68% interval are then taken to be the final lag and its uncertainty.

Table 4 lists the ICCF lags for H$\beta$ and He II $\lambda$4686 measured against the 1158 Å, 5100 Å, and V-band continua. The lag between the 5100 Å and 1158 Å continua is also given. For comparison with Paper I, which presents the UV emission-line lags against the 1367 Å continuum, we also include H$\beta$ and He II lags measured against this continuum. Distributions of $\tau_{\mathrm{cen}}$ values from the Monte Carlo bootstrap analysis using the 1158 Å continuum are shown in the top panels of Figure 8.

We also computed H$\beta$ and He II $\lambda$4686 lags for time periods T1 and T2, and found the T1 lags to be consistently longer by about 2$\sigma$. For H$\beta$, the T1 and T2 lags bracket the full-campaign lag, while for He II $\lambda$4686, the full-campaign lag is shorter than both T1 and T2 lags. For comparison, Paper I found the Ly$\alpha$, Si IV, C IV, and He II $\lambda$1640 lags to be longer for T2 than for T1, which is the opposite of what we find for H$\beta$.

To illustrate the effects of spectral decomposition, we also include in Table 4 the ICCF lags for the H$\beta$ light curve where the fluxes were measured using the straight-line continuum-subtraction method and without spectral decomposition. We calculated lags for both the MDM-only and the multi-site H$\beta$ light curves, which were truncated at THJD = 6828.75 to exclude the last 10 epochs of MDM data (see Section 2). The lags measured with and without using spectral decomposition are consistent to within 1$\sigma$.

In addition to the ICCF method, we computed the emission-line lags using the JAVELIN suite of Python codes (Zu et al., 2011). We used JAVELIN to linearly detrend the light curves and model the AGN continuum variability as a damped random walk process (DRW; Kelly et al., 2009; Zu et al., 2013). JAVELIN explicitly models the emission-line light curves as smoothed, scaled, and lagged versions of the continuum light curve. Since the decorrelation of the line and continuum light curves during the latter half of the campaign clearly violates these assumptions, it is of interest to examine the consequences for the JAVELIN models. For these models, we simultaneously fit the H$\beta$ and He II light curves using either the 1158 Å, 5100 Å, or V-band continuum light curve. The AGN STORM light curves are too short to accurately determine the DRW damping timescale, so this value was fixed to $\tau_{\mathrm{DRW}}\approx 164$ days as derived by Zu et al. (2011) from fits to the 13-year light curve of NGC 5548 (Peterson et al., 2002). The precise value of $\tau_{\mathrm{DRW}}$ is not critical for the algorithm to work provided that it is approximately correct.

The top three panels of Figure 9 show the results of using JAVELIN directly on the observed data. Despite the long DRW time scale, the light-curve models show rapid fluctuations. The algorithm tries to match the suppressed flux in the line light curve to the continuum light curve, and because of the small uncertainties, it strongly prefers lag times for which the continuum and line light curves have minimal temporal overlap. This caused the posterior lag distribution to have multiple narrow peaks corresponding to lags that best desynchronize the light curves.

We can attempt to compensate for this problem by increasing the uncertainties to encompass the amplitude of the decorrelation. This requires scaling up the full-campaign light-curve uncertainties by factors of 5 and 3 for the continuum and line light curves, respectively. The lower panels in Figure 9 show the JAVELIN results from fitting to the light curves with scaled errors, where the broader uncertainties allow the algorithm to construct smooth light-curve models. Since there is more statistical weight from fitting both lines simultaneously, the models track the line light curves best and show a smooth systematic offset for the continuum where the line and continuum light curves are decorrelated. When computing the half-campaign lags, JAVELIN favors a smaller line flux scale factor for T2 to account for the suppressed line fluxes, which begins near the epoch separating T1 and T2. We therefore did not need to scale the flux errors by as much as for the full campaign to account for the decorrelation, and used an error scaling factor of 3 for both continuum and line light curves. The resulting JAVELIN lags, shown in Table 4, are consistent with those measured using the ICCF method. Similar to the ICCF lags, the JAVELIN lags for the full-campaign light curves are also shorter than those for T1 and T2. However, since the JAVELIN assumptions of the relationship between the line and continuum light curves are not valid for this campaign and the flux errors were scaled for the sole purpose of producing convergent solutions, we do not use the JAVELIN lags for subsequent analysis.

We examined the velocity-resolved emission-line response by dividing the H$\beta$ line profile into bins with velocity width of 500 km s*-1* and set zero velocity using the peak of the narrow H$\beta$ component in the mean spectrum. We constructed light curves for each velocity bin separately, and Figure 11 shows the light curves for 1500 km s*-1* velocity bins across the H$\beta$ line profile (black), with the velocity at the center of each bin shown in the top left of each panel. There were six epochs of spectra222THJD = [6689.0, 6699.0, 6757.8, 6758.9, 6796.8, 6797.8] that produced outlying H$\beta$ fluxes and significantly higher-than-average flux uncertainties for individual velocity bins, so we removed them from the velocity-resolved light curves in order to improve the lag measurements.

We determined the ICCF lag for each of these binned light curves with respect to both the UV and optical continua, and show these lags as a function of line-of-sight velocity in the top left panel of Figure 10. The second left panel shows the velocity-resolved H$\beta$–UV lag for T1 and T2, and the maximum cross-correlation coefficients ($r_{\mathrm{max}}$) are shown in the right-side panels. The bottom panel of Figure 10 shows the MDM full-campaign mean spectrum for reference, and Table 5 lists the velocity-resolved lags. Emission-line variations have lower amplitudes during T2 compared with T1, which led to poorly-constrained ICCF lags with large uncertainties for velocity bins with low line flux. We therefore reduced the upper limit of the CCF lag range from 40 days to 20 days when computing the T2 lags.

The velocity-resolved H$\beta$–UV lag for the full campaign is shortest ($\tau_{\mathrm{cen}}\approx 2$ days) in the line wings where $v_{\mathrm{r}}\approx\pm 7000$ km s*-1*. The lag increases as $v_{\mathrm{r}}$ approaches zero from both sides of the lag profile, reaches local maxima of $\tau_{\mathrm{cen}}\approx 10$ days at about $v_{\mathrm{r}}\approx\pm 3000$ km s*-1*, and then steadily decreases until it reaches a local minimum of $\tau_{\mathrm{cen}}\approx 4$ days near the line profile center. The lag profile measured against the optical continuum has a similar shape, but with all lags $\sim 2$–3 days shorter, as we would expect from the $\sim 2$-day lag between the two continua (Table 4). A similar double-peaked lag profile is also observed for Ly$\alpha$ (see Paper I). The T1 lag profile closely resembles that of the full campaign, but the T2 lag profile shows a slightly different structure. The bins where the T1 and T2 lags are most discrepant are also where the T2 light curves are least correlated with the continuum ($r_{\mathrm{max}}<0.4$). However, even excluding these outliers, there are still discernible differences between the T1 and T2 lag profiles. Additionally, the T2 $r_{\mathrm{max}}$ values are lower than those of T1 in every velocity bin, which clearly demonstrates that the line and continuum light curves are less correlated in T2 than in T1.

The shape of the velocity-resolved lag profile can provide qualitative information about the kinematics of the line-emitting gas (e.g., Kollatschny, 2003; Bentz et al., 2009b; Denney et al., 2010; Barth et al., 2011; Du et al., 2016a). In simple models of the BLR (Ulrich et al., 1984; Gaskell, 1988; Welsh & Horne, 1991; Horne et al., 2004; Goad et al., 2012; Gaskell & Goosmann, 2013; Grier et al., 2013), pure infall motion would lead to longer lags on the blue side of the line profile, and for outflow, the most redshifted gas would have the longest lag. For gas in Keplerian orbits, the shortest lags would be in the line wings, since gas with higher $v_{\mathrm{r}}$ is closer to the central black hole. Gas with very low $v_{\mathrm{r}}$ could have a wide range of lags, and a spherical or flat disk distribution of BLR clouds in Keplerian motion could lead to a double-peaked velocity-resolved lag profile if the ionizing source is emitting anisotropically (Welsh & Horne, 1991; Goad & Wanders, 1996; Horne et al., 2004).

Previous studies of the UV and optical lines in NGC 5548 have inferred either Keplerian orbits (Horne et al., 1991; Wanders et al., 1995; Denney et al., 2009; Bentz et al., 2010b) or infalling motion (Crenshaw & Blackwell, 1990; Done & Krolik, 1996; Welsh et al., 2007; Pancoast et al., 2014; Gaskell & Goosmann, 2016) for the BLR gas. From our data, the shape of the H$\beta$ velocity-resolved lag profile suggests a BLR dominated by Keplerian motion. The discrepancy between the T1 and T2 lag profiles may suggest a change in the distribution or dynamics of the BLR gas, though such changes typically occur on timescales much longer than our campaign. More detailed interpretation requires comparison with transfer functions generated for various dynamical models of the BLR, which will be the subject of future work in this series (Pancoast et al., in prep.).

Figure 11 also shows modified versions of the 1158 Å continuum light curve in red. These light curves were shifted in time by the average T1 lag of the bins incorporated in each H$\beta$ light curve, which are shown in the bottom right of each panel. The fluxes were roughly scaled and shifted to match the first half of the H$\beta$ light curves. Comparing the line and continuum light curves, it is evident that the H$\beta$ response in T2 is heavily dependent on the line-of-sight velocity.

4.2. Anomalous Emission-Line Response to Continuum

Our data show that previous assumptions about the relationship between the continuum and emission-line light curves—namely that the emission-line light curves are smoothed, scaled, and time-shifted versions of the continuum light curve—are not always valid, as we observed all the UV and optical emission-line light curves decorrelating from the UV continuum about halfway through the monitoring period. Paper IV examined this effect for the UV emission lines by measuring changes in the emission-line equivalent width (EW),

[TABLE]

and the responsivity ($\eta_{\mathrm{eff}}$), which is the power-law index that relates the driving continuum flux to the responding emission-line fluxes,

[TABLE]

In the case of no line response, $\eta_{\mathrm{eff}}=0$, and if the line responds linearly to continuum variations (i.e. the transfer function is a $\delta$ function), then $\eta_{\mathrm{eff}}=1$. The emission-line EW and the continuum flux are related via the Baldwin relation (Baldwin, 1977), which is described by

[TABLE]

where the choice for $F_{\mathrm{cont}}$ is assumed to be a reasonable proxy for the ionizing continuum. Thus, $\beta$ is also known as the slope of the Baldwin relation.

Following the same procedures as in Paper IV, we compute the responsivity $\eta_{\mathrm{eff}}$ and EW for the portion of the H$\beta$ light curve that correlates with the UV continuum, then examine how these values change in different segments of the light curves. The values $F_{\mathrm{cont}}$ and $F_{\mathrm{line}}$ refer to the continuum and emission-line fluxes after removing non-variable components such as host galaxy and narrow-line flux contributions, and after correcting for the mean time delay between the continuum and line light curves (Pogge & Peterson, 1992; Gilbert & Peterson, 2003; Goad et al., 2004). There is very little host galaxy flux in the 1158 Å continuum, which is dominated by the variable AGN. For the line fluxes, we took H$\beta$ fluxes measured after linear continuum removal (without spectral decomposition) and subtracted a constant narrow H$\beta$ flux measured from the MDM mean spectrum fit to remove the non-variable line component. To correct for the emission-line time delay, we shifted the H$\beta$ light curve by 8 days, which corresponds to the lag for the portion of the line light curve closely correlated with the UV continuum. Figure 12(a) shows the 1158 Å continuum light curve (black) with the time-shifted and flux-scaled H$\beta$ light curve, which has been truncated at the beginning to match the first epoch of continuum observations. To show the H$\beta$ light curve’s general behavior toward the end of the campaign, we have shown here the full 143 epochs of H$\beta$ flux measurements instead of the 133-epoch light curve we used in the H$\beta$ lag analysis. However, since the last 10 epochs of spectra suffer from inconsistent spectral flux calibration (see Section 2), we do not use these points in calculating $\eta_{\mathrm{eff}}$ or $\beta$.

We divided the H$\beta$ light curve into five segments. The first segment corresponds to the period when the line light curve closely follows the continuum light curve (blue points); the second and third segments correspond to periods when the line light curve decouples from the continuum (cyan points) and remains in a state of depressed flux (red points); the last two segments correspond to the line light curve recovering from the depressed state (magenta points) and correlating once again with the continuum light curve (green points). The epochs that divide these segments are THJD = [6743, 6772, 6812, 6827].

Figure 12(d) and 12(e) show the H$\beta$ broad-line flux and EW as a function of the 1158 Å continuum flux density determined from the HST epoch closest to each MDM epoch. The red lines represent linear least-squares fits to Equations 8 and 9 using only the blue points. We found that the cyan, red, and magenta points—corresponding to when the light curves are not well-correlated—lie well below the best fits of Equations 8 and 9 to the blue points. Furthermore, the epochs during the anomaly (red points) are characterized by an $\eta_{\mathrm{eff}}$ value consistent with zero ($\eta_{\mathrm{eff}}$ $=0.02\pm 0.03$, black dotted line in Figure 12d), which shows that the emission-line strength remained constant independent of continuum strength. Similar results were found for the C IV, Ly$\alpha$, He II(+O III]) and Si IV(+O IV]) emission lines in Paper IV. Table 6 summarizes the $\eta_{\mathrm{eff}}$ and $\beta$ values for all broad emission lines measured during this campaign. Note that the values from Paper IV were computed using epochs both before and after the anomaly (blue and green points), whereas our values were calculated using only epochs before the anomaly (blue points).

We measured the amplitude of the anomaly by comparing the observed H$\beta$ light curve to a simulated light curve that represents what the line response would have been without the anomaly. The simulated line fluxes ($F_{\mathrm{sim}}$) were calculated using the observed continuum fluxes and Equation 8, where $\eta_{\mathrm{eff}}$ and $A$ are chosen to be the best-fit values for the period when the light curves are well-correlated (blue points). This simulated light curve is shown in Figure 12(b) with gray dots, and comparing it with the observed data (colored points) shows that the H$\beta$ line had lower-than-expected variability amplitude and mean flux during T2. If we define the fractional flux loss as $F_{\mathrm{lost}}=(F_{\mathrm{sim}}-F_{\mathrm{obs}})/F_{\mathrm{sim}}$ (Figure 12c), then this implies an H$\beta$ flux deficit of $\sim 6\%$ during the anomaly (red points), which is close to the deficit for Ly$\alpha$ but much smaller than that of the other UV emission lines, as shown in Table 6.

Table 7 summarizes the H$\beta$ responsivities and AGN optical continuum flux densities measured by Goad et al. (2004, GKK04) for every year of the 13-year monitoring campaign carried out by the AGN Watch consortium (Peterson et al., 2002)333The continuum flux densities from Peterson et al. (2002) have been updated by Bentz et al. (2013) and Kilerci Eser et al. (2015)., along with the values of $\eta_{\mathrm{eff}}=0.59\pm 0.01$ and $\langle F_{\mathrm{5100,\mathrm{AGN}}}\rangle=(7.44\pm 0.50)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}~{}\mathrm{\AA}^{-1}$ calculated from this dataset. Both the responsivity and optical continuum flux density for this campaign are close to those measured from the 1998 data. Using the mean H$\beta$ flux without narrow-line contributions of $\langle F_{\mathrm{H\beta}}\rangle=(690\pm 2.40)\times 10^{-15}~{}\mathrm{erg~{}s}^{-1}~{}\mathrm{cm}^{-2}$, we find a mean EW for this campaign of $92.75\pm 2.45$ Å, which is lower than values measured by Goad & Korista (2014) for previous campaigns.

The H$\beta$ light curve appears to decorrelate from the UV continuum light curve at a somewhat earlier time (THJD $\approx 6742$) than C IV (THJD $\approx 6765$ as found in Paper IV). If we assume that the H$\beta$ light curve decorrelates and re-correlates with the continuum light curve at the same times as C IV and use the same dividing epochs as those used in Paper IV (THJD = [6766, 6777, 6814, 6830]), then $\eta_{\mathrm{eff}}=0.13\pm 0.01$ and $\beta=-0.88\pm 0.01$ for the blue points.

While He II $\lambda$4686 also shows anomalous behavior during the T2 period (Fig. 7), its broad line component is very weak and the fitted profile is poorly constrained in the spectral decomposition process. The He II light curve is thus noisier than that of H$\beta$, and the $F_{\mathrm{line}}$ and EW values are poorly correlated with $F_{\mathrm{cont}}$. We therefore do not perform a detailed analysis of the responsivity for He II.

4.3. Line Width and $M_{\mathrm{BH}}$ Estimate

The black hole mass in NGC 5548 has been estimated for several previous RM campaigns, including the AGN Watch consortium (Peterson et al., 2002), Bentz et al. (2007), Denney et al. (2010), and Bentz et al. (2009b). The AGN Watch group determined the BH mass using data from each of the program’s monitoring years, and subsequent campaigns each produced an independent BH mass value. Here we compute the BH mass using data from this campaign and compare it with previous results.

We measured the line widths from the H$\beta$ mean and rms spectra as shown in Figure 5. The narrow H$\beta$ line essentially disappears in the rms spectrum but is still present in the mean spectrum. We removed this emission component by subtracting a Gaussian fit to the H$\beta$ narrow line in the mean MDM spectrum, then linearly interpolated over the narrow H$\beta$ residuals and the [O III] $\lambda$4959 and $\lambda$5007 residuals.

Two emission-line width values are typically measured in RM: the FWHM and the line dispersion, defined by

[TABLE]

where $S_{i}$ is the flux density at wavelength bin $\lambda_{i}$ and $\lambda_{0}$ is the flux-weighted centroid wavelength of the line profile. We measured $\sigma_{\mathrm{line}}$ and FWHM for the mean and rms spectra using the line profile within the rest-frame wavelength range 4669.8$-$5063.0 Å. We treated the mean profile as double-peaked and followed the procedures described by Peterson et al. (2004) to measure its FWHM. From each of the two peaks at 4827.1 Å and 4886.1 Å, we traced the line profile outward until the flux reached $0.5F_{\mathrm{max}}$, then traced the line profile inward from the continuum until the flux again reached $0.5F_{\mathrm{max}}$. The two wavelengths at $0.5F_{\mathrm{max}}$—which generally agree well for a smooth line profile—are then averaged to obtain the wavelength at half maximum on each side of the profile. The rms profile is more complicated because it has more than two peaks and the troughs between them can reach well below half of the peak fluxes. We therefore identified a single maximum flux $F_{\mathrm{max}}$ and traced the profile from the continuum toward the center on both sides until the flux reached $0.5F_{\mathrm{max}}$. The separation between the two wavelengths at $0.5F_{\mathrm{max}}$ is taken to be the FWHM of the rms spectrum.

The H$\beta$ line widths and their uncertainties were determined using Monte Carlo bootstrap analysis. With $n$ total spectra, we randomly selected $n$ spectra from the dataset with replacement, constructed mean and rms line profiles, and measured the line dispersion and FWHM. The median and standard deviation of $10^{4}$ bootstrap realizations are used for the $\sigma_{\mathrm{line}}$ and FWHM and their estimated uncertainties.

There are additional systematic uncertainties in the line widths from using different Fe II templates in the spectral decomposition. We repeated the bootstrap analysis after performing spectral decompositions using each of the Boroson & Green (1992), Véron-Cetty et al. (2004) and Kovačević et al. (2010) Fe II templates, then took the standard deviation of the H$\beta$ widths from using the different templates as the systematic uncertainty for the line width. This systematic error dominates the error budget for all $\sigma_{\mathrm{line}}$ measurements and the FWHM of the mean spectrum, but is comparable to the uncertainty from bootstrapping analysis for the FWHM of the rms spectrum. We added this systematic uncertainty in quadrature to the statistical uncertainty from the Kovačević et al. (2010) line width to obtain the final H$\beta$ line width uncertainty.

The line widths are also affected by instrumental broadening due to the use of a wide slit. The observed line width is the quadratic sum of the intrinsic and instrumental line widths ($\sigma_{\mathrm{observed}}^{2}=\sigma_{\mathrm{intrinsic}}^{2}+\sigma_{\mathrm{instrumental}}^{2}$ or similarly for FWHM). To calculate the instrumental broadening for a 5″ slit, we followed the methods described by Bentz et al. (2009b) and compare the [O III] $\lambda$5007 line width measured from our observations to the width measured from a higher-resolution observation taken using a narrow ($\sim 2$″) slit, which represents the intrinsic line width. The FWHM of the [O III] $\lambda$5007 model from the MDM rest-frame mean spectrum is 9.79 Å (572 km s*-1*), while Whittle (1992) found a FWHM of 410 km s*-1* using a 2″ slit. This implies an instrumental broadening of FWHM${}_{\mathrm{instrument}}=$ 399 km s*-1*, corresponding to $\sigma_{\mathrm{instrument}}=$ 170 km s*-1* for a Gaussian model of the [O III] line profile. We subtract this instrumental width in quadrature from the H$\beta$ line width measurements to obtain the intrinsic H$\beta$ line widths, which are listed in Table 8. The large uncertainties on the FWHM measurements from the rms spectrum are a result of the jagged shape of the rms profile.

We use the H$\beta$ line dispersion measured from the rms spectrum as the velocity dispersion $\Delta V$, as is common practice in RM (Peterson et al. 2004, but also see Collin et al. 2006 and Mejía-Restrepo et al. 2016 for comparison of FWHM and sigma as indicators of BLR virial velocity for BH mass estimation), and calculate the virial product, defined as $\mathrm{VP}=c\tau\Delta V^{2}/G$. The f factor in Equation 2, which incorporates the geometry and kinematics of the BLR, is generally unknown for any individual AGN, so a single value $\langle f\rangle$ is often used to represent the average normalization for all AGNs. This value is usually taken to be the scale factor that puts the sample of RM virial products onto the same $M_{\mathrm{BH}}-\sigma_{\star}$ relation as nearby inactive galaxies (see Kormendy & Ho 2013 for a discussion of the uncertainties in the $M_{\mathrm{BH}}-\sigma_{\star}$ relation), and can vary depending on the sample of AGNs used in the fit as well as the fitting method (Onken et al., 2004; Watson et al., 2007; Woo et al., 2010; Graham et al., 2011; Park et al., 2012a; Grier et al., 2013; Ho & Kim, 2015). We adopt a value of $\langle f\rangle=4.47$ as calculated by Woo et al. (2015). Since $\langle f\rangle$ is calibrated using H$\beta$ lags measured against the optical continuum, we used the H$\beta$–5100 Å lag to calculate the VP and BH mass.

Table 9 lists the virial products and the inferred BH masses for the full, T1, and T2 segments. Our $M_{\mathrm{BH}}$ uncertainties do not include uncertainties in $\langle f\rangle$ (Woo et al., 2015) or scatter in the distribution of $f$ for different AGNs (e.g. Ho & Kim, 2015). The T1 and T2 $M_{\mathrm{BH}}$ estimates are consistent even though the two lags are different by more than 1$\sigma$. Compared to recent measurements, the full-campaign virial product ($1.49^{+0.49}_{-0.49}\times 10^{7}$ $M_{\odot}$) is entirely consistent with the values of $1.50^{+0.37}_{-0.51}\times 10^{7}$ $M_{\odot}$ and $1.38^{+0.51}_{-0.41}\times 10^{7}$ $M_{\odot}$ obtained by Bentz et al. (2009b) and Lu et al. (2016), respectively. Our BH mass measurement ($6.66^{+2.17}_{-2.17}\times 10^{7}$ $M_{\odot}$) is also consistent to 1$\sigma$ with the dynamical mass of $3.24^{+2.26}_{-0.90}\times 10^{7}$ $M_{\odot}$ as determined by Pancoast et al. (2014, P14). P14 use the V-band continuum as the ionizing source but do not use the f-factor, which compensates for the difference in lag between using the UV and optical continua (see next section). This will lead to an underestimate of $R_{\mathrm{BLR}}$ and hence a proportional underestimate of the BH mass (Eq. 2). If we scale the P14 mass by the ratio between the H$\beta$–1158 Å and the H$\beta$–V-band lags from Table 4 ($\tau_{\mathrm{H\beta-UV}}$/$\tau_{\mathrm{H\beta}-V}=1.64\pm 0.12$), then the two $M_{\mathrm{BH}}$ measurements become much more congruent.

5. Discussion

5.1. Implications of UV and Optical H$\beta$ Lags

Ground-based RM campaigns have traditionally used the optical continuum light curve—by necessity—to determine emission-line lags, even though the far-UV continuum is a better proxy for the ionizing source. Lags relative to the UV continuum ($\tau_{\mathrm{H\beta-UV}}$) thus should yield more accurate estimates of the BLR characteristic radius than lags measured relative to the optical continuum. Our H$\beta$–UV lag ($\tau_{\mathrm{H\beta-UV}}$ $=6.23^{+0.39}_{-0.44}$ days) is $\sim 2$ days longer than the H$\beta$–optical lag ($\tau_{\mathrm{H\beta-opt}}$ $=4.17^{+0.36}_{-0.36}$ days). Given that past measurements of the H$\beta$-optical lag for this object range from $\sim 4$ to $\sim 25$ days (Bentz et al., 2013, and references therein), and assuming that $\tau_{\mathrm{H\beta-opt}}$ is always $\sim 2$ days longer than $\tau_{\mathrm{H\beta-UV}}$, then the BLR characteristic radius estimated from previous campaigns using only optical data is biased low by $10-50\%$. The difference in these H$\beta$ lags is also consistent with the optical-to-UV continuum lag of $2.24^{+0.24}_{-0.24}$ days found in Paper III.

Since virial estimates of $M_{\mathrm{BH}}$ scale with $R_{\mathrm{BLR}}$, it may seem that this difference between $\tau_{\mathrm{H\beta-UV}}$ and $\tau_{\mathrm{H\beta-opt}}$ will change the BH mass estimate for NGC 5548 and other reverberation-mapped AGNs. However, the virial product—not the BH mass—is the quantity that is directly affected, and the normalization factor f is still needed to scale the virial product to a calibrated BH mass (Eq. 2). If the ratio of $\tau_{\mathrm{H\beta-UV}}/\tau_{\mathrm{H\beta-opt}}$ is the same for all AGNs, then all RM virial products would be scaled up by a constant value, so simply changing the value of f would remove this bias and leave the RM black hole masses unchanged. Even if the lag ratio is not constant for all AGNs, its effect on the BH mass scale is still small because the largest source of $M_{\mathrm{BH}}$ uncertainty in RM is the calibration uncertainty for f at $\sim 0.12$ dex (Woo et al., 2015). Furthermore, for NGC 5548, the H$\beta$ lag measured during this campaign is relatively short compared to its historical values (see Section 5.3). If the lag had been longer, the ratio of $\tau_{\mathrm{H\beta-UV}}/\tau_{\mathrm{H\beta-opt}}$ would be closer to unity and the bias in the BLR characteristic size and BH mass estimate would be much smaller.

While the discrepancy between $\tau_{\mathrm{H\beta-UV}}$ and $\tau_{\mathrm{H\beta-opt}}$ may not change $M_{\mathrm{BH}}$ measurements that use the f-factor, dynamical models that directly infer $M_{\mathrm{BH}}$ and BLR characteristics (e.g., Pancoast et al., 2011; Li et al., 2013) using only optical continuum data are affected since they do not depend on this virial normalization. The BLR characteristic size for each AGN inferred using optical data alone would thus be biased by a factor that depends on the value of $\tau_{\mathrm{H\beta-UV}}/\tau_{\mathrm{H\beta-opt}}$ for that particular object. Changes in the inferred $R_{\mathrm{BLR}}$ could also impact single-epoch $M_{\mathrm{BH}}$ estimates, which rely on the empirical relation between the BLR characteristic radius and the AGN continuum luminosity ($R_{\mathrm{BLR}}\propto L_{\mathrm{AGN}}^{\alpha}$ with $\alpha\approx\nicefrac{{1}}{{2}}$, e.g., Laor, 1998; Wandel et al., 1999; Kaspi et al., 2000; McLure & Jarvis, 2002; Kaspi et al., 2005; Vestergaard & Peterson, 2006; Bentz et al., 2006, 2009a, 2013). If $\tau_{\mathrm{H\beta-UV}}/\tau_{\mathrm{H\beta-opt}}$ correlates with AGN luminosity, then the expected value of $\alpha$ would change, and if this ratio is different for all AGNs but is uncorrelated with any other AGN properties, then this would introduce additional scatter to the scaling relation.

We can examine the expected scaling of $\tau_{\mathrm{H\beta-UV}}/\tau_{\mathrm{H\beta-opt}}$ with respect to $M_{\mathrm{BH}}$ and $L_{\mathrm{AGN}}$ by using simple disk and BLR ionization models. Assuming $\tau_{\mathrm{H\beta-UV}}$ is the sum of $\tau_{\mathrm{H\beta-opt}}$ and the optical–UV inter-band continuum lag $\tau_{\mathrm{opt-UV}}$, then

[TABLE]

For a standard thin disk, the characteristic size scale of the disk region emitting at wavelength $\lambda$ scales as

[TABLE]

where $l=L/L_{\mathrm{Edd}}$. For a simple photoionization equilibrium model, the BLR characteristic radius scales as

[TABLE]

If the inter-band continuum lags are due to light-travel time across the accretion disk, then

[TABLE]

The lag ratio is thus weakly dependent on both $M_{\mathrm{BH}}$ and accretion rate, where even a factor of $10^{3}$ increase in the black hole mass or accretion rate will change the ratio by only a factor of three. Empirically, Bentz et al. (2013) found a low scatter of 0.13 dex around the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation for 41 AGNs over four orders of magnitude in luminosity, which suggests this effect is indeed small for most AGNs. However, it is important to directly examine the potential consequences of this effect by obtaining more simultaneous observations of the UV/optical continua and the broad emission lines for AGNs over a wide range of luminosities.

5.2. Anomalous Emission-Line Light Curve Behavior

Despite apparent differences in the C IV and H$\beta$ decorrelation start times and flux deficits during T2, it is likely that the cause of the anomalous light-curve behavior is the same for the UV and optical emission lines. The line response during the anomaly is also heavily dependent on the line-of-sight velocity for H$\beta$ (Fig. 11) and the UV emission lines (Goad et al, in prep.). Paper IV suggests two scenarios that could produce the anomaly: (1) a temporary obscuration of the ionizing source from parts of the BLR by a moving veil of gas between the accretion disk and BLR, or (2) a temporary change in spectral energy distribution of the ionizing source. Future papers will investigate in detail the timing and magnitude of the anomaly for all UV and optical lines using the full multi-wavelength dataset from this campaign.

The bottom panel of Figure 7 shows that the H$\beta$ light curve decorrelation was also clearly detected using only the optical data. If our campaign had lasted for only the duration of T2, we still would have been able to measure the H$\beta$ lags with fairly high precision ($\tau_{\mathrm{H\beta-UV}}=5.99^{+0.71}_{-0.75}$ days and $\tau_{\mathrm{H\beta-opt}}=3.10^{+0.77}_{-0.80}$ days), but the lag signal would be contaminated by other unknown factors and would lead to a somewhat biased estimate of the BLR characteristic radius. Depending on how common this decorrelation behavior is for H$\beta$, this effect could contribute to additional scatter in the single-object $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relations for NGC 5548 and other AGNs, which can account for about half of the observed scatter in the global $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation for the entire sample of reverberation-mapped AGNs (Kilerci Eser et al., 2015).

The high cadence and long duration of this campaign have both been crucial in detecting this decorrelation phenomenon. Horne et al. (2004) found that a campaign duration of at least three times the maximum BLR light-crossing time is needed to recover high-fidelity velocity-delay maps from reverberation mapping data. Given the H$\beta$–UV lag of $\sim 6$ days for NGC 5548 during our monitoring period, the BLR characteristic light-crossing time is $\sim 12$ days and the minimum campaign length needed to recover velocity-delay maps would correspond to $\sim 40$ days, which would not have allowed us to see this decorrelation. In order to detect and characterize these anomalous behaviors in the emission-line light curves, RM campaigns must be much longer than the minimum requirement for obtaining velocity-delay maps.

Finally, while there are no previously published results documenting similar emission-line light curve behavior, it is possible that this decorrelation phenomenon was indeed observed in other AGNs in previous RM campaigns, but were not recognized as such because the campaign had relatively low cadence and/or short duration. RM programs designed to study large numbers of sources with lower cadence would also not be able to detect these decorrelations. This further highlights the importance of high-cadence, long-duration and high-SNR multi-wavelength reverberation datasets in order to determine the prevalence of this phenomenon.

5.3. Comparison to Previous Campaigns: the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ Relation

We compare the NGC 5548 H$\beta$–optical lags and optical continuum luminosity ($L_{\mathrm{5100}}$) from this campaign to those from previous RM campaigns targeting this object, as compiled by Kilerci Eser et al. (2015). The total 5100 Å flux and the AGN continuum flux densities for the full, T1, and T2 segments are listed in columns 6 and 7 of Table 9, respectively. We applied a Galactic extinction correction of $E(B-V)=0.017$ mag (Schlegel et al., 1998; Schlafly & Finkbeiner, 2011) and used a luminosity distance of 75 Mpc in converting fluxes to luminosities. Figure 13 shows the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation for NGC 5548, where the uncertainties are from absolute photometric calibration using [O III] $\lambda 5007$ and do not include the luminosity distance uncertainty. The lags from this campaign have smaller uncertainties due to the high cadence and long duration of the ground-based monitoring. Denney et al. (2010) monitored NGC 5548 as part of a multi-object RM campaign and found $\tau_{\mathrm{H\beta-opt}}=12.40^{+2.74}_{-3.85}$ days. However, the light curves were dominated by a large, long-term trend, and after detrending the light curves with a 3rd order polynomial, the H$\beta$ lag becomes $5.07^{+2.46}_{-2.37}$ days (open circle in Fig. 13).

The surprising result is that, given the AGN luminosity during this campaign, the H$\beta$ lags are nearly five times shorter than expected based on past measurements. NGC 5548 was also monitored in 2012 (AGN12, De Rosa et al., in prep.) and had an H$\beta$ lag of $4.1\pm 0.9$ days and a luminosity of $\mathrm{log}(\lambda L_{\mathrm{5100}})=43.22\pm 0.1$, whereas Lu et al. (2016) found the H$\beta$ lag to be $7.2^{+1.33}_{-0.35}$ days and the AGN luminosity to be $\mathrm{log}(\lambda L_{\mathrm{5100}})=43.21\pm 0.1$ during their 2015 campaign, as shown in Figure 13. Both of these lags are also shorter than expected based on past results, though the point from Lu et al. (2016) suggests that the AGN may be returning to its previously measured $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation.

In numerical simulations of the emission-line response to continuum variations for a model BLR with fixed radial extent, Goad & Korista (2014) found that if the characteristic continuum variability timescale ($\tau_{\mathrm{char}}$) is smaller than the BLR light-crossing time, then there exists a strong correlation between $\tau_{\mathrm{char}}$, the line responsivity ($\eta_{\mathrm{eff}}$), and measured lag, such that both $\eta_{\mathrm{eff}}$ and the lag decrease as $\tau_{\mathrm{char}}$ decreases (their Figure 9). This is because short timescale continuum variations only probe the inner parts of the BLR, so the measured lag and responsivity would be biased low. If we use the FWHM of the continuum light curve auto-correlation function as a crude proxy for the characteristic variability time scale, then $\tau_{\mathrm{char}}\sim 10$ days, which is significantly shorter than the values measured for this source for the 1989 IUE campaign and the 1993 HST campaign, and could explain the shorter than expected H$\beta$ lags.

Regardless of the physical causes of the short lags, these results suggest that the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation is more complex than previously realized. Bentz et al. (2013) found a tight correlation between $R_{\mathrm{BLR}}$ and $L_{\mathrm{AGN}}$ for a sample of $\sim 40$ reverberation mapped AGNs, but recent studies by Du et al. (2015) and Du et al. (2016b) have shown that many AGNs with very high accretion rates tend to have considerably shorter H$\beta$ lags compared to low-accretion AGNs with similar luminosities. Now, we have shown that even for a single AGN with low accretion rate ($\lambda_{\mathrm{Edd}}=0.021$ for NGC 5548; Vasudevan et al., 2010), the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation does not always follow a tight power law and that more complex physical processes may contribute significantly to the scatter. In order to further investigate the single-object $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation, repeated monitoring campaigns for individual AGNs are needed to track the behavior of each object over a range of timescales and luminosity states. This will, in turn, help improve our understanding of the global $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation.

6. Summary

We present the results of an optical spectroscopic monitoring program in 2014 targeting the galaxy NGC 5548 as part of the AGN STORM project. Our campaign spanned six months and observed the AGN with an almost daily cadence. Our main findings are as follows.

(1) We determined H$\beta$ and He II $\lambda$4686 emission-line lags relative to the far-UV and optical continua, and found that the lag measured against the UV continuum is $\sim 2$ days longer than that measured against the optical continuum, consistent with the lag between the UV and optical continua. Given that past measurements of the H$\beta$ lag against the optical continuum for this object range from $\sim 4$ to $\sim 25$ days and assuming that this 2-day lag difference is constant over time, then the characteristic size of the BLR inferred from previous data is biased low by 10%–50%. Depending on how the ratio of UV and optical H$\beta$ lags scales with other AGN properties, the RM black hole mass scale and the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation may be affected, which would, in turn, impact single-epoch $M_{\mathrm{BH}}$ estimates for high-redshift AGNs.

(2) We measured velocity-resolved lags for the broad H$\beta$ line and found a double-peaked lag profile as a function of line-of-sight velocity, with shorter lags in the high-velocity wings. The overall shape of the lag profile is qualitatively similar to those of Keplerian models (e.g., Horne et al., 2004), and is very similar to what is found for Ly$\alpha$ (De Rosa et al., 2015).

(3) Both the H$\beta$ and He II $\lambda$4686 emission lines exhibit significant changes in their response to UV continuum variations halfway through our monitoring campaign. The line light curves decorrelate from that of the continuum and remain in a suppressed state until near the end of the campaign. The same anomalous behavior is observed for all the UV emission lines (De Rosa et al., 2015; Goad et al., 2016). Further investigation into the simultaneous UV and optical line responses during this campaign may elucidate the cause of this anomaly. Depending on how frequently this phenomenon occurs in the AGN population as a whole, this effect could contribute to the scatter in both single-object and global $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relations. This type of anomalous line behavior is likely only detectable with monitoring campaigns that have a combination of high cadence, long duration, and high data quality.

(4) Given the optical luminosity of NGC 5548 during our campaign, the H$\beta$ lag measured against the optical continuum is a factor of five shorter than the expected value based on the $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation for NGC 5548 from past monitoring campaigns. Our results, combined with other recent H$\beta$ lag measurements, suggest that this object does not follow a simple power-law $R_{\mathrm{BLR}}-L_{\mathrm{AGN}}$ relation at all times.

We thank the staffs at the various observatories used to obtain the data in this paper. Support for HST program number GO-13330 was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. L.P. and A.J.B. have been supported by National Science Foundation (NSF) grant AST-1412693. M.M.F., G.D.R., B.M.P., C.J.G., and R.W.P. are grateful for the support of NSF grant AST-1008882 to The Ohio State University. M.C. Bentz gratefully acknowledges support through NSF CAREER grant AST-1253702 to Georgia State University. A.V.F.’s group at UC Berkeley is grateful for financial assistance from NSF grant AST-1211916, the TABASGO Foundation, and the Christopher R. Redlich Fund. C.S.K. is supported by NSF grant AST-1515876. V.N.B. gratefully acknowledges assistance from a NSF Research at Undergraduate Institutions (RUI) grant AST-1312296. M.C. Bottorff acknowledges HHMI for support through an undergraduate science education grant to Southwestern University. K.D.D. is supported by an NSF Fellowship awarded under grant AST-1302093. M.E. thanks the members of the Center for Relativistic Astrophysics at Georgia Tech and the Department of Astronomy at the University of Washington, where he was based during the observing campaign, for their warm hospitality. R.E. gratefully acknowledges support from NASA under awards NNX13AC26G, NNX13AC63G, and NNX13AE99G. J.M.G. gratefully acknowledges support from NASA under award NNH13CH61C. P.B.H. is supported by NSERC. K.H. acknowledges support from the UK Science and Technology Facilities Council through grant ST/M001296/1. TW-SH is supported by the DOE Computational Science Graduate Fellowship, grant number DE-FG02- 97ER25308. M.I. acknowledges support from the Creative Initiative program, No. 2008-0060544, of the National Research Foundation of Korea (NRFK) funded by the Korean government (MSIP). M.D.J. acknowledges NSF grant AST-0618209. SRON is financially supported by NWO, the Netherlands Organization for Scientific Research. B.C.K. is partially supported by the UC Center for Galaxy Evolution. C.S.K. acknowledges the support of NSF grant AST-1009756. D.C.L. acknowledges support from NSF grants AST-1009571 and AST-1210311. P.L. acknowledges support from Fondecyt grant No. 1161184. A.P. acknowledges support from a NSF graduate fellowship and a UCSB Dean’s Fellowship. J.S.S. acknowledges CNPq, National Council for Scientific and Technological Development (Brazil) for partial support and The Ohio State University for warm hospitality. N.T. acknowledges support from CONICYT PAI/82140055. T.T. has been supported by NSF grant AST-1412315. T.T. and B.C.K. acknowledge support from the Packard Foundation in the form of a Packard Research Fellowship to T.T. Also, T.T. thanks the American Academy in Rome and the Observatory of Monteporzio Catone for kind hospitality. The Dark Cosmology Centre is funded by the Danish National Research Foundation. M.V. gratefully acknowledges support from the Danish Council for Independent Research via grant no. DFF 4002-00275. J.-H.W. acknowledges support by the National Research Foundation of Korea (NRF) grant funded by the Korean government (No. 2010-0027910). This work is based partly on observations obtained with the Apache Point Observatory 3.5-m telescope, which is owned and operated by the Astrophysical Research Consortium. Research at Lick Observatory is partially supported by a generous gift from Google. This research has made use of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. We thank the anonymous referee for the thorough review and helpful comments that helped to improve the clarity of this manuscript.