The first spectroscopic dust reverberation programme on active galactic nuclei: the torus in NGC 5548

TL;DR

This study presents the first spectroscopic dust reverberation campaign on AGN, revealing detailed dust properties and dynamics in NGC 5548, challenging previous assumptions about dust temperature and sublimation near the AGN.

Contribution

It introduces a novel spectroscopic approach to measure dust properties and response times in AGN, providing new insights into dust composition and the structure of the dusty torus.

Findings

Dust responds to flux changes with a ~70 light-day lag.

Dust temperature is around 1450 K, indicating carbon dominance.

The dusty torus is more like a 'dusty wall' than a thin shell.

Abstract

We have recently initiated the first spectroscopic dust reverberation programme on active galactic nuclei (AGN) in the near-infrared. Spectroscopy enables measurement of dust properties, such as flux, temperature and covering factor, with higher precision than photometry. In particular, it enables measurement of both luminosity-based dust radii and dust response times. Here we report results from a one-year campaign on NGC 5548. The hot dust responds to changes in the irradiating flux with a lag time of ~70 light-days, similar to what was previously found in photometric reverberation campaigns. The mean and rms spectra are similar, implying that the same dust component dominates both the emission and the variations. The dust lag time is consistent with the luminosity-based dust radius only if we assume a wavelength-independent dust emissivity-law, i.e. a blackbody, which is appropriate for grains of large sizes (of a few microns). For such grains the dust temperature is ~1450 K. Therefore, silicate grains have most likely evaporated and carbon is the main chemical component. But the hot dust is not close to its sublimation temperature, contrary to popular belief. This is further supported by our observation of temperature variations largely consistent with a heating/cooling process. Therefore, the inner dust-free region is enlarged and the dusty torus rather a "dusty wall", whose inner radius is expected to be luminosity-invariant. The dust-destruction mechanism that enlarges the dust-free region seems to partly affect also the dusty region. We observe a cyclical decrease in dust mass with implied dust reformation times of ~5-6 months.

Tables

Table 1. Table 1: IRTF Journal of observations

Observation	exposure	airmass	aperture	PA	continuum S/N			telluric standard star		seeing	cloud	correction
Date	(s)		(arcsec2)	(∘)	$J$	$H$	$K$	name	airmass	(arcsec)	condition	factor
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)	(12)	(13)
2016 Aug 2	16$\times$120	1.180	0.3$\times$3.2	96	80	141	220	HD 131951	1.382	0.7	clear	$0.85\pm 0.02$
2016 Aug 11	16$\times$120	1.209	0.3$\times$3.6	92	66	118	187	HD 131951	1.425	0.6	photom.	$1.03\pm 0.03$
2016 Aug 21	5$\times$120	1.544	0.3$\times$4.0	85	25	48	80	HD 131951	1.962	0.6	clear	$0.71\pm 0.06$
2016 Dec 21	16$\times$120	2.184	0.3$\times$4.2	281	45	82	133	HD 131951	1.765	0.4	photom.	$0.78\pm 0.03$
2017 Jan 6	16$\times$120	1.355	0.3$\times$3.6	274	46	84	144	HD 131951	1.178	0.9	clear	$1.02\pm 0.03$
2017 Jan 20	16$\times$120	1.216	0.3$\times$4.2	269	48	84	142	HD 131951	1.109	0.6	cirrus	$0.81\pm 0.03$
2017 Jan 24	21$\times$120	1.240	0.3$\times$5.0	270	39	72	125	HD 131951	1.078	1.2	clear	$0.82\pm 0.04$
2017 Feb 5	16$\times$120	1.064	0.3$\times$3.4	260	60	108	187	HD 131951	1.015	0.7	photom.	$0.67\pm 0.02$
2017 Feb 15	12$\times$120	1.006	0.3$\times$4.0	248	45	83	142	HD 131951	1.006	0.8	cirrus	$0.73\pm 0.02$
2017 Feb 24	16$\times$120	1.005	0.3$\times$3.2	214	64	109	183	HD 131951	1.022	0.6	cirrus	$1.15\pm 0.04$
2017 Mar 17	16$\times$120	1.027	0.3$\times$3.6	250	48	93	165	HD 131951	1.005	0.7	photom.	$1.00$
2017 Mar 22	18$\times$120	1.030	0.3$\times$3.4	254	71	133	224	HD 131951	1.005	0.5	photom.	$0.77\pm 0.02$
2017 May 9	16$\times$120	1.041	0.3$\times$3.6	256	52	94	157	HD 131951	1.009	0.6	clear	$0.89\pm 0.03$
2017 May 25	16$\times$120	1.042	0.3$\times$3.2	256	55	99	171	HD 131951	1.007	0.4	photom.	$1.11\pm 0.04$
2017 Jun 10	15$\times$120	1.055	0.3$\times$3.4	258	52	95	162	HD 121996	1.054	0.8	photom.	$0.95\pm 0.04$
2017 Jun 20	22$\times$120	1.049	0.3$\times$3.8	257	54	101	173	HD 121996	1.043	0.7	cirrus	$0.80\pm 0.03$
2017 Jun 28	24$\times$120	1.027	0.3$\times$3.4	251	70	125	211	HD 121996	1.025	0.6	clear	$0.78\pm 0.03$
2017 Jul 3	20$\times$120	1.005	0.3$\times$3.4	178	80	136	217	HD 121996	1.001	0.5	photom.	$0.98\pm 0.04$

Table 2. Table 2: Optical and near-IR nuclear and host galaxy fluxes from GROND photometry

Observation	$g^{\prime }$			$r^{\prime }$
Date	${\lambda }_{\mathrm{eff}}=4587$ Å			${\lambda }_{\mathrm{eff}}=6220$ Å
	nuclear flux	PSF	host galaxy flux	nuclear flux	PSF	host galaxy flux
	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)
2017 Mar 11	(6.73$\pm$0.18)e$-$15	1.12	(2.04$\pm$0.07)e$-$14	(4.79$\pm$0.06)e$-$15	1.09	(2.23$\pm$0.04)e$-$14
2017 Mar 18	(5.73$\pm$0.15)e$-$15	1.23	(2.15$\pm$0.07)e$-$14	(4.91$\pm$0.06)e$-$15	1.19	(2.07$\pm$0.04)e$-$14
Observation	$i^{\prime }$			$z^{\prime }$
Date	${\lambda }_{\mathrm{eff}}=7641$ Å			${\lambda }_{\mathrm{eff}}=8999$ Å
	nuclear flux	PSF	host galaxy flux	nuclear flux	PSF	host galaxy flux
	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)
2017 Mar 11	(3.01$\pm$0.01)e$-$15	0.92	(2.70$\pm$0.01)e$-$14	(2.96$\pm$0.07)e$-$15	0.98	(2.85$\pm$0.07)e$-$14
2017 Mar 18	(2.55$\pm$0.01)e$-$15	1.06	(2.71$\pm$0.01)e$-$14	(2.24$\pm$0.05)e$-$15	0.98	(2.16$\pm$0.05)e$-$14
Observation	$J$			$H$
Date	${\lambda }_{\mathrm{eff}}=12399$ Å			${\lambda }_{\mathrm{eff}}=16468$ Å
	nuclear flux	PSF	host galaxy flux	nuclear flux	PSF	host galaxy flux
	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)
2017 Mar 11	(2.30$\pm$0.01)e$-$15	1.45	(1.57$\pm$0.01)e$-$14	(2.12$\pm$0.02)e$-$15	1.55	(9.72$\pm$0.08)e$-$15
2017 Mar 18	(2.17$\pm$0.01)e$-$15	1.58	(1.57$\pm$0.01)e$-$14	(1.95$\pm$0.02)e$-$15	1.87	(1.14$\pm$0.01)e$-$14
Observation	$K$
Date	${\lambda }_{\mathrm{eff}}=21706$ Å
	nuclear flux	PSF	host galaxy flux
	(erg/s/cm2/Å)	(arcsec)	(erg/s/cm2/Å)
2017 Mar 11	(1.97$\pm$0.06)e$-$15	1.71	(3.04$\pm$0.16)e$-$15
2017 Mar 18	(1.93$\pm$0.06)e$-$15	1.99	(2.81$\pm$0.13)e$-$15

Table 3. Table 3: Physical parameters for the calculation of luminosity-weighted dust radii

Observation	accretion	blackbody			silicate dust			carbon dust
Date	disc	($\beta =0$)			($\beta =-1$)			($\beta =-2$)
	log $L_{\mathrm{uv}}$	$T_d$	log $L_d$	$R_{d,\mathrm{lum}}$	$T_d$	log $L_d$	$R_{d,\mathrm{lum}}$	$T_d$	log $L_d$	$R_{d,\mathrm{lum}}$
	(erg/s)	(K)	(erg/s)	(lt-days)	(K)	(erg/s)	(lt-days)	(K)	(erg/s)	(lt-days)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)	(10)	(11)
2016 Aug 2	44.44	1432$\pm$12	44.11	59	1212$\pm$9	43.12	564	1052$\pm$6	44.04	367
2016 Aug 11	44.44	1467$\pm$14	44.10	56	1240$\pm$10	43.23	539	1075$\pm$7	44.00	351
2016 Aug 21	44.55	1458$\pm$29	44.15	64	1231$\pm$21	43.14	620	1067$\pm$15	44.02	404
2016 Dec 21	44.46	1444$\pm$20	43.97	59	1219$\pm$14	43.10	570	1056$\pm$10	43.76	372
2017 Jan 6	44.43	1430$\pm$16	44.00	58	1207$\pm$12	43.09	562	1046$\pm$9	43.93	367
2017 Jan 20	44.49	1432$\pm$20	43.97	62	1211$\pm$14	43.09	598	1050$\pm$11	43.94	390
2017 Jan 24	44.50	1466$\pm$25	43.93	60	1232$\pm$17	43.07	585	1063$\pm$13	43.97	385
2017 Feb 5	44.45	1390$\pm$13	44.03	63	1180$\pm$10	43.10	602	1026$\pm$8	43.61	390
2017 Feb 15	44.41	1411$\pm$17	43.98	58	1196$\pm$12	43.08	559	1039$\pm$9	43.68	363
2017 Feb 24	44.48	1396$\pm$15	43.95	64	1188$\pm$11	43.04	615	1035$\pm$8	43.67	397
2017 Mar 17	44.36	1478$\pm$17	43.98	50	1245$\pm$12	43.10	487	1076$\pm$9	43.58	320
2017 Mar 22	44.36	1453$\pm$14	44.03	52	1230$\pm$9	43.13	499	1067$\pm$7	43.91	325
2017 May 9	44.54	1471$\pm$18	43.96	62	1238$\pm$13	43.04	606	1070$\pm$9	43.91	398
2017 May 25	44.60	1439$\pm$14	44.04	70	1215$\pm$10	43.10	675	1053$\pm$7	43.70	440
2017 Jun 10	44.60	1515$\pm$19	43.99	63	1270$\pm$13	43.13	617	1094$\pm$10	43.98	408
2017 Jun 20	44.54	1477$\pm$17	44.02	62	1245$\pm$11	43.11	600	1077$\pm$8	43.94	392
2017 Jun 28	44.63	1474$\pm$13	44.02	69	1243$\pm$9	43.12	667	1076$\pm$7	43.60	436
2017 Jul 3	44.63	1474$\pm$13	44.08	69	1248$\pm$9	43.17	662	1083$\pm$7	43.62	430

Table 4. Table 4: Near-IR spectral light-curves corrsponding to the accretion disc

Observation	HJD	8700$-$9000 Å		9730$-$9790 Å
Date	-2,400,000	original	corrected⋆	original	corrected⋆
		(erg/s/cm2/Å)	(erg/s/cm2/Å)	(erg/s/cm2/Å)	(erg/s/cm2/Å)
2016 Aug 2	57602.76	(2.87$\pm$0.14)e$-$15	(2.45$\pm$0.13)e$-$15	(2.49$\pm$0.07)e$-$15	(2.12$\pm$0.08)e$-$15
2016 Aug 11	57611.74	(2.42$\pm$0.11)e$-$15	(2.50$\pm$0.13)e$-$15	(2.10$\pm$0.08)e$-$15	(2.16$\pm$0.10)e$-$15
2016 Aug 21	57621.76	(3.96$\pm$0.48)e$-$15	(2.82$\pm$0.40)e$-$15	(3.97$\pm$0.45)e$-$15	(2.83$\pm$0.39)e$-$15
2016 Dec 21	57744.08	(3.24$\pm$0.19)e$-$15	(2.54$\pm$0.17)e$-$15	(2.78$\pm$0.10)e$-$15	(2.18$\pm$0.11)e$-$15
2017 Jan 6	57760.10	(2.44$\pm$0.15)e$-$15	(2.49$\pm$0.17)e$-$15	(2.02$\pm$0.11)e$-$15	(2.06$\pm$0.13)e$-$15
2017 Jan 20	57774.09	(3.19$\pm$0.21)e$-$15	(2.60$\pm$0.19)e$-$15	(2.84$\pm$0.11)e$-$15	(2.31$\pm$0.11)e$-$15
2017 Jan 24	57778.08	(3.31$\pm$0.26)e$-$15	(2.71$\pm$0.25)e$-$15	(2.82$\pm$0.16)e$-$15	(2.31$\pm$0.18)e$-$15
2017 Feb 5	57790.09	(3.96$\pm$0.20)e$-$15	(2.64$\pm$0.15)e$-$15	(3.17$\pm$0.09)e$-$15	(2.11$\pm$0.08)e$-$15
2017 Feb 15	57800.11	(3.27$\pm$0.20)e$-$15	(2.40$\pm$0.17)e$-$15	(2.77$\pm$0.11)e$-$15	(2.03$\pm$0.10)e$-$15
2017 Feb 24	57809.10	(2.32$\pm$0.10)e$-$15	(2.66$\pm$0.14)e$-$15	(1.90$\pm$0.05)e$-$15	(2.18$\pm$0.09)e$-$15
2017 Mar 17	57830.01	(2.28$\pm$0.13)e$-$15	(2.28$\pm$0.13)e$-$15	(1.87$\pm$0.08)e$-$15	(1.87$\pm$0.08)e$-$15
2017 Mar 22	57834.99	(2.81$\pm$0.11)e$-$15	(2.17$\pm$0.10)e$-$15	(2.43$\pm$0.07)e$-$15	(1.88$\pm$0.07)e$-$15
2017 May 9	57882.85	(3.24$\pm$0.17)e$-$15	(2.89$\pm$0.18)e$-$15	(2.70$\pm$0.09)e$-$15	(2.41$\pm$0.11)e$-$15
2017 May 25	57898.81	(2.96$\pm$0.18)e$-$15	(3.28$\pm$0.23)e$-$15	(2.40$\pm$0.08)e$-$15	(2.66$\pm$0.13)e$-$15
2017 Jun 10	57914.76	(3.51$\pm$0.20)e$-$15	(3.32$\pm$0.22)e$-$15	(2.77$\pm$0.08)e$-$15	(2.62$\pm$0.12)e$-$15
2017 Jun 20	57924.73	(3.56$\pm$0.24)e$-$15	(2.86$\pm$0.22)e$-$15	(3.08$\pm$0.11)e$-$15	(2.48$\pm$0.13)e$-$15
2017 Jun 28	57932.72	(4.32$\pm$0.18)e$-$15	(3.36$\pm$0.18)e$-$15	(3.50$\pm$0.09)e$-$15	(2.73$\pm$0.12)e$-$15
2017 Jul 3	57937.75	(3.47$\pm$0.15)e$-$15	(3.39$\pm$0.20)e$-$15	(2.91$\pm$0.05)e$-$15	(2.84$\pm$0.13)e$-$15

Table 5. Table 5: Near-IR spectral light-curves corresponding to the hot dust

Observation	HJD	1.55$-$1.60 $\mu$m			2.34$-$2.39 $\mu$m
Date	-2,400,000	original	corrected⋆	corrected	original	corrected⋆	corrected
		(erg/s/cm2/Å)	(erg/s/cm2/Å)	(subtracted)†	(erg/s/cm2/Å)	(erg/s/cm2/Å)	(subtracted)†
				(erg/s/cm2/Å)			(erg/s/cm2/Å)
2016 Aug 2	57602.76	(2.76$\pm$0.04)e$-$15	(2.35$\pm$0.07)e$-$15	(1.66$\pm$0.06)e$-$15	(2.46$\pm$0.04)e$-$15	(2.10$\pm$0.07)e$-$15	(1.83$\pm$0.06)e$-$15
2016 Aug 11	57611.74	(2.32$\pm$0.03)e$-$15	(2.39$\pm$0.07)e$-$15	(1.70$\pm$0.06)e$-$15	(2.03$\pm$0.03)e$-$15	(2.09$\pm$0.07)e$-$15	(1.82$\pm$0.06)e$-$15
2016 Aug 21	57621.76	(3.81$\pm$0.12)e$-$15	(2.71$\pm$0.23)e$-$15	(1.88$\pm$0.17)e$-$15	(3.32$\pm$0.14)e$-$15	(2.36$\pm$0.21)e$-$15	(2.03$\pm$0.19)e$-$15
2016 Dec 21	57744.08	(2.49$\pm$0.04)e$-$15	(1.96$\pm$0.07)e$-$15	(1.24$\pm$0.05)e$-$15	(2.06$\pm$0.04)e$-$15	(1.61$\pm$0.06)e$-$15	(1.33$\pm$0.06)e$-$15
2017 Jan 6	57760.10	(1.92$\pm$0.04)e$-$15	(1.96$\pm$0.07)e$-$15	(1.28$\pm$0.06)e$-$15	(1.68$\pm$0.02)e$-$15	(1.72$\pm$0.06)e$-$15	(1.45$\pm$0.05)e$-$15
2017 Jan 20	57774.09	(2.41$\pm$0.04)e$-$15	(1.96$\pm$0.07)e$-$15	(1.22$\pm$0.05)e$-$15	(2.03$\pm$0.03)e$-$15	(1.65$\pm$0.06)e$-$15	(1.36$\pm$0.05)e$-$15
2017 Jan 24	57778.08	(2.35$\pm$0.05)e$-$15	(1.93$\pm$0.10)e$-$15	(1.17$\pm$0.07)e$-$15	(1.84$\pm$0.03)e$-$15	(1.51$\pm$0.08)e$-$15	(1.21$\pm$0.07)e$-$15
2017 Feb 5	57790.09	(2.95$\pm$0.05)e$-$15	(1.96$\pm$0.06)e$-$15	(1.26$\pm$0.05)e$-$15	(2.73$\pm$0.03)e$-$15	(1.82$\pm$0.05)e$-$15	(1.54$\pm$0.05)e$-$15
2017 Feb 15	57800.11	(2.52$\pm$0.05)e$-$15	(1.85$\pm$0.07)e$-$15	(1.19$\pm$0.05)e$-$15	(2.23$\pm$0.03)e$-$15	(1.64$\pm$0.06)e$-$15	(1.37$\pm$0.05)e$-$15
2017 Feb 24	57809.10	(1.56$\pm$0.03)e$-$15	(1.80$\pm$0.06)e$-$15	(1.06$\pm$0.04)e$-$15	(1.37$\pm$0.02)e$-$15	(1.57$\pm$0.06)e$-$15	(1.28$\pm$0.05)e$-$15
2017 Mar 17	57830.01	(1.93$\pm$0.03)e$-$15	(1.93$\pm$0.07)e$-$15	(1.32$\pm$0.05)e$-$15	(1.63$\pm$0.02)e$-$15	(1.63$\pm$0.05)e$-$15	(1.40$\pm$0.05)e$-$15
2017 Mar 22	57834.99	(2.60$\pm$0.04)e$-$15	(2.01$\pm$0.06)e$-$15	(1.40$\pm$0.05)e$-$15	(2.30$\pm$0.04)e$-$15	(1.78$\pm$0.05)e$-$15	(1.54$\pm$0.05)e$-$15
2017 May 9	57882.85	(2.33$\pm$0.04)e$-$15	(2.08$\pm$0.07)e$-$15	(1.26$\pm$0.05)e$-$15	(1.86$\pm$0.03)e$-$15	(1.66$\pm$0.06)e$-$15	(1.34$\pm$0.05)e$-$15
2017 May 25	57898.81	(2.09$\pm$0.03)e$-$15	(2.32$\pm$0.09)e$-$15	(1.43$\pm$0.06)e$-$15	(1.77$\pm$0.03)e$-$15	(1.96$\pm$0.07)e$-$15	(1.61$\pm$0.06)e$-$15
2017 Jun 10	57914.76	(2.47$\pm$0.04)e$-$15	(2.34$\pm$0.10)e$-$15	(1.44$\pm$0.07)e$-$15	(1.89$\pm$0.03)e$-$15	(1.79$\pm$0.07)e$-$15	(1.44$\pm$0.06)e$-$15
2017 Jun 20	57924.73	(2.83$\pm$0.04)e$-$15	(2.28$\pm$0.09)e$-$15	(1.46$\pm$0.07)e$-$15	(2.29$\pm$0.05)e$-$15	(1.84$\pm$0.08)e$-$15	(1.52$\pm$0.07)e$-$15
2017 Jun 28	57932.72	(3.05$\pm$0.03)e$-$15	(2.37$\pm$0.09)e$-$15	(1.44$\pm$0.06)e$-$15	(2.43$\pm$0.06)e$-$15	(1.89$\pm$0.08)e$-$15	(1.53$\pm$0.07)e$-$15
2017 Jul 3	57937.75	(2.63$\pm$0.03)e$-$15	(2.57$\pm$0.11)e$-$15	(1.64$\pm$0.08)e$-$15	(2.18$\pm$0.14)e$-$15	(2.13$\pm$0.16)e$-$15	(1.77$\pm$0.15)e$-$15

Table 6. Table 6: Results from MCMCRev  assuming the response function of a standard accretion disc

Wavelength	band	${\tau }_0$	$⟨\tau ⟩$	${\chi }_{\nu }^2$
(Å)		(lt-days)	(lt-days)
23600	$K$ (spectrum)	12$\pm$2	55$\pm$9	1.59
15800	$H$ (spectrum)	7$\pm$1	35$\pm$7	1.60
8850	$z$ (spectrum)	3$\pm$6	17$\pm$3	0.80
8700	$z_s$ (photometry)	3$\pm$6	16$\pm$3	1.41
5448	$V$ (photometry)	2$\pm$3	10$\pm$2	1.14

Equations

$$F_{\mathrm{d}}(\lambda )=\pi B_{\lambda }(T)Q_{\lambda }(a),$$

$$\frac{\pi a^2}{4\pi R_{\mathrm{d},\mathrm{lum}}^2}\int L_{\lambda }{Q}_{\lambda }^{\mathrm{abs}}d\lambda =4\pi a^2\int F_{\mathrm{d}}(\lambda )d\lambda ,$$

$$\frac{L_{\mathrm{uv}}}{4\pi R_{\mathrm{d},\mathrm{lum}}^2}=4\sigma T^4⟨Q^{\mathrm{em}}⟩,$$

$$F_{\mathrm{uv}}(t)$$

$$F_{\nu }(t)$$

$$⟨R_{\mathrm{d},\mathrm{rev}}⟩=c⟨\tau ⟩=\frac{\int c\tau \Psi (\tau )d\tau }{\int \Psi (\tau )d\tau }$$

Full text

††pagerange: The first spectroscopic dust reverberation programme on active galactic nuclei: the torus in NGC 5548–References††pubyear: 2019

1 Introduction

Some active galactic nuclei (AGN) exhibit in their optical spectra both broad ($\sim 1-5\%$ the speed of light) and narrow emission lines (type 1 AGN) and some only narrow lines (type 2 AGN). The existence of these two apparently distinct classes of AGN has been explained within unified schemes by orientation effects: a geometrically and optically thick, dusty torus (or warped disc) located outside the accretion disc obscures the broad emission line region (BLR) for some lines of sight (see reviews by Lawrence, 1987; Antonucci, 1993; Urry & Padovani, 1995; Netzer, 2015). The strongest evidence that such a torus exists comes from observations of broad emission lines in the polarized, scattered light of numerous type 2 AGN (e.g., Antonucci & Miller, 1985; Cohen et al., 1999; Lumsden et al., 2001; Tran, 2003) and of strong infrared (IR) continuum emission at wavelengths $\lambda>1~{}\mu$m in most AGN. The latter is expected since the dusty torus will absorb a significant fraction of the accretion disc radiation and reradiate it at longer wavelengths. However, the precise location and geometry of the dusty torus remain unknown, mainly because it cannot be spatially resolved. Mid-IR imaging at high angular resolution (Packham et al., 2005; Ramos Almeida et al., 2011), near-IR and mid-IR interferometric observations (Tristram et al., 2009; Pott et al., 2010; Kishimoto et al., 2011; Burtscher et al., 2013; Kishimoto et al., 2013) and sub-mm interferometric observations (García-Burillo et al., 2016; Imanishi et al., 2018) of a few of the brightest AGN have so far imposed mainly upper limits on its extent (of $\sim 0.1-6$ pc).

The most promising technique to constrain the size of the dusty torus when it cannot be spatially resolved is reverberation, i.e., measuring the time with which the dust responds to changes in the irradiating flux. For $\sim 20$ AGN, the radius of the innermost part of the torus, where dust grains are expected to be close to their sublimation temperature, has been successfully determined through coordinated, long-term optical and near-IR photometric campaigns (Clavel et al., 1989; Sitko et al., 1993; Nelson, 1996; Oknyanskij & Horne, 2001; Glass, 2004; Minezaki et al., 2004; Suganuma et al., 2006; Koshida et al., 2014; Vazquez et al., 2015). Interestingly, the hot dust radius is variable in most AGN that have multiple determinations. When photometry can be performed simultaneously at several near-IR wavelengths, a crude estimate of the dust temperature is also available. Its variability behaviour in addition to that of the dust lag times can reveal the mechanism by which dust is destroyed and reforms and can also constrain the dust grain type and size (Barvainis, 1992; Koshida et al., 2009; Schnülle et al., 2013, 2015). Dust lag times are found to be systematically smaller than dust radii measured with near-IR interferometry or estimated from the AGN bolometric luminosity and the dust sublimation temperature (Oknyanskij & Horne, 2001; Kishimoto et al., 2007; Nenkova et al., 2008). Kawaguchi & Mori (2010, 2011) ascribe this to a dust geometry that is bowl-shaped rather than spherical due to the anisotropy of the accretion disc emission. Such a geometry would place the hottest dust further in than the bulk of the emission, thus reducing the observed dust reverberation times relative to the dust radius estimated from the continuum radiation.

There are two approaches that have been taken for the dynamics of the dusty torus. The hydrostatic scenario depicts the torus as a structure populated by clouds whose motion is dominated by the gravitational potential of the central super-massive black hole (Krolik & Begelman, 1988; Pier & Krolik, 1992; Krolik, 2007; Nenkova et al., 2008). The vertical motion needed to sustain the clouds in a hydostatic structure with a height-to-radius ratio of roughly unity could be provided by radiation pressure. The other scenario advocates the outflow of clouds embedded in a hydromagnetic disc wind (Königl & Kartje, 1994; Elitzur & Shlosman, 2006; Keating et al., 2012). In this approach, the torus is that region in the wind wherein the clouds are dusty.

With the advent of efficient near-IR cross-dispersed spectrographs (mounted in particular at 8 m class telescopes) spectroscopic and not only photometric dust reverberation campaigns are now feasible in reasonable amounts of observing time. We have recently initiated the first such programme. The advantages of spectroscopy over photometry are enormous: (i) a large portion of the hot dust spectral energy distribution (SED) can be monitored at once and so the temporal evolution of the dust temperature determined with high accuracy; (ii) the contribution of the accretion disc to the near-IR continuum flux can be reliably separated from that of the hot dust, with the former being a large source of uncertainty in photometry campaigns; (iii) with accurate dust temperatures and accretion disc luminosities in hand, luminosity-based dust radii can be determined simultaneously with the dust lag times and temporal changes in dust covering factor studied; and (iv) emission lines believed to be associated with the dusty torus can also be sampled. For example, the coronal (i.e. forbidden high-ionisation) emission line region, which is thought to lie betweeen the broad and narrow emission line region (NELR), although itself dust-free (Ferguson et al., 1997), is assumed to mark the onset of the dusty torus. Determining how the line flux and profile shape vary in relation to the variable hot dust emission can explain how this region is related to the torus, e.g., through an X-ray heated wind launched from it (Pier & Voit, 1995). The hot dust presumably also influences the BLR by limiting its extent (Netzer & Laor, 1993; Czerny & Hryniewicz, 2011; Mor & Netzer, 2012). Evidence for a limited BLR is provided by flat-topped (Landt et al., 2014) or double-peaked (Storchi-Bergmann et al., 2017) line profiles as well as velocity-delay maps (Horne et al., in preparation). Then, if the BLR is indeed dust- rather than radiation-bounded, we expect the widths of the profile top to increase as the hot dust radius decreases. In this case, the BLR might have a bowl-shaped geometry like the dust (Goad et al., 2012).

Here we present our results from a one-year long spectroscopic dust reverberation campaign on the nearby, well-known AGN NGC 5548 conducted between 2016 August and 2017 July. In this study, we concentrate on the variability of the hot dust and will discuss the variability of the coronal lines and broad emission lines elsewhere. The paper is structured as follows. After we briefly discuss our science target in Section 2, we give details of the observations, data reduction and measurements in Section 3. In Sections 4, we derive response-weighted dust radii. In Section 5, we constrain the location, dynamics and geometry of the hot dust torus. Finally, in Section 6, we present a short summary and our conclusions.

2 The science target

We have selected the type-1 AGN NGC 5548 as our science target for several reasons: (i) it has a short and variable dust response time ($\sim 40-80$ days; Koshida et al., 2014), which means that in one year we can measure several dust lags and so can investigate why they vary; (ii) it has strong and variable near-IR coronal lines (Landt et al., 2015); (iii) its broad emission line profiles have clearly discernible and therefore easily separable broad- and narrow-line components, and (iv) multiple optical reverberation campaigns have determined its black hole mass ($M_{\rm BH}=5\times 10^{7}~{}M_{\odot}$; Peterson et al., 2004; Bentz et al., 2009; Bentz & Katz, 2015). Furthermore, two recent ambitious programs, namely, the X-ray spectroscopic monitoring campaign in summer 2013 for $\sim 1.5$ months with XMM-Newton (Kaastra et al., 2014; Mehdipour et al., 2015) and the ultraviolet (UV) spectroscopic monitoring campaign with the Hubble Space Telescope (HST) during the first half of 2014 (De Rosa et al., 2015; Edelson et al., 2015; Fausnaugh et al., 2016; Goad et al., 2016; Pei et al., 2017; Starkey et al., 2017; Mathur et al., 2017), both supported by extensive space- and ground-based multiwavelength observations, make it one of the best-documented AGN yet.

NGC 5548 (J2000 sky coordinates R.A. $14^{h}17^{m}59.5^{s}$, Decl. $+25^{\circ}08^{\prime}12^{\prime\prime}$) is observable from both the northern and southern hemispheres. It is at a low redshift ($z=0.01718$), which means that a cross-dispersed near-IR spectrum samples a large portion of the hot dust component and several strong near-IR coronal lines. The AGN, which is hosted by a spiral galaxy, is relatively bright in the near-IR (2MASS $J=11.8$ mag, $K_{s}=10.1$ mag; Skrutskie et al., 2006), requiring only modest exposure times on a 4 m class telescope. We adopt here cosmological parameters $H_{0}=70$ km s*-1* Mpc*-1*, $\Omega_{\rm M}=0.3$, and $\Omega_{\Lambda}=0.7$, which give a luminosity distance to NGC 5548 of 74.6 Mpc and an angular scale at the galaxy of 349 pc per arcsec.

3 The observations

3.1 The near-IR spectroscopy

We observed the source NGC 5548 between 2016 August and 2017 July (semesters 2016B and 2017A) with the recently refurbished SpeX spectrograph (Rayner et al., 2003) at the NASA Infrared Telescope Facility (IRTF), a 3 m telescope on Maunakea, Hawaii. Our approved proposal requested a cadence of about a week during the period of a year starting 2016 August 1. Excluding the 3.5-month period (Sep to mid-Dec) when the source is unobservable. We were scheduled 24 observing windows of 2-2.5 hours each between 2016 August and 2017 July with an observing gap of about a month (mostly during April) due to engineering time for the spectrograph. Of the scheduled observing windows, we lost six due to weather and engineering issues, resulting in a total of 18 near-IR spectra with an average cadence of about ten days. We list the journal of observations in Table 1.

We used the short cross-dispersed mode (SXD, $0.7-2.55~{}\mu$m) equipped with the $0.3^{\prime\prime}\times 15^{\prime\prime}$ slit, which we oriented at the parallactic angle. This set-up gives an average spectral resolution of $R=2000$ or full width at half maximum (FWHM) $\sim 150$ km s*-1*, which is sufficiently high to study line profiles and to clearly discern the narrow and broad line components for the permitted transitions. The narrow slit also minimizes the flux contamination from the host galaxy. The on-source exposure time was usually $16\times 120$ s, however, when there was time available we obtained additional frames. We note that the exposure time for the run on 2016 August 21 was only $5\times 120$ s, since we lost most data due to telescope problems. The chosen exposure time ensured that we obtained spectra with a high signal-to-noise ($S/N$) ratio in order to reliably measure emission line profiles. The observations were done in the usual ABBA nodding pattern. Since the source is extended in the near-IR, we nodded off onto a blank patch of sky for the background subtraction.

After the science target, we observed the nearby (in position and airmass) A0 V star HD 131951 that has accurate optical magnitudes ($B=5.870$ mag, $V=5.901$ mag). We used this standard star to correct our science spectrum for telluric absorption and for flux calibration. For the observations in 2017 June and July we used the A0 V star HD 121996, which is of a similar brightness ($B=5.743$ mag, $V=5.748$ mag). Flats and arcs were taken mostly after the science target.

We reduced the data using Spextool (version 4.1), an Interactive Data Language (IDL)-based software package developed for SpeX users (Cushing et al., 2004). Spextool carries out all the procedures necessary to produce fully reduced spectra. This includes preparation of calibration frames, processing and extraction of spectra from science frames, wavelength calibration of spectra, telluric correction and flux-calibration of spectra, and merging of the different orders into a single, continuous spectrum. We used for both the science target and the telluric standard star an optimally weighted extraction (Horne, 1986), which requires the fitting and subtracting of a local background. Based on the Galactic hydrogen column densities published by Dickey & Lockman (1990), the Galactic extinction towards the source NGC 5548 is negligible and we did not correct for it. In Fig. 1, we show the spectrum from our run on 2017 February 5 as a representative example.

3.2 The complementary photometry

3.2.1 The GROND observations

We have obtained simultaneous optical and near-IR photometry in seven bands with the Gamma-Ray Burst Optical and Near-Infrared Detector (GROND; Greiner et al., 2008) mounted on the MPG 2.2 m telescope, La Silla, Chile. The observations were taken on 2017 March 11 and 18 with the latter date only one day after our near-IR spectroscopy observation from 2017 March 17. The optical channels have a field-of-view (FOV) of $5.4^{\prime}\times 5.4^{\prime}$ and a plate scale of $0.158^{\prime\prime}$ per pixel, whereas both the FOV and plate scale of the near-IR channels is larger ($10^{\prime}\times 10^{\prime}$ and $0.6^{\prime\prime}$ per pixel, respectively). The filter passbands are relatively wide ($\sim 1000-3000$ Å). The seeing on both nights was $\sim 0.8^{\prime\prime}-0.9^{\prime\prime}$, with clear skies and thin clouds on 2017 March 11 and 18, respectively.

Images were reduced with the IRAF-based pipelines devoloped by R. Decarli and G. De Rosa (Morganson et al., 2012), modified and extended for our purpose. After bias and dark subtraction and flat-fielding of each individual image, a scaled median sky image constructed from separate sky exposures was subtracted from each frame. The dither step size was $18^{\prime\prime}$ for both the science and sky images, with five dither positions obtained for each. After sky subtraction, the individual frames were realigned and combined into the final image. In Table 2, we list the optical and near-IR nuclear fluxes with their $1~{}\sigma$ errors, which were obtained with GALFIT (Peng et al., 2002), a software package that models the object’s surface brightness profile with a point spread function (PSF) and a host galaxy component. We modelled the host galaxy with only a bulge, since the $S/N$ of the images is not high enough to warrant including also a disc component. The surface brightness profile of the bulge was approximated with a Sérsic profile with Sérsic index $n=3.0$ and effective radius $R_{\rm e}=12.6^{\prime\prime}$. In addition, GALFIT simultaneously models the PSF fluxes of the reference stars, which we used to flux-calibrate the optical and near-IR images based on the known stellar fluxes as listed in the SDSS (Abazajian et al., 2009) and 2MASS (Skrutskie et al., 2006).

3.2.2 The LCOGT observations

We monitored NGC 5548 with the 1 m robotic telescope network of the Las Cumbres Observatory (LCOGT; Brown et al., 2013) as part of the 2014 AGN Key Project almost daily in $V$ between 2016 July 16 and August 8 and in six bands ($V$, $u$, $g$, $r$, $i$ and $z_{s}$) between 2016 December 12 and 2017 September 9. We make extensive use here of the light-curves in the $V$ and $z_{s}$ filters, which have a wavelength width of 840 Å and 1040 Å around their central wavelength of 5448 Å and 8700 Å, respectively. The Sinistro cameras used for imaging have a FOV of $26.5^{\prime}\times 26.5^{\prime}$ and a plate scale of $0.389^{\prime\prime}$ per pixel. The frames were first processed by LCOGT’s BANZAI pipeline in the usual way (bias and dark subtraction, flat-fielding and image correction) and were subsequently analysed with AstroImageJ (Collins et al., 2017). The total flux of the science target was extracted in a circular aperture with a radius of 15 pixels, from which we subtracted a background flux determined from a circular region with a minimum and maximum radius of 25 pixels and 35 pixels, respectively. Comparison stars were selected from the same FOV, whereby we omitted bright, saturated stars and set the number of stars such that the sum of their flux is $\sim 8-10$ times the flux of the science target. In Fig. 2 (top panels), we show the $i$ and $z_{s}$ light-curves, which overlap in wavelength with the IRTF spectrum.

3.3 The absolute spectral flux scale

In order to obtain meaningful light-curves and estimates of phyiscal parameters such as, e.g. luminosities, it is crucial to know the absolute flux scale of the spectra as accurately as possible. In optical spectroscopic reverberation studies, this absolute spectrophotometry is usually obtained by using the flux of the strong, narrow emission line [O III] $\lambda 5007$, which is typically assumed to stay constant during the campaign and to have a value as measured in the best-weather spectra. The assumption of constancy is usually justified for forbidden line transitions, since the emitting gas is located at pc-to-kpc distances from the ionising source and has a number density low enough for recombination time-scales to be large. Then, in order to ensure that the same amount of narrow line flux is enclosed in all spectra, a wide slit (of a few arcsec) is used, which is kept at the same position angle (e.g. Bentz et al., 2008; Denney et al., 2010). The choice of a high-dispersion grating makes up for the lost spectral resolution and, due to its relatively narrow wavelength range, atmospheric dispersion effects can be ignored. The standard star used for the initial flux calibration is usually observed on the same night but not necessarily close in time.

The situation is very different in our case. Our near-IR spectra have a relatively large wavelength range, which we need to study the hot dust SED. Therefore, we observed with the slit at the parallactic angle in order to mitigate flux losses due to atmospheric dispersion. We used the narrowest slit offered ($0.3^{\prime\prime}$) for two reasons: to get the spectral resolution necessary to study line shapes and to minimize the contribution from the host galaxy to the continuum flux. However, the changing slit position angle and the narrow slit size mean that we have not necessarily enclosed the same amount of narrow line flux in all spectra, in particular, if the line-emitting region is elongated rather than circular or has considerable small-scale structure like filaments and knots. On the other hand, the telluric standard star, which we use for the initial flux calibration, was observed very close in time and airmass to the science target. Therefore, if slit losses due to seeing for the two were the same, the spectral flux scale is expected to be accurate even though the enclosed narrow line flux might vary with slit position angle.

We can best assess if our initial spectral flux scale needs to be adjusted to reflect an absolute flux scale by using a strong, non-variable, forbidden narrow emission line in the near-IR as observed with integral field unit (IFU) spectroscopy. Such observations are not subject to slit losses since they perform imaging spectroscopy and so can be used to measure the line flux expected in a given spectral aperture oriented at a certain angle. The strongest forbidden narrow emission line in the near-IR is [S III] $\lambda 9531$, but there are currently no IFU observations of this line for NGC 5548. Instead, Schönell et al. (2017) observed NGC 5548 in 2012 at the Gemini North 8 m telescope with the Near-Infrared Integral Field Spectrograph (NIFS; McGregor et al., 2003) and the Adaptive Optics system ALTAIR in the $J$ band ($1.14-1.36~{}\mu$m). These observations, obtained under a seeing of $0.28^{\prime\prime}$, which is similar to our IRTF slit width, and with a FOV of $3^{\prime\prime}\times 3^{\prime\prime}$, cover the weak, forbidden narrow emission line [Fe II] $1.2567~{}\mu$m. They show that the [Fe II] emission region is slightly extended (on subarcsec scales) and has a smooth but elongated structure (see their Fig. 4). For NGC 4151, there are near-IR IFU observations covering both the [S III] $\lambda 9531$ and [Fe II] $1.2567~{}\mu$m lines (Storchi-Bergmann et al., 2009). These observations show that the extended morphologies of the two lines are very similar and so we will assume in the following that any variation we measure for the [Fe II] line with changing spectral aperture and/or orientation angle holds also for the [S III] line. In any case, narrow-line reverberation results for NGC 5548 show that the [O III] emitting region is 1-3 pc in extent with a density of $\sim 10^{5}$ cm*-3* (Peterson et al., 2013). Since photoionisation models indicate that all the high-ionisation lines arise in this same component (Kraemer et al., 1998), we expect the [S III] emission to be dominated by a nuclear component.

We have re-reduced the NIFS $J$ band datacube for NGC 5548 and have produced spectra summed within a spectral aperture of $0.3^{\prime\prime}\times 1.2^{\prime\prime}$ centered on the nucleus and oriented at the same position angle as the IRTF SpeX spectral aperture (see Fig. 3 for an example). The chosen aperture length ensures that we enclose all the [Fe II] flux in any given direction (Schönell et al., 2017) and avoids the noisy outer regions of the NIFS FOV. We have then used the PrepSpec routine developed by K. Horne (described in Shen et al., 2015, 2016) to determine the [Fe II] flux variation in the extracted NIFS spectra. PrepSpec models the profiles of all emission lines, both broad and narrow, and the total continuum. A time-dependent scaling factor is included to match the profiles of selected narrow emission lines in all the spectra. As van Groningen & Wanders (1992) and Fausnaugh (2017) showed, matching the line profiles rather than directly measuring the line flux gives a much improved photometric alignment (usually \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}1-5\%). The line profiles are modelled as splines to keep them smooth, but are otherwise determined by the data. The median photometric correction factor is held at one, appropriate for a majority of the spectra being photometric. Using only the [Fe II] $1.2567~{}\mu$m line, we find a relatively small spread around the median for the different NIFS spectra, with an average photometric correction factor of only $\sim 2\%$. This is not unexpected since our slit orientation was very similar for most observations. This then means that we expect the [Fe II] line flux, and so the [S III] line flux, to be constant in our IRTF spectra. We have then used PrepSpec to determine the spectral photometric correction factors based on the strong [S III] $\lambda 9531$ line.

In one case, we can compare the initial spectral flux scale directly to GROND photometry results (Fig. 4). The GROND observation from 2017 March 18 is only one day after our near-IR spectrum from 2017 March 17, and so no significant variability is expected. The GROND nuclear fluxes in the $i^{\prime}$, $z^{\prime}$, $J$, $H$ and $K$ filters, which overlap in wavelength with our near-IR spectrum, are within $\sim 4\%$, 2%, 12%, 0.5% and 5%, respectively, of the spectral fluxes. Note that the much larger value for the $J$ band is due to two strong, broad emission lines being included in the filter passband. This high similarity between the GROND nuclear fluxes and the near-IR spectral fluxes not only confirms that our initial flux scale for this spectrum, which was observed under photometric weather conditions, is accurate at the few per cent level, but also strongly supports our assertion below that the IRTF SpeX spectra are barely contaminated by host galaxy flux (see Section 3.4). We have then assumed a photometric correction factor of unity for this spectrum and scaled the photometric correction factors from PrepSpec for all other spectra relative to it (instead of relative to the median). The results are listed in Table 1, column (13). The photometric correction factors lie in the range of $\sim 2-30\%$, with an average of $\sim 15\%$.

In Fig. 2, bottom panels, we compare the IRTF spectral continuum flux in the wavelength range of $8700-9000$ Å, both original and corrected, with the $z$ band light-curve from the LCOGT photometry. We have corrected the $z$ band light-curve for a constant host galaxy flux contribution, which we estimated using our host galaxy modelling of the GROND $z^{\prime}$ band image from 2017 March 18 and scaled it to match the near-IR spectral flux of the observation from 2017 March 17. Whereas large discrepancies between spectral and photometry points are evident in about half the cases for the original data, the two data sets align very well once we have applied the photometric correction factors based on the [S III] $\lambda 9531$ line.

Finally, we have also determined the photometric correction factors inlcuding only the three strongest, permitted narrow emission lines, namely, He I $1.08~{}\mu$m, Pa$\beta$, and Pa$\alpha$. Permitted transitions can be produced in gas of a higher number density and so can be potentially variable on the time-scales of our monitoring campaign. We obtain results similar to those using only the [S III] $\lambda 9531$ line, which indicates that these lines did not strongly vary.

3.4 The spectral continuum components

In this study, we intend to measure dust lag times and determine luminosity-based dust radii in order to compare the two. For the first we need to assemble the light-curves of both the irradiating (accretion disc) flux and hot dust emission, for the latter we need to measure the dust temperature and estimate the UV/optical accretion disc luminosity responsible for heating the dust. All of this is possible using our cross-dispersed near-IR spectra. Their large wavelength range covers about half of the extent of the hot dust SED, which means that we determine the dust temperature with high precision and measure the dust flux near its peak, where it is least contaminated by other emission components. Then, as Landt et al. (2011a, b) showed, the accretion disc still dominates the continuum flux at $\sim 1~{}\mu$m. So with the increased wavelength coverage in the blue of the refurbished SpeX (as short as $\sim 0.7$ $\mu$m) we sample a considerable part of this component, simultaneously with the hot dust SED.

Our near-IR spectral aperture also includes flux from the host galaxy. This flux, however, constitutes in this low-redshift source the entire local background and so is largely removed from the final spectrum during the data reduction process (see Section 3.1). This assertion is supported by the close similarity between the GROND PSF fluxes and near-IR spectral fluxes for the observation of 2017 March 17 (see Section 3.3 and Fig. 4). As a further test, we have, on the one hand, estimated for this observation the host galaxy flux enclosed in the near-IR spectral aperture using the correlation of Landt et al. (2011a) between spectral aperture and host galaxy flux based on Hubble Space Telescope images and, on the other hand, repeated the data reduction process without subtracting a local background. Both approaches give a similar result, indicating that the spectral flux in the $J$ band, where the host galaxy spectrum has its maximum (see e.g. Fig. 6 of Landt et al., 2011a), would be $\sim 20-30\%$ higher if the host galaxy flux was still included in the final spectrum.

In order to determine the dust temperature and luminosity and to estimate the UV/optical accretion disc luminosity, we have decomposed the spectral continuum into its two components. For this purpose, we have first approximated the rest-frame wavelength range of $0.7-1~{}\mu$m with an accretion disc spectrum, which we have subsequently subtracted from the total spectrum. We have then fitted the resultant hot dust spectrum at wavelengths $>1~{}\mu$m with a blackbody, representing emission by large dust grains, and with two modified blackbodies, approximating the emissivity of sub-micron silicate and carbon dust grains. Fig. 5 shows as an example the observation from 2017 March 17 for which the quasi-simultaneous GROND photometry extends the near-IR spectral range into the optical. Table 3 lists the relevant physical parameters extracted from the spectral decomposition.

3.4.1 The accretion disc spectrum

We calculated the accretion disc spectrum for a steady, geometrically thin, optically thick accretion disc, in which case the emitted flux is independent of viscosity and each element of the disc face radiates roughly as a blackbody with a characteristic temperature depending only on the mass of the black hole, $M_{\rm BH}$, the accretion rate, $\dot{M}$, and the radius of the innermost stable orbit (e.g., Peterson, 1997; Frank et al., 2002). We adopted the Schwarzschild geometry (non-rotating black hole) and for this the innermost stable orbit is at $r_{\rm in}=6\,r_{\rm g}$, where $r_{\rm g}=GM_{\rm BH}/c^{2}$ is the gravitational radius, with $G$ the gravitational constant and $c$ the speed of light. Furthermore, we have assumed that the disc is viewed face-on and that it extends to $r_{\rm out}=10^{4}r_{\rm g}$. With these assumptions, only the mass and accretion rate of the black hole are additionally required in order to constrain the accretion disc spectrum. The black hole mass of NGC 5548 has been derived by numerous optical reverberation campaigns and is estimated to be $M_{\rm BH}=(5.2\pm 0.2)\times 10^{7}~{}M_{\odot}$ (Peterson et al., 2004; Bentz et al., 2009; Bentz & Katz, 2015), assuming a geometrical scaling factor of $f=4.3$ (Grier et al., 2013) to convert the measured virial product to black hole mass. This black hole mass is uncertain by $\sim 0.4$ dex due to the intrinsic scatter of $f$ values for individual AGN, but is consistent with the value of $M_{\rm BH}=(3.9^{+2.9}_{-1.5})\times 10^{7}~{}M_{\odot}$ derived by Pancoast et al. (2014) based on dynamical modelling and the preliminary value of $\sim(7\pm 2)\times 10^{7}~{}M_{\odot}$ (Horne et al., in preparation) derived from the velocity-delay map of the 2014 HST campaign data, both of which do not require a scaling factor. The accretion rate can be obtained directly from an approximation of the accretion disc spectrum to the data. In Table 3, column (2), we list the resultant total accretion disc luminosity, $L_{\rm uv}$, which we have used in Section 3.4.2 to calculate the luminosity-based dust radius. We note that our approach attributes the observed changes in accretion disc luminosity to changes in accretion rate, although on the time-scales of our monitoring campaign the accretion disc variability is almost certainly due to reprocessing of the variable high-energy radiation (see also Section 4). However, the total accretion disc spectrum in the low energy range covered by our near-IR spectrum does not change substantially if reprocessed flux is added, and by omitting this component we expect to underestimate the total accretion disc luminosity by only $\sim 10-20\%$ (Gardner & Done, 2017, see their Fig. 3).

The spectrum of a standard accretion disc such as the one assumed here is observed in AGN for the variable (e.g. Collier et al., 1999) and polarized components of the optical/near-IR continuum (Kishimoto et al., 2008) and also for the total flux in AGN with negigible contribution from host galaxy starlight (e.g. Koratkar & Blaes, 1999; Landt et al., 2011b). Therefore, we are confident that our spectral decomposition reliably isolates the hot dust component. However, the 2014 HST reverberation campaign on NGC 5548 showed that, whereas the mean accretion disc spectrum as inferred from the light-curves by Starkey et al. (2017) is compatible with the standard accretion disc spectrum, albeit assuming a relatively low accretion rate, the mean observed UV/optical AGN continuum fluxes (corrected for host galaxy and line emission) presented by Fausnaugh et al. (2016) lie well below the extrapolation of the standard accretion disc spectrum to shorter wavelengths (Fig. 6). This spectral difference is also evident in our data. In Fig. 7, we show the mean and variable (rms) spectra for our campaign, which we have calculated following Peterson et al. (2004). Whereas the mean spectrum at rest-wavelengths \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}1~{}\mum is redder than the spectrum of a standard accretion disc, this model approximates well the rms spectrum in this wavlength range. In this respect, see also the example in Fig. 5, which shows the GROND $r^{\prime}$ band photometry to agree well with the calculated accretion disc spectrum but the observed $g^{\prime}$ band flux to lie well below it. The flux difference between the mean spectrum of Fausnaugh et al. (2016) and the accretion disc spectrum of Starkey et al. (2017) indicates the presence of a ‘constant red component’. This component has roughly the spectrum of a blackbody with a temperature of $T\sim 10^{4}$ K and, therefore, could be the reprocessor predicted by Gardner & Done (2017).

The presence of the ‘constant red component’ means that our accretion disc spectrum overestimates the UV/optical luminosity available to heat the dust, and so the dust radii based on it. Since in the two interpretations of the data from the 2014 HST reverberation campaign the total accretion disc luminosities differ by a factor of $\sim 10$, we estimate that we may have overestimated the hot dust radius by up to a factor of $\sim 3$. The dominance of the ‘constant red component’ at the shortest wavelengths of our near-IR spectrum could also mean that it is its spectral shape rather than that of the accretion disc that determines the resultant dust spectrum. If we assume this, more flux at short wavelengths is attributed to dust than in the case of the accretion disc, thus broadening the dust spectrum (Fig. 5). However, this scenario is unlikely. Since the shape of the variable (rms) and mean spectrum are similar at wavelengths \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}1~{}\mum (Fig. 7), only the variable components (accretion disc and hot dust) contribute to the near-IR flux.

The spectrum of the isolated hot dust component (Section 3.4.2) strongly depends on the assumed accretion disc spectrum at long wavelengths, which in turn is sensitive to the chosen outer radius of the accretion disc. We have assumed a very large outer radius, which gives a maximum contribution from the accretion disc at near-IR wavelengths. In order for the accretion disc to dominate the total flux at $\sim 1~{}\mu$m, as established by the near-IR radius-luminosity relationship (Landt et al., 2011a), the outer radius can be as small as $r_{\rm out}\sim 1500\,r_{\rm g}\sim 4.4$ light-days. Choosing instead this smaller value, does not alter significantly the resultant hot dust spectrum. Accretion discs are assumed to become unstable and fragment if their extent considerably exceeds the self-gravity radius. In our case, this value is reached at $r_{\rm out}\sim 500\,r_{\rm g}\sim 1.5$ light-days. If we assume this small outer radius, the accretion disc spectrum at near-IR wavelengths is very similar to that of the ‘constant red component’, thus resulting in a similarly broad dust spectrum (Fig. 5). We further consider this case below and show that such a broad dust spectrum is unlikely, since it implies a large temperature range and thus that the hottest dust is much more variable than the coldest dust. This would change the spectral shape of the variable (rms) spectrum relative to that of the mean spectrum, which is contrary to what we observe (Fig. 7). An unexpectedly large accretion disc radius for NGC 5548 was also inferred by Edelson et al. (2015) based on the 2014 HST reverberation campaign data. They found a lag time of $\tau=0.35$ light-days at 1367 Å, which translates to a radius of r_{\rm out}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}5 light-days for wavelengths \lambda\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}1~{}\mum using the standard relationship of $\tau\propto\lambda^{4/3}$.

3.4.2 The hot dust spectrum

Assuming that the hot dust is optically thick to the incident UV/optical radiation but optically thin to its own thermal radiation, its emission depends on the temperature, $T$, and the wavelength, $\lambda$, as:

[TABLE]

where $B_{\lambda}(T)$ is the Planck function and $Q_{\lambda}(a)$ is the emission efficiency of the dust at IR wavelengths, which depends on the grain size $a$. Assuming radiative equilibrium, the dependence of the dust temperature on the luminosity of the irradiating source can then be derived by equating the flux absorbed by a spherical dust grain to that radiated by it:

[TABLE]

where $R_{\rm d,lum}$ is the luminosity-weighted dust radius and $Q_{\lambda}^{\rm abs}$ is the absorption efficiency of the dust in the UV/optical. We assume that $Q_{\lambda}^{\rm abs}=1$, i.e. that the absorption cross section equals the geometrical cross section (e.g Rees et al., 1969) and Planck-average to get:

[TABLE]

where $\sigma$ is the Stefan-Boltzmann constant and $\langle Q^{\rm em}\rangle$ is the Planck-averaged value of $Q_{\lambda}(a)$. We consider here only a single grain size, whereas in general the dust will have a grain size distribution. However, it can be shown that this is a reasonable approximation for our limited near-IR wavelength range if it samples a single temperature, since the emission is dominated by the largest and hottest grains.

We have subtracted from the total near-IR spectrum that of the accretion disc and have fitted eq. 1 to the resulting dust spectrum at wavelengths $>1~{}\mu$m (Fig. 5). For this we have assumed that the IR emissivity can be approximated with a power-law of the form $Q_{\lambda}(a)\propto\lambda^{\beta}$. We have used the extinction efficiencies recently calculated by B. Draine and collaborators to determine the power-law slope $\beta$ in the wavelength range of interest ($\lambda=1-3~{}\mu$m) for different grain sizes ($a=0.001-10~{}\mu$m) for (turbostratic and single-crystal) graphite (Draine, 2016) and (amorphous) silicates (Draine & Hensley, 2017). The results are shown in Fig. 8. Both dust types span a relatively large range in $\beta$, however, the assumption of $\beta=-2$ and $\beta=-1$ is appropriate for carbon and silicate dust, respectively, for small grain sizes of a\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.1~{}\mum. The blackbody case ($\beta=0$) is reached for grain sizes of a\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}0.4~{}\mum and \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}2~{}\mum for carbon and silicate dust, respectively. We have performed the fits with the C routine MPFIT (version 1.3a; Markwardt, 2009), which solves the least-squares problem with the Levenberg-Marquardt technique, and have fitted for the temperature and flux scaling assuming three emissivity laws ($\beta=0$, $-1$ and $-2$). We have included in the fit only the continuum part of the near-IR spectrum, i.e., we have excluded emission lines, and have rebinned it to $\Delta\log\nu=0.01$ Hz. The errors were calculated in the usual way as the standard error on the binned flux average. We note that, since about half the hot dust spectrum is constrained by the data, its temperature and flux are well defined. We obtain average temperatures of $\langle T\rangle=1450\pm 7$ K, $1225\pm 5$ K and $1061\pm 4$ K for an emissivity law with $\beta=0$, $-1$ and $-2$, respectively.

We have then calculated luminosity-weighted dust radii using the best-fit dust temperatures and approximated $L_{\rm uv}$ with the accretion disc luminosity. For the Planck-averaged emission efficiencies we used in the case of $\beta=-1$ a value of $\langle Q^{\rm em}\rangle=0.0210$ appropriate for silicates of $T=1259$ K and $a=0.1~{}\mu$m (Laor & Draine, 1993) and in the case of $\beta=-2$ a value of $\langle Q^{\rm em}\rangle=0.0875$ appropriate for graphite of $T=1000$ K and $a=0.1~{}\mu$m (Draine, 2016)111The calculated data are provided by B. Draine in machine readable form at http://www.astro.princeton.edu/~draine/dust/. We obtain average luminosity-weighted dust radii of $\langle R_{\rm d,lum}\rangle=61\pm 1$ light-days, $385\pm 8$ light-days and $590\pm 12$ light-days in the case of a blackbody, and small-grain carbon and silicate dust, respectively. Since for the latter we have considered only solutions for $a=0.1~{}\mu$m, the resulting hot dust radius will be larger by up to a factor of $\sim 10$ for grain sizes as small as $a=0.001~{}\mu$m. In addition, we might have overestimated in all three cases the total accretion disc luminosity by up to a factor of $\sim 10$ (see Section 3.4.1), resulting in an overestimation of the hot dust radius by a factor of $\sim 3$.

The unknown outer radius of the accretion disc introduces an uncertainty in the shape of the hot dust spectrum. As discussed in Section 3.4.1, assuming that the disc cuts off already at the self-gravity radius, gives a much broader dust spectrum. Such a broad spectrum is much wider than a single blackbody (Fig. 5). Fitting the sum of two blackbodies to this broad dust spectrum for the observation from 2017 March 17, we obtain best-fit temperatures of $T\sim 1800$ K and $\sim 960$ K, which results in dust radii of $R_{\rm d,lum}\sim 35$ light-days and $\sim 120$ light-days for the two components. This large difference in radius is expected to affect the observed variability, which should be dominated by the fastest variations and so by the hottest component. Therefore, in this scenario, we expect the spectral shape of the variable (rms) spectrum to be dominated by the hottest dust component and so to be significantly different from that of the mean spectrum. This is contrary to what we observe; the spectral shape at wavelengths \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}1~{}\mum in the rms spectrum is similar to that in the mean spectrum (Fig. 7), indicating that the hot dust region spans a very narrow range in radius thus mimicking a single dust component. A narrow range in dust radius is also suggested by our reverberation results, which do not show any significant difference between the $H$ and $K$ bands (Section 4.3).

4 Response-weighted dust radii

It is currently assumed that the variability of the hot dust is driven by the variability of the accretion disc, since this provides the UV/optical photons necessary to heat the dust grains. In turn, the emission variability of the accretion disc flux on time-scales less than a year is attributed to reprocessing of high-energy radiation from an intrinsically variable source located above the accretion disc and illuminating it (the lamp-post model). Assuming that the dust absorbs and re-emits instantaneously the incident UV/optical photon flux, we can estimate the dust radius based on light-travel time arguments alone as $R_{\rm d,rev}=\tau\,c$, where $\tau$ is the reverberation lag time, i.e. the time between the emission of the accretion disc flux and the re-emission of this flux by the dust, with the lag in the observer’s frame increased by a factor of $1+z$. In practice, the dust will have an extended geometry and so an emissivity distribution that has a lag distribution (also referred to as delay map, transfer function or response function) associated with it.

The difficulty lies in deriving the delay map from the observed light-curves. The echo light-curve is just a smoothed and delayed version of the driving light-curve, i.e. the driving light-curve shifted by a lag and convolved with the geometry-dependent lag distribution (Blandford & McKee, 1982; Peterson, 1993). With very high-quality data, i.e. extremely well-sampled light-curves, the observed transfer function can be recovered by Fourier transforming in time, and so knowledge about the geometry and physical conditions of the reprocessor gained. If the data are of relatively poor-quality, i.e. irregularly and sparsely sampled light-curves, the observed transfer function can be recovered by other methods, e.g. the maximum entropy method (Horne, 1994) and Markov chain Monte Carlo (MCMC) sampling of model parameters (Pancoast et al., 2011; Starkey et al., 2016), or just the mean lag time estimated using, e.g., cross-correlation methods (Gaskell & Peterson, 1987; White & Peterson, 1994; Zu et al., 2011).

Koshida et al. (2014) used cross-correlation analysis to estimate the mean hot dust response time in NGC 5548 for six monitoring periods and found values varying in the range of $\sim 40-80$ days. This result is interesting since it may imply that the hot dust distribution is not constant with time but changes, possibly due to dust grain destruction and reformation processes (Koshida et al., 2009). The length of our near-IR monitoring campaign is 336 days and so spans $\sim 4-8$ times the mean reverberation lag time. However, our time sampling is not frequent enough to allow us to test for lag variability. In the following, we first derive the light-curves corresponding to both the accretion disc and hot dust (Section 4.1). We then use three different models to infer the lag distribution for the observed hot dust light-curve (Section 4.2), which we find to converge to a very similar result (Section 4.3).

4.1 The observed spectral light-curves

We have assembled the light-curves of both the accretion disc (Table 4) and the hot dust component (Table 5) from fluxes measured in the near-IR spectrum. Ideally, we would like to measure the accretion disc flux at the shortest wavelengths, where it is least contaminated by the hot dust component. However, the $S/N$ ratio of the near-IR spectra decreases significantly towards short wavelengths and is relatively low in particular in the region $<0.8~{}\mu$m ($S/N\sim 5$). Therefore, we have chosen to measure the accretion disc flux in the 60 Å wide rest-frame wavelength region of $\lambda=9730-9790$ Å, which lies between the two broad hydrogen emission lines Pa$\epsilon$ and Pa$\delta$ and is known to be line-free, and also in the wider (300 Å) rest-frame wavelength region of $\lambda=8700-9000$ Å, which does not appear to be contaminated by significant emission-line flux in NGC 5548 (see Fig. 1). We have measured the average flux in the selected wavelength ranges and have calculated error bars in the usual way as the standard error on the flux average. We note that these errors are only a few per cent and so on average much smaller than the applied photometric correction factors.

We have measured the flux of the hot dust in the near-IR spectrum after the subtraction of the accretion disc component. This is important not only because the accretion disc flux dominates almost the entire $J$ band, but also because the accretion disc contribution, although it decreases with increasing wavelength, remains non-negligible even at the long-wavelength end of the $K$ band, where the hot dust spectrum peaks (e.g. $\sim 15\%$ of the total flux at $2.4~{}\mu$m in the observation from 2017 March 17; see Table 5). Therefore, we have chosen to measure the flux of the hot dust in two line-free, 500 Å wide rest-frame wavelength regions, namely, $\lambda=2.34-2.39~{}\mu$m, which is at the red end of the $K$ band but close to the peak of the blackbody spectrum, and $\lambda=1.55-1.60~{}\mu$m, which is approximately in the middle of the $H$ band.

During this campaign, NGC 5548 was in a relatively low emission state (on average a factor of $\sim 2$ lower than during the 2014 HST campaign), but was highly variable (Fig. 2). The continuum fluxes at the shortest wavelengths corresponding to the accretion disc had maximum variations of $\sim 50-60\%$ and the light-curves show several clearly isolated peaks and dips. This high variability of the observed light-curves is a prerequisite for recovering lag times from the data.

4.2 The MEMEcho and MCMCRev formalisms

MEMEcho (Horne et al., 1991; Horne, 1994) is a formalism based on the maximum-entropy method and infers the lag distribution given a driver and a reprocessor light-curve. In short, MEMEcho fits a linearised echo model:

[TABLE]

where both the observed driver light-curve, $F_{\rm uv}(t)$, and the observed reprocessor light-curve, $F_{\nu}(t)$, are assumed to be composed of a constant and a variable component. The total flux of the reprocessor is then modelled as the sum of the constant (mean) reprocessor flux, $\overline{F}_{\nu}$, and the convolution of the response function, $\Psi_{\nu}(\tau)$, with the variable part of the driver light-curve, $\Delta F_{\rm uv}(t)$. MEMEcho recovers from the data the union of the three positive functions, $\Psi_{\nu}(\tau)$, $F_{\rm uv}(t)$ and $\overline{F}_{\nu}$. The total entropy of the image is the sum of the entropies of these three subimages and of the many different images that satisfy the data constraint, the image with maximum entropy is selected.

MCMCRev (similar to CREAM; Starkey et al., 2016) uses the MCMC technique to sample posterior parameter distributions for a model response function given the light-curves for the driver and the reprocessor. We have modelled the response function with a log-normal distribution and with a temperature-radius relationship for a standard accretion disc seen face-on ($T\propto R^{-3/4}$). The latter is motivated by recent models of the dusty torus as a radiatively accelerated outflow launched from the outer regions of the accretion disc (Czerny & Hryniewicz, 2011; Czerny et al., 2017). In effect, the temperature-radius relationship translates into a wavelength-lag relationship.

4.3 Reverberation results

We have run MEMEcho with both the $8700-9000$ Å spectroscopic light-curve and the $V$ band photometric light-curve as the driver. The delay distributions for both the $H$ and $K$ dust light-curves have a clear peak, which is at $\sim 40-45$ days and at $\sim 70-80$ days using the spectroscopic and photometric light-curves as the driver, respectively. The width of the MEMEcho delay map is much narrower in the former case, most likely due to the spectroscopic light-curve being too smooth.

For an MCMC fit, posteriors on the model parameters can be interpreted as measurements with error bars. In this case, the centroid of the response function gives the average response-weighted dust radius, which is defined as:

[TABLE]

Using the $V$ band light-curve as the driver and modelling the response function with a log-normal distribution, we get similar mean values for the dust delay map in the $H$ and $K$ band of $\langle R_{\rm d,rev}\rangle=73\pm 7$ light-days and $79\pm 8$ light-days, respectively (Fig. 9). But sparse time sampling and relatively large error bars on the infrared light-curve data limit the information that can be obtained on the delay distribution. In this situation, the mean delay can be constrained fairly well, but not the detailed shape. With this caveat, the log-normal delay maps for the $H$ and $K$ bands peak at $\sim 20$ and 35 days, respectively. We note also that evidence for a lack of prompt response in $H$ and $K$ occurs just after the $V$ band flux has dropped sharply between days 7810 and 7825. In response to this the $8700-9000$ Å flux has dropped by day 7830, but the $H$ and $K$ fluxes remain high at days 7810 and 7835.

Assuming instead for the response function the temperature-radius relationship of a standard accretion disc ($\tau\propto\lambda^{4/3}$) does not change these values significantly (Fig. 10). We get average values of $\langle R_{\rm d,rev}\rangle=35\pm 7$ light-days and $55\pm 9$ light-days for the $H$ and $K$ band, respectively (Table 6). However, the results show consistency with the accretion disc temperature-radius relationship for both the average and peak values. Moreover, the delay distribution is strongly skewed towards low values with the peak at $\tau_{0}=7\pm 1$ light-days and $12\pm 2$ light-days for the $H$ and $K$ band, respectively. We discuss this further in Section 5.3.

5 The torus in NGC 5548

There are three ways to measure or estimate the radius of the hot dust in AGN: (i) through the time delay in the response of the dust emission to the central heating source; (ii) assuming radiative equilibrium for known temperature and grain properties, and (iii) measuring the distance of the dust emission by geometrical means, e.g. with near-IR interferometry or high spatial-resolution imaging. In general, it is found that interferometric dust radii are a factor of $\sim 2$ larger than reverberation dust radii (Kishimoto et al., 2009, 2011; Koshida et al., 2014). Since the former are effectively luminosity-weighted dust radii, this result can be understood if one assumes that response-weighted radii are dominated by the fastest variations, which are expected from the smallest radii, whereas luminosity-weighted radii are dominated by the emission from grains with the highest emissivity, i.e. the hottest and largest grains. The highest temperature of the largest grains will, in any case, be lower than the highest temperature of the smallest grains, due to their lower heat capacity, and so their emission will come from further out. Alternatively, Kawaguchi & Mori (2010, 2011) ascribe the small reverberation radii to a dust geometry that is concave (or bowl-shaped) rather than spherical due to the anisotropy of the accretion disc emission; since the UV/optical accretion disc luminosity is reduced for larger viewing angles (e.g. Netzer, 1987; Hubeny et al., 2000), so is the dust sublimation radius (on average by a factor of $\sim 2-3$). The interferometric radii would then be dominated by the region containing the bulk of the dust, which in this geometry is located further out. Comparisons between reverberation dust radii and luminosity-weighted dust radii derived from SEDs have yielded so far ambiguous results (Mor & Netzer, 2012; Landt et al., 2014; Koshida et al., 2014), most likely due to the unknown grain properties and a poorly constrained dust temperature and irradiating luminosity. In this study, we combine for the first time the simultaneous measurement of a reverberation radius and of a well-constrained luminosity-based radius and show that it can constrain the grain size (Section 5.1), the chemical composition (Section 5.2) and the origin (Section 5.3) of the dust.

5.1 The location

Our response-weighted radius of $\sim 70$ light-days (Section 4.3) is only $\sim 15\%$ larger than the average luminosity-weighted radius for dust of large grains, but much smaller than the average luminosity-weighted radius for dust of small grains (by a factor \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}6 and \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}8 for carbon and silicate dust, respectively; Section 3.4.2). Taking into account the uncertainty of the irradiating luminosity in the calculation of the luminosity-weighted radius gives that in the blackbody case the reverberation radius is larger than the average luminosity-based radius by a factor of $\sim 3$ and in the case of small-grain carbon and silicate dust smaller by a factor \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}2 and \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}3, respectively. Clearly, the smallest discrepancy between the response- and luminosity-weighted radius is for dust of large grains and for an irradiating luminosity as estimated.

Given that we have almost certainly overestimated the UV/optical accretion disc luminosity (see Section 3.4.1), the similarity between the reverberation and luminosity-based dust radius in the case of a blackbody emissivity law is intrigueing. A reason for this could be that the dust is also heated by X-rays. For small grains, the absorption efficiency drops sharply for frequencies larger than the UV, but for large grains it continues to be high. Mehdipour et al. (2015) have presented a well-sampled, multiwavelength SED of NGC 5548, which extends to hard X-ray frequencies and was corrected for host galaxy emission but included the ‘constant red component’ (see their Fig. 6, right panel). Scaling this SED to our $8700-9000$ Å continuum flux gives an irradiating luminosity $\sim 40-70\%$ higher than our (over)estimated accretion disc luminosity, which translates into an increase in dust radius of $\sim 20-30\%$ and so an average radius for the blackbody case of $\sim 70-80$ light-days. This plausible match with the reverberation radius indicates that X-rays are important in heating the dust and constrains the grain size further to a\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}2~{}\mum for both graphite and silicates (Draine, 2016; Draine & Hensley, 2017).

It was previously proposed that the dust composition in the circumnuclear regions of AGN is dominated by large grains (of a few $\mu$m). Laor & Draine (1993) found that small grains in an optically thin dust configuration are likely to be destroyed on short timescales in an AGN environment, whereas very large grains are not. Large grains in AGN were also postulated by Maiolino et al. (2001) to explain the observed lack of prominent silicate absorption features in the mid-IR spectra of Seyfert 2s and the lack of a strong 2175 Å (carbon) feature in the UV spectra of reddened Seyfert 1s, since they make the extinction curve flatter and featureless. Furthermore, the high-density environment at the centre of an AGN is expected to be a natural catalyst for the formation of large grains by coagulation. On the other hand, almost certainly not all dust is concentrated in large grains as shown by IR polarisation studies, which are sensitive to the smallest dust grains and give constraints of a\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.1~{}\mum (e.g. Lopez-Rodriguez et al., 2013, 2017).

5.2 The dynamics

It is usually assumed that the inner edge of the torus is set by thermal sublimation processes, which implies that dust can exist as far in as is physically possible without being destroyed via direct heating by the central irradiating source, and that the dust radius changes mainly in response to the AGN luminosity (the ‘receding torus model’; Lawrence, 1991). We can now test this assumption with our measured dust temperatures and their variation with time, since, as we have shown in Section 3.4.1, the most variable dust component also dominates the emission. If the dust is close to sublimation, the observed temperature can constrain the chemical composition of the dust, since the expected sublimation temperatures are very different for carbon and silicate dust ($T_{\rm sub}\sim 1800-2000$ K and $\sim 1300-1500$ K, respectively; Salpeter, 1977; Lodders, 2003). If the hot dust is much cooler than its sublimation temperature, this implies that the dust-free inner region in AGN is enlarged and that the inner edge of the torus is set by a different process and so maybe independent of the irradiating luminosity. In this respect we note that if the central luminosity is powerful enough to heat the dust to sublimation, the dust temperature will generally be close to sublimation and there will not be an enlarged dust-free inner region “evacuated” by a previous high-state in luminosity stable on time-scales longer than grain reformation times (i.e. a few months). In addition, since $L_{\rm uv}\propto T^{4}$ (eq. 3), only a strong increase in luminosity (by factors of a few) and lasting for several months (i.e. much longer than typical grain reformation times) could have produced the current radius at which the dust temperature is now well below sublimation. This required variability behaviour is not usually observed for radio-quiet AGN. In particular for NGC 5548, large UV flux variations occur only over shorter periods of weeks (e.g. De Rosa et al., 2015).

We have found a temperature range of $T\sim 1390-1515$ K for the blackbody case222We will consider here and in the following only the blackbody case, since we have shown in (Section 5.1) that the nuclear dust in NGC 5548 is composed mainly of large grains. (Section 3.4.2), which is consistent with silicates at sublimation and is $\sim 300-500$ K below the sublimation temperatures of carbon dust. Since the growth and survival of dust particles requires low enough grain temperatures (well below the condensation limit), this result indicates that the hot dust consists mainly of carbon, with silicates having sublimated. This conclusion is also supported by previous near-IR spectra of NGC 5548, where blackbody dust temperatures as high as $T\sim 1700$ K were found (Landt et al., 2015). The fact that we have measured a dust reverberation signal in this campaign also points to the dust being well below sublimation so that further heating is possible without destroying the dust. Both temperature and mass changes (e.g., due to grain destruction) can introduce variations in the dust flux. However, assuming instantanous heating, the change in dust flux due to a change in temperature according to eq. 3 is expected to dominate the variations since grain reformation times are relatively large (of the order of a few months). Then, if the dust was too close to sublimation, only a change in dust mass could introduce a significant change in dust flux. But such a flux change will not resemble a reverberation signal since it will not be the delayed and smoothed version of the driver light-curve. Since sublimation timescales are short and sublimation is a runaway process once the grain starts to evaporate, an increase in irradiating luminosity will quickly decrease the dust flux and this will stay low and unchanged until grains have reformed and the flux can increase again. But grain reformation not only takes a relatively long time, it also requires the irradiating luminosity to stay continuously low for that (prolonged) period.

In Fig. 11 we show the temporal behaviour of the dust flux and temperature. The two are very similar, i.e. the change in dust flux during our campaign is mainly due to a change in temperature, and the temperature variations are consistent with direct heating by the accretion disc. We get a $1\sigma$ dispersion around the mean $8700-9000$ Å continuum flux of $\sim 13\%$ ($\langle f_{\lambda}\rangle=(2.74\pm 0.38)\times 10^{-15}$ erg s*-1* cm*-2* Å*-1*), based on which we expect a dispersion around the mean temperature of $\sim 3\%$. The observed dispersion is $\sim 2\%$ ($\langle T\rangle=1450\pm 32$ K). Based on the maximum variation in the $8700-9000$ Å continuum flux of $\sim 56\%$, we expect a difference between the minimum and maximum temperature of $\sim 12\%$, which translates to $\sim 170$ K, similar to the observed value of $\sim 130$ K.

If it is not sublimation, then what sets the location of the inner edge of the torus? Or, put differently, why is the hot dust located much further out and so the dust-free inner region enlarged? There are two possibilities; either there was initially carbon dust in the inner regions permitted by sublimation but it has since been destroyed by localized processes or this inner region has always been dust-free. The latter possibility is unlikely if gas is to be expected in this region, since where there is gas there usually is dust. Furthermore, the gas densities in the inner region are expected to be higher than further out, which favours grain formation and growth. Therefore, most likely, dust destruction has taken place. If so, we would expect this process, which should be observable as a change in dust mass, to also affect the dusty region to some degree, since its efficiency would drop only gradually with radius.

In Fig. 12, we show the dust flux versus the dust temperature. Also plotted are heating tracks that show the expected flux increase due to a temperature increase alone (i.e. for a constant dust mass). It is clear that the data cannot be fit by a single such heating track, i.e. the spread in flux is larger than what is expected from a change in temperature alone, and so a changing dust mass has to be assumed additionally. In Fig. 11 (bottom panel) we show the normalised dust mass relative to the start of the campaign and connected in time. As expected for a grain destruction and reformation scenario, the dust mass is the highest at the beginning of the campaign. This is followed by the largest mass loss (by $\sim 40\%$, significant at the $6\sigma$ level), from which the dust mass recovers close to its initial value after $\sim 5.5$ months, followed by a renewed mass loss (again by $\sim 40\%$, significant at the $7\sigma$ level), from which the dust starts to recover only at the end of our campaign, i.e. again after $\sim 5$ months. Our estimate of a grain reformation time of $\sim 5-6$ months is within the usual theoretical range and of the same order as that estimated for Fairall 9 by Barvainis (1992) (a few months) and for NGC 4151 by Koshida et al. (2009) (about a year). In this respect, we note that if the hot dust was close to its sublimation temperature, we would also expect a decrease in dust mass due to grain evaporation, but this would not necessarily be observable as a decrease in near-IR flux. For a large distribution in grain size, sublimation will affect mainly the grains with the smallest sizes, since they have the smallest heat capacity, but these grains do not dominate the emission. If the grain size distribution is relatively narrow, as we find in our case, all grains will be similarly affected by sublimation. But, since sublimation timescales are short (Baskin & Laor, 2018), this loss of dust mass loss will not only be drastic but sudden. Finally, we note that our observation of a dust mass change that is not accompanied by large variations in dust temperature further indicates that the enlarged dust-free inner region was not “evacuated” by a previous high-state in luminosity and is now being replenished, but that it is rather kept dust-free by on-going in-situ processes. We will return to this point in Section 5.3.

Temperature estimates in photometric dust reverberation campaigns are rare since they require simultaneous observations in several filters and in any case have a high uncertainty (Clavel et al., 1989; Glass, 2004). Temperature variations in such campaigns have only recently been studied. Schnülle et al. (2013, 2015) performed a 5-month photometric monitoring campaign on NGC 4151 with GROND covering the $J$, $H$ and $K$ bands. They also found large blackbody temperature variations (of $\sim 200$ K) and changes in covering factor (both a decrease and increase).

Our result that the hot dust is not close to sublimation predicts that the location of the inner edge of the torus is largely independent of the irradiating luminosity and so rather a “torus wall”. This is confirmed for NGC 4151 by the observed variability of both the interferometric and reverberation dust radius. Pott et al. (2010) presented a comprehensive analysis of all then available near-IR interferometric observations and concluded that the hot dust radius is largely luminosity-invariant. Koshida et al. (2009) measured the reverberation radius for eight separate epochs over a six-year period and found it to vary largely independent of luminosity, making dust destruction and reformation a necessary part in the explanation of the observations. The length of our near-IR monitoring campaign covers in principle the mean reverberation lag several times, but, unfortunately, our time sampling is not frequent enough to allow us to test for lag variability.

Finally, we note that, as discussed by Baskin & Laor (2018), the dust sublimation temperature also depends weakly on the gas number density. The sublimation temperatures assumed here for graphite and silicate are relevant for the high-density gas expected in the torus (of $n_{e}\sim 10^{7}-10^{9}$ cm*-3*). If we assume that graphite was at sublimation at $T\sim 1450$ K, it would imply a density of only $n_{e}\sim 10^{4}$ cm*-3* (see their eq. (32)). But such rarified gas at a distance of $\sim 70$ light-days would be very hot if it was exposed to the total bolometric luminosity and the dust grains would be destroyed by thermal sputtering. An interpretation of our results in terms of silicate dust being close to sublimation would imply that it is carbon-poor. However, amorphous carbon as the main constituent of the hot dust is required by dust-driven wind models since it dominates the overall dust opacity and so provides the driving force of the wind. We will expand on this point in the next Section 5.3. In any case, it lies at hand to associate the enlarged dust-free inner region with the coronal line region (Pier & Voit, 1995) and we plan to explore this possibility in the future.

5.3 The geometry

In recent years, an understanding of the tight connection between the extent of the BLR and the innermost part of the dusty torus has started to emerge (Czerny & Hryniewicz, 2011; Goad et al., 2012; Landt et al., 2014). This connection is made in particular within models for radiatively accelerated dusty outflows launched from the outer regions of the accretion disc (e.g. Elitzur & Shlosman, 2006; Czerny et al., 2016, 2017). In this case, the large opacity of carbon dust leads to an inflated disc structure of which the (cold) back side naturally forms a hot dusty torus (Baskin & Laor, 2018). But this hot dust, will not necessarily be externally illuminated, since it lies largely in the plane of the accretion disc, and so it will not reverberate in the classical sense. In addition to dust located within a toroidal structure and to that condensing out in the accretion disc material, a third component could be present in the polar rather than the equatorial region if the polar dust observed at mid-IR wavelengths (Landt et al., 2010; Hönig et al., 2013; Fuller et al., 2019) extends close into the nucleus. Our results can constrain which of these three hot dust components is observed.

We start by using our simultaneous measurement of the temperature and reverberation radius together with its emissivity law to “calibrate” the temperature-radius relationship for the hot dust. We show this relationship in Fig. 13 (red solid line) for the observation from 2017 March 17. This temperature-radius relationship, which is that of a blackbody and so of the form $T\propto R^{-0.5}$, is that of dust exposed at all times to the total bolometric luminosity. The bolometric luminosty for this observation is $\log L_{\rm bol}\sim 44.6$ erg s*-1*, if we scale the multiwavelength SED of NGC 5548 presented by Mehdipour et al. (2015) to our $8700-9000$ Å continuum flux. We then compare the temperature-radius relationship of the hot dust to that for the accretion disc, which is of the standard form $T\propto R^{-0.75}$, and so steeper than that for the hot dust (black solid line). We have considered both the accretion disc luminosity that approximates the observed continuum flux as well as the lowest possible luminosity, with the latter lower than the former by a factor of $\sim 10$, due to the presence of the ‘constant red component’ (Section 3.4.1). From these curves, one can immediately see that dust is expected to condense out in the accretion disc matrial at radii much smaller than those permitted by the total bolometric luminosity. Then, if this dust was to be raised off the accretion disc to considerable heights, e.g., by a dusty outflow, it would immediately be destroyed, as discussed in Czerny et al. (2017). But, if this dust exists, it would emit thermal near-IR radiation and we should be able to observe it. If it is carbon dust, we expect it in our case at a radius of $\sim 3-5$ light-days if it is close to sublimation ($T\sim 1900$ K) and at a radius of $\sim 4-9$ light-days if it is at our observed temperature ($T\sim 1450$ K). Assuming for the response function the temperature-radius relationship of a standard accretion disc, we obtained in Section 4.3 average delay times much larger than this, but peak values of $\sim 7\pm 1$ and $12\pm 2$ light-days for the $H$ and $K$ band, respectively (vertical black dashed lines in Fig. 13). Therefore, within the errors, this reverberation signal is consistent with dust in the accretion disc. But if so, this implies relatively large accretion disc outer radii of at least $R_{\rm out}\sim 1700-3000\,r_{\rm g}$.

Our observation of a changing dust mass implies that on short time-scales (of about a year) the mass loss rate is not balanced by the mass inflow rate, i.e. the inner edge of the torus is not in a steady-state. This is expected in wind scenarios for the torus (Königl & Kartje, 1994; Elitzur & Shlosman, 2006; Keating et al., 2012). But even in the hydrostatic picture, radiation pressure is expected to dominate over gravity at the inner edge and so clouds in this (small) region will not be able to maintain bound orbits (Krolik & Begelman, 1986). The observed loss in dust mass could also be due to dust destruction (of the largest grains) by in-situ processes. Possible dust destruction processes that could happen locally, e.g., in fast shocks, are grain sputtering caused by collisions with high-energy ($>10$ eV) atoms or ions and grain fragmentation (or even vaporization) through grain-grain collisions, which are expected to almost completely remove large grains ($a\geq 0.1~{}\mu$m) (Jones, 2004). Fast shock waves are naturally expected in the outer accretion disc through the dissipation of turbulence, dynamical instabilities or spiral shocks. In the case of NGC 5548, the fast ($\sim 1000-5000$ km s*-1*) and long-lived outflow from locations just outside the BLR detected by Kaastra et al. (2014) during the 2013 HST/XMM-Newton monitoring campaign could also be a source for fast shocks.

Finally, if there was hot polar dust present on subarcsecond scales (Hönig et al., 2013), in addition to the dusty torus, we would have a situation similar to that described in Kawaguchi & Mori (2010) (see their Fig. 2) and Baskin & Laor (2018) (see their Fig. 13); due to the anisotropy of the irradiating (accretion) disc emission ($L_{\rm bol}\propto\cos\theta$, with $\theta$ the angle between the accretion disc rotation axis and the location of the dust), the dust in the polar regions will be located further out than dust in the equatorial plane if both have the same temperature. How much further out the polar dust will be depends on the scale height of the torus, i.e. on its average $\theta$. If we assume for this the usual value of $\theta\sim 60-70^{\circ}$, the polar dust would be $\sim 50\%$ further away than the torus dust. But this radius difference will increase fast as the average $\theta$ of the torus increases further, as would be expected, e.g. for the rather flat dusty structure of Baskin & Laor (2018). Clearly, if the polar and torus hot dust have the same temperature, the polar dust is not expected to dominate the variations. If we assume instead that the polar dust is hotter than the hot torus dust, it will have to dominate both the emission and the variations, since in our campaign we observe the variable (rms) spectrum and the mean spectrum to have the same spectral shape at wavelengths \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}1~{}\mum (Section 3.4.1). This would then imply that the hot dust in the torus is at much lower temperatures ($T<<1000$ K) in order for it to not dominate the variations (if not the emission). The emission of such dust peaks at wavelengths well outside our near-IR spectral range ($\lambda>>3~{}\mu$m) and so we would not have observed it. But, in this case, it would be difficult to explain why there is hot polar dust but not hot torus dust (at $T\sim 1450$ K).

6 Summary and conclusions

We have conducted the first spectroscopic dust reverberation programme on the well-known AGN NGC 5548. Spectroscopy, unlike photometry, can simultaneously determine several dust properties, such as its flux, temperature and covering factor. In particular, it can measure both response-weighted and luminosity-weighted dust radii. Our main results can be summarised as follows.

(i) We present the first near-IR variable (rms) spectrum for an AGN and find that it shows the usual spectral inflection at $\sim 1~{}\mu$m, indicating that both the accretion disc and the hot dust contribute to the variability. Most importantly, we find that the spectral shape at wavelengths \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}1~{}\mum is similar in the variable (rms) and mean spectrum, which means that the same dust component dominates both the emission and the variations. Based on this, we constrain the hot dust properties.

(ii) We measure the response time in the near-IR using several techniques and find them to converge to a value of $\sim 70$ light-days. This is at the high end of the range of $\sim 40-80$ light-days measured previously by photometric reverberation campaigns. Assuming thermal equilibrium for optically thin dust, we derive the luminosity-based dust radius using our measurement of the temperature and an estimate of the irradiating luminosity, both obtained from the near-IR spectrum. We find that the two dust radii are in excellent agreement if we assume a wavelength-independent dust emissivity-law, i.e. a blackbody, which is appropriate for grains of relatively large sizes (of a few $\mu$m).

(iii) We measure an average dust temperature of $T\sim 1450$ K, which indicates that silicates have sublimated ($T_{\rm sub}\sim 1300-1500$ K) and graphite dominates the chemical composition. However, the measured temperature is several 100 K below the sublimation temperature of astrophysical carbon ($T_{\rm sub}\sim 1800-2000$ K) and we observe temperature and flux variations for the hot dust consistent with it experiencing heating and cooling by the central irradiating source. Since the observations also require the dust mass to change, which we find is destroyed by up to a maximum of $\sim 40\%$ and reforms on time scales of $\sim 5-6$ months during our campaign, this would imply that the inner dust-free region is enlarged by in-situ dust destruction proceses. This then means that the inner radius of the dusty torus is not set by sublimation and is expected to be luminosity-invariant, i.e. it is a “dusty wall”.

(iv) We have also considered the response function of a reprocessing accretion disc, which predicts $\tau\propto\lambda^{4/3}$. We find that, whereas the mean reverberation radii in the near-IR are consistent with $\sim 70$ light-days, the delay map is strongly skewed toward small values and peaks at radii consistent with dust in the accretion disc. The existence of this dust, which implies that the accretion disc extends out to relatively large radii (of $R_{\rm out}\sim 1700-3000\,r_{\rm g}$), is a prerequisite for the recent models of Czerny et al. (2017) and Baskin & Laor (2018). These models explain both the broad emission line region and the dusty torus as part of the same dusty outflow launched by radiation pressure from the outer regions of the accretion disc. It would be desirable to follow up on this result with a near-IR reverberation campaign of a higher cadence than ours.

In the future, we plan to use our near-IR spectroscopy of NGC 5548 to also study the coronal emission lines and the profile variations of the broad emission lines, both of which are expected to further constrain torus models. In general, our spectroscopic dust reverberation programme currently includes four more sources (Mrk 876, Mrk 110, Mrk 509 and 3C 273) and the results from these campaigns are forthcoming.

Acknowledgments

H.L. thanks Nancy Levenson for supporting this project in its initial stages and Chris Done for the many enlightening conversions on the physics of accretion discs. H.L., M.J.W. and D.K. acknowledge the Science and Technology Facilities Council (STFC) through grant ST/P000541/1 for support. C.P. acknowledges support from the National Science Foundation (NSF) grant no. 1616828. G.J.F. acknowledges support by NSF (1816537), NASA (ATP 17-ATP17-0141), and STScI (HST-AR- 15018). K.H. acknowledges support from STFC grant ST/R000824/1.