Optical variability power spectrum analysis of blazar sources on intranight timescales
DDraft version January 5, 2021
Typeset using L A TEX twocolumn style in AASTeX63
Optical variability power spectrum analysis of blazar sources on intranight timescales
Arti Goyal Astronomical Observatory of the Jagiellonian University, Orla 171, 30-244 Krakow, Poland (Received; Revised; Accepted)
Submitted to ApJABSTRACTWe report the first results of a systematic investigation to characterize blazar variability powerspectral densities (PSDs) at optical frequencies using densely sampled (5–15 minutes integration time),high photometric accuracy ( (cid:46) ∼
15 minutes. Our sample consists of 14 optically bright blazars, including nine BLLacertae objects (BL Lacs) and five flat-spectrum radio quasars (FSRQs) which have shown statisticallysignificant variability during 29 monitoring sessions. We model the intranight PSDs as simple power–laws and derive the best-fit slope along with uncertainty using the ‘power spectral response’ method.Our main results are the following: (1) on 19 out of 29 monitoring sessions, the intranight PSDs showan acceptable fit to simple power-laws at the rejection confidence ≤ ∼ ≥
3) noise stochastic processes; (3) the average PSD slopes for the BL Lacs and FSRQsare indistinguishable from one another; (4) the normalization of intranight PSDs for individual blazarsources which were monitored on more than one occasion turns out to be consistent with one anotherwith a few exceptions. The average PSD slope, 2.9 ± σ uncertainty) is steeper than the red-noisetype character of variability found on longer timescales (many decades to days), indicative of a cutoff inthe variability spectrum on timescales around a few days at the synchrotron frequencies of the emissionspectrum. Keywords:
Galaxies:active–galaxies:jets–acceleration of particles–radiation mechanisms: non-thermal INTRODUCTIONIntense emission and rapid flux variability provideimportant clues to understand the underlying physi-cal process operating in the relativistic, magnetized jetsof blazar sources. This subset of active galactic nu-clei (AGN) is comprised of BL Lacertae objects (BLLacs) and the high optical polarization flat-spectrum ra-dio quasars (FSRQs), with jets launched from a super-massive black hole (SMBH)–accretion disk systems(fora recent review see, Hovatta & Lindfors 2019). Theflux variability is observed at all frequencies of the elec-tromagnetic spectrum on both long-term (decades to (cid:39) day; Ulrich et al. 1997; Aller et al. 1999; Ghoshet al. 2000; Gopal-Krishna et al. 2011; Goyal et al. 2012;
Corresponding author: Arti [email protected]
Gupta et al. 2016; Gaur et al. 2019; Abdalla et al. 2017;Abdo et al. 2010a) and intranight timescales ( ≤ day;Aharonian et al. 2007; Albert et al. 2007; Goyal et al.2013a; Bachev 2015; Ackermann et al. 2016; Nalewa-jko 2017; Zhu et al. 2018; Shukla et al. 2018). Thetwo-component broadband spectral energy distributionis nonthermal radiation arising within the relativistic jet(Ghisellini & Tavecchio 2008). Within the leptonic sce-nario, the particle (electron and positron) pairs, acceler-ated to GeV/TeV energies, produce synchrotron radia-tion in the presence of magnetic field at lower frequencies(radio–to–optical/X-rays) and inverse Comptonizationof the seed photons (same or thermal from the accre-tion disk) by the synchrotron emitting particles produceemission at higher frequencies (X-rays–to–TeV γ − rays).Alternatively, direct synchrotron radiation by the pro-tons accelerated to PeV/EeV energies or the emissionfrom secondaries can give rise to high energy radia- a r X i v : . [ a s t r o - ph . H E ] J a n A. Goyal tion within the hadronic scenario (e.g., Blandford et al.2019). Well–defined flare emission has often been at-tributed to particle acceleration mechanisms related toshocks in the jet (Spada et al. 2001; Marscher et al. 2008)while turbulence can mimic observed fluctuations onlong-term as well as small timescales (Marscher 2014).Annihilation of magnetic field lines at the reconnectionsites within the jet plasma can also impart energy tothe particles (Giannios 2013; Sironi et al. 2015); thisscenario is supported by the recent detection of minute-like variability at GeV energies for the blazar 3C 279(Shukla & Mannheim 2020). On the other hand, fluxchanges on the intranight timescales have also been as-sociated with variable Doppler boosting factors relatedto changes in the viewing angle of the emitting plasma(e.g., Gopal-Krishna & Wiita 1992), although in such ascenario frequency–independent variability is expected(see, in this context, Pasierb et al. 2020).The noise-like appearance of blazar light curves haveprompted efforts to investigate variability power spec-tral densities (PSDs) which is a distribution of vari-ability amplitudes over different Fourier frequencies(=timescale − ). The blazar PSDs are mostly repre-sented by power-law shapes defined as P( ν k ) ∝ ν − βk where β (=1–3) is the slope and ν k is the temporalfrequency which indicate that variability is a corre-lated colored–noise type stochastic processes (see, Goyal2020, and references therein). Specifically, β (cid:39) (cid:39) (cid:38) β (cid:39) uncorrelated , white–noise type stochastic process (Press 1978; Schroeder1991). For a colored noise–type stationary stochasticprocess, one expects the slope of PSDs to change to 0on longer timescales to preserve the finite variance of theprocess, leading to a relaxation timescale beyond whichthe variations should be generated due to uncorrelatedprocesses. Moreover, it also means that different randomrealizations of the process will have different statisticalmoments (e.g., mean, sigma) due to statistical fluctua-tion of the process itself and not due to the change ofnature of the process which indicates that the process isweakly non-stationary (Vaughan et al. 2003). Fluctua-tions resulting from such stochastic processes obey cer-tain probability distributions, so the light curves tendto produce predictable PSDs. PSD slope and normal-ization, as well as the breaks, are of particular interestas they carry information about the parameters of thestochastic process and the ‘characteristic timescales’ inthe system which can be related to physical parame-ters shaping the variability, such as the size of the emis-sion zone or the particle cooling timescales (Sobolewska et al. 2014; Finke & Becker 2014; Chen et al. 2016). Thenoise-like appearance of light curves has been modeledwhere the aggregate flux arises from many cells behinda shock (Calafut & Wiita 2015; Pollack et al. 2016, seealso, Marscher (2014) who models the flux and polar-ization light curves but does not provide PSDs) againstthe emission from well-defined flares which could be at-tributed to single emission zones (Hughes et al. 1985;Abdo et al. 2010b). The models of Calafut & Wi-ita (2015) and (Pollack et al. 2016) compute the lightcurves and the PSDs, which are shaped by the combina-tion bulk Lorentz factor fluctuations and the turbulencewithin the jets. In their hydrodynamic simulations of2D jets, the changes in bulk Lorentz factor produce PSDslopes in the range 2.1 to 2.9 while the turbulence pro-duces PSD slopes in the range 1.7 to 2.3, respectively(Pollack et al. 2016). The model of O’Riordan et al.(2017), on the other hand, hypothesizes that the turbu-lence in the magnetically arrested disk (MAD) shapesthe variability of synchrotron and IC emission compo-nents from the jet, with a cutoff of variability power atthe timescales governed by the light crossing time of theevent horizon of the SMBH.Unlike the long–term variability timescales where theslopes of PSDs of multiwavelength variability have beenestimated for large samples of blazar sources (in particu-lar, β ∼ β ∼ γ − rays, Max-Moerbeck et al. 2014; Park & Trippe 2017; Nilsson et al.2018; Meyer et al. 2019; Goyal 2020), such studies re-main scarce on intranight timescales. This is due to thefact that it requires continuous pointing of an observingfacility to a single target for many hours which is usu-ally not feasible due to scheduling constraints, weather,limited photon sensitivity (at high energies), etc. Theintranight PSD slopes at GeV and TeV energies exhib-ited β ∼ − Fermi-
LATand the High Energy Stereoscopic System data, butonly when the blazars were in a flaring state (Aharo-nian et al. 2007; Ackermann et al. 2016). Zhang et al.(2019) obtained the X-ray PSD slopes of 1.5, 3.1 and1.4 using the 40–180 ks long
Suzaku observations. Inanother study, Bhattacharyya et al. (2020) obtained in-tranight X-ray PSD slopes equal to 2.7, 2.6, 1.9 and, 2.7for the Mrk 421 and 2.2, 2.8, and 2.9, respectively, for thePKS 2155 −
304 using the 30–90 ks long XMM–
Newton observations. Goyal et al. (2017) obtained the opticalintranight PSD slopes in the range 1.5–4.1 for five mon-itoring sessions of the BL Lac PKS 0735+178. Wehrleet al. (2013) and Wehrle et al. (2019) obtained the PSDslopes of blazar sources using the
Kepler –satellite in therange 1.2–3.8 on timescales ranging from half a day to ptical intranight variability PSDs Table 1.
Sample properties.IAU name RA(J2000) Dec(J2000) SED z V–mag M BH Reference for M BH (h m s) (d (cid:48) (cid:48)(cid:48) ) (M (cid:12) )(1) (2) (3) (4) (5) (6) (7) (8)0109+224 01 12 05.824 +22 44 38.78 BL Lac c f c f × Raiteri et al. (2007)0420 −
014 04 23 15.800 −
01 20 33.06 FSRQ c f × Liang & Liu (2003)0716+714 07 21 53.448 +71 20 36.36 BL Lac c f × Liang & Liu (2003)0806+315 08 09 13.440 +31 22 22.90 BL Lac d g e h × Wu et al. (2002)0851+202 a
08 54 48.874 +20 06 30.64 BL Lac c f × Liang & Liu (2003)1011+496 10 15 04.139 +49 26 00.70 BL Lac c f × Wu et al. (2002)1156+295 11 59 31.833 +29 14 43.82 FSRQ c f × Liang & Liu (2003)1216 −
010 12 18 34.929 −
01 19 54.34 BL Lac e i c f × Liang & Liu (2003)1253 − b
12 56 11.166 −
05 47 21.52 FSRQ c f × Sbarrato et al. (2012)1510 −
089 15 12 50.532 −
09 05 59.82 FSRQ d j × Sbarrato et al. (2012)1553+113 15 55 43.044 +11 11 24.36 BL Lac c f Note — (1) the name of the blazar following the IAU convention. a also known as OJ 287; b also known as3C 279; (2) right ascension; (3) declination; (4) SED classification. c Healey et al. (2008); d V´eron-Cetty &V´eron (2006); e Plotkin et al. (2008); (5) spectroscopic redshift. f Healey et al. (2008); g Falco et al. (1998); h Bade et al. (1998); i Dunlop et al. (1989); j Thompson et al. (1990). (6) typical optical V–band magnitude(V´eron-Cetty & V´eron 2010); (7) mass of the SMBH; (8) reference for the mass of the SMBH. few months. Recently, Raiteri et al. (2020) obtained thePSD slope ∼ TESS )–2 min integration time light curve forthe blazar 0716+714 at variability timescales between amonth and few minutes.In this respect, a few optical observatories with 1–2 m class optical telescopes and fitted with CCDs havebeen devoted to blazar/AGN monitoring programs since1990 (see, for a review, Gopal-Krishna & Wiita 2018).The goal of the present study is to characterize the in-tranight variability of a large sample of blazars using theARIES monitoring program which was carried out be-tween 1998 to 2010, the results of which are presented inGoyal et al. (2013a). The paper is organized as follows.Sample selection is given in Section 2 while Section 3provides the details on the analysis method, in partic-ular, the derivation of power spectral densities and theestimation of best-fit PSD shapes using extensive nu-merical simulations of light curves. Section 4 providesthe main results while a discussion and conclusions aregiven in Section 5. SAMPLE The blazar light curves studied here are obtained fromthe samples of Goyal et al. (2013a) who studied theintra-night variability properties of different types of ac-tive galactic nuclei (AGN) using 262 intranight lightcurves. The AGNs monitored belong to radio–quietquasar, radio–intermediate quasar, radio–loud quasarand blazar types, including BL Lacs and the FSRQs.Goyal et al. blazar sample consists of 24 sources mon-itored on 85 sessions. The details of data gather-ing, reduction procedure, generation of differential lightcurves, and the statistical tests used to infer intranightvariability are given in Goyal et al. (2013a) which webriefly describe here. On each monitoring session, con-tinuous CCD observations ( > Image reduction and analysis facility(IRAF) software and the instrumental magnitudes ofthe target blazar and the comparison stars on the sameCCD chip were derived using aperture photometry. Therelative instrumental magnitude of the blazar was com-puted against one steady ‘star–star’ pair, thereby pro-
A. Goyal
15 16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.7 4.3 p b b
17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level (b)
Raw per.Log-binned per.Best-fit =3.9 2.8 p b b
17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.7 3.2 p b b
14 15 16 17 18 19 20 21 22
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level (b)
Raw per.Log-binned per.Best-fit =4.0 2.5 p b b
16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.7 2.8 p b b Figure 1.
The intra-night variability power spectrum of blazar sources obtained in this study. Panel (a) presents the lightcurve on a linear scale (see text). Panel (b) presents the derived power spectrum down to the Nyquist sampling frequency of the(mean) observed data. The dashed line shows the ‘raw’ periodogram while the blue triangles and red circles give ‘logarithmicallybinned power spectrum’ and the best-fit power spectrum, respectively. The error on the best-fit PSD slope corresponds to a98% confidence limit. The dashed horizontal line corresponds to the statistical noise floor level due to measurement noise. Panel(c) shows the probability curve as a function of the input power spectrum slope. The source name and date of monitoring arepresented at the top of each panel. ptical intranight variability PSDs
18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level (b)
Raw per.Log-binned per.Best-fit =2.7 2.9 p b b
20 21 22
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.5 2.3 p b b
16 17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.0 2.4 p b b
16 17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.8 3.8 p b b
19 20 21 22 23 24 25
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.8 2.8 p b b Figure 1. (continued) ducing two differential light curves (DLCs) for a blazar on a given monitoring session. The differential photom-
A. Goyal
16 17 18 19 20 21
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.9 2.8 p b b
15 16 17 18 19 20 21 22
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.2 2.4 p b b
14 15 16 17 18 19 20
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.9 2.9 p b b
16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.1 2.9 p b b
16 17 18 19 20 21 22
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.3 2.8 p b b Figure 1. (continued) etry technique is widely used in variability studies as it counters the spurious AGN variability occurring due to ptical intranight variability PSDs
14 15 16 17 18 19 20 21
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.5 2.9 p b b
14 15 16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =1.4 3.1 p b b
15 16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.6 4.0 p b b
15 16 17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =4.0 2.8 p b b
17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.6 2.0 p b b Figure 1. (continued) varying atmospheric conditions (variable seeing or pres- ence of thin clouds during the monitoring session) as
A. Goyal
16 17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.3 1.8 p b b
19 20 21 22 23 24 25
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =1.9 1.7 p b b
17 18 19 20 21 22 23 24
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3-2 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.3 1.7 p b b
14 15 16 17 18 19 20
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3-2-1 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.4 2.2 p b b
17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.3 2.8 p b b Figure 1. (continued) any changes seen in blazar light should be accompanied by the same changes in comparison starlight, thereby ptical intranight variability PSDs
17 18 19 20 21 22
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =4.0 2.8 p b b
16 17 18 19 20
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =2.4 2.0 p b b
16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.9 2.3 p b b
16 17 18 19 20 21 22 23
UT (hr) F l u x den s i t y ( a . u . ) log ( k /d -1 ) -8-7-6-5-4-3 l og ( P ( k ) / (r m s / m ean ) d ) Noise Floor Level
Raw per.Log-binned per.Best-fit =3.1 2.7 p b b Figure 1. (continued) keeping the difference unaltered. Next, we used F − testto infer the statistical significance of variability at sig-nificance levels, α = 0.01 and 0.05, corresponding to p value > > F − test statistic for a DLC turned out to be more thanthe critical value at α =0.05(0.01), the DLCs were as-signed a ‘probable variable (confirmed variable)’ status.If the test statistic turned out to be smaller than thecritical value at α =0.05 for the DLC, it was assigned a‘non-variable’ status. In the present study, all the mon- itoring sessions where both DLCs showed a ‘confirmedvariable’ status are used. We also included two mon-itoring sessions where one of the two DLCs showed a‘probable variable’ status while the other showed a ‘con-firmed variable’ status. The monitoring sessions wherethe blazar DLCs showed ‘non-variable’ status are notused in this analysis as the PSDs of ‘non-variable’ lightcurves will be consistent with β ∼
0, resulting from fluc-tuations arising from measurement errors. The abovecriterion reduced the sample to 15 blazars and 34 in-0
A. Goyal tranight light curves. We further exclude the BL Lac ob-ject PKS 0735+178 from the sample as the PSD analysisof its five intranight light curves has been reported inGoyal et al. (2017). Therefore, the current sample con-sists of 14 blazars which have shown statistically signif-icant variability on 29 intranight monitoring sessions .Table 1 lists the basic properties of the blazars.Next, we converted the differential magnitudes ofblazars from logarithmic scale to linear scale using therelation F obs =F × − . × m obs , where F (=1) is the ar-bitrary zero point magnitude flux, m obs is the differen-tial blazar magnitude relative to one comparison starand F obs is the corresponding differential flux density asit contains the contribution from the steady compari-son star flux. The errors in were derived using stan-dard error propagation (Bevington & Robinson 2003)and scaled by a factor 1.54 to account for the underes-timation of the photometric errors by the IRAF (see,for details, Goyal et al. 2013b, and references theirin).We note that the differential blazar flux density canbe scaled to proper fluxes using the appropriate F forthe R–band and the apparent magnitudes of compar-ison star used. Figure 1 (panel a) shows flux densityintranight light curves. PSD ANALYSIS3.1.
Derivation of PSDs: discrete Fourier transform
Since the aim of the study is to obtain reliable shapesof PSDs, we subject the light curves Fourier transforma-tion using the discrete Fourier transform (DFT) method(see, for details, Goyal 2020, and references therein).The fractional rms-squared-normalized periodogram isgiven as the squared modulus of its DFT for the evenlysampled light curve f ( t i ), observed at discrete times t i and consisting of N data points and the total monitoringduration T , P ( ν k ) = 2 Tµ N (cid:40)(cid:34) N (cid:88) i =1 f ( t i ) cos(2 πν k t i ) (cid:35) + (cid:34) N (cid:88) i =1 f ( t i ) sin(2 πν k t i ) (cid:35) (cid:41) , (1)where µ is the mean of the light curve and is subtractedfrom the flux values, f ( t i ). The DFT is computed forevenly spaced frequencies ranging from the total dura-tion of the light curve down to the Nyquist sampling fre-quency ( ν Nyq ) of the observed data. Specifically, the fre-quencies corresponding to ν k = k/T with k = 1 , ..., N/ Light curves can be obtained upon request. ν Nyq = N/ T , and T = N ( t k − t ) / ( N −
1) are consid-ered. The normalized periodogram as defined in Eq. 1corresponds to total excess variance when integratedover positive frequencies. The constant noise floor levelfrom measurement uncertainties is given as (e.g., Isobeet al. 2015; Vaughan et al. 2003)P stat = 2 T µ N σ . (2)where, σ stat = (cid:80) j = Nj =1 ∆ f ( t j ) /N is the mean variance ofthe measurement uncertainties on the flux values ∆ f ( t j )in the observed light curve at times t j , with N denotingthe number of data points in the original light curve.The intranight light curves are roughly evenly sampledbut the application of the DFT method requires stricteven sampling of the time series, otherwise, the ‘spectralwindow function’ corresponding to the sampling timesgives a non-zero response in the Fourier-domain, result-ing in false powers in the periodograms (see, AppendixA of Goyal 2020; Deeming 1975). Therefore, in orderto perform the DFT, we obtained the regular samplingonly by linearly interpolating between the two consecu-tive observed data points with an interpolation intervalof 1 minute which is roughly 5–15 times smaller thanthe original sampling interval (see column 7 of Table 2).Even though the choice of interpolation interval is arbi-trary, we note that it cannot be longer than the meansampling interval. We tested our procedure by also usingan interpolation interval about half of the original sam-pling interval which did not change the results. We referthe reader to Goyal et al. (2017) and Max-Moerbecket al. (2014) for a discussion on the distortions intro-duced in the PSDs due to the discrete sampling andthe finite duration of the light curve, known as ‘red-noise leakage’ and ‘aliasing’ respectively. To minimizethe effects of red-noise leak, the PSDs are generatedusing the ‘Hanning’ window function (e.g., Press et al.1992; Max-Moerbeck et al. 2014). Aliasing, on the otherhand, contributes an equal amount of power (around theNyquist frequency) to the periodograms (Uttley et al.2002), hence will not distort the shape of PSDs.The periodogram obtained using equation (1), knownas the ‘raw’ periodogram, provides a noisy estimate ofthe spectral power as it consists of independently dis-tributed χ variables with two degrees of freedom (DOF)(Timmer & Koenig 1995; Papadakis & Lawrence 1993;Vaughan et al. 2003). Therefore, a number of PSD esti-mates should be averaged in order to obtain a reliable es-timate of the spectral power. The periodograms fallingwithin a factor of 1.6 in frequency range are averagedwith the representative frequency taken as the geomet-ric mean of each bin (Isobe et al. 2015; Goyal et al. 2017; ptical intranight variability PSDs Table 2.
Summary of the observations and the PSD analysis.IAU name Date of obs. Tel. T obs N obs ψ T mean log (P stat ) log ( ν k ) range β ± err p ∗ β (hr) (min) ( r ms m ean ) d (d − )(1) (2) (3) (4) (5) (6) (7) (8) (9) (10) (11)0109+224 2005 Oct 29 ST 7.1 36 3.98 11.9 − ± − ± − ± − ± −
014 2003 Nov 19 ST 6.7 38 2.22 10.55 − ± − ± − ± − ± − ± − ± − ± − ± − ± − ± − ± − ± − ± − ± −
010 2002 Mar 16 ST 8.2 22 14.27 22.35 − ± − ± − ± −
055 2006 Jan 26 ST 4.7 21 3.23 13.56 − ± − ± − ± −
089 2009 May 1 ST 6.0 25 7.73 14.45 − ± − ± − ± − ± − ± Note — (1) name of the blazar following the IAU convention; (2) the date of observations; (3) the telescope facilityused. ST = 1 m Sampurnanand Telescope of Aryabhatta Research Institute of Observational Sciences, India; IGO= 2 m IUCAA-Girawali Observatory of Inter-University Centre of Astronomy and Astrophysics, India. (4) theduration of the observed light curve; (5) number of data points in the light curve; (6) peak-to-peak variabilityamplitude (Eq. 9 of Goyal et al. 2013a); (7) the mean sampling interval for the observed light curve (light curveduration/number of data points); (8) the noise level in PSD due to the measurement uncertainty; (9) the temporalfrequency range covered by the binned logarithmic power spectra; (10) the best-fit power-law slope of the PSDalong with the corresponding errors representing 98% confidence limit (see Section 3.2); (11) corresponding p β . ∗ power law model is considered as a bad-fit if p β ≤ ≥
90% (Section 3.2). A. Goyal
Goyal 2020). Except for the first bin, this choice of bin-ning factor provides at least two periodograms in eachfrequency bin.Since the observed power-spectrum is related to the‘true’ power spectrum by P ( ν k ) = P true ( ν k ) χ for anoise-like process (Papadakis & Lawrence 1993; Tim-mer & Koenig 1995; Vaughan et al. 2003). The trans-formation to log-log space, offsets the the observed pe-riodograms aslog [ P ( ν k )] = log [( P true ( ν k )] + log (cid:104) χ (cid:105) . (3)This offset is the expectation value of χ distributionwith 2 DOF in log-log space and is equal to − Estimation of the spectral shape: PSRESP method
Since the aim of the present study is to derive shapesof intranight PSDs, we use the ‘power spectral response’(PSRESP) method (e.g., Uttley et al. 2002; Chatterjeeet al. 2008; Max-Moerbeck et al. 2014; Isobe et al. 2015;Meyer et al. 2019; Goyal 2020) which further mitigatesthe deleterious effects of red-noise leak and aliasing. Inthis method, an (input) PSD model is tested against theobserved PSD. The estimation of best-fit model param-eters and their uncertainties is performed by varyingthe model parameters. To achieve this, a large num-ber of light curves are generated with a known underly-ing power-spectral shape using Monte Carlo (MC) sim-ulations. Rebinning of the light curve to mimic thesampling pattern and interpolation is performed for theDFT application. The DFT of such light curve givesthe distorted PSD due to effects mentioned above. Av-eraging large number of such PSDs gives the mean ofthe distorted model (input) power spectrum. The stan-dard deviation around the mean gives errors on the mod-eled (input) power spectrum. The goodness of fit of themodel is estimated by computing two functions, similarto χ , defined as χ = ν k = ν max (cid:88) ν k = ν min [log P sim ( ν k ) − log P obs ( ν k )] ∆log P sim ( ν k ) (4)and χ , i = ν k = ν max (cid:88) ν k = ν min [log P sim ( ν k ) − log P sim , i ( ν k )] ∆log P sim ( ν k ) , (5)where log P obs and log P sim , i are the observed andthe simulated log-binned periodograms, respectively. log P sim and ∆log P sim are the mean and the standarddeviation obtained by averaging a large number of PSDs; k represents the number of frequencies in the log-binnedpower spectrum (ranging from ν min to ν max ), while i runs over the number of simulated light curves for agiven β .Here χ determines the minimum χ for the modelcompared to the data and the χ values determine thegoodness of the fit corresponding to the χ . We notethat χ and χ are not the same as a standard χ dis-tribution because log P obs ( ν k )’s are not normally dis-tributed variables since the number of power spectrumestimates averaged in each frequency bin are small (Pa-padakis & Lawrence 1993). Therefore, a reliable good-ness of fit is computed using the distribution of χ values. For this, the χ values are sorted in ascendingorder. The probability or p β , that a given model can berejected is then given by the percentile of χ distribu-tion above which χ is found to be greater than χ for a given β (success fraction; Chatterjee et al. 2008). Alarge value of p β represents a good–fit in the sense thata large fraction of random realizations of the model (in-put) power spectrum are able to recover the shape ofthe intrinsic PSD. Therefore, this analysis essentiallyuses the MC approach toward a frequentist estimationof the quality of the model compared to the data. Thisis a well-known approach to estimate the goodness of fitwhen the fitting statistic is not well understood (see, fordetails, Press et al. 1992).In this study, the light curve simulations are per-formed using the method of Emmanoulopoulos et al.(2013) which preserves the probability density function(PDF) of the flux distribution as well as the underlyingpower spectral shape. In addition to assuming the powerspectral shape, the method requires supplying a value ofmean and standard deviation ( σ ) of the flux values toreproduce the flux distribution and match the variance(Meyer et al. 2019). We have assumed single power-lawPSDs with a given β (to reproduce the PSD shape) andsupplied mean and σ of the logarithmically transformedflux values which is found to be an adequate represen-tation of flux distribution on shorter ( (cid:46) days) timescalesfor a few cases (e.g., H. E. S. S. Collaboration et al. 2010;Kushwaha & Pal 2020). For this purpose, the meanand the σ are computed by fitting a Gaussian functionto the flux distribution. Finally, the measurement errorsin the simulated flux values were incorporated by addinga Gaussian random variable with mean 0 and standarddeviation equal to the mean error of the measurementuncertainties on the observed flux values (Meyer et al.2019; Goyal 2020). In such a manner, 1,000 light curvesare simulated in the β range 0.1 to 4.0, with a step ptical intranight variability PSDs p β value and the uncertainty isgiven as 2.354 σ of the p β curve where σ is the standarddeviation of fitted Gaussian. This gives roughly a 98%confidence limit on the best-fit PSD slope.Details on the intranight light curves used for the anal-ysis and the derived PSDs, along with the best-fit PSDsand the maximum p β , are summarized in Table 2. Fig-ure 1 presents the analyzed light curves (panel a), thecorresponding best-fit PSD (panel b) and the probabil-ity distribution curves ( p β as a function of β ; panel c)for the duration of the light curve down to the meansampling intervals for the given light curve. In our anal-ysis, we have not subtracted the constant noise floorlevel (shown by the dashed horizontal lines in panel bof the figure), as some of the data points are below thislevel. The PSRESP method also allows us to computethe rejection confidence for the input PSD shape; themaximum probability lower than 10% means that therejection confidence (1– p β value) is higher than 90% forthe (input) PSD model. In our analysis, we use p β < β is drawn from a Gaussian distri-bution of mean and standard deviation equal to β anderror/2.354 (note that Table 2 reports errors equal to2.354 σ ). The mean β is computed. These two stepsare repeated 500 times. The error is given as the stan-dard deviation of the distribution of mean β values. Inaddition, we provide the ν k P ( ν k ) vs. ν k curves in Fig-ure 3 for blazars observed on more than one occasion, tocompare the ‘square’ fractional variability on timescalesprobed by our analysis (Goyal 2020). RESULTSIn this study, we have derived the optical intranightvariability PSDs of the blazar sources, covering the tem-poral frequency range from 10 . day − to 10 . day − (timescales corresponding to 7.4 hours and ∼
15 min-utes). The intranight light curves showed modest in-tranight variability with peak–to–peak variability am-plitudes, > N u m be r pe r b i n Total( =2.9 0.3)BL Lacs( =3.1 0.3)FSRQs( =2.6 0.4)
Figure 2.
Histograms of the best-fit PSD slopes derived forthe entire blazar sample (cyan line; 10 sources and 19 moni-toring sessions), BL Lacs (red line; seven sources and 13 mon-itoring sessions), and FSRQs (blue line; three sources and sixmonitoring sessions), respectively. The sample mean alongwith 1 σ uncertainty estimated using the bootstrap methodfor different groups is given in parentheses (see Section 3.2).
1. Out of the 29 intranight light curves analyzed inthe present study, the PSDs shows an acceptablefit to the single power-law spectral shapes for 19monitoring sessions (see Section 3.2). The maxi-mum p β is higher than 10% and reaches as high as100% for these sessions (column 10; Table 2; panelc of Figure 1).2. For these 19 acceptable PSD fits, the simplepower-law slopes range from 1.4 to 4.0 (albeit witha large scatter); consistent with a statistical char-acters of red ( β ∼
2) and black ( β ≥
3) noise stochas-tic processes (Table 2, Figure 2). The mean β turns out to be 2.9 ± σ uncertainty) for blazarsources.3. The computed mean value PSD slopes for the BLLac objects (seven sources and 13 light curves)and FSRQs (three sources and six light curves)are 3.1 ± ± σ uncertainty (Figure 2).4. The PSD slopes for a few sources whose intranightPSDs show an acceptable fit to single power-law onmultiple occasions are consistent with each other(column 9; Table 2).5. The normalization of the PSDs for the sourcesmonitored on different epochs turns out to be con-sistent with each another within 1 σ uncertainty forthe blazars 1156+295, 1219+285. However, oneorder of magnitude change is noted in the normal-ization of PSDs between 2003 November 19 and4 A. Goyal log ( k /d -1 ) -6-5-4-3-2 l og ( k P ( k ) / (r m s / m ean ) ) log ( k /d -1 ) -6-5-4-3-2 l og ( k P ( k ) / (r m s / m ean ) ) log ( k /d -1 ) -5-4-3-2-1 l og ( k P ( k ) / (r m s / m ean ) ) log ( k /d -1 ) -5-4-3-2 l og ( k P ( k ) / (r m s / m ean ) ) log ( k /d -1 ) -5-4-3-2-1 l og ( k P ( k ) / (r m s / m ean ) ) log ( k /d -1 ) -6-5-4-3-2 l og ( k P ( k ) / (r m s / m ean ) ) Figure 3. ν k P( ν k ) PSDs for individual blazars which showed acceptable fit in the analysis. The lines show the log-binnedperiodograms and the filled symbols show mean and standard deviation of best-fit PSDs given by the PSRESP method fordifferent epochs. − − − − β > −
014 on 2009 October 25, 0806+524 on 2005 Febru-ary 4, 1011+496 on 2010 February 19, and 1156+295 on2012 April 1. A possible cause of this could be the lim-ited number of data points in the studied light curves (20–67; column 4 of Table 2). Aleksi´c et al. (2015) stud-ied the effects of changing the number of data points inthe light curve and the estimation of the best-fit PSDslope (and uncertainty) using the long-term multiwave-length light curves having ≥
30 data points for the blazarMrk 421. They note that the location of the maximumin the probability distribution curve does not changenoticeably for different binning factors but the width,shape, and amplitude change significantly; however, itis unclear if the broadening of the probability distribu-tion curve and hence the estimation of the uncertainty inthe best-fit PSD slope is related to the gradual increaseof binning factors (i.e., a decrease of a number of datapoints) for different light curves (see, Figure 4 of Aleksi´cet al. 2015). Also, we note that the reported uncertain-ties on the best-fit X-ray intra-night PSD slopes usingthe PSRESP method also show large scatter, despitehaving >
300 data points in the examined light curves forthe AGNs Mrk 421, PKS 2155 − DISCUSSION AND CONCLUSIONSWe report the first systematic study to characterizethe intranight variability PSD properties comprising of ptical intranight variability PSDs ∼
15 minutes.All the analyzed light curves were of duration ≥ ∼ (cid:46) β for the entire sample is ∼ β ’s forthe BL Lacs and the FSRQs subclasses (due to the smallnumber of sources and the intranight light curves ana-lyzed) give ∼ β =1.4 to 2.1).The majority of obtained slopes could be reconciled ifthe intranight fluctuations are driven by changes in thebulk Lorentz factors of the jet, provided that the turbu-lence is dominant on smaller than few minutes timescales( β ∼ Kepler- satellite data withlong duration ( >
75 days), nearly uniformly sampledlight curves with sampling intervals 30 minutes (long–cadence) and 1 minute (short cadence data for the blazarOJ 287), respectively. First, the PSD slopes obtained fortheir sample using the long–cadence data range between β ∼ ν k ∼ -6.5 Hz and ∼ -5.0 Hz with the slope tend-ing to white noise at higher frequencies (timescales ≤ Kepler ’s long– cadence data ( ∼ β ∼ ν k = − − ν k ≥ − . Hz. Ourintranight PSD slopes for the OJ 287, obtained on threeseparate occasions, are 3.8, 2.9, and 3.1, respectively,consistent with their result on overlapping variabilityfrequencies (Table 2).In Goyal et al. (2017), Goyal et al. (2018), and Goyal(2020), based on PSD analyses of long-term variabilityusing decade–long GHz–band radio–to–TeV γ − ray lightcurves of a few selected blazar sources, we hypothe-sized that the broadband emission is generated in anextended yet highly turbulent jet. The variability ap-peared to be driven by a single stochastic process atsynchrotron frequencies but seemed to require the lin-ear superposition of two stochastic processes at IC fre-quencies with relaxation timescales ≥ ∼ days, respectively. Stochastic fluctuations in the lo-cal jet conditions (e.g., velocity fluctuations in the jetplasma or magnetic field variations) lead to energy dis-sipation over all spatial scales. The radiative responseof the accelerated particles is delayed with respect tothe input perturbations and this forms the red–noisesegment of the PSD at synchrotron frequencies. At ICfrequencies, however, due to inhomogeneities in the lo-cal photon population available for upscattering, the ad-ditional relaxation timescale of about ∼ one day, i.e.,the light crossing time of the emission region, can re-sult in a jet with Doppler boosting factor, 30, formingthe pink-noise segment of the PSD. The steeper thanred–noise PSD slopes on intranight timescales obtainedin this analysis against the strict red–noise character oflong–term variability at optical frequencies ( β ∼
2; Chat-terjee et al. 2008; Goyal et al. 2017; Nilsson et al. 2018;Goyal 2020), indicate a cutoff of variability power ontimescales around ∼ days. We note that such a cutoff ofvariability power on timescales ∼ days have been notedin the X–ray PSD of the blazar Mrk 421 for which the X-ray emission, although it originates in the non-thermaljet, it is believed that the variability process is driven byaccretion disk processes (Chatterjee et al. 2018). More-over, our conclusion is only tentative, as joint analysis offull variability spectrum using long-term and intranightdata, covering many orders of frequencies without gaps,is needed to reach robust conclusions. The normaliza-tion of PSDs for a few sources which were monitoredon multiple occasions turns out to be consistent withone another within 1 σ uncertainty with a few excep-tions. For the blazars 0420 − A. Goyal between different epochs (Figure 3). This indicates ahint of non-stationarity of the variability process on in-tranight timescales (similar conclusions are obtained forthe intranight X–ray variability of the blazar Mrk 421for which the intranight light curves are modeled as anon-stationary stochastic process; Bhattacharyya et al.2020).At this point, we note that the duty cycle of in-tranight blazar variability at optical frequencies is foundto be ∼
40% when monitored for a duration > > daysto years; Stalin et al. 2004; Sagar et al. 2004; Gopal-Krishna et al. 2011; Goyal et al. 2012) which exhibitsa red–noise character down to a few days timescales(see above). However, ∼
60% of the monitoring ses-sions, these sources turned out to be non-variable atshort timescales meaning that small-scale flux variabil-ity, if present, is below the measurement uncertainties.This would imply, occurring intermittently, a cutoff ofvariability power on timescales longer than ∼ ∼ solar mass SMBH; Begelman et al.2008). Using the black hole masses of blazars studiedhere (column 7; Table 1), the light crossing time of theevent horizon ranges from ∼ ∼ − ∼
24 seconds–12 minutes in the observer’sframe assuming typical bulk Lorentz factor, 10, for thejet plasma (Lister et al. 2016). Such timescales are notcovered by us, given the typical sampling intervals ∼ > Facilities:
NED, ST:1.0m, IGO:2.0m ptical intranight variability PSDs